hsquared v0.9.0

hsquared v0.9.0 is an experimental, evidence-bounded R release. It preserves the documented 0.9.0 scientific contract and does not promote any capability beyond that boundary.\n\nThis is a GitHub release only and has not been submitted to CRAN. Validation evidence and known limitations remain documented in the project status pages.

hsquared 0.1.0

hsquared 0.1.0

New features

• The default hsquared() call now fits the v0.1 univariate Gaussian animal model y ~ fixed + animal(1 | id, pedigree = ped) by REML (average-information) through the HSquared.jl engine, returning heritability, variance components, breeding values, fixed effects, fitted values, residuals, and diagnostics. Fitting requires a local Julia, the JuliaCall package, and an HSquared.jl checkout; without them the default call errors with install guidance, and hs_control(engine = "validate") validates the contract without fitting. ML is not implemented — REML = FALSE is rejected on the fit path. The fit is validated by known-truth recovery, the published gryphon REML anchor (within the signed-off comparator band), and sommer agreement (#6, #7). These engine-recovery checks run locally through the R-to-Julia bridge (a local Julia and HSquared.jl checkout are required); public CI exercises the equivalent pure-R REML reference and skip-guards the live-engine tests.
• animal() is now exported as an inert formula marker, and hsquared() parses the narrow v0.1 formula contract animal(1 | id, pedigree = ped) (#4, #6).
• EBV() and BLUP() now alias breeding_values() for hsquared_fit objects, and accuracy() derives square-root reliability when reliability estimates are present (#5).
• fitted() and residuals() now work for hsquared_fit objects that contain fitted-value predictions and response values (#5).
• Added fit_diagnostics() for hsquared_fit objects so users and developers can inspect engine, method, target, convergence, optimizer status, iterations, log-likelihood, and dense-validation-path metadata without refitting or implying production support (#5, #6).
• coef() now aliases fixed-effect extraction for hsquared_fit objects, and nobs() reports the result observation count or response-payload fallback when available (#5).
• hs_control() has an engine = "julia" option that exposes advanced engine targets (supplied-variance Henderson MME, the opt-in sparse REML optimizer, explicit ai_reml control, repeatability/two-effect/genomic/SNP-BLUP targets, and the opt-in multivariate target) through a sibling HSquared.jl checkout via JuliaCall. The default engine = "fit" already fits the v0.1 model via ai_reml; non-Gaussian fitted-model support is still planned (#6).
• hs_control() now recognizes an experimental, opt-in engine_control = list(target = "sparse_reml", initial = ..., iterations = ...) path that surfaces the Julia-owned HSquared.fit_sparse_reml() REML-only sparse optimizer through the local bridge. It is not the default hsquared() behaviour, not variance-component estimation in the public R interface, and not a production sparse fitting, AI-REML, or ASReml-parity claim (#6, #7).
• Experimental, opt-in repeatability (permanent-environment) model. hsquared() now parses animal(1 | id, pedigree = ped) + permanent(1 | id) and fits it through engine_control = list(target = "repeatability"), surfacing the Julia-owned HSquared.fit_repeatability_reml() REML-only optimizer. It returns three variance components (animal, permanent, residual), repeatability(), heritability(), breeding values, and permanent_effects(). The default engine = "fit" path stays single-effect and rejects permanent() with a pointer to the opt-in target. This is experimental and REML-only; the additive (σ²a) and permanent-environment (σ²pe) variances are identifiable only with repeated records per individual. It is not the default, not ML, and not yet a production or comparator-validated claim.
• Experimental, opt-in common-environment (two-effect) model. hsquared() now parses animal(1 | id, pedigree = ped) + common_env(1 | group) and fits it through engine_control = list(target = "two_effect"), surfacing the Julia-owned HSquared.fit_two_effect_reml() REML-only optimizer (additive genetic effect + an IID common-environment effect, e.g. litter or cage). It returns three variance components (animal, common_env, residual), heritability(), breeding values, and common_env_effects(). The default engine = "fit" path stays single-effect. This is experimental and REML-only; not the default, not production or comparator-validated; correlated direct–maternal (2×2 G) effects remain planned.
• Experimental, opt-in maternal-genetic two-effect model. hsquared() now parses animal(1 | id, pedigree = ped) + maternal_genetic(1 | dam) and fits it through engine_control = list(target = "two_effect"): a direct additive genetic effect plus a maternal genetic effect expressed through the dam, both carrying the pedigree relationship (A₂ = pedigree A). It returns three variance components (animal, maternal_genetic, residual), heritability(), breeding values, and maternal_effects(); the dams must be animals in the animal() pedigree. Experimental and REML-only; the correlated direct–maternal (2×2 G) model remains planned.
• Experimental, opt-in genomic GREML model. hsquared() now parses genomic(1 | id, Ginv = Ginv) — a primary genomic effect with a user-supplied genomic relationship inverse — and fits it through engine_control = list(target = "genomic"), surfacing HSquared.fit_ai_reml() on a Ginv-based spec (genomic REML). It returns the genomic and residual variance components, genomic heritability(), and genomic breeding values (GEBVs). You can either supply a precomputed Ginv, or pass a raw marker matrix (genomic(1 | id, markers = M)) and let the engine build the genomic relationship and its inverse. The default engine = "fit" path is the pedigree animal model; genomic() is opt-in. Experimental and REML-only; comparator validation and weighted/standardized-marker variants remain planned.
• Experimental, opt-in single-step model. hsquared() now parses single_step(1 | id, Hinv = Hinv) — a primary effect with a user-supplied single-step relationship inverse — and fits it through engine_control = list(target = "single_step"), reusing the same fit_ai_reml-on-a-supplied-inverse path as genomic(). Experimental and REML-only; building Hinv from a pedigree and genomic relationship (single-step HBLUP construction) and comparator validation remain planned.
• Experimental, opt-in SNP-BLUP / RR-BLUP marker-effect model. hsquared() now fits genomic(1 | id, markers = M) through engine_control = list(target = "snp_blup", variance_components = c(sigma_g2 = ..., sigma_e2 = ...)), surfacing the Julia-owned HSquared.fit_snp_blup(). At the supplied genomic/residual variances it estimates per-marker effects — extractable with marker_effects() — together with per-individual genomic breeding values and fixed effects (the engine centres the markers, VanRaden method 1). marker_variance_explained() reports a descriptive fitted-marker contribution table for this path, computed from effect squared times centered marker variance and normalized across markers; it is not a marker-scan, p-value, or QTL claim. This is a supplied-variance solve (it does not estimate the variance components), opt-in, and not the default; it mirrors the twin V2-SNPBLUP gate (the GBLUP↔SNP-BLUP genomic-breeding-value equivalence). REML estimation of the marker variance, weighted/Bayesian marker priors, and comparator parity remain planned.
• Experimental, opt-in multivariate Gaussian animal model. hsquared() now parses cbind(trait1, trait2, ...) ~ fixed + animal(1 | id, pedigree = ped) and fits it through engine_control = list(target = "multivariate"), surfacing the Julia-owned HSquared.fit_multivariate_reml() REML-only dense estimator. It returns G/R covariance matrices, genetic and residual correlations, per-trait heritability, and cross-trait EBVs through genetic_covariance(), residual_covariance(), genetic_correlation(), residual_correlation(), heritability(), and breeding_values(). Missing trait cells are accepted as NA and marshal to Julia NaN; rank-deficient fixed-effect designs are rejected up front; non-converged fits do not expose logLik()/AIC(). This is opt-in, dense validation-scale, and partial; t>=2 known-truth recovery, external comparator parity, and long-format/structured covariance grammar remain planned.
• hsquared() now fences the reserved multivariate engine_control$genetic_structure and engine_control$rank fields: "unstructured" is accepted for the current opt-in multivariate bridge, while "diagonal", "lowrank", "factor_analytic", and rank error as planned rather than being silently ignored.
• formula_status() and structured-covariance error messages now list the full planned covariance vocabulary (cov = us(), cov = diag(), cov = lowrank(K), and cov = fa(K)) while keeping it planned, not fitted.
• Multivariate cbind() responses now require unique, non-empty trait names before fitting, so G/R matrices, per-trait h², EBVs, comparator files, and future wide/long response paths share an unambiguous trait-order contract (#10).
• G_matrix() and R_matrix() now alias genetic_covariance() and residual_covariance() for hsquared_fit objects, giving multivariate users the familiar G/R matrix names without changing the underlying extractor contract.
• Added a "Reading G matrices" pkgdown article that explains the current G/R matrix extractors, genetic and residual correlations, per-trait h2, cross-trait EBVs, and the boundaries around P_matrix(), factor-analytic loadings, and selection-response claims.
• Added a "Genomic prediction" pkgdown article that separates the current opt-in supplied-Ginv, marker-built GREML, SNP-BLUP, and supplied-Hinv single-step paths from planned genomic construction, APY, GWAS/QTL/eQTL, and production-comparator work.
• Added a "QTL, GWAS, and eQTL status" pkgdown article that explains the current reserved scan vocabulary, live SNP-BLUP marker effects / descriptive marker variance shares, scale caveats, and ...