Skip to contents

Run an advanced Laplace-consistency check

Source: R/check-consistency.R
gllvmTMB_check_consistency.Rd

gllvmTMB_check_consistency() is an advanced validation helper for the Laplace approximation, not a first-use diagnostic. Start with check_gllvmTMB() or gllvmTMB_diagnose(). Use this function when a fit looks numerically acceptable but you need a direct check of whether the approximate marginal score is centred.

Usage

gllvmTMB_check_consistency(fit, n_sim = 100L, seed = NULL, estimate = FALSE)

Arguments

fit

A fit returned by gllvmTMB().

n_sim

Integer number of simulate-evaluate replicates. Default 100L. Cost is roughly n_sim joint likelihood evaluations; budget 5-30 seconds on a Tier-1 fixture and proportionally longer for larger fits.

seed

Optional integer seed for the simulation RNG. NULL (default) draws a random seed.

estimate

Logical. When TRUE, also refits the model on each simulated dataset and reports the joint-score p_value + bias. Default FALSE (much faster; marginal-score check only).

Value

An object of class gllvmTMB_check_consistency with components:

$marginal_p_value

Joint chi-squared p-value across all marginal-score components. NA if TMB's information matrix could not be inverted (often the case on tiny fixtures).

$marginal_bias

Named numeric vector: per-parameter bias (mean of the simulated score, normalised by its SE). Large absolute values flag specific parameters where the Laplace approximation is biased.

$joint_p_value

Joint-score p-value. NA unless estimate = TRUE. Slower because each replicate refits.

$flagged_parameters

Character names of parameters whose marginal bias exceeds |bias| > 0.5 (a heuristic; finite-n_sim noise can produce small spurious biases, but |0.5| is generally beyond noise for n_sim = 100).

$n_sim

The n_sim actually run.

$threshold

The marginal-bias threshold used.

$diagnostics

Character vector of one or more of: "no_centring_warning" (the finite simulation did not detect non-centring; this is not an accuracy certificate), "marginal_score_non_centred", "joint_score_non_centred", "information_matrix_singular" (TMB couldn't invert it; interpret with caution), "marginal_p_value_unavailable" (TMB returned NA without a captured warning; usually happens on tiny / weakly-identified fixtures).

$raw

The full TMB::checkConsistency() return value in case the user wants the per-replicate gradient matrix.

$call

match.call() of the invocation.

Details

It simulates n_sim datasets from the fitted model under the joint parameter vector at the MLE, then tests whether the approximate marginal score function is centred at zero across the simulations. A non-centred score is a sign that the Laplace approximation is unreliable for this fit – typically because the random-effects posterior is far from Gaussian (saturated binomial / Beta, sparse counts, weakly identified random effects) or the data don't constrain the random effects well.

Scope: local Laplace-consistency diagnostics for a fitted model. This is a diagnostic signal rather than a substitute for bootstrap, profile, or external Bayesian calibration. Broader calibration evidence is not yet available.

This is a complementary signal to sanity_multi() (which checks Hessian definiteness, gradient magnitude, convergence flags) and to a target-explicit known-DGP study of rank-selection behaviour. Slower than the first-line checks – n_sim likelihood evaluations – but the only diagnostic that targets the Laplace approximation itself rather than the optimisation outcome or the parameter identification.

A large marginal-score p-value is an absence of evidence against centring in this finite simulation; it is not a certificate that the Laplace approximation is accurate. If the result is small (p_value <= 0.05) or large flagged-parameter list), consider:

  • Refitting with a richer random-effects structure (so the conditional posterior is more Gaussian).

  • Switching to a parametric bootstrap CI (which does not depend on the Laplace approximation being a good fit).

  • Validating against tmbstan::tmbstan(fit$tmb_obj) (the audit's recommended Bayesian-comparison path).

See also

sanity_multi() (structural / convergence checks), check_auto_residual() (mixed-family safeguard), gllvmTMB_diagnose() (holistic fit summary), TMB::checkConsistency() (the underlying TMB call).

Examples

if (FALSE) { # \dontrun{
fit <- gllvmTMB(value ~ 0 + trait +
                latent(0 + trait | site, d = 1),
                data  = sim$data,
                trait = "trait",
                unit  = "site")
res <- gllvmTMB_check_consistency(fit, n_sim = 50L, seed = 1)
res
res$marginal_bias
res$flagged_parameters
} # }