%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% this script produces autocorrelation (AC) plots
% to assess the efficiency of the posterior simulator
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Using posterior draws based on efficient (2 block) algorithm
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

load c:\klaus\AAEC6564\mlab\worksp\mod2application betamat sig2mat;

% Example: normal regression model with independent priors.
% Here we have 3 slope coefficients and 1 variance, so let's create a
% figure with a 2 x 2 set of subplots

figure(1);
% focus on 2nd equation coeff's and elements of E
subplot(2,2,1);
d=50; %we're onlyn interested into lags up to 50
int=betamat(end-2,:)';
m=length(int);
j=1;
for j=1:d
    int1=int(1:m-j,1);
    int2=int(j+1:m,1);
    int3=corrcoef(int1,int2);
    lagcorr(j,1)=int3(1,2);
end

xgrid=linspace(1,d,d)';
bar(xgrid,lagcorr,'y');
xlabel('Lag education');
ylabel('Correlation');

subplot(2,2,2);
int=betamat(end-1,:)';
m=length(int);
j=1;
for j=1:d
    int1=int(1:m-j,1);
    int2=int(j+1:m,1);
    int3=corrcoef(int1,int2);
    lagcorr(j,1)=int3(1,2);
end

xgrid=linspace(1,d,d)';
bar(xgrid,lagcorr,'y');
xlabel('Lag KL6');
ylabel('Correlation');

subplot(2,2,3);
int=betamat(end,:)';
m=length(int);
j=1;
for j=1:d
    int1=int(1:m-j,1);
    int2=int(j+1:m,1);
    int3=corrcoef(int1,int2);
    lagcorr(j,1)=int3(1,2);
end

xgrid=linspace(1,d,d)';
bar(xgrid,lagcorr,'y');
xlabel('Lag K618');
ylabel('Correlation');

subplot(2,2,4);
int=sig2mat';
m=length(int);
j=1;
for j=1:d
    int1=int(1:m-j,1);
    int2=int(j+1:m,1);
    int3=corrcoef(int1,int2);
    lagcorr(j,1)=int3(1,2);
end

xgrid=linspace(1,d,d)';
bar(xgrid,lagcorr,'y');
xlabel('Lag sig2');
ylabel('Correlation');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Using posterior draws based on inefficient (5 block) algorithm
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

load c:\klaus\AAEC6564\mlab\worksp\mod2blocking beta1mat beta2mat beta3mat sig2mat;

figure(2);
% focus on 2nd equation coeff's and elements of E
subplot(2,2,1);
d=50; %we're onlyn interested into lags up to 50
int=beta2mat(end,:)';
m=length(int);
j=1;
for j=1:d
    int1=int(1:m-j,1);
    int2=int(j+1:m,1);
    int3=corrcoef(int1,int2);
    lagcorr(j,1)=int3(1,2);
end

xgrid=linspace(1,d,d)';
bar(xgrid,lagcorr,'y');
xlabel('Lag education');
ylabel('Correlation');

subplot(2,2,2);
int=beta3mat(end-1,:)';
m=length(int);
j=1;
for j=1:d
    int1=int(1:m-j,1);
    int2=int(j+1:m,1);
    int3=corrcoef(int1,int2);
    lagcorr(j,1)=int3(1,2);
end

xgrid=linspace(1,d,d)';
bar(xgrid,lagcorr,'y');
xlabel('Lag KL6');
ylabel('Correlation');

subplot(2,2,3);
int=beta3mat(end,:)';
m=length(int);
j=1;
for j=1:d
    int1=int(1:m-j,1);
    int2=int(j+1:m,1);
    int3=corrcoef(int1,int2);
    lagcorr(j,1)=int3(1,2);
end

xgrid=linspace(1,d,d)';
bar(xgrid,lagcorr,'y');
xlabel('Lag K618');
ylabel('Correlation');

subplot(2,2,4);
int=sig2mat';
m=length(int);
j=1;
for j=1:d
    int1=int(1:m-j,1);
    int2=int(j+1:m,1);
    int3=corrcoef(int1,int2);
    lagcorr(j,1)=int3(1,2);
end

xgrid=linspace(1,d,d)';
bar(xgrid,lagcorr,'y');
xlabel('Lag sig2');
ylabel('Correlation');

% Some really bad cases:
%%%%%%%%%%%%%%%%%%%%%%%%

figure(3);
% focus on 2nd equation coeff's and elements of E
subplot(2,1,1);
d=50; %we're onlyn interested into lags up to 50
int=beta1mat(end,:)';
m=length(int);
j=1;
for j=1:d
    int1=int(1:m-j,1);
    int2=int(j+1:m,1);
    int3=corrcoef(int1,int2);
    lagcorr(j,1)=int3(1,2);
end

xgrid=linspace(1,d,d)';
bar(xgrid,lagcorr,'y');
xlabel('Lag age');
ylabel('Correlation');

subplot(2,1,2);
d=50; %we're onlyn interested into lags up to 50
int=beta2mat(end-1,:)';
m=length(int);
j=1;
for j=1:d
    int1=int(1:m-j,1);
    int2=int(j+1:m,1);
    int3=corrcoef(int1,int2);
    lagcorr(j,1)=int3(1,2);
end

xgrid=linspace(1,d,d)';
bar(xgrid,lagcorr,'y');
xlabel('Lag age2');
ylabel('Correlation');
