' ==============================================================================================
'                           Granger and Newbold (1974): "Spurious Regressions in Econometrics"
'
'     1.-  Generate two independent random walks  yt = muy +  yt-1 + eyt  and xt = mux +  xt-1 + ext
'     2.-  Estimate by OLS: 
'                   yt  = a + b xt + et
'     3.-  Compute the t-statistics, the DW test and the R2 and save them in the matrices "t", "DW", and "R2", respectively
'     4.-  Iterate 1000 times
'     5.- Calculate the rejection frequency of the significance test of b=0
'     6.-  Order the t-statistics, the DW and the R2
'     7.-  Compute the corresponding percentiles: palpha, pbeta, pDW, pR2
'           By columns:
'             col 1 : 5th percentile
'             col 2 : 50th percentile
'             col 3 : 95th percentile
'
' Other things to do:
' HAC (Newey-West) : equation eq1.ls(cov=hac) yt c xt  or 
' equation eq1.ls(cov=hac, bw=neweywest) yt c xt
' Regression in forst difference: equation eq1.ls(cov=hac) d(yt) c d(xt) 
'==============================================================================================

create u 1 1000

rndseed 12345678

' Generate Data and Statistics

scalar rep=1000

vector(rep) alpha = 0
vector(rep) beta = 0
vector(rep) talpha = 0
vector(rep) tbeta = 0
vector(rep) DW = 0
vector(rep) R2 = 0

scalar muy = 0
scalar mux = 0

for !i = 1 to 1000

     smpl @first @first
     series yt = nrnd
     smpl @first+1 1000
     series yt =  muy + yt(-1)+ nrnd 

     smpl @first @first
     series xt = nrnd
     smpl @first+1 1000
     series xt =  mux + xt(-1)+ nrnd 
     
     equation eq1.ls yt c xt

     smpl @first @last

     alpha(!i)=eq1.@coefs(1)
     beta(!i)=eq1.@coefs(2)
     talpha(!i)=eq1.@tstats(1)
     tbeta(!i)=eq1.@tstats(2)
     DW(!i)=eq1.@dw
     R2(!i)=eq1.@r2 

next

smpl @first @last

' Rejection frequency Ho: beta=0

scalar freqbeta0=0

for !i = 1 to 1000
     if  abs(tbeta(!i)) > 1.96 then
     freqbeta0=freqbeta0+(1/1000)
     endif
next

smpl @first @last

' Sort vectors of stats

vector alphasort=@sort(alpha,"a")
vector betasort=@sort(beta,"a")

vector talphasort=@sort(talpha,"a")
vector tbetasort=@sort(tbeta,"a")
vector DWsort=@sort(DW,"a")
vector R2asort=@sort(R2,"a")

' Percentile

matrix(1,3) palpha = 0
matrix(1,3) pbeta = 0
matrix(1,3) pDW = 0
matrix(1,3) pR2 = 0

' 5th, 50th, and 95th percentile t-alpha
palpha(1,1) = talphasort(50)
palpha(1,2) = talphasort(500)
palpha(1,3) = talphasort(950)
' 5th, 50th, and 95th percentile t-beta
pbeta(1,1) = tbetasort(50)
pbeta(1,2) = tbetasort(500)
pbeta(1,3) = tbetasort(950)
' 5th, 50th, 95th, percentile DW
pDW (1,1) = DWsort(50)
pDW (1,2) = DWsort(500)
pDW (1,3) = DWsort(950)
' 5th, 50th, 95th percentile R2
pR2 (1,1) = R2asort(50)
pR2 (1,2) = R2asort(500)
pR2 (1,3) = R2asort(950)