]> code.communitydata.science - ml_measurement_error_public.git/blobdiff - simulations/simulation_base.R
make first simulation with precise accuracies and explained variances
[ml_measurement_error_public.git] / simulations / simulation_base.R
index 26916abb7dee3f0a509648cb2ef3ba6ecc64cffd..345d14e34a092e5e3239e7c1bc92153a70d3f011 100644 (file)
@@ -82,26 +82,26 @@ run_simulation <-  function(df, result){
                           ))
 
     ## What if we can't observe k -- most realistic scenario. We can't include all the ML features in a model.
-    amelia.out.nok <- amelia(df, m=200, p2s=0, idvars=c("x","k","w_pred"), noms=noms)
-    mod.amelia.nok <- zelig(y~x.obs+g, model='ls', data=amelia.out.nok$imputations, cite=FALSE)
-    (coefse <- combine_coef_se(mod.amelia.nok, messages=FALSE))
-
-    est.x.mi <- coefse['x.obs','Estimate']
-    est.x.se <- coefse['x.obs','Std.Error']
-    result <- append(result,
-                     list(Bxy.est.amelia.nok = est.x.mi,
-                          Bxy.ci.upper.amelia.nok = est.x.mi + 1.96 * est.x.se,
-                          Bxy.ci.lower.amelia.nok = est.x.mi - 1.96 * est.x.se
-                          ))
-
-    est.g.mi <- coefse['g','Estimate']
-    est.g.se <- coefse['g','Std.Error']
-
-    result <- append(result,
-                     list(Bgy.est.amelia.nok = est.g.mi,
-                          Bgy.ci.upper.amelia.nok = est.g.mi + 1.96 * est.g.se,
-                          Bgy.ci.lower.amelia.nok = est.g.mi - 1.96 * est.g.se
-                          ))
+    ## amelia.out.nok <- amelia(df, m=200, p2s=0, idvars=c("x","w_pred"), noms=noms)
+    ## mod.amelia.nok <- zelig(y~x.obs+g, model='ls', data=amelia.out.nok$imputations, cite=FALSE)
+    ## (coefse <- combine_coef_se(mod.amelia.nok, messages=FALSE))
+
+    ## est.x.mi <- coefse['x.obs','Estimate']
+    ## est.x.se <- coefse['x.obs','Std.Error']
+    ## result <- append(result,
+    ##                  list(Bxy.est.amelia.nok = est.x.mi,
+    ##                       Bxy.ci.upper.amelia.nok = est.x.mi + 1.96 * est.x.se,
+    ##                       Bxy.ci.lower.amelia.nok = est.x.mi - 1.96 * est.x.se
+    ##                       ))
+
+    ## est.g.mi <- coefse['g','Estimate']
+    ## est.g.se <- coefse['g','Std.Error']
+
+    ## result <- append(result,
+    ##                  list(Bgy.est.amelia.nok = est.g.mi,
+    ##                       Bgy.ci.upper.amelia.nok = est.g.mi + 1.96 * est.g.se,
+    ##                       Bgy.ci.lower.amelia.nok = est.g.mi - 1.96 * est.g.se
+    ##                       ))
 
     N <- nrow(df)
     m <- nrow(df[!is.na(x.obs)])
@@ -148,8 +148,8 @@ run_simulation <-  function(df, result){
                      )
 
 ##    clean up memory
-    rm(list=c("df","y","x","g","w","v","train","p","amelia.out.k","amelia.out.nok", "mod.calibrated.mle","gmm.res","mod.amelia.k","mod.amelia.nok", "model.true","model.naive","model.feasible"))
+##    rm(list=c("df","y","x","g","w","v","train","p","amelia.out.k","amelia.out.nok", "mod.calibrated.mle","gmm.res","mod.amelia.k","mod.amelia.nok", "model.true","model.naive","model.feasible"))
     
-    gc()
+##    gc()
     return(result)
 }

Community Data Science Collective || Want to submit a patch?