knitr::opts_chunk$set(echo = TRUE)
coverage <- c() for(i in 1:200){ n <- 1500 W1 <- runif(n, -1 , 1) W2 <- runif(n, -1 , 1)# rbinom(n, 1 , plogis(W1)) W <- cbind(W1, W2) A <- rbinom(n, 1 , plogis(0.5*(W1 + W2 ))) Y <- rpois(n, plogis( W1 + W2 ) * exp(A*(1+ 0.5*(W1 + W2)))) quantile(Y) LRR <- (1 + W1 + W2) true_coefs <- c(1, 0.5, 0.5) R <- as.numeric(1:n %in% c(which(Y>0), which(Y==0)[rbinom(sum(Y==0), size = 1, prob = 0.5)==1])) R <- rep(1,n) pR0 <- mean(R[Y==0]) pR1 <- mean(R[Y>0]) weights <- R / ifelse(Y>0, pR1, pR0) sl3_Learner_EYAW <- Lrnr_hal9001$new(family = "poisson", smoothness_orders =1, num_knots = c(10,5), max_degree = 2) sl3_Learner_pA1W <- Lrnr_gam$new() fit_msm <- npRRWorkingModel(~ 1 + W1 + W2, W, A, Y, weights, EY1W, EY0W, pA1W, sl3_Learner_EYAW, sl3_Learner_pA1W, outcome_type = "continuous") coverage <- cbind(coverage, fit_msm$coefficients$ci_right >= true_coefs & fit_msm$coefficients$ci_left <= true_coefs) } rowMeans(coverage)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.