cal 1976 2 4
allocate 1998:4
* PROGRAM TO COMPUTE R/S Test Statistic
* as described in Andrew Lo (1991) ECTA
* MODIFIED TO COMPUTE THE ANDREW LO (1991) 
* PROGRAM by O. Mikhail
* http://www.bus.ucf.edu/omikhail

*declare series acov
*declare vector[series] autocov


* this program in RATS will read the raw data on UNEMPLOYMENT
* then will do RANGE/SCALE TEST.
* Q == Quarterly   92 obs
* UE0? == total uemployment                       D980712
* UE1? == total uemployment goods sector          D980733
* UE4? == total uemployment manufacturing sector  D980737
* UE5? == total uemployment services sector       D980739

* RATS PROGRAM TO COMPUTE R/S Test Statistic
* as described in Andrew Lo (1991) ECTA
* Written by O. Mikhail

open data c:\oss\data\hurst\ueQc.prn
data(format=free,org=obs) 1976:2 1998:4 ue0Q ue1Q ue4Q ue5Q

print / ue0Q ue1Q ue4Q ue5Q


labels ue0Q ue1Q ue4Q ue5Q
# 'Total UE (Q)' $
'UE GOODS  (Q)' $
'UE MANUF. (Q)' $
'UE SERVICE (Q)'


correlate(covariances,number=90) ue0Q * * autocov
correlate(number=90) ue0Q * * autocor

print / autocov autocor
source(noecho) c:\oss\data\hurst\rs.src

* i will loop to get the weighted series to modify the R/S
* note that autocov is the autocovariance series 
* computed outside the procedure
* sacov is the sum of the weighted autocovariances
* q will loop to produce the modified R/S under different q

compute [real] sacov = 0.0
compute q = 0
do q = 1,8
  compute i = 0
  do i=1,q
  compute [real] acov = (1.0 - (i/(q+1.0))) * autocov(i) 
  * display 'i' i ',acov' acov ',autocov(i)' autocov(i)
  * compute the sum of the weighted autocovariances
  compute [real] sacov = sacov + acov
  dis '****'
  display 'i=' i ',q=' q ',acov=' acov ',sacov=' sacov
  dis '****'
  end do i
dis ' '
dis 'For q =' q 'The sum of the weighted Autocovariances =' sacov
dis 'Now doing Modified R/S with q=' q
dis ' '
@rs ue0Q 1976:2 1998:4 sacov q
compute [real] sacov = 0.0
end do q
table


* HERE HOW TO CHOOSE Q ?
* FROM LO (1991) ECTA , page 1302
stats(noprint) ue0q
compute n = %nobs
compute [real] toto = autocor(2)
compute [real] to2 = 1.0 - (toto)**2.0
compute [real] kn = ((3.0 * n / 2.0 )**(1.0/3.0))*((2.0 * toto/(to2))**(2.0/3.0))
dis ' '
dis 'Number of Obs=' n ' ,Autocorr(2)=' autocor(2)
dis ' '
display 'kn' kn
dis ' '
* choice of q = [kn] where [] is greatest integer les than or equal to kn

* Do it only once for q =[kn]
  compute [real] sacov = 0.0
* Change kn from real to integer by using the function fix
  compute q = fix(kn)
display ' *** only once for q=' q
  compute i = 0
  do i=1,q
  compute [real] acov = (1.0 - (i/(q+1.0))) * autocov(i) 
  * display 'i' i ',acov' acov ',autocov(i)' autocov(i)
  * compute the sum of the weighted autocovariances
  compute [real] sacov = sacov + acov
  dis '****'
  display 'i=' i ',q=' q ',acov=' acov ',sacov=' sacov
  dis '****'
  end do i
dis ' '
dis 'For q =' q 'The sum of the weighted Autocovariances =' sacov
dis 'Now doing Modified R/S with q=' q
dis ' '
@rs ue0Q 1976:2 1998:4 sacov q

* BYE BYE NOW
* MAY THE FORCE BE WITH YOU
* O. Mikhail http://www.bus.ucf.edu/omikhail
end