%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Example: Weekly Themostat Sales
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%




%%%% Import Weekly Themostat Sales data to matlab

load Thermostat_Sales.mat;

%%%%  Compute the sample size of data

n=length(Thermostat_Sales);

%%%%  change data to z_t, get the ordinary least square estimate for beta

z=Thermostat_Sales';

X=[ones(n,1) (1:n)'];

beta=inv(X'*X)*X'*z;

%%% computing residuals for model
 
Error=z-X*beta;

%%%%  Compute the auto-correlation of residuals from lag=1 to 6 and show it
%%%%  on the screen 

for i=1:6
    r(i)=auto_corr(Error,i);
end

[1:6; r]

%%%% Compute SSE from alpha=0.02 to 0.3

alpha=0.02:0.01:0.3;

m=length(alpha);

sse=zeros(m,1);

r_error=zeros(m,6);

error_mean=zeros(m,1);

error_std=zeros(m,1);

% sse: SSE(alpha), 
% t_error: auto_correlation of one step forecast error,
% error_mean: mean of forecast error. 
% error_std: standard deviation of the estimate of the mean of forecast error

for i=1:m
   
    [sse(i), r_error(i,:), error_mean(i), error_std(i)]=SSE1(z,alpha(i));
   
end

% show the value of SSE on screen according to the order of alpha 


[alpha' sse]

%%% Show the autocorrelation value, mean and standard deviation of residuals
%%% on screen when alpha=0.14

[1:6; r_error(13,:)]

[error_mean(13) error_std(13)]

%%% Plot a figure to show SSE with alpha, and select alpha which minimize
%%% SSE. 

plot(alpha,sse);

title('Plot of SSE(\alpha), for double exponential smoothing---thermostat sales');

xlabel('\alpha');

ylabel('SSE(\alpha)');