[Gmm-commits] r170 - in pkg/momentfit: . R man vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu May 21 21:16:48 CEST 2020
Author: chaussep
Date: 2020-05-21 21:16:47 +0200 (Thu, 21 May 2020)
New Revision: 170
Added:
pkg/momentfit/man/setCoef-methods.Rd
Modified:
pkg/momentfit/DESCRIPTION
pkg/momentfit/NAMESPACE
pkg/momentfit/R/allClasses.R
pkg/momentfit/R/momentModel-methods.R
pkg/momentfit/R/rModel-methods.R
pkg/momentfit/R/rsysMomentModel-methods.R
pkg/momentfit/R/sysMomentModel-methods.R
pkg/momentfit/man/Dresiduals-methods.Rd
pkg/momentfit/man/coef-methods.Rd
pkg/momentfit/man/evalDMoment-methods.Rd
pkg/momentfit/man/evalGmm-methods.Rd
pkg/momentfit/man/evalMoment-methods.Rd
pkg/momentfit/man/getRestrict-methods.Rd
pkg/momentfit/man/model.matrix-methods.Rd
pkg/momentfit/man/modelDims-methods.Rd
pkg/momentfit/man/print-methods.Rd
pkg/momentfit/man/printRestrict-methods.Rd
pkg/momentfit/man/residuals-methods.Rd
pkg/momentfit/man/restModel-methods.Rd
pkg/momentfit/man/subsetting.Rd
pkg/momentfit/vignettes/gmmS4.Rnw
pkg/momentfit/vignettes/gmmS4.pdf
Log:
start to add nonlinear systems of equations
Modified: pkg/momentfit/DESCRIPTION
===================================================================
--- pkg/momentfit/DESCRIPTION 2020-05-08 13:57:42 UTC (rev 169)
+++ pkg/momentfit/DESCRIPTION 2020-05-21 19:16:47 UTC (rev 170)
@@ -1,6 +1,6 @@
Package: momentfit
Version: 0.1-1
-Date: 2020-02-20
+Date: 2020-05-21
Title: Methods of Moments
Author: Pierre Chausse <pchausse at uwaterloo.ca>
Maintainer: Pierre Chausse <pchausse at uwaterloo.ca>
Modified: pkg/momentfit/NAMESPACE
===================================================================
--- pkg/momentfit/NAMESPACE 2020-05-08 13:57:42 UTC (rev 169)
+++ pkg/momentfit/NAMESPACE 2020-05-21 19:16:47 UTC (rev 170)
@@ -42,7 +42,7 @@
printRestrict, restModel, getRestrict, gmm4, sysMomentModel, ThreeSLS,
rhoET, rhoEL, rhoEEL, rhoHD, Wu_lam, EEL_lam, REEL_lam, getLambda,
solveGel, rhoETEL, rhoETHD, ETXX_lam, gelFit, evalGel, getImpProb,
- evalGelObj, momFct, gel4)
+ evalGelObj, momFct, gel4, setCoef)
### S3 methods ###
Modified: pkg/momentfit/R/allClasses.R
===================================================================
--- pkg/momentfit/R/allClasses.R 2020-05-08 13:57:42 UTC (rev 169)
+++ pkg/momentfit/R/allClasses.R 2020-05-21 19:16:47 UTC (rev 170)
@@ -320,6 +320,51 @@
centeredVcov=from at centeredVcov, data=dat)
})
+setAs("slinearModel", "snonlinearModel",
+ function(from) {
+ spec <- modelDims(from)
+ X <- model.matrix(from)
+ theta0 <- rep(1,sum(spec$k))
+ names(theta0) <- paste("theta", 1:sum(spec$k), sep="")
+ eqNames <- paste("Eqn", 1:length(X), sep="")
+ xk <- c(0,cumsum(from at k))
+ theta0 <- lapply(1:length(X), function(i) theta0[(1+xk[i]):(xk[i+1])])
+ parNames <- lapply(theta0, names)
+ rhs <- lapply(1:length(X), function(i){
+ n <- paste("*", colnames(X[[i]]), sep="")
+ n[n=="*(Intercept)"] <- ""
+ n <- paste(names(theta0[[i]]), n, sep="")
+ parse(text=paste(n, collapse="+", sep=""))
+ })
+ lhs <- lapply(1:length(X), function(i)
+ parse(text=from at modelT[[i]][[2]]))
+ varNames <- lapply(1:length(lhs), function(i) {
+ v1 <- all.vars(lhs[[i]])
+ v2 <- all.vars(rhs[[i]])
+ v2 <- v2[!(v2%in%names(theta0[[i]]))]
+ c(v1,v2)})
+
+ Y <- do.call(cbind, modelResponse(from))
+ colnames(Y) <- sapply(lhs, all.vars)
+ X <- do.call(cbind, X)
+ X <- X[,!duplicated(colnames(X))]
+ X <- X[,colnames(X)!="(Intercept)"]
+ Z <- do.call(cbind, model.matrix(from, type="instruments"))
+ Z <- Z[,!duplicated(colnames(Z))]
+ Z <- Z[,colnames(Z) != "(Intercept)"]
+ dat <- cbind(X, Y[,!(colnames(Y) %in% colnames(X))])
+ dat <- cbind(dat, Z[,!(colnames(Z)%in%colnames(dat))])
+ new("snonlinearModel", data=as.data.frame(dat), instT=from at instT,
+ vcov=from at vcov, theta0=theta0, n=spec$n, q=spec$q,k=spec$k,
+ parNames=parNames, momNames=from at momNames, fRHS=rhs,
+ fLHS=lhs, eqnNames=eqNames, vcovOptions=from at vcovOptions,
+ centeredVcov=from at centeredVcov, sameMom=from at sameMom,
+ SUR=from at SUR, varNames=varNames, isEndo=from at isEndo,
+ omit=from at omit, survOptions=from at survOptions,
+ sSpec=from at sSpec, smooth=from at smooth)
+ })
+
+
setAs("sysMomentWeights", "momentWeights",
function(from) {
w <- quadra(from)
Modified: pkg/momentfit/R/momentModel-methods.R
===================================================================
--- pkg/momentfit/R/momentModel-methods.R 2020-05-08 13:57:42 UTC (rev 169)
+++ pkg/momentfit/R/momentModel-methods.R 2020-05-21 19:16:47 UTC (rev 170)
@@ -64,7 +64,33 @@
setMethod("show", "momentModel", function(object) print(object))
+### setCoef #######
+### A generic for setting and validating the format of the coefficient
+setGeneric("setCoef", function(model, ...) standardGeneric("setCoef"))
+
+setMethod("setCoef", "momentModel",
+ function(model, theta) {
+ spec <- modelDims(model)
+ if (length(theta) != length(spec$parNames))
+ stop(paste("Wrong number of coefficients (should be ",
+ spec$k, ")", sep=""))
+ if (!is.null(names(theta)))
+ {
+ chk <- !(names(theta)%in%spec$parNames)
+ if (any(chk))
+ {
+ mess <- paste("The following theta names are invalid: ",
+ paste(names(theta)[chk], collapse=", ", sep=""),
+ sep="")
+ stop(mess)
+ }
+ theta <- theta[match(spec$parNames, names(theta))]
+ } else {
+ names(theta) <- spec$parNames
+ }
+ theta})
+
##### coef ########
### For this, it only attach the names to theta
@@ -142,14 +168,7 @@
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")
+ theta <- setCoef(object, theta)
varList <- c(as.list(theta), as.list(object at modelF))
if (!is.null(res$fLHS))
{
@@ -181,7 +200,7 @@
setMethod("evalMoment", signature("functionModel"),
function(object, theta) {
- theta <- coef(object, theta)
+ theta <- setCoef(object, theta)
gt <- object at fct(theta, object at X)
if (!is.null(sub <- attr(object at X, "subset")))
gt <- gt[,sub]
@@ -193,14 +212,7 @@
setMethod("evalMoment", signature("formulaModel"),
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")
+ theta <- setCoef(object, theta)
varList <- c(as.list(theta), as.list(object at modelF))
gt <- sapply(1:res$q, function(i) {
if (!is.null(res$fLHS[[i]]))
@@ -237,19 +249,11 @@
setMethod("Dresiduals", signature("nonlinearModel"),
function(object, theta) {
- theta <- coef(object, theta)
+ theta <- setCoef(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))
De <- numeric()
- for (i in nt)
+ for (i in names(theta))
{
if (!is.null(res$fLHS))
d <- eval(D(res$fLHS, i), varList)
@@ -257,6 +261,7 @@
d <- 0
De <- cbind(De, d-matrix(eval(D(res$fRHS, i), varList),res$n,1))
}
+ colnames(De) <- object at parNames
De
})
@@ -349,6 +354,7 @@
function(object, theta, impProb=NULL, lambda=NULL)
{
spec <- modelDims(object)
+ theta <- setCoef(object, theta)
if (object at smooth && !is.null(object at dfct))
{
object at dfct <- NULL
@@ -405,18 +411,10 @@
setMethod("evalDMoment", signature("formulaModel"),
function(object, theta, impProb=NULL, lambda=NULL)
{
- theta <- coef(object, theta)
+ theta <- setCoef(object, theta)
spec <- modelDims(object)
- nt <- names(theta)
- nt0 <- names(spec$theta0)
if (is.null(impProb))
impProb <- 1/spec$n
- 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))
if (!is.null(lambda))
{
@@ -466,7 +464,7 @@
nG <- list(spec$momNames, spec$parNames)
}
G <- numeric()
- for (i in nt)
+ for (i in names(theta))
{
lhs <- f(i, lambda, spec$fLHS)
rhs <- f(i, lambda, spec$fRHS)
@@ -1062,16 +1060,7 @@
Call <- try(match.call(call=sys.call(sys.parent())), silent=TRUE)
if (inherits(Call,"try-error"))
Call <- NULL
- if (!is.null(names(theta)))
- {
- if (!all(names(theta) %in% spec$parNames))
- stop("You provided a named theta with wrong names")
- theta <- theta[match(spec$parNames, names(theta))]
- } else {
- if (class(model) %in% c("formulaModel","nonlinearModel"))
- stop("To evaluate nonlinear models, theta must be named")
- names(theta) <- spec$parNames
- }
+ theta <- setCoef(model, theta)
if (is.null(wObj))
wObj <- evalWeights(model, theta)
new("gmmfit", theta=theta, convergence=NULL, convIter=NULL,
@@ -1333,16 +1322,7 @@
if (inherits(Call,"try-error"))
Call <- NULL
spec <- modelDims(model)
- if (!is.null(names(theta)))
- {
- if (!all(names(theta) %in% spec$parNames))
- stop("You provided a named theta with wrong names")
- theta <- theta[match(spec$parNames, names(theta))]
- } else {
- if (class(model) %in% c("formulaGel","nonlinearGel", "formulaGel"))
- stop("To evaluate nonlinear models, theta must be named")
- names(theta) <- spec$parNames
- }
+ theta <- setCoef(model, theta)
type <- paste("Eval-", gelType, sep="")
if (is.null(lambda))
{
Modified: pkg/momentfit/R/rModel-methods.R
===================================================================
--- pkg/momentfit/R/rModel-methods.R 2020-05-08 13:57:42 UTC (rev 169)
+++ pkg/momentfit/R/rModel-methods.R 2020-05-21 19:16:47 UTC (rev 170)
@@ -358,7 +358,7 @@
n <- length(hypothesis)
k <- length(cnames)
# an attempt to rename all special variable names (from transformed I() e.g. or
- # interection :.
+ # interaction :.
newN <- paste("theta", 1:k, sep="")
tmp <- cnames
hasI <- grepl("I(", cnames, fixed=TRUE)
@@ -613,6 +613,7 @@
setMethod("getRestrict", "rlinearModel",
function(object, theta) {
+ theta <- setCoef(as(object, "linearModel"), theta)
R <- c(object at cstLHS%*%theta)
cst <- .printHypothesis(object at cstLHS, object at cstRHS, object at parNames)
list(dR=object at cstLHS, R=R, q=object at cstRHS, hypo=cst,
@@ -623,6 +624,7 @@
function(object, theta) {
dR <-numeric()
R <- numeric()
+ theta <- setCoef(as(object, "nonlinearModel"), theta)
for (r in object at R)
{
dlhs <- sapply(object at parNames, function(pn)
@@ -685,18 +687,7 @@
function(object, theta)
{
spec <- modelDims(object)
- if (length(theta)>0)
- {
- if (is.null(names(theta)))
- {
- if (length(theta)!=length(spec$parNames))
- stop("Wrong number of coefficients")
- names(theta) <- spec$parNames
- } else {
- if (!all(names(theta)%in%spec$parNames))
- stop("theta has wrong names")
- }
- }
+ theta <- setCoef(object, theta)
theta2 <- rep(0,object at k)
names(theta2) <- object at parNames
theta2[names(theta)] <- theta
Modified: pkg/momentfit/R/rsysMomentModel-methods.R
===================================================================
--- pkg/momentfit/R/rsysMomentModel-methods.R 2020-05-08 13:57:42 UTC (rev 169)
+++ pkg/momentfit/R/rsysMomentModel-methods.R 2020-05-21 19:16:47 UTC (rev 170)
@@ -1,6 +1,69 @@
######## Methods for restricted System of Equations
######################################################
+## Hidden function to arrange the restrictions
+
+.imposeSNLRestrict <- function(R, object)
+{
+ spec <- modelDims(object)
+ allPar <- do.call("c", c(list(),spec$parNames))
+ theta0 <- do.call("c", c(list(),spec$theta0))
+ fRHS <- object at fRHS
+ fLHS <- object at fLHS
+ chk <- sapply(R, function(r) all(all.vars(r) %in% allPar))
+ chk2 <- sapply(R, function(r) sum(sapply(spec$parNames,
+ function(ni) any(all.vars(r) %in% ni))))
+ crossRest <- any(chk2>1)
+ if (!all(chk))
+ stop("Wrong coefficient names in some of the restrictions")
+ rest <- sapply(R, function(r) as.character(r[[2]]))
+ 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(allPar, function(pn)
+ eval(D(r[[2]], pn), as.list(theta0)))
+ rhs <- sapply(allPar, function(pn)
+ eval(D(r[[3]], pn), as.list(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 <- lapply(fRHS, function(ri) as.character(ri))
+ lhs <- lapply(fLHS, function(ri) ifelse(is.null(ri),NULL,as.character(ri)))
+ for (r in R)
+ {
+ rhs <- gsub(as.character(r[2]), paste("(", as.character(r[3]),
+ ")", sep=""), rhs)
+ if (!is.null(lhs))
+ lhs <- gsub(as.character(r[2]),
+ paste("(", as.character(r[3]),
+ ")", sep=""), lhs)
+ }
+ rhs <- lapply(rhs, function(ri) parse(text=ri))
+ lhs <- lapply(lhs, function(ri) parse(text=ri))
+ k <- sum(spec$k)-length(R)
+ parNames <- if (is.list(spec$parNames))
+ lapply(spec$parNames, function(pi) pi[!(pi %in% rest)])
+ else
+ spec$parNames[!(spec$parNames %in% rest)]
+ theta0 <- if (is.list(spec$theta0))
+ lapply(spec$theta0, function(ti) ti[!(names(ti) %in% rest)])
+ else
+ spec$theta0[!(names(spec$theta0) %in% rest)]
+ if (length(rhs) == 1)
+ {
+ rhs <- rhs[[1]]
+ lhs <- lhs[[1]]
+ }
+ list(rhs=rhs, lhs=lhs, parNames=parNames,
+ theta0=theta0, k=k, crossEquRest=crossRest)
+}
+
+
## Constructor of restricted models
setMethod("restModel", signature("slinearModel"),
@@ -7,9 +70,12 @@
function(object, R, rhs=NULL)
{
parNames <- paste(rep(object at eqnNames, object at k), ".",
- do.call("c", object at parNames), sep = "")
+ do.call("c", object at parNames), sep = "")
+
if (is.character(R))
{
+ parNames <- gsub("\\(Intercept\\)", "Intercept", parNames)
+ R <- gsub("\\(Intercept\\)", "Intercept", R)
res <- .makeHypothesis(parNames, R, rhs)
R <- res$R
rhs <- res$rhs
@@ -75,16 +141,46 @@
cstSpec=res, object)
})
+
+setMethod("restModel", signature("snonlinearModel"),
+ 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(!inherits(R,"formula"))
+ stop("R must be a formula or a list of formulas")
+ R <- list(R)
+ } else {
+ chk <- sapply(R, function(r) inherits(r,"formula"))
+ if (!all(chk))
+ stop("R must be a formula, a list of formulas or a vector of characters")
+ }
+ }
+ res <- .imposeSNLRestrict(R, object)
+ cstSpec <- list(newParNames = res$parNames,
+ originParNames=object at parNames,
+ k=res$k, theta0=res$theta0, fRHS=res$rhs, fLHS=res$lhs,
+ crossEquRest=res$crossEquRest)
+ new("rsnonlinearModel", R=R, cstSpec=cstSpec, object)
+ })
+
## coef
setMethod("coef","rslinearModel",
function (object, theta)
- {
+ {
spec <- modelDims(object)
- if (!is.list(theta))
- stop("theta must be a list")
- if (length(theta) != length(spec$eqnNames))
- stop("Wrong number of coefficients")
+ theta <- setCoef(object, theta)
if (!object at cstSpec$crossEquRest)
{
tet <- lapply(1:length(theta), function(i)
@@ -97,6 +193,24 @@
theta
})
+setMethod("coef", "rsnonlinearModel",
+ function(object, theta)
+ {
+ spec <- modelDims(object)
+ theta <- setCoef(object, theta)
+ names(theta) <- NULL
+ theta <- do.call("c", theta)
+ theta2 <- rep(0,sum(object at k))
+ names(theta2) <- do.call("c", object at parNames)
+ theta2[names(theta)] <- theta
+ chk <- sapply(object at R, function(r) is.numeric(r[[3]]))
+ for (r in object at R[chk])
+ theta2[as.character(r[[2]])] <- r[[3]]
+ for (r in object at R[!chk])
+ theta2[as.character(r[[2]])] <- eval(r[[3]], as.list(theta2))
+ .tetReshape(theta2, spec$eqnNames, object at parNames)
+ })
+
## modelDims
setMethod("modelDims", "rslinearModel",
@@ -107,6 +221,19 @@
isEndo=res$isEndo)
})
+setMethod("modelDims", "rsnonlinearModel",
+ function(object) {
+ cst <- object at cstSpec
+ k <- cst$k
+ parNames <- cst$newParNames
+ theta0 <- cst$theta0
+ fRHS <- cst$fRHS
+ fLHS <- cst$fLHS
+ list(k=k, q=object at q, n=object at n, parNames=parNames,
+ momNames=object at momNames, theta0=theta0,
+ fRHS=fRHS, fLHS=fLHS, eqnNames=object at eqnNames)
+ })
+
## printRestrict
setMethod("printRestrict", "rslinearModel",
@@ -119,6 +246,37 @@
cat("\t", cst[i], "\n")
})
+setMethod("printRestrict", "rsnonlinearModel",
+ function(object){
+ cat("Constraints:\n")
+ parNames <- object at parNames
+ eqn <- object at eqnNames
+ chk <- lapply(object at R, function(ri)
+ which(sapply(parNames, function(pi) any(all.vars(ri) %in% pi))))
+ for (i in unique(chk))
+ {
+ chk2 <- which(sapply(chk, function(ci) isTRUE(all.equal(ci,i))))
+ nr <- length(i)
+ if (nr>1)
+ {
+ cat("\tInvolving equations: ")
+ if (nr>2)
+ sep <- c(rep(", ", nr-2), " and ", "")
+ else
+ sep <- c(" and ", "")
+ cat(paste(eqn[i], sep, collapse="", sep=""), "\n")
+ } else {
+ cat("\tInvolving equation: ", eqn[i], "\n", sep="")
+ }
+ for (ri in object at R[chk2])
+ {
+ cat("\t\t")
+ print(ri)
+ }
+ }
+ cat("\n")
+ })
+
## print
setMethod("print", "rslinearModel",
@@ -142,6 +300,14 @@
}
})
+setMethod("print", "rsnonlinearModel",
+ function (x, ...) {
+ print(as(x, "snonlinearModel"))
+ printRestrict(x)
+ cat("\n")
+ } )
+
+
## getRestrict
setMethod("getRestrict", "sysModel",
@@ -152,6 +318,7 @@
setMethod("getRestrict", "rslinearModel",
function(object, theta) {
+ theta <- setCoef(as(object, "slinearModel"), theta)
theta <- do.call("c", theta)
R <- c(object at cstLHS%*%theta)
parNames <- paste(rep(object at eqnNames, object at k), ".",
@@ -161,6 +328,35 @@
orig.R=object at cstLHS, orig.rhs=object at cstRHS)
})
+setMethod("getRestrict", "rsnonlinearModel",
+ function(object, theta) {
+ theta <- setCoef(as(object, "snonlinearModel"), theta)
+ names(theta) <- NULL
+ theta <- do.call("c", theta)
+ dR <-numeric()
+ R <- numeric()
+ parNames <- names(theta)
+ for (r in object at R)
+ {
+ dlhs <- sapply(parNames, function(pn)
+ eval(D(r[[2]], pn), as.list(theta)))
+ drhs <- sapply(parNames, function(pn)
+ eval(D(r[[3]], pn), as.list(theta)))
+ dR <- rbind(dR, dlhs-drhs)
+ lhs <- eval(r[[2]], as.list(theta))
+ rhs <- eval(r[[3]], as.list(theta))
+ R <- c(R, lhs-rhs)
+ }
+ if (any(is.na(c(R,dR))) || any(!is.finite(c(dR,R))))
+ stop("Some values in R or dR at theta are either infinite or NAN")
+ if (qr(dR)$rank < length(R))
+ stop("The matrix of derivatives of the constraints is not full rank")
+ hypo <- sapply(object at R, function(r) capture.output(print(r)))
+ list(dR=dR, R=R, q=rep(0, nrow(dR)), hypo=hypo,
+ orig.R=object at R, orig.rhs=NULL)
+ })
+
+
### Subset
setMethod("[", c("rslinearModel", "numeric", "missing"),
@@ -189,6 +385,26 @@
restModel(m, R, rhs)
})
+
+setMethod("[", c("rsnonlinearModel", "numeric", "missing"),
+ function(x, i, j)
+ {
+ i <- as.integer(i)
+ if (x at cstSpec$crossEquRest)
+ stop("Cannot select an equation when the system has cross-equation restrictions")
+ if (length(i) > 1)
+ {
+ warning("You can only select one equation, the first element of i is used")
+ i <- i[1]
+ }
+ m <- as(x, "snonlinearModel")[i]
+ chk <- sapply(x at R, function(ri) any(all.vars(ri) %in% m at parNames))
+ R <- x at R[chk]
+ if (length(R)==0)
+ return(m)
+ restModel(m, R)
+ })
+
### model.matrix
setMethod("model.matrix", "rslinearModel",
@@ -218,20 +434,23 @@
## evalMoment
-setMethod("evalMoment", "rslinearModel",
+setMethod("evalMoment", "rsysModel",
function (object, theta)
{
theta <- coef(object, theta)
- evalMoment(as(object, "slinearModel"), theta)
+ evalMoment(as(object, "sysModel"), theta)
})
+
## evalDMoment
setMethod("evalDMoment","rslinearModel",
function (object, theta)
{
+ chk <- ifelse(is.null(object at cstSpec$crossEquRest), FALSE,
+ object at cstSpec$crossEquRest)
neqn <- length(object at eqnNames)
- if (!object at cstSpec$crossEquRest)
+ if (!chk)
dgt <- lapply(1:neqn, function(i) evalDMoment(object[i]))
else
dgt <- list(neqn*evalDMoment(as(object, "rlinearModel")))
@@ -239,15 +458,86 @@
dgt
})
+setMethod("evalDMoment", signature("rsnonlinearModel"),
+ function(object, theta, impProb=NULL)
+ {
+ De <- Dresiduals(object, theta)
+ Z <- model.matrix(as(object, "snonlinearModel"), "instrument")
+ spec <- modelDims(object)
+ if (is.null(impProb))
+ impProb <- 1/spec$n
+ sG <- lapply(1:length(spec$eqnNames), function(eq) {
+ G <- apply(De[[eq]],2, function(x)
+ {
+ tmp <- Z[[eq]]*x
+ if (object at smooth)
+ tmp <- stats::kernapply(tmp, object at sSpec@w)
+ colSums(tmp*impProb)
+ })
+ G <- as.matrix(G)
+ dimnames(G) <- list(spec$momNames[[eq]], colnames(De[[eq]]))
+ G
+ })
+ names(sG) <- spec$eqnNames
+ sG
+ })
+
## residuals
-setMethod("residuals", "rslinearModel",
+setMethod("residuals", "rsysModel",
function (object, theta)
{
theta <- coef(object, theta)
- residuals(as(object, "slinearModel"), theta)
+ residuals(as(object, "sysModel"), theta)
})
+## model.matrix
+
+setMethod("model.matrix", "rsnonlinearModel",
+ function(object, type = c("regressors", "instruments"))
+ {
+ type <- match.arg(type)
+ model.matrix(as(object, "snonlinearModel"), type)
+ })
+
+
+## Dresiduals
+
+setMethod("Dresiduals", signature("rsnonlinearModel"),
+ function(object, theta) {
+ spec <- modelDims(object)
+ theta <- setCoef(object, theta)
+ if (!object at cstSpec$crossEquRest)
+ {
+ De <- Dresiduals(as(object, "snonlinearModel"), coef(object, theta))
+ for (i in 1:length(De))
+ De[[i]] <- De[[i]][,colnames(De[[i]]) %in% spec$parNames[[i]]]
+ return(De)
+ }
+ neqn <- length(object at eqnNames)
+ fLHS <- spec$fLHS
+ fRHS <- spec$fRHS
+ names(theta) <- NULL
+ theta <- do.call("c", theta)
+ nt <- names(theta)
+ varList <- c(as.list(theta), as.list(object at data))
+ sDe <- lapply(1:neqn, function(eq){
+ De <- numeric()
+ for (i in nt)
+ {
+ if (!is.null(fLHS[[eq]]))
+ d <- eval(D(fLHS[[eq]], i), varList)
+ else
+ d <- 0
+ De <- cbind(De, d-matrix(eval(D(fRHS[[eq]], i), varList),spec$n,1))
+ }
+ colnames(De) <- nt
+ De
+ })
+ names(sDe) <- spec$eqnNames
+ sDe
+ })
+
## solveGmm
setMethod("solveGmm", c("rslinearModel","sysMomentWeights"),
Modified: pkg/momentfit/R/sysMomentModel-methods.R
===================================================================
--- pkg/momentfit/R/sysMomentModel-methods.R 2020-05-08 13:57:42 UTC (rev 169)
+++ pkg/momentfit/R/sysMomentModel-methods.R 2020-05-21 19:16:47 UTC (rev 170)
@@ -48,6 +48,72 @@
fRHS=object at fRHS, fLHS=object at fLHS, eqnNames=object at eqnNames,
isEndo=object at isEndo)
})
+
+## setCoef
+## Used to validate and format the coefficient into a named list
+
+setMethod("setCoef", signature("sysModel"),
+ function(model, theta) {
+ spec <- modelDims(model)
+ k <- sum(spec$k)
+ if (!is.list(theta))
+ {
+ if (length(theta) != k)
+ stop(paste("Wrong number of coefficients (should be ",
+ k, ")", sep=""))
+ if (!is.null(names(theta)))
+ {
+ parNames <- do.call("c", spec$parNames)
+ if (any(duplicated(names(theta))))
+ stop("Cannot have more than one coefficient with the same name")
+ chk <- !(names(theta) %in% parNames)
+ if (any(chk))
+ {
+ mes <- "The following coefficients have invalid name: "
+ mes <- paste(mes, paste(names(theta)[chk], collapse=", ",
+ sep=""), sep="")
+ stop(mes)
+ }
+ theta <- theta[match(parNames, names(theta))]
+ }
+ theta <- .tetReshape(theta, spec$eqnNames, spec$parNames)
+ } else {
+ if (length(theta) != length(spec$eqnNames))
+ stop("Wrong number of equations")
+ if (is.null(names(theta)))
+ {
+ names(theta) <- spec$eqnNames
+ } else {
+ if(!all(names(theta) %in% spec$eqnNames))
+ stop("Wrong equation names")
+ theta <- theta[match(spec$eqnNames,names(theta))]
+ }
+ chk <- sapply(1:length(theta), function(i)
+ length(theta[[i]]) == length(spec$parNames[[i]]))
+ if (!all(chk))
+ {
+ mes <- "Wrong number of coefficients in the following equation: "
+ mes <- paste(mes, paste(names(theta)[!chk], collapse=", ", sep=""),
+ sep="")
+ stop(mes)
+ }
+ theta <- lapply(1:length(theta), function(i) {
+ ti <- theta[[i]]
+ if (is.null(names(ti)))
+ {
+ names(ti) <- spec$parNames[[i]]
+ } else {
+ if (!all(names(ti) %in% spec$parNames[[i]]))
+ stop(paste("Wrong coefficient names in equation ",
+ names(theta)[i], sep=""))
+ ti <- ti[match(spec$parNames[[i]], names(ti))]
+ }
+ ti})
+ names(theta) <- spec$eqnNames
+ }
+ theta
+ })
+
## Subsetting '['
setMethod("[", c("sysModel", "missing", "list"),
@@ -478,15 +544,32 @@
setMethod("evalGmmObj", signature("sysModel", "list", "sysMomentWeights"),
function(object, theta, wObj, ...)
- {
- gt <- evalMoment(object, theta)
- gt <- lapply(gt, function(g) colMeans(g))
- gt <- do.call("c", gt)
- n <- object at n
- obj <- quadra(wObj, gt)
- n*obj
- })
+ {
+ gt <- evalMoment(object, theta)
+ gt <- lapply(gt, function(g) colMeans(g))
+ gt <- do.call("c", gt)
+ n <- object at n
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/gmm -r 170
More information about the Gmm-commits
mailing list