Add extract.iv_robust as clone of extract.lm_robust; add some tests; c…

…loses #202

NAMESPACE

@@ -32,6 +32,8 @@ S3method(vcov,lm_robust)
 export(commarobust)
 export(declaration_to_condition_pr_mat)
 export(difference_in_means)
+export(extract.iv_robust)
+export(extract.lm_robust)
 export(gen_pr_matrix_cluster)
 export(horvitz_thompson)
 export(iv_robust)

R/helper_extract.R

@@ -1,12 +1,13 @@
 # This code modified from
 # https://github.com/leifeld/texreg/blob/master/R/extract.R (no LICENSE)
 #' Extract model data for \pkg{texreg} package
+#' @rdname extract.lm_robust
 #'
-#' @description Prepares an \code{"lm_robust"} object for the \pkg{texreg}
+#' @description Prepares a \code{"lm_robust"} or \code{"iv_robust"} object for the \pkg{texreg}
 #' package. This is largely a clone of the \code{\link[texreg]{extract.lm}}
 #' method.
 #'
-#' @param model an object of class \code{\link{lm_robust}}
+#' @param model an object of class \code{\link{lm_robust}} or \code{"iv_robust"}
 #' @param include.ci logical. Defaults to TRUE
 #' @param include.rsquared logical. Defaults to TRUE
 #' @param include.adjrs logical. Defaults to TRUE
@@ -15,7 +16,7 @@
 #' @param include.rmse logical. Defaults to TRUE
 #' @param ... unused
 #'
-extract.lm_robust <- function(model,
+extract.robust_default <- function(model,
                               include.ci = TRUE,
                               include.rsquared = TRUE,
                               include.adjrs = TRUE,
@@ -84,3 +85,14 @@ extract.lm_robust <- function(model,
   )
   return(tr)
 }
+
+#' @rdname extract.lm_robust
+#'
+#' @export
+extract.lm_robust <- extract.robust_default
+
+#' @rdname extract.lm_robust
+#'
+#' @export
+extract.iv_robust <- extract.robust_default
+

R/zzz.R

@@ -9,6 +9,10 @@
       signature = className("lm_robust", pkgname),
       definition = extract.lm_robust
     )
+    setMethod("extract",
+              signature = className("iv_robust", pkgname),
+              definition = extract.iv_robust
+    )
   }
   invisible()
 }

tests/testthat/run-stata-iv-models.do

@@ -9,32 +9,26 @@ file open outf using stata-iv-ests.txt, write r
 
 ivregress 2sls mpg (hp am = wt gear), small
 mat V=e(V)
-scalar F=e(F)
-file write outf _n "classical" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (F)
+file write outf _n "classical" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (e(F)) _tab (e(r2)) _tab (e(r2_a)) _tab (e(rmse))
 
 ivregress 2sls mpg (hp am = wt gear), small rob
 mat V=e(V)
-scalar F=e(F)
-file write outf _n "rob" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (F)
+file write outf _n "rob" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (e(F)) _tab (e(r2)) _tab (e(r2_a)) _tab (e(rmse))
 
 ivregress 2sls mpg (hp am = wt gear), small vce(cluster cyl)
 mat V=e(V)
-scalar F=e(F)
-file write outf _n "cl" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (F)
+file write outf _n "cl" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (e(F)) _tab (e(r2)) _tab (e(r2_a)) _tab (e(rmse))
 
 ivregress 2sls mpg (hp am = wt gear) [aweight = w], small
 mat V=e(V)
-scalar F=e(F)
-file write outf _n "classical_w" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (F)
+file write outf _n "classical_w" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (e(F)) _tab (e(r2)) _tab (e(r2_a)) _tab (e(rmse))
 
 ivregress 2sls mpg (hp am = wt gear) [aweight = w], small rob
 mat V=e(V)
-scalar F=e(F)
-file write outf _n "rob_w" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (F)
+file write outf _n "rob_w" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (e(F)) _tab (e(r2)) _tab (e(r2_a)) _tab (e(rmse))
 
 ivregress 2sls mpg (hp am = wt gear) [aweight = w], small vce(cluster cyl)
 mat V=e(V)
-scalar F=e(F)
-file write outf _n "cl_w" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (F)
+file write outf _n "cl_w" _tab (V[1,1])  _tab (V[2,2])  _tab (V[3,3]) _tab (e(F)) _tab (e(r2)) _tab (e(r2_a)) _tab (e(rmse))
 
 file close outf

tests/testthat/stata-iv-ests.txt

@@ -1,7 +1,7 @@
 
-classical	.0001602	2.2793852	5.2396632	44.808289
-rob	.00015781	2.2003453	4.7873499	41.366981
-cl	.00004096	.7950331	.31612954	108.37591
-classical_w	.00015649	2.3424049	5.2400158	46.736059
-rob_w	.00014704	2.0834808	4.5677404	45.410829
-cl_w	.00006651	.64892154	.54670656	81.595912
+classical	.0001602	2.2793852	5.2396632	44.808289	.75605596	.73923223	3.0776867
+rob	.00015781	2.2003453	4.7873499	41.366981	.75605596	.73923223	3.0776867
+cl	.00004096	.7950331	.31612954	108.37591	.75605596	.73923223	3.0776867
+classical_w	.00015649	2.3424049	5.2400158	46.736059	.76415188	.74788649	3.0971208
+rob_w	.00014704	2.0834808	4.5677404	45.410829	.76415188	.74788649	3.0971208
+cl_w	.00006651	.64892154	.54670656	81.595912	.76415188	.74788649	3.0971208

tests/testthat/test-stata-output.R

@@ -60,26 +60,26 @@ test_that("iv_robust matches stata", {
 
   stata_ests <- read.table(
     "stata-iv-ests.txt",
-    col.names = c("model", "v1", "v2", "v3", "fstat"),
+    col.names = c("model", "v1", "v2", "v3", "fstat", "r2", "r2_a", "rmse"),
     stringsAsFactors = FALSE
   )
 
   mtcars$w <- mtcars$drat / 5
 
-  estimatr_mat <- matrix(NA, 6, 4)
+  estimatr_mat <- matrix(NA, 6, 7)
   iv_c <- iv_robust(mpg ~ hp + am | wt + gear, data = mtcars, se_type = "classical")
-  estimatr_mat[1, ] <- c(iv_c$std.error ^ 2, iv_c$fstatistic[1])
+  estimatr_mat[1, ] <- c(iv_c$std.error ^ 2, iv_c$fstatistic[1], iv_c$r.squared, iv_c$adj.r.squared, sqrt(iv_c$res_var))
   iv_hc1 <- iv_robust(mpg ~ hp + am | wt + gear, data = mtcars, se_type = "HC1")
-  estimatr_mat[2, ] <- c(iv_hc1$std.error ^ 2, iv_hc1$fstatistic[1])
+  estimatr_mat[2, ] <- c(iv_hc1$std.error ^ 2, iv_hc1$fstatistic[1], iv_hc1$r.squared, iv_hc1$adj.r.squared, sqrt(iv_hc1$res_var))
   iv_stata <- iv_robust(mpg ~ hp + am | wt + gear, clusters = cyl, data = mtcars, se_type = "stata")
-  estimatr_mat[3, ] <- c(iv_stata$std.error ^ 2, iv_stata$fstatistic[1])
+  estimatr_mat[3, ] <- c(iv_stata$std.error ^ 2, iv_stata$fstatistic[1], iv_stata$r.squared, iv_stata$adj.r.squared, sqrt(iv_stata$res_var))
 
   iv_c_w <- iv_robust(mpg ~ hp + am | wt + gear, data = mtcars, weights = w, se_type = "classical")
-  estimatr_mat[4, ] <- c(iv_c_w$std.error ^ 2, iv_c_w$fstatistic[1])
+  estimatr_mat[4, ] <- c(iv_c_w$std.error ^ 2, iv_c_w$fstatistic[1], iv_c_w$r.squared, iv_c_w$adj.r.squared, sqrt(iv_c_w$res_var))
   iv_hc1_w <- iv_robust(mpg ~ hp + am | wt + gear, data = mtcars, weights = w, se_type = "HC1")
-  estimatr_mat[5, ] <- c(iv_hc1_w$std.error ^ 2, iv_hc1_w$fstatistic[1])
+  estimatr_mat[5, ] <- c(iv_hc1_w$std.error ^ 2, iv_hc1_w$fstatistic[1], iv_hc1_w$r.squared, iv_hc1_w$adj.r.squared, sqrt(iv_hc1_w$res_var))
   iv_stata_w <- iv_robust(mpg ~ hp + am | wt + gear, clusters = cyl, weights = w, data = mtcars, se_type = "stata")
-  estimatr_mat[6, ] <- c(iv_stata_w$std.error ^ 2, iv_stata_w$fstatistic[1])
+  estimatr_mat[6, ] <- c(iv_stata_w$std.error ^ 2, iv_stata_w$fstatistic[1], iv_stata_w$r.squared, iv_stata_w$adj.r.squared, sqrt(iv_stata_w$res_var))
 
   expect_true(
     max(abs(estimatr_mat[, 1] - as.numeric(stata_ests[, 4]))) < 2e-05
@@ -89,4 +89,24 @@ test_that("iv_robust matches stata", {
     max(abs(estimatr_mat[, 4] - as.numeric(stata_ests[, 5]))) < 3e-05
   )
 
+  # Note, RMSE is different for stata with weights than ivreg or iv_robust
+  expect_true(
+    max(abs(estimatr_mat[, 5:6] - stata_ests[, 6:7])) < 4e-08
+  )
+
+  ivrego_w <- AER::ivreg(mpg ~ hp + am | wt + gear, data = mtcars, weights = w)
+
+  expect_equal(
+    ivrego_w$sigma,
+    sqrt(iv_c_w$res_var)
+  )
+  expect_equal(
+    ivrego_w$sigma,
+    sqrt(iv_hc1_w$res_var)
+  )
+  expect_equal(
+    ivrego_w$sigma,
+    sqrt(iv_stata_w$res_var)
+  )
+
 })