update real data examples code and rerun project.

[ml_measurement_error_public.git] / simulations / simulation_base.R

diff --git a/simulations/simulation_base.R b/simulations/simulation_base.R
index 08b11ec9595a49a34553b3b1eaaaa3fb1463e27f..e715edfaf61b9b1f423001808187ad66828a43b0 100644 (file)
--- a/simulations/simulation_base.R
+++ b/simulations/simulation_base.R
@@ -180,26 +180,35 @@ run_simulation_depvar <- function(df, result, outcome_formula=y~x+z, proxy_formu
     
 
     # amelia says use normal distribution for binary variables.
     
 
     # amelia says use normal distribution for binary variables.

-
-    amelia.out.k <- amelia(df, m=200, p2s=0, idvars=c('y','ystar','w'))
-    mod.amelia.k <- zelig(y.obs~x+z, model='ls', data=amelia.out.k$imputations, cite=FALSE)
-    (coefse <- combine_coef_se(mod.amelia.k, messages=FALSE))
-    est.x.mi <- coefse['x','Estimate']
-    est.x.se <- coefse['x','Std.Error']
-    result <- append(result,
-                     list(Bxy.est.amelia.full = est.x.mi,
+    amelia_result <- list(Bxy.est.amelia.full = NA,
+                          Bxy.ci.upper.amelia.full = NA,
+                          Bxy.ci.lower.amelia.full = NA,
+                          Bzy.est.amelia.full = NA,
+                          Bzy.ci.upper.amelia.full = NA,
+                          Bzy.ci.lower.amelia.full = NA
+                          )
+
+    tryCatch({
+        amelia.out.k <- amelia(df, m=200, p2s=0, idvars=c('y','ystar','w'))
+        mod.amelia.k <- zelig(y.obs~x+z, model='ls', data=amelia.out.k$imputations, cite=FALSE)
+        (coefse <- combine_coef_se(mod.amelia.k, messages=FALSE))
+        est.x.mi <- coefse['x','Estimate']
+        est.x.se <- coefse['x','Std.Error']
+
+        est.z.mi <- coefse['z','Estimate']
+        est.z.se <- coefse['z','Std.Error']
+        amelia_result <- list(Bxy.est.amelia.full = est.x.mi,

                           Bxy.ci.upper.amelia.full = est.x.mi + 1.96 * est.x.se,
                           Bxy.ci.upper.amelia.full = est.x.mi + 1.96 * est.x.se,

-                          Bxy.ci.lower.amelia.full = est.x.mi - 1.96 * est.x.se
-                          ))
-
-    est.z.mi <- coefse['z','Estimate']
-    est.z.se <- coefse['z','Std.Error']
-
-    result <- append(result,
-                     list(Bzy.est.amelia.full = est.z.mi,
+                          Bxy.ci.lower.amelia.full = est.x.mi - 1.96 * est.x.se,
+                          Bzy.est.amelia.full = est.z.mi,

                           Bzy.ci.upper.amelia.full = est.z.mi + 1.96 * est.z.se,
                           Bzy.ci.lower.amelia.full = est.z.mi - 1.96 * est.z.se
                           Bzy.ci.upper.amelia.full = est.z.mi + 1.96 * est.z.se,
                           Bzy.ci.lower.amelia.full = est.z.mi - 1.96 * est.z.se

-                          ))
+                          )
+    },
+    error = function(e){
+    result[['error']] <- e}
+    )
+    result <- append(result,amelia_result)

 
     return(result)
 
 
     return(result)