Fit the Manuscript-style Regime Integration pGLS — regime_integration

Reproduce the representative pGLS test used in the post-hoc integration analysis: log regime rate is modeled as a function of log post-hoc mean variance and Fisher-Z transformed mean absolute trait correlation on a collapsed regime phylogeny.

Usage

regime_integration_pgls(
  summary_data,
  search = NULL,
  tree = NULL,
  model = "BM",
  min_tips = NULL,
  ...
)

Arguments

summary_data: A data frame from summarize_regime_covariances() or a manuscript-compatible vars_cors table with columns rate, vars, corrs, and State. Regime IDs must be unique and non-empty.
search: A bifrost_search object containing the mapped regime tree.
tree: Optional SIMMAP-style mapped tree. Ignored when search is supplied.
model: Evolutionary model passed to phylolm::phylolm(). Defaults to "BM", matching the manuscript.
min_tips: Optional minimum tip count for downstream inclusion when summary_data contains a tip_count column. Downstream summaries are retained when tip_count >= min_tips; this differs from the strict PCA manuscript filter used by regime_correlation_pca().
...: Additional arguments passed to phylolm::phylolm().

Value

A phylolm fit.

Details

The collapse step follows the manuscript implementation: a monophyletic regime is represented by one collapsed tip, whereas all tips assigned to a nonmonophyletic regime are removed. When removals occur, the function emits one warning listing every dropped regime ID. The collapse stops if regime relabeling would create duplicated output tip labels. Standardization is computed across all summary rows surviving the optional min_tips filter before rows absent from the collapsed regime phylogeny are dropped. This intentional ordering matches manuscript preprocessing.

Examples

if (requireNamespace("phylolm", quietly = TRUE)) {
  tree <- ape::read.tree(text = paste0(
    "((((a:1,b:1):1,(c:1,d:1):1):1,",
    "((e:1,f:1):1,(g:1,h:1):1):1):1,(i:1,j:1):3);"
  ))
  tree <- phytools::paintSubTree(
    tree, node = ape::Ntip(tree) + 1L, state = "root"
  )
  tip_pairs <- list(
    r1 = c("a", "b"), r2 = c("c", "d"), r3 = c("e", "f"),
    r4 = c("g", "h"), r5 = c("i", "j")
  )
  for (regime in names(tip_pairs)) {
    tree <- phytools::paintSubTree(
      tree,
      node = ape::getMRCA(tree, tip_pairs[[regime]]),
      state = regime
    )
  }
  summary_data <- data.frame(
    regime = names(tip_pairs),
    rate = c(0.8, 1.1, 1.7, 2.2, 3.4),
    mean_variance = c(0.6, 0.9, 1.3, 2.0, 2.7),
    mean_abs_correlation = c(0.10, 0.18, 0.25, 0.32, 0.41)
  )
  fit <- regime_integration_pgls(summary_data, tree = tree)
  stats::coef(fit)
}
#>                      (Intercept)                 scale(log(vars)) 
#>                     -0.003447776                     -0.002379431 
#> scale(fisher_z_transform(corrs)) 
#>                      0.987224597