[Gmm-commits] r128 - in pkg/gmm4: . R man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Tue Sep 11 22:35:04 CEST 2018


Author: chaussep
Date: 2018-09-11 22:35:04 +0200 (Tue, 11 Sep 2018)
New Revision: 128

Added:
   pkg/gmm4/man/formulaGmm-class.Rd
   pkg/gmm4/man/rformulaGmm-class.Rd
Modified:
   pkg/gmm4/NAMESPACE
   pkg/gmm4/R/allClasses.R
   pkg/gmm4/R/gmm4.R
   pkg/gmm4/R/gmmData.R
   pkg/gmm4/R/gmmModel.R
   pkg/gmm4/R/gmmModels-methods.R
   pkg/gmm4/R/gmmfit-methods.R
   pkg/gmm4/R/rGmmModel-methods.R
   pkg/gmm4/R/summaryGmm-methods.R
   pkg/gmm4/man/coef-methods.Rd
   pkg/gmm4/man/evalDMoment-methods.Rd
   pkg/gmm4/man/evalMoment-methods.Rd
   pkg/gmm4/man/getRestrict-methods.Rd
   pkg/gmm4/man/gmmFit-methods.Rd
   pkg/gmm4/man/gmmModel.Rd
   pkg/gmm4/man/modelDims-methods.Rd
   pkg/gmm4/man/momentStrength-methods.Rd
   pkg/gmm4/man/print-methods.Rd
   pkg/gmm4/man/printRestrict-methods.Rd
   pkg/gmm4/man/restGmmModel-methods.Rd
   pkg/gmm4/man/solveGmm-methods.Rd
   pkg/gmm4/man/subsetting.Rd
Log:
added formula based moment conditions and all methods

Modified: pkg/gmm4/NAMESPACE
===================================================================
--- pkg/gmm4/NAMESPACE	2018-06-28 21:10:11 UTC (rev 127)
+++ pkg/gmm4/NAMESPACE	2018-09-11 20:35:04 UTC (rev 128)
@@ -8,7 +8,7 @@
            "model.response", "na.omit", "terms", "residuals",
            "D", "numericDeriv", "sd", "optim", "lm", "pf", "coef", "update",
            "fitted", "lm.fit", "pchisq", "pnorm", "printCoefmat", "anova",
-           "model.frame", "reformulate", "formula")
+           "model.frame", "reformulate", "formula", "nlminb")
 importFrom("sandwich", "vcovHAC", "estfun","kernHAC",
            "bread","bwAndrews","bwNeweyWest","weightsAndrews",
            "weightsLumley", "vcovHC")
@@ -18,7 +18,8 @@
               "specTest", "summaryGmm", "rlinearGmm", "hypothesisTest",
               "numericORcharacter", "tsls", "rnonlinearGmm", "rfunctionGmm",
               "slinearGmm", "snonlinearGmm", "sysGmmModels",
-              "sgmmfit","stsls", "rslinearGmm", "rsnonlinearGmm", "rsysGmmModels")
+              "sgmmfit","stsls", "rslinearGmm", "rsnonlinearGmm", "rsysGmmModels",
+              "formulaGmm","rfunctionGmm")
 exportMethods(residuals, print, show, vcovHAC, coef, vcov, bread, summary, update,
               model.matrix, hypothesisTest, "[", merge, subset)
 

Modified: pkg/gmm4/R/allClasses.R
===================================================================
--- pkg/gmm4/R/allClasses.R	2018-06-28 21:10:11 UTC (rev 127)
+++ pkg/gmm4/R/allClasses.R	2018-09-11 20:35:04 UTC (rev 128)
@@ -44,9 +44,21 @@
                                        isEndo="logical"),
          prototype(vcov="MDS", kernel="Quadratic Spectral", bw="Andrews", prewhite=1L,
                    ar.method="ols", approx="AR(1)", tol=1e-7, dfct=NULL))
+setClass("formulaGmm", representation(modelF="data.frame", 
+                                      vcov="character",theta0="numeric",
+                                      n="integer", q="integer",k="integer",
+                                      parNames="character", momNames="character",
+                                      fRHS="list", fLHS="list",
+                                      kernel="character", bw="numericORcharacter",
+                                      prewhite="integer", ar.method="character",
+                                      approx="character", tol="numeric",
+                                      centeredVcov="logical", varNames="character",
+                                      isEndo="logical", isMDE="logical"),
+         prototype(vcov="MDS", kernel="Quadratic Spectral", bw="Andrews", prewhite=1L,
+                   ar.method="ols", approx="AR(1)", tol=1e-7))
 setClassUnion("regGmm", c("linearGmm", "nonlinearGmm"))
-setClassUnion("allNLGmm", c("nonlinearGmm", "functionGmm"))
-setClassUnion("gmmModels", c("linearGmm", "nonlinearGmm", "functionGmm"))
+setClassUnion("allNLGmm", c("nonlinearGmm", "functionGmm", "formulaGmm"))
+setClassUnion("gmmModels", c("linearGmm", "nonlinearGmm", "functionGmm", "formulaGmm"))
 
 ## gmmWeights
 
@@ -84,8 +96,12 @@
 setClass("rfunctionGmm", representation(R="list", cstSpec="list"),
          contains="functionGmm")
 
-setClassUnion("rgmmModels", c("rlinearGmm", "rnonlinearGmm", "rfunctionGmm"))
+setClass("rformulaGmm", representation(R="list", cstSpec="list"),
+         contains="formulaGmm")
 
+setClassUnion("rgmmModels", c("rlinearGmm", "rnonlinearGmm", "rfunctionGmm",
+                              "rformulaGmm"))
+
 ## hypothesisTest
 
 setClass("hypothesisTest", representation(test="numeric", hypothesis="character",

Modified: pkg/gmm4/R/gmm4.R
===================================================================
--- pkg/gmm4/R/gmm4.R	2018-06-28 21:10:11 UTC (rev 127)
+++ pkg/gmm4/R/gmm4.R	2018-09-11 20:35:04 UTC (rev 128)
@@ -1,8 +1,6 @@
 ################### the main gmm functions ###################
 ########## These functions ar to avoid having to builf model objects
 
-
-
 gmm4 <- function (g, x, tet0 = NULL, grad = NULL, 
                   type = c("twostep", "iter", "cue", "onestep"),
                   vcov = c("MDS", "HAC", "iid", "TrueFixed"),
@@ -30,10 +28,18 @@
         }
     if (is.list(g))
         {
-            model <- sysGmmModel(g=g, h=x, tet0=tet0, vcov=vcov,
-                                 kernel=kernel, crit=crit, bw=bw, prewhite=prewhite,
-                                 ar.method=ar.method, approx=approx, tol=kerntol,
-                                 centeredVcov=centeredVcov, data=data)
+            ## Formula of sysGMM? Need to find a better way.
+            model <- NULL
+            if (is.null(x) & !is.null(tet0))
+                model <- try(gmmModel(g=g, x=x, tet0=tet0, grad=grad, vcov=vcov,
+                                      kernel=kernel, crit=crit, bw=bw, prewhite=prewhite,
+                                      ar.method=ar.method, approx=approx, tol=kerntol,
+                                      centeredVcov=centeredVcov, data=data), silent=TRUE)
+            if (is.null(model) || class(model)=="try-error")
+                model <- sysGmmModel(g=g, h=x, tet0=tet0, vcov=vcov,
+                                     kernel=kernel, crit=crit, bw=bw, prewhite=prewhite,
+                                     ar.method=ar.method, approx=approx, tol=kerntol,
+                                     centeredVcov=centeredVcov, data=data)
         } else {
             model <- gmmModel(g=g, x=x, tet0=tet0, grad=grad, vcov=vcov,
                               kernel=kernel, crit=crit, bw=bw, prewhite=prewhite,

Modified: pkg/gmm4/R/gmmData.R
===================================================================
--- pkg/gmm4/R/gmmData.R	2018-06-28 21:10:11 UTC (rev 127)
+++ pkg/gmm4/R/gmmData.R	2018-09-11 20:35:04 UTC (rev 128)
@@ -106,6 +106,36 @@
              parNames=parNames, isEndo=isEndo, varNames=parNames)
     }
 
+
+
+.formGmmData <- function(formula, tet0, data)
+    {
+        res <- lapply(formula, function(f) .nlGmmData(f, ~1, tet0, data))
+        fRHS <- lapply(res, function(r) r$fRHS)
+        fLHS <- lapply(res, function(r) r$fLHS)
+        parNames <- res[[1]]$parNames
+        varNames <- do.call("c", lapply(res, function(r) r$varNames))
+        varNames <- unique(varNames)       
+        chkLHS <- sapply(fLHS, function(r) any(all.vars(r) %in% names(tet0)))
+        chkRHS <- sapply(fRHS, function(r) any(all.vars(r) %in% names(tet0)))
+        isMDE <- all(chkLHS) |  all(chkRHS)        
+        modelF <- sapply(varNames, function(n) data[[n]])
+        modelF <- as.data.frame(modelF)        
+        k <- length(tet0)
+        q <- length(formula)
+        if (is.null(names(formula)))
+            momNames <- paste("Mom_", 1:q, sep="")
+        else
+            momNames <- names(formula)
+        isEndo <- rep(FALSE, length(varNames))
+        n <- nrow(modelF)
+        list(modelF=modelF,  fRHS=fRHS, fLHS=fLHS, n=n, k=k, q=q,
+             momNames=momNames, parNames=parNames, varNames=varNames, isEndo=isEndo,
+             isMDE=isMDE)
+    }
+
+
+
 .nlGmmData <- function(formula, h, tet0, data)
     {
         varNames <- all.vars(formula)
@@ -178,17 +208,17 @@
              momNames=momNames, parNames=parNames, varNames=varNames, isEndo=isEndo)
     }
 
-.fGmmData <- function(g, x, theta0)
+.fGmmData <- function(g, x, thet0)
     {
-        mom <- try(g(theta0, x))
-        k <- length(theta0)        
-        if (is.null(names(theta0)))
+        mom <- try(g(thet0, x))
+        k <- length(thet0)        
+        if (is.null(names(thet0)))
             parNames <- paste("tet", 1:k, sep="")
         else
-            parNames <- names(theta0)
+            parNames <- names(thet0)
         if (any(class(mom)=="try-error"))
             {
-                msg <- paste("Cannot evaluate the moments at theta0\n",
+                msg <- paste("Cannot evaluate the moments at thet0\n",
                              attr(mom,"conditon"))
                 stop(msg)
             } else {

Modified: pkg/gmm4/R/gmmModel.R
===================================================================
--- pkg/gmm4/R/gmmModel.R	2018-06-28 21:10:11 UTC (rev 127)
+++ pkg/gmm4/R/gmmModel.R	2018-09-11 20:35:04 UTC (rev 128)
@@ -3,10 +3,10 @@
 
 ##################  Constructor for the gmmModels Classes  #####################
 
-gmmModel <- function(g, x, tet0=NULL,grad=NULL,
+gmmModel <- function(g, x=NULL, tet0=NULL,grad=NULL,
                      vcov = c("HAC", "MDS", "iid"),
                      kernel = c("Quadratic Spectral",  "Truncated", "Bartlett", "Parzen",
-                         "Tukey-Hanning"), crit = 1e-06,
+                          "Tukey-Hanning"), crit = 1e-06,
                      bw = "Andrews", prewhite = 1L, ar.method = "ols", approx = "AR(1)", 
                      tol = 1e-07, centeredVcov = TRUE, data=parent.frame())
     {
@@ -15,7 +15,7 @@
         if (is.numeric(bw))
             names(bw) <- "Fixed"
         if (!is.list(data) && !is.environment(data)) 
-            stop("'data' must be a list or an environment")
+            stop("'data' must be a list or an environment")    
         if (any(class(g)=="formula"))
             {
                 chk <- names(tet0) %in% all.vars(g)
@@ -58,7 +58,7 @@
                                       momNames=model$momNames, varNames=model$varNames,
                                       isEndo=model$isEndo)
                     }
-            } else {
+            } else if (class(g)=="function") {
                 model <- .fGmmData(g, x, tet0)
                 gmodel <- new("functionGmm", X=x, fct=g,
                               theta0=tet0, vcov=vcov, kernel=kernel, bw=bw,
@@ -68,8 +68,24 @@
                               n=model$n, parNames=model$parNames,
                               momNames=model$momNames, varNames=model$varNames,
                               isEndo=model$isEndo)
+            } else {
+                if (!is.null(x))
+                    stop("For formula GMM, x must be NULL. The moments are only defined as a list of formulas")
+                if (class(g) != "list")
+                    stop("For formula GMM, g must be a list of formulas")
+                if (any(sapply(g, function(gi) class(gi)) != "formula"))
+                    stop("For formula GMM, g must be a list of formulas")
+                model <- .formGmmData(g, tet0, data)
+                gmodel <- new("formulaGmm", modelF=model$modelF, 
+                              vcov=vcov, theta0=tet0,fRHS=model$fRHS,
+                              fLHS=model$fLHS, kernel=kernel, bw=bw,
+                              prewhite=as.integer(prewhite),
+                              ar.method=ar.method, approx=approx, tol=tol,
+                              centeredVcov = centeredVcov, k=model$k, q=model$q,
+                              n=model$n, parNames=model$parNames,
+                              momNames=model$momNames, varNames=model$varNames,
+                              isEndo=model$isEndo, isMDE=model$isMDE)
             }
         gmodel
     }
 
-

Modified: pkg/gmm4/R/gmmModels-methods.R
===================================================================
--- pkg/gmm4/R/gmmModels-methods.R	2018-06-28 21:10:11 UTC (rev 127)
+++ pkg/gmm4/R/gmmModels-methods.R	2018-09-11 20:35:04 UTC (rev 128)
@@ -138,7 +138,39 @@
               spec$fct(theta, object at X)
           })
 
+setMethod("evalMoment", signature("formulaGmm"),
+          function(object, theta) {
+              res <- modelDims(object)
+              nt <- names(theta)
+              nt0 <- names(res$theta0)
+              if (length(theta) != length(nt0))
+                  stop("The length of theta is not equal to the number of parameters")
+              if (is.null(nt))
+                  stop("theta must be a named vector")
+              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))
+              sapply(1:res$q, function(i) {
+                  if (!is.null(res$fLHS[[i]]))
+                      {
+                          lhs <- try(eval(res$fLHS[[i]], varList))
+                          if (any(class(lhs)=="try-error"))
+                              stop("Cannot evaluate the LHS")
+                      } else {
+                          lhs <- 0
+                      }
+                  if (!is.null(res$fRHS[[i]]))
+                      {
+                          rhs <- try(eval(res$fRHS[[i]], varList))
+                          if (any(class(lhs)=="try-error"))
+                              stop("Cannot evaluate the RHS")
+                      } else {
+                          lhs <- 0
+                      }
+                  c(lhs-rhs)})
+          })
 
+
 ################ evalDresiduals ##########################
 
 setGeneric("Dresiduals", function(object, theta, ...) standardGeneric("Dresiduals"))
@@ -196,6 +228,14 @@
                    fct=object at fct, dfct=object at dfct, isEndo=object at isEndo)
           })
 
+setMethod("modelDims", "formulaGmm",
+          function(object) {
+              list(k=object at k, q=object at q, n=object at n, parNames=object at parNames,
+                   momNames=object at momNames, theta0=object at theta0,
+                   fRHS=object at fRHS, fLHS=object at fLHS, isEndo=object at isEndo,
+                   isMDE=object at isMDE)
+          })
+
 ################ evalDMoment ##########################
 
 setGeneric("evalDMoment", function(object, ...) standardGeneric("evalDMoment"))
@@ -233,6 +273,37 @@
               G
               })
 
+setMethod("evalDMoment", signature("formulaGmm"),
+          function(object, theta) {
+              res <- modelDims(object)
+              nt <- names(theta)
+              nt0 <- names(res$theta0)
+              if (length(theta) != length(nt0))
+                  stop("The length of theta is not equal to the number of parameters")
+              if (is.null(nt))
+                  stop("theta must be a named vector")
+              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))
+              G <- numeric()
+              for (i in nt)
+                  {
+                      lhs <- sapply(1:res$q, function(j) {
+                          if (!is.null(res$fLHS[[j]]))
+                              d <- mean(eval(D(res$fLHS[[j]], i), varList))
+                          else
+                              d <- 0
+                          c(d)})
+                      rhs <- sapply(1:res$q, function(j) {
+                          if (!is.null(res$fRHS[[j]]))
+                              d <- mean(eval(D(res$fRHS[[j]], i), varList))
+                          else
+                              d <- 0
+                          c(d)})
+                      G <- cbind(G, lhs-rhs)
+                  }
+              G
+          })
 
 ###########   estfun :  Don't like it ###############
 
@@ -296,7 +367,6 @@
                   }
               w})
 
-
 ################### weights Object and methods: Is it too much??? #################
 
 
@@ -387,23 +457,28 @@
           })
 
 setMethod("solveGmm", signature("allNLGmm", "gmmWeights"),
-          function(object, wObj, theta0=NULL, ...)
-          {
-              if (is.null(theta0))
-                  theta0 <- modelDims(object)$theta0
+          function(object, wObj, theta0=NULL, algo=c("optim","nlminb"), ...)
+              {
+                  algo <- match.arg(algo)
+                  if (is.null(theta0))
+                      theta0 <- modelDims(object)$theta0
                   g <- function(theta, wObj, object)
                       evalObjective(object, theta, wObj)
                   dg <- function(theta, wObj, object)
-                  {
-                      gt <- evalMoment(object, theta)
-                      n <- nrow(gt)
-                      gt <- colMeans(gt)
-                      G <- evalDMoment(object, theta)
-                      obj <- 2*n*quadra(wObj, G, gt)
-                      obj
-                  }
-                  res <- optim(par=theta0, fn=g, gr=dg, method="BFGS", object=object,
-                               wObj=wObj, ...)
+                      {
+                          gt <- evalMoment(object, theta)
+                          n <- nrow(gt)
+                          gt <- colMeans(gt)
+                          G <- evalDMoment(object, theta)
+                          obj <- 2*n*quadra(wObj, G, gt)
+                          obj
+                      }
+                  if (algo == "optim")
+                      res <- optim(par=theta0, fn=g, gr=dg, method="BFGS", object=object,
+                                   wObj=wObj, ...)
+                  else
+                      res <- nlminb(start=theta0, objective=g, gradient=dg,
+                                    object=object, wObj=wObj, ...)
                   theta <- res$par
                   names(theta) <- modelDims(object)$parNames
                   list(theta=theta, convergence=res$convergence)
@@ -426,6 +501,12 @@
               list(strength=NULL, mess=NULL)
           })
 
+setMethod("momentStrength", signature("formulaGmm"), 
+          function(object, theta=NULL, ...)
+          {
+              list(strength=NULL, mess=NULL)
+          })
+
 setMethod("momentStrength", signature("linearGmm"), 
           function(object, theta, vcovType=c("OLS","HC","HAC")){
               spec <- modelDims(object)
@@ -513,6 +594,25 @@
                x at q <- length(momNames)
                x at momNames <- momNames
                x
+          })
+
+setMethod("[", c("formulaGmm", "numeric", "missing"),
+          function(x, i, j){
+              i <- unique(as.integer(i))
+              spec <- modelDims(x)
+              q <- spec$q
+              if (!all(abs(i) %in% (1:q))) 
+                  stop("SubMoment must be between 1 and q")
+               if (length(i)==q)
+                  return(x)
+               momNames <- x at momNames[i]
+               if (length(momNames)<spec$k)
+                   stop("The model is under-identified")
+              x at fRHS <- x at fRHS[i]
+              x at fLHS <- x at fLHS[i]
+              x at q <- length(momNames)
+              x at momNames <- momNames
+              x
            })
 
 setMethod("[", c("gmmModels", "missing", "missing"),
@@ -539,10 +639,40 @@
               x at n <- nrow(x at X)
               x})
 
+setMethod("subset", "formulaGmm",
+          function(x, i) {
+              x at modelF <- x at modelF[i,,drop=FALSE]
+              x at n <- nrow(x at modelF)
+              x})
+
 ## gmmFit
 
 setGeneric("gmmFit", function(object, ...) standardGeneric("gmmFit"))
 
+setMethod("gmmFit", signature("formulaGmm"), valueClass="gmmfit", 
+          definition = function(object, type=c("twostep", "iter","cue", "onestep"),
+              itertol=1e-7, initW=c("ident", "tsls"), weights="optimal", 
+              itermaxit=100, efficientWeights=FALSE, start=NULL, ...)
+              {
+                  if (object at isMDE && object at centeredVcov)
+                      {
+                          if (is.character(weights) && weights == "optimal")
+                              {
+                                  spec <- modelDims(object)
+                                  wObj <- evalWeights(object, spec$theta0, "optimal")
+                                  met <- getMethod("gmmFit", "gmmModels")
+                                  res <- met(object, weights=wObj, efficientWeights=TRUE,
+                                             ...)
+                                  res at type <- "mde"
+                                  return(res)
+                              } else {
+                                  callNextMethod()
+                              }
+                      } else {
+                          callNextMethod()
+                      }
+              })
+
 setMethod("gmmFit", signature("gmmModels"), valueClass="gmmfit", 
           definition = function(object, type=c("twostep", "iter","cue", "onestep"),
               itertol=1e-7, initW=c("ident", "tsls"), weights="optimal", 
@@ -563,7 +693,10 @@
                   spec <- modelDims(object)
                   if (spec$q==spec$k)
                       {
-                          weights <- "ident"
+                          # This allow to weight the moments in case of
+                          # large scale difference.
+                          if (!is.matrix(weights) && class(weights)!="gmmWeights")
+                              weights <- "ident"
                           type <- "onestep"
                       } else if (type == "onestep" && !is.matrix(weights)) {
                           weights <- "ident"
@@ -573,21 +706,21 @@
                           type <- "onestep"
                       }
                   if (type == "onestep")
-                  {
-                      if (class(weights)=="gmmWeights")
-                          wObj <- weights
-                      else
-                          wObj <- evalWeights(object, w=weights)
-                      res <- solveGmm(object, wObj, start, ...)
-                      convergence <- res$convergence
-                      efficientGmm <- ifelse(is.character(weights), FALSE,
+                      {
+                          if (class(weights)=="gmmWeights")
+                              wObj <- weights
+                          else
+                              wObj <- evalWeights(object, w=weights)
+                          res <- solveGmm(object, wObj, start, ...)
+                          convergence <- res$convergence
+                          efficientGmm <- ifelse(is.character(weights), FALSE,
                                              efficientWeights)
-                      ans <- new("gmmfit", theta=res$theta,
-                                 convergence=convergence, convIter=NULL, type=type,
-                                 wObj=wObj, model=object, call=Call, niter=i,
-                                 efficientGmm=efficientGmm)
-                      return(ans)
-                  }
+                          ans <- new("gmmfit", theta=res$theta,
+                                     convergence=convergence, convIter=NULL, type=type,
+                                     wObj=wObj, model=object, call=Call, niter=i,
+                                     efficientGmm=efficientGmm)
+                          return(ans)
+                      }
                   if (class(object) == "linearGmm")
                       {
                           if (object at vcov == "iid")
@@ -603,7 +736,7 @@
                           theta0 <- coef(tsls(object))
                       } else {
                           wObj <- evalWeights(object, NULL, "ident")
-                          theta0 <- solveGmm(object, wObj, start, ...)$theta                        
+                          theta0 <- solveGmm(object, wObj, start, ...)$theta
                       }
                   bw <- object at bw
                   if (type != "cue")
@@ -705,7 +838,7 @@
                           stop("You provided a named theta with wrong names")
                       theta <- theta[match(spec$parNames, names(theta))]
                   } else {
-                      if (class(object) == "nonlinearGmm")
+                      if (class(object) %in% c("formulaGmm","nonlinearGmm"))
                           stop("To evaluate nonlinear models, theta must be named")
                       names(theta) <- spec$parNames
                   }

Modified: pkg/gmm4/R/gmmfit-methods.R
===================================================================
--- pkg/gmm4/R/gmmfit-methods.R	2018-06-28 21:10:11 UTC (rev 127)
+++ pkg/gmm4/R/gmmfit-methods.R	2018-09-11 20:35:04 UTC (rev 128)
@@ -17,7 +17,8 @@
               ntype <- matrix(c("Two-Step GMM", "Iterated GMM", "CUE",
                                 "One-Step GMM with fixed weights","Two-Stage Least Squares",
                                 "Evaluated at a fixed Theta; No estimation",
-                                "twostep","iter","cue","onestep","tsls", "eval"),
+                                "One-Step Efficient M.D.E.",
+                                "twostep","iter","cue","onestep","tsls", "eval","mde"),
                               ncol=2)
               type <- ntype[match(x at type, ntype[,2]),1]
               spec <- modelDims(x at model)

Modified: pkg/gmm4/R/rGmmModel-methods.R
===================================================================
--- pkg/gmm4/R/rGmmModel-methods.R	2018-06-28 21:10:11 UTC (rev 127)
+++ pkg/gmm4/R/rGmmModel-methods.R	2018-09-11 20:35:04 UTC (rev 128)
@@ -182,6 +182,56 @@
         list(fct=fct, dfct=dfct, parNames=parNames, theta0=theta0, k=k)
     }
 
+.imposeFORMRestrict <- function(R, object)
+    {
+        chk <- sapply(R, function(r) all(all.vars(r) %in% object at parNames))
+        if (!all(chk))
+            stop("Wrong coefficient names in some of the restrictions")
+        rest <- sapply(R, function(r) as.character(r[[2]]))
+        if (any(duplicated(rest)))
+            stop("LHS of R must not have duplicated variables")
+        if (!all(sapply(rest, function(x) length(x)==1)))
+            stop("LHS of R formulas must contain only one coefficient")
+        dR <-numeric()
+        for (r in R)
+            {
+                lhs <- sapply(object at parNames, function(pn)
+                    eval(D(r[[2]], pn), as.list(object at theta0)))
+                rhs <- sapply(object at parNames, function(pn)
+                    eval(D(r[[3]], pn), as.list(object at theta0)))
+                dR <- rbind(dR, lhs-rhs)
+            }
+        if (any(is.na(dR)) || any(!is.finite(dR)))
+            stop("The derivative of the constraints at theta0 is either infinite or NAN")
+        if (qr(dR)$rank < length(R))
+            stop("The matrix of derivatives of the constraints is not full rank")
+        rhs <- list()
+        lhs <- list()
+        for (i in 1:length(object at fRHS))
+            {
+                rhs[[i]] <- as.character(object at fRHS[[i]])
+                if (!is.null(object at fLHS[[i]]))
+                    lhs[[i]] <- as.character(object at fLHS[[i]])
+                else
+                    lhs[[i]] <- NULL      
+                for (r in R)
+                    {
+                        rhs[[i]] <- gsub(as.character(r[2]), paste("(", as.character(r[3]),
+                                                                   ")", sep=""), rhs[[i]])
+                        if (!is.null(lhs[[i]]))
+                            lhs[[i]] <- gsub(as.character(r[2]),
+                                             paste("(", as.character(r[3]),
+                                                   ")", sep=""), lhs[[i]])
+                    }
+                rhs[[i]] <- parse(text=rhs[[i]])
+                lhs[[i]] <- parse(text=lhs[[i]])
+            }
+        k <- object at k-length(R)
+        parNames <- object at parNames[!(object at parNames %in% rest)]
+        theta0 <- object at theta0[!(object at parNames %in% rest)]        
+        list(rhs=rhs, lhs=lhs, parNames=parNames, theta0=theta0, k=k)
+    }
+
 ################## model.matrix and modelResponse #################
 ### I did not make model.response as generic because it is not
 ### a method in stats and I want different arguments
@@ -229,6 +279,14 @@
                    momNames=object at momNames, isEndo=res$isEndo)
           })
 
+setMethod("modelDims", "rformulaGmm",
+          function(object) {
+              res <- object at cstSpec
+              list(k=res$k, q=object at q, n=object at n, parNames=res$newParNames,
+                   momNames=object at momNames, theta0=res$theta0,
+                   fRHS=res$fRHS, fLHS=res$fLHS)
+          })
+
 setMethod("modelDims", "rnonlinearGmm",
           function(object) {
               res <- object at cstSpec
@@ -302,7 +360,6 @@
     hyp
 }
 
-
 .makeHypothesis <- function (cnames, hypothesis, rhs = NULL) 
 {
     l <- list()
@@ -394,6 +451,16 @@
                   }
           })
 
+setMethod("printRestrict", "rformulaGmm",
+          function(object){
+              cat("Constraints:\n")
+              for (i in 1:length(object at R))
+                  {
+                      cat("\t")
+                      print(object at R[[i]])
+                  }
+          })
+
 setMethod("printRestrict", "rfunctionGmm",
           function(object){
               cat("Constraints:\n")
@@ -411,6 +478,13 @@
               printRestrict(x)
           })
 
+setMethod("print", "rformulaGmm",
+          function(x)
+          {
+              callNextMethod()
+              printRestrict(x)
+          })
+
 setMethod("print", "rnonlinearGmm",
           function(x)
           {
@@ -512,6 +586,35 @@
               new("rfunctionGmm", R=R, cstSpec=cstSpec, object)
           })
 
+setMethod("restGmmModel", signature("formulaGmm"),
+          function(object, R, rhs=NULL) {
+              if (!is.null(rhs))
+                  warning("rhs is ignored for nonlinear models")
+              if (is.character(R))
+                  {
+                      R2 <- list()
+                      R <- gsub("=", "~", R, fixed=TRUE)
+                      for (r in R)
+                          R2 <- c(R2, as.formula(r, .GlobalEnv))
+                      R <- R2
+                  } else {
+                      if (!is.list(R))
+                          {
+                              if(class(R) != "formula")
+                                  stop("R must be a formula or a list of formulas")
+                              R <- list(R)
+                          } else {
+                              chk <- sapply(R, function(r) class(r)=="formula")
+                              if (!all(chk))
+                                  stop("R must be a formula, a list of formulas or a vector of characters")
+                          }
+                  }
+              res <- .imposeFORMRestrict(R, object)
+              cstSpec <- list(newParNames = res$parNames,
+                              originParNames=object at parNames,
+                              k=res$k, theta0=res$theta0, fRHS=res$rhs, fLHS=res$lhs)
+              new("rformulaGmm", R=R, cstSpec=cstSpec, object)
+          })
 
 ### Get the restriction matrices
 
@@ -550,6 +653,12 @@
                    orig.R=object at R, orig.rhs=NULL)
           })
 
+setMethod("getRestrict", "rformulaGmm",
+          function(object, theta) {
+              getMethod("getRestrict", "rnonlinearGmm")(object, theta)
+          })
+
+
 setMethod("getRestrict", "rfunctionGmm",
           function(object, theta){
               getMethod("getRestrict", "rnonlinearGmm")(object, theta)
@@ -607,7 +716,12 @@
               getMethod("coef","rnonlinearGmm")(object, theta)
           )
 
+setMethod("coef", "rformulaGmm",
+          function(object, theta)
+              getMethod("coef","rnonlinearGmm")(object, theta)
+          )
 
+
 ## Subsetting '['
 
 setMethod("[", c("rfunctionGmm", "numeric", "missing"),
@@ -658,7 +772,26 @@
                   }
           })
 
+setMethod("gmmFit", signature("rformulaGmm"), valueClass="gmmfit", 
+          definition = function(object, type=c("twostep", "iter","cue", "onestep"),
+              itertol=1e-7, initW=c("ident", "tsls"), weights="optimal", 
+              itermaxit=100, efficientWeights=FALSE, start=NULL, ...) {
+              cst <- object at cstSpec
+              if (cst$k==0)
+                  {
+                      theta <- coef(object, numeric())
+                      object <- as(object, "formulaGmm")                      
+                      if (class(weights)=="gmmWeights")
+                          wObj <- weights
+                      else
+                          wObj <- evalWeights(object, theta=theta, w=weights)
+                      return(evalGmm(object, theta, wObj))
+                  } else {
+                      callNextMethod()
+                  }
+          })
 
+
 ### momentStrength
 ### For now, there is no measure of moment strength in restricted models
 ### Have to figure out how to identify exluded instruments after

Modified: pkg/gmm4/R/summaryGmm-methods.R
===================================================================
--- pkg/gmm4/R/summaryGmm-methods.R	2018-06-28 21:10:11 UTC (rev 127)
+++ pkg/gmm4/R/summaryGmm-methods.R	2018-09-11 20:35:04 UTC (rev 128)
@@ -5,9 +5,12 @@
           function(x, digits=5, ...)
           {
               print(x at model)
-              ntype <- matrix(c("Two-Step GMM", "Iterated GMM", "CUE", 
-                                "One-Step GMM with fixed weights", "Two-Stage Least Squares", 
-                                "twostep", "iter", "cue", "onestep", "tsls"), ncol = 2)
+              ntype <- matrix(c("Two-Step GMM", "Iterated GMM", "CUE",
+                                "One-Step GMM with fixed weights","Two-Stage Least Squares",
+                                "Evaluated at a fixed Theta; No estimation",
+                                "One-Step Efficient M.D.E.",
+                                "twostep","iter","cue","onestep","tsls", "eval","mde"),
+                              ncol=2)              
               type <- ntype[match(x at type, ntype[, 2]), 1]
               spec <- modelDims(x at model)
               if (spec$q == spec$k) 
[TRUNCATED]

To get the complete diff run:
    svnlook diff /svnroot/gmm -r 128


More information about the Gmm-commits mailing list