Extract a cross-lineage Gamma block — extract

extract_Gamma() slices the shared covariance matrix returned by extract_Sigma() to give the row-lineage by column-lineage block Gamma_shape = Lambda_row Lambda_col^T. In a fitted cross-lineage kernel model, level is the named kernel_*() tier, row_traits are the host traits, and col_traits are the partner traits. This is a point-estimate extractor for the shared covariance block of a fitted dense-kernel model, covered by the package's validation tests at the fixed-rho scale: when the fitted kernel tier was built from make_cross_kernel(), scale = "effect" returns Gamma_effect = rho * Gamma_shape using the fixed rho recorded on that kernel. Uncertainty for Gamma is not yet reported by this helper; use bootstrap workflows when interval estimates are required. Use profile_cross_rho() to compare a fixed rho grid; in-engine estimation of this kernel rho and profile intervals remain planned work. (This is a different quantity from the source-strength rho that kernel_latent(..., rho = NULL) estimates.)

Usage

extract_Gamma(fit, level, row_traits, col_traits, scale = c("shape", "effect"))

Arguments

fit: A fitted gllvmTMB_multi object.
level: Character scalar naming the covariance tier. For kernel_*() fits this is the name argument supplied in the formula.
row_traits, col_traits: Character vectors of trait names defining the rows and columns of the returned block.
scale: Character scalar. "shape" (default) returns Gamma_shape = Lambda_row Lambda_col^T. "effect" returns rho * Gamma_shape for fixed cross-lineage kernels built by make_cross_kernel(). The current engine does not estimate rho.

Value

A numeric matrix with rows row_traits and columns col_traits.

Details

rho is part of the supplied K matrix, not a fitted parameter in the current engine. To profile it, use profile_cross_rho() to rebuild K_star over a small rho grid, refit the same formula, and compare logLik() values. scale = "effect" is therefore a fixed-kernel transformation, not an estimate of rho. Treat Gamma from a single association matrix W as data-condition sensitive: the package's recovery tests include a sparse-versus-dense W check, and sparse or poorly replicated host-partner links should be reported as weaker evidence rather than as a precise coevolution estimate.

For fixed multi-kernel fits, extract_Gamma() also inspects the fitted kernel-similarity diagnostics. If the requested component participates in a high-overlap kernel pair, the returned block is still available, but a warning reminds the caller that component-specific separation is weak evidence.

Examples

set.seed(1)
## A small host-partner system: 12 host and 12 partner species. Each
## lineage has its own similarity matrix. Here similarity decays with
## distance along a species ordering, a stand-in for a phylogenetic
## correlation matrix such as ape::vcv(tree, corr = TRUE).
host <- paste0("H", 1:12)
partner <- paste0("P", 1:12)
A_H <- 0.6^abs(outer(1:12, 1:12, "-"))
A_P <- A_H
dimnames(A_H) <- list(host, host)
dimnames(A_P) <- list(partner, partner)
## W records which hosts interact with which partners: host i interacts
## with partners i - 1, i and i + 1.
W <- 1 * (abs(outer(1:12, 1:12, "-")) <= 1)
dimnames(W) <- list(host, partner)
## One kernel over all 24 species; the fixed rho sets how strongly
## interacting hosts and partners are linked.
K <- make_cross_kernel(A_H, A_P, W, rho = 0.9)

## Simulate one latent axis shared by all four traits. Host traits are
## measured on hosts only and partner traits on partners only: four
## replicate measurements per species, each with noise (SD 0.3).
lambda <- c(h_size = 0.9, h_defence = 0.6, p_size = 0.8, p_attack = -0.5)
z <- drop(t(chol(K)) %*% rnorm(24))
names(z) <- rownames(K)
dat <- expand.grid(rep = 1:4, trait = names(lambda),
                   species = rownames(K), stringsAsFactors = FALSE)
dat <- dat[startsWith(dat$species, "H") == startsWith(dat$trait, "h_"), ]
dat$value <- lambda[dat$trait] * z[dat$species] +
  rnorm(nrow(dat), sd = 0.3)
dat$trait <- factor(dat$trait, levels = names(lambda))
dat$species <- factor(dat$species, levels = rownames(K))
## Each record (one species, one replicate) is a sampling unit
dat$obs <- factor(paste(dat$species, dat$rep))

## One latent axis whose species scores are correlated as in K. The
## term is named "cross"; the functions below refer to it by name.
fit <- gllvmTMB(
  value ~ 0 + trait + kernel_latent(species, K = K, d = 1, name = "cross"),
  data = dat, trait = "trait", unit = "obs", cluster = "species"
)
## Gamma, the host-by-partner block of the fitted trait covariance: how
## each host trait covaries with each partner trait. It is a point
## estimate; this helper reports no uncertainty.
extract_Gamma(fit, level = "cross",
              row_traits = c("h_size", "h_defence"),
              col_traits = c("p_size", "p_attack"))
#>              p_size   p_attack
#> h_size    0.7668029 -0.4140754
#> h_defence 0.5121792 -0.2765780
## The simulated truth, for comparison
outer(lambda[1:2], lambda[3:4])
#>           p_size p_attack
#> h_size      0.72    -0.45
#> h_defence   0.48    -0.30
## The same block multiplied by the kernel's fixed rho (0.9)
extract_Gamma(fit, level = "cross",
              row_traits = c("h_size", "h_defence"),
              col_traits = c("p_size", "p_attack"), scale = "effect")
#>              p_size   p_attack
#> h_size    0.6901226 -0.3726679
#> h_defence 0.4609613 -0.2489202