Confidence intervals for fitted model parameters

confint() returns confidence intervals for a fitted drmTMB model. Wald intervals are fast and are returned for fixed-effect coefficients and direct response-scale parameter targets by default. Direct Wald targets include constant residual-scale, random-effect standard-deviation, random-effect correlation, and constant residual-correlation rows when the fitted TMB parameter and TMB::sdreport() covariance are available. Correlation Wald intervals are computed on the fitted TMB correlation-link scale, equivalent to a guarded Fisher z/atanh transform, and then returned on the correlation scale. Bootstrap intervals simulate and refit direct targets. For positive scale and SD targets, percentile endpoints are taken on the fitted log scale before back-transforming to the response scale. Profile-likelihood intervals are slower because nuisance parameters are re-optimized; this first public profile path supports explicit fixed-effect, constant distributional-scale, random-effect standard-deviation, random-effect correlation, bivariate phylogenetic q=2 mean-mean correlation, block-diagonal bivariate phylogenetic mu1/mu2 and sigma1/sigma2 correlations, and constant residual-correlation targets. For predictor-dependent scale, residual-correlation, or currently supported corpair() formulae, supply newdata with parm = "sigma", parm = "rho12", or the fitted corpair(...) dpar to profile the fitted response-scale value for each supplied row.

Usage

# S3 method for class 'drmTMB'
confint(
  object,
  parm = NULL,
  level = 0.95,
  method = c("wald", "profile", "bootstrap"),
  newdata = NULL,
  trace = FALSE,
  profile_precision = c("default", "fast"),
  profile_maxit = NULL,
  profile_engine = c("auto", "endpoint", "tmbprofile"),
  profile_endpoint_max_eval = NULL,
  R = 199L,
  seed = NULL,
  parallel = c("none", "multicore"),
  workers = NULL,
  refit_control = NULL,
  sd_boundary = 1e-04,
  rho_boundary = 0.98,
  small_sample_df = c("location", "none", "group"),
  bias_correct = c("location", "none", "group"),
  ...
)

Arguments

object: A drmTMB fit.
parm: Optional character or integer vector selecting interval targets. NULL selects all direct Wald-ready targets for Wald intervals. Profile intervals require explicit target names or target-set shortcuts. Supported shortcuts are "fixed_effects", "random_effects", "variance_components", and "correlations".
level: Confidence level.
method: Interval method: "wald", "profile", or "bootstrap". If newdata is supplied and method is omitted, method = "profile" is used.
newdata: Optional data frame for response-scale profile intervals for predictor-dependent sigma, sigma1, sigma2, rho12, or fitted corpair() values. Each row is profiled separately by profiling its fixed-effect linear predictor and then transforming the interval to the response scale.
trace: Logical; passed to TMB::tmbprofile() when the tmbprofile profile engine is used.
profile_precision: Profile-control shortcut. "default" leaves TMB::tmbprofile() controls unchanged. "fast" supplies ystep = 0.5 and ytol = 2 unless the caller supplies those controls in ..., giving a quicker first-pass profile for long variance-component or correlation targets.
profile_maxit: Optional positive whole number passed to TMB::tmbprofile() as maxit when method = "profile". Use this as a per-target adaptive-step budget for long or exploratory profile runs.
profile_engine: Profile engine for direct fitted-object targets. "auto" uses a scalar endpoint solver for direct scale, SD, and correlation targets when no TMB::tmbprofile() controls are supplied, and otherwise uses TMB::tmbprofile(). "endpoint" requires the scalar endpoint solver, while "tmbprofile" preserves the previous full-profile route for comparison and debugging.
profile_endpoint_max_eval: Optional positive whole number limiting constrained endpoint evaluations per endpoint side when the scalar endpoint engine is used. This is a diagnostic escape hatch for long variance-component or correlation profiles; when the budget is reached the row is returned with conf.status = "profile_failed" and missing endpoints.
R: Number of parametric-bootstrap refits when method = "bootstrap".
seed: Optional seed for bootstrap simulation.
parallel: Profile or bootstrap backend: "none" or Unix "multicore". For profile intervals, targets or newdata rows are split across workers. For bootstrap intervals, refits are split across workers.
workers: Requested profile or bootstrap workers. If NULL and parallel = "multicore", drmTMB uses about half the detected CPU cores. Multicore execution is capped at the number of jobs and at 10 workers.
refit_control: Optional drm_control() object used for bootstrap refits. The default skips TMB::sdreport() and drops the TMB object because bootstrap intervals use refit point estimates.
sd_boundary, rho_boundary: Boundary thresholds used by method = "wald" to flag intervals where the symmetric Wald interval is unreliable: a variance-component standard deviation within sd_boundary of zero, or a correlation within rho_boundary of +/-1, is returned with conf.status = "wald_at_boundary" and a warning pointing to method = "profile". Defaults match check_drm() (1e-4, 0.98).
small_sample_df: Small-sample reference distribution for the method = "wald" interval. "location" (the default) references a t-quantile with df = g - 1 – a group-based, Satterthwaite-style between-group degrees of freedom – for each location-axis (mu, mu1, mu2) structured random-effect SD target (phylo, spatial, animal, relmat) with a resolvable group count g, and keeps the ordinary normal quantile z = qnorm((1 + level) / 2) for every other target. This widens only the location-axis structured-RE SD intervals, which under-cover at small g under the normal quantile; the t-quantile lifts their coverage and converges back to z as g grows. Dispersion (sigma, sigma1, sigma2) structured SD targets are deliberately left at the normal quantile because they already over-cover; non-structured and fixed-effect targets are never widened; and a plain labelled covariance block such as (1 + x | p | id) is also left at the normal quantile by default, because the correction magnitude was calibrated for structured blocks only. "none" is the opt-out: it uses the normal quantile for every target and is byte-identical to the pre-default behaviour. "group" widens every resolvable SD target – structured and labelled-covariance, location and dispersion axis alike; use it only when you deliberately want to t-reference the dispersion axis or a labelled covariance block too. Note that the t-quantile only widens the interval; it does not move its centre. The residual small-g gap is ML shrinkage bias in the variance-component point estimate (a job for REML), not a quantile problem, so the centre shift is handled by bias_correct.
bias_correct: Small-sample point-estimate bias correction for the method = "wald" interval. "location" (the default) adds a simulation-calibrated shift log(g / (g - 1)) to the log-scale point estimate – before back-transforming – of each location-axis (mu, mu1, mu2) structured random-effect SD target with a resolvable group count g (the same targets small_sample_df = "location" widens), and leaves every other target unshifted. The ML estimate of a structured-RE variance-component SD is biased low on the log scale by about log(g / (g - 1)) – a downward shrinkage that REML, or simply a larger g, removes asymptotically; this shift moves the interval centre up to counter that shrinkage. Used together with small_sample_df = "location", the centre shift and the t(df = g - 1) width act independently, giving exp((log(sigma_hat) + log(g / (g - 1))) +/- qt(p, g - 1) * se_log) and lifting small-g coverage of location-axis structured-RE SD targets to nominal. "none" is the opt-out: it leaves every centre at the ML estimate and is byte-identical to the pre-default behaviour. "group" shifts every resolvable SD target – structured and labelled-covariance, location and dispersion axis alike. The default deliberately excludes the dispersion axis (sigma SD intervals already over-cover, so the upward shift would push them further conservative) and plain labelled covariance blocks (whose correction magnitude is not yet simulation-calibrated).
...: Additional arguments passed to TMB::tmbprofile() when method = "profile" and the tmbprofile profile engine is used. drmTMB supplies the profiled obj, name, lincomb, and trace arguments internally; set the profile target with parm.

Value

A data frame with columns parm, level, lower, upper, scale, transformation, tmb_parameter, index, method, and profile.engine, conf.status, profile.boundary, and profile.message. Successful rows currently use conf.status = "wald", "profile", or "bootstrap". Failed numeric profile rows use "profile_failed" with missing endpoints; profile rows mark intervals that land near a lower SD boundary or correlation boundary. Bootstrap interval results carry a "bootstrap.diagnostics" attribute with one diagnostic row per refit and target, including refit convergence, target availability, draw use, and the refit message.

Details

Target names follow the profile target namespace. For fixed effects, use names such as "fixef:mu:x", "fixef:sigma:(Intercept)", or "fixef:rho12:w". Compact coefficient labels from summary(fit), such as "mu:x", are also accepted. Random-effect SD intervals are reported on the SD scale, and random-effect correlation intervals are reported on the correlation scale. For bivariate Gaussian fits with constant residual correlation, parm = "rho12" profiles the residual correlation and reports the interval on the response correlation scale. For fits with constant sigma, sigma1, or sigma2, parm = "sigma" and friends report response-scale intervals.

The fastest routine route is confint(fit), which uses Wald intervals for fixed effects and direct response-scale targets. For long phylogenetic, spatial, animal-model, or relatedness fits, profile only the needed variance-component or correlation rows with the default profile_engine = "auto" first; direct scalar scale, SD, and correlation targets use the endpoint engine when no full-profile controls are supplied. Use profile_engine = "tmbprofile" or profile_precision = "fast" when you want the previous full-curve TMB::tmbprofile() route for comparison, diagnostics, or control tuning.

References

The small-sample corrections applied by small_sample_df and bias_correct (whether the default "location" scope or the broader "group" scope) are motivated by established mixed-model theory but are simulation-calibrated, not derived: the log(g / (g - 1)) SD-scale centre shift is about twice the leading-order REML SD correction (0.5 * log(g / (g - 1))), matching the larger ML shrinkage measured on these structured/bivariate cells (their effective df is well below g - 1). The magnitude's authority is the per-model-class simulation in docs/design/219-structured-re-small-sample-bias-correction.md, not a single source. Relevant references:

Restricted maximum likelihood and the variance-component bias that motivates the centre shift (REML debiases the variance by g/(g-1); the shift here is on the SD log scale and ~2x the leading-order REML SD term):

Patterson, H. D., & Thompson, R. (1971). Recovery of inter-block information when block sizes are unequal. Biometrika, 58(3), 545-554. doi:10.1093/biomet/58.3.545

Harville, D. A. (1977). Maximum likelihood approaches to variance component estimation and to related problems. Journal of the American Statistical Association, 72(358), 320-338. doi:10.1080/01621459.1977.10480998

Searle, S. R., Casella, G., & McCulloch, C. E. (1992). Variance Components. New York: Wiley. doi:10.1002/9780470316856

Analytic and penalized first-order bias reduction of maximum likelihood estimates (the general framework an additive log-scale shift instantiates):

Cox, D. R., & Snell, E. J. (1968). A general definition of residuals. Journal of the Royal Statistical Society, Series B, 30(2), 248-265. doi:10.1111/j.2517-6161.1968.tb00724.x

Firth, D. (1993). Bias reduction of maximum likelihood estimates. Biometrika, 80(1), 27-38. doi:10.1093/biomet/80.1.27

The t-quantile / between-group effective degrees of freedom used by small_sample_df:

Satterthwaite, F. E. (1946). An approximate distribution of estimates of variance components. Biometrics Bulletin, 2(6), 110-114. doi:10.2307/3002019

Kenward, M. G., & Roger, J. H. (1997). Small sample inference for fixed effects from restricted maximum likelihood. Biometrics, 53(3), 983-997. doi:10.2307/2533558

Boundary regime (a variance component at or near zero), where neither the t-width nor the centre shift restores nominal coverage:

Self, S. G., & Liang, K.-Y. (1987). Asymptotic properties of maximum likelihood estimators and likelihood ratio tests under nonstandard conditions. Journal of the American Statistical Association, 82(398), 605-610. doi:10.1080/01621459.1987.10478472

Stram, D. O., & Lee, J. W. (1994). Variance components testing in the longitudinal mixed effects model. Biometrics, 50(4), 1171-1177. doi:10.2307/2533455

Parametric-bootstrap bias correction of mixed-model variance components, and the general delicacy of bootstrap bias estimation (a single-level parametric bootstrap does not recover the centre bias for these targets at small g):

Kubokawa, T., & Nagashima, B. (2012). Parametric bootstrap methods for bias correction in linear mixed models. Journal of Multivariate Analysis, 106, 1-16. doi:10.1016/j.jmva.2012.01.011

Efron, B., & Tibshirani, R. J. (1993). An Introduction to the Bootstrap. New York: Chapman & Hall.

Small-sample variance-component (repeatability) interval coverage in ecology and evolution:

Wolak, M. E., Fairbairn, D. J., & Paulsen, Y. R. (2012). Guidelines for estimating repeatability. Methods in Ecology and Evolution, 3(1), 129-137. doi:10.1111/j.2041-210X.2011.00125.x

Examples

dat <- data.frame(y = c(0.2, 0.5, 1.1, 1.4), x = c(-1, 0, 1, 2))
fit <- drmTMB(bf(y ~ x, sigma ~ 1), data = dat)
confint(fit)
#>                      parm level       lower      upper    scale
#> 1    fixef:mu:(Intercept)  0.95  0.51798631  0.6620137     link
#> 2              fixef:mu:x  0.95  0.36120108  0.4787990     link
#> 3 fixef:sigma:(Intercept)  0.95 -3.39479048 -2.0088865     link
#> 4                   sigma  0.95  0.03354758  0.1341380 response
#>     transformation tmb_parameter index method profile.engine conf.status
#> 1 linear_predictor       beta_mu     1   wald           <NA>        wald
#> 2 linear_predictor       beta_mu     2   wald           <NA>        wald
#> 3 linear_predictor    beta_sigma     1   wald           <NA>        wald
#> 4              exp    beta_sigma     1   wald           <NA>        wald
#>   profile.boundary profile.message
#> 1               NA            <NA>
#> 2               NA            <NA>
#> 3               NA            <NA>
#> 4               NA            <NA>
confint(fit, parm = "variance_components")
#>    parm level      lower    upper    scale transformation tmb_parameter index
#> 1 sigma  0.95 0.03354758 0.134138 response            exp    beta_sigma     1
#>   method profile.engine conf.status profile.boundary profile.message
#> 1   wald           <NA>        wald               NA            <NA>
confint(fit, parm = "sigma", method = "profile")
#>    parm level      lower     upper    scale transformation tmb_parameter index
#> 1 sigma  0.95 0.03817528 0.1645079 response            exp    beta_sigma     1
#>    method profile.engine conf.status profile.boundary profile.message
#> 1 profile       endpoint     profile            FALSE              ok
# Use the full-profile engine when you need the older tmbprofile route:
# confint(fit, parm = "sigma", method = "profile",
#   profile_engine = "tmbprofile", profile_precision = "fast")
# Direct-target parametric bootstrap is available when refit cost is worth it:
# confint(fit, parm = "sigma", method = "bootstrap", R = 99)
# Bootstrap intervals for positive scale and SD targets use link-scale
# percentiles before back-transforming to the response scale.