Holt Winters Python Implementation

gistfile1.txt

def holtwinters(y, alpha, beta, gamma, c, debug=True):
    """
    y - time series data.
    alpha , beta, gamma - exponential smoothing coefficients 
                                      for level, trend, seasonal components.
    c -  extrapolated future data points.
          4 quarterly
          7 weekly.
          12 monthly
 
 
    The length of y must be a an integer multiple  (> 2) of c.
    """
    #Compute initial b and intercept using the first two complete c periods.
    ylen =len(y)
    if ylen % c !=0:
        return None
    fc =float(c)
    ybar2 =sum([y[i] for i in range(c, 2 * c)])/ fc
    ybar1 =sum([y[i] for i in range(c)]) / fc
    b0 =(ybar2 - ybar1) / fc
    if debug: print "b0 = ", b0
 
    #Compute for the level estimate a0 using b0 above.
    tbar  =sum(i for i in range(1, c+1)) / fc
    print tbar
    a0 =ybar1  - b0 * tbar
    if debug: print "a0 = ", a0
 
    #Compute for initial indices
    I =[y[i] / (a0 + (i+1) * b0) for i in range(0, ylen)]
    if debug: print "Initial indices = ", I
 
    S=[0] * (ylen+ c)
    for i in range(c):
        S[i] =(I[i] + I[i+c]) / 2.0
 
    #Normalize so S[i] for i in [0, c)  will add to c.
    tS =c / sum([S[i] for i in range(c)])
    for i in range(c):
        S[i] *=tS
        if debug: print "S[",i,"]=", S[i]
 
    # Holt - winters proper ...
    if debug: print "Use Holt Winters formulae"
    F =[0] * (ylen+ c)   
 
    At =a0
    Bt =b0
    for i in range(ylen):
        Atm1 =At
        Btm1 =Bt
        At =alpha * y[i] / S[i] + (1.0-alpha) * (Atm1 + Btm1)
        Bt =beta * (At - Atm1) + (1- beta) * Btm1
        S[i+c] =gamma * y[i] / At + (1.0 - gamma) * S[i]
        F[i]=(a0 + b0 * (i+1)) * S[i]        
        print "i=", i+1, "y=", y[i], "S=", S[i], "Atm1=", Atm1, "Btm1=",Btm1, "At=", At, "Bt=", Bt, "S[i+c]=", S[i+c], "F=", F[i]
        print i,y[i],  F[i]
    #Forecast for next c periods:
    for m in range(c):
        print "forecast:", (At + Bt* (m+1))* S[ylen + m]
 
# the time-series data.
y =[146, 96, 59, 133, 192, 127, 79, 186, 272, 155, 98, 219]
 
holtwinters(y, 0.2, 0.1, 0.05, 4)