# Helping in understaing new Pesaran and Timmermann test (and R code validation)

I am currently working with Pesaran & Timmermann test version from year 2009. Since I could not find any R package that contains function to calculate it (rugarch has 1992 version in DACTest). I recently posted thread in cross validated with similar question where kind user mlofton advised me to look into this list. Since cross-posting is not acceptable practice after over week of hard research I managed to develop R code for newPT_test. It seems to be correct, but since I had lots of doubts maybe someone will clarify my concerns.

Here is the formula I used (Pönkä 2017):

$$PT=(T-1)(S^{-1}_{yy,w}S_{y\hat{y},w}S^{-1}_{\hat{y}\hat{y},w}S_{\hat{y}y,w}) \sim {\chi}^2_1,$$

$$S_{yy,w}=(T-1)^{-1}Y'M_wY,$$ $$S_{\hat{y}\hat{y},w}=(T-1)^{-1}\hat{Y}'M_w\hat{Y},$$ $$S_{y\hat{y},w}=(T-1)^{-1}\hat{Y}'M_w\hat{Y},$$ $$S_{y\hat{y},w}=(T-1)^{-1}Y'M_w\hat{Y},$$ $$M_w=I_{T-1}-W(W'W)^{-1}W',$$ $$W=(\tau_{T-1},Y_{-1},\hat{Y}_{-1}),$$

"...and $$Y=(y_2,...,y_T)'$$",$$\hat{Y}=(\hat{y}_2,...,\hat{y}_T)$$,$$Y_{-1}=(y_1,...,y_{T-1})$$,$$\hat{Y}=(\hat{y}_1,...,\hat{y}_{T-1})$$ and $$\tau_t$$ is as $$(T-1)\times1$$ vector of ones..." From other source I learned that $$I_t$$ is identity matrix.

My main issue was notation. I understand that $$y_t$$ represents ones and zeroes depending on whether change between actual value and its last observation was positive ($$1$$) or negative ($$0$$)? Same for forecast $$\hat{y}_t$$. I am not sure how to understand vectors: $$Y=(y_2,...,y_t)'$$ and $$Y'$$ in $$S$$ elements. I know that $$'$$ at the end of $$(y_2,...,y_t)'$$ means it is column vector, so $$Y'$$ indicates that now its row vector now?. Also it starts from $$y_2$$ which means I remove first observation from binary series while in $$Y_{-1}=(y_1,...,y_{T-1})$$ I remove last? That makes $$\tau_{t-1}$$ a $$(T-1)\times3$$ matrix right?

Here is R code:

nwPT_test=function(actual,forecast){

yt=actual #assign actual to yt to make code shorter
xt=forecast #assign forecast to xt...

delta_yt=as.matrix(cbind(ifelse(yt-lag(yt)>0,1,0)[-1])) #calc change of yt (delta)
delta_xt=as.matrix(cbind(ifelse(xt-lag(xt)>0,1,0)[-1])) #calc change of xt(delta)
nT=length(delta_yt) #number of Time periods
Yt=cbind(delta_yt[-1]) #Yt=(y2,...,yT)
Xt=cbind(delta_xt[-1]) #Xt=(x2,...,xT)
Yt2=as.vector(rbind(delta_yt[-nT])) #Yt2=(y1,...,yT-1)
Xt2=as.vector(rbind(delta_xt[-nT])) #Xt2=(x1,...,xT-1)
teta=rep(1,nT-1) #T-1 vector of ones
I=diag(nT-1) #Identity matrix
W=cbind(teta,Yt2,Xt2) #W matrix as in formula

Mw=I-(W%*%((t(W)%*%W)^(-1))%*%t(W))#calcualting Mw as in formula

#calculating elements S as in formula
Syy.w=((nT-1)^(-1))*t(Yt)%*%Mw%*%Yt
Sxx.w=((nT-1)^(-1))*t(Xt)%*%Mw%*%Xt
Sxy.w=((nT-1)^(-1))*t(Xt)%*%Mw%*%Yt
Syx.w=((nT-1)^(-1))*t(Yt)%*%Mw%*%Xt

PT=(nT-1)*(Syy.w^(-1)*Syx.w*Sxx.w^(-1)*Sxy.w)#finally calculating PT

p.value=1-pchisq(PT,df=1)#calculating p-value
#some code to make it looks nicer
summary=c(PT,p.value)
names(summary)=c("PT statistic","p.value")
summary
}


I managed to solve most issues from my previous thread. Here is the link:: https://stats.stackexchange.com/questions/372776/confusions-about-pesaran-timmermann-test-2009-version?noredirect=1#comment701286_372776 It might help understand the way I approached it. So far my code seems to give correct results.

Please answer only questions from current thread. Hope I updated my question enough to not be treated as double post... And thanks again to mlofton for informing me about this great SE subsite.