options crt;  ? 6e5
in kopcke;
freq q;
smpl 52:1,86:4;
x1 = y/ce;
x2 = y(-1)/ce;  ? (a)
x2n = -x2;
?
smpl 56:1,79:4;
ar1 ie c x1(4,13,none) x2n(4,13,none) kelag;  ? (b)
rename @logl loglu;
?
ar1 ie c x1(4,13,none) kelag;  ? (c.1)
set lrtest = 2*(loglu-@logl);  ? (c.2)
cdf(chisq,df=4) lrtest;
?
?  these do loops estimate all PDLs up to 5th degree for up to 13 lags
?
smpl 56:1,86:4;
const minsbic 999999 nl 0 nd 0;
supres @coef;  ? turn off all output
do nlags=1,13;
 set maxd = nlags - pos(nlags-6);  ? maxd = min(nlags,6)
 do ndeg=1,maxd;
   ar1(silent) ie c x1(ndeg,nlags,none) x2n(ndeg,nlags,none) kelag;
   ? SBIC includes C, RHO, KELAG, free params of PDL
   set @sbic = -@logl + (3+2*ndeg)*log(@nob)/2;
   if (@sbic < minsbic); then; do;
     set minsbic = @sbic;
     set nl = nlags;
     set nd = ndeg;
   enddo;
 enddo;
enddo;
supres;
print minsbic nl nd;  ? (d) -- min of @SBIC
?  print selected model
ar1 ie c x1(nd,nl,none) x2n(nd,nl,none) kelag;