[Gmm-commits] r231 - in pkg/momentfit: . R vignettes

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Fri Mar 22 22:10:17 CET 2024


Author: chaussep
Date: 2024-03-22 22:10:17 +0100 (Fri, 22 Mar 2024)
New Revision: 231

Modified:
   pkg/momentfit/DESCRIPTION
   pkg/momentfit/R/gelfit-methods.R
   pkg/momentfit/R/weak.R
   pkg/momentfit/vignettes/weak.Rmd
   pkg/momentfit/vignettes/weak.pdf
Log:
worked on weak test and fixed abig with smoothed gel

Modified: pkg/momentfit/DESCRIPTION
===================================================================
--- pkg/momentfit/DESCRIPTION	2024-03-06 20:53:34 UTC (rev 230)
+++ pkg/momentfit/DESCRIPTION	2024-03-22 21:10:17 UTC (rev 231)
@@ -1,6 +1,6 @@
 Package: momentfit
 Version: 0.6
-Date: 2024-01-24
+Date: 2024-03-21
 Title: Methods of Moments
 Author: Pierre Chausse <pchausse at uwaterloo.ca>
 Maintainer: Pierre Chausse <pchausse at uwaterloo.ca>

Modified: pkg/momentfit/R/gelfit-methods.R
===================================================================
--- pkg/momentfit/R/gelfit-methods.R	2024-03-06 20:53:34 UTC (rev 230)
+++ pkg/momentfit/R/gelfit-methods.R	2024-03-22 21:10:17 UTC (rev 231)
@@ -553,7 +553,8 @@
               }
               if (type %in% c("All","J"))
               {
-                  J <- sum(lm.fit(gt, rep(1,n))$fitted.values)
+                  J <- sum(lm.fit(gt, rep(1,NROW(gt)))$fitted.values)/
+                      object at model@sSpec at k[1]^2
                   test <- c(test, J)
                   names(test)[length(test)] <- " J: "                  
               }

Modified: pkg/momentfit/R/weak.R
===================================================================
--- pkg/momentfit/R/weak.R	2024-03-06 20:53:34 UTC (rev 230)
+++ pkg/momentfit/R/weak.R	2024-03-22 21:10:17 UTC (rev 231)
@@ -667,11 +667,14 @@
     eps <- .Machine$double.eps
     m <- nrow(A)
     n <- ncol(A)
-    if (m != n) {
+    if (m != n)
+    {
         stop("Argument 'A' must be a square matrix.")
+    } else if (n == 1 && A <= 0) {
+        return(as.matrix(eps))
+    } else if (n == 1) {
+        return(A)
     }
-    else if (n == 1 && A <= 0) 
-        return(as.matrix(eps))
     B <- (A + t(A))/2
     svdB <- svd(B)
     H <- svdB$v %*% diag(svdB$d) %*% t(svdB$v)
@@ -859,7 +862,7 @@
         RR <- qr.R(res)
         diagRR <- sign(diag(RR))
         if (any(diagRR < 0))
-            Q <- Q%*%diag(diagRR)
+            Q <- Q%*%diag(diagRR, length(diagRR))
         list(Q=Q, R=RR)
     }
 
@@ -946,7 +949,6 @@
         Qp <- Q
         Q <- gamma*Qp + 1
         Cval = (gamma*Qp*Cval + F)/Q
-
     }
     if (itr >= opts$mxitr)
         out$msg = "exceed max iteration"

Modified: pkg/momentfit/vignettes/weak.Rmd
===================================================================
--- pkg/momentfit/vignettes/weak.Rmd	2024-03-06 20:53:34 UTC (rev 230)
+++ pkg/momentfit/vignettes/weak.Rmd	2024-03-22 21:10:17 UTC (rev 231)
@@ -538,7 +538,7 @@
 to `TRUE`.
 
 
-## Testing for weak instrument: @montiel-olea-pflueger13
+## Testing for weak instruments: @montiel-olea-pflueger13
 
 In most applied economic studies, it is unrealistic to assume that the
 errors are conditionally homoskedastic. When the errors are
@@ -676,7 +676,7 @@
 ```{r}
 w2
 ```
-	
+
 This is what the function `getMOPW` computes. For now, it is not
 exported, so we need to run it using `momentfit:::`.  The function
 returns the different $\hat{W}$'s separately for convenience. Here we
@@ -717,7 +717,7 @@
 argument `simplified` to `FALSE`.
 
 ```{r}
-MOPtest(mod2, simplified=FALSE, estMethod="LIML")
+MOPtest(mod5, simplified=FALSE, estMethod="LIML")
 ```
 
 The test is the same as above, but the critical value is smaller,
@@ -726,12 +726,9 @@
 small. We can compare the test with the one for TSLS:
 
 ```{r}
-MOPtest(mod2, simplified=FALSE, estMethod="TSLS")
+MOPtest(mod5,  simplified=FALSE, estMethod="TSLS")
 ```
 
-It seems that we are rejecting the weak instrument hypothesis for TSLS
-and not for LIML. This is surprising. We'll need to investigate.
-
 As mentioned by the authors, the efficient F is the robust F when the
 model is just identified. The model `mod4` created above is
 just-identified, but we need to change the argument `vcov` to `"MDS"`:
@@ -858,8 +855,44 @@
 m <- momentModel(g, h, data=dat, vcov="MDS")
 ``` 
 
+```{r, echo=FALSE}
+## n: sample size
+## N: number of endogenous
+## K: number of Instruments
+## Nx: number of exogenous regressors
 
+genDat <- function(n=200, N=2, K=5, Nx=2)
+{
+    beta <- rep(1,N) 
+    DL <- 0.08
+    X <- qr.R(qr(matrix(rnorm(K^2),K,K)))
+    L0 <- t(X[,1:N])
+    Pi <- sqrt(K)*DL^0.5*L0
+    Pi <- t(Pi)
+    By <- rnorm(Nx+1)
+    BY <- matrix(rnorm((Nx+1)*N), Nx+1, N)
+    u <- rnorm(n)
+    v <- matrix(rnorm(n*N), n, N)
+    X <- cbind(matrix(rnorm(n*Nx), n, Nx), 1)
+    Z <- matrix(rnorm(n*K), n, K)
+    Y <- Z%*%Pi+X%*%BY+v
+    y <- c(Y%*%beta+X%*%By+u)
+    X <- X[,-ncol(X),drop=FALSE]    
+    colnames(Y) <- paste("Y", 1:ncol(Y), sep="")
+    colnames(X) <- paste("X", 1:ncol(X), sep="")
+    colnames(Z) <- paste("Z", 1:ncol(Z), sep="")
+    dat <- as.data.frame(cbind(Y,X,Z))
+    dat$y <- y
+    dat
+}
+dat <- genDat()
+m <- momentModel(y~Y1+X1+X2, ~Z1+Z2+Z3+Z4+Z5+X1+X2, data=dat,
+                 vcov="MDS")
+LewMertest(m, simplified=FALSE, npoints=100)
+MOPtest(m, estMethod="TSLS", simplified=FALSE)
+```
 
+
 # References
 
 

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



More information about the Gmm-commits mailing list