So now that I have the initial values for the three elements – level, trend and seasonality, I am onto the most challenging part of Holt Winters algorithm. That is to calculate $latex \alpha$, $latex \beta$, and $latex \gamma$ that yield the least RMSE. With my limited knowledge of programming (if anyone who happens to stumble upon this blog knows a better way, please kindly suggest in the comment section), I decided to recursively find the best $latex \alpha$, $latex \beta$, and $latex \gamma$. From the get go, it does not seem very efficient. Because even if I do not want any more granularity below nearest one percent, I will have to iterate at least 1.030301 million times ($latex 101^3$). But if the speed is not an issue (though it was for us), this does the job.

***You can get more efficiency by using multi-threading and the ladder and jar algorithm but it still is not in the practical range of speed for web service. Seth Hoenig, our super genius developer, converted my Python code into C and that’s when the speed became so fast that it’s almost instantaneous.

In order to calculate $latex \alpha$, $latex \beta$, and $latex \gamma$, I first have to set up algorithms that calculate subsequent values of level, trend and seasonality and smoothed values.

- Level: $latex s_t = \alpha \frac{x_t}{c_{t-L}} + (1 – \alpha)(s_{t-1} + b_{t – 1})$
- Trend: $latex b_t = \beta (s_t – s_{t-1}) + (1- \beta) b_{t-1}$
- Seasonality: $latex c_t = \gamma \frac{x_t}{s_t} + (1 – \gamma) c_{t-L}$
- Smoothed values: $latex F_{t + m} = (s_t + mb_t) c_t$

Function to calculate level, trend, seasonality and smoothed values

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 | def calculate_smoothed_values(data, cycle, level_init_value, trend_init_value, season_list, alpha, beta, gamma): smoothed_value_list = [] previous_level = level_init_value previous_trend = trend_init_value previous_seasonality = season_list[-1] new_season_list = [0] * len(season_list) for i in xrange(len(data) - cycle): if i < len(season_list): new_season_list[i] = season_list[i] level = (alpha * data[i + cycle] / new_season_list[i]) + ((1 - alpha)*(previous_level + previous_trend)) trend = (beta * (level - previous_level)) + ((1 - beta) * previous_trend) seasonality = (gamma * data[i + cycle] / level) + ((1 - gamma) * previous_seasonality) smoothed_value = (level + trend) * seasonality previous_level = level previous_trend = trend previous_seasonality = seasonality new_season_list.append(seasonality) smoothed_value_list.append(smoothed_value) return smoothed_value_list |

Then I call the above function in an recursive function with different alpha, beta and gamma combos and calculate RMSE for each case. (rmse function that’s being used below is a simple function that calculates root mean squared errors)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 | def least_rmse_combo(data, cycle, level_init_value, trend_init_value, season_init_list): ''' Try different permutation of alpha, beta and gamma to get smoothed data of Level, Trend and Seasonality with the least RMSE ''' least_rmse = float("inf") lrc = [0, 0, 0] for a in xrange(101): alpha = a/100.0 for b in xrange(101): beta = b / 100.0 for g in xrange(101): gamma = g / 100.0 smoothed_values = calculate_smoothed_values(data, cycle, level_init_value, trend_init_value, season_init_list, alpha, beta, gamma) rmse_value = rmse(smoothed_values, data[cycle:]) if rmse_value < least_rmse: lrc = [alpha, beta, gamma] least_rmse = rmse_value return lrc |

Exponential smoothing. (n.d.). In Wikipedia. Retrieved Mar 30, 2013, from

http://en.wikipedia.org/wiki/Exponential_smoothing