Skip to contents

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 ||.

Usage

phylo_coef(formula, tree = NULL, vcv = NULL, rho = NULL)

Arguments

formula

A coefficient-basis bar expression such as 1 + x | trait, 0 + x | trait, or 1 + x || trait.

tree

An ape::phylo tree whose tip labels match the response columns. Mutually exclusive with vcv.

vcv

A labelled positive-definite covariance matrix, or a labelled sparse precision matrix, for the response columns. Mutually exclusive with tree.

rho

NULL to estimate an interior phylogenetic mixture, or one numeric value in [0, 1] to fix it.

Value

A formula marker; never evaluated directly.

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."
#>