Documentation of PDO_looks


Global Index (short | long) | Local contents | Local Index (short | long)


Help text

  Look at partial autocorrelation function

Cross-Reference Information

This script calls

Listing of script PDO_looks


clear
cd /home/disk/tao/dvimont/Data/Indices
load pdo.asc;
yr = pdo(:,1); pdo = pdo(:,2:13);
pdo2 = pdo'; pdo2 = pdo2(:);
pdo2(pdo2 == -9999) = NaN;
pdo2 = annave(pdo2);
yr = 1900:1/12:2002; yr = yr(1:1224);
pdo3 = repmat(NaN, [100, 5]);
for i = 1:100;
  ind = 12*(i-1);
  pdo3(i,1) = mean2(pdo2(ind+[11:22]));
  pdo3(i,2) = mean2(pdo2(ind+[12:14]));
  pdo3(i,3) = mean2(pdo2(ind+[15:17]));
  pdo3(i,4) = mean2(pdo2(ind+[18:20]));
  pdo3(i,5) = mean2(pdo2(ind+[21:23]));  
end
yr3 = 1901:2000;
figure_tall(1); clf;
subplot(4,1,1);
  p1 = plot(yr, pdo2, 'k');
  axis([1900 2005 -3 3]);
subplot(4,1,2);
  p2 = plot(yr3, pdo3, 'k');
  axis([1900 2005 -3 3]);
  set(p2(1), 'linewidth', 2);
  set(p2(2), 'color', 'b');
  set(p2(4), 'color', 'r');
a = NaN*ones(5, 1);
for i = 1:5;
  a(i) = corr(pdo3(:,i), pdo3(:,i), 1);
end
z = randn(500, 6);
b = zeros(500, 5);
for j = 1:5;
for i = 2:500;
  b(i,j) = a(1)*b(i-1,j)+0.3*z(i,j)+0.7*z(i,6);
end
end
for blah = 1:4;
tim = 151:250;
ctann = b(tim,blah);

lags = 1:20; nlag = length(lags);
a = NaN*ones(size(lags));
for i = 1:nlag;
  a(i) = corr(ctann, ctann, lags(i));
end

alpha = zeros(nlag);
alpha(1,1) = a(1);
for p = 2:(nlag);
  k = 1:(p-1);
  ind2 = p-k;
  alpha(p,p) = ( a(p) - sum(alpha(p-1,k).*a(ind2)) ) / ...
        ( 1 - sum(alpha(p-1,ind2).*a(ind2)) );
  for k = 1:(p-1);
    alpha(p,k) = alpha(p-1,k) - alpha(p,p)*alpha(p-1,p-k);
  end
end

ntim = length(ctann);
ci = 1.96*sqrt(1/ntim);

%  Plot result
subplot(4,2,4+blah);
h1 = bar([lags], [diag(alpha)]);
set(h1, 'facecolor', 0.7*[1 1 1]);
axis([-0.5 20.5 -0.3 0.8]);
hold on;
  h2 = hline([-1*ci 1*ci], '--k');
hold off;
  set(gca, 'XTick', 0:2:10, 'YTick', [-1:.25:1]);
  set(gca, 'fontsize', 8);
%  xl = xlabel('lag (years)');
%  set(xl, 'fontsize', 9);
end