[Gmm-commits] r163 - in pkg: causalGel/R causalGel/vignettes gmm4 gmm4/R gmm4/man gmm4/vignettes

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Fri Dec 6 23:44:23 CET 2019


Author: chaussep
Date: 2019-12-06 23:44:22 +0100 (Fri, 06 Dec 2019)
New Revision: 163

Added:
   pkg/gmm4/man/momFct-methods.Rd
Modified:
   pkg/causalGel/R/causalGel.R
   pkg/causalGel/R/causalMethods.R
   pkg/causalGel/vignettes/causalGel.Rnw
   pkg/causalGel/vignettes/causalGel.pdf
   pkg/gmm4/NAMESPACE
   pkg/gmm4/R/gelModels-methods.R
   pkg/gmm4/R/gelfit-methods.R
   pkg/gmm4/R/gmmModels-methods.R
   pkg/gmm4/man/evalDMoment-methods.Rd
   pkg/gmm4/man/vcov-methods.Rd
   pkg/gmm4/vignettes/gelS4.Rnw
   pkg/gmm4/vignettes/gelS4.pdf
Log:
added robust to misspecification se

Modified: pkg/causalGel/R/causalGel.R
===================================================================
--- pkg/causalGel/R/causalGel.R	2019-12-05 22:35:49 UTC (rev 162)
+++ pkg/causalGel/R/causalGel.R	2019-12-06 22:44:22 UTC (rev 163)
@@ -96,7 +96,7 @@
     } else {
         theta0 <- NULL
     }
-    fit <- modelFit(object=model, initTheta=initTheta, theta0=theta0,
+    fit <- modelFit(model=model, initTheta=initTheta, theta0=theta0,
                     lambda0=lambda0, vcov=getVcov, coefSlv=coefSlv,
                     lamSlv=lamSlv, tControl=tControl, lControl=lControl)
     fit at call <- Call

Modified: pkg/causalGel/R/causalMethods.R
===================================================================
--- pkg/causalGel/R/causalMethods.R	2019-12-05 22:35:49 UTC (rev 162)
+++ pkg/causalGel/R/causalMethods.R	2019-12-06 22:44:22 UTC (rev 163)
@@ -93,12 +93,12 @@
 ## modelFit
 
 setMethod("modelFit", signature("causalGel"), valueClass="causalGelfit", 
-          definition = function(object, gelType=NULL, rhoFct=NULL,
+          definition = function(model, gelType=NULL, rhoFct=NULL,
                                 initTheta=c("gmm", "modelTheta0"), theta0=NULL,
                                 lambda0=NULL, vcov=FALSE, ...)
           {
               Call <- try(match.call(call=sys.call(sys.parent())), silent=TRUE)
-              if (class(Call)=="try-error")
+              if (inherits(Call,"try-error"))
                   Call <- NULL              
               res <- callNextMethod()
               res at call <- Call

Modified: pkg/causalGel/vignettes/causalGel.Rnw
===================================================================
--- pkg/causalGel/vignettes/causalGel.Rnw	2019-12-05 22:35:49 UTC (rev 162)
+++ pkg/causalGel/vignettes/causalGel.Rnw	2019-12-06 22:44:22 UTC (rev 163)
@@ -685,10 +685,41 @@
 and $\Ex(\mu(X)|Z=1)$. We can see that the moments are well balanced,
 at least up to six decimals.
 
-\section{The aceGEL function}
+\section{The causalGEL function}
 
+The function allows to estimate the causal effect without having to go
+through the step of creating the model. The different arguments are a
+mixture of the arguments of causalModel(), the \textit{solveGel}
+method of the ``gmm4'' package, and the \textit{modelFit}
+method. The average causal effect of the training program, assuming
+random assignment and using \EL, can be obtained as follows:
 
+<<>>=
+data(nsw)
+nsw$re78 <- nsw$re78/1000
+nsw$re75 <- nsw$re75/1000
+fit1 <- causalGEL(re78~treat, ~age+ed+black+hisp+re75, nsw, gelType="EL",
+                 momType="uncondBal")
+@ 
 
+Similarly, the \ACE and \ACT can be computed as 
+
+<<>>=
+fit2 <- causalGEL(re78~treat, ~age+ed+black+hisp+re75, nsw, gelType="EL",
+                  momType="ACE")
+fit3 <- causalGEL(re78~treat, ~age+ed+black+hisp+re75, nsw, gelType="EL",
+                  momType="ACT")
+@ 
+
+The results are presented in Table \ref{t2}.
+
+<<echo=FALSE, results='asis'>>=
+texreg(list(fit1,fit2,fit3), digits=4, label="t2", 
+       custom.model.names=c("ACE(rand.)","ACE(non-random)", "ACT"),
+       caption="Causal Effect for a Training Program", ci.force=TRUE,
+       fontsize='footnotesize')
+@ 
+
 \bibliography{causal}
 
 \appendix

Modified: pkg/causalGel/vignettes/causalGel.pdf
===================================================================
(Binary files differ)

Modified: pkg/gmm4/NAMESPACE
===================================================================
--- pkg/gmm4/NAMESPACE	2019-12-05 22:35:49 UTC (rev 162)
+++ pkg/gmm4/NAMESPACE	2019-12-06 22:44:22 UTC (rev 163)
@@ -9,7 +9,7 @@
 
 importFrom("grDevices", rgb, col2rgb)
 
-importFrom("utils", capture.output)
+importFrom("utils", capture.output, head, tail)
 
 importFrom("stats", "ar", "as.formula", "model.matrix","vcov",
            "model.response", "na.omit", "terms", "residuals",
@@ -41,7 +41,7 @@
        tsls, modelFit, meatGmm, specTest, gmm4, restModel, modelResponse, DWH,
        modelDims, printRestrict, getRestrict, sysGmmModel, ThreeSLS, gelModel,
        rhoET, rhoEL, rhoEEL, rhoHD, Wu_lam, EEL_lam, REEL_lam, getLambda, 
-       solveGel, getImpProb, rhoETEL, rhoETHD, gel4)
+       solveGel, getImpProb, rhoETEL, rhoETHD, gel4, momFct)
  
 ###  S3 methods ###
 

Modified: pkg/gmm4/R/gelModels-methods.R
===================================================================
--- pkg/gmm4/R/gelModels-methods.R	2019-12-05 22:35:49 UTC (rev 162)
+++ pkg/gmm4/R/gelModels-methods.R	2019-12-06 22:44:22 UTC (rev 163)
@@ -48,35 +48,58 @@
 
 ################ evalDMoment ##########################
 
-setMethod("evalDMoment", "gelModels", function(object, theta, impProb=NULL)
-    {
-        if (object at vcov != "HAC")
-            {
-                evalDMoment(as(object, "gmmModels"), theta, impProb)
-            } else {
-                f <- function(theta, object, impProb)
-                    {
-                        gt <- evalMoment(object, theta)
-                        if (is.null(impProb))
-                            colMeans(gt)
-                        else
-                            colSums(gt*impProb)
-                    }
-                env <- new.env()
-                assign("theta", theta, envir = env)
-                assign("object", object, envir = env)
-                assign("f", f, envir = env)
-                assign("impProb", impProb, envir=env)
-                G <- numericDeriv(quote(f(theta, object, impProb)), "theta", 
-                                  env)
-                G <- attr(G, "gradient")
-                spec <- modelDims(object)
-                if (!is.matrix(G))
-                        G <- matrix(G,  spec$q, spec$k)
-                dimnames(G) <- list(spec$momNames, spec$parNames)
-                G
-            }
-    })
+setMethod("evalDMoment", "gelModels",
+          function(object, theta, impProb=NULL, lambda=NULL)
+          {
+              spec <- modelDims(object)
+              if (object at vcov != "HAC")
+              {
+                  evalDMoment(as(object, "gmmModels"), theta, impProb, lambda)
+              } else {
+                  if (!is.null(lambda))
+                  {
+                      if (length(lambda) != spec$q)
+                          stop("The length of lambda must be equal to the number of conditions")
+                      if (is.null(impProb))
+                      impProb <- 1/spec$n              
+                      f_lam <- function(theta, object, lambda)
+                      {
+                          gt <- evalMoment(object, theta)
+                          c(gt%*%lambda)
+                      }
+                      env <- new.env()
+                      assign("theta", theta, envir=env)
+                      assign("object", object, envir=env)
+                      assign("lambda", lambda, envir=env)
+                      assign("f_lam", f_lam, envir=env)
+                      G <- numericDeriv(quote(f_lam(theta, object, lambda)), "theta", env)
+                      G <- attr(G, "gradient")*c(impProb)
+                      colnames(G) <- spec$parNames
+                      rownames(G) <- NULL
+                      return(G)
+                  } 
+                  f <- function(theta, object, impProb)
+                  {
+                      gt <- evalMoment(object, theta)
+                      if (is.null(impProb))
+                          colMeans(gt)
+                      else
+                          colSums(gt*impProb)
+                  }
+                  env <- new.env()
+                  assign("theta", theta, envir = env)
+                  assign("object", object, envir = env)
+                  assign("f", f, envir = env)
+                  assign("impProb", impProb, envir=env)
+                  G <- numericDeriv(quote(f(theta, object, impProb)), "theta", 
+                                    env)
+                  G <- attr(G, "gradient")
+                  if (!is.matrix(G))
+                      G <- matrix(G,  spec$q, spec$k)
+                  dimnames(G) <- list(spec$momNames, spec$parNames)
+                  G
+              }
+          })
 
 ################ momentVcov  ##########################
 

Modified: pkg/gmm4/R/gelfit-methods.R
===================================================================
--- pkg/gmm4/R/gelfit-methods.R	2019-12-05 22:35:49 UTC (rev 162)
+++ pkg/gmm4/R/gelfit-methods.R	2019-12-06 22:44:22 UTC (rev 163)
@@ -200,8 +200,22 @@
 ## vcov
 
 setMethod("vcov", "gelfit",
-          function(object, withImpProb=FALSE, tol=1e-10) {
+          function(object, withImpProb=FALSE, tol=1e-10, robToMiss=FALSE) {
               spec <- modelDims(object at model)
+              if (robToMiss)
+              {
+                  eta <- c(coef(object), object at lambda)
+                  names(eta) <- NULL
+                  mod <- gmmModel(g=momFct, x=object, theta0=eta, vcov="MDS")
+                  fit <- evalModel(mod, theta=eta)
+                  v <- vcov(fit)
+                  spec <- modelDims(object at model)
+                  Sigma <- v[1:spec$k, 1:spec$k]
+                  dimnames(Sigma) <- list(spec$parNames, spec$parNames)
+                  SigmaLam <- v[(spec$k+1):nrow(v), (spec$k+1):ncol(v)]
+                  dimnames(SigmaLam) <- list(spec$momNames, spec$momNames)
+                  return(list(vcov_par = Sigma, vcov_lambda = SigmaLam))
+              }              
               q <- spec$q
               gt <- evalMoment(object at model, object at theta)
               n <- nrow(gt)
@@ -563,7 +577,6 @@
               arg <- list(...)
               ev <- new.env(parent.frame())
               theta0 <- arg$theta0
-              
               if (object at call[[1]] == "modelFit")
               {
                   model <- if(is.null(newModel))
@@ -636,4 +649,25 @@
                   call
           })
 
-   
+
+### This method is for specific moment functions
+
+setGeneric("momFct", function(eta, object, ...) standardGeneric("momFct"))
+
+## That moment function is used to rewrite GEL models into
+## GMM just identified models. It is useful to compute robust-to-misspecification s.e.
+
+setMethod("momFct", signature("numeric", "gelfit"),
+          function(eta, object) {
+              spec <- modelDims(object at model)
+              if (length(eta) != (spec$k+spec$q))
+                  stop("eta must include theta and lambda")
+              object at theta <- head(eta, spec$k)
+              object at lambda <- tail(eta, spec$q)
+              pt <- getImpProb(object, FALSE, FALSE)$pt*spec$n
+              gt <- evalMoment(object at model, object at theta)*pt
+              Gtl <- evalDMoment(object at model, object at theta, pt, object at lambda)
+              cbind(Gtl, gt)
+          })
+
+

Modified: pkg/gmm4/R/gmmModels-methods.R
===================================================================
--- pkg/gmm4/R/gmmModels-methods.R	2019-12-05 22:35:49 UTC (rev 162)
+++ pkg/gmm4/R/gmmModels-methods.R	2019-12-06 22:44:22 UTC (rev 163)
@@ -254,22 +254,56 @@
 setGeneric("evalDMoment", function(object, ...) standardGeneric("evalDMoment"))
 
 setMethod("evalDMoment", signature("regGmm"),
-          function(object, theta, impProb=NULL) {
+          function(object, theta, impProb=NULL, lambda=NULL)
+          {
               De <- Dresiduals(object, theta)
               Z <- model.matrix(object, "instrument")
+              spec <- modelDims(object)
               if (is.null(impProb))
                   impProb <- 1/nrow(Z)
-              G <- apply(De,2, function(x) colSums(Z*x*impProb))
-              spec <- modelDims(object)
-              if (!is.matrix(G))
+              if (!is.null(lambda))
+              {
+                  if (length(lambda) != spec$q)
+                      stop("The length of lambda must be equal to the number of conditions")
+                  Z <- c(Z%*%lambda)
+                  G <- De*Z*impProb
+                  colnames(G) <- spec$parNames
+                  rownames(G) <- NULL
+              } else {
+                  G <- apply(De,2, function(x) colSums(Z*x*impProb))
+                  if (!is.matrix(G))
                       G <- matrix(G,  spec$q, spec$k)
-              dimnames(G) <- list(spec$momNames, spec$parNames)
+                  dimnames(G) <- list(spec$momNames, spec$parNames)
+              }
               G
           })
 
 setMethod("evalDMoment", signature("functionGmm"),
-          function(object, theta, impProb=NULL) {
+          function(object, theta, impProb=NULL, lambda=NULL)
+          {
               spec <- modelDims(object)
+              if (!is.null(lambda))
+              {
+                  if (length(lambda) != spec$q)
+                      stop("The length of lambda must be equal to the number of conditions")
+                  if (is.null(impProb))
+                      impProb <- 1/spec$n              
+                  f_lam <- function(theta, object, lambda)
+                  {
+                      gt <- evalMoment(object, theta)
+                      c(gt%*%lambda)
+                  }
+                  env <- new.env()
+                  assign("theta", theta, envir=env)
+                  assign("object", object, envir=env)
+                  assign("lambda", lambda, envir=env)
+                  assign("f_lam", f_lam, envir=env)
+                  G <- numericDeriv(quote(f_lam(theta, object, lambda)), "theta", env)
+                  G <- attr(G, "gradient")*c(impProb)
+                  colnames(G) <- spec$parNames
+                  rownames(G) <- NULL
+                  return(G)
+              } 
               if (is.null(object at dfct))
               {
                   f <- function(theta, object, impProb)
@@ -294,14 +328,15 @@
                   G <- matrix(G,  spec$q, spec$k)
               dimnames(G) <- list(spec$momNames, spec$parNames)
               G
-              })
+          })
 
+
 setMethod("evalDMoment", signature("formulaGmm"),
-          function(object, theta, impProb=NULL) {
-              res <- modelDims(object)
+          function(object, theta, impProb=NULL, lambda=NULL)
+          {
+              spec <- modelDims(object)              
               nt <- names(theta)
-              nt0 <- names(res$theta0)
-              spec <- modelDims(object)
+              nt0 <- names(spec$theta0)
               if (is.null(impProb))
                   impProb <- 1/spec$n
               if (length(theta) != length(nt0))
@@ -311,41 +346,61 @@
               if (!all(nt%in%nt0 & nt0%in%nt))
                   stop("names in theta dont match parameter names")
               varList <- c(as.list(theta), as.list(object at modelF))
+              if (!is.null(lambda))
+              {
+                  if (length(lambda) != spec$q)
+                      stop("The length of lambda must be equal to the number of conditions")
+                  f <- function(theta, lambda, eq)
+                  {
+                      Gi <- sapply(1:spec$q, function(j)
+                      {
+                          if (!is.null(eq[[j]]))
+                          {
+                              tmp <- eval(D(eq[[j]], i), varList)*lambda[j]
+                              if (length(tmp)==1)
+                                  tmp <- rep(tmp, spec$n)
+                              tmp <- c(tmp*impProb)
+                          } else {
+                              tmp <- numeric(spec$n)
+                          }
+                          tmp
+                      })
+                      rowSums(Gi)
+                  }
+                  nG <- list(NULL, spec$parNames)
+              } else {
+                  f <- function(theta, lambda, eq)
+                  {
+                      Gi <- sapply(1:spec$q, function(j)
+                      {                          
+                          if (!is.null(eq[[j]]))
+                          {
+                              tmp <- eval(D(eq[[j]], i), varList)
+                              if (length(tmp)>1)
+                                  tmp <- sum(tmp*impProb)
+                          } else {
+                              tmp <- 0
+                          }
+                          c(tmp)
+                      })
+                      Gi
+                  }
+                  nG <- list(spec$momNames, spec$parNames)
+              }                  
               G <- numeric()
               for (i in nt)
-                  {
-                      lhs <- sapply(1:res$q, function(j) {
-                          if (!is.null(res$fLHS[[j]]))
-                              {
-                                  tmp <- eval(D(res$fLHS[[j]], i), varList)
-                                  if (length(tmp)>1)
-                                      d <- sum(tmp*impProb)
-                                  else
-                                      d <- tmp
-                              } else {
-                                  d <- 0
-                              }
-                          c(d)})
-                      rhs <- sapply(1:res$q, function(j) {
-                          if (!is.null(res$fRHS[[j]]))
-                              {
-                                  tmp <- eval(D(res$fRHS[[j]], i), varList)
-                                  if (length(tmp)>1)
-                                      d <- sum(tmp*impProb)
-                                  else
-                                      d <- tmp
-                              } else {
-                                  d <- 0
-                              }
-                          c(d)})
-                      G <- cbind(G, lhs-rhs)
-                  }
+              {
+                  lhs <- f(i, lambda, spec$fLHS)
+                  rhs <- f(i, lambda, spec$fRHS)
+                  G <- cbind(G, lhs-rhs)
+              }
               if (!is.matrix(G))
                   G <- matrix(G,  spec$q, spec$k)
-              dimnames(G) <- list(spec$momNames, spec$parNames)
+              dimnames(G) <- nG 
               G
           })
 
+
 ###########   estfun :  Don't like it ###############
 
 estfun.gmmFct <- function(x, ...) x

Modified: pkg/gmm4/man/evalDMoment-methods.Rd
===================================================================
--- pkg/gmm4/man/evalDMoment-methods.Rd	2019-12-05 22:35:49 UTC (rev 162)
+++ pkg/gmm4/man/evalDMoment-methods.Rd	2019-12-06 22:44:22 UTC (rev 163)
@@ -13,6 +13,36 @@
 It computes the matrix of derivatives of the sample moments with respect
 to the coefficients.
 }
+
+\usage{
+\S4method{evalDMoment}{functionGmm}(object, theta, impProb=NULL,
+lambda=NULL)
+
+\S4method{evalDMoment}{formulaGmm}(object, theta, impProb=NULL,
+lambda=NULL)
+
+\S4method{evalDMoment}{regGmm}(object, theta, impProb=NULL,
+lambda=NULL)
+
+\S4method{evalDMoment}{sysGmmModels}(object, theta)
+
+\S4method{evalDMoment}{rslinearGmm}(object, theta)
+
+\S4method{evalDMoment}{gelModels}(object, theta, impProb=NULL,
+lambda=NULL)
+}
+
+\arguments{
+  \item{object}{An model object}
+  \item{theta}{A numerical vector of coefficients} 
+  \item{impProb}{If a vector of implied probablities is provided, the
+    sample means are computed using them. If not provided, the means are
+    computed using the uniform weight}
+  \item{lambda}{A vector of Lagrange multipliers associated with the
+  moment conditions. Its length must therefore match the number of
+  conditions. See details below.}
+}
+  
 \section{Methods}{
 \describe{
 
@@ -34,6 +64,22 @@
 \item{\code{signature(object = "rslinearGmm")}}{
 }
 }}
+
+\details{
+Without the argument \code{lambda}, the method returns a \eqn{q \times
+  k} matrix, where \eqn{k} is the number of coefficients, and \eqn{q} is
+the number of moment conditions. That matrix is the derivative of the
+sample mean of the moments with respect to the coefficient. 
+
+If \code{lambda} is provided, the method returns an \eqn{n \times k}
+  matrix, where \eqn{n} is the sample size. The ith row is
+  \eqn{G_i'\lambda}, where $G_i$ is the derivative of the moment
+  function evaluated at the ith observation. For now, this option is
+  used to compute robust-to-misspecified standard errors of GEL
+  estimators.
+}
+
+
 \examples{
 data(simData)
 theta <- c(1,1)

Added: pkg/gmm4/man/momFct-methods.Rd
===================================================================
--- pkg/gmm4/man/momFct-methods.Rd	                        (rev 0)
+++ pkg/gmm4/man/momFct-methods.Rd	2019-12-06 22:44:22 UTC (rev 163)
@@ -0,0 +1,30 @@
+\name{momFct-methods}
+\docType{methods}
+
+\alias{momFct-methods}
+\alias{momFct}
+\alias{momFct,numeric,gelfit-method}
+
+\title{Methods for Function \code{momFct} in Package \pkg{gmm4}}
+\description{
+ The methods computes the moment matrix. It is use to create special
+ moment functions
+}
+
+\usage{
+\S4method{momFct}{numeric,gelfit}(eta, object)
+}
+
+\arguments{
+  \item{eta}{A vector that includes the coefficient and the Lagrange multipliers}
+  \item{object}{An object of class \code{"gmmfit"}}
+}
+
+\section{Methods}{
+\describe{
+
+  \item{\code{signature(eta = "numeric", object = "gelfit")}}{
+}
+}}
+\keyword{methods}
+\keyword{misspecified}

Modified: pkg/gmm4/man/vcov-methods.Rd
===================================================================
--- pkg/gmm4/man/vcov-methods.Rd	2019-12-05 22:35:49 UTC (rev 162)
+++ pkg/gmm4/man/vcov-methods.Rd	2019-12-06 22:44:22 UTC (rev 163)
@@ -18,7 +18,8 @@
 
 \S4method{vcov}{tsls}(object, sandwich=TRUE, df.adj=FALSE)
 
-\S4method{vcov}{gelfit}(object, withImpProb=FALSE, tol=1e-10)
+\S4method{vcov}{gelfit}(object, withImpProb=FALSE, tol=1e-10,
+                        robToMiss=FALSE)
 }
 \arguments{
   \item{object}{An object of class \code{"gmmfit"}}
@@ -42,6 +43,8 @@
   \item{withImpProb}{Should we compute the moments with the implied
     probabilities}
   \item{tol}{Any diagonal less than \code{"tol"} is set to tol}
+  \item{robToMiss}{Should we compute a covariance matrix that is robust
+    to misspecification?}
 }
 \section{Methods}{
 \describe{

Modified: pkg/gmm4/vignettes/gelS4.Rnw
===================================================================
--- pkg/gmm4/vignettes/gelS4.Rnw	2019-12-05 22:35:49 UTC (rev 162)
+++ pkg/gmm4/vignettes/gelS4.Rnw	2019-12-06 22:44:22 UTC (rev 163)
@@ -653,7 +653,8 @@
 \item \textit{getImpProb}: It computes the vector of implied
   probabilities $\hat{p}_i$'s. By default, they are defined as:
   \[
-  \hat{p}_i = \frac{-\rho'(g_i(\hat{\theta}))/n}{\sum_{j=1}^n -\rho'(g_j(\hat{\theta}))/n}
+  \hat{p}_i = \frac{-\rho'(\hat{\lambda}'g_i(\hat{\theta}))/n}{\sum_{j=1}^n -
+    \rho'(\hat{\lambda}'g_j(\hat{\theta}))/n}
   \]
   The options are
   \begin{itemize}
@@ -901,6 +902,91 @@
 \textit{restModel} is applied to the model in fit3.
 \end{itemize}
 
+\subsection{Robust to misspecification standard errors}\label{sec:robtomiss}
+
+The covariance matrix of $\hat{\theta}$ that we described above
+assumes that the model is correctly specified. In that case, the
+constraints associated with the moment conditions are asymptotically
+not binding, which implies that $\hat{\lambda}=O_p(n^{-1/2})$. In
+misspecified models, some moment conditions are not satisfied and the
+associated Lagrange multipliers do not converge to 0. Under some
+regularity conditions \citep{schennach07}, the GEL estimator is still
+root-n consistent for a pseudo true value. However, inference using
+the above expression for the covariance matrices is not valid.
+
+One way to obtain valid inference is to rewrite the GEL model into a
+just-identified GMM one. If we collect the first order conditions of
+the GEL problem for $\theta$ and $\lambda$ into one vector, we get a
+GMM model with moment condition $E(\psi_i(\eta))=0$, where
+
+\[
+\psi_i(\eta) = \begin{pmatrix}
+  \rho'(\lambda'g_i(\theta))G_i(\theta)'\lambda\\
+  \rho'(\lambda'g_i(\theta))g_i(\theta)\\
+  \end{pmatrix}\,,
+\]
+
+and $\eta=\{\theta',\lambda'\}'$. If we define 
+\[
+\hat{\Omega}=\frac{1}{n}\sum_{i=1}^n \psi_i(\hat{\eta})\psi_i(\hat{\eta})'
+\]
+and 
+\[
+\hat{\Gamma}=\frac{1}{n}\sum_{i=1}^n \nabla_\eta \psi_i(\hat{\eta})\,,
+\]
+where $\nabla_\eta$ is the gradiant operator, then
+$\hat{\Gamma}^{-1}\hat{\Omega}\hat{\Gamma}^{-1}$ is a consistent
+estimator of the variance $\sqrt{n}(\hat{\eta}-\eta^*)$, where
+$\eta^*$ is the pseudo true value. The \textit{momFct} method for
+numeric vector and ``gelfit'' object computes the $\psi_i(\eta)$
+matrix. We could therefore obtain a robust-to-misspecification
+covariance matrix using the following steps:
+
+<<>>=
+mod1 <- gelModel(y~x1+x2, ~x2+z1+z2+z3, gelType="EL", data=simData)
+fit1 <- modelFit(mod1)
+@ 
+
+We create the $\eta$ vector:
+
+<<>>=
+eta <- c(coef(fit1), fit1 at lambda)
+names(eta) <- NULL
+@ 
+
+We create a GMM model using the \textit{momFct} method:
+
+<<>>=
+mod2 <- gmmModel(momFct, fit1, theta0=eta, vcov="MDS")
+@ 
+
+After, we just need to evaluate the model and use the \textit{vcov} method:
+
+<<>>=
+fit2 <- evalModel(mod2, eta)
+v <- vcov(fit2)[1:3,1:3]
+sqrt(diag(v))
+@ 
+
+Which differs a little from the non-robust standard errors:
+
+<<>>=
+sqrt(diag(vcov(fit1)$vcov_par))
+@ 
+
+A quicker way is to use set the option ``robToMiss'' to TRUE:
+
+<<>>=
+sqrt(diag(vcov(fit1, robToMiss=TRUE)$vcov_par))
+@ 
+
+Since the ... of the \textit{summary} method is passed to
+\textit{vcov}, we can obtain a robust-to-misspecification summary table:
+
+<<>>=
+summary(fit1, robToMiss=TRUE)@coef
+@ 
+
 \section{The \textit{evalModel} Method} \label{sec:gelmodels-evalmodel}
 
 This method create a ``gelfit'' object without estimation. It is

Modified: pkg/gmm4/vignettes/gelS4.pdf
===================================================================
(Binary files differ)



More information about the Gmm-commits mailing list