Learn R Programming

HOIFCar (version 1.1.1)

fit.JASA: Covariate Adjustment via JAckknife Score-based Adjustment (JASA) for Generalized Linear Working Models

Description

Implements the Jackknife Score-Based Adjustment (JASA) method and its calibration for covariate adjustment in simple randomized experiments where covariate dimension p may be large relatively to sample size n. Handles Continuous, Binary, and Poisson outcomes.

Usage

fit.JASA(
  Y,
  X,
  A,
  family = "gaussian",
  pi1 = NULL,
  is.parallel = FALSE,
  core_num = 4,
  opt_obj = c("beta", "mu")[1]
)

Value

List with components:

tau_vec

Named vector of average treatment effect estimates (JASA, JASA-cal)

tau0_vec

Treatment effect estimates on the control group

tau1_vec

Treatment effect estimates on the treatment group

var_tau_vec

Variance estimates for average treatment effect estimates

var_tau0_vec

Variance estimates for treatment effect estimates on the control group

var_tau1_vec

Variance estimates for treatment effect estimates on the treatment group

y_hat_mat

Matrix of predicted outcomes (columns: control, treatment) using varied Leave-One-Out strategy.

obj_value_mat

Matrix of objective values from optimization (columns: control, treatment)

Arguments

Y

Numeric vector of outcome values.

X

Matrix of centered baseline covariates (may include intercept). Dimensions: n x p.

A

Binary treatment vector (0 = control, 1 = treatment). Assignment is assumed to follow simple randomization.

family

GLM family specification: "gaussian", "binomial", or "poisson". Default: "gaussian".

pi1

The assignment probability for the simple randomization. If NULL (default), the empirical assignment probability is used.

is.parallel

Boolean for parallelization. Default: FALSE.

core_num

Number of cores for parallel computing (default is 4 if is.parallel = TRUE).

opt_obj

Ways to optimization: 'beta' (GLM coefficients) or 'mu' (expected outcomes). Default: 'beta'.

Examples

Run this code
generate_data_SR <- function(n, family, pi1, p_n_ratio = 0.05, seed = 123){
  set.seed(seed)
  alpha0 <- 0.15
  p0 <- ceiling(round(n * alpha0))
  beta0_full <- 1/(1:p0)^(1/4)*(-1)^c(1:p0)
  Sigma_true <- matrix(0, nrow = p0, ncol = p0)
  for (i in 1:p0) {
    for (j in 1:p0) {
      Sigma_true[i, j] <- 0.1 ** (abs(i - j))
    }
  }

  if(family != 'poisson'){
    X <- mvtnorm::rmvt(n, sigma = Sigma_true, df = 3)
  }else{
    X0 <- mvtnorm::rmvt(n, sigma = Sigma_true, df = 3)
    X <- pmin(pmax(X0, -3), 3)
    rm(X0)
  }

  beta <- beta0_full / norm(beta0_full,type='2')
  lp0 <- X %*% beta

  delta_X <- 1 - 1/2 * pmin(X[, 1]^2, 5) + 1/4 * X[, 1:10] %*% beta[1:10]
  lp1 <- lp0 + delta_X


  if (family == 'binomial') {
    r0 <- plogis(2 * lp0)
    r1 <- plogis(2 * lp1)
    Y1 <- rbinom(n, size=1, prob=r1)
    Y0 <- rbinom(n, size=1, prob=r0)
  }else if(family == 'poisson'){
    # quantile(lp1);quantile(lp0)
    lp1_tran <- pmin(lp1, 4)
    lp0_tran <- pmin(lp0, 4)
    r1 <- exp(lp1_tran)
    r0 <- exp(lp0_tran)

    Y1 <- rpois(n,r1)
    Y0 <- rpois(n,r0)
  }else if(family == 'gaussian'){
    r1 <- lp1;
    r0 <- lp0
    Y1 <- r1 + rnorm(n)
    Y0 <- r0 + rnorm(n)
  }

  A <- rbinom(n, size=1, prob=pi1)
  Y <- A * Y1 + (1 - A) * Y0

  p <- ceiling(round(n * p_n_ratio))
  if(p > ncol(X)){
    if(family != 'poisson'){
      X_noise <- rmvt(n, sigma = diag(p - ncol(X)), df = 3)
    }else{
      X0_noise <- rmvt(n, sigma = diag(p - ncol(X)), df = 3)
      X_noise <- pmin(pmax(X0_noise, -3), 3)
    }
    X_obs <- cbind(X, X_noise)
  }else{
    X_obs <- X[, 1:p, drop = FALSE]
  }

  data_ls <- list(
    X = X_obs, Y = Y, A = A,
    Y1 = Y1, Y0 = Y0,
    r1 = r1, r0 = r0
  )
  return(data_ls)
}


n <- 400; pi1 <- 1/3

family <- 'gaussian'; p_n_ratio <- 0.05
data_ls <- generate_data_SR(n, family, pi1, p_n_ratio)
X <- data_ls$X;
A <- data_ls$A
Y <- data_ls$Y

if (FALSE) {
Xc <- scale(X, scale = FALSE)
Xc_aug <- cbind(1, Xc)
result.jasa.ls <- fit.JASA(Y, Xc_aug, A, family, opt_obj = 'beta')
result.jasa.ls$tau_vec
result.jasa.ls$var_tau_vec


family <- 'poisson'; p_n_ratio <- 0.05
data_ls <- generate_data_SR(n, family, pi1, p_n_ratio)
X <- data_ls$X;
A <- data_ls$A
Y <- data_ls$Y

Xc <- scale(X, scale = FALSE)
Xc_aug <- cbind(1, Xc)
result.jasa.ls <- fit.JASA(Y, Xc_aug, A, family, opt_obj = 'mu')
result.jasa.ls$tau_vec
result.jasa.ls$var_tau_vec


family <- 'binomial'; p_n_ratio <- 0.05
data_ls <- generate_data_SR(n, family, pi1, p_n_ratio)
X <- data_ls$X;
A <- data_ls$A
Y <- data_ls$Y

Xc <- scale(X, scale = FALSE)
Xc_aug <- cbind(1, Xc)
result.jasa.ls <- fit.JASA(Y, Xc_aug, A, family, opt_obj = 'beta')
result.jasa.ls$tau_vec
result.jasa.ls$var_tau_vec
}

Run the code above in your browser using DataLab