我想对AR(2)的phi2, 也就是y_{t-2}的系数构造一个95%的likelihood ratio confidence interval，最大似然估计采用exact MLE，具体用kalman filter来做
model： (y_t - mu)=phi1*(y_{t-1} - mu) + phi2*(y{t-2} - mu) + epsilon_t, epsilon_t ~ N(0, var)，所以，参数theta=( mu, var, phi1, phi2)
因为采用exact MLE, 所以必须限制| phi1 |<2，我想固定住phi1, 然后用 qnewton来求likelihood的最大值，并且求出相应的phi2的估计，如果这个最大值超出-341.368，那么phi2就在95%置信区间内，但是程序提示error G0003 : Indexing a matrix as a vector，这个应该是 Kalman Filter 里面用列向量 phi1 来定义矩阵F出了问题，但我不知道为什么会有问题。。。
我已经花了2个小时来修改程序，可是还是搞不定，希望板上有大侠能帮忙解答哈。。。感激不尽！
程序提示错误如下(刚开始没看到第一行的错误提示，汗。。。） ：
H:\coursework\S1_cousework\Applied Ecometrix\hw3\Q2\ii\LR CI(58) :
error G0003 : Indexing a matrix as a vector
Currently active call: lik_fcn [59] H:\coursework\S1_cousework\Applied Ecometrix\hw3\Q2\ii\LR CI
Stack trace:
lik_fcn called from d:\src\qnewton.src, line 178
_bfgs called from d:\src\qnewton.src, line 98
qnewton called from H:\coursework\S1_cousework\Applied Ecometrix\hw3\Q2\ii\LR CI, line 27


程序如下 （data见附件）：
new;
cls;
library pgraph;
format /m1 /rd 9,3;

load data[256,1]=gdp.txt;

yy=100*ln(data);


dy=yy[2:rows(yy)]-yy[1:rows(yy)-1];
T=rows(dy); @sample size@
n=1;
phi1=seqa(-2+0.00001,0.01,399);
s=rows(phi1);
phi2_mat={};
lik_mat={};

do until n > s;
@MAXIMUM LIKELIHOOD ESTIMATION@
prmtr_in=zeros(3,1); @starting values for parameters@
prmtr_in[1]=meanc(dy);@Use the sample mean as starting value for mu@
prmtr_in[2]=ln((stdc(dy))^2); @Use sample variance of dy as starting value for var_e@

alpha=0.05; @level at which to set Type I error@
{xout,fout,gout,cout}=qnewton(&lik_fcn,prmtr_in);

prm_fnl=trans(xout); @final constrained point estimates@
lik=-fout; @final likelihood value@

if lik ge -341.368;
phi2_mat=phi2_mat|prm_fnl[3];
lik_mat=lik_mat|lik;
endif;

n=n+1;
endo;

end;

@========================================================================@
@========================================================================@
proc LIK_FCN(PRMTR1); /* Kalman Filter*/

local prmtr,mu,var,phi1,phi2,theta1,theta2,f,q,h,x_ll,p_ll,
llikv,j,x_tl,p_tl,eta_tl,f_tl,x_tt,p_tt,vP_ll,A,z,R,y;

prmtr=trans(prmtr1); @constrain parameters@

@PARAMETERS, AR(2), phi1 fixed@
mu=prmtr[1]; @unconditional mean@
var=prmtr[2]; @standard deviation of error term@
phi2=prmtr[3]; @ar coefficient@


@SET UP STATE SPACE FORM@

F= phi1[n]~phi2|
1~0;

Q= var~0|
0~0;

H= 1~0;

z=ones(T,1); @Let "z" denote a vector of ones to pick up a constant term in the observation equation@
A= mu; @coefficient on the constant term@

R=0; @no error in the observation equation@

y=dy; @Let "y" denote the name of a generic observation variable or vector, even though it is real GDP growth in this case@


@INITIAL STATE ESTIMATES@
X_LL= zeros(rows(F),1); @unconditional expectation of state vector@
vP_LL = inv( (eye(rows(F)^2)-F.*.F) )*vec(Q);
P_LL=reshape(vP_LL,rows(F),rows(F));

llikv = 0;
j=1;
do until j>T;

@KALMAN FILTER@
X_TL = F * X_LL;
P_TL = F* P_LL * F' + Q;

ETA_TL = y[j] - H * X_TL - A*z[j];
F_TL = H * P_TL * H' + R; @lambda_t@

X_TT = X_TL + P_TL*H'*INV(F_TL) * ETA_TL;
P_TT = P_TL - P_TL * H'*INV(F_TL) * H * P_TL;


@LIKELIHOOD FUNCTION@
llikv=llikv - 0.5*(ln(2*PI*F_TL) + (ETA_TL^2)/F_TL);

@UPDATE for NEXT ITERATION@
X_LL=X_TT;
P_LL=P_TT;

j=j+1;
endo;

retp(-llikv); @return negative of likelihood because qnewton is a minimization routine@
endp;

@========================================================================@


@PROCEDURE TO CONSTRAIN PARAMETERS@
proc trans(c0);

local c1,aaa,ccc;

c1 = c0;
c1[2]=exp(c0[2]); @positive variance@
aaa=phi1[n]./2;
ccc=(1-abs(aaa))*c0[3]./(1+abs(c0[3]))+abs(aaa)-aaa^2;

c1[3]=-1* (aaa^2+ccc);
retp(c1);
endp;