? =========================== GMMVAR ============================

proc gmmvar nob neq nparam z dulist difu a0 msel rvar trace ndf ;
?
?   Procedure to compute the robust form of the estimated covariance
?   of GMM estimated parameters, given the instruments and the 
?   derivatives of the equations with respect to parameters. 
?   The notation is that of Hansen (1984) and Davidson and McKinnon. 
?
?   This procedure must be run immediately after a GMM estimation.
?   The sample is assumed to be continuous form 1 to nob. 
?   The input arguments are the following:
?       nob      number of observations.
?       neq      numbr of equations.
?       nparam   number of estimated parameters.
?       z        the list of instruments used in GMM estimation (in order).
?       dulist   the list of variables into which the GMM residuals will 
?                be placed (corresponds to u). 
?       difu     the list of variables which contain the derivatives of 
?                the equations u w.r.t. the parameters; neq*nparam long.
?                The order is (eq 1, param 1), (eq 2, param 1), (eq 3, 
?                param 1), ...., (eq 1, param2), and so forth...
?       a0       The weighting matrix (covoc option) that was used by 
?                GMM; this matrix does not have to be optimal.
?       msel     A diagonal noc by noc matrix containing the mask that
?                was used in estimation to zero out some of the moments.
?   The output is 
?       rvar     An nparam by nparam symmetric matrix rvar that contains
?                the estimated robust variance-covariance matrix of the 
?                estimated params.
?       trace    The chi-squared value for the test of overidentifying
?                restrictions, equal to the optimal criterion function.
?       ndf      The actual degrees of freedom for the chi-squared test.
?
  length z ninst ;                          ? The number of instruments.
  set noc = neq*ninst ;                     ? The number of OCs.
  mform (type=sym,nrow=noc) covoc ;
  mat covoc = 0 ;
  smpl 1 nob ;
  unmake @res dulist ;
  do i = 1 to nob ;
    smpl i i ;
    mmake zmat z ; mmake umat dulist ;
    mat covoc = covoc + (umat'#zmat')*(umat#zmat) ;
  enddo ;
  mform (type=sym,nrow=noc) covoc ;
  mat covoc = msel'covoc*msel ;             ? Estimate of phi, the cov(OC).
  mat a0inv = (msel'a0*msel)" ;             ? The A matrix used in estimation.
  mat apa = a0inv'covoc*a0inv ;
  mform (nrow=noc,ncol=nparam) dfsum ;
  mat dfsum = 0 ;
  do i = 1 to nob ;
    smpl i i ;
    mmake zmat z ; mmake dumat difu ;
    mform (nrow=neq,ncol=nparam) dumat ;
    mat dfsum = dfsum + dumat#zmat' ;
  enddo ;

  mat dfsum = dfsum/nob ;
  mat dapad = dfsum'apa*dfsum ;
  mat dadinv = ((dfsum'a0inv*dfsum)")/nob ;
  print dadinv ;                            ? The non-robust v-cov estimate.
  mform (type=sym,nrow=nparam) dapad ;
  mat rvar = dadinv*dapad*dadinv' ;         ? The robust v-cov estimate.
  print rvar ;
  smpl 1 nob ;
?
?       Compute the chi-squared test of overidentifying restrictions.
?  
  mat ndf = rank(covoc)-nparam ;
  mat trace = nob*nob*@oc'covoc"@oc ;
  cdf(chisq,df=ndf) trace ;
  tstats @coef rvar ;
endproc gmmvar ;
stop ; end ;