Fit response-column-specific random intercepts and/or slopes with source
covariance K_rho = rho * K + (1 - rho) * diag(diag(K)). Supply a numeric
rho in [0, 1] to fix the mixture or use rho = NULL (the default) to
estimate one interior value. The fitted covariance across the coefficient
basis is full for | and diagonal for ||.
Arguments
- formula
A coefficient-basis bar expression such as
1 + x | trait,0 + x | trait, or1 + x || trait.- tree
An
ape::phylotree whose tip labels match the response columns. Mutually exclusive withvcv.- vcv
A labelled positive-definite covariance matrix, or a labelled sparse precision matrix, for the response columns. Mutually exclusive with
tree.- rho
NULLto estimate an interior phylogenetic mixture, or one numeric value in[0, 1]to fix it.
Details
This point-model route is covered for Gaussian multivariate data in long or
traits(...) wide form, with a labelled positive-definite tree covariance
source and bare numeric row predictors. Numeric rho and one estimated
interior rho are supported. Interval inference and non-Gaussian coefficient
models remain unavailable; animal_coef(), kernel_coef(), and
spatial_coef() cover their bounded public source regimes.
Existing phylo_slope() remains current and warning-free.
For exact compatibility with the released slope engine, a no-intercept
dense-vcv fit with rho = 1 uses the existing phylo_slope()
conditioning seam, K + 1e-8 I. Tree sources use their released sparse
precision. Interior fixed rho, estimated rho, and intercept-bearing
rho = 1 fits use the raw covariance-scale mixture shown above.
Estimating rho requires genuine between-column correlation contrast in
the standardized source; a diagonal source cannot identify the mixture and
is rejected. Supply exactly one of tree or vcv, never both.
Examples
set.seed(2)
dat <- expand.grid(unit = factor(1:12), trait = factor(paste0("sp", 1:4)))
dat$x <- rnorm(12)[dat$unit]
dat$value <- rnorm(nrow(dat))
K <- diag(4); dimnames(K) <- list(levels(dat$trait), levels(dat$trait))
fit <- gllvmTMB(value ~ 1 + phylo_coef(0 + x | trait, vcv = K, rho = 0.5),
data = dat, trait = "trait", unit = "unit", family = gaussian(),
control = gllvmTMBcontrol(se = FALSE), silent = TRUE)
extract_Sigma(fit, level = "column_coef")
#> $Sigma
#> x
#> x 0.322712
#>
#> $R
#> x
#> x 1
#>
#> $level
#> [1] "column_coef"
#>
#> $part
#> [1] "dep"
#>
#> $basis
#> [1] "x"
#>
#> $source
#> $source$type
#> [1] "phylo"
#>
#> $source$grouping
#> [1] "trait"
#>
#> $source$labels
#> [1] "sp1" "sp2" "sp3" "sp4"
#>
#>
#> $rho
#> [1] 0.5
#>
#> $rho_status
#> [1] "fixed"
#>
#> $K_rho
#> sp1 sp2 sp3 sp4
#> sp1 1 0 0 0
#> sp2 0 1 0 0
#> sp3 0 0 1 0
#> sp4 0 0 0 1
#>
#> $note
#> [1] "Response-column coefficient covariance; the response-column source supplies the other Kronecker factor."
#>
