X-Git-Url: https://code.communitydata.science/ml_measurement_error_public.git/blobdiff_plain/003733f22f42b435315803fd5f47d483c712d72d..588bdd7ed74cf8fe8fd0f15df58a6a40c26ebae5:/simulations/01_two_covariates.R

diff --git a/simulations/01_two_covariates.R b/simulations/01_two_covariates.R
index 419403d..7b8e12e 100644
--- a/simulations/01_two_covariates.R
+++ b/simulations/01_two_covariates.R
@@ -50,8 +50,8 @@ simulate_data <- function(N, m, B0=0, Bxy=0.2, Bgy=-0.2, Bgx=0.2, y_explained_va
     }
 
     df <- df[,w_pred:=x]
-
     df <- df[sample(1:N,(1-prediction_accuracy)*N),w_pred:=(w_pred-1)**2]
+    w <- predict(glm(x ~ w_pred,data=df,family=binomial(link='logit')),type='response')
     df <- df[,':='(w=w, w_pred = w_pred)]
     return(df)
 }
@@ -61,15 +61,20 @@ parser <- add_argument(parser, "--N", default=500, help="number of observations
 parser <- add_argument(parser, "--m", default=100, help="m the number of ground truth observations")
 parser <- add_argument(parser, "--seed", default=4321, help='seed for the rng')
 parser <- add_argument(parser, "--outfile", help='output file', default='example_1.feather')
+parser <- add_argument(parser, "--y_explained_variance", help='what proportion of the variance of y can be explained?', default=0.005)
+parser <- add_argument(parser, "--gx_explained_variance", help='what proportion of the variance of x can be explained by g?', default=0.15)
+parser <- add_argument(parser, "--prediction_accuracy", help='how accurate is the predictive model?', default=0.73)
+
 args <- parse_args(parser)
 
 B0 <- 0
 Bxy <- 0.2
 Bgy <- -0.2
-Bgx <- 0.5
+Bgx <- 0.4
+
+df <- simulate_data(args$N, args$m, B0, Bxy, Bgy, Bgx, seed=args$seed, y_explained_variance = args$y_explained_variance, gx_explained_variance = args$gx_explained_variance, prediction_accuracy=args$prediction_accuracy)
 
-df <- simulate_data(args$N, args$m, B0, Bxy, Bgy, Bgx, seed=args$seed, y_explained_variance = 0.025, gx_explained_variance = 0.15)
-result <- list('N'=args$N,'m'=args$m,'B0'=B0,'Bxy'=Bxy,'Bgy'=Bgy, 'Bgx'=Bgx, 'seed'=args$seed)
+result <- list('N'=args$N,'m'=args$m,'B0'=B0,'Bxy'=Bxy,'Bgy'=Bgy, 'Bgx'=Bgx, 'seed'=args$seed, 'y_explained_variance' = args$y_explained_variance, 'gx_explained_variance' = args$gx_explained_variance, "prediction_accuracy"=args$prediction_accuracy)
 outline <- run_simulation(df, result)
 
 outfile_lock <- lock(paste0(args$outfile, '_lock'),exclusive=TRUE)