好像网上有
/* ** 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;
扫码加好友，拉您进群
全部版块
我的主页
收藏
计量都是挺难的.
我刚刚开始起步。。。
