*******************************************************************************
* PROCEDURE GARCH.SRC
*
*  GARCH will estimate (via MLE) a number of ARCH-type models,
*    with up to an ARMA(13,13) for the conditional mean equation,
*    up twenty-five exogenous variables, and up to an xARCH(4,4)
*    for the conditional variance with up to thirteen exogenous/
*    user-input variables in the conditional variance.
*
*  Syntax:
*     @GARCH(options)  series  start  end  *nresids  *cvar
*     #<suppl. card>   (OPTIONAL exogenous series for the
*                         mean equation w/ EXOG option)
*     #<suppl. card>   (OPTIONAL exogenous series for the
*                         variance equation w/ VEXOG)
*
*     NOTE: If both EXOG and VEXOG are used, the suppl.
*            cards must be in the above order.
*
*           Also, this version includes non-negativity constraints.
*
*
*  Parameters:
*     depvar     is the series for the dependant variable
*     start end  are the start and end points for the depvar;
*                  defaults to the entire length of the series
*     nresids    is the series to which the normalized
*                  residuals are saved
*     cvar       is the series to which the estimated
*                  conditional variance is saved.
*
*  Options:
*
*    MOD   =      ARCH/[GARCH]/IGARCH/EGARCH/LGARCH/LIGARCH
*                 The type of ARCH model to be estimated
*                 The ARCH(q) model has q lagged squared
*                     residuals in the conditional variance
*                 The GARCH(p,q) model adds p lagged
*                     conditional variance terms to the ARCH(p)
*                 The IGARCH(p,q) model forces the coefficients
*                     in the conditional variance to sum to one
*                 The EGARCH(p,q) is the exponential GARCH
*                     model of Nelson (1991)
*                 The LGARCH(p,q) adds a leverage term to account
*                     for asymmetry/skewness to a GARCH(p,q) model
*                 The LIGARCH(p,q) is the LGARCH model with the IGARCH
*                     constraint
*
*    DIST  =      [NORMAL]/TDIST/GED Likelihood based on Normal, Student-t,
*                   or Generalized Error Distribution
*
*    DET   =      NONE/[CONSTANT]  The deterministic components for the mean.
*
*    AR    =      [0] Number of autoregressive terms in mean, 0 to 13
*
*    MA    =      [0] Number of moving average terms in mean, 0 to 13
*
*    P     =      [0] Number lagged variance terms in conditional
*                     variance, 0 to 4
*
*    Q     =      [1] Number of lagged squared residuals in
*                     conditional variance, 0 or 1 to 4
*
*    EXOG/[NOEXOG] Use EXOG to alert the proc of an additional suppl.
*                     card for exogenous variables additively entering
*                     the mean equation.  UP TO TWENTY-FIVE EXOGENOUS
*                     VARIABLES CAN BE USED
*    VEXOG/[NOVEXOG] Use VEXOG to alert the proc of an additional suppl.
*                     card for exogenous variables additively entering
*                     the conditional variance equation.  UP TO
*                     THIRTEEN EXOGENOUS VARIABLES CAN BE USED
*
*    AIM/[NOAIM]  Add an ARCH-IN-MEAN term: a function of the conditional
*                     variance term enters the mean equation
*    TERM  =      [ROOT]/NONE
*
*                 If ARCH-in-Mean is selected (AIM), the conditional
*                    variance term (NONE) or its square root (ROOT) appears
*                    in the mean equation
*
*    SCLAGS=      [24] Number of lags to test for serial correlation in the
*                    normalized residuals and squared normalized residuals
*
*    SCSPAN=      [4]  Span for the Ljung-Box test for serial
*                    correlation
*
*    TRACE/[NOTRACE] Show the results for each iteration
*
*    ITER = [100] The maximum number of iterations for
*                    maximizing the log-likelihood function
*
*    CC QC1 QC2 QC3 QC4 PC1 PC2 PC3 PC4 LC LEVC MC VC1...VC13 = [0.0]
*       [cc defaults to %seesq of OLS regression if possible, else 0.5]
*                 Initial guesses for the conditional variance parameters.
*                    QCx are on lagged squared residuals,
*                    PCx are on lagged cond. var. terms,
*                    CC is the constant in variance,
*                    LC is the asymmtery term for EGARCH,
*                    LEVC is for the leverage GARCH.
*                    MC is the ARCH-in-Mean term
*                    VC1 to VC13 are the exogenous variables
*                      in the conditional variance equation
*
*    [SIMPLEX]/NOSIMPLEX Do a few iterations maximizing
*                    the log-likelihood function with
*                    the simplex algorithm before handing
*                    it off to BHHH.  If you are reasonably
*                    confident with your inputted initial guess,
*                    turn it off, as it's kind of slow.
*
*    SITER = [5]   Number of iterations for the simplex algorithm
*
*    [GRAPH]/NOGRAPH Graph (or not) the square root of the estimated
*                    conditional variance
*
*    [INTERACTIVE]/NOINTERACTIVE
*                 Menu-driven or not
*
*    NOTE:  The menu driven part allows one to input the
*           deterministic components, AR length, GARCH model,
*           and lags (p,q) in the GARCH model using menus.
*           The other options must be set, if different from their
*           default values, in the @GARCH(...) statement.
*
*
*
*******************************************************************************
*
* GARCH WILL ESTIMATE THE MEAN EQUATION:
*
* y(t) =  a0 + a1*y(t-1) + ... + a13*y(t-13) + b1*x1(t) + ... + b25*x25(t)
*              |--    AR Lags 1 to 13   --|    |-  up to 25 exog. vars. -|
*         + m*H(t)**.5 + u(t) - ma1*u(t-1) - ... - ma13*u(t-13)
*            ^                  |--     MA Lags 1 to 13     --|
*            |
*      ARCH-in-Mean Term
*
*     u(t) = SQRT(H(t))*e(t), where e(t) is mean zero, unit variance
*
* AND THE EQUATION FOR THE CONDITIONAL VARIANCE OF u(t):
*
*  ARCH(q)      q = 1 to 4
*    H(t) = c + q1*u(t-1)**2 + ... + q4*u(t-4)**2
*
*  GARCH(p,q)   p = 0 to 4, q = 1 to 4
*    H(t) = c + q1*u(t-1)**2 + ... + q4*u(t-4)**2 + p1*h(t-1) + ... + p4*h(t-4)
*
*  IGARCH(p,q)  p = 1 to 4, q = 1 to 4
*    same as GARCH(p,q), but it forces p1+...+p4+q1+...+q4 = 1
*
*  LGARCH(p,q)  p = 0 to 4, q = 1 to 4
*    adds a leverage term, lev*I(t-1)*u(t-1)**2, to the GARCH(p,q) model,
*    where I(t) = 1 if u(t) < 0, and 0 otherwise, to allow for asymmetry
*    in the conditional variance
*
*  LIGARCH(p,q) p = 1 to 4, q = 1 to 4
*    adds the leverage term to the IGARCH model
*
*  EGARCH(p,q)  p = 0 to 4, q = 1 to 4
*    allows for asymmetry in the conditional variance and allows for
*    negative coefficients, since it is exponetiated.
*    H(t) = exp(c + q1*G(t-1) + ... + q4*G(t-4) + p1*log(H(t-1)) + ...
*           + p4*log(H(t-4))
*       where
*    G(t) = |u(t)/SQRT(H(t))| - SQRT(2.0/%PI) - l1*u(t)/SQRT(H(t))
*       for the noraml distribution
*    G(t) = |u(t)/SQRT(H(t)) - 2*SQRT(D-2)/((D-1)*GAMMA((D+1)/2.)
*        /(GAMMA(D/2.)*GAMMA(0.5)) - l1*u(t)/SQRT(H(t))
*       for the t-distibution, and
*    G(t) = |u(t)/SQRT(H(t))| - LAMBDA*(2^(1/D))*GAMMA(2/D)/GAMMA(1/D) $
*           - l1*u(t)/SQRT(H(t))
*       for the generalized error distribution
*       where
*           LAMBDA = SQRT((2^(-2/D))*GAMMA(1/D)/GAMMA(3/D))
*       and  GAMMA(.) is the gamma function
*
*
* The conditional variance equation, with the VEXOG option,
*   allows one to add up to THIRTEEN exogenous varibales to the conditional
*   mean equation.  With the VEXOG option, the conditional variance estimated
*   is
*
*   H(t) = GARCH-MODEL-CHOSEN(t) + v1*vx1(1)+...+v13*vx13(t)
*
*       or, for the EGARCH model,
*
*   H(t) = EXP(log(EGARCH-MODEL-CHOSEN(t)) + v1*vx1(1)+...+v13*vx13(t))
*
*
* Estimation is done by maximum likelihood, assuming a
*   normal, a Student t, or a generalized error distribution
*   for the normalized residuals.  The GED was originally
*   used with the EGARCH model (Nelson, 1991).
*
* For all models except the EGARCH, the p and q coefficients are
*   forced to be positive.  If you see a negative coefficient
*   on q1,q2,q3,q4,p1,p2,p3,p4 for non-EGARCH models, the model thinks
*   it's positive, as abs(px) and abs(qx) is estimated.
*
* A GED is equivalent to a normal distribution for D = 2
*   (a paramter to be estimated), and has excess kurtosis, i.e.,
*   fatter tails for D < 2.  For the Student t, it is equivalent
*   to a normal distribution as D -> oo, but gets very close for
*   D > 30 or so.
*
*
* After estimation, the procedure will then test the normalized
*   residuals (i.e., the residuals divided by the square root of
*   the conditional variance) for serial correlation and ARCH,
*   and will test the squared normalized residuals for serial
*   correlation.
*
* The procedure then tests for sign and size bias effects:
*
*   Let m(t) = 1 if the normalized residual is negative, 0 otherwise
*       and let p(t) = 1 if it is positive, 0 otherwise.  Let
*       sm(t) = p(t)*e(t), where e(t) is the normalized residual, and
*       sp(t) = p(t)*e(t).
*   After regressing the squared normalized residuals on a constant,
*       m,sm,sp, the sign test tests if the coefficient on m(t) is
*       zero, the negative size tests the coefficient on sm(t), and
*       the positive size tests that of sp(t).  The joint tests tests
*       all three simultaneously.
*
* The estimated coefficients are saved, in the order shown as output,
*   to a vector called %COEFFS.
*
*
*******************************************************************************
*
*
*  Revision History:
*     January 2008:
*       - Changed INQUIRE command on main GARCH Procedure to properly
*           handle cases where the user does not supply explicit start and
*           end values.
*
*     Minor updates
*       - October, 2001: Added "D" to list of LOCAL REALs
*
*     Version 6, July-August, 1999:
*       - Added the ability to model the mean as an ARMA
*           process, rather than just an AR process -- up
*           to an ARMA(13,13)
*       - Tom Maycock reduced the amount of code considerably,
*           making the program much more efficient (THANKS TOM!)
*     Version 5, November-December, 1997:
*       - Added the GED to the choices for distributions
*       - Allowed up to 25 exogenous variables (formerly 13)
*       - Allowed up to 13 exogenous variables in the variance
*       - Allowed up to FOUR lags in the ARCH equations (formerly 3)
*     Version 4, October, 1997:
*       - Added the EXOG option to allow for exogenous variables in
*           the mean -- up to 13 of them
*       - Allowed up to an AR(13) in the mean (formerly AR(6))
*       - Added the ARCH-in-Mean specification
*       - Added the LIGARCH(p,q) model option
*       - Added the ability to choose the number of SIMPLEX iterations
*     Version 3, April, 1997:
*       - Added the [INTERACTIVE]/NOINTERACTIVE option, for a
*          _menu-driven_ version.
*       - Added the choice of having or not having a constant in the
*          mean equation.
*       - Fixed two mistakes in the leveraged-GARCH model pointed out
*          to me by Professor Trevor at Macquarie University.
*     Version 2, March 1996
*       - Added the LGARCH model
*       - Fixed several bugs
*
*  References:
*       Bera, A.K., and M.L. Higgins (1993): "ARCH Models: Properties,
*        Estimation and Testing," Journal of Economic Surveys, 7:4, pp.
*        305-366.
*       Engle, R.F. (1982): "Autoregressive Conditional
*        Heteroskedasticity with Estimates of the Variance
*        of United Kingdom Inflation,"  Econometrica
*        50, pp. 987-1007.
*       Hamilton, J.D. (1994): TIME SERIES ANALYSIS, Princeton
*        University Press, Chapter 21.
*       Nelson, D.B. (1991): "Conditional Heteroskedasticity
*        in Asset Returns:  A New Approach," Econometrica
*        52, pp. 347-370.
*       Trevor, R. (1994) GARCH.UV Procedure from the Estima BBS
*
*
*    A huge thanks to Rob Trevor of Macquarie University for doing
*       the hard part -- setting up the log-likelihood functions for
*       the EGARCH model and the models with a t-distribution;
*       I just did bookkeeping, more or less.  I would also like
*       to thank him for noticing a few mistakes in a previous
*       version of GARCH.SRC
*    A huge thanks to Tom Maycock of Estima for helpful suggestions
*       and for finding an elegant way to handle exogenous variables which
*       are lags (and explaining it to me).  He also greatly reduced
*       the length of the procedure with much more effient programming
*       than I am capable of, and found a couple of errors.
*    I would also like to thank Christos Tsoukalas for noticing a couple
*       errors; Nikos Nomikos for supplying a formula and references for
*       part of the EGARCH model when a t-distribution is used; Pete
*       Summers for discovering a problem with dates and exogenous variables.
*
*    Any errors, of course, are my fault.  If you find any, email me
*       at nmorin@frb.gov
*
*      --Norm, v6, 9/ 9/99
*              v5 12/ 2/97
*              v4 10/31/97
*              v3  4/ 1/97
*              v2  3/17/96
*              v1  6/25/94
*
******************************************************************************

******************************************************************************
*
PROCEDURE LB series start end

DECLARE real lbstats
TYPE SERIES     series
TYPE INTEGER     start end

OPTION INTEGER    LAGS 24
OPTION INTEGER DFC     0
OPTION INTEGER SPAN    0
OPTION SWITCH PRINT 1

LOCAL INTEGER     startl endl spanl count nobs
LOCAL SERIES    corr corr2 csum

INQUIRE(series=series) startl>>start endl>>end
COMPUTE nobs = endl - startl + 1

IF ((span.eq.0).or.(span.ge.lags))
   COMPUTE spanl = lags
ELSE
   COMPUTE spanl = span

CORRELATE(noprint,number=lags+1) series startl endl corr

SET TREND = t
SET corr2 1 lags+1 = corr**2/(nobs-trend)
SET csum 1 lags = 0

ACCUMULATE corr2 2 lags+1 csum 1

IF PRINT.eq.1 {
DISPLAY ' '
DISPLAY 'The Ljung-Box Q-Test for Serial Correlation in' %LABEL(series)
}

IF (spanl.lt.lags) {
  COMPUTE count = spanl
  WHILE (count.lt.lags) {
    IF (count.gt.dfc) {
      DISPLAY 'LB(' @-1 count @-1 ')'  @10 'Test Statistic:' $
	@-1 ####.#### nobs*(nobs+2)*csum(count) $
	@40 'Significance Level:' @-1 ####.#####  $
	%CHISQR(nobs*(nobs+2)*csum(count),count-dfc)
      COMPUTE count = count+spanl
      }
    ELSE
      COMPUTE count = count+spanl
    }
   }

COMPUTE lbstats = nobs*(nobs+2)*csum(lags)
COMPUTE lbsig = %CHISQR(nobs*(nobs+2)*csum(lags),lags-dfc)

IF PRINT.eq.1 {
DISPLAY 'LB(' @-1 lags @-1 ')'  @10 'Test Statistic:' @-1 ####.#### $
	   nobs*(nobs+2)*csum(lags) $
	  @40 'Significance Level:' @-1 ####.#####  $
	  %CHISQR(nobs*(nobs+2)*csum(lags),lags-dfc)
}


END
*
******************************************************************************

******************************************************************************
*
PROCEDURE ARCHTEST depvar start end

TYPE SERIES         depvar
TYPE INTEGRER       start end

OPTION INTEGER      lags    4
OPTION INTEGER      dfc     0
OPTION INTEGER      span    0

LOCAL SERIES        ressqr
LOCAL INTEGER       startl endl spanl lagnum num

IF ((span.eq.0).or.(span.gt.lags))
   COMPUTE spanl = lags
ELSE
   COMPUTE spanl = span

INQUIRE(series=depvar)  startl>>start  endl>>end

SET ressqr = depvar**2

DISPLAY ' '
DISPLAY 'F-Test of no ARCH vs. ARCH in' %LABEL(depvar)

IF (spanl.eq.lags) {
    LINREG(noprint) ressqr startl+lags endl
    #constant ressqr{1 to lags}
    EXCLUDE(noprint)
    # ressqr{1 to lags}
   DISPLAY 'ARCH(' @-1 lags @-1 ')'  @10 'Test Statistic:' $
	   @-1 ####.#### %CDSTAT  @40 'Significance Level:' $
	  @-1 ####.##### %SIGNIF
    }

ELSE  {
    COMPUTE lagnum = 0
    COMPUTE num = FIX(lags/spanl)
    DO jj = 1,num
	COMPUTE lagnum = lagnum + spanl
	LINREG(noprint,dfc=dfc) ressqr startl+lagnum endl
	#constant ressqr{1 to lagnum}
	EXCLUDE(noprint)
	# ressqr{1 to lagnum}
	DISPLAY 'ARCH(' @-1 lagnum @-1 ')'  @10 'Test Statistic:' $
	 @-1 ####.#### %CDSTAT  @40 'Significance Level:' $
	 @-1 ####.##### %SIGNIF
    END DO

    IF (FLOAT(num).lt.(FLOAT(lags)/spanl)) {
	LINREG(noprint,dfc=dfc) ressqr startl+lags endl
	#constant ressqr{1 to lags}
	EXCLUDE(noprint)
	# ressqr{1 to lags}
	DISPLAY 'ARCH(' @-1 lags @-1 ')'  @10 'Test Statistic:' $
	 @-1 ####.#### %CDSTAT @40 'Significance Level:' $
	 @-1 ####.##### %SIGNIF

	}

   }


END
*
******************************************************************************

******************************************************************************
*
PROCEDURE JARQUEBERA depvar start end

TYPE SERIES         depvar
TYPE INTEGER         start end

LOCAL INTEGER startl endl

COMPUTE jbstat = 0., jbsignif = 0.

INQUIRE(SERIES=depvar) startl>>start endl>>end

STATISTICS(noprint) depvar startl endl

COMPUTE jbstat = (%nobs/6.)*%skewness**2 + (%nobs/24.)*%kurtosis**2

COMPUTE jbsignif = %CHISQR(jbstat,2)
DISPLAY ' '
DISPLAY 'The Jarque-Bera Normality Test, ChiSqr(2), for' %LABEL(depvar)
DISPLAY @10 'Test Statistic:' @-1 ####.#### jbstat  $
	@40 'Significance Level:' @-1 ####.##### jbsignif


END
*
******************************************************************************

******************************************************************************
******************************************************************************
*
PROCEDURE GARCH series start end nresids1 cvar1

TYPE     SERIES      series *nresids1 *cvar1
TYPE     INTEGER    start end

OPTION    CHOICE     MOD      2  ARCH GARCH IGARCH EGARCH LGARCH LIGARCH
OPTION    CHOICE     DIST     1  NORMAL TDIST GED
OPTION    INTEGER    AR       0
OPTION    INTEGER    MA       0
OPTION    INTEGER    P        0
OPTION    INTEGER    Q        1
OPTION    CHOICE     DET      2  NONE CONSTANT
OPTION    SWITCH     TRACE    0
OPTION    INTEGER    SCLAGS  24
OPTION    INTEGER    SCSPAN   4
OPTION    INTEGER    ITER   100
OPTION    INTEGER    SITER    5
OPTION    SWITCH     EXOG     0
OPTION    SWITCH     VEXOG    0
OPTION    SWITCH     SIMPLEX  1
OPTION    SWITCH     INTERACTIVE  1
OPTION    SWITCH     GRAPH    1
OPTION    SWITCH     AIM 0
OPTION    CHOICE     TERM     1  ROOT NONE
OPTION    REAL cc    0.0
OPTION    REAL qc1   0.0
OPTION    REAL qc2   0.0
OPTION    REAL qc3   0.0
OPTION    REAL qc4   0.0
OPTION    REAL pc1   0.0
OPTION    REAL pc2   0.0
OPTION    REAL pc3   0.0
OPTION    REAL pc4   0.0
OPTION    REAL lc    0.0
OPTION    REAL levc  0.0
OPTION    REAL mc    0.0
OPTION    REAL vc1   0.0
OPTION    REAL vc2   0.0
OPTION    REAL vc3   0.0
OPTION    REAL vc4   0.0
OPTION    REAL vc5   0.0
OPTION    REAL vc6   0.0
OPTION    REAL vc7   0.0
OPTION    REAL vc8   0.0
OPTION    REAL vc9   0.0
OPTION    REAL vc10  0.0
OPTION    REAL vc11  0.0
OPTION    REAL vc12  0.0
OPTION    REAL vc13  0.0

LOCAL REAL c d
LOCAL SERIES U H res indic lcvar GG
LOCAL INTEGER sclag spanl ars qs ps maxlag iters dfc deter
LOCAL INTEGER n inp aimt nx nvx traceit distr ii jj kk i j k
LOCAL INDEX XLIST VXLIST
LOCAL RECT XMAT VXMAT
LOCAL VECT[series] XSERIES VXSERIES
LOCAL VECT[labels] XLABEL VXLABEL
LOCAL FRML logl lambda yd regres resgresx garchvar g garchvarx
LOCAL FRML resaim regresxx regstub

******************************************************************************


*** Old way of setting range--doesn't work with current versions of RATS:
*INQUIRE(series=series) startl endl
*COMPUTE startl = %imax(start,startl)
*COMPUTE endl = %imin(end,endl)

*** Better way:
INQUIRE(series=series) startl>>start endl>>end


IF (EXOG.eq.1) {
   ENTER(varying) xlist

   INQUIRE(regressorlist)  STARTL2 ENDL2
   #xlist
   COMPUTE startl = %imax(startl,startl2)
   COMPUTE endl = %imin(endl,endl2)

   MAKE XMAT startl endl
   # xlist
   COMPUTE nx = %cols(xmat)
   DIM XSERIES(nx)
   DO i = 1,nx
      SET XSERIES(i) startl endl =  XMAT(t-startl+1,i)
   END DO

   DIM xlabel(nx)

   * Assign labels:
   SET depvar startl startl = 1.0
   EQUATION dummyeq depvar
   # xlist
   COMPUTE eqnarray = %eqntable(dummyeq)

   DO i = 1,nx
     IF eqnarray(2,i)>0
       COMPUTE xlabel(i) = $
	 %label( [series] eqnarray(1,i) )+'{'+eqnarray(2,i)+'}'
     ELSE
       COMPUTE xlabel(i) = %label([series] eqnarray(1,i))
   END DO
   DO i = 1,nx
     lables xseries(i)
     #xlabel(i)
   END DO
 }

ELSE {
   COMPUTE nx = 0
   DIM XSERIES(1)
   SET XSERIES(1) startl endl = 0.
  }

IF (VEXOG.eq.1) {
   ENTER(varying) vxlist

   INQUIRE(regressorlist)  STARTL2 ENDL2
   #vxlist
   COMPUTE startl = %imax(startl,startl2)
   COMPUTE endl = %imin(endl,endl2)

   MAKE VXMAT startl endl
   # vxlist
   COMPUTE nvx = %cols(vxmat)
   DIM VXSERIES(nvx)
   DO i = 1,nvx
      SET VXSERIES(i) startl endl =  VXMAT(t-startl+1,i)
   END DO

   DIM vxlabel(nvx)

   * Assign labels:
   SET depvar startl startl = 1.0
   EQUATION dummyeq depvar
   # vxlist
   COMPUTE eqnarray = %eqntable(dummyeq)

   DO i = 1,nvx
     IF eqnarray(2,i)>0
       COMPUTE vxlabel(i) = $
	 %label( [series] eqnarray(1,i) )+'{'+eqnarray(2,i)+'}'
     ELSE
       COMPUTE vxlabel(i) = %label([series] eqnarray(1,i))
   END DO
   DO i = 1,nvx
     lables vxseries(i)
     #vxlabel(i)
   END DO
 }

ELSE
   COMPUTE nvx = 0


DISPLAY ' '

IF %DEFINED(iter).and.(iter.gt.0)
   COMPUTE iters = iter
ELSE
   COMPUTE iters = 100
IF %DEFINED(siter).and.(siter.gt.0)
   COMPUTE siters = siter
ELSE
   COMPUTE siters = 5

COMPUTE nobs = endl-startl+1
IF sclags.gt.0
   COMPUTE sclag = sclags
ELSE
   COMPUTE sclag = %IMIN(36,FIX(%MIN(nobs/4.,3*SQRT(nobs))))
IF (scspan.gt.sclag).or.(scspan.le.0)
   COMPUTE spanl = sclag
ELSE
   COMPUTE spanl = scspan

COMPUTE MODN = MOD
COMPUTE DETER = DET
COMPUTE DISTR = DIST
COMPUTE AIMT  = AIM
COMPUTE TRACEIT = TRACE
COMPUTE miter = iters
COMPUTE simp = SIMPLEX
COMPUTE simpiter = siter

IF interactive.eq.0 {
    IF (ar.le.0)
	COMPUTE ars = 0
    ELSE IF (ar.gt.13)
	COMPUTE ars = 13
    ELSE
	COMPUTE ars = ar

    IF (ma.le.0)
	COMPUTE mas = 0
    ELSE IF (ma.gt.13)
	COMPUTE mas = 13
    ELSE
	COMPUTE mas = ma

    IF (p.le.0)
	COMPUTE ps = 0
    ELSE IF (p.gt.4)
	COMPUTE ps = 4
    ELSE
	COMPUTE ps = p

    IF (q.le.0)
	COMPUTE qs = 0
    ELSE IF (q.gt.4)
	COMPUTE qs = 4
    ELSE
	COMPUTE qs = q
 }

ELSE {
  LOOP ***************
   MENU 'Mean and Variance Equations'


    CHOICE 'Deterministic Components'
     {
     MENU  'Deterministic Components'
	 CHOICE 'None'
	    COMPUTE DETER = 1
	 CHOICE 'Constant'
	    COMPUTE DETER = 2
     END MENU
     }


    CHOICE 'Choose the AR lags for the mean'
     {
      DISPLAY(store=qstring) 'Number of AR lags'
	Query(prompt=qstring) ars
      IF (ars.lt.0).or.(ars.gt.13) {
	 DISPLAY 'Lag lenth must be from 0 to 13'
	 DISPLAY(store=qstring) 'Number of AR lags'
	 Query(prompt=qstring) ars
	}
      }

    CHOICE 'Choose the MA lags for the mean'
     {
      DISPLAY(store=qstring) 'Number of MA lags'
	Query(prompt=qstring) mas
      IF (mas.lt.0).or.(mas.gt.13) {
	 DISPLAY 'Lag lenth must be from 0 to 13'
	 DISPLAY(store=qstring) 'Number of MA lags'
	 Query(prompt=qstring) mas
	}
      }

    CHOICE 'Choose the ARCH Model'
     {
      MENU 'ARCH Model'
	 CHOICE 'ARCH(q)'
	    COMPUTE MODN = 1
	 CHOICE 'GARCH(p,q)'
	    COMPUTE MODN = 2
	 CHOICE 'IGARCH(p,q)'
	    COMPUTE MODN = 3
	 CHOICE 'EGARCH(p,q)'
	    COMPUTE MODN = 4
	 CHOICE 'Leverage GARCH(p,q)'
	    COMPUTE MODN = 5
	 CHOICE 'Leverage IGARCH(p,q)'
	    COMPUTE MODN = 6
      END MENU
      IF mod.ne.5 {
	 MENU 'Choose q -- lagged squared residuals'
	     CHOICE 'q = 1'
		 COMPUTE qs = 1
	     CHOICE 'q = 2'
		 COMPUTE qs = 2
	     CHOICE 'q = 3'
		 COMPUTE qs = 3
	     CHOICE 'q = 4'
		 COMPUTE qs = 4
	 END MENU
	 }
      ELSE {
	 MENU 'Choose q -- lagged squared residuals'
	     CHOICE 'q = 0'
		 COMPUTE qs = 0
	     CHOICE 'q = 1'
		 COMPUTE qs = 1
	     CHOICE 'q = 2'
		 COMPUTE qs = 2
	     CHOICE 'q = 3'
		 COMPUTE qs = 3
	     CHOICE 'q = 4'
		 COMPUTE qs = 4
	 END MENU
	 }
      IF MODN.gt.1
	{
	 MENU 'Choose p -- lagged variance terms'
	     CHOICE 'p = 1'
		 COMPUTE ps = 1
	     CHOICE 'p = 2'
		 COMPUTE ps = 2
	     CHOICE 'p = 3'
		 COMPUTE ps = 3
	     CHOICE 'p = 4'
		 COMPUTE ps = 4
	     CHOICE 'p = 0'
		 COMPUTE ps = 0
	 END MENU
	}
       }

    CHOICE 'ARCH-in-Mean'
      {
       MENU 'Add an ARCH-in-mean term'
	 CHOICE 'Yes'
	   COMPUTE aimt = 1
	 CHOICE 'No'
	   COMPUTE aimt = 0
       END MENU
      }

    CHOICE 'Innovation Distribution'
     {
     MENU  'Innovation Distribution'
	 CHOICE 'Normal'
	    COMPUTE DISTR = 1
	 CHOICE 'Student-t'
	    COMPUTE DISTR = 2
	 CHOICE 'Generalized Error Distribution'
	    COMPUTE DISTR = 3
     END MENU
     }

    CHOICE 'Trace Iterations'
      {
      MENU 'Trace Iterations'
	CHOICE 'Yes'
	  COMPUTE TRACEIT = 1
	CHOICE 'No'
	  COMPUTE TRACEIT = 0
      END MENU

      }

    CHOICE 'Maximization Options'
      {
      LOOP
       MENU 'Maximization Options'
	 CHOICE 'Maximum Iterations'
	   {
	   QUERY(prompt='Maximum iterations') miter
	   COMPUTE miter = %imax(miter,1)
	   }
	 CHOICE 'Simplex Algorithm'
	   {
	   MENU 'Start off with simplex algorithm'
	     CHOICE 'Yes'
	       COMPUTE simp = 1
	     CHOICE 'No'
	       COMPUTE simp = 0
	   END MENU
	   IF SIMP.eq.1 {
	     QUERY(prompt='Simplex iterations') simpiter
	     COMPUTE simpiter = %imax(miter,1)
	     }
	   }
	 CHOICE 'Continue'
	    BREAK
       END MENU
      END LOOP
       }


    CHOICE 'Continue'
	BREAK


   END MENU
  END LOOP ***********


} ; **** IF INTERACTIVE


IF (MOD.ne.5)
 COMPUTE qs = %IMAX(qs,1)

COMPUTE a0=a1=a2=a3=a4=a5=a6=a7=0.
COMPUTE a8=a9=a10=a11=a12=a13=0.
COMPUTE ma1=ma2=ma3=ma4=ma5=ma6=0.
COMPUTE ma7=ma8=ma9=ma10=ma11=ma12=ma13=0.
COMPUTE b1=b2=b3=b4=b5=b6=b7=b8=b9=b10=b11=b12=b13=0.
COMPUTE b14=b15=b16=b17=b18=b19=b20=b21=b22=b23=b24=b25=0.
COMPUTE q1=qc1,q2=qc2,q3=qc3,q4=qc4
COMPUTE p1=pc1,p2=pc2,p3=pc3,p4=pc4
COMPUTE lev=levc,l1=lc,m=mc
COMPUTE c   = cc
COMPUTE vx1 = vc1
COMPUTE vx2 = vc2
COMPUTE vx3 = vc3
COMPUTE vx4 = vc4
COMPUTE vx5 = vc5
COMPUTE vx6 = vc6
COMPUTE vx7 = vc7
COMPUTE vx8 = vc8
COMPUTE vx9 = vc9
COMPUTE vx10 = vc10
COMPUTE vx11 = vc11
COMPUTE vx12 = vc12
COMPUTE vx13 = vc13

************************************************************************


SET res startl endl = 0.

NONLIN

* IF exogenous variables:
IF EXOG.eq.1 {
 * IF any AR terms:
  IF (ars>=1.and.deter==2) {
    LINREG(noprint,def=resideq) series startl+ars endl res
    #constant series{1 to ars} xseries
    FRML(equation=resideq,names='a',addparms) regstub
  }
  ELSE IF (ars>=1.and.deter==1) {
    LINREG(noprint,def=resideq) series startl+ars endl res
    # series{1 to ars} xseries
    FRML(equation=resideq,names='a',addparms) regstub
  }
  ELSE IF (ars==0.and.DETER.eq.2) {
    LINREG(noprint,def=resideq) series startl+ars endl res
    #constant xseries
    FRML(equation=resideq,names='a',addparms) regstub
  }
  ELSE IF (ars==0.and.DETER.eq.1) {
    LINREG(noprint,def=resideq) series startl+ars endl res
    #xseries
    FRML(equation=resideq,names='a',addparms) regstub
  }
}
IF (EXOG.le.0) {
 * IF any AR terms:
  IF (ars>=1.and.DETER.eq.2) {
    LINREG(noprint,def=resideq) series startl+ars endl res
    #constant series{1 to ars}
    FRML(equation=resideq,names='a',addparms) regstub
  }
  ELSE IF (ars>=1.and.DETER.eq.1) {
    LINREG(noprint,def=resideq) series startl+ars endl res
    # series{1 to ars}
    FRML(equation=resideq,names='a',addparms) regstub
  }
  ELSE IF (ars==0.and.DETER.eq.2) {
    LINREG(noprint,def=resideq) series startl+ars endl res
    # constant
    FRML(equation=resideq,names='a',addparms) regstub
  }
  ELSE IF (ars==0.and.DETER.eq.1) {
    LINREG(noprint,def=resideq) series startl+ars endl res
    # xseries
    FRML regstub = 0.0
  }
}

FRML regres = series - regstub

IF cc<>0.0
 compute c = cc
ELSE IF (%valid(%seesq))
  COMPUTE c = %seesq
ELSE
  COMPUTE c = 0.5

*******************************************************************************

SET U startl endl = res

IF mas.le.0
   FRML regresxx = regres
IF mas>=1 {
   linreg(noprint,def=maeq) series
   # u{1 to mas}
   FRML(equation=maeq,addparms,names='ma') regresxxstub
   FRML regresxx = regresxxstub + regres
   }

*******************************************************************************

COMPUTE dfc = ars

SET indic startl endl = 0.
    DO t = startl+ars, endl
	IF (U(t).lt.0)
	COMPUTE indic(t) = 1.
    END DO
SET H startl endl = %seesq

STATISTICS(noprint) series
DISPLAY ' '
DISPLAY '********************************************************************'

DISPLAY '*   Estimating an ARMA(' @-1 ars @-1 ',' @-1 mas @-1 ') for the mean of' $
    %LABEL(series) @67 '*'


IF EXOG.eq.1 {
    DISPLAY '*     w/' nx 'exogenous variables' @67 '*'
    DO i = 1,nx
       DISPLAY '*        ' XLABEL(i) @67 '*'
    END DO
}

IF VEXOG.eq.1 {
    DISPLAY '*     w/' nvx $
       'exogenous variables in the conditional variance' @67 '*'
    DO i = 1,nvx
       DISPLAY '*        ' VXLABEL(i) @67 '*'
    END DO
}

nonlin(add) c

IF (MODN.eq.1) {

    DISPLAY '*     with an ARCH(' @-1 qs @-1 $
	 ') model for the conditional variance' @67 '*'

    IF (qs.eq.1)
       FRML garchvar = c + abs(q1)*U{1}**2
    ELSE IF (qs.eq.2)
       FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2
    ELSE IF (qs.eq.3)
       FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 + abs(q3)*U{3}**2
    ELSE IF (qs.eq.4)
       FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 + abs(q3)*U{3}**2 $
	    + abs(q4)*U{4}**2

    IF (qs.eq.1) {
	NONLIN(ADD) q1
	IF q1.le.0; COMPUTE q1 = .2
	}
    ELSE IF (qs.eq.2) {
	NONLIN(ADD) q1 q2
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	}
    ELSE IF (qs.eq.3) {
	NONLIN(ADD) q1 q2 q3
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	}
    ELSE IF (qs.eq.4) {
	NONLIN(ADD) q1 q2 q3 q4
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	IF q4.le.0; COMPUTE q4 = .1
	}

    }


ELSE IF (MODN.eq.2) {

    DISPLAY '*     w/ a GARCH(' @-1 ps @-1 ','   $
	 @-1 qs @-1 ') model for the conditional variance' @67 '*'

    IF (qs.eq.1) {
       IF (ps.eq.0)
	  FRML garchvar = c + abs(q1)*U{1}**2
       ELSE IF (ps.eq.1) {
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(p1)*H{1}
          }
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(p1)*H{1} + abs(p2)*H{2}
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(p1)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3}
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(p1)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3} + abs(p4)*H{4}
       }
    ELSE IF (qs.eq.2) {
       IF (ps.eq.0)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2
       ELSE IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 + abs(p1)*H{1}
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 + abs(p1)*H{1} $
	       + abs(p2)*H{2}
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 + abs(p1)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3}
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 + abs(p1)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3} + abs(p4)*H{4}
       }
    ELSE IF (qs.eq.3) {
       IF (ps.eq.0)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2
       ELSE IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(p1)*H{1}
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(p1)*H{1} + abs(p2)*H{2}
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(p1)*H{1} + abs(p2)*H{2} + abs(p3)*H{3}
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(p1)*H{1} + abs(p2)*H{2} + abs(p3)*H{3} $
	       + abs(p4)*H{4}
       }
    ELSE IF (qs.eq.4) {
       IF (ps.eq.0)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2
       ELSE IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 + abs(p1)*H{1}
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 + abs(p1)*H{1} $
	       + abs(p2)*H{2}
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 + abs(p1)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3}
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 + abs(p1)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3} + abs(p4)*H{4}
       }


    IF (qs.eq.1) {
	NONLIN(ADD) q1
	IF q1.le.0; COMPUTE q1 = .2
	}
    ELSE IF (qs.eq.2) {
	NONLIN(ADD) q1 q2
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	}
    ELSE IF (qs.eq.3) {
	NONLIN(ADD) q1 q2 q3
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	}
    ELSE IF (qs.eq.4) {
	NONLIN(ADD) q1 q2 q3 q4
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	IF q4.le.0; COMPUTE q4 = .1
	}
    IF (ps.eq.1) {
	NONLIN(ADD) p1
	IF p1.le.0; COMPUTE p1 = .2
	}
    ELSE IF (ps.eq.2) {
	NONLIN(ADD) p1 p2
	IF p1.le.0; COMPUTE p1 = .15
	IF p2.le.0; COMPUTE p2 = .1
	}
    ELSE IF (ps.eq.3) {
	NONLIN(ADD) p1 p2 p3
	IF p1.le.0; COMPUTE p1 = .15
	IF p2.le.0; COMPUTE p2 = .1
	IF p3.le.0; COMPUTE p3 = .1
	}
    ELSE IF (ps.eq.4) {
	NONLIN(ADD) p1 p2 p3 p4
	IF p1.le.0; COMPUTE p1 = .15
	IF p2.le.0; COMPUTE p2 = .1
	IF p3.le.0; COMPUTE p3 = .1
	IF p4.le.0; COMPUTE p4 = .1
	}

    }


ELSE IF (MODN.eq.3) {

    DISPLAY '*     w/ an IGARCH(' @-1 ps @-1 ','   $
	 @-1 qs @-1 ') model for the conditional variance' @67 '*'

    IF ps.eq.0 {
       DISPLAY
       DISPLAY ' For the IGARCH model, p cannot be 0; it will be set to 1'
       DISPLAY
       COMPUTE ps = 1
       }

    IF (qs.eq.1) {
       IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1}
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2}
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3}
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3} + abs(p4)*H{4}
       }
    ELSE IF (qs.eq.2) {
       IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1}
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2}
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2} + abs(p3)*H{3}
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2} + abs(p3)*H{3} $
	       + abs(p4)*H{4}
       }
    ELSE IF (qs.eq.3) {
       IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1}
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2}
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3}
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3} + abs(p4)*H{4}
       }
    ELSE IF (qs.eq.4) {
       IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1}
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2}
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3}
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3} + abs(p4)*H{4}
       }


    IF (qs.eq.1) {
	NONLIN(ADD) q1
	IF q1.le.0; COMPUTE q1 = .2
	}
    ELSE IF (qs.eq.2) {
	NONLIN(ADD) q1 q2
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	}
    ELSE IF (qs.eq.3) {
	NONLIN(ADD) q1 q2 q3
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	}
    ELSE IF (qs.eq.4) {
	NONLIN(ADD) q1 q2 q3 q4
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	IF q4.le.0; COMPUTE q4 = .1
	}
    IF (ps.eq.2) {
	NONLIN(ADD) p2
	IF p2.le.0; COMPUTE p2 = .1
	}
    ELSE IF (ps.eq.3) {
	NONLIN(ADD) p2 p3
	IF p2.le.0; COMPUTE p2 = .1
	IF p3.le.0; COMPUTE p3 = .1
	}
    ELSE IF (ps.eq.4) {
	NONLIN(ADD) p2 p3 p4
	IF p2.le.0; COMPUTE p2 = .1
	IF p3.le.0; COMPUTE p3 = .1
	IF p4.le.0; COMPUTE p4 = .1

	}

    }


ELSE IF (MODN.eq.4) {

    DISPLAY '*     w/ an EGARCH(' @-1 ps @-1 ','   $
	 @-1 qs @-1 ') model for the conditional variance' @67 '*'

    IF DISTR.eq.3 {
      COMPUTE D = 2.
      FRML LAMBDA = EXP((-1./D)*LOG(2) + 0.5*%LNGAMMA(1./D) $
	   - 0.5*%LNGAMMA(3./D))
      FRML G = ABS(U(T)/SQRT(H(T)))  - $
	   LAMBDA*2**(1./D)*EXP(%LNGAMMA(2./D)-%LNGAMMA(1./D)) - $
	   l1*U(T)/SQRT(H(T))
      }
    ELSE IF DISTR.eq.2 {
      COMPUTE D = 50.
      FRML G = ABS(U(T)/SQRT(H(T))) - $
	   2*SQRT(D-2)/((D-1)*EXP(%LNGAMMA((D+1)/2.)-%LNGAMMA(D/2.) - $
	   %LNGAMMA(0.5))) - l1*U(T)/SQRT(H(T))
     }
    ELSE
      FRML G = ABS(U(T)/SQRT(H(T))) - SQRT(2.0/%PI) - l1*U(T)/SQRT(H(T))

    FRML garchvar = $
       EXP( c + p1*LOG(H{1}) + p2*LOG(H{2}) + p3*LOG(H{3}) + p4*LOG(H{4}) $
	  + q1*G{1} + q2*G{2} + q3*G{3} +q4*G{4})

    IF (qs.eq.1) {
       IF (ps.eq.0)
	   FRML garchvar = exp(c + q1*G{1})
       ELSE IF (ps.eq.1)
	   FRML garchvar = exp(c + q1*G{1} + p1*log(H{1}))
       ELSE IF (ps.eq.2)
	   FRML garchvar = exp(c + q1*G{1} + p1*log(H{1}) + p2*log(H{2}))
       ELSE IF (ps.eq.3)
	   FRML garchvar = exp(c + q1*G{1} + p1*log(H{1}) $
		+ p2*log(H{2}) + p3*log(H{3}))
       ELSE IF (ps.eq.4)
	   FRML garchvar = exp(c + q1*G{1} + p1*log(H{1}) $
		+ p2*log(H{2}) + p3*log(H{3}) + p4*log(H{4}))
       }
    ELSE IF (qs.eq.2) {
       IF (ps.eq.0)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2})
       ELSE IF (ps.eq.1)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} $
		+ p1*log(H{1}))
       ELSE IF (ps.eq.2)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} $
		+ p1*log(H{1}) + p2*log(H{2}))
       ELSE IF (ps.eq.3)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} $
		+ p1*log(H{1}) + p2*log(H{2}) + p3*log(H{3}))
       ELSE IF (ps.eq.4)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} $
		+ p1*log(H{1}) + p2*log(H{2}) + p3*log(H{3}) $
		+ p4*log(H{4}))
       }
    ELSE IF (qs.eq.3) {
       IF (ps.eq.0)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} + q3*G{3})
       ELSE IF (ps.eq.1)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} + q3*G{3} $
		+ p1*log(H{1}))
       ELSE IF (ps.eq.2)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} + q3*G{3} $
		+ p1*log(H{1}) + p2*log(H{2}))
       ELSE IF (ps.eq.3)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} + q3*G{3} $
		+ p1*log(H{1}) + p2*log(H{2}) + p3*log(H{3}))
       ELSE IF (ps.eq.4)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} + q3*G{3} $
		+ p1*log(H{1}) + p2*log(H{2}) + p3*log(H{3}) $
		+ p4*log(H{4}))
       }
    ELSE IF (qs.eq.4) {
       IF (ps.eq.0)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} + q3*G{3} $
		+ q4*G{4})
       ELSE IF (ps.eq.1)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} + q3*G{3} $
		+ q4*G{4} + p1*log(H{1}))
       ELSE IF (ps.eq.2)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} + q3*G{3} $
		+ q4*G{4} + p1*log(H{1}) + p2*log(H{2}))
       ELSE IF (ps.eq.3)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} + q3*G{3} $
		+ q4*G{4} + p1*log(H{1}) + p2*log(H{2}) + $
		+ p3*log(H{3}))
       ELSE IF (ps.eq.4)
	   FRML garchvar = exp(c + q1*G{1} + q2*G{2} + q3*G{3} $
	     + q4*G{4} + p1*log(H{1}) + p2*log(H{2}) + p3*log(H{3}) $
	     + p4*log(H{4}))
       }


    IF (qs.eq.1) {
	NONLIN(ADD) q1
	IF q1.le.0; COMPUTE q1 = .2
	}
    ELSE IF (qs.eq.2) {
	NONLIN(ADD) q1 q2
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	}
    ELSE IF (qs.eq.3) {
	NONLIN(ADD) q1 q2 q3
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	}
    ELSE IF (qs.eq.4) {
	NONLIN(ADD) q1 q2 q3 q4
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	IF q4.le.0; COMPUTE q4 = .1
	}
    IF (ps.eq.1) {
	NONLIN(ADD) p1
	IF p1.le.0; COMPUTE p1 = .2
	}
    ELSE IF (ps.eq.2) {
	NONLIN(ADD) p1 p2
	IF p1.le.0; COMPUTE p1 = .15
	IF p2.le.0; COMPUTE p2 = .1
	}
    ELSE IF (ps.eq.3) {
	NONLIN(ADD) p1 p2 p3
	IF p1.le.0; COMPUTE p1 = .15
	IF p2.le.0; COMPUTE p2 = .1
	IF p3.le.0; COMPUTE p3 = .1
	}
    ELSE IF (ps.eq.4) {
	NONLIN(ADD) p1 p2 p3 p4
	IF p1.le.0; COMPUTE p1 = .15
	IF p2.le.0; COMPUTE p2 = .1
	IF p3.le.0; COMPUTE p3 = .1
	IF p4.le.0; COMPUTE p4 = .1
	}

    NONLIN(ADD) l1
    }


ELSE IF (MODN.eq.5)  {

    DISPLAY '*     w/ a Leveraged GARCH(' @-1 ps @-1 ','  $
	 @-1 qs @-1 ') model for the conditional variance' @67 '*'

    IF (qs.eq.0) {
       IF (ps.eq.0)
	  FRML garchvar = c + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.1)
	  FRML garchvar = c + abs(p1)*H{1} $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(p1)*H{1} $
	       + abs(p2)*H{2} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(p1)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(p1)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3} + abs(p4)*H{4} + lev*indic{1}*U{1}**2
       }
    ELSE IF (qs.eq.1) {
       IF (ps.eq.0)
	  FRML garchvar = c + abs(q1)*U{1}**2 $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(p1)*H{1} $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(p1)*H{1} $
	       + abs(p2)*H{2} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(p1)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(p1)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3} + abs(p4)*H{4} + lev*indic{1}*U{1}**2
       }
    ELSE IF (qs.eq.2) {
       IF (ps.eq.0)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(p1)*H{1} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(p1)*H{1} + abs(p2)*H{2} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(p1)*H{1} + abs(p2)*H{2} + abs(p3)*H{3} $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(p1)*H{1} + abs(p2)*H{2} + abs(p3)*H{3} $
	       + abs(p4)*H{4} + lev*indic{1}*U{1}**2
       }
    ELSE IF (qs.eq.3) {
       IF (ps.eq.0)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(p1)*H{1} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(p1)*H{1} + abs(p2)*H{2} $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(p1)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(p1)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3} + abs(p4)*H{4} + lev*indic{1}*U{1}**2
       }
    ELSE IF (qs.eq.4) {
       IF (ps.eq.0)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 + abs(p1)*H{1} $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 + abs(p1)*H{1} $
	       + abs(p2)*H{2} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 + abs(p1)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 + abs(p1)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3} + abs(p4)*H{4} $
	       + lev*indic{1}*U{1}**2
       }


    IF (qs.eq.1) {
	NONLIN(ADD) q1
	IF q1.le.0; COMPUTE q1 = .2
	}
    ELSE IF (qs.eq.2) {
	NONLIN(ADD) q1 q2
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	}
    ELSE IF (qs.eq.3) {
	NONLIN(ADD) q1 q2 q3
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	}
    ELSE IF (qs.eq.4) {
	NONLIN(ADD) q1 q2 q3 q4
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	IF q4.le.0; COMPUTE q4 = .1
	}
    IF (ps.eq.1) {
	NONLIN(ADD) p1
	IF p1.le.0; COMPUTE p1 = .2
	}
    ELSE IF (ps.eq.2) {
	NONLIN(ADD) p1 p2
	IF p1.le.0; COMPUTE p1 = .15
	IF p2.le.0; COMPUTE p2 = .1
	}
    ELSE IF (ps.eq.3) {
	NONLIN(ADD) p1 p2 p3
	IF p1.le.0; COMPUTE p1 = .15
	IF p2.le.0; COMPUTE p2 = .1
	IF p3.le.0; COMPUTE p3 = .1
	}
    ELSE IF (ps.eq.4) {
	NONLIN(ADD) p1 p2 p3 p4
	IF p1.le.0; COMPUTE p1 = .15
	IF p2.le.0; COMPUTE p2 = .1
	IF p3.le.0; COMPUTE p3 = .1
	IF p4.le.0; COMPUTE p4 = .1
	}

    IF (levc.eq.0).and.(qs.eq.0)
	COMPUTE lev = .1

    NONLIN(ADD) lev
  }


ELSE IF (MODN.eq.6) {

    DISPLAY '*     w/ a LIGARCH(' @-1 ps @-1 ','   $
	 @-1 qs @-1 ') model for the conditional variance' @67 '*'

    IF ps.eq.0 {
       DISPLAY
       DISPLAY ' For the IGARCH model, p cannot be 0; it will be set to 1'
       DISPLAY
       COMPUTE ps = 1
       }

    IF (qs.eq.1) {
       IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3} + abs(p4)*H{4} $
	       + lev*indic{1}*U{1}**2
       }
    ELSE IF (qs.eq.2) {
       IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2} $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2} + abs(p3)*H{3} $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2} + abs(p3)*H{3} $
	       + abs(p4)*H{4} + lev*indic{1}*U{1}**2
       }
    ELSE IF (qs.eq.3) {
       IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2}  + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} $
	       + abs(p2)*H{2} + abs(p3)*H{3} + abs(p4)*H{4} $
	       + lev*indic{1}*U{1}**2
       }
    ELSE IF (qs.eq.4) {
       IF (ps.eq.1)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.2)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2}  $
	       + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.3)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3} + lev*indic{1}*U{1}**2
       ELSE IF (ps.eq.4)
	  FRML garchvar = c + abs(q1)*U{1}**2 + abs(q2)*U{2}**2 $
	       + abs(q3)*U{3}**2 + abs(q4)*U{4}**2 $
	       + (1-q1-q2-q3-q4-p2-p3-p4)*H{1} + abs(p2)*H{2} $
	       + abs(p3)*H{3} + abs(p4)*H{4} + lev*indic{1}*U{1}**2
       }

    IF (qs.eq.1) {
	NONLIN(ADD) q1
	IF q1.le.0; COMPUTE q1 = .2
	}
    ELSE IF (qs.eq.2) {
	NONLIN(ADD) q1 q2
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	}
    ELSE IF (qs.eq.3) {
	NONLIN(ADD) q1 q2 q3
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	}
    ELSE IF (qs.eq.4) {
	NONLIN(ADD) q1 q2 q3 q4
	IF q1.le.0; COMPUTE q1 = .15
	IF q2.le.0; COMPUTE q2 = .1
	IF q3.le.0; COMPUTE q3 = .1
	IF q4.le.0; COMPUTE q4 = .1
	}
    IF (ps.eq.2) {
	NONLIN(ADD) p2
	IF p2.le.0; COMPUTE p2 = .1
	}
    ELSE IF (ps.eq.3) {
	NONLIN(ADD) p2 p3
	IF p2.le.0; COMPUTE p2 = .1
	IF p3.le.0; COMPUTE p3 = .1
	}
    ELSE IF (ps.eq.4) {
	NONLIN(ADD) p2 p3 p4
	IF p2.le.0; COMPUTE p2 = .1
	IF p3.le.0; COMPUTE p3 = .1
	IF p4.le.0; COMPUTE p4 = .1
	}
    NONLIN(ADD) lev
    }

*******************************************************************************

IF (VEXOG.eq.1) {

  IF MOD.ne.4 {

    IF nvx.eq.1 {
      NONLIN(ADD) VX1
      FRML garchvarx = garchvar + vx1*vxseries(1)
      }
    ELSE IF nvx.eq.2 {
      NONLIN(ADD) VX1 VX2
      FRML garchvarx = garchvar + vx1*vxseries(1) + vx2*vxseries(2)
      }
    ELSE IF nvx.eq.3 {
      NONLIN(ADD) VX1 VX2 VX3
      FRML garchvarx = garchvar + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3)
      }
    ELSE IF nvx.eq.4 {
      NONLIN(ADD) VX1 VX2 VX3 VX4
      FRML garchvarx = garchvar + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4)
      }
    ELSE IF nvx.eq.5 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5
      FRML garchvarx = garchvar + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5)
      }
    ELSE IF nvx.eq.6 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6
      FRML garchvarx = garchvar + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6)
      }
    ELSE IF nvx.eq.7 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7
      FRML garchvarx = garchvar + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7)
      }
    ELSE IF nvx.eq.8 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8
      FRML garchvarx = garchvar + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8)
      }
    ELSE IF nvx.eq.9 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8 VX9
      FRML garchvarx = garchvar + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8) $
	   + vx9*vxseries(9)
      }
    ELSE IF nvx.eq.10 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8 VX9 VX10
      FRML garchvarx = garchvar + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8) $
	   + vx9*vxseries(9) + vx10*vxseries(10)
      }
    ELSE IF nvx.eq.11 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8 VX9 VX10 VX11
      FRML garchvarx = garchvar + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8) $
	   + vx9*vxseries(9) + vx10*vxseries(10) + vx11*vxseries(11)
      }
    ELSE IF nvx.eq.12 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8 VX9 VX10 VX11 VX12
      FRML garchvarx =  garchvar + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8) $
	   + vx9*vxseries(9) + vx10*vxseries(10) + vx11*vxseries(11) $
	   + vx12*vxseries(12)
      }
    ELSE IF nvx.eq.13 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8 VX9 VX10 VX11 VX12 VX13
      FRML garchvarx = garchvar + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8) $
	   + vx9*vxseries(9) + vx10*vxseries(10) + vx11*vxseries(11) $
	   + vx12*vxseries(12) + vx13*vxseries(13)
      }

    }

  ELSE {
    IF nvx.eq.1 {
      NONLIN(ADD) VX1
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1))
      }
    ELSE IF nvx.eq.2 {
      NONLIN(ADD) VX1 VX2
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2))
      }
    ELSE IF nvx.eq.3 {
      NONLIN(ADD) VX1 VX2 VX3
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3))
      }
    ELSE IF nvx.eq.4 {
      NONLIN(ADD) VX1 VX2 VX3 VX4
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4))
      }
    ELSE IF nvx.eq.5 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5))
      }
    ELSE IF nvx.eq.6 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6))
      }
    ELSE IF nvx.eq.7 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7))
      }
    ELSE IF nvx.eq.8 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8))
      }
    ELSE IF nvx.eq.9 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8 VX9
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8) $
	   + vx9*vxseries(9))
      }
    ELSE IF nvx.eq.10 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8 VX9 VX10
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8) $
	   + vx9*vxseries(9) + vx10*vxseries(10))
      }
    ELSE IF nvx.eq.11 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8 VX9 VX10 VX11
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8) $
	   + vx9*vxseries(9) + vx10*vxseries(10) + vx11*vxseries(11))
      }
    ELSE IF nvx.eq.12 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8 VX9 VX10 VX11 VX12
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8) $
	   + vx9*vxseries(9) + vx10*vxseries(10) + vx11*vxseries(11) $
	   + vx12*vxseries(12))
      }
    ELSE IF nvx.eq.13 {
      NONLIN(ADD) VX1 VX2 VX3 VX4 VX5 VX6 VX7 VX8 VX9 VX10 VX11 VX12 VX13
      FRML garchvarx = exp(log(garchvar) + vx1*vxseries(1) + vx2*vxseries(2) $
	   + vx3*vxseries(3) + vx4*vxseries(4) + vx5*vxseries(5) $
	   + vx6*vxseries(6) + vx7*vxseries(7) + vx8*vxseries(8) $
	   + vx9*vxseries(9) + vx10*vxseries(10) + vx11*vxseries(11) $
	   + vx12*vxseries(12) + vx13*vxseries(13))
      }
    }

  }

ELSE
  FRML garchvarx = garchvar



*******************************************************************************


IF AIMT.eq.1 {
  DISPLAY '*     and an ARCH-in-Mean term in the conditional mean ' @67 '*'
  }
DISPLAY '********************************************************************'
DISPLAY ' '
DISPLAY
DISPLAY
IF DET.eq.2
  DISPLAY '  A1 is the constant in the mean equation'
IF ars.gt.0
  DISPLAY '  The A2->An coefficients refer to the AR equation for the mean'
IF EXOG.eq.1
  DISPLAY '  The B coefficients are on the exogenous variables'
IF AIMT.eq.1
  DISPLAY '  M is the coefficient on the ARCH-in-Mean term'
IF MAS.gt.0
   DISPLAY '  The MA coefficients refer to the moving average terms'
DISPLAY '  C is the constant in the conditional variance equation'
IF qs.gt.0
   DISPLAY '  The Q coefficients refer to the lagged squared residuals'
IF ps.gt.0
   DISPLAY '  The P coefficients refer to the lagged conditional variance'
DISPLAY '  The numbers in the A,MA,Q,P coefficients refer to the lag'
IF (MOD.eq.5).or.(MOD.eq.6)
  DISPLAY '  LEV is the coefficient on the leverage term'
IF VEXOG.eq.1
  DISPLAY $
 '  The VX coefficients are on the exogenous variables in the variance'
IF MOD.eq.4
  DISPLAY '  L is the asymmetry term'
IF DISTR.eq.2
  DISPLAY '  D is the degrees of freedom for the t-distribution'
IF DISTR.eq.3
  DISPLAY '  D is the scale parameter for the generalized error distribution'
DISPLAY


*******************************************************************************

IF (DISTR.eq.1) {

   IF AIMT.eq.1 {
     NONLIN(ADD) M
	IF TERM.eq.1
	   FRML resaim = regresxx-M*H**.5
	ELSE IF TERM.eq.2
	   FRML resaim = regresxx-M*H

     FRML logl = H(T) = garchvarx(T), U(T) = resaim(T), $
       indic(T) = %IF(resaim(T)<0,1,0), $
       -.5*(log(2*%PI)+log(H)+U**2/H)
	}

   ELSE {
       FRML logl = H(T) = garchvarx(T), U(T) = regresxx(T), $
	 indic(T) = %if(regresxx(T)<0,1,0), $
	 -.5*(log(2*%PI)+log(H)+U**2/H)
      }
  }


ELSE IF (DISTR.eq.2) {

   DISPLAY ' '
   DISPLAY '******************************************'
   DISPLAY '*  The distribution used is a Student-t  *'
   DISPLAY '******************************************'
   COMPUTE D = 50.0

   IF AIMT.eq.1 {
     NONLIN(ADD) M
     IF TERM.eq.1 ;
	  FRML resaim = regresxx-M*H**.5
     ELSE IF TERM.eq.2
	  FRML resaim = regresxx-M*H

     FRML logl = (H(T) = garchvarx(T)), U(T) = resaim(T), $
       indic(T) = %if(resaim(T)<0,1,0), $
       %LNGAMMA(0.5*(D+1.0)) - %LNGAMMA(0.5*D) - $
	0.5*LOG(D-2) - 0.5*LOG(H) - 0.5*log(%PI) -$
	0.5*(D+1.0)*LOG( 1.0+(U**2/(H*(D-2))) )
       }

    ELSE {

       FRML logl = (H(T)=garchvarx(T)),(U(T)=regresxx(T)),$
	  indic(T) = %if(regresxx(T)<0,1,0), $
	  %LNGAMMA(0.5*(D+1.0)) - %LNGAMMA(0.5*D) - $
	  0.5*LOG(D-2) - 0.5*LOG(H) - 0.5*log(%PI) - $
	  0.5*(D+1.0)*LOG( 1.0+(U**2/(H*(D-2))) )
       }
   }


ELSE IF DISTR.eq.3 {

   DISPLAY ' '
   DISPLAY '***************************************'
   DISPLAY '*  The distribution used is a G.E.D.  *'
   DISPLAY '***************************************'
   COMPUTE D = 2.0

   IF MODN.ne.4 {
      FRML LAMBDA = EXP((-1./D)*LOG(2) + 0.5*%LNGAMMA(1./D) $
	   - 0.5*%LNGAMMA(3./D))
    }

   IF AIMT.eq.1 {
     NONLIN(ADD) M
     IF TERM.eq.1
	  FRML resaim = regresxx-M*H**.5
     ELSE IF TERM.eq.2
	  FRML resaim = regresxx-M*H

     FRML logl = (H(T) = garchvarx(T)), U(T) = resaim(T), $
	 indic(T) = %if(resaim(T)<0,1,0), $
	 LOG(D) - 0.5*ABS(U/(LAMBDA*SQRT(H)))**D - $
	 (LOG(LAMBDA) + (1+1./D)*LOG(2) + %LNGAMMA(1./D)) - $
	    .5*LOG(H)
      }


   ELSE {

       FRML logl = (H(T)=garchvarx(T)),(U(T)=regresxx(T)),$
	  indic(T) = %if(regresxx(T)<0,1,0), $
	  LOG(D) - 0.5*ABS(U/(LAMBDA*SQRT(H)))**D - $
	  (LOG(LAMBDA) + (1+1./D)*LOG(2) + %LNGAMMA(1./D)) - $
	     .5*LOG(H)

      }
   }


COMPUTE maxlag = %imax(ps,qs)+%imax(ars,mas)

IF (SIMP.eq.1) {
MAXIMIZE(method=simplex,print,iter=simpiter) $
    logl startl+maxlag endl
    }
IF DISTR.gt.1 {
    NONLIN(ADD) D
    }
NLPAR(subiter=250,cvcrit=0.000001)
MAXIMIZE(method=bhhh,iter=miter,recursive,trace=traceit) $
    logl startl+maxlag endl

COMPUTE %COEFFS = %beta

IF (%converged==0)
 {
  DISPLAY '********************************************************'
  DISPLAY ' The model did not converge after ' %iters 'iterations.'
  DISPLAY ' You should re-execute the GARCH procedure using the '
  DISPLAY ' ITERS option to specify a larger number of iterations,'
  DISPLAY ' or input new initial values (or ones different from'
  DISPLAY ' the default values).'
  DISPLAY '********************************************************'
  DISPLAY
 }


IF (MODN.eq.3).or.(MODN.eq.6) {
    COMPUTE p1 = 1-q1-q2-q3-q4-p2-p3-p4
     DISPLAY 'The IGARCH constraint implies'
    DISPLAY '  P1 =' #.########## p1
    }


*******************************************************************************

SET cvar startl endl = c/(1-p1-p2-p3-p4)
SET lcvar startl endl = c/(1-p1-p2-p3-p4)
SET nresids startl endl = 0.
SET residz startl endl = 0.

DO t = startl+maxlag,endl
    COMPUTE residz(t) = U(t)
    COMPUTE cvar(t) = H(t)
END DO

SET nresids startl+maxlag endl = residz/sqrt(cvar)
SET sqnresids startl+maxlag endl = nresids**2

DISPLAY ' '
DISPLAY ' '
DISPLAY $
 '* Analyzing the normalized residuals for serial correlation and normality *'
DISPLAY ' '

COMPUTE spanl2 = spanl
IF (spanl.lt.(dfc+1)) {
   COMPUTE spanl = %IMIN(2*spanl,sclag)
    }

STATISTICS(noprint) nresids startl+maxlag endl
DISPLAY @10  'Skewness' @40  'Kurtosis'
DISP @10 ##.##### %skewness @40 ##.##### %kurtosis+3

*IF (DISTR.eq.2) {
*    DISPLAY ' '
*    DISPLAY 'The Kurtosis implied by a Student t distribution'
*    DISPLAY '    with ' @-1  ##.## D @-1 ' degrees of freedom is: ' @-1 $
*          ##.#### 3+6./(D-4)
*        }

@LB(lags=sclag,span=spanl,dfc=dfc) nresids startl+maxlag endl
@JARQUEBERA nresids startl+maxlag endl
@ARCHTEST(lags=sclag,span=spanl2,dfc=dfc+nx) nresids startl+maxlag endl

DISPLAY ' '
DISPLAY ' '
DISPLAY $
  '* Analyzing the squared normalized residuals for serial correlation *'

@LB(lags=sclag,span=spanl,dfc=dfc) sqnresids startl+maxlag endl
** @LB(lags=sclag,span=spanl) sqnresids startl+maxlag endl


*******************************************************************************

DISPLAY ' '
DISPLAY ' '
DISPLAY '* Sign and Size Bias Tests *'
DISPLAY ' '
SET splus startl endl = 0.; SET sminus startl endl = 0.
SET plus startl endl = 0.; SET minus startl endl = 0.
DO t = startl+maxlag,endl
   IF (nresids(t).lt.0) {
	COMPUTE splus(t) = 0.
	COMPUTE sminus(t) = 1.
	}
   ELSE {
	COMPUTE splus(t) = 1.
	COMPUTE sminus(t) = 0.
	}
   COMPUTE plus(t) = splus(t)*residz(t)
   COMPUTE minus(t) = sminus(t)*residz(t)
END DO


LINREG(noprint,dfc=dfc+nx) sqnresids startl+maxlag endl
#constant sminus{1} minus{1} plus{1}
COMPUTE sz1 = %beta(2)/sqrt(%seesq*%xx(2,2))
COMPUTE sgn1 = %beta(3)/sqrt(%seesq*%xx(3,3))
COMPUTE sgn2 = %beta(4)/sqrt(%seesq*%xx(4,4))

DISPLAY 'The Sign Bias Test'
DISPLAY  @10 'Test Statistic:' @-1 ####.#### sz1 $
    @40 'Significance Level:'  @-1 ####.#####  %TTEST(sz1,%NDF)
DISPLAY 'The Negative Size Bias Test'
DISPLAY  @10 'Test Statistic:' @-1 ####.#### sgn1 $
    @40 'Significance Level:'  @-1 ####.#####  %TTEST(sgn1,%NDF)
DISPLAY 'The Positive Size Bias Test'
DISPLAY  @10 'Test Statistic:' @-1 ####.#### sgn2 $
    @40 'Significance Level:'  @-1 ####.#####  %TTEST(sgn2,%NDF)
DISPLAY 'The Joint Test for the Three Effects'
EXCLUDE(noprint) ; # sminus{1} minus{1} plus{1}
DISPLAY @10 'Test Statistic:' @-1 ####.#### %CDSTAT $
    @40 'Significance Level:'  @-1 ####.##### %SIGNIF


DISPLAY ' '
DISPLAY '********************************************************************'


SET lcvar = cvar**.5

IF GRAPH.eq.1 {
   GRAPH(header='Square Root of Estimated Conditional Variance') 1
   #lcvar startl+maxlag+1 *
   }

IF %DEFINED(cvar1)
  SET cvar1 = cvar
IF %DEFINED(nresids1)
  SET nresids1 = nresids


END
*******************************************************************************




