PROCEDURE VHOLS series start end 

***********************************************************************
***
***    MAY 17, 1996
***
***    UNIX VERSION
***
***    To use interactively (i.e., menues):
*** 
***    @VHOLS series start end 
***    # the dependent variables
***
***        NOTE THAT YOU ANSWERS QUESTIONS BY TYPING THE ITEM NUMBER!!!
***
***
***    To use noninteractively:
***
***    @VARHAC(nointeractive,det= ,crit= ,maxlag= ) series start end 
***    # the dependent variables
***
***    where:
***        
***    det  =  "none" or "constant"    --  deterministic component in the VAR
***    crit =  "AIC", "SIC" or "Kbar"  --  model selection criteria
***        (Kbar means use Kbar as the lag in each eq. in the VAR)
***    maxlag =  maxlag to consider for the model selection criteria
***
***
***    For interactive and noninteractive use, NOPRINT will 
***    suppress printed output
***
***
***    The Spectral density estimator is saved to the symmetric matrix SDM  
***
***
***		REFERENCES
***
*** 	DenHaan, Wouter J. and Andrew T. Levin, Inferences from parametric
***		and non-parametric covariance matrix estimation procedures, 
***		NBER technical working paper #195.
***
*** 	DenHaan, Wouter J. and Andrew T. Levin, A practioner's guide to
***		robust covariance matrix estimation, NBER technical working paper #197.
***
***
***    EXAMPLES
***
***    source(noecho) VHOLS2.SRC
***
***    @vhols(nointeractive,det=constant,crit=aic,maxlag=4) y 59:1 93:4
***    #constant x1 x2 x3 x4
***
***
***    or
***
***    @vhols(nointeractive,noprint,det=none,crit=aic,maxlag=4) y 59:1 93:4
***    #constant x1 x2 x3 x4
***
***
***
***    !!!  ALL VARIABLES IN THE VARHAC PROCEDURE ARE LOCAL.  HOWEVER,
***    THE VARHAC PROCEDURE ESTIMATES A VAR SYSTEM, AND THE PARAMETERS
***    ESTIMATED BY THIS SYSTEM ARE NOT.
***
***
*********************************************************************

***
***
*** To use noninteractively:
***
***   @VARHAC(nointeractive,det= ,crit= ,maxlag= ) start end
***   # the dependent variables   
*** 
***   NOTE:  Enter the deteministic components for the primary
***          regression on the supplemental card with the other
***          explanatory variables
***
***   where:
***        
***    det   =  none constant --  deterministic component in the VAR
***    crit   =  AIC SIC Kbar  --  model selction criteria
***        (Kbar means use Kbar as the lag in each eq. in the VAR)
***    maxlag =  maxlag to consider for the model selection criteria
***
***
*** IN BOTH CASES:  NOPRINT will supress all printed output except the
***                     final coefficient estimates
***
***
*** The Spectral density estimator is saved to the symmetric matrix SDM  
***
*** References:
***    Den Haan, Wouter J., and Andrew Levin, 1994:  "Inferences from 
***        Parametric and Non-parametric Covariance Matrix 
***        Estimation Procedures"  UCSD Discussion Paper, 94-30.
***
*********************************************************************

DECL SYMM VCOVMAT SDM

TYPE INTEGER  start end
TYPE SERIES   series
 
OPTION CHOICE    DET 2  NONE CONSTANT 
OPTION INTEGER   MAXLAG 4
OPTION CHOICE    CRIT 1  AIC SIC Kbar
OPTION SWITCH    INTERACTIVE 1
OPTION SWITCH    PRINT 1

LOCAL INDEX    varlist levsx ccmat  laglevsx  depvarlist Ydex Udex
LOCAL INTEGER  lag vars counter startl endl nprint num deter 
LOCAL INTEGER  critmax startl2 endl2 vdeter 
LOCAL REAL     nin nin2
LOCAL SERIES   AICVEC SICVEC  
LOCAL RECT     lev laglev f11 F FP cmat AIC SIC U Y YLAGS Umat DISPMAT
LOCAL VECT[series]  DEPVAR RESIDVEC RESIDS
LOCAL SYMM     VCOV VCORR
LOCAL VECT     bb bb2 STDEV RSSMAT AICMINENT SICMINENT 
LOCAL VECT     AICMIN SICMIN SEvec Tvec SIGvec
LOCAL VECT[rect]    outmatrix VMAMAT COEFFMAT
LOCAL VECT[labels]  vlabel vlabel2 
LOCAL VECT[INDEX]   laglevsx2
LOCAL VECT[VECT]    COEFFVEC 

LOCAL INTEGER  Kbar Kmax ii jj kk i j k nvars msc 
  
ENTER(varying) varlist
COMPUTE vars = %rows(varlist)

INQUIRE(regressorlist)  startl>>start  endl>>end
#varlist

INQUIRE(series=series)  startl2>>start  endl2>>end
 
COMPUTE startl = %IMAX(startl,startl2)
COMPUTE endl = %IMIN(endl,endl2)

DIM vlabel(vars)
DO ii = 1,vars
   COMPUTE vlabel(ii) = %label([series] varlist(ii)) 
END DO

IF PRINT.eq.1 {
 DISP; DISP  
 DISP '*******************************************************'
 DISP '* The VARHAC Estimator of the Spectral Density Matrix *'
 DISP '*     by Norman Morin                                 *'
 DISP '*******************************************************' 
 DISP; DISP
}

DIMENSION levsx(0) 
DIMENSION laglevsx(0)
DIMENSION ccmat(0)

IF INTERACTIVE.eq.0 {
   COMPUTE Kbar = maxlag 
   COMPUTE deter = det
   COMPUT MSC = crit
 }

ELSE {
   DISP
   MENU 'Deterministic components in the VAR'
      CHOICE 'None'
          COMPUTE deter = 1
      CHOICE 'Constant'
          COMPUTE deter = 2
   END MENU
   
   DISP
   DISPLAY(store=qstring) $
     'Input the maximum lag to consider for AIC/Schwarz' 
   QUERY(prompt=qstring,status=check) Kbar 

   DISP
   MENU 'Model Selection Criterion'
      CHOICE 'AIC'
          COMPUTE MSC = 1
      CHOICE 'Schwarz'
          COMPUTE MSC = 2
      CHOICE 'The maximum lag previously input'
          COMPUTE MSC = 3
   END MENU
   DISP

   }

DO ii = 1,vars
    ENTER(varying) levsx
    # levsx varlist(ii) 
END DO

MAKE(trans) lev startl endl
# levsx


LINREG(noprint) series startl endl residseries
#varlist
COMPUTE betavec = %beta
 

DIM Umat(1,endl-startl+1)
DO jj = startl,endl
   COMPUTE Umat(1,jj-startl+1) = RESIDSERIES(jj)
END DO

DIMENSION depvar(vars)
DO ii = 1,vars
   SET depvar(ii) startl endl = $
         Umat(1,t-startl+1)*lev(ii,t-startl+1) 
END DO 


DIM DEPVARLIST(vars)

DO jj = 1,vars
   COMPUTE depvarlist(ii) = depvar(ii)
END DO

DO i = 1,vars
    ENTER(varying) laglevsx 
    # laglevsx depvarlist(i){1 to Kbar}
END DO

MAKE(trans) laglev startl endl
# laglevsx

DIM vlabel2(vars) 
DO ii = 1,vars 
   IF ii.eq.1 {
      IF VLABEL(1).eq.'Constant'
        COMPUTE Vlabel2(1) = 'U'
      ELSE
        COMPUTE vlabel2(1) = 'U' + VLABEL(1) 
      }
   ELSE
      COMPUTE vlabel2(ii) = 'U' + VLABEL(ii) 
END DO
DO ii = 1,vars
   LABEL depvar(ii)
   #Vlabel2(ii)
END DO
 
DIMENSION residvec(vars) 
DIMENSION coeffvec(vars)

IF MSC.lt.3 {
 IF (deter.eq.1)
    {
    DIM aic(vars,Kbar)
    COMPUTE aic = %CONST(0.) 
    DIM sic(vars,Kbar)
    COMPUTE sic = %CONST(0.) 
    DO lagnum = 1,Kbar
        DIM RESIDVEC(vars)
        
        SYSTEM 1 to vars
        LAGS 1 to lagnum  
        VARIABLES depvar
        END(SYSTEM)
        ESTIMATE(noftest,noprint) startl+Kbar endl RESIDVEC(1)
     
        DO jj = 1,vars
          set rr1 startl+Kbar endl = RESIDVEC(jj)**2
          acc rr1 startl+Kbar endl rr2
          COMPUTE aic(jj,lagnum) = $
                log(rr2(endl)/%nobs) + 2.*%nreg/%nobs
          COMPUTE sic(jj,lagnum) = $
                log(rr2(endl)/%nobs) + log(%nobs)*(1.*%nreg/%nobs) 
        END DO

    END DO 
    }
ELSE 
   {
    DIM aic(vars,Kbar)
    COMPUTE aic = %CONST(0.) 
    DIM sic(vars,Kbar)
    COMPUTE sic = %CONST(0.) 
    DO lagnum = 1,Kbar

        DIM RESIDVEC(vars)
        SYSTEM 1 to vars
        LAGS 1 to lagnum  
        VARIABLES depvar
        DETERMINISTIC CONSTANT
        END(SYSTEM)
        ESTIMATE(noftest,noprint) startl+Kbar endl RESIDVEC(1)
     
        DO jj = 1,vars
           set rr1 startl+Kbar endl = RESIDVEC(jj)**2
           acc rr1 startl+Kbar endl rr2 
           COMPUTE aic(jj,lagnum) = $
               log(rr2(endl)/%nobs) + 2.*%nreg/%nobs 
           COMPUTE sic(jj,lagnum) = $
               log(rr2(endl)/%nobs) + log(%nobs)*(1.*%nreg/%nobs) 
        END DO

    END DO 
    }
   
SET AICVEC 1 Kbar = 0.; SET SICVEC 1 Kbar = 0.
DIM AICMIN(vars) SICMIN(vars) AICMINENT(vars) SICMINENT(vars)
 
DO jj = 1,vars
   DO ii = 1,Kbar
      COMPUTE aicvec(ii) = aic(jj,ii)
      COMPUTE sicvec(ii) = sic(jj,ii)
   END DO
   EXTREMUM(noprint) aicvec 1 Kbar
   COMPUTE aicminent(jj) = %minent
   COMPUTE aicmin(jj) = %minimum 
   EXTREMUM(noprint) sicvec  1 Kbar
   COMPUTE sicminent(jj) = %minent
   COMPUTE sicmin(jj) = %minimum
END DO  

 IF MSC.eq.1 {
    COMPUTE critmax = FIX(AICMINENT(1))
    DO jj = 1,vars
       IF critmax.lt.(FIX(AICMINENT(jj)))
          COMPUTE critmax = FIX(AICMINENT(jj))
    END DO
    }
 ELSE {
    COMPUTE critmax = FIX(SICMINENT(1))
    DO jj = 1,vars
       IF critmax.lt.(FIX(SICMINENT(jj)))
          COMPUTE critmax = FIX(SICMINENT(jj))
    END DO
    }

IF PRINT.eq.1 { 
 DISP '  The maximum lag length considered is' Kbar
 DISP  
 WRITE(noskip,format='(6X8A9)') VLABEL2 
 DISP ' '
 DISP 'SIC'
 WRITE(noskip,format='(8F9.0)') SICMINENT 
 DISP 'AIC'
 WRITE(noskip,format='(8F9.0)') AICMINENT  
 DISP
 IF MSC.eq.1
    DISP '  The AIC is being used as the model selection criterion'
 ELSE
    DISP '  The SIC is being used as the model selection criterion'
 DISP
 }
}

ELSE {
  IF PRINT.eq.1 {
   DISP 'Using' Kbar 'lags in each equation in the VAR'
   DISP
   }
   DIM AICMIN(vars) SICMIN(vars) AICMINENT(vars) SICMINENT(vars)
   COMPUTE AICMINENT = %const(kbar)
   COMPUTE SICMINENT = %const(kbar) 
   COMPUTE critmax = kbar
  }


*****************************************************************
 
 DIM laglevsx2(vars) RESIDS(vars) COEFFVEC(vars)
 DO jj = 1,vars
    DIM COEFFVEC(jj)(FIX(AICMINENT(jj)))
    DIM LAGLEVSX2(jj)(0)
    DO i = 1,vars
        IF (MSC.eq.1)
           ENTER(varying) laglevsx2(jj) 
           # laglevsx2(jj) depvar(i){1 to FIX(AICMINENT(jj))} ;  
        ELSE
           ENTER(varying) laglevsx2(jj) 
           # laglevsx2(jj) depvar(i){1 to FIX(SICMINENT(jj))}
    END DO
    IF (DETER.eq.1) {
      LINREG(noprint) depvar(jj) startl+Kbar endl RESIDS(jj)
      #laglevsx2(jj)  
      COMPUTE COEFFVEC(jj) = %beta
      }
    ELSE {
      LINREG(noprint) depvar(jj) startl+Kbar endl RESIDS(jj)
      #laglevsx2(jj) constant
      COMPUTE COEFFVEC(jj) = %beta
      }
 END DO
  
 DIM YDEX(0) UDEX(0)

COMPUTE nin =  1./(endl-(startl+Kbar)+1.) 
 
DO jj = 1,vars
   ENTER(varying) UDEX
   #UDEX resids(jj)
END DO
CMOMENT startl+Kbar endl
#UDEX 
COMPUTE VCOVMAT = nin*%cmom

DIM COEFFMAT(critmax)
DO jj = 1,critmax
   DIM COEFFMAT(jj)(vars,vars)
   COMPUTE COEFFMAT(jj) = %CONST(0.)
END DO


IF MSC.eq.1 {
  DO ii = 1,vars                   ; ********  EQ #
     DO kk = 1,FIX(AICMINENT(ii))  ; ********  # LAGS
       DO jj = 1,vars             ; ********  # of VARIABLE in iith EQ
         COMPUTE COEFFMAT(kk)(ii,jj) = $
            coeffvec(ii)(kk+(jj-1)*FIX(AICMINENT(ii)) )
       END DO
    END DO 
  END DO   
  }
ELSE {
  DO ii = 1,vars                   ; ********  EQ #
     DO kk = 1,FIX(SICMINENT(ii))  ; ********  # LAGS
       DO jj = 1,vars             ; ********  # of VARIABLE in iith EQ
         COMPUTE COEFFMAT(kk)(ii,jj) = $
             coeffvec(ii)(kk+(jj-1)*FIX(SICMINENT(ii)) )
       END DO
    END DO 
  END DO   
  }

COMPUTE SUMMAT = %IDENTITY(vars)
COMPUTE SUMMAT2 = %IDENTITY(vars)

DO jj = 1,critmax
   COMPUTE SUMMAT = SUMMAT - COEFFMAT(jj)
   COMPUTE SUMMAT2 = SUMMAT2 - tr(COEFFMAT(jj))
END DO
COMPUTE SDM = INV(SUMMAT)*VCOVMAT*INV(SUMMAT2)

IF PRINT.eq.1 {
 DISP ' ' 
 DISP 'The Estimated Spectral Density Matrix at the Zero Frequency'
 DISP ' '
 WRITE(noskip,format='(7X8A14)') VLABEL2
 WRITE(noskip,format='(8g14.5)') SDM
}

 
CMOMENT startl endl
#varlist 

COMPUTE nin2 = 1./(endl-startl+1)
COMPUTE Dmat = nin2*%CMOM

COMPUTE Vmat = inv(Dmat)*SDM*inv(tr(Dmat))

DIM SEvec(vars) Tvec(vars) DISPmat(ii,5) SIGvec(vars)
DO ii = 1,vars
  COMPUTE SEvec(ii) =  Sqrt(Vmat(ii,ii)*nin2)
  COMPUTE Tvec(ii)  =  betavec(ii)/SEvec(ii)
  COMPUTE SIGvec(ii) = %TTEST(Tvec(ii),endl-startl+1-vars)
  COMPUTE DISPMAT(ii,2) = betavec(ii),DISPMAT(ii,3) = SEvec(ii), $
          DISPMAT(ii,4) = Tvec(ii), DISPMAT(ii,1) = ii, $
          DISPMAT(ii,5) = SIGvec(ii)
END DO 

DISP ; DISP 
DISP $
   '         Coefficients     Std. Errors    T-Statistics    Significance'
WRITE(noskip,format='(1f3.0,8g16.5)') DISPMAT
 
END 
