qgcomp_with_penalization.R

# simple example of using Bayes to penalize estimates in quantile g-computation

library(qgcomp)
library(MASS)
n= 50
Sigma <- matrix(c(10,3,3,2),2,2)
Sigma

set.seed(1232)
X = MASS::mvrnorm(n, mu=rep(0,2), Sigma)
qX = qgcomp::quantize(data.frame(X), c("X1", "X2"))$data
z = runif(n)
py = -5 + cbind(as.matrix(qX), z) %*% c(1.0, 1.0, -.5)
y_binom = rbinom(n, 1, plogis(py))


#"vanilla" qgcomp
qgcomp.noboot(y_binom~., expnms=c("X1", "X2"), data=data.frame(y_binom,X,z=z), family=binomial())

# gcomp with penalization using Normal priors on regression coefficients
# help via: ?arm::bayesglm
qgcomp.noboot(y_binom~., expnms=c("X1", "X2"), data=data.frame(y_binom,X,z=z), family=binomial(), 
              bayes=TRUE, 
              prior.mean.for.intercept = 0,
              prior.scale.for.intercept = 100.0,
              prior.mean = rep(0,3),      # prior mean for coefficients
              prior.scale = c(1, 1, 2.0) # prior standard devaition (e.g. stronger prior for mixture vs. covariate)
)

# true value of psi = 2.0