Skip to contents

Reduced-rank phylogenetic latent factors: phylo_latent(species, d = K)

Source: R/brms-sugar.R
phylo_latent.Rd

Canonical name for the reduced-rank phylogenetic random effect. Formerly phylo_rr(species, d = K) – same engine, new name.

Usage

phylo_latent(
  species,
  d = 1,
  tree = NULL,
  vcv = NULL,
  A = NULL,
  Ainv = NULL,
  unique = FALSE,
  rho = 1
)

Arguments

species

Unquoted column name for the species factor.

d

Integer; number of phylogenetic latent factors.

tree

An ape::phylo object. Canonical. Use this if you have a tree.

vcv

A tip-only phylogenetic correlation matrix (n_species x n_species). Legacy alias of A =.

A

Tip-level relatedness matrix (n_species x n_species) – alias of vcv =, aligned with the animal_* family's argument naming. Supply one of tree, vcv, or A / Ainv.

Ainv

Precision matrix (inverse of A). Sparse inputs are preserved for the sparse precision route.

unique

Logical; TRUE auto-includes the phylo-structured diagonal trait-specific \(\boldsymbol\Psi_{phy}\) companion, folding the shared rank-K loadings and the diagonal companion into a single term (\(\boldsymbol\Sigma_{phy} = \boldsymbol\Lambda \boldsymbol\Lambda^\top \otimes \mathbf{A} + \boldsymbol\Psi_{phy} \otimes \mathbf{A}\)). The default FALSE preserves the loadings-only / rotation-invariant subset.

rho

Source strength: a number in [0,1] fixes rho * K + (1-rho) * diag(diag(K)) on the legacy-resolved source scale. Omitted or explicit 1 preserves the existing model. NULL estimates strength for one Gaussian structured trait-intercept block with complete replicated multivariate observations and no competing ordinary covariance. Estimated latent terms require rank one and at least four traits. The same strength applies to the entire latent-plus-Psi covariance. This parameter is not a variance-share summary. Fixed attenuation and the admitted Gaussian estimator have implementation checks; recovery is regime-specific.

Value

A formula marker; never evaluated.

Details

Two phylogeny inputs: tree = (canonical) and vcv = (legacy)

Pass the phylogeny inside the keyword via one of two arguments:

  • tree = phylo (recommended when the tree is available) – the full ape::phylo object. The package constructs the sparse phylogenetic precision using its Hadfield–Nakagawa implementation.

  • vcv = Cphy – a tip-level n_species x n_species correlation matrix for analyses that begin from a supplied covariance matrix.

These inputs encode the same tip-level covariance target when they are constructed from the same tree and aligned identically. Numerical agreement still depends on labels, scaling, and fitting health; the function does not estimate or report ancestral states.

The two routes differ in one respect worth knowing. The vcv = route adds a fixed 1e-8 ridge to the supplied matrix before inverting it; the tree = route builds the sparse precision analytically from branch lengths and adds no ridge. The two therefore agree with each other only to roughly 1e-5 in log-density, and less on large or poorly conditioned trees.

Species are matched by FACTOR LEVEL, not by tree tip order

The species column is aligned to the phylogeny through the order of its factor levels: levels(data$species) is looked up against the tip labels, and whatever order those levels happen to be in is the order the latent scores are indexed by. Tree tip order is not consulted.

If the factor's levels do not correspond to the species you intend, the model still fits and issues no warning. It simply fits a different phylogeny. Set the levels explicitly:

data$species <- factor(data$species, levels = tree$tip.label)

Every example below does this. It is not decoration.

See the phylogenetic covariance article for the benchmark.

References

Hadfield JD, Nakagawa S (2010). General quantitative genetic methods for comparative biology: phylogenies, taxonomies and multi-trait models for continuous and categorical characters. J. Evol. Biol. 23: 494-508. doi:10.1111/j.1420-9101.2009.01915.x

See also

phylo_scalar(), phylo_indep(), phylo_dep(), phylo_rr() (deprecated alias).

Examples

if (FALSE) { # \dontrun{
  tree <- ape::rcoal(20); tree$tip.label <- paste0("sp", seq_len(20))
  sim <- simulate_site_trait(
    n_sites = 1, n_species = 20, n_traits = 4,
    mean_species_per_site = 20,
    Cphy = ape::vcv(tree, corr = TRUE),
    sigma2_phy = rep(0.3, 4), seed = 1
  )
  sim$data$species <- factor(sim$data$species, levels = tree$tip.label)
  fit <- gllvmTMB(
    value ~ 0 + trait + phylo_latent(species, d = 2),
    data       = sim$data,
    trait      = "trait",
    unit       = "species",
    cluster    = "species",
    phylo_tree = tree
  )
} # }