#time series, generate ar5 and add intercept, nonstat-stat data, choice of lambda, c, urbic,main program table on setup 2 design 3 p=5 t=100 ite=1000 rej=matrix(0,ite,1) for (jjj in 1:ite) { #zy=numeric(t+p) zy=matrix(0,t+p,1) e=rnorm(t+p) zy[1]=e[1] zy[2]=e[2] zy[3]=e[3] zy[4]=e[4] zy[5]=e[5] for (i in (p+1):(t+p)){zy[i]=0.15*zy[i-1]+0.65*zy[i-2]+0*zy[i-3]+0*zy[i-4]+0*zy[i-5]+e[i]} ones=rep(1,t+p) int=0.5*ones int=as.matrix(int) zy=as.matrix(zy) tt=0.1*seq(1,t+p,1) tt=as.matrix(tt) y=zy+int+tt y=as.matrix(y) #data is generated #dy=matrix(0,t+p-1,1) dy=y[(p+1):(t+p),]-y[(p):(t+p-1),] yl=y[(p):(t+p-1),] dyl=matrix(0,t,(p-1)) dyl=cbind(y[(p):(t+p-1),]-y[(p-1):(t+p-2),],y[(p-1):(t+p-2),]-y[(p-2):(t+p-3),],y[(p-2):(t+p-3),]-y[(p-3):(t+p-4),],y[(p-3):(t+p-4),]-y[(p-4):(t+p-5),]) ones=as.matrix(ones) x=cbind(yl,rep(1,100),seq(1,100,1),dyl) x=as.matrix(x) rhoi=-0.025 ini=0.48 tti=0.02 zeta1i=-0.62 zeta2i=-0.02 zeta3i=-0.62 zeta4i=-0.02 zetai=rbind(zeta1i,zeta2i,zeta3i,zeta4i) zetai=as.matrix(zetai) g1=1/4 g2=1/2 g3=1/3 g4=1/2 lam1=rbind(0.1,1,10,100) lam1=as.matrix(lam1) c1=seq(15,35,5) c1=as.matrix(c1) noc=5 #for penalty 4 s1=colSums(abs(zetai)^g2) s22=g2*colSums(abs(zetai)^(g2-2)) #done nb1=matrix(0,4,noc) nb=matrix(0,p+1,noc) #4 is 4 lambda choice bic=matrix(0,4,noc) for (ii in 1:4) { for (ji in 1:(noc)) { #now penalty terms and function brid=function(bet){ ob=t((dy-x%*%bet))%*%(dy-x%*%bet) pen1=(lam1[ii]*abs(rhoi)^g1)+(0.5*(bet[1]-rhoi)^2)*(lam1[ii]*g1*abs(rhoi)^(g1-2)) pen2=(lam1[ii]*abs(ini)^g2)+(0.5*(bet[2]-ini)^2)*(lam1[ii]*g2*abs(ini)^(g2-2)) pen3=(lam1[ii]*abs(tti)^g3)+(0.5*(bet[3]-tti)^2)*(lam1[ii]*g3*abs(tti)^(g3-2)) pen4=lam1[ii]*s1+(0.5*t((bet[4:7]-zetai))%*%(bet[4:7]-zetai)*lam1[ii]*s22) fval=ob+pen1+pen2+pen3+pen4} bet0=c(0,0,0,0,0,0,0) brid(bet0) a=optim(bet0,brid,method="BFGS",hessian=TRUE) bhat=a\$par bhat1=t(bhat) #end function and optimization #start truncation bhat1=matrix(bhat1,7,1) if (bhat1[1,]>-c1[ji,]/t) {nb1[ii,ji]=0} else {nb1[ii,ji]=bhat1[1,]} if (abs(bhat1[2,])