[Gmm-commits] r76 - in pkg: gmm/R gmm/vignettes gmmExtra gmmExtra/R

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Fri Feb 6 23:50:15 CET 2015


Author: chaussep
Date: 2015-02-06 23:50:15 +0100 (Fri, 06 Feb 2015)
New Revision: 76

Modified:
   pkg/gmm/R/getModel.R
   pkg/gmm/R/gmm.R
   pkg/gmm/R/momentEstim.R
   pkg/gmm/vignettes/gmm_with_R.rnw
   pkg/gmmExtra/DESCRIPTION
   pkg/gmmExtra/NAMESPACE
   pkg/gmmExtra/R/KConfig.R
   pkg/gmmExtra/R/bootGmm.R
Log:
Fixed bugs from last change and worked in gmmExtra

Modified: pkg/gmm/R/getModel.R
===================================================================
--- pkg/gmm/R/getModel.R	2015-02-03 21:28:39 UTC (rev 75)
+++ pkg/gmm/R/getModel.R	2015-02-06 22:50:15 UTC (rev 76)
@@ -78,7 +78,7 @@
 	obj$x$k <- obj$x$k-nrow(object$eqConst)
 	if (obj$x$k<=0)
 		stop("Nothing to estimate")
-	} else {
+        } else {
 	attr(obj$x,"eqConst") <- list(unConstg = obj$g, unConstgradv = obj$gradv, eqConst = object$eqConst)   	   	
 	obj$g <- function(tet, dat)
 		{
@@ -106,8 +106,7 @@
   obj$type <- paste(obj$type,"(with equality constraints)",sep=" ")	
   mess <- paste(rownames(object$eqConst), " = " , object$eqConst[,2], "\n",collapse="")
   mess <- paste("#### Equality constraints ####\n",mess,"##############################\n\n",sep="")
-  obj$specMod <- mess  
-                                      
+  obj$specMod <- mess                                   
   return(obj)
   }
 

Modified: pkg/gmm/R/gmm.R
===================================================================
--- pkg/gmm/R/gmm.R	2015-02-03 21:28:39 UTC (rev 75)
+++ pkg/gmm/R/gmm.R	2015-02-06 22:50:15 UTC (rev 76)
@@ -167,27 +167,33 @@
   hm <- as.matrix(x[,(ny+k+1):(ny+k+nh)])
   includeExo <- which(colnames(xm)%in%colnames(hm))
   if (!is.null(type))
-  	{
+  	{            
   	if(type=="2sls")
-	  	{
+	  	{                   
                 if (length(includeExo) > 0)
-                    {
+                    {                        
                     endo <- xm[, -includeExo, drop = FALSE]
                     endoName <- colnames(endo)
-                    restsls <- lm(endo~hm-1)
-                    fsls <- xm
-                    fsls[, -includeExo] <- restsls$fitted
+                    if (ncol(endo) != 0)
+                        {
+                            restsls <- lm(endo~hm-1)
+                            fsls <- xm
+                            fsls[, -includeExo] <- restsls$fitted
+                        } else {
+                            fsls <- xm
+                            restsls <- NULL
+                        }
                 } else {
                     restsls <- lm(xm~hm-1)
                     fsls <- restsls$fitted
                     endoName <- colnames(xm)
-                }
+                }                
   	     	par <- lm.fit(as.matrix(fsls), ym)$coefficients
 		if (ny == 1)
-		{
+		{                                    
   	     	e2sls <- ym-xm%*%par
  	     	v2sls <- lm.fit(as.matrix(hm), e2sls)$fitted
-  	     	value <- sum(v2sls^2)/sum(e2sls^2)
+  	     	value <- sum(v2sls^2)/sum(e2sls^2)                
   	     }
   	     else
   	     {
@@ -200,7 +206,7 @@
 	  	}
   	}
   else
-  	{
+  	{            
   if (ny>1)
   	{
      if (inv) 
@@ -239,9 +245,12 @@
   if (!is.null(type))
      {    
      if (type == "2sls")
-     res$firstStageReg <- restsls	
-     res$fsRes <- summary(restsls)
-     attr(res$fsRes, "Endo") <- endoName
+     res$firstStageReg <- restsls
+     if (!is.null(restsls))
+         {
+             res$fsRes <- summary(restsls)
+             attr(res$fsRes, "Endo") <- endoName
+         }
      }
   return(res)
   }

Modified: pkg/gmm/R/momentEstim.R
===================================================================
--- pkg/gmm/R/momentEstim.R	2015-02-03 21:28:39 UTC (rev 75)
+++ pkg/gmm/R/momentEstim.R	2015-02-06 22:50:15 UTC (rev 76)
@@ -197,9 +197,9 @@
     z = list(coefficients = res2$par, objective = res2$value, dat = dat, k = k, k2 = k2, n = n, q = q, df = df, df.residual = (n-k))
     }
   else
-    {
+    {  
     if (P$vcov == "iid")
-    	{
+   	{            
       res1 <- .tetlin(dat, w, P$gradv, P$g, type="2sls")
       initTheta <- res1$par
       gmat <- g(res1$par, dat)
@@ -209,7 +209,7 @@
       res2$fsRes <- res1$fsRes
       }
     if (P$vcov == "HAC")
-      {
+      {          
       res1 <- .tetlin(dat, w, P$gradv, P$g, type="2sls")
       initTheta <- res1$par
       gmat <- g(res1$par, dat)
@@ -235,7 +235,6 @@
   z$g <- P$g
   z$WSpec <- P$WSpec
   z$w0 <- w
-
   names(z$coefficients) <- P$namesCoef
   colnames(z$gt) <- P$namesgt
  

Modified: pkg/gmm/vignettes/gmm_with_R.rnw
===================================================================
--- pkg/gmm/vignettes/gmm_with_R.rnw	2015-02-03 21:28:39 UTC (rev 75)
+++ pkg/gmm/vignettes/gmm_with_R.rnw	2015-02-06 22:50:15 UTC (rev 76)
@@ -318,6 +318,10 @@
 For this sub-sample, the hypothesis that the return follows a stable distribution is  rejected. The normality assumption can be analyzed  by testing $H_0:\alpha=2,~\beta=0$ using \code{linearHypothesis} from the \pkg{car} package:
 <<>>=
 library(car)
+# the following line is to go around a bug with car version 2.0.24
+# It is fixed in version 2.0.25 but not on CRAN yet
+# It does not affect the result of the test
+res3$df.residual=1500-4
 linearHypothesis(res3,cbind(diag(2),c(0,0),c(0,0)),c(2,0))
 @
 It is clearly rejected. The result is even stronger if the whole sample is used. 

Modified: pkg/gmmExtra/DESCRIPTION
===================================================================
--- pkg/gmmExtra/DESCRIPTION	2015-02-03 21:28:39 UTC (rev 75)
+++ pkg/gmmExtra/DESCRIPTION	2015-02-06 22:50:15 UTC (rev 76)
@@ -5,7 +5,7 @@
 Author: Pierre Chausse <pchausse at uwaterloo.ca>
 Maintainer: Pierre Chausse <pchausse at uwaterloo.ca>
 Description: Tools for GMM such as additional tests or robust confidence regions. They only apply to gmm class object for the gmm package.
-Depends: R (>= 2.10.0), gmm (>= 1.4), parallel, car
+Depends: R (>= 2.10.0), gmm (>= 1.5), parallel, car
 Suggests: mvtnorm, stabledist, MASS, timeDate, timeSeries
 Imports: stats
 License: GPL (>= 2)

Modified: pkg/gmmExtra/NAMESPACE
===================================================================
--- pkg/gmmExtra/NAMESPACE	2015-02-03 21:28:39 UTC (rev 75)
+++ pkg/gmmExtra/NAMESPACE	2015-02-06 22:50:15 UTC (rev 76)
@@ -1,5 +1,7 @@
 import(stats)
-importFrom(car, linearHypothesis)
+importFrom(car, linearHypothesis, makeHypothesis, printHypothesis)
+importFrom(parallel, mclapply)
+importFrom(gmm, gmm, specTest, gmmWithConst, KTest)
 
 export(KConfid, bootGmm, plot.bootGmm, summary.bootGmm, bootJ, linearHypothesis.bootGmm, print.bootGmm)
 

Modified: pkg/gmmExtra/R/KConfig.R
===================================================================
--- pkg/gmmExtra/R/KConfig.R	2015-02-03 21:28:39 UTC (rev 75)
+++ pkg/gmmExtra/R/KConfig.R	2015-02-06 22:50:15 UTC (rev 76)
@@ -1,6 +1,9 @@
-KConfid <- function(obj, which, type = c("K", "KJ"), alpha = 0.05, alphaJ = 0.01, n = 4)
+KConfid <- function(obj, which, type = c("K", "KJ"), alpha = 0.05, alphaJ = 0.01, n = 4,
+                    mc.cores=2)
 	{
-	theApply <- mclapply
+        theApply <- mclapply
+        if (Sys.info()[["sysname"]] == "Windows")
+            mc.cores <- 1
 
 	type <- match.arg(type)
 	if ( (obj$df == 0) & (type == "KJ"))
@@ -60,19 +63,23 @@
 		X0 <- obj$coef[which]
 		getBoth <- function(d)
 			{
-			res <- try(uniroot(getUniInt,c(X0,X0+d*step), obj = obj, which = which, type = type, 
-					alpha = alpha, alphaJ = alphaJ, alphaK = alphaK, tol=1e-4), silent=TRUE)
+			res <- try(uniroot(getUniInt,c(X0,X0+d*step),
+                                           obj = obj, which = which, type = type, 
+                                           alpha = alpha, alphaJ = alphaJ, alphaK = alphaK,
+                                           tol=1e-4), silent=TRUE)
 			if (class(res) == "try-error")
 				return(NA)
 			else
 				return(res$root)
 			}
-		res <- theApply(c(-1,1), getBoth)
+		res <- theApply(c(-1,1), getBoth, mc.cores=mc.cores)
 		c(res[[1]],res[[2]])
 	} else {
 		step <- 5*sqrt(diag(vcov(obj)))[which]
-		sol <- .getCircle(obj$coef[which][1],obj$coef[which][2],gForMulti,n , step, obj=obj, 
-					which=which, type=type, alpha=alpha, alphaJ=alphaJ, alphaK=alphaK)	
+		sol <- .getCircle(obj$coef[which][1],obj$coef[which][2],
+                                  gForMulti,n , step, obj=obj, 
+                                  which=which, type=type, alpha=alpha,
+                                  alphaJ=alphaJ, alphaK=alphaK)	
 		colnames(sol) <- names(obj$coef)[which]
 		sol
 		}
@@ -103,8 +110,12 @@
 	# get the four points of the cross
 	selectf <- list(f,f,f2,f2)
 	selectd <- c(1,-1,1,-1)
-	xy12 <- theApply(1:2, function(i) try(uniroot(selectf[[i]],c(x0,x0+selectd[i]*b[1]),y0=y0,tol=tol)$root,silent=TRUE))
-	xy34 <- theApply(3:4, function(i) try(uniroot(selectf[[i]],c(y0,y0+selectd[i]*b[2]),y0=x0,tol=tol)$root,silent=TRUE))
+	xy12 <- theApply(1:2, function(i)
+            try(uniroot(selectf[[i]],c(x0,x0+selectd[i]*b[1]),y0=y0,tol=tol)$root,
+                silent=TRUE))
+	xy34 <- theApply(3:4, function(i)
+            try(uniroot(selectf[[i]],c(y0,y0+selectd[i]*b[2]),y0=x0,tol=tol)$root,
+                silent=TRUE))
 
 	x1 <- xy12[[1]]
 	x2 <- xy12[[2]]
@@ -115,7 +126,8 @@
 		{
 		yi <- y1*lambda + y0*(1-lambda)
 		xi <- x0*lambda + x*(1-lambda)
-		xt <- try(uniroot(f,c(x0,x0+dir*b[1]),y0=y0, x0=x0, xi=xi, yi=yi,tol=tol)$root,silent=TRUE)
+		xt <- try(uniroot(f,c(x0,x0+dir*b[1]),y0=y0, x0=x0,
+                                  xi=xi, yi=yi,tol=tol)$root,silent=TRUE)
 		xt <- f3(xt)
 		yt <- y0 + (yi-y0)/(xi-x0)*(xt-x0)
 		return(c(xt,yt))

Modified: pkg/gmmExtra/R/bootGmm.R
===================================================================
--- pkg/gmmExtra/R/bootGmm.R	2015-02-03 21:28:39 UTC (rev 75)
+++ pkg/gmmExtra/R/bootGmm.R	2015-02-06 22:50:15 UTC (rev 76)
@@ -137,6 +137,8 @@
 
 bootGmm <- function(obj, N, seed = NULL, niter=8, trace=TRUE)
 	{
+        if (Sys.info()[["sysname"]] == "Windows")
+            niter <- 8
 	if(!is.null(seed))
 		set.seed(seed)
 	argCall <- obj$allArg
@@ -304,7 +306,7 @@
 
 	# The following is borrowed from the car::linearHypothesis.default
 	if (is.character(hypothesis.matrix)) {
-		L <- car:::makeHypothesis(names(attr(obj,"gmm")$coefficients), hypothesis.matrix, rhs = NULL)
+		L <- makeHypothesis(names(attr(obj,"gmm")$coefficients), hypothesis.matrix, rhs = NULL)
 		if (is.null(dim(L))) L <- t(L)
 		rhs <- L[, NCOL(L)]
 		L <- L[, -NCOL(L), drop = FALSE]
@@ -315,7 +317,7 @@
 				else hypothesis.matrix
 		if (is.null(rhs)) rhs <- rep(0,nrow(L))
 	}
-	hyp <- car:::printHypothesis(L, rhs, names(attr(obj,"gmm")$coefficients))
+	hyp <- printHypothesis(L, rhs, names(attr(obj,"gmm")$coefficients))
 	
 	##############
 	rhsB <- L%*%attr(obj,"gmm")$coefficients



More information about the Gmm-commits mailing list