/* Henderson's 23-term Moving Average */
new; cls;
let data[50,1]=
      60
      61
      60
      62
      60
      63
      64
      65
      69
      70
      69
      66
      64
      64
      62
      64
      65
      55
      53
      57
      60
      58
      63
      63
      63
      68
      69
      69
      68
      69
      70
      67
      66
      67
      68
      69
      69
      77
      80
      80
      78
      76
      76
      72
      71
      70
      73
      71
      78
      80
;
data=rev(data);               /* original data in reverse order */
t=seqa(1,1,rows(data));

library pgraph;
graphset;
_pltype=6;
_plegctl=1;
_plegstr="actual\000Henderson 23";
title("Henderson 23");
xy(t,data~henderson23(data));




proc henderson23(x);
   local w,t,k,i,m,w1,w2,w3,w4,w5,w6,w7,w8,w9,w10,w11;
/* Symmetric part in the middle */
   w={0.148,0.138,0.122,0.097,0.068,0.039,0.013,-0.005,-0.015,-0.016,-0.011,-0.004};
   t=rows(x);
   k=rows(w)-1;
   w=rev(w[2:k+1])|w;
   m=zeros(rows(x),cols(x));
   i=k+1;
   do while i<=t-k;
      m[i,.]=w'x[i-k:i+k,.];
      i=i+1;
   endo;
/* Asymmetric part at end points */
   w1={-0.077,-0.064,-0.049,-0.028,0.002,0.039,0.084,0.133,0.1819999999999999,0.227,0.2629999999999999,0.2879999999999999}; 
   w2={-0.04599999999999999,-0.041,-0.035,-0.02399999999999999,-0.004,0.025,0.06099999999999999,0.101,0.1409999999999999,0.1759999999999999,0.203,0.2189999999999999,0.224};
   w3={-0.02199999999999999,-0.025,-0.025,-0.019,-0.005,0.01799999999999999,0.049,0.08200000000000001,0.116,0.1459999999999999,0.166,0.177,0.1759999999999999,0.166};
   w4={-0.008,-0.014,-0.01799999999999999,-0.015,0.004,0.015,0.042,0.07299999999999999,0.103,0.129,0.1469999999999999,0.154,0.15,0.134,0.112};
   w5={-0.001,-0.008,-0.013,-0.01199999999999999,-0.003,0.015,0.04009999999999999,0.06809999999999999,0.0981,0.1211,0.1370999999999999,0.1421,0.1350999999999999,0.1191,0.09510000000000001,0.0661};
   w6={0.003,-0.006,-0.01099999999999999,-0.01099999999999999,-0.002,0.015,0.039,0.067,0.095,0.1189999999999999,0.134,0.139,0.131,0.114,0.08799999999999998,0.05899999999999999,0.027};
   w7={0.002,-0.006,-0.01199999999999999,-0.01099999999999999,-0.003,0.015,0.039,0.068,0.09599999999999999,0.1179999999999999,0.134,0.138,0.132,0.114,0.08899999999999998,0.05899999999999999,0.027,0.001};
   w8={0.001,-0.007,-0.013,-0.01099999999999999,-0.003,0.015,0.039,0.068,0.09599999999999999,0.1199999999999999,0.135,0.14,0.133,0.116,0.09,0.06,0.031,0.005,-0.015};
   w9={-0.002,-0.007,-0.013,-0.013,-0.003,0.014,0.039,0.068,0.09699999999999999,0.1199999999999999,0.137,0.14,0.136,0.1179999999999999,0.09399999999999998,0.064,0.034,0.008,-0.01,-0.021};
   w10={-0.003,-0.01,-0.015,-0.014,-0.005,0.014,0.04,0.069,0.09699999999999999,0.1219999999999999,0.138,0.1429999999999999,0.137,0.1199999999999999,0.095,0.067,0.037,0.01099999999999999,-0.007,-0.017,-0.019};
   w11={-0.004,-0.01099999999999999,-0.016,-0.015,-0.005,0.013,0.039,0.068,0.09699999999999999,0.1219999999999999,0.138,0.1439999999999999,0.138,0.1219999999999999,0.09699999999999999,0.068,0.039,0.013,-0.005,-0.015,-0.016,-0.01099999999999999};
   m[1,.]=rev(w1)'x[1:k+1,.]; m[t,.]=w1'x[t-k:t,.];
   m[2,.]=rev(w2)'x[1:k+2,.]; m[t-1,.]=w2'x[t-k-1:t,.];
   m[3,.]=rev(w3)'x[1:k+3,.]; m[t-2,.]=w3'x[t-k-2:t,.];
   m[4,.]=rev(w4)'x[1:k+4,.]; m[t-3,.]=w4'x[t-k-3:t,.];
   m[5,.]=rev(w5)'x[1:k+5,.]; m[t-4,.]=w5'x[t-k-4:t,.];
   m[6,.]=rev(w6)'x[1:k+6,.]; m[t-5,.]=w6'x[t-k-5:t,.];
   m[7,.]=rev(w7)'x[1:k+7,.]; m[t-6,.]=w7'x[t-k-6:t,.];
   m[8,.]=rev(w8)'x[1:k+8,.]; m[t-7,.]=w8'x[t-k-7:t,.];
   m[9,.]=rev(w9)'x[1:k+9,.]; m[t-8,.]=w9'x[t-k-8:t,.];
   m[10,.]=rev(w10)'x[1:k+10,.]; m[t-9,.]=w10'x[t-k-9:t,.];
   m[11,.]=rev(w11)'x[1:k+11,.]; m[t-10,.]=w11'x[t-k-10:t,.];
   retp(m);
endp;