/*
** This procedure performs Tsay's (1998) multivariate
** TAR specification test
** and was written by Ming Chien Lo on August 27, 1999.
**
** Last updated: August 28, 1999.
**
** Input: y T x k endogenous data series
** x T x q exogenous data series
** z T x 1 threshold data series
** p order of the autoregressive process for y
** q order of the autoregressive process for x
** d lag index for threshold variable z{t-d}
** m0 starting value for the resursive regression
**
** Output: Cd resulting LR-statistics
**
*/
proc tsaymul(y,x,z,p,q,d,m0);
local h, t, k, v, nr, ztd, yt, ytl, i, yl, kk, xtl, vv, xt, gs, y_, x_, Vm, phim, phim1, etham1,
ym1, xm1, Km1, x1, psi, w, Cd, S0, S1, xl, em1, invXX;
h=maxc(p|q|d); @ starting index of the arranged data set @
t=rows(y); @ total no. of observations @
k=cols(y);
v=cols(x);
nr=p*k+q*v+1; @ no. of RHS regressors @
ztd=z[h+1-d:t-d]; @ potential threshold @
yt=y[h+1:t,.]; @ truncated y series @
ytl=zeros(rows(y),p*k); @ form and do the same thing for y{t-p} matrix @
i=1;
do until i>k;
yl=shiftr(y[.,i]',seqa(1,1,p),0);
kk=(i-1)*p+1;
ytl[.,kk:kk+p-1]=zeros(p,p)|trimr(yl',p,0);
i=i+1;
endo;
ytl=ytl[h+1:t,.];
xt=ytl;
if x NE 0;
xtl=zeros(rows(x),q*v); @ form and do the same thing for x{t-q} matrix @
i=1;
do until i>v;
xl=shiftr(x[.,i]',seqa(1,1,q),0);
vv=(i-1)*q+1;
xtl[.,vv:vv+q-1]=zeros(q,q)|trimr(xl',q,0);
i=i+1;
endo;
xtl=xtl[h+1:t,.];
xt=ytl~xtl; @ truncated independent variable matrix @
endif;
gs=sortind(ztd); @ index for sorted z{t-d} @
y=yt[gs,.]; @ arranged Y according to sorted y{t-d} @
x=xt[gs,.]; @ arranged X according to sorted y{t-d} @
x=ones(t-h,1)~x;
@ derive first phi_m for the first m cases according to m0 @
etham1=zeros(k,t-h-m0);
y_=y[1:m0,.];
x_=x[1:m0,.];
Vm=invpd(moment(x_,0));
phim=Vm*x_'y_;
@ recursive regression @
i=1;
do until i>t-h-m0;
ym1=y[m0+i,.]'; @ yt(m+1)+d @
xm1=x[m0+i,.]'; @ xt(m+1)+d @
@ compute predictive error @
em1=ym1-phim'xm1;
etham1[.,i]=em1/sqrt(1+xm1'Vm*xm1);
@ recursion @
Km1=Vm*xm1/(1+xm1'*Vm*xm1);
phim1=phim+Km1*em1';
x_=x[1:m0+i,.];
Vm=invpd(moment(x_,0)); @ Vm for next loop @
phim=phim1;
i=i+1;
endo;
@ calculate F-statistics @
x1=x[m0+1:t-h,.];
invXX=invpd(moment(x1,0));
psi=invXX*x1'etham1';
w=etham1'-x1*psi;
w=w';
S0=(1/(t-h-m0))*(etham1*etham1');
S1=(1/(t-h-m0))*(w*w');
Cd=(t-h-m0-nr)*(ln(det(S0))-ln(det(S1)));
retp(Cd);
endp;
扫码加好友，拉您进群
全部版块
我的主页
收藏
