[Gmm-commits] r152 - in pkg/gmm4: . R man vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Fri Nov 15 20:56:25 CET 2019
Author: chaussep
Date: 2019-11-15 20:56:25 +0100 (Fri, 15 Nov 2019)
New Revision: 152
Added:
pkg/gmm4/R/gel4.R
pkg/gmm4/man/gel4.Rd
pkg/gmm4/man/mconfint-class.Rd
pkg/gmm4/man/plot-methods.Rd
Modified:
pkg/gmm4/DESCRIPTION
pkg/gmm4/NAMESPACE
pkg/gmm4/R/allClasses.R
pkg/gmm4/R/gel.R
pkg/gmm4/R/gelModels-methods.R
pkg/gmm4/R/gelfit-methods.R
pkg/gmm4/R/rGelModel-methods.R
pkg/gmm4/man/confint-methods.Rd
pkg/gmm4/man/gelModel.Rd
pkg/gmm4/man/modelFit-methods.Rd
pkg/gmm4/man/print-methods.Rd
pkg/gmm4/man/show-methods.Rd
pkg/gmm4/vignettes/empir.bib
pkg/gmm4/vignettes/gelS4.Rnw
pkg/gmm4/vignettes/gelS4.pdf
Log:
worked on the vignettes and added some features for GEL
Modified: pkg/gmm4/DESCRIPTION
===================================================================
--- pkg/gmm4/DESCRIPTION 2019-11-06 22:42:07 UTC (rev 151)
+++ pkg/gmm4/DESCRIPTION 2019-11-15 19:56:25 UTC (rev 152)
@@ -1,12 +1,12 @@
Package: gmm4
Version: 0.1-0
-Date: 2019-11-01
+Date: 2019-11-15
Title: S4 Generalized Method of Moments
Author: Pierre Chausse <pchausse at uwaterloo.ca>
Maintainer: Pierre Chausse <pchausse at uwaterloo.ca>
Description: This is a complete restructured version of the 'gmm' package (Chausse 2010; <doi:10.18637/jss.v034.i11>) using 'S4' only type of classes and methods. It provides tools for estimating single equations and system of equations using the Generalized Method of Moments (Hansen 1982; <doi:10.2307/1912775>). It is in a very early stage and suggestions are welcome. See the vignette for more details.
Depends: R (>= 3.0.0), sandwich
-Imports: stats, methods
+Imports: stats, methods, parallel
Suggests: lmtest, knitr, texreg
Collate: 'allClasses.R' 'validity.R' 'gmmData.R' 'gmmModels-methods.R'
'gmmfit-methods.R' 'tsls-methods.R' 'specTest-methods.R'
@@ -14,7 +14,7 @@
'rGmmModel-methods.R' 'hypothesisTest-methods.R'
'sysGmmModel.R' 'sysGmmModels-methods.R' 'rsysGmmModels-methods.R'
'sgmmfit-methods.R' 'gmm4.R' 'gel.R' 'gelModels-methods.R'
- 'rGelModel-methods.R' 'gelfit-methods.R'
+ 'rGelModel-methods.R' 'gelfit-methods.R' 'gel4.R'
License: GPL (>= 2)
NeedsCompilation: yes
VignetteBuilder: knitr
Modified: pkg/gmm4/NAMESPACE
===================================================================
--- pkg/gmm4/NAMESPACE 2019-11-06 22:42:07 UTC (rev 151)
+++ pkg/gmm4/NAMESPACE 2019-11-15 19:56:25 UTC (rev 152)
@@ -3,6 +3,10 @@
"getClassDef", "selectMethod", "callNextMethod", "as", "setAs",
"getMethod")
+importFrom("parallel", mclapply)
+
+importFrom("graphics", plot, polygon, grid)
+
importFrom("utils", capture.output)
importFrom("stats", "ar", "as.formula", "model.matrix","vcov",
@@ -25,10 +29,10 @@
"sgmmfit","stsls", "rslinearGmm", "rsnonlinearGmm", "rsysGmmModels",
"formulaGmm","rfunctionGmm", "gelfit", "summaryGel", "confint",
"rlinearGel", "nonlinearGel", "rfunctionGel", "rformulaGel",
- "rgelModels","callORNULL")
+ "rgelModels","callORNULL", "mconfint")
exportMethods(residuals, print, show, vcovHAC, coef, vcov, bread, summary, update,
model.matrix, hypothesisTest, "[", merge, subset, confint, gmmToGel,
- kernapply)
+ kernapply, plot)
export(gmmModel, evalMoment, Dresiduals, evalDMoment, momentVcov, estfun.gmmFct,
evalWeights, quadra, evalObjective, solveGmm, momentStrength,evalModel,
@@ -35,7 +39,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)
+ solveGel, getImpProb, rhoETEL, rhoETHD, gel4)
### S3 methods ###
Modified: pkg/gmm4/R/allClasses.R
===================================================================
--- pkg/gmm4/R/allClasses.R 2019-11-06 22:42:07 UTC (rev 151)
+++ pkg/gmm4/R/allClasses.R 2019-11-15 19:56:25 UTC (rev 152)
@@ -88,6 +88,10 @@
setClass("confint", representation(interval = "matrix", type="character",
level="numeric"))
+
+setClass("mconfint",
+ representation(areaPoints="matrix", type="character", level="numeric"))
+
## summaryGmm
setClass("summaryGmm", representation(coef="matrix", specTest = "specTest",
Modified: pkg/gmm4/R/gel.R
===================================================================
--- pkg/gmm4/R/gel.R 2019-11-06 22:42:07 UTC (rev 151)
+++ pkg/gmm4/R/gel.R 2019-11-15 19:56:25 UTC (rev 152)
@@ -139,70 +139,75 @@
getLambda <- function (gmat, lambda0=NULL, gelType=NULL, rhoFct=NULL,
tol = 1e-07, maxiter = 100, k = 1, method="BFGS",
algo = c("nlminb", "optim", "Wu"), control = list())
+{
+ if (!is.null(gelType))
{
- algo <- match.arg(algo)
- gmat <- as.matrix(gmat)
- if (is.null(lambda0))
- lambda0 <- rep(0, ncol(gmat))
- if (is.null(rhoFct))
+ if (length(algo) == 3 & gelType == "EL")
+ algo <- "Wu"
+ }
+ algo <- match.arg(algo)
+ gmat <- as.matrix(gmat)
+ if (is.null(lambda0))
+ lambda0 <- rep(0, ncol(gmat))
+ if (is.null(rhoFct))
+ {
+ if (is.null(gelType))
+ stop("Without rhoFct, gelType must be given")
+ rhoFct <- get(paste("rho",gelType,sep=""))
+ } else {
+ gelType <- "Other"
+ }
+ if (algo == "Wu" & gelType != "EL")
+ stop("Wu (2005) algo to compute Lambda is for EL only")
+ if (algo == "Wu")
+ return(Wu_lam(gmat, tol, maxiter, k))
+ if (gelType == "EEL")
+ return(EEL_lam(gmat, k))
+ if (gelType == "REEL")
+ return(REEL_lam(gmat, NULL, maxiter, k))
+ if (gelType %in% c("ETEL", "ETHD"))
+ return(ETXX_lam(gmat, lambda0, k, gelType, algo, method, control))
+
+ fct <- function(l, X, rhoFct, k) {
+ r0 <- rhoFct(X, l, derive = 0, k = k)
+ -mean(r0)
+ }
+ Dfct <- function(l, X, rhoFct, k)
+ {
+ r1 <- rhoFct(X, l, derive = 1, k = k)
+ -colMeans(r1 * X)
+ }
+ DDfct <- function(l, X, rhoFct, k)
+ {
+ r2 <- rhoFct(X, l, derive = 2, k = k)
+ -crossprod(X * r2, X)/nrow(X)
+ }
+ if (algo == "optim") {
+ if (gelType == "EL")
{
- if (is.null(gelType))
- stop("Without rhoFct, gelType must be given")
- rhoFct <- get(paste("rho",gelType,sep=""))
+ ci <- -rep(1, nrow(gmat))
+ res <- constrOptim(lambda0, fct, Dfct, -gmat, ci, control = control,
+ X = gmat, rhoFct = rhoFct, k = k)
+ } else if (gelType == "HD") {
+ ci <- -rep(1, nrow(gmat))
+ res <- constrOptim(lambda0, fct, Dfct, -gmat, ci, control = control,
+ X = gmat, rhoFct = rhoFct, k = k)
} else {
- gelType <- "Other"
+ res <- optim(lambda0, fct, gr = Dfct, X = gmat, rhoFct = rhoFct,
+ k = k, method = method, control = control)
}
- if (algo == "Wu" & gelType != "EL")
- stop("Wu (2005) algo to compute Lambda is for EL only")
- if (algo == "Wu")
- return(Wu_lam(gmat, tol, maxiter, k))
- if (gelType == "EEL")
- return(EEL_lam(gmat, k))
- if (gelType == "REEL")
- return(REEL_lam(gmat, NULL, maxiter, k))
- if (gelType %in% c("ETEL", "ETHD"))
- return(ETXX_lam(gmat, lambda0, k, gelType, algo, method, control))
-
- fct <- function(l, X, rhoFct, k) {
- r0 <- rhoFct(X, l, derive = 0, k = k)
- -mean(r0)
- }
- Dfct <- function(l, X, rhoFct, k)
- {
- r1 <- rhoFct(X, l, derive = 1, k = k)
- -colMeans(r1 * X)
- }
- DDfct <- function(l, X, rhoFct, k)
- {
- r2 <- rhoFct(X, l, derive = 2, k = k)
- -crossprod(X * r2, X)/nrow(X)
- }
- if (algo == "optim") {
- if (gelType == "EL")
- {
- ci <- -rep(1, nrow(gmat))
- res <- constrOptim(lambda0, fct, Dfct, -gmat, ci, control = control,
- X = gmat, rhoFct = rhoFct, k = k)
- } else if (gelType == "HD") {
- ci <- -rep(1, nrow(gmat))
- res <- constrOptim(lambda0, fct, Dfct, -gmat, ci, control = control,
- X = gmat, rhoFct = rhoFct, k = k)
- } else {
- res <- optim(lambda0, fct, gr = Dfct, X = gmat, rhoFct = rhoFct,
- k = k, method = method, control = control)
- }
- } else {
- res <- nlminb(lambda0, fct, gradient = Dfct, hessian = DDfct,
- X = gmat, rhoFct = rhoFct, k = k, control = control)
- }
- lambda0 <- res$par
- if (algo == "optim")
- conv <- list(convergence = res$convergence, counts = res$counts,
- message = res$message)
- else
- conv <- list(convergence = res$convergence, counts = res$evaluations,
- message = res$message)
- return(list(lambda = lambda0, convergence = conv,
- obj= mean(rhoFct(gmat,lambda0,0,k))))
+ } else {
+ res <- nlminb(lambda0, fct, gradient = Dfct, hessian = DDfct,
+ X = gmat, rhoFct = rhoFct, k = k, control = control)
}
+ lambda0 <- res$par
+ if (algo == "optim")
+ conv <- list(convergence = res$convergence, counts = res$counts,
+ message = res$message)
+ else
+ conv <- list(convergence = res$convergence, counts = res$evaluations,
+ message = res$message)
+ return(list(lambda = lambda0, convergence = conv,
+ obj= mean(rhoFct(gmat,lambda0,0,k))))
+}
Added: pkg/gmm4/R/gel4.R
===================================================================
--- pkg/gmm4/R/gel4.R (rev 0)
+++ pkg/gmm4/R/gel4.R 2019-11-15 19:56:25 UTC (rev 152)
@@ -0,0 +1,43 @@
+################### the main gel function ###################
+########## The functions is to avoid having to build model objects
+
+gel4 <- function (g, x=NULL, theta0=NULL,lambda0=NULL, getVcov=FALSE,
+ gelType = c("EL","ET","EEL","HD", "REEL","ETEL","ETHD"),
+ vcov = c("MDS","iid","HAC"), grad=NULL,
+ vcovOptions=list(), centeredVcov = TRUE,
+ cstLHS=NULL, cstRHS=NULL, lamSlv=NULL,
+ rhoFct=NULL, initTheta=c("gmm", "theta0"),
+ data = parent.frame(),
+ coefSlv=c("optim","nlminb","constrOptim"),
+ lControl=list(), tControl=list())
+{
+ Call <- match.call()
+ vcov <- match.arg(vcov)
+ gelType <- match.arg(gelType)
+ coefSlv <- match.arg(coefSlv)
+ initTheta <- match.arg(initTheta)
+ if (initTheta=="theta0" & is.null(theta0))
+ stop("theta0 is required when initTheta='theta0'")
+ model <- gelModel(g, x, gelType, rhoFct, theta0, grad, vcov, vcovOptions,
+ centeredVcov, data=data)
+ if (initTheta == "theta0")
+ {
+ initTheta <- "modelTheta0"
+ names(theta0) = model at parNames
+ } else {
+ theta0 <- NULL
+ }
+ if (!is.null(cstLHS))
+ {
+ model <- restModel(model, cstLHS, cstRHS)
+ spec <- modelDims(model)
+ if (!is.null(theta0))
+ theta0 <- theta0[(names(theta0) %in% spec at parNames)]
+ }
+ fit <- modelFit(object=model, initTheta=initTheta, theta0=theta0,
+ lambda0=lambda0, vcov=getVcov, coefSlv=coefSlv,
+ lamSlv=lamSlv, tControl=tControl, lControl=lControl)
+ fit at call <- Call
+ fit
+}
+
Modified: pkg/gmm4/R/gelModels-methods.R
===================================================================
--- pkg/gmm4/R/gelModels-methods.R 2019-11-06 22:42:07 UTC (rev 151)
+++ pkg/gmm4/R/gelModels-methods.R 2019-11-15 19:56:25 UTC (rev 152)
@@ -115,49 +115,112 @@
function(object, theta0=NULL, lambda0=NULL, lamSlv=NULL,
coefSlv=c("optim","nlminb","constrOptim"),
lControl=list(), tControl=list())
+ {
+ coefSlv <- match.arg(coefSlv)
+ f <- function(theta, model, lambda0, slv, lcont,returnL=FALSE)
{
- coefSlv <- match.arg(coefSlv)
- f <- function(theta, model, lambda0, slv, lcont,returnL=FALSE)
- {
- gt <- evalMoment(model, theta)
- gelt <- model at gelType
- k <- model at wSpec$k
- args <- c(list(gmat=gt, lambda0=lambda0, gelType=gelt$name,
- rhoFct=gelt$fct), lcont, k=k[1]/k[2])
- res <- do.call(slv, args)
- if (returnL)
- return(res)
- res$obj
- }
- if (is.null(lambda0))
- lambda0 <- rep(0, modelDims(object)$q)
- if (is.null(theta0))
- {
- if (!("theta0"%in%slotNames(object)))
- stop("Theta0 must be provided")
- theta0 <- modelDims(object)$theta0
- }
- if (is.null(lamSlv))
- lamSlv <- getLambda
- if (coefSlv == "nlminb")
- args <- c(list(start=theta0, objective=f,
- model=object, lambda0=lambda0,
- slv=lamSlv, lcont=lControl), tControl)
- else
- args <- c(list(par=theta0, fn=f, model=object, lambda0=lambda0,
- slv=lamSlv, lcont=lControl), tControl)
- res <- do.call(get(coefSlv), args)
- resl <- f(res$par, object, lambda0, lamSlv, lControl, TRUE)
- names(resl$lambda) <- modelDims(object)$momNames
- theta <- res$par
- names(theta) <- modelDims(object)$parNames
- list(theta=theta, convergence=res$convergence,
- lambda=resl$lambda, lconvergence=resl$convergence)
+ names(theta) <- modelDims(model)$parNames
+ gt <- evalMoment(model, theta)
+ gelt <- model at gelType
+ k <- model at wSpec$k
+ args <- c(list(gmat=gt, lambda0=lambda0, gelType=gelt$name,
+ rhoFct=gelt$fct), lcont, k=k[1]/k[2])
+ res <- do.call(slv, args)
+ if (returnL)
+ return(res)
+ res$obj
+ }
+ if (is.null(lambda0))
+ lambda0 <- rep(0, modelDims(object)$q)
+ if (is.null(theta0))
+ {
+ if (!("theta0"%in%slotNames(object)))
+ stop("Theta0 must be provided")
+ theta0 <- modelDims(object)$theta0
+ }
+ if (is.null(lamSlv))
+ lamSlv <- getLambda
+ if (coefSlv == "nlminb")
+ args <- c(list(start=theta0, objective=f,
+ model=object, lambda0=lambda0,
+ slv=lamSlv, lcont=lControl), tControl)
+ else
+ args <- c(list(par=theta0, fn=f, model=object, lambda0=lambda0,
+ slv=lamSlv, lcont=lControl), tControl)
+
+ res <- do.call(get(coefSlv), args)
+ resl <- f(res$par, object, lambda0, lamSlv, lControl, TRUE)
+ names(resl$lambda) <- modelDims(object)$momNames
+ theta <- res$par
+ names(theta) <- modelDims(object)$parNames
+ list(theta=theta, convergence=res$convergence,
+ lambda=resl$lambda, lconvergence=resl$convergence)
})
######################### modelFit #########################
+setMethod("modelFit", signature("linearGel"), valueClass="gelfit",
+ definition = function(object, 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")
+ Call <- NULL
+ met <- getMethod("modelFit","gelModels")
+ obj <- met(object, gelType, rhoFct, initTheta, theta0,
+ lambda0, vcov, ...)
+ obj at call <- Call
+ obj
+ })
+setMethod("modelFit", signature("nonlinearGel"), valueClass="gelfit",
+ definition = function(object, 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")
+ Call <- NULL
+ met <- getMethod("modelFit","gelModels")
+ obj <- met(object, gelType, rhoFct, initTheta, theta0,
+ lambda0, vcov, ...)
+ obj at call <- Call
+ obj
+ })
+
+setMethod("modelFit", signature("formulaGel"), valueClass="gelfit",
+ definition = function(object, 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")
+ Call <- NULL
+ met <- getMethod("modelFit","gelModels")
+ obj <- met(object, gelType, rhoFct, initTheta, theta0,
+ lambda0, vcov, ...)
+ obj at call <- Call
+ obj
+ })
+
+setMethod("modelFit", signature("functionGel"), valueClass="gelfit",
+ definition = function(object, 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")
+ Call <- NULL
+ met <- getMethod("modelFit","gelModels")
+ obj <- met(object, gelType, rhoFct, initTheta, theta0,
+ lambda0, vcov, ...)
+ obj at call <- Call
+ obj
+ })
+
+
+
setMethod("modelFit", signature("gelModels"), valueClass="gelfit",
definition = function(object, gelType=NULL, rhoFct=NULL,
initTheta=c("gmm", "modelTheta0"), theta0=NULL,
Modified: pkg/gmm4/R/gelfit-methods.R
===================================================================
--- pkg/gmm4/R/gelfit-methods.R 2019-11-06 22:42:07 UTC (rev 151)
+++ pkg/gmm4/R/gelfit-methods.R 2019-11-15 19:56:25 UTC (rev 152)
@@ -41,10 +41,10 @@
{
tControl <- list(method="Brent", lower=rang[1], upper=rang[2])
fit2 <- suppressWarnings(update(fit, newModel=model, tControl=tControl,
- start.tet=coef(fit)[-which]))
+ theta0=coef(fit)[-which]))
} else {
fit2 <- suppressWarnings(update(fit, newModel=model,
- start.tet=coef(fit)[-which]))
+ theta0=coef(fit)[-which]))
}
test <- specTest(fit2, type=type, ...)@test[1] - test0
if (is.null(corr))
@@ -115,11 +115,15 @@
setMethod("getImpProb", "gelfit",
function(object, posProb=FALSE, normalize=TRUE) {
- rhoFct <- object at model@gelType
+ rhoFct <- object at model@gelType
if (is.null(rhoFct$fct))
+ {
+ if (rhoFct$name %in% c("ETEL","ETHD"))
+ rhoFct$name <- "ET"
rhoFct <- get(paste("rho", rhoFct$name, sep=""))
- else
+ } else {
rhoFct <- rhoFct$fct
+ }
gt <- evalMoment(object at model, object at theta)
k <- object at model@wSpec$k
pt <- -rhoFct(gt, object at lambda, 1, k[1]/k[2])/nrow(gt)
@@ -269,7 +273,7 @@
}
dimnames(ans) <- list(nlam,
c((1 - level)/2, 0.5 + level/2))
- return(new("confint", interval=matrix,
+ return(new("confint", interval=ans,
type=ntest, level=level))
}
if (type == "Wald")
@@ -297,6 +301,168 @@
new("confint", interval=ans, type=ntest, level=level)
})
+setMethod("confint", "numeric",
+ function (object, parm, level = 0.95, gelType="EL",
+ type = c("Wald", "invLR", "invLM", "invJ"),
+ fact = 3, vcov="iid", ...)
+ {
+ Call <- try(match.call(call=sys.call(sys.parent())), silent=TRUE)
+ if (class(Call)=="try-error")
+ Call <- NULL
+ type <- match.arg(type)
+ object <- as.data.frame(object)
+ if (!is.null(Call))
+ names(object) <- as.character(Call)[2]
+ g <- as.formula(paste(names(object),"~1",sep=""))
+ n <- nrow(object)
+ s <- sd(object[[1]], na.rm=TRUE)/sqrt(n)
+ m <- mean(object[[1]], na.rm=TRUE)
+ mod <- gelModel(g, ~1, gelType=gelType, vcov=vcov, data=object)
+ fit <- modelFit(mod, tControl=list(method="Brent",lower=m-s,upper=m+s),
+ ...)
+ ans <- confint(fit, 1, level=level, type=type, fact=fact)
+ rownames(ans at interval) <- names(object)
+ ans
+ })
+
+
+.lineInterval <- function(obj, slope, which=1:2, ...)
+ {
+ if (length(which)!=2)
+ stop("You must select 2 coefficients")
+ if (length(coef(obj))<2)
+ stop("You need at least two coefficients")
+ mu <- coef(obj)[which]
+ b <- mu[2]-slope*mu[1]
+ cst <- paste(names(mu)[2],"=",slope,"*",names(mu)[1],"+",b)
+ mod <- restModel(obj at model, cst)
+
+ if (length(coef(obj)==2))
+ cont <- list(method="Brent", lower=mu[1]-1, upper=mu[1]+1)
+ else
+ cont <- list()
+ fit <- modelFit(mod, tControl=cont)
+ ci <- confint(fit, names(mu)[1], ...)
+ low <- coef(mod, ci at interval[1,1])
+ up <- coef(mod, ci at interval[1,2])
+ ans <- rbind(low, up)
+ colnames(ans) <- names(mu)
+ rownames(ans) <- c("low","up")
+ ans
+ }
+
+setMethod("confint", "data.frame",
+ function (object, parm, level = 0.95, gelType="EL",
+ type = c("Wald", "invLR", "invLM", "invJ"),
+ fact = 3, vcov="iid", n=10, cores=4, ...)
+ {
+ type <- match.arg(type)
+ if (Sys.info()["sysname"] == "Windows")
+ cores <- 1L
+ n <- floor(n/2)
+ if (missing(parm))
+ parm <- 1
+ if (length(parm) == 1)
+ {
+ ans <- confint(object[[parm]], level=level, gelType=gelType,
+ type=type, fact=fact, vcov=vcov, ...)
+ rownames(ans at interval) <- names(object[parm])
+ return(ans)
+ }
+ if (length(parm) != 2)
+ stop("You can only select 2 variable from your data.frame")
+ object <- object[parm]
+ parNames <- paste("mu_", names(object), sep="")
+ g <- list(as.formula(paste(names(object)[1], "~", parNames[1], sep="")),
+ as.formula(paste(names(object)[2], "~", parNames[2], sep="")))
+ theta0 <- colMeans(object)
+ names(theta0) <- parNames
+ mod <- gelModel(g, gelType=gelType, vcov=vcov, data=object,
+ theta0=theta0)
+ fit <- modelFit(mod, ...)
+ mu <- coef(fit)
+ c1 <- confint(fit,1:2,level, type=type, fact=fact)@interval
+ delta <- seq(1,0,length.out=n+2)[-c(1,n+2)]
+ p1 <- c(mu[1],c1[2,2])
+ p2 <- c(c1[1,2], mu[2])
+ p3 <- c(mu[1],c1[2,1])
+ p4 <- c(c1[1,1],mu[2])
+ slU <- sapply(delta, function(d) {
+ p <- (p1*d + p2*(1-d))
+ (p[2]-mu[2])/(p[1]-mu[1])})
+ slD <- sapply(delta, function(d) {
+ p <- (p2*d + p3*(1-d))
+ (p[2]-mu[2])/(p[1]-mu[1])})
+ slope <- c(slU, slD)
+ resU <- mclapply(slU, function(s) .lineInterval(fit, s, ...),
+ mc.cores=8)
+ resD <- mclapply(slD, function(s) .lineInterval(fit, s, ...),
+ mc.cores=8)
+ Q1 <- cbind(p3, sapply(resU, function(l) l[1,]))
+ Q2 <- cbind(p4, sapply(resD, function(l) l[1,]))
+ Q3 <- cbind(p1, sapply(resU, function(l) l[2,]))
+ Q4 <- cbind(p2, sapply(resD, function(l) l[2,]))
+ ans <- rbind(t(Q1),t(Q2), t(Q3), t(Q4))
+ colnames(ans) <- names(mu)
+ rownames(ans) <- NULL
+ type <- paste("2D ", type, " confidence region", sep="")
+ new("mconfint", areaPoints=ans, type=type, level=level)
+ })
+
+setMethod("confint", "ANY",
+ function (object, parm, level = 0.95, ...)
+ {
+ stats::confint(object, parm, level, ...)
+ })
+
+
+### print, show and plot for mconfint
+
+setMethod("print", "mconfint",
+ function(x, digits=4, ...)
+ {
+ cat(x at type, "\n")
+ cat(paste(rep("*", nchar(x at type)), collapse="", sep=""), "\n", sep="")
+ cat("Level: ", x at level, "\n", sep="")
+ cat("Number of points: ", nrow(x at areaPoints), "\n", sep="")
+ cat("Variables:\n")
+ v <- colnames(x at areaPoints)
+ r <- formatC(apply(x at areaPoints, 2, range), digits=digits, ...)
+ cat("\tRange for ", v[1], ": [", r[1,1], ", ", r[2,1], "]\n", sep="")
+ cat("\tRange for ", v[2], ": [", r[1,2], ", ", r[2,2], "]\n", sep="")
+ })
+
+setMethod("show", "mconfint", function(object) print(object))
+
+setGeneric("plot")
+
+setMethod("plot", "mconfint", function(x, y, main=NULL, xlab=NULL, ylab=NULL,
+ pch=21, bg=1, ...)
+{
+ v <- colnames(x at areaPoints)
+ if (is.null(main))
+ main <- paste(100*x at level, "% confidence region for ", v[1], " and ", v[2],
+ sep="")
+ if (is.null(xlab))
+ xlab <- v[1]
+ if (is.null(ylab))
+ ylab <- v[2]
+ plot(x at areaPoints, xlab=xlab, ylab=ylab, main=main, pch=pch, bg=bg)
+ grid()
+ polygon(x at areaPoints[,1], x at areaPoints[,2], ...)
+ invisible()
+})
+
+setMethod("plot", "ANY", function(x, y, ...)
+ if(!missing(y))
+ graphics::plot(x, y, ...)
+ else
+ graphics::plot(x, ...)
+ )
+
+
+
+
## specTest
setMethod("specTest", signature("gelfit", "missing"),
Modified: pkg/gmm4/R/rGelModel-methods.R
===================================================================
--- pkg/gmm4/R/rGelModel-methods.R 2019-11-06 22:42:07 UTC (rev 151)
+++ pkg/gmm4/R/rGelModel-methods.R 2019-11-15 19:56:25 UTC (rev 152)
@@ -96,7 +96,6 @@
## coef
-
setMethod("coef", "rgelModels",
function(object, theta)
{
@@ -116,6 +115,62 @@
## modelFit
+setMethod("modelFit", signature("rlinearGel"), valueClass="gelfit",
+ definition = function(object, 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")
+ Call <- NULL
+ met <- getMethod("modelFit", "rgelModels")
+ obj <- met(object, gelType, rhoFct, initTheta, theta0, lambda0, vcov, ...)
+ obj at call <- Call
+ obj
+ })
+
+setMethod("modelFit", signature("rnonlinearGel"), valueClass="gelfit",
+ definition = function(object, 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")
+ Call <- NULL
+ met <- getMethod("modelFit", "rgelModels")
+ obj <- met(object, gelType, rhoFct, initTheta, theta0, lambda0, vcov, ...)
+ obj at call <- Call
+ obj
+ })
+
+setMethod("modelFit", signature("rformulaGel"), valueClass="gelfit",
+ definition = function(object, 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")
+ Call <- NULL
+ met <- getMethod("modelFit", "rgelModels")
+ obj <- met(object, gelType, rhoFct, initTheta, theta0, lambda0, vcov, ...)
+ obj at call <- Call
+ obj
+ })
+
+setMethod("modelFit", signature("rfunctionGel"), valueClass="gelfit",
+ definition = function(object, 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")
+ Call <- NULL
+ met <- getMethod("modelFit", "rgelModels")
+ obj <- met(object, gelType, rhoFct, initTheta, theta0, lambda0, vcov, ...)
+ obj at call <- Call
+ obj
+ })
+
setMethod("modelFit", signature("rgelModels"), valueClass="gelfit",
definition = function(object, gelType=NULL, rhoFct=NULL,
Modified: pkg/gmm4/man/confint-methods.Rd
===================================================================
--- pkg/gmm4/man/confint-methods.Rd 2019-11-06 22:42:07 UTC (rev 151)
+++ pkg/gmm4/man/confint-methods.Rd 2019-11-15 19:56:25 UTC (rev 152)
@@ -4,11 +4,66 @@
\alias{confint,ANY-method}
\alias{confint,gelfit-method}
\alias{confint,gmmfit-method}
+\alias{confint,numeric-method}
+\alias{confint,data.frame-method}
+\alias{confint,ANY-method}
\title{ ~~ Methods for Function \code{confint} in Package \pkg{stats} ~~}
\description{
Method to contruct confidence intervals for objects of class
\code{"gmmfit"} and \code{"gelfit"}.
}
+
+\usage{
+\S4method{confint}{gmmfit}(object, parm, level = 0.95, vcov=NULL, \dots)
+
+\S4method{confint}{gelfit}(object, parm, level = 0.95, lambda = FALSE,
+ type = c("Wald", "invLR", "invLM", "invJ"),
+ fact = 3, corr = NULL, vcov=NULL, \dots)
+
+\S4method{confint}{data.frame}(object, parm, level = 0.95, gelType="EL",
+ type = c("Wald", "invLR", "invLM", "invJ"),
+ fact = 3, vcov="iid", n=10,
+ cores=4, \dots)
+
+\S4method{confint}{ANY}(object, parm, level = 0.95, \dots)
+}
+
+\arguments{
+
+ \item{object}{Object of class \code{"gmmfit"},
+ \code{"gelfit"}, \code{"numeric"} or \code{"data.frame"}.}
+
+ \item{parm}{Vector of integers or characters for selecting the
+ elements for which the intervals should be computed.}
+
+ \item{level}{The confidence level.}
+
+ \item{lambda}{Should be compute intervals for the Lagrange
+ multipliers?}
+
+ \item{type}{The type of confidence intervals. The default is the Wald
+ interval, and the others are computed by inverting the LR, LM or J
+ specification test.}
+
+ \item{fact}{For the inversion of the specification tests,
+ \code{\link{uniroot}} searches within \code{fact} standard error of
+ the coefficient estimates}
+
+ \item{corr}{Correction to apply to the specification tests}
+
+ \item{vcov}{For Wald intervals, an optional covariance matrix can be
+ provided}
+
+ \item{cores}{The number of cores for \code{\link{mclapply}}. It is set
+ to 1 for Windows OS.}
+
+ \item{gelType}{Type of GEL confidence interval.}
+
+ \item{n}{Number of lines to construct the confidence region.}
+
+ \item{\dots}{Other arguments to pass to \code{\link{modelFit}}}
+}
+
\section{Methods}{
\describe{
@@ -23,6 +78,19 @@
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/gmm -r 152
More information about the Gmm-commits
mailing list