% Time Series. Econ232C. Spring 2006. Homework 1
clear all;
%close all;
clc;
%define parameters
N=100;  %cross-section size (50,100,500)
T=500;  %time horizon
ro=0.7; %AR(1) parameter
L=10;   %Lag number
sigmax=10; %error variance
%construct data
Xjt=sigmax*randn(N,T);  %generate X
err=randn(1,T);
Zt=zeros(1,T);
for i=2:T       %generate Z as AR(1)
    Zt(i)=ro*Zt(i-1)+err(i);
end
Yjt=Xjt+Kron(Zt,ones(N,1)); %construct Y
%compute ACF
h=Yjt(1,:);
gamma0=(h*h')/T;
for j=1:L
    gamma(j)=( h(1,(j+1):T) * (h(1,1:(T-j)))')/(T-j);
    ster(j)=2/sqrt(T);
end
ACF=gamma/gamma0;   
figure(1);
subplot(1,2,1);
stem_hadles=stem(1:L,ACF,'fill');
hold on;
plot_hadles=plot(1:L,ster,'r',1:L,-ster,'r');
title('ACF')
hold off;

%compute PACF
h=Yjt(1,:);
for j=1:L
    Y=h(1,(j+1):T)';
    X=zeros((T-j),j);
    for k=1:j
        X(:,k)=h(1,(j+1-k):(T-k))';
    end
    beta=inv(X'*X)*(X'*Y);
    alpha(j)=beta(j);
    ster(j)=2/sqrt(T);
end
PACF=alpha;
subplot(1,2,2);
stem_hadles=stem(1:L,PACF,'fill');
hold on;
plot_hadles=plot(1:L,ster,'r',1:L,-ster,'r');
title('PACF')
hold off;

SY=ones(1,N)*Yjt;
%compute ACF
h=SY;
gamma0=(h*h')/T;
for j=1:L
    gamma(j)=( h(1,(j+1):T) * (h(1,1:(T-j)))')/(T-j);
    ster(j)=2/sqrt(T);
end
ACF=gamma/gamma0;   
figure(2);
subplot(1,2,1);
stem_hadles=stem(1:L,ACF,'fill');
hold on;
plot_hadles=plot(1:L,ster,'r',1:L,-ster,'r');
title('ACF')
hold off;

%compute PACF
h=SY(1,:);
for j=1:L
    Y=h(1,(j+1):T)';
    X=zeros((T-j),j);
    for k=1:j
        X(:,k)=h(1,(j+1-k):(T-k))';
    end
    beta=inv(X'*X)*(X'*Y);
    alpha(j)=beta(j);
    ster(j)=2/sqrt(T);
end
PACF=alpha;
subplot(1,2,2);
stem_hadles=stem(1:L,PACF,'fill');
hold on;
plot_hadles=plot(1:L,ster,'r',1:L,-ster,'r');
title('PACF')
hold off;
