Pedigree PEV, reliability, and trajectory uncertainty

The pedigree controls accept prediction_error = "none" (the default), "selected_inverse", or "selected_solve". Exact uncertainty is available for the single-trait, incomplete-multitrait, and longitudinal pedigree interfaces when a compatible accepted sparse-direct factor is retained. selected_inverse remains CHOLMOD-only. selected_solve is qualified for PARDISO and may be used with CHOLMOD for parity.

The example assumes prepared single-trait records and a supported PARDISO build. Omitted starts use the existing initializer.

control <- gs_pedigree_reml_control(
  linear_solver = "sparse_direct",
  direct_backend = "pardiso",
  prediction_error = "selected_solve",
  prediction_error_animals = c("animal-17", "animal-3"),
  maximum_prediction_error_target_systems = 64L,
  prediction_error_rhs_batch_size = 8L,
  maximum_prediction_error_entries = 4096L,
  maximum_prediction_error_memory_bytes = 8 * 1024^3)
fit <- gs_pedigree_fit(records, control = control)
stopifnot(isTRUE(fit$status$ok), isTRUE(fit$optimizer_converged),
  isTRUE(fit$prediction_error$status$ok),
  identical(fit$prediction_error$availability, "available"))
fit$prediction_error[c("prediction_error_variance", "reliability")]

The result is conditional plug-in uncertainty from the random-effect block of the inverse of the complete MME. Single-trait reliability uses the pedigree kernel’s animal-specific \(A_{ii}=1+F_i\). Multitrait output contains full within-animal covariance blocks. Longitudinal output contains coefficient blocks and PEV, SEP, prior variance, reliability, and accuracy for observed and prediction basis rows in caller order.

For a single-trait animal effect, the reported quantities satisfy

\[ \mathrm{Rel}_i = 1 - \frac{\mathrm{PEV}_i}{\hat\sigma_a^2 A_{ii}}, \qquad \mathrm{SEP}_i = \sqrt{\mathrm{PEV}_i}, \qquad \mathrm{Accuracy}_i = \sqrt{\mathrm{Rel}_i}. \]

These calculations condition on the supplied or estimated covariance parameters; they do not include uncertainty from estimating those parameters. Accuracy here is a model-based quantity, not an observed correlation with known breeding values or a cross-validation result. Only roundoff-sized violations are adjusted and reported; materially invalid PEV or reliability values produce an error status.

The data-frame gfit() pedigree route accepts the corresponding named pedigree controls and retains uncertainty in diagnostics(fit)$native$prediction_error. The separate compiled route selected by gs_fit_control() currently requires prediction_error = "none". The persistent and component interfaces have their own uncertainty contracts; these pedigree settings do not apply to every gs_fit object.

For selected_solve, let T be the canonical coefficient targets implied by the requested animals and let E_T contain their coordinate columns. The native calculation solves \(CX=E_T\) in bounded column-major batches and extracts the target rows, giving exactly \(C^{-1}[T,T]\). Targets are sorted for the internal solve and mapped back to the caller’s requested animal and trait order. Duplicate or invalid targets are rejected.

This is not selected inversion. It uses ordinary documented retained-factor multi-RHS solves, and its work grows with the number \(K=\lvert T\rvert\) of target coefficient columns. Preflight reports K systems, \(\lceil K/\texttt{batch\_size}\rceil\) batches, \(K^2\) returned inverse entries, RHS workspace, returned-block workspace, projected peak workspace, configured system/entry/memory limits, factor validity, and direct/support threads. An over-limit request fails before the first target solve. Neither a dense MME nor a full dense inverse is built.

Both exact routes perform no additional symbolic analysis, numerical factorization, response solve, likelihood, trace, score, or AI evaluation. Target solves are counted separately from response, trace, and AI solves. PCG, an incompatible direct backend, or an invalid retained factor fails explicitly; there is no coordinate, stochastic, backend, or prediction-error-method fallback. CHOLMOD selected inversion remains the scalable exact route for large requested inverse closures, while PARDISO selected_solve is intended for resource-bounded target sets.