Posterior-tree prediction sensitivity

Run pigauto's full imputation pipeline on each of T posterior phylogenies, generating m_per_tree datasets per tree, either as point predictions or stochastic completions. This returns T * m_per_tree datasets. Each completed dataset is conditional on a specific posterior tree (recorded in mi$tree_index).

Usage

multi_impute_trees(
  traits,
  trees,
  m_per_tree = 1L,
  species_col = NULL,
  trait_types = NULL,
  multi_proportion_groups = NULL,
  log_transform = TRUE,
  missing_frac = 0.25,
  covariates = NULL,
  epochs = 2000L,
  verbose = TRUE,
  seed = NULL,
  share_gnn = TRUE,
  reference_tree = NULL,
  draws_method = c("mc_dropout", "conformal"),
  gnn = FALSE,
  ...,
  ordinal_levels = NULL
)

# S3 method for class 'pigauto_mi_trees'
print(x, ...)

Arguments

traits: data.frame. Same format as multi_impute() and impute().
trees: list of phylo objects (class multiPhylo or plain list). For the default shared-reference workflow, supply trees with a common tip set, including the reference tree if supplied. With share_gnn = FALSE, tip sets may differ. Input species absent from the modeled tree, or from the reference tree in the shared workflow, remain in returned datasets but are not imputed there. Modeled tree tips without trait rows are represented internally as all-missing rows. Posterior samples from BirdTree.org (Jetz et al. 2012) are an example; the bundled trees300 dataset provides 50 posterior trees for avonet300.
m_per_tree: integer. Number of completed datasets returned per tree (default 1). Total datasets = length(trees) * m_per_tree. With draws_method = "mc_dropout" and the default of one, the prediction is deterministic conditional on each fitted model; this setting does not estimate within-tree variation.
species_col: character or NULL. See impute().
trait_types: named character vector overriding auto-detected trait types. Required for "proportion" and "zi_count". See impute() / preprocess_traits(). Default NULL (auto-detect).
multi_proportion_groups: named list declaring compositional trait groups (rows summing to 1), forwarded to impute() / preprocess_traits(). Default NULL.
log_transform: logical. Log-transform every continuous trait whose observed values are all positive (default TRUE), then back-transform the imputations. Set it to FALSE when traits are already on a log scale (for example log body mass): otherwise they are logged a second time and imputed on a log(log) scale, which can bias analyses done on the log scale. With a mix of raw and already-logged traits, log the raw ones yourself and set FALSE. In a real-data check (PanTHERIA, 4,027 mammal species, four already-logged traits, posterior draws), the default shifted Rubin-pooled regression slopes by up to 3.9 standard errors of the complete-data slope, in a direction that depended on which values were missing; with FALSE the largest shift was 1.8.
missing_frac: numeric. Fraction held out for validation/test during training (default 0.25).
covariates: data.frame or matrix of environmental covariates (fully observed), passed through to impute(). Numeric and integer columns are standardized; factor and ordered columns are one-hot encoded. Follow the row-alignment rules in impute() and provide covariates for modeled species or observations, including tree tips without trait rows. Default NULL (no covariates).
epochs: integer. Maximum GNN training epochs per tree (default 2000).
verbose: logical. Print progress (default TRUE).
seed: optional integer. Base random seed. How it is used depends on share_gnn and draws_method; see the "Reproducibility" section. The default NULL uses the current RNG stream.
share_gnn: logical. With gnn = TRUE, TRUE (the default) trains one GNN on a reference tree and reuses it while recomputing the baseline for each posterior tree. Set FALSE to refit the GNN on each tree. With the default gnn = FALSE, no GNN is trained in either mode. See the "Share-GNN" section below.
reference_tree: optional phylo used for the reference fit when share_gnn = TRUE. Default NULL selects the maximum-clade-credibility tree via phangorn::maxCladeCred(trees). If phangorn is not installed, falls back to trees[[1]] with a warning.
draws_method: character. How missing cells are completed. One of "mc_dropout" (default) or "conformal". The methods differ in whether they add stochastic variation; see the section below. A single mc_dropout completion is a deterministic prediction, while "conformal" runs its sampling step even when m_per_tree = 1.
gnn: logical. Passed through to impute() / fit_pigauto() for every tree. When TRUE, the usual GNN correction is trained (once, on the reference tree, when share_gnn = TRUE). When FALSE, no GNN is used anywhere in the pipeline: baseline-only fits throughout (the default since 0.11.0.9001; see impute()). With share_gnn = TRUE, each tree's baseline is fit without held-out splits (splits = NULL), matching the production baseline_full semantics fit_pigauto() uses.
...: additional arguments forwarded to fit_pigauto() via impute().
ordinal_levels: named list giving the scientific order of ordinal levels, coerced to ordered factors before fitting. See impute() / read_traits(). Default NULL.
x: A posterior-tree prediction-sensitivity object.

Value

An object of class "pigauto_mi_trees", inheriting from "pigauto_mi", with components:

datasets: List of T * m_per_tree data.frames. Observed cells are preserved. Missing cells are filled with tree-conditional predictions or draws where available. A cell marked for imputation can remain NA if no prediction is available. These datasets are for sensitivity diagnostics.
m: Total number of datasets (T * m_per_tree).
n_trees: Number of posterior trees used.
m_per_tree: Number of prediction or completion datasets per tree.
tree_index: Integer vector of length m; element i gives the tree index (1..T) for dataset i.
pooled_point: Single data.frame averaging numeric trait columns across all T * m_per_tree datasets; an NA in a numeric trait cell propagates to its average. Factor and ordered-factor columns are taken from the first dataset. Compare the datasets directly to describe nonnumeric prediction sensitivity. For reporting, not inference.
se: On the shared-GNN path, the arithmetic average of the per-tree prediction uncertainty summaries (see predict.pigauto_fit()). Descriptive reporting only, not a pooled inferential SE. NULL when unavailable, including the separate-fit path.
imputed_mask: Logical matrix aligned to the input rows and trait columns; TRUE marks an originally missing target in a row included in the fit, not necessarily a successful fill. A marked cell can remain NA if no prediction was available. When species_col is supplied, its column is excluded.
share_gnn: Logical; TRUE if the shared-GNN path was used.
draws_method: Character; "mc_dropout" or "conformal", echoing the argument this call used.
mi_workflow: "pigauto_tree_sensitivity_diagnostic", recording that these are prediction-sensitivity completions and cannot be passed to with_imputations() or pool_mi().
fit: Single pigauto_fit fitted on the reference tree when share_gnn = TRUE; NULL otherwise.
fits: List of T pigauto_fit objects (one per tree) when share_gnn = FALSE; NULL when share_gnn = TRUE.
reference_tree: The reference phylo used for the shared fit when share_gnn = TRUE; NULL otherwise.
trees: The input posterior trees.
species_col: Passed-through species column name.
lambda_per_trait_by_tree: Returned only when share_gnn = TRUE. A list indexed in the same order as trees; element i is the named lambda_per_trait vector returned by fit_baseline() for trees[[i]]. Values can be fixed rather than estimated. Categorical one-vs-rest lambda estimates are not included in this vector, so a categorical column's default value of 1 is not an estimated class-specific lambda.

The print method returns x invisibly.

Details

This is an experimental prediction-sensitivity check. Use it to see whether point predictions change across the trees you supply. It does not provide inferential uncertainty for downstream model estimates. See vignette("tree-uncertainty") for a step-by-step example.

Each dataset reflects the tree assigned to it. Comparing datasets across trees shows prediction sensitivity to the supplied tree sample. When multiple draws are requested per tree, they can also vary within a tree, depending on the selected draw method. These results are diagnostic, not analysis-ready imputations. The downstream fitting and pooling helpers reject them, so do not use them for inference. This workflow does not combine tree-to-tree variation with imputation variation into an inferential estimate.

For each tree the function runs preprocessing, baseline fitting and prediction. A GNN is also fitted only when gnn = TRUE; with share_gnn = TRUE, that GNN is trained once on a reference tree. Topologies and branch lengths vary across trees, so the phylogenetic baseline covariance differs for each tree.

The returned datasets may be compared descriptively to assess sensitivity of point imputations to the tree sample. No calibrated downstream standard error or Rubin-pooling claim is made for this path.

Computation time. The following historical planning estimates apply to gnn = TRUE, with one shared GNN fit plus T baseline passes when share_gnn = TRUE. They are not timings for the current GNN-off default; estimate the intended tree and sample setting before a campaign.

Species n	1 fit	T = 50 share_gnn=TRUE	T = 50 share_gnn=FALSE
300	~30-60 s	~3-5 min	25-50 min
5,000	~5-10 min	~10-20 min	4-8 hr
10,000	~20-40 min	~30-60 min	17-33 hr

When to use this

Choose the prediction-diagnostic function based on how many trees you have:

One tree (single published phylogeny, single time-calibrated tree): use multi_impute() and select a draw method explicitly when needed.
Multiple posterior trees (BirdTree samples, BEAST posterior, etc.): use multi_impute_trees() only to examine prediction sensitivity. Its results cannot be pooled for downstream inference.

Reproducibility

Pigauto routes the base seed separately to tree fits and sampling steps, so the offset depends on the fitting mode. Here, $t$ is a tree's position in the supplied list. With share_gnn = FALSE, the separate fit for tree $t$ receives seed + t - 1. With share_gnn = TRUE, the reference-tree fit receives seed, and each tree's baseline refit receives seed + t. Stochastic predictions reuse the stored seed from the reference fit; deterministic single-pass predictions do not set a new seed. If draws_method = "conformal", the tree-specific completion sampling step receives seed + t - 1; each completion within a tree also receives its own seed offset. With seed = NULL, these steps use the current random-number stream.

When gnn = TRUE, share_gnn = TRUE, neural weights, the reference graph, spectral features, and calibrated blend weights are reused from the reference fit (MCC by default). The baseline is recomputed for each posterior tree and passed into prediction again. Because baseline values also enter the neural inputs, shared weights do not imply a constant neural prediction across trees. No direction or magnitude of tree-variance bias has been established for this approximation.

Set share_gnn = FALSE with gnn = TRUE to fit a separate GNN for each tree. With the default gnn = FALSE, there is no GNN channel.

How the completion methods vary (and their limitations)

With draws_method = "mc_dropout" (the default), one completion per tree is a deterministic prediction conditional on that tree's fitted model. When m_per_tree > 1, pigauto adds independent Gaussian perturbations to latent cells, scaled by their baseline standard deviations. If gnn = TRUE, dropout can add further variation. These perturbations are not joint Brownian-motion posterior draws. Variation can be zero when a source has no uncertainty or its contribution is weighted out; decoding can also map different latent values to the same category.

With draws_method = "conformal", missing-value sampling occurs even when m_per_tree = 1. For continuous, count, ordinal and proportion traits, it uses conformal widths where available and reported uncertainty. When no finite conformal score is available for these traits, the route uses the reported pred$se as the fallback scale input. For log-transformed continuous traits, log1p-transformed counts, and logit-transformed proportions, pred$se is reported on the original scale, but the draw is made on the transformed scale. Pigauto does not fully convert pred$se between these scales, so diagnostic draw spread may be too wide or too narrow. Fitted point predictions are unchanged. These tree-based conformal completions use global trait scores; they do not carry per-cell near/far Mondrian scores into the draws, even when prediction intervals use Mondrian calibration. Binary and categorical traits are sampled from predicted probabilities; zero-inflated counts combine a Bernoulli gate with a magnitude draw; compositions use latent standard deviations. A fixed seed makes the draws reproducible. Neither method validates downstream inference.

When share_gnn = TRUE, conformal scores and calibrated blend weights come from the reference fit and are reused across trees. With share_gnn = FALSE, each tree uses scores and weights from its own fit.

Safety floor with a shared fit

When share_gnn = TRUE with safety_floor = TRUE, the grand-mean baseline mean_baseline_per_col and the three calibrated weights (r_cal_bm, r_cal_gnn, r_cal_mean) are computed ONCE on the reference tree and reused across all posterior trees. They are properties of the reference fit and are not recalibrated for each topology. Each posterior tree recomputes the baseline. The three weights stay fixed, while the GNN delta can change with its per-tree baseline inputs. This preserves prediction sensitivity to the per-tree baseline without re-calibrating the safety floor; it does not establish downstream inferential validity.

Printing

Prints tree, row, trait, and targeted-cell counts and the diagnostic-only scope. A targeted cell can remain missing when the prediction is unavailable. The draw-method label does not imply dropout ran with GNN off or with one completion per tree. The ... argument is currently unused.

References

Nakagawa S, de Villemereuil P (2019). "A general method for simultaneously accounting for phylogenetic and species sampling uncertainty via Rubin's rules in comparative analysis." Systematic Biology 68(4): 632-641.

Jetz W, Thomas GH, Joy JB, Hartmann K, Mooers AO (2012). "The global diversity of birds in space and time." Nature 491(7424): 444-448.

Examples

# \donttest{
data(avonet300, tree300)
tree <- ape::keep.tip(tree300, tree300$tip.label[seq_len(30L)])
trees <- structure(list(tree, tree), class = "multiPhylo")
df <- avonet300[match(tree$tip.label, avonet300$Species_Key),
                c("Mass", "Wing.Length"), drop = FALSE]
rownames(df) <- tree$tip.label
df$Mass[seq_len(3L)] <- NA_real_
mi <- multi_impute_trees(df, trees, m_per_tree = 1L, epochs = 5L,
                         verbose = FALSE)
#> Warning: Few tree-sensitivity datasets (T * m_per_tree = 2). Descriptive sensitivity summaries may be unstable.
#> Warning: Small validation set for 2 trait(s): Mass (n=3), Wing.Length (n=4). Calibrated gate and conformal scores will be noisy for these trait(s). 95% split-conformal coverage is NOT achievable for 2 trait(s) with fewer than 19 validation cells (Mass (n=3), Wing.Length (n=4)): the achievable ceiling is n_val / (n_val + 1), which only reaches 0.95 at n_val >= 19. See `?fit_pigauto` under 'Calibration at small n' for smoothing options.
print(mi)
#> pigauto experimental posterior-tree prediction sensitivity
#>   Trees     : 2 posterior phylogenies
#>   Per tree  : 1 dataset(s) (mc_dropout method; single deterministic prediction)
#>   Total     : 2 completed datasets
#>   Rows      : 30
#>   Traits    : 2 -- Mass, Wing.Length
#>   Cells     : 3 targeted / 60 total (5.0%)
#> 
#>   Access datasets:           mi$datasets[[i]]
#>   Tree index for dataset i: mi$tree_index[i]
#>   Downstream inference:     unsupported for this path
mass_by_tree <- vapply(mi$datasets, function(dat) dat$Mass,
                       numeric(nrow(df)))
apply(mass_by_tree, 1L, stats::sd)
#>  [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
# }