Skip to contents

Standard errors of random-effect blocks without a scores accessor

Source: R/re-uncertainty.R
getREsd.Rd

getLV() exposes standard errors for the ordinary latent-variable score blocks (z_B / z_W). Several other random-effect families fitted by gllvmTMB() have the identical machinery available – a marginal variance for every random-effect coordinate, already computed by the fit's TMB sdreport() – but no accessor exposed it. getREsd() closes that gap for the blocks whose point-estimate reshape convention is already established elsewhere in the package (see block below).

Usage

getREsd(
  fit,
  block = c("diag_unit", "diag_unit_obs", "diag_species", "phylo", "re_int", "equalto")
)

Arguments

fit

A fitted multivariate model returned by gllvmTMB().

block

Which random-effect block to read. One of:

  • "diag_unit" – the unit-level (site) diagonal random effect from an indep(0 + trait | <unit>) / diag() term (s_B in the TMB template). Returns an n_traits x n_sites matrix, matching the orientation predict.gllvmTMB_multi() already uses when it reads this block's point estimate.

  • "diag_unit_obs" – the within-unit diagonal random effect from an indep(0 + trait | <unit_obs>) term (s_W). Returns an n_traits x n_site_species matrix in the same orientation.

  • "diag_species" – the species-level diagonal random effect from an indep(0 + trait | <species>) term (q_sp). Returns an n_traits x n_species matrix.

  • "phylo" – the propto() per-species phylogenetic random effect (p_phy). Returns an n_species x n_traits matrix, matching the orientation predict.gllvmTMB_multi() already uses when it reads this block's point estimate.

  • "re_int" – the (1 | group) random intercept(s) (u_re_int). Returns a named list, one numeric vector per bar term (named by that term's grouping-column name), because multiple (1 | group) terms are packed into a single flat parameter vector.

  • "equalto" – the equalto() known-covariance random effect (e_eq). Returns a plain numeric vector, one entry per observation row (the block's native shape; there is no unit/trait grid to reshape into).

Value

As described per block above. Standard errors are always NA (with a warning) when the fit's Hessian is not positive-definite.

What this is, and is not

Every number returned here is a Wald / sdreport() marginal standard error, propagated through the delta method the same way getLV()'s se = TRUE already is – not a resampled, profiled, or simulation-based quantity. Coverage of intervals built from it has not been measured for these blocks; do not describe or advertise it as exact, and do not construct a CI from it and call the CI "certified" without separately measuring its coverage. diag.cov.random reflects fixed-effect-uncertainty -propagated marginal variance for each coordinate individually; it does NOT give cross-block or random-vs-fixed covariance (that needs sdreport(getJointPrecision = TRUE), which this fit's single production sdreport() call does not request). For eligible Gaussian fits, the ordinary cell effects are integrated out analytically. Their standard errors combine the propagated uncertainty of their reconstructed conditional means with their conditional variances. This preserves the same first-order uncertainty convention. Their removed coordinates are not present in the reduced tape's raw joint precision; cross-covariances involving these cells also require conditional reconstruction, rather than only getJointPrecision = TRUE.

Because of that fixed-effect-uncertainty propagation, a near-saturated or otherwise degenerate design (e.g. an indep() diagonal term at a grouping where every group has one observation, with its residual variance driven to the boundary) can return SEs one or more orders of magnitude larger than the block's own conditional (GMRF) variance would suggest – in that regime the random effect deterministically tracks the fixed-effect residual, so its delta-method SE inherits the fixed effect's own uncertainty almost entirely. This is correct Wald behaviour, not a bug, but it means the magnitude is not always "the block's own noise" in an intuitive sense; treat a surprisingly large SE as a cue to check whether the corresponding variance component sits at or near a boundary.

Blocks not covered here (s_B_slope, s_W_slope, r_c2 / diagonal cluster2, g_phy, g_phy_diag, the SPDE spatial fields omega_spde*, the dense-kernel blocks g_kernel*, every augmented random-slope block, and the missing-predictor latent blocks x_mis / u_mi_group / g_x) are, in principle, reachable the same way – the underlying sdreport() call already covers their marginal variances too – but either their point-estimate reshape convention is not yet established anywhere else in the package (so getting the row/column order right cannot be cross-checked against existing code), or their dimensions are not carried as top-level fields on the fit object. Adding them safely needs that groundwork first; this function raises a clear error naming the block rather than silently mis-reshaping.

See also

getLV() for the ordinary latent-score standard errors this function generalises.

Examples

if (FALSE) { # \dontrun{
getREsd(fit, block = "re_int")
getREsd(fit, block = "diag_unit")
} # }