Calculate Marginal Effects for Hierarchical Models — marginal

Computes the posterior distribution of average marginal effects by varying a target covariate over a grid.

Usage

marginal_effects(object, ...)

# S3 method for class 'rhierMnlRwMixture'
marginal_effects(object, z_values, Z, burn = 0, verbose = TRUE, ...)

Arguments

object: A fitted rhierMnlRwMixture object.
...: Additional arguments (currently ignored).
z_values: A numeric matrix where each row represents a specific grid point for evaluation. Columns correspond to the covariates in the model's original Z matrix. Non-NA values in a row fix the corresponding covariate to that value for the grid point. NA values indicate that the covariate should take its value from the corresponding row in the base Z matrix provided.
Z: A numeric matrix representing the base population or context over which the effects are averaged. It must have the same number of columns as required by the model (nz). For each grid point (row) in z_values, a counterfactual matrix is constructed based on Z, and the average effect across the rows of this counterfactual matrix is calculated.
burn: Integer, the number of initial MCMC draws to discard.
verbose: Logical, whether to print progress messages.

Value

The result depends on the method implementation.

A list object of class "marginal_effects" containing:

z_values: The user-provided matrix of grid points used for evaluation.
avg_betabar_draws: A list where each element [[i]] is an ncoef x ndraws_use matrix representing the posterior draws of the average betabar (DeltaZ + mu), evaluated at the i-th grid point (row of z_values) and averaged over the rows of the base Z matrix.
call: The matched call to the function.
burn: The number of burn-in draws discarded.

Details

This is a generic function. The main implementation for this package is marginal_effects.rhierMnlRwMixture.

Examples

# --- Simulate Data (requires bayesm.HART package) ---
if (requireNamespace("bayesm.HART", quietly = TRUE)) {
  sim_data <- bayesm.HART::sim_hier_mnl(nlgt = 50, nT = 10, p = 3, nz = 3, 
                                    nXa = 1, nXd = 0, seed = 123)
  Data <- list(p = sim_data$p, lgtdata = sim_data$lgtdata, Z = sim_data$Z)
  ncoef <- sim_data$true_values$dimensions$ncoef
  # --- Fit Model (minimal run for example) ---
  # Note: Use much larger R and keep for real analysis!
  Prior <- list(ncomp = 1)
  Mcmc <- list(R = 100, keep = 1)
  # Use try() to avoid errors stopping the example build if model fails
  fit <- try(bayesm.HART::rhierMnlRwMixture(Data = Data, Prior = Prior, Mcmc = Mcmc, 
                                        r_verbose = FALSE), silent = TRUE)
  if (!inherits(fit, "try-error")) {
    # --- Define Grid for Marginal Effects (vary Z1) ---
    target_z_index <- 1
    grid_z1 <- seq(min(Data$Z[, target_z_index]), 
                   max(Data$Z[, target_z_index]), 
                   length.out = 5)
    z_grid <- matrix(NA, nrow = length(grid_z1), ncol = ncol(Data$Z))
    z_grid[, target_z_index] <- grid_z1
    # --- Calculate Marginal Effects ---
    mfx_result <- marginal_effects(fit, 
                                     z_values = z_grid, 
                                     Z = Data$Z, 
                                     burn = 20, # Discard first 20 draws
                                     verbose = FALSE)
    print(names(mfx_result))
    print(dim(mfx_result$avg_betabar_draws[[1]])) # Check dimensions
    # --- Summarize Effects (see example for summary.marginal_effects) ---

    # --- Plot Effects (see example for plot.summary.marginal_effects) ---
  }
}
#>  MCMC Iteration (est time to end - min) 
#>  100 (0.0)
#>  Total Time Elapsed: 0.00 
#> [1] "z_values"          "avg_betabar_draws" "call"             
#> [4] "burn"             
#> [1]  3 80