Random-regression (reaction-norm) extractors

Usage

rr_covariance(object, ...)

random_coefficients(object, ...)

rr_genetic_variance(object, at = NULL, n = 25L, ...)

rr_heritability(object, at = NULL, n = 25L, ...)

rr_correlation(object, at = NULL, n = 25L, ...)

rr_eigenfunctions(object, at = NULL, n = 25L, ...)

Arguments

object: A random-regression hsquared_fit object.
...: Reserved for future arguments.
at: Covariate values on the original scale at which to evaluate the trajectory. NULL (the default) uses an evenly spaced grid over the fitted covariate range. A supplied value outside the fitted covariate range errors (see the "Out-of-range at" section below).
n: Number of grid points used when at = NULL.

Value

rr_covariance() returns a numeric matrix; random_coefficients() returns a data frame; the trajectory extractors return a data frame with a covariate column and the evaluated values.

Details

These extractors summarize an opt-in, experimental random-regression (reaction-norm) fit (target = "random_regression"), fitted with an animal(rr(covariate, order = k) | id, pedigree = ped) term. The model estimates a k x k genetic covariance matrix K_g among an animal's normalized-Legendre random-regression coefficients plus a single homogeneous residual variance.

rr_covariance() returns the estimated k x k coefficient genetic covariance matrix K_g.
random_coefficients() returns the per-animal predicted Legendre coefficients (long format: id, coefficient, value).
rr_genetic_variance() returns the additive genetic variance trajectory v_g(t) = phi(t)' K_g phi(t) across covariate points.
rr_heritability() returns the heritability trajectory h^2(t) = v_g(t) / (v_g(t) + sigma_e^2). Because the residual is homogeneous and there is no permanent-environment term yet, this can OVERSTATE h^2(t) for repeated-records designs (test-day, growth curves).
rr_correlation() returns the genetic correlation surface among the covariate points.
rr_eigenfunctions() returns the eigen-decomposition of the coefficient covariance K_g evaluated as functions of the covariate: the eigenvalues (additive genetic variance carried by each axis), the proportion of genetic variance explained, the sign-canonicalized eigen-coefficients, and the eigenfunctions psi_j(t) = phi(t)' v_j over the covariate grid. These are rotation-invariant functionals of K_g (mirroring the engine's rr_eigenfunctions); no raw, rotation-arbitrary loadings are returned.

The trajectories are computed in R from the estimated K_g and the recorded covariate standardization range; at is supplied on the ORIGINAL covariate scale (defaulting to a grid over the fitted range) and re-standardized to [-1, 1] internally, matching the Julia engine's basis convention.

Response scale

Rescale the response to about unit variance before fitting. On an unscaled trait the dense Nelder-Mead search can report converged = TRUE below the REML optimum and return a wrong K_g (off-diagonals stuck near 0). Dividing y by sd(y) recovers the optimum; multiply fitted variances by sd(y)^2 to return to the original scale. A converged = TRUE flag is not enough on the raw scale.

Out-of-range `at`

at must lie inside the fitted covariate range (object$result$random_regression$lower to $upper); a value outside that range errors with a message naming the fitted range, instead of silently returning the nearest endpoint's value. This matches the Julia engine (HSquared.legendre_basis), which refuses a standardized covariate outside [-1, 1] rather than clamping it (#213).

Examples

if (FALSE) {
fit_rr <- hsquared(
  weight ~ sex + animal(rr(age, order = 2) | id, pedigree = ped),
  data = long_records,
  family = gaussian(),
  REML = TRUE,
  control = hs_control(
    engine = "julia",
    engine_control = list(target = "random_regression")
  )
)

rr_covariance(fit_rr)
random_coefficients(fit_rr)
rr_genetic_variance(fit_rr)
rr_heritability(fit_rr)
rr_correlation(fit_rr, at = c(1, 3, 5))
rr_eigenfunctions(fit_rr)
}