Introduction to causaldef

causaldef implements Le Cam deficiency theory to provide decision-theoretic diagnostics for causal inference.

Unlike standard sensitivity analysis tools that output abstract bias parameters (like E-values or partial R²), causaldef quantifies the information gap (\(\delta\)) between the data you have throughout your observational study and the data you wish you had (a perfect randomized trial).

The Core Workflow

  1. Specify: Define your causal problem using causal_spec().
  2. Estimate: Calculate a deficiency proxy \(\delta\) for your adjustment strategy using estimate_deficiency().
  3. Diagnose: Test assumptions with nc_diagnostic() or confounding_frontier().
  4. Decide: Compute safety bounds for decision making using policy_regret_bound() with a pre-specified adjustment method.

Example: Gene Perturbation Analysis

We’ll use a dataset of gene expression outcomes from a CRISPR knockout experiment. This case study illustrates how to detect unobserved confounding using negative controls.

library(causaldef)
library(ggplot2) # Ensure ggplot2 is loaded

# DAG Helper (using ggplot2 only)
plot_dag <- function(coords, edges, title = NULL) {
  edges_df <- merge(edges, coords, by.x = "from", by.y = "name")
  colnames(edges_df)[c(3,4)] <- c("x_start", "y_start")
  edges_df <- merge(edges_df, coords, by.x = "to", by.y = "name")
  colnames(edges_df)[c(5,6)] <- c("x_end", "y_end")
  
  ggplot2::ggplot(coords, ggplot2::aes(x = x, y = y)) +
    ggplot2::geom_segment(data = edges_df, ggplot2::aes(x = x_start, y = y_start, xend = x_end, yend = y_end),
                 arrow = ggplot2::arrow(length = ggplot2::unit(0.3, "cm"), type = "closed"), 
                 color = "gray40", size = 1, alpha = 0.8) +
    ggplot2::geom_point(size = 14, color = "white", fill = "#8FBC8F", shape = 21, stroke = 1.5) + # DarkSeaGreen
    ggplot2::geom_text(ggplot2::aes(label = name), fontface = "bold", size = 2.5, color = "white") +
    ggplot2::ggtitle(title) + 
    ggplot2::theme_void(base_size = 14) +
    ggplot2::theme(plot.title = ggplot2::element_text(hjust = 0.5, face = "bold", margin = ggplot2::margin(b = 10))) +
    ggplot2::coord_fixed()
}

data(gene_perturbation)

head(gene_perturbation)
#>   sample_id batch library_size knockout_status target_expression
#> 1        S1     3       948877        Knockout             8.643
#> 2        S2     4       984900         Control            11.545
#> 3        S3     4       868332         Control            12.325
#> 4        S4     2      1052105         Control            12.941
#> 5        S5     2       880271         Control             8.289
#> 6        S6     2      1001246        Knockout             9.400
#>   housekeeping_gene
#> 1             6.594
#> 2             9.424
#> 3             8.568
#> 4             9.224
#> 5             6.380
#> 6             6.706

# Visualize the Structural Assumption
coords <- data.frame(
  name = c("Unobserved", "Knockout", "TargetExpr", "Housekp"),
  x = c(0, -1.5, 1.5, 0),
  y = c(1, 0, 0, -1)
)
edges <- data.frame(
  from = c("Unobserved", "Unobserved", "Unobserved", "Knockout"),
  to =   c("Knockout", "TargetExpr", "Housekp", "TargetExpr")
)
plot_dag(coords, edges, title = "Gene Perturbation DAG")

1. Specification

We are interested in the effect of knockout_status on target_expression, adjusting for technical confounders batch and library_size.

2. Deficiency Estimation

How close can we get to a randomized experiment using these covariates? The following analysis compares the distributional distance (\(\delta\)) of the unadjusted observational data against the reweighted data using Inverse Probability Weighting (IPTW). Compare the deficiency values to determine the information gap.

Interpreting the PS-TV Deficiency Proxy (\(\delta\)):

The current proxy \(\delta\) quantifies how far the reweighted propensity-score distribution is from the target population:

What to do with these results:

  1. Choose the method with the lowest \(\delta\) (here, IPTW is better than unadjusted)
  2. Check the confidence interval (if bootstrap was used) to assess uncertainty
  3. If \(\delta\) is large, consider: (a) adding more covariates, (b) using a different adjustment method, or (c) acknowledging the limitation

3. Diagnose with Negative Controls

Is the adjustment sufficient? We use the housekeeping gene (which shouldn’t change) to test for residual confounding.

Interpreting Negative Control Results:

The negative control diagnostic tests whether your adjustment removes ALL confounding:

If the test fails:

  1. Add more covariates that might capture hidden confounding
  2. Try different adjustment methods (e.g., AIPW instead of IPTW)
  3. Consider sensitivity analysis to quantify the impact of unmeasured confounding
  4. Acknowledge the limitation in your conclusions

4. Decision Making

Suppose we need to decide whether to proceed with this gene target. What is the risk of being wrong due to remaining confounding?

Using policy_regret_bound(), we can compute two decision-theoretic quantities from a pre-specified deficiency proxy:

Interpreting Policy Regret Bounds:

The bounds answer: “How much regret could confounding add to my decision, and what regret is unavoidable?”

When results comes from estimate_deficiency(), these are plug-in bounds based on the PS-TV proxy rather than an identified exact deficiency.

Decision Rule:

Example: If your gene target could save 5 units of benefit, but the transfer penalty is 0.5 units (and the minimax safety floor is 0.25 units), you have roughly a 10:1 (and 20:1) benefit-to-risk ratio — often acceptable.

5. Effect Estimation

Finally, once we are confident in our adjustment strategy (low deficiency, passed diagnostics), we can estimate the causal effect.

Interpreting the Effect Estimate:

Reporting Guidelines:

  1. Report the effect with the method used: “Using IPTW adjustment, the estimated ATE is…”
  2. Reference the proxy: “The adjustment achieved \(\delta\) = 0.XX, indicating [excellent/good/moderate] approximation to an RCT under the PS-TV diagnostic”
  3. Mention diagnostic results: “The negative control diagnostic [passed/failed], suggesting [adequate/inadequate] control of confounding”
  4. State the bounds: “Transfer penalty is XX; minimax floor is YY on the chosen utility scale”

Complete Statement Example:
“Using IPTW adjustment (deficiency \(\delta\) = 0.01), we estimate the knockout reduces target expression by 2.3 units (95% CI: [1.8, 2.8]). The negative control diagnostic passed (p = 0.24). On a 0–10 utility scale (M = 10), the transfer penalty is 0.10 units and the minimax safety floor is 0.05 units, well below our risk tolerance threshold.”