compare: Compare two Blimp models

In rblimp: Integration of 'Blimp' Software into R

Compare two Blimp models

Description

Compares two Bayesian models by calculating the proportion of posterior samples where the comparison model's parameters exceed (or fall below) the reference model's summary statistic. This is useful for model comparison and assessing incremental variance explained (e.g., R-squared differences).

Usage

compare(
  model0,
  model,
  use = "mean",
  greaterThan = TRUE,
  suffixes = c("R2: Coefficients", "R2: Level-2 Random Intercepts",
    "R2: Level-2 Random Slopes", "R2: Level-3 Random Slopes",
    "R2: Level-3 Random Intercepts", "R2: Residual Variation",
    "R2: Level-1 Residual Variation")
)

Arguments

model0	A blimp_obj. The baseline or simpler model used as the reference point.
model	A blimp_obj. The comparison model (typically more complex) to evaluate.
use	Summary statistic to use as the cutpoint from model0. Options: Character: "mean" or "median" Numeric < 1: Quantile proportion (e.g., 0.5 for median) Numeric >= 1: Percentile (e.g., 50 for median) Function: Custom function applied to model0 iterations List: Multiple summary statistics
greaterThan	Logical. If TRUE (default), calculates the proportion of model iterations greater than the model0 cutpoint. If FALSE, calculates proportion less than.
suffixes	Character vector of parameter name suffixes to compare. Defaults to all R-squared values (coefficients, random effects, residual variation).

Details

The comparison works by:

1. Computing a summary statistic (e.g., mean) from model0's posterior samples

2. Calculating what proportion of model's posterior samples exceed this value

3. Reporting this proportion for each parameter matching the specified suffixes

Value

A blimp_cp object containing a matrix of comparison proportions.

Note

Due to R restrictions, lists of functions will not give useful printed names.

Examples

# Generate data
mydata <- rblimp_sim(
    c(
        'x1 ~ normal(0, 1)',
        'x2 ~ normal(0, 1)',
        'y ~ normal(10 + 0.5*x1 + 0.3*x2, 1)'
    ),
    n = 200,
    seed = 123
)

# Fit baseline model (x1 only)
model0 <- rblimp(
    'y ~ x1',
    mydata,
    seed = 123,
    burn = 1000,
    iter = 1000
)

# Fit comparison model (x1 + x2)
model1 <- rblimp(
    'y ~ x1 x2',
    mydata,
    seed = 123,
    burn = 1000,
    iter = 1000
)

# Compare models - proportion of model1 R-squared > mean(model0 R-squared)
compare(model0, model1)

rblimp documentation built on May 18, 2026, 9:07 a.m.