Standalone YAML model specification

The experimental user interface combines short model formulas with named random effects, reusable relationships, bases and covariance matrices. Trait and basis channels are generated automatically. Users describe a model once; the compiler generates its explicit numerical representation and makes it available for inspection before fitting. The terminology is intended for plant, human and animal data: an individual may be a person, organism, line, cultivar or other unit identified by the relationship model. Effect names are chosen by the user.

gsuite.model/v1-draft is implemented in the standalone gsuite run model.yml frontend with an updated greml example driver. The examples include tiny synthetic CSV and PLINK inputs. Supplied-covariance solving and deterministic sparse-direct REML for unstructured/FA covariance and serial parameters are supported. Packed EM-REML supports unstructured genetic and independent residual covariance with explicit exact/stochastic trace controls. Selected ordinary packed genetic predictions can request exact or explicitly budgeted stochastic uncertainty conditional on supplied covariances. This schema is the sole model interface. This draft is an experimental development interface, without an industrial-scale parser qualification.

Complete examples

Example filenames below refer to models/ in the supplied standalone distribution. The YAML files and their synthetic inputs are not downloads from this website.

File	Model	Current implementation boundary
Multitrait pedigree (models/multitrait-pedigree.yml)	Two illustrative human traits, fixed cohort and age, correlated individual genetic effects	Compiled-MME supplied-covariance solve; small dense-equation comparison.
Multitrait single step (models/multitrait-single-step.yml)	Five traits, pedigree/genomic individual effects, site, family and individual environment	Packed conditional PCG; shifted-pedigree and Qnn-reuse preconditioners; small dense-H comparison.
Long-format single step (models/multitrait-single-step-long.yml)	The same five-trait model from ID/trait/value phenotype rows	Grouped residual occasions retain trait covariance and missing traits; shared native solve.
Direct/maternal single step (models/direct-maternal-single-step.yml)	Correlated direct and maternal genetic channels, separate individual/mother incidence and maternal environment	One shared H and joint channel covariance; packed all-record prediction with separate channel outputs.
Selected genetic uncertainty (models/multitrait-single-step-uncertainty.yml)	Two traits and explicit individual/trait targets	Exact conditional PCG PEV and reliability; retained factors, bounded work and dense-inverse comparison.
Stochastic genetic uncertainty (models/multitrait-single-step-stochastic-uncertainty.yml)	All genotyped individuals and explicit probe settings	Variance-reduced PCG diagonal estimation, MCSE and per-target qualification.
Multitrait random regression (models/multitrait-random-regression.yml)	Repeated-measure curves with genetic and individual-environment coefficients and site-visit effects	Joint compiled-MME solve; six-channel curve effects, correlated coefficients and prediction grid checked against dense equations.
External covariance and grouped residuals (models/multitrait-external-covariance.yml)	Supplied individual relationship, group-specific trait residual covariance and fixed interaction	Natural, labelled covariance CSV; compiled direct or PCG solve.
Nested regression (models/multitrait-nested-regression.yml)	Legendre genetic curves, individual-by-visit environment and period-specific fixed curves	Explicit interaction expansion, composite level IDs and curve prediction.
Discrete AR1 (models/longitudinal-ar1.yml) / continuous exponential (models/longitudinal-exponential.yml)	Univariate repeated measurements with missing visits and actual position gaps	Native structured residuals; compiled direct or PCG solve.
Packed AR1 (models/single-step-ar1.yml) / packed exponential (models/single-step-exponential.yml)	PLINK/pedigree ssSNPBLUP with supplied serial residual parameters	Exact sparse Markov innovations, all-record conditional PCG and actual observed time gaps.
Exact ssGBLUP (models/multitrait-ssgblup.yml)	Supplied genomic covariance, pedigree and non-genotyped individuals	Implicit H precision through compiled PCG; no full H matrix.
Raw-genotype APY ssGBLUP (models/multitrait-ssgblup-apy-plink.yml)	PLINK panel, ordered core-ID CSV and explicit genomic scaling/ridge	Tiled native preparation, file-backed cross products and implicit H PCG; small exact references, industrial scaling unqualified.
APY ssGBLUP (models/multitrait-ssgblup-apy.yml)	Explicit core covariance, noncore/core covariance and marginal diagonal	Native APY genomic preparation followed by implicit H PCG; approximation labelled in diagnostics.
Metafounders (models/multitrait-metafounders.yml) / metafounder ssGBLUP (models/multitrait-ssgblup-metafounders.yml)	Supplied base Gamma and explicit unknown-parent assignments	Compiled pedigree or supplied-genomic H; latent base levels retained.
Genetic groups (models/multitrait-genetic-groups.yml) / group-mean APY ssGBLUP (models/multitrait-ssgblup-genetic-groups.yml)	Independent random group means and inherited contributions	Compiled a + F gamma model; coefficients and combined individual predictions.
Packed metafounders (models/multitrait-packed-metafounders.yml) / packed genetic groups (models/multitrait-packed-genetic-groups.yml)	PLINK marker products with supplied base covariance or independent group means	Full or ancestral Qnn; metafounders retain their correlated base block, group inheritance uses implicit sweeps.

From the checkout, with an updated native executable:

python tools/build/hpc_single_step/gsuite.py describe \
  tools/examples/model_yaml/multitrait-single-step.yml
python tools/build/hpc_single_step/gsuite.py run \
  tools/examples/model_yaml/multitrait-single-step.yml \
  --exe /path/to/gsuite_model_driver --output /project/results/small-single-step

The final output directory must be fresh; managed scratch is collected there, including prepared inputs/, run diagnostics and resource reports. run uses output.directory when --output is omitted. Relative input and output paths in YAML resolve beside that YAML file. validate MODEL checks inputs; describe MODEL prints expanded model settings without fitting. Compact models take execution settings only from YAML; solver/thread command-line overrides are not exposed. Use a matching frontend and rebuilt native driver. The standalone source distribution and complete workflow are described in the GenomeDK bundle guide. For installation and a complete workflow start with the standalone user guide.

The canonical model architecture owns the compiled representation. The model catalog owns capability and qualification claims. A new input keyword does not establish a supported model or a production-scale performance claim.

Compact model declarations

The random-regression example’s model section is:

model:
  data: observations
  fixed: "factor(site_period) + curve"
  random:
    genetic:
      group: individual
      relationship: A
      basis: curve
      covariance: genetic_curve
    individual_environment:
      group: individual
      basis: curve
      covariance: individual_curve
    site_visit: {group: site_visit, covariance: site_visit}
  residual: residual

Effect names are stable model identities, not reserved biological terms. genetic and individual_environment use the same individual-ID binding but remain independent components. The first uses relationship A; the second defaults to identity. Each applies to all declared traits. Define a basis once:

bases:
  curve: {column: time, degree: 2, center: 155, scale: 150}

With traits [trait1, trait2], each curve effect expands to six channels: trait1.b0, trait1.b1, trait1.b2, trait2.b0, trait2.b1, trait2.b2. Site-visit has two intercept channels. The user does not write these channel lists or an incidence matrix. Output request names such as genetic_values are also user defined; a breeding application may choose breeding_values.

Specifying all input files

All input files are declared under inputs in the same YAML file as the model. Logical names on the left of columns and responses are the names used by formulas, grouping, bases and trait ordering. Values on the right are actual source headers. The pedigree example contains:

traits: [trait1, trait2]
inputs:
  observations:
    type: table
    file: data/phenotypes.csv
    format: csv
    id: record_id
    columns: {individual: sample_id, cohort: cohort, age: age_years}
    responses: {trait1: height_cm, trait2: bmi}
  pedigree:
    type: pedigree
    file: data/pedigree.csv
    columns: {individual: sample_id, parent1: parent_1, parent2: parent_2}
relationships:
  A: {type: pedigree, input: pedigree}
model:
  data: observations
  # fixed, random and residual declarations follow here.

For example, the phenotype header is record_id,sample_id,cohort,age_years,height_cm,bmi and the pedigree header is sample_id,parent_1,parent_2. A record ID identifies a measurement, visit or plot. An individual ID identifies the unit carrying genetic or other individual effects; it can appear on several records. Equality of individual IDs provides the declared record-to-relationship link, not an automatic table join.

parent1 and parent2 are neutral parent-role bindings. Actual columns can be named father, mother, sire, dam or otherwise; users map those headers explicitly. These labels do not select a different pedigree algorithm, infer sex or request probabilistic parentage.

The single-step example additionally declares the genotype panel in inputs:

inputs:
  genotypes:
    type: plink
    prefix: data/genotypes
    id: IID

This names data/genotypes.bed, .bim and .fam. FAM IID matches pedigree individual IDs exactly. The H relationship refers to genotypes: genotypes, where the second name is the input identifier. No separate data configuration or command-line file list is required. All paths resolve against this YAML file’s directory.

The curve example likewise declares prediction_grid as a table input, with columns: {individual: sample_id, time: time_days}. Its prediction request uses data: prediction_grid. Covariance matrices may stay inline; for external matrices declare a type: covariance input and reference it from covariances. Thus prediction and matrix input paths can also remain in the one input block.

Long-format phenotypes

Choose either layout: wide (the default) or layout: long on the observation table referenced by model.data. For long input, replace the wide trait-to-column response mapping with two source-column bindings:

traits: [trait1, trait2]
inputs:
  observations:
    type: table
    file: data/phenotypes_long.csv
    format: csv                    # tsv is also supported
    layout: long
    id: visit_id                   # Shared residual occasion, repeated by trait
    columns: {individual: sample_id, time: day, site: site_id}
    responses: {trait: measured_trait, value: measurement}
    missing_values: ['', NA, '.']
# Relationships, model terms, covariances, execution and outputs are the same
# as for wide input. Covariance rows/columns follow traits: [trait1, trait2].

For example:

visit_id,sample_id,day,site_id,measured_trait,measurement
v001,003,0,S1,trait2,3.2
v002,003,2,S1,trait1,1.8
v001,003,0,S1,trait1,1.4
v002,003,2,S1,trait2,NA

responses.trait and responses.value name source headers; values in the trait column must exactly match declared trait names. The ID identifies a residual occasion, such as a visit or plot, globally across the file. It repeats across traits; an ID/trait pair must occur at most once, even if its value is missing. Use a different ID for each visit of an individual. Rows can be interleaved, and traits need not appear in declaration order. The frontend retains ID order of first appearance and orders scalars within an ID by declared trait order.

All logical columns bindings must have identical source strings within an ID, including individual, time, effect groups and fixed covariates. Conflicting values are rejected rather than inferred or averaged. Trait/value source columns must be distinct from the ID and mapped columns. Without columns, other source headers become logical columns with the same names. As with wide input, explicit bindings are recommended.

Absent ID/trait pairs and configured missing-value tokens both mean an unobserved scalar. Unknown/empty trait labels, duplicate pairs, blank IDs and nonfinite or nonnumeric observed responses are rejected. An entirely missing occasion is omitted from fitting and counted in omitted_records; an occasion absent from the file cannot be counted. A prediction request with data: observations emits each ID once for every requested output trait, including explicitly supplied all-missing occasions. Other helper inputs and separate prediction grids remain wide tables with unique IDs where declared.

The two observed traits at v001 use the joint residual covariance submatrix; they are not treated as independent records. Distinct visit IDs are independent under a simple residual covariance, while declared AR1/exponential group and position columns connect visits for a serial model. Fixed formulas, random-effect channels, relationships, natural covariance matrices and solver settings do not change between layouts. This supports different formulas by trait with shared mapped column values; within-occasion trait-specific column values are not an additional incidence syntax.

describe and model.json report phenotype_layout, physical phenotype_rows, retained residual records, omitted_records and observed_scalars separately. Long input is grouped in Python memory, with no dense population-sized pivot or incidence matrix. It does not add streaming input or an application memory cap. The complete five-trait example (models/multitrait-single-step-long.yml) includes the matching CSV, PLINK inputs, covariances, thread settings and outputs.

Focused Windows checks compare equivalent wide/long models with independently assembled dense equations for multivariate pedigree direct/PCG solves, conditional packed ssSNPBLUP, quadratic random regression and packed AR1/exponential residuals. The largest coefficient error is 1.10e-10 (compiled PCG); the other cases are below 6.3e-14. Native acceptance retains both requested residual gates. Tests include interleaved rows, reversed trait order, absent/explicit missing traits, all-missing occasions, CSV/TSV and per-occasion predictions. Duplicate pairs and inconsistent mapped columns are rejected. All 60 focused YAML tests pass, including twenty-two runnable examples and equivalent coefficients for the complete five-trait wide/long pair. These are small correctness checks, not industrial parser-memory or throughput qualification. Validation uses NumPy 2.3.5 and PyYAML 6.0.3; native model/solver and PEV precision work are unchanged.

Sections and defaults

Section	Contents and defaults
schema, name	Required version literal and descriptive model identifier
traits	Required ordered trait names; determines automatic covariance ordering
inputs	Named table, pedigree, plink or covariance inputs with file/prefix and explicit column/ID bindings; CSV is the table default
relationships	Named relationship definitions referring to inputs or other relationships
bases	Optional reusable power, Legendre or Chebyshev polynomials; includes constant, generated columns b0 through bd
covariances	Named covariance matrices inline, labelled inline, or by {input: name}
model.data	Required reference to the table input supplying observations and response bindings
model.fixed	Common fixed formula, or a complete map of trait names to fixed formulas; implicit intercept
model.random	Named effect mappings; group and covariance required; relationship defaults to identity; no basis means an intercept; traits default to all
model.residual	Covariance name for independent observation units, or a typed grouped/AR1/exponential declaration
analysis	Explicit task: solve or supported reml; packed EM adds explicit traces and optional likelihood controls
execution	Solver settings separate from the model; applicable thread counts default to 1
output	Fresh result directory and named prediction requests; default tables are fixed effects, random effects and diagnostics

Defaults are recorded in the expanded model. All defaults are versioned with the schema. Relationship choice, covariance values, polynomial center/scale and genomic-base adjustment must not be inferred from effect names or fitted automatically.

Object names match [A-Za-z][A-Za-z0-9_]*; generated channel labels may contain periods. References are case sensitive and identity is reserved. Unknown or duplicate keys, unresolved/cyclic references, nonfinite values and unsupported combinations are errors. Use one plain YAML document without executable tags, merge keys, anchors or environment substitution. Multiple phenotype tables are not implicitly merged; model.data selects one, already aligned table for a fit.

Fixed formulas and regression bases

The formula grammar supports +, :, *, 1, 0, numeric column names, factor(column) and named bases. A numeric column gives a slope; a factor uses treatment coding with the first observed level as reference. Levels are sorted by UTF-8 byte order independently of host locale. No column type is guessed from numeric-looking IDs. Basis/column name collisions are errors.

An intercept is included unless 0 is present. A full basis containing the constant b0 and an implicit intercept generates that same constant once, reported as (Intercept); this structural deduplication does not remove other linearly dependent columns. With 0, a requested basis retains its constant as curve.b0. With 0, factors use all observed level indicators. Combining that full indicator set with a basis constant is rank deficient and is rejected. Other rank deficiencies are validation errors.

A formula string applies separately to every trait. A trait-to-formula map must name every declared trait and allows different fixed effects; neither spelling constrains fixed coefficients equal across traits. a:b multiplies expanded columns; a*b expands to a + b + a:b, with repeated terms deduplicated. factor(site)*curve fits site-specific fixed curves using a reference curve and differences. factor(site):curve alone contains only the interaction contrasts under an intercept, so use * when all site curves are intended. Products may have more than two terms. Parenthesized subformulas, subtraction, arbitrary functions and general R/Python expression evaluation are unsupported.

A polynomial basis requires a numeric column, nonnegative degree d, finite center and positive scale. It generates 1,x,...,x^d for x=(value-center)/scale by default (family: power). family: legendre uses the ordinary Legendre polynomials, beginning 1,x,(3*x*x-1)/2; family: chebyshev uses first-kind polynomials, beginning 1,x,2*x*x-1. There is no sample-dependent orthogonalization or additional normalization. Center/scale are explicit, including for orthogonal-polynomial families; values outside [-1,1] are evaluated as extrapolations. The same family, center, scale and channel order apply to fitting and prediction.

Random effects and automatic expansion

model.random is a map from effect names to declarations. group names the logical column containing level IDs; it is always explicit, even when it matches the effect name. For identity relationships, group: 'individual:visit' builds the joint domain and fits a separate coefficient vector for every observed individual/visit combination. With basis: curve this is a nested random regression. Composite IDs are JSON arrays of the constituent string values; delimiter characters in source IDs cannot merge distinct combinations. Prediction tables bind each constituent column. A named relationship requires a single grouping column and retains its complete declared level domain, including unobserved pedigree individuals. Identity uses the sorted domain of levels on rows with at least one observed response.

Without basis, an effect has one channel per selected trait. With a basis, it has the trait-by-basis product, in root trait order then basis-column order. Optional traits: [trait1] restricts an effect to a subset; selection does not reorder the root trait list. Coefficients are level-major then channel-major. For effect j, the covariance is R_j tensor K_j, where R is its relationship and K its supplied channel covariance. All effects are independent of each other and residuals. Sharing a covariance reference ties covariance parameters, not random-effect realizations or cross-effect covariance.

An advanced channels list may replace automatic traits/basis expansion. Each entry explicitly names a channel, trait and optional basis column, for example {name: trait1_slope, trait: trait1, basis: curve.b1}. It is mutually exclusive with effect-level traits and basis. This preserves explicit trait-specific basis layouts without requiring them for routine models. The effect-level group supplies the default incidence column. Each explicit channel may override it with its own group. All channels retain the same relationship level domain and joint covariance. For example:

covariances:
  genetic:
    - [1.0, -0.18]
    - [-0.18, 0.5]
model:
  random:
    genetic:
      group: individual
      relationship: H
      channels:
        - {name: direct, trait: weight}
        - {name: maternal, trait: weight, group: mother}
      covariance: genetic

This declares direct[individual] + maternal[mother] in each weight record, with covariance H tensor genetic in channel order [direct, maternal]. It is one correlated effect, rather than two independent random declarations. The names direct, maternal, individual and mother are user-chosen. Each root trait can have several channels; multitrait roles use an explicit joint channel covariance. Channel basis weights compose with the incidence override. describe reports every channel’s trait, grouping and basis binding.

The observation input must bind mother to an explicit source column. No parent lookup, file join, unknown-parent substitution or missing-incidence inference is performed. Referenced IDs must exist in the relationship; blank or unknown IDs are errors when that channel’s response is observed. An unknown-parent token such as 0 is not silently substituted or omitted. Metafounder levels remain latent. Named relationships require one column per channel; composite groups remain available for identity effects. Independent random genetic-group means still require their shared incidence mapping. Multiple contributors within one channel and per-level covariance transforms remain outside this interface.

Packed models with distinct genetic incidence bindings require the native capability packed_multi_individual_incidence. They use explicit all-record, fixed-covariance prediction preparation with the same bounded record cache. They do not supply owner-filtered training metadata or perform marker-weight learning. Existing owner-filtered training checks remain strict. Both Qnn factor modes and both packed PCG preconditioners remain selectable.

Focused Windows validation covers eight independent dense-equation references: an ordinary pedigree model, all four packed factor/preconditioner combinations, a pedigree-blended genomic covariance, and compiled/packed basis-weighted incidence. The two-trait fixture has negative direct/maternal covariance, inbreeding, reordered pedigree/genotype IDs, mixed genotyping, a missing call, missing responses, repeated observations and maternal environment. Maximum primary-coefficient difference is below 4.6e-13; selected-channel and mapped record predictions also agree. Both residual gates retain their requested thresholds. The driver solves internally ten times more tightly for this all-record path, matching the existing group/metafounder policy.

All 53 focused YAML tests pass, including nineteen runnable examples and invalid-ID/channel/old-executable rejection. Four native marker-training regressions pass, covering the new cached/streamed all-record actions, cancellation/provider failures and unchanged ownership/held-out-response checks. Receipts are build/validation/yaml-correlated-incidence-regression.log and build/validation/yaml-correlated-incidence-native-tests.log. These are small exact-reference and regression checks, not industrial/HPC qualification, new uncertainty qualification or a full R-package check. Those observations precede the independent-pilot precision extension below.

Covariance input and estimation

The examples keep their complete matrices inline and omit repeated labels:

covariances:
  genetic:
    - [100, 12]
    - [12, 4]
  residual:
    - [150, 8]
    - [8, 6]

An unlabelled matrix follows the generated channel order; residuals follow root trait order. Every dimension is checked. An external matrix is declared in the same input block as other files:

inputs:
  genetic_matrix: {type: covariance, file: data/genetic_covariance.csv}
covariances:
  genetic: {input: genetic_matrix}

The matrix CSV has a first header coordinate, named numeric columns, and matching row labels in that first column. Complete row and column labels must match the generated coordinates; no automatic matrix reordering occurs. A mismatch reports the required order. An inline {labels: [...], values: [...]} mapping provides the same check. type: covariance always means a covariance matrix, never a precision matrix.

Matrix inputs are variances and covariances, never inverses or packed Cholesky parameters. Matrices must be finite, symmetric and positive definite for the initial adapter. Referencing a matrix from several effects requires matching channel coordinates. Semidefinite boundaries and zero-variance components need an explicit later contract.

Factor-analytic genetic covariance

A random-effect covariance can instead supply natural factor loadings and specific variances. Greml uses K = loadings * loadings' + diag(specific_variances), with independent unit-variance factors. For example:

covariances:
  genetic:
    type: factor_analytic
    rank: 2
    labels: [north, south, east, west, central]
    factors: [shared1, shared2]
    loadings:
      - [0.9, 0.0]
      - [0.5, 0.8]
      - [-0.4, 0.3]
      - [0.7, -0.2]
      - [0.2, 0.6]
    specific_variances: [0.6, 0.4, 0.7, 0.5, 0.8]

Rows follow the effect’s generated channel order; columns follow factor order. Optional labels must match that channel order exactly. Optional factors defaults to factor1, factor2, etc. Loadings have response units, and specific_variances are variances, not standard deviations or inverse values. The upper triangle of the leading rank-by-rank loading block is exactly zero, its diagonal pivots are positive, and all specific variances are strictly positive. Free remaining loadings may have either sign. Rank must be positive and smaller than the number of channels. The identified parameter count q*r - r*(r-1)/2 + q cannot exceed q*(q+1)/2. Greml additionally checks the actual covariance Jacobian rank during native preparation/evaluation; a rank-deficient state is rejected without rotation, clipping or fallback. describe checks shapes, natural values and structural constraints; it does not run that native identification check.

The three-trait example (models/multitrait-factor-analytic.yml) uses FA1 in a pedigree model. The five-environment example (models/multienvironment-factor-analytic.yml) uses FA2 with packed ssSNPBLUP. It represents one response measured across environments through the existing multivariate channels:

traits: [north, south, east, west, central]
inputs:
  observations:
    type: table
    file: data/environment_yield.csv
    layout: long
    id: record_id
    columns: {individual: sample_id, cohort: cohort, replicate: replicate}
    responses: {trait: environment, value: yield}

The source environment column selects a declared channel; yield contains the measured value. Fixed formulas are evaluated separately for each channel, so their intercepts represent environment-specific means. Genetic values of the same individual in different environments share the FA covariance. A diagonal residual covariance makes residuals from different environments independent. Reusing an ID groups paired measurements with matching mapped columns; separate plot IDs can instead carry one observed environment each, including plot-specific covariates. Joint trait/environment outcomes can be listed explicitly as distinct channel names; the adapter does not construct an implicit Cartesian product or a separable trait-by-environment model. The output trait column contains the declared environment/channel labels.

Compiled models use greml’s public FA random-group and natural parameter APIs. Packed prediction validates the same native FA layout and implied covariance, then reuses its existing dense channel precision and PCG route. Native capability factor_analytic_covariance and compact format 8 are required. run.json retains the specification, and covariance_models.json reports labels, factors, loadings, specific variances and the implied covariance. Predictions remain genetic values per individual and trait/environment. The genetic coefficient count is unchanged; this parameterization alone does not establish a faster or smaller fixed-covariance solve. Existing native FA validation includes derivative/Jacobian work and needs separate scaling qualification for many channels.

The interface accepts supplied FA parameters for prediction or starting FA parameters for sparse-direct REML. FA is accepted for random-effect covariances, not residual components. Uncertainty requests for FA YAML models are rejected pending separate validation. Small independent dense references cover FA1/FA2, missing responses, environment-labelled long input, partial genotyping and both packed pedigree factors/preconditioners. These checks do not qualify industrial scale.

analysis: {task: solve} fixes all supplied covariance values. A future task: reml additionally requires estimate: [genetic, residual]; the named matrices use supplied values as starts while other matrices stay fixed. Requests must match the native model/backend/estimation capability. A zero off-diagonal entry is simply a supplied zero or a starting zero, not an implicit constraint that fixes that covariance during estimation. Partial parameter constraints require a separate syntax.

Input binding and relationship semantics

A type: table input has file, optional format: csv or tsv, and column bindings. Observation/prediction tables require an explicit unique id source column. IDs are strings with leading zeros preserved. A columns mapping binds logical model names to source headers; if omitted, headers themselves are the logical names. When supplied, only the listed logical names, declared responses and record ID are available to the model. Referenced but unmapped fields, conflicting logical names and absent headers are errors. Unused source columns may remain in the input table. No files are joined or renamed on disk.

An observation input additionally needs responses, binding each ordered trait to a distinct numeric source column. Wide layout is the initial default; missing response tokens default to "", "NA", ".". Required grouping, factor, basis and covariate values must exist on observed rows. All-missing response rows are omitted and counted. Partially missing traits use the principal residual covariance submatrix, not imputation. Other missing-value tokens can be declared on the input; they apply to responses only.

A wide row or one grouped long record ID is one residual unit by default. model.residual: residual means independent occasions with the named correlated trait covariance. Repeated individual records share their declared random effects.

Grouped covariance uses complete trait covariance matrices, with exact observed group names as keys (quote numeric-looking keys):

model:
  residual:
    type: grouped
    group: site
    covariances: {S1: residual_S1, S2: residual_S2}

Both covariance names are declared in covariances, in root trait order. Rows remain independent conditional on random effects. Missing traits select the principal covariance submatrix for that row’s group, then invert it; they never select a submatrix of the full precision. The group map must match the observed groups exactly. This is supported by compiled MME and ssSNPBLUP.

Univariate repeated measurements may instead use:

model:
  residual:
    type: ar1
    group: individual
    position: time
    variance: 1.7
    correlation: 0.4

Within each group, Cov(e_i,e_j) = variance * correlation^abs(time_i-time_j). Positions must be distinct integers within +/- (2^53-1) within a group and correlation lies in (-1,1). Input order need not be temporal. Missing responses do not renumber visits. For continuous positions, use type: exponential and rate: 0.7 instead of correlation; covariance is variance * exp(-rate*distance) with positive variance and rate. Distinct groups are independent. These serial structures use one supplied variance and correlation/rate and one trait. They support both compiled MME and packed ssSNPBLUP. Within-series changes in the marginal trait covariance remain outside this adapter. The compiled structured residual implementation constructs local dense work per series and has separate long-series memory/performance limits.

Packed models require the native capability packed_serial_residuals. Greml prepares exact stationary Markov innovations once, transforming both the response and the full sparse incidence row. For the first observed record in each group the residual variance is the supplied variance. Subsequent rows subtract phi times the preceding observed row, with innovation variance variance * (1-phi^2). Here phi = correlation^gap for AR1, including its sign at odd/even integer gaps, or exp(-rate*gap) for exponential covariance. Actual observed time gaps are retained after missing responses are omitted. Singleton groups use the supplied marginal variance. No visits are imputed, renumbered, cut into independent chunks or dropped to meet a block-size cap.

The native adapter sorts only group/position metadata, retains original record order, and produces scalar innovation blocks with coalesced sparse incidence. Metadata uses O(records) storage and preparation reads at most two source rows per innovation; no dense per-series covariance or inverse is built. Prepared records then reuse the existing bounded record cache and PCG operator. This selects all-record fixed-covariance prediction, without owner-filtered training or marker-weight updates. Both Qnn factor modes and preconditioners are available. No application RAM ceiling or timeout is introduced; record entry limits and caller-specified component allowances retain their scope. Covariance fitting and new PEV precision qualification are not added.

For multivariate serial residuals, replace variance with a named natural trait covariance, in root trait order:

covariances:
  residual:
    - [2.0, 0.3]
    - [0.3, 1.5]
model:
  residual:
    type: ar1
    group: individual
    position: time
    covariance: residual
    correlation: -0.4

Within a group, Cov(e_it, e_js) = residual[t,s] * correlation^abs(i-j); AR1 positions must be exactly represented integers. For continuous time, use type: exponential and positive rate, giving the multiplier exp(-rate * abs(i-j)). Supply exactly one of variance or covariance; multiple traits require covariance. Distinct groups remain independent. The packed AR1 example (models/multitrait-single-step-ar1.yml) and compiled exponential example (models/multitrait-longitudinal-exponential.yml) include input files, matrices, execution settings and predictions.

Both prediction routes require the new native capability multivariate_serial_residuals. Greml prepares exact q-state conditional innovations, updating the full trait residual state at the actual observed gaps. At each visit it uses only the observed trait subset. This marginalizes missing traits across visits; independently whitening each observed subset or filling missing responses with zero would give a different model. Entirely missing visits supply no measurement; prediction grids can retain their IDs. Source rows may be unsorted and use wide or long layout.

Greml transforms both responses and complete sparse incidence, using public gsolve decompositions of matrices no larger than the trait count. Prepared innovations are cached once in unit/time order with identity residual precision. Compiled supplied-covariance prediction consumes that cache; packed PCG reuses its existing all-record cache and pedigree preconditioners. No population-sized or dense time-series covariance is formed. Metadata is O(records times traits); covariance state is O(traits squared). Sparse design history can grow with distinct coefficient columns in an incomplete series, and transformed record size can therefore grow. Existing record adapter shape limits still apply; no fixed application RAM ceiling or series-length cap is added. The original scalar two-row innovations path remains available.

The prepared residual log determinant is reported for auditing; this does not expose REML estimation of serial parameters. Small independent dense references check missing traits, irregular gaps, singleton series, wide/long layouts and both packed factors/preconditioners. Industrial preparation cost and memory for these multivariate models remain unqualified. Uncertainty requests with these models are rejected pending separate validation. Existing PEV qualification limits remain in force.

Focused packed serial checks compare eight factor/preconditioner configurations against independent dense GLS/BLUP equations for mixed-genotyping data with shuffled records, missing visits, actual gaps, singleton series and fixed time effects. Maximum coefficient difference is below 6.3e-14; both residual gates remain unchanged. Additional positive/zero AR1 checks include agreement with the independent-residual model at zero correlation. All 56 focused YAML tests pass, including twenty-one runnable examples and old-driver/parameter rejection. The native prepared-record test also passes dense serial matrix/RHS/diagonal oracles, a 300-record series, near-boundary stability, failure-atomic reads and the existing checksum/corruption regressions. Receipts are build/validation/yaml-packed-serial-regression.log and build/validation/yaml-packed-serial-record-tests.log. This is local focused Windows validation, not Linux/HPC throughput, broad model-combination or full R-package qualification. Existing PEV qualification limits remain in force.

Supplied relationships use the same covariance input mechanism:

inputs:
  relationship_matrix: {type: covariance, file: data/relationship.csv}
relationships:
  K: {type: covariance, input: relationship_matrix}
model:
  random:
    genetic: {group: individual, relationship: K, covariance: genetic}

Here the CSV coordinates are individual/line/sample IDs, not trait labels. The first header is coordinate; row and column IDs must match in order and be unique. Leading-zero IDs remain strings. Every record must map to the declared domain; unobserved levels remain in the model and may be predicted. Several supplied relationship domains may coexist, also alongside the one supported pedigree domain. The native adapter factors the supplied covariance and prepares precision through the public greml relationship constructor; users never provide its inverse. Preparation is dense, with O(n^2) storage and O(n^3) arithmetic per effect; it is not a scalable implicit relationship or APY selector. There is no fixed application RAM ceiling. Sparse relationship input and additional relationship effects alongside ssSNPBLUP remain unsupported.

A type: pedigree input has file and explicit columns bindings for individual, parent1, parent2. Unknown-parent tokens default to "", "0", "NA"; unlisted_parents defaults to error. Duplicate IDs, self-parenting, cycles and unlisted known parents fail validation. A pedigree relationship uses {type: pedigree, input: pedigree_name} and retains supplied external individual order. Probabilistic parents are not exposed. Metafounders require the explicit declaration below; ordinary parent columns do not imply them.

Metafounders and genetic groups

Metafounders define the pedigree covariance base. Declare a natural Gamma covariance and a table assigning each unknown parent slot to a named base:

inputs:
  base_covariance: {type: covariance, file: metafounder_gamma.csv}
  base_assignments:
    type: table
    file: metafounder_assignments.csv
    columns: {individual: sample_id, parent1_group: parent_1_group, parent2_group: parent_2_group}
relationships:
  A:
    type: pedigree
    input: pedigree
    metafounders: {covariance: base_covariance, assignments: base_assignments}

Gamma uses the labelled covariance CSV format: first header coordinate, then ordered metafounder IDs, with matching row labels. It must be finite, exactly symmetric and positive definite, with each diagonal strictly between zero and two. No inversion, clipping or jitter is requested from the user. Base IDs must be disjoint from ordinary pedigree IDs. For example:

coordinate,MF1,MF2
MF1,0.4,0.12
MF2,0.12,0.7

The assignment table contains every ordinary pedigree individual exactly once, joined by ID in any row order. Unknown parent slots require a declared group; known parent slots require an empty cell. A known parent remains an ordinary pedigree ID in the pedigree CSV. Do not insert metafounders into those rows. This miniature pedigree has two founders and two offspring of both founders:

sample_id,parent_1_group,parent_2_group
002,MF2,MF2
004,,
001,MF1,MF1
003,,

Native Pedigree retains the correlated Gamma base and appends latent metafounder levels to the ordinary external order. They have genetic coefficients but cannot carry observed phenotypes or genotypes. The same A can enter H: {type: single_step, pedigree: A, genomic: G}, with supplied effective G or APY blocks; the user owns compatibility of that genomic base and Gamma. No extra base transformation or Gamma estimation is inferred. Reuse the same explicit A declaration across effects in a metafounder model.

Genetic groups instead add explicit independent random mean coefficients: y = X b + Z(a + F gamma) + other + e. The ordinary A or supplied-genomic H prior acts on a. The native pedigree inheritance operator constructs F by propagating half of each parent’s group contribution. Declare this extra effect:

model:
  random:
    genetic: {group: individual, relationship: A, covariance: genetic}
    ancestral_means:
      type: genetic_groups
      group: individual
      pedigree: A
      assignments: group_assignments
      groups: [G1, G2]
      covariance: group_mean
covariances:
  group_mean:
    - [0.8, 0.15]
    - [0.15, 0.5]
output:
  predictions:
    genetic_values: {effects: [genetic, ancestral_means], levels: individuals}

group_assignments is a table input with the same column bindings and coverage rules as above, using G1/G2 labels. groups fixes the coefficient order. Each group level has the same natural trait/basis covariance group_mean, and different group levels are independent: Var(gamma) = I_groups tensor K_group. basis, traits and explicit channels work as for other random effects. The default is mode: random. Compiled models also accept mode: fixed with an explicit reference group and no covariance prior. Combined metafounder/group models retain separate covariance and mean declarations. Correlated between-group priors are not inferred; see the extended adapters below.

Both options support compiled pedigree direct/PCG, supplied-genomic ssGBLUP/APY PCG, and packed ssSNPBLUP. Packed metafounders use supplied or explicitly estimated Gamma and either pedigree_factorization: full_qnn or ancestral_sweeps. The ancestral factor retains every metafounder, including bases without genotyped descendants, and the complete supplied Gamma precision block. Ordinary peripheral individuals use exact innovation sweeps. Either shifted or Qnn-reuse preconditioning is available. The metafounder factors use exact symmetric diagonal equilibration (S Qfactor S) with matching RHS/solution maps for the selected full or ancestral factor. Original Qnn prior actions and the statistical equations remain unscaled. This balances high-degree latent base rows without jitter or relaxed original-equation checks. The extra scale vector and solve staging are accounted for in native resource allowances. Packed genetic groups require one group-mean effect, declared last, sharing the genetic effect’s pedigree, individual column and ordered channels. They support full Qnn or ordinary ancestral sweeps and either preconditioner. Matching curve channels are supported; different group/genetic channel maps remain a compiled-route feature. Combined metafounder/group models are still rejected. R interfaces remain unchanged.

On the compiled route the group adapter materializes F once through native PedigreeGroupProportions, then expands weighted record incidence and releases the dense F after compilation. This costs O(n times groups) setup storage plus expanded incidence; it has no fixed application memory ceiling and is not an industrial-scale qualification. The packed route instead composes public GeneticGroupDataOperator with the cached original records and repeated natural channel covariance. F and F’ remain pedigree sweeps: no n-by-g matrix or expanded group-record incidence is built. Exact group diagonals for preconditioning are prepared by one operator probe per group/channel coordinate, reusing the unchanged base diagonal. This adds setup work proportional to the number of group coefficients and O(primary dimension plus n times channels plus groups times channels) scratch. The original RHS is transformed consistently and reporting uses the same F. Native capabilities packed_metafounders / packed_genetic_groups are checked before result creation. For these new routes PCG solves to an internal target one order tighter than the requested tolerance; both requested transformed and original-equation residual gates remain mandatory and unchanged.

random_effects.csv reports latent metafounder coefficients or gamma group coefficients in their declared order. A genetic-group prediction reports the individual contribution F gamma, so selecting both effects reports their sum. levels: individuals excludes latent metafounders; levels: all includes them for a metafounder genetic effect and uses ordinary individuals for group contributions. Combining effects requires matching prediction domains. Curve predictions use an explicit grid as before. Diagnostics include supplied base/group counts and the materialized F byte count. These routes are checked against independent small tabular/dense equations, not industrial benchmarks.

A type: plink input has prefix and id: IID. Every FAM IID must be unique and match a pedigree individual ID when linked by a single-step relationship. For other sample identities, an optional id_map may reference a table input whose columns bind FID, IID and individual. It requires one row per FAM pair, unique target individuals and exact pedigree matches. This explicit map replaces direct IID matching; no ID harmonization is guessed. Genotype calls count BIM allele 1, with mean imputation of missing calls. All panel samples and markers are used; QC and selection need separate typed extensions.

A marker-based single-step relationship references a named pedigree relationship and a PLINK input. The packed ssSNPBLUP model is specified as:

relationships:
  H:
    type: single_step
    pedigree: A
    genotypes: genotypes
    genomic:
      adjustment: pedigree_moments
      polygenic_weight: 0.25
      independent_ridge: 0.0

Centering defaults to selected-sample allele means and normalization to sum_2p1mp, with unit marker weights. G0=s ZZ'/sum(2p(1-p)), where positive genomic scale: s defaults to 1. Scaling is a covariance-model choice, separate from execution. Adjustment is explicitly none or pedigree_moments. The latter forms G_adjusted=b G0+a11' matching the all-entry mean and mean diagonal to A22. The model is G*=(1-w)G_adjusted+w A22+rI, with H^-1=A^-1+embed(G*^-1-A22^-1). The genomic block must state the adjustment, w=polygenic_weight and r=independent_ridge; these scientific choices have no hidden defaults. The first adapter supports positive w<1 with r=0, or w=0 with r>0. It rejects other combinations until qualified. The existing marker adjustment additionally requires positive b and nonnegative a. No automatic ridge or change of base repairs a failed model.

Supplied genomic covariance and APY ssGBLUP

Exact ssGBLUP accepts the effective genomic covariance directly:

inputs:
  genomic_matrix: {type: covariance, file: data/genomic_covariance.csv}
relationships:
  A: {type: pedigree, input: pedigree}
  G: {type: covariance, input: genomic_matrix}
  H: {type: single_step, pedigree: A, genomic: G}
model:
  random:
    genetic: {group: individual, relationship: H, covariance: genetic}
execution:
  formulation: ss_gblup
  solver: pcg
  sparse_backend: pardiso
  preconditioner: block_jacobi
  threads: {factor: 1, relationship: 8}

The matrix CSV follows the supplied relationship format above. Its coordinates are the genotyped subset, with unique IDs present in the pedigree; its order need not match pedigree or phenotype order. H retains the full pedigree domain. The native operator applies H^-1=A^-1+embed(G^-1-A22^-1) with exact marginal pedigree precision, not a submatrix of A^-1. This compiled route currently uses PCG only, and never constructs full H. Other random effects can use identity, supplied covariance or the same pedigree input. Correlated channels and polynomial regression follow the ordinary compiled-model contract.

The supplied G is the effective covariance, already prepared in the intended pedigree/genomic base. This spelling applies no centering, scaling, moment adjustment, pedigree blend or ridge. To compare against ssSNPBLUP, supply the same effective G* defined by the marker model; these scientific choices are not repeated. The route does not construct G from PLINK. Raw genotype-to-G preparation and additional genomic-base options need a later interface.

APY replaces G’s declaration with explicit covariance blocks:

inputs:
  core_matrix: {type: covariance, file: data/genomic_core.csv}
  cross_matrix: {type: covariance, file: data/genomic_noncore_core.csv}
  diagonal:
    type: table
    file: data/genomic_noncore_diagonal.csv
    columns: {individual: sample_id, variance: genomic_variance}
relationships:
  G:
    type: apy
    core_covariance: core_matrix
    cross_covariance: cross_matrix
    noncore_diagonal: diagonal
    minimum_conditional_variance: 0.0
  H: {type: single_step, pedigree: A, genomic: G}

Core selection is explicit: the ordered row/column IDs of core_matrix are the chosen core. The cross CSV has coordinate followed by the same ordered core columns, with distinct noncore row IDs disjoint from the core. Diagonal table rows must match that noncore order exactly. Values are the original marginal variances G_ii, not conditional variances or precision entries. All blocks must come from one intended effective genomic covariance.

Native preparation computes P=Gnc Gcc^-1 and D_i=G_ii-(P Gcn)_ii. D must be finite and strictly above the explicit nonnegative minimum_conditional_variance (default zero). Failure is rejected, without clipping or added jitter. The APY model replaces conditional noncore/noncore covariance by diag(D), so its implied noncore covariance is P Gcc P' + diag(D). It generally differs from exact G. The resulting H is exact for that supplied APY covariance; convergence does not remove the APY approximation. An all-core declaration is allowed, with an empty noncore cross block and a header-only diagonal table, and recovers exact G. APY currently enters YAML models through the single-step H relationship.

Diagnostics record genomic_approximation: apy, core/noncore counts, minimum prepared conditional variance, operator applications, preconditioner storage, pedigree factor statistics, phase timings and peak native RSS. model.json and run.json retain the explicit core/noncore IDs, genomic subset order, input hashes and covariance-model declarations. They do not report a guessed approximation error; compare with an exact reference where feasible.

These are supplied-matrix/block integrations, not qualified production-scale alternatives to the packed marker route. Exact G preparation uses O(ng^2) storage and O(ng^3) work per effect. Native APY preparation stores O(k^{2+rk+r) numeric/precision data and includes O(k^3+rk}2) work. This frontend currently loads CSV blocks into Python memory, including temporary object/numeric copies; the alternative raw PLINK APY spelling instead reuses native streamed preparation with a retained core factor and file-backed cross products, avoiding full G or genomic precision storage. Complete H covariance metadata streams omitted-individual covariance calculations, which can dominate large-population setup. Fixed-covariance compiled prediction now defaults to execution.relationship_metadata: prediction_only, which skips that unused covariance metadata. The exact precision diagonal, determinant and preconditioner remain unchanged. complete requests the previous metadata work; this control applies only to compiled ssGBLUP. Compact format 7 requires fresh preparation with the updated driver. It does not retain a full H or A22 matrix, but bounded storage alone does not establish fast setup. Keep the existing packed ssSNPBLUP route for the current industrial benchmark; assess these new routes on manageable data before a separate HPC scaling study.

Execution and reporting

Use one thread count for all applicable operations:

execution:
  formulation: ss_snp
  solver: pcg
  sparse_backend: pardiso
  preconditioner: qnn_reuse
  pedigree_factorization: ancestral_sweeps
  threads: 16

Use threads: {default: 16, decode: 4} to override selected operations. Packed keys are marker_products, marker_setup, decode, pedigree, preconditioner and vectors. Compiled models use factor and relationship; raw-genotype APY also accepts decode and genomic_products. Omitted threads or omitted mapping entries without a default use one. These request native phase concurrency; Python input preparation and dependency-ordered work can remain serial. Inapplicable keys are errors. No command-line thread override or compatibility alias is provided.

execution retains explicit solver/backend choices, single-step formulation and preconditioner, linear convergence tolerances and phase thread counts. Compiled PCG supports automatic, diagonal and block_jacobi preconditioners. The default is block_jacobi for supplied-genomic ssGBLUP and automatic for other compiled PCG models. threads.relationship controls retained genomic/ pedigree sparse factor backend work; threads.factor controls the compiled direct-factor request. Both default to one. Dense covariance/APY preparation also includes serial native arithmetic; requested threads do not imply every setup phase is parallel. Explicit CHOLMOD requests require one relationship factor thread; use PARDISO for threaded sparse factors. Backend requests do not silently fall back. The native subset diagonal workspace allowance is unlimited in this adapter, with the existing bounded RHS batching retained. For ssSNPBLUP, execution.pedigree_factorization selects full_qnn (default) or ancestral_sweeps. The latter retains all individual effects and records, factors the non-genotyped ancestral core and uses exact pedigree sweeps for the remaining individuals. Use preconditioner: qnn_reuse to share that complete solve and avoid a separate full shifted factor. This is an execution choice; supplied covariance matrices, genomic scaling and model terms are unchanged. The option is rejected for ordinary compiled-MME routes. The native executable must advertise ancestral_pedigree support before preparation. The ancestral reduction admits ordinary and supplied-metafounder pedigrees. Metafounder models additionally require the native ancestral_metafounders capability; every correlated base remains in the factor, with all original coefficients reconstructed. Existing thread controls apply to the native phases, but dependency-ordered pedigree sweeps can contain serial work. The single-step toy example explicitly uses a tighter 1e-12 PCG tolerance to pass both transformed and original-equation checks. Omitted applicable thread counts are 1; inapplicable options and unsupported route combinations are errors. REML convergence is separate from linear convergence. There is no fixed application RAM ceiling or timeout. Optional future budgets must distinguish factor/workspace allowances from process RSS.

output.directory must be fresh. Default tables are fixed_effects, random_effects and diagnostics; an explicit tables list may select them. Named output.predictions requests specify contributing effects. levels: all reports existing level values for intercept-only effects. Curve requests name a table input through data, with an explicit unique ID and the grouping/basis bindings required by the requested effects. Traits default to all. Prediction inputs do not need response columns. Requests report only the named effects, without adding fixed effects or clipping time. New effect levels and unsupported reporting requests fail before fitting. For one selected effect, optional channels: [direct] restricts the sum to named channels in that effect. Channel names must be unique and present; the default trait set is then inferred from those channels. Explicit trait requests must be represented by the selected channels. For example:

output:
  directory: results/direct-maternal
  predictions:
    direct_values: {effects: [genetic], channels: [direct], levels: individuals}
    maternal_values: {effects: [genetic], channels: [maternal], levels: individuals}
    record_contributions: {effects: [genetic], data: observations}

Level predictions evaluate the selected coefficients at the same named level; they do not substitute that individual’s parent. Data predictions instead use each selected channel’s grouping column, so the final request reports direct[individual] + maternal[mother]. Omitting channels sums all applicable channels. A selected maternal-only prediction table needs the mother binding, but does not need an unused direct-individual column. Basis channels require data predictions; selecting only intercept channels permits level predictions. Separate channel selection currently requires exactly one effect per request; multi-effect data predictions can still sum the full requested effects. Selected ordinary packed genetic values can additionally request uncertainty:

output:
  directory: results/my-fit
  predictions:
    genetic_values: {effects: [genetic], levels: individuals}
  uncertainty:
    method: exact
    effect: genetic
    individuals: ["004", "001"]
    traits: [trait1, trait2]
    maximum_targets: 4

This writes genetic_uncertainty.csv, ordered by the requested individuals and traits. Columns are individual, trait, prediction, prior_variance, pev, prediction_error_standard_error, reliability, method, prior_mcse, pev_mcse, reliability_mcse and qualified. Exact extraction has zero MCSE. PEV includes uncertainty from jointly estimated fixed and other random effects; the standard error is sqrt(PEV). Reliability is 1 - PEV / (Hii * Ktt) for the supplied effective genomic/pedigree covariance and trait variance. It is conditional on the supplied covariance parameters, not uncertainty in estimated variance components or a claim of empirical calibration.

method: exact extracts selected diagonal entries within the numerical solve tolerances. It performs one PCG inverse query per individual/trait pair, reusing the fitted conditional operator, preconditioner and pedigree factor. Both transformed and original augmented residual gates must pass for every target. No dense inverse or stochastic fallback is used. maximum_targets is a required positive solve-work cap, not a RAM ceiling; the request fails before preparation if its number of individual/trait pairs exceeds that cap. All IDs must be unique quoted strings present in the pedigree. Trait order is preserved and defaults to the genetic effect’s declared channels.

Prior variances reuse native exact genomic diagonal preparation and bounded H unit-column extraction. Non-genotyped targets require additional inner PCG work for their prior Hii; these solves and iterations are reported separately. This stage can be expensive for many targets even when the prediction fit is fast. Native uncertainty stage time, target-solve iterations and both maximum residual-to-threshold ratios appear in diagnostics. Estimates are not clipped; failed solves or invalid PEV/reliability bounds fail the requested evaluation.

The initial YAML scope is ordinary-pedigree packed ss_snp PCG with an intercept-only genetic effect with one channel per genetic trait, either pedigree-factor route and either preconditioner, including supplied ridge or pedigree-blended genomic covariance. Metafounders, genetic groups, genetic regression trajectories, compiled-MME uncertainty and arbitrary contrasts are rejected. The native executable must advertise packed_uncertainty_exact before result creation. See the complete example (models/multitrait-single-step-uncertainty.yml).

For a large target set, the alternative is an explicit stochastic request:

output:
  directory: results/genotyped-uncertainty
  uncertainty:
    method: stochastic
    effect: genetic
    individuals: genotyped
    traits: [trait1, trait2]
    probes: 32
    deflation_rank: 16
    seed: 4401
    maximum_solves: 48
    maximum_reliability_mcse: 0.02
    prior: {method: exact, maximum_individuals: 0}

The native estimator uses a fixed Nyström sketch followed by independent Rademacher residual probes. deflation_rank + probes posterior PCG queries must fit maximum_solves; each query reuses fitted factors and passes both equation checks. Rank zero selects plain diagonal probing. At least two residual probes are required, rank is at most the target dimension and 2048, and all seeds are explicit unsigned 64-bit integers. There is no automatic probe increase or exact fallback. These are work budgets, not memory ceilings.

individuals accepts explicit quoted IDs, genotyped (mapped FAM order), or all (pedigree order). Trait order is retained. For genotype-only requests, Hii is exact and requires no unit-column query. The exact prior option’s maximum_individuals caps additional non-genotyped prior diagonals. To probe non-genotyped prior Hii instead, specify:

    prior:
      method: stochastic
      probes: 32
      deflation_rank: 8
      seed: 4402
      maximum_relative_mcse: 0.01

Prior sketch/probe counts bound H applications separately from the posterior maximum_solves; H applications perform inner pedigree PCG work, whose solves and iterations are reported. Genotyped diagonals remain exact. Prior and posterior seeds must differ and cannot differ by the sketch-stream constant 0x9e3779b97f4a7c15. This prevents reuse of either random stream.

PEV MCSE is estimated from residual probes conditional on the sketch. Reliability MCSE uses a first-order delta approximation combining independent PEV and prior MCSE. Neither is a confidence-interval or calibration guarantee; PCG error, covariance-parameter estimation and model misspecification are outside MCSE. Nonpositive/nonfinite prior estimates fail uncertainty. Finite PEV outside [0,prior] remains an unqualified candidate, without clipping; negative PEV has a blank prediction-error standard error.

A row is qualified only when PEV/reliability satisfy their bounds, reliability MCSE meets maximum_reliability_mcse and prior relative MCSE meets its target. run.json retains prediction accepted separately from uncertainty_accepted and uncertainty_status (passed, unqualified or failed). Numerical stochastic failure preserves validated predictions but publishes no uncertainty CSV. Precision/bounds failure retains the candidate CSV and per-row flags. The driver must advertise packed_uncertainty_stochastic. This mode has the same ordinary-pedigree/intercept-only scope as exact extraction. See the complete example (models/multitrait-single-step-stochastic-uncertainty.yml).

Stochastic requests can set estimator: information instead of the default estimator: pev. This probes the explained covariance I = reported_prior - reported_posterior, reusing the same PCG queries, sketch/probe budget and residual checks. It reconstructs PEV = prior - diag(I) and reliability = diag(I)/prior. This can reduce stochastic error for weakly informed targets; the benefit depends on the covariance structure. For genotyped-only targets the additional prior actions apply G and the trait covariance directly, with no pedigree solves or dense G storage. Mixed targets use H actions and incur separately reported inner pedigree solves. These prior actions are additional to maximum_solves, which caps posterior queries. The driver must advertise packed_uncertainty_information.

The CSV includes estimator and information_mcse (zero in PEV/exact mode). For information probing, pev_mcse = hypot(prior_mcse, information_mcse) and reliability_mcse = hypot(information_mcse/prior, reliability*prior_mcse/prior). The latter accounts for the shared estimated prior in reconstructed PEV and its denominator. Prior-diagonal probes use an independent seed. Scope, qualification thresholds and the no-clipping policy remain the same. See the complete information example (models/multitrait-single-step-information-uncertainty.yml).

For stochastic estimator: pev, optional probing: trait masks residual probes to each requested trait in turn; unstructured remains the default. The Nyström sketch is unchanged. probes counts all masked posterior queries, so it must be divisible by the number of requested traits and supply at least two complete groups. At four traits, deflation_rank: 4, probes: 28 and maximum_solves: 32 use four sketch queries plus seven independent groups of four masked queries. MCSE uses seven samples per target, not 28. Ordinary probing with the same controls uses 28 samples per target. Removing cross-trait noise can compensate for this smaller sample count, depending on the model. Diagnostics label the structure and report independent groups and color count. The native executable must advertise packed_uncertainty_trait_probing. This choice is not available for exact extraction or information probing.

The describe operation prints the expanded effects and channels, covariance order, fixed terms, level/record counts and route decisions before factorization. For the curve example it shows 6 genetic channels, 6 individual environmental channels and 2 site-visit channels, with a 2-by-2 residual matrix. Generating that metadata must not materialize a population-sized incidence matrix. Primary scientific coefficients and marker/pedigree-auxiliary solver coordinates are counted separately in native diagnostics after fitting.

Diagnostics retain the resolved model/defaults, identities, named results, phase timings, process RSS, and prediction/original-equation residual checks. run.json is always retained as the execution audit. Selected tables produce fixed_effects.csv (trait,term,estimate), random_effects.csv (effect,level,channel,trait,estimate) and diagnostics.json. Each named prediction writes <name>.csv, with level,trait,prediction for level requests or record_id,trait,prediction for a prediction grid. Grid IDs are retained verbatim. Native phase RSS covers the native process; it does not include Python preparation. End-to-end Python preparation memory is not yet profiled.

Preparation and implementation limits

Ordinary pedigree and joint random-regression models lower through greml’s public compiled-MME facade; single-step models reuse its packed conditional PCG driver. No new numerical solver or R selector is introduced. The focused checks in tools/validation/test_hpc_model_yaml.py cover independent dense equations, missing traits, reordered pedigree IDs, covariance labels, curve predictions, checksums and CLI execution. These checks do not establish performance at industrial scale or full-package qualification.

The frontend spools CSV/TSV input tables, normalized observations, long-input grouping, residual assignments and record prediction grids to a disposable SQLite snapshot. By default this uses the operating system temporary directory (TMPDIR on Linux); optional execution.preparation.directory selects an existing directory, with relative paths resolved beside the YAML file:

execution:
  # Other required execution settings remain here.
  preparation:
    directory: /scratch/myproject

Each CLI invocation owns a fresh temporary child directory and removes it on normal completion or handled failure; an externally killed process may leave that child behind. The parent directory and source files are preserved. Disk use grows with input rows, and preparation makes multiple sequential passes over the snapshot. No application RAM or disk ceiling is imposed. Pedigree and effect-level IDs/indexes, PLINK sample/marker metadata, supplied dense relationships and covariance matrices remain in memory. Rank checking uses a fixed-term Gram matrix, quadratic in fixed coefficients per trait. Prepared record prediction targets live in checksum-protected prediction_*.jsonl files rather than an embedded population-sized JSON array; reporting reads these files one row at a time. The frontend writes sparse record references, without materializing a population-sized incidence matrix. The native text adapter has a maximum of 256 responses per record, 65,536 incidence entries per record and 64 genetic channels on the packed marker route; these are adapter shape restrictions, not RAM ceilings. Packed ssSNPBLUP requires one H effect, declared first, with identity relationships for the other effects. Ordinary compiled models and supplied-genomic ssGBLUP can reuse one pedigree input across multiple random effects. Supplied-genomic ssGBLUP currently allows one H effect; genomic and pedigree references must be explicit and valid. The current compiled-MME facade requires at least one fixed coefficient across the model; an entirely empty fixed design is rejected before preparation. The packed marker route supports an empty fixed design. Covariance matrices must be positive definite. New prediction levels, unsupported genomic combinations, unsupported REML structures and uncertainty outside the ordinary packed genetic scope above are rejected rather than inferred. YAML duplicate keys, anchors, aliases, custom tags and unknown keys are rejected. The draft schema has runtime validation; a separate formal machine-readable schema has not been published.

Supported routes

Features in this table do not imply that every cross-combination is supported. Invalid combinations are rejected before fitting. See the complete example catalog for runnable files rather than copying isolated fragments.

Route	Prediction	Covariance estimation	Likelihood
Compiled pedigree/identity/dense relationships	Direct or PCG; regression, correlated/weighted incidence, groups/metafounders, FA and serial/general residual blocks	Sparse-direct unstructured/FA/serial REML with named fixed coordinates; operator relationships excluded	Exact absolute restricted likelihood and line search
Compiled exact/APY ssGBLUP	PCG; supplied effective G/blocks or streamed raw PLINK APY; pedigree individuals retained	Not exposed for operator H relationships	Not an estimation route
Packed ssSNPBLUP	PCG; ancestral/full Qnn, reuse/shifted preconditioning, separate/weighted incidence, supplied FA/serial covariance, groups/metafounders	Ordinary pedigree, one genetic channel per trait, unstructured genetic/independent residual matrices; whole-block fixing; fixed additional identity-effect covariances	Optional exact or SLQ start-to-final change; zero global determinant work when omitted

Wide/long CSV/TSV observations lower to the same joint residual occasions. Sparse external natural covariance uses retained sparse factors; factor fill and exact diagonal work remain model-dependent. Raw-genotype APY uses explicit core IDs and tiled/file-backed products; industrial throughput is unqualified. Ancestral metafounder reduction retains every correlated base. Packed fixed/combined groups prepare O(individuals × groups) inheritance columns; ordinary random groups remain implicit. Packed Gamma pseudo-EM is conditional on an explicit absolute marker base, separate from joint covariance REML.

Selected exact/stochastic PEV is conditional on supplied covariances and has its own work/precision qualification, separate from prediction acceptance. General PEC/trajectory uncertainty, probabilistic pedigrees, general level-specific or structured marker covariance, scalar interval constraints, shared parameter tying and zero-variance boundary fitting are not exposed. Industrial REML precision/runtime, industrial APY preparation, new large-data combinations and complete package qualification remain separate work.

Estimation and additional model declarations

Compiled models use transport format 8 with the matching native driver. The declarations below use the same input/model/analysis/execution sections. Existing R interfaces and PEV qualification limits are unchanged. Complete examples include small CSV/PLINK inputs; small reference checks do not establish industrial performance.

REML starts, constraints and acceptance

Random-intercept REML example (models/reml-random-intercept.yml):

covariances:
  genetic:
    - [2.0]
  residual:
    - [2.0]
analysis:
  task: reml
  estimate: [genetic, residual]
  optimizer:
    maximum_iterations: 100
    maximum_line_search_reductions: 24
    maximum_step: 1.0
    score_tolerance: 1.0e-7
execution:
  solver: sparse_direct
  sparse_backend: cholmod
  threads: {factor: 1}

The natural covariance matrices are starting values. estimate names the free whole covariance blocks; omitted blocks are held exactly at their supplied values. Omitting estimate estimates every referenced block. At least one must be free. Log-Cholesky coordinates enforce positive definite covariance, and line search accepts nondecreasing restricted likelihood. maximum_step bounds the transformed-coordinate step. Convergence requires stationarity of the free transformed scores, not small changes alone. A fixed block’s score need not be zero. Named optimizer coordinates may also be fixed, as detailed below. Unstructured covariance entries are not independently fixed by that control: it fixes the named Cholesky element. Zero-variance boundary fits, cross-component equality tying and user-defined scalar interval bounds are not exposed.

REML supports ordinary/explicit-base compiled models, distinct covariance names for each component and sparse-direct solving, using PARDISO or CHOLMOD. FA loadings/specifics and scalar or multivariate AR1/exponential parameters are estimated through native derivatives. Operator relationship REML is not exposed by this compiled adapter. Packed EM-REML is a separate route below. Exact derivative work can be substantial; this is a moderate-problem REML route, not an industrial marker REML claim.

variance_components.csv labels starting and fitted covariance entries and fixed/free status. reml_history.csv retains likelihood/score/step diagnostics. convergence.json reports convergence and covariance/BLUP availability. Its structured likelihood block (also written as likelihood.json) reports absolute likelihood for compiled REML and optional start-to-final likelihood change for packed REML, with separate status/precision diagnostics. Nonconverged fits return failure and retain their diagnostics; they do not publish qualified predictions. Linear residual checks still apply at the accepted final covariance.

Packed EM-REML and optional likelihood

The fast packed route uses the same covariance starts and estimate names, with independent covariance-fit and likelihood work controls:

analysis:
  task: reml
  estimate: [genetic, residual]
  optimizer: {maximum_iterations: 200, score_tolerance: 1.0e-6, step_fraction: 1.0}
  traces:
    method: stochastic
    probes: 64
    seed: 4401
    maximum_solves: 10000
    maximum_score_mcse: 0.001
  likelihood:
    method: slq
    maximum_operator_products: 6144

step_fraction is the EM interpolation fraction in (0,1]; the compiled optimizer’s maximum_step is inapplicable. Whole blocks can be fixed; packed coordinate constraints are not yet supported. Exact trace references use traces: {method: exact, maximum_solves: 10000}. These budgets are examples, not production precision recommendations. The ordinary-pedigree/unstructured scope, score/MCSE gates, independent final bank and solve accounting are specified in packed REML.

Omit likelihood to do no global determinant work. method: exact requires maximum_dimension for bounded small matrix references. SLQ requires an operator-product cap and accepts probes, seed, lanczos_steps, check_steps, maximum_mcse and maximum_depth_change. It evaluates both start/final states at both depths, using paired samples. Likelihood outputs describe a change relative to the supplied start; they are not absolute likelihoods comparable across different models or starts. Sampling MCSE and depth sensitivity are separate; the latter is not an error bound. This optional diagnostic does not control the fit or reject qualified predictions if its precision targets fail. The complete exact (models/multitrait-packed-reml.yml) and SLQ (models/multitrait-packed-reml-slq.yml) examples include all files and controls.

Fixed group means, explicit phantom codes and combined bases

Fixed groups (models/multitrait-fixed-genetic-groups.yml) use the group declaration family, with no random prior:

model:
  random:
    genetic: {group: individual, relationship: A, covariance: genetic}
    ancestral_means:
      type: genetic_groups
      mode: fixed
      group: individual
      pedigree: A
      assignments: group_assignments
      groups: [G1, G2]
      reference: G1

This adds inherited, unpenalized group contrasts to the fixed design for each selected trait. The reference group’s mean is zero. Full fixed-design rank is checked natively; a reference declaration does not guarantee identifiability in every dataset. Estimates appear in fixed_effects.csv; named predictions can sum genetic and ancestral_means as usual.

Explicit phantom group-parent codes (models/multitrait-phantom-groups.yml) can replace the separate assignment table:

inputs:
  pedigree:
    type: pedigree
    file: data/pedigree_phantom.csv
    columns: {individual: sample_id, parent1: parent_1, parent2: parent_2}
    phantom_groups: {PG1: G1, PG2: G2}

Each code designates a missing-parent group slot and is absorbed into the native inherited proportion design. It does not create an individual, a shared biological parent, or a latent phantom coefficient. Every unknown parent slot needs an explicit group, either through codes or an assignment table. If both describe a slot they must agree. Unlisted biological parent IDs still cause an error. A group-mean model must explicitly consume these codes.

Combined bases (models/multitrait-combined-base.yml) keep supplied metafounder covariance and fixed or independent random inherited group means separate: a has the metafounder pedigree covariance, while F is computed from the original ordinary pedigree and missing-parent assignments. Latent metafounders have no observed group-mean incidence. Group and genetic random coefficients are independent. No covariance/base conversion is inferred. The compiled adapter materializes ordinary pedigree group proportions once; the existing packed random-group route keeps its implicit inheritance sweeps. Fixed groups and combined bases currently use compiled pedigree or supplied ssGBLUP/APY models; their packed adapters are rejected explicitly.

Gamma estimation (models/multitrait-estimated-gamma.yml) uses the existing native pseudo-EM estimator during packed preparation:

relationships:
  A:
    type: pedigree
    input: pedigree
    metafounders:
      covariance: Gamma_start
      assignments: base_assignments
      estimation: {method: pseudo_em, maximum_iterations: 2000}
  H:
    type: single_step
    pedigree: A
    genotypes: genotypes
    genomic:
      adjustment: none
      scale: 1.0
      polygenic_weight: 0.0
      independent_ridge: 0.05
      centering: absolute_half
      normalization: markers_over_two

This explicitly selects allele-frequency one-half centering (dosage center one) and normalization M/2. Pseudo-EM uses the raw absolute-base marker covariance; the declared independent ridge/blend belongs to the subsequent evaluation model. The supplied Gamma is a start. Failed or nonconverged estimation blocks prediction. Prepared Gamma is fingerprinted and used to rebuild the prediction pedigree. metafounders.csv and run.json.gamma_estimation report the estimate, iterations and fixed-point residual. Pseudo-EM is not exact REML, does not certify biological calibration, and prediction is conditional on the estimate.

Weighted multiple-individual incidence

Complete example (models/multitrait-weighted-incidence.yml):

model:
  random:
    genetic:
      group: individual
      relationship: A
      covariance: genetic
      incidence:
        - {individual: individual, weight: 0.6}
        - {individual: contributor, weight: contribution}

individual names a mapped logical ID column. weight is a finite literal or a mapped numeric column. This yields 0.6*a[individual] + contribution*a[contributor]. Weights are not normalized; negative contrasts are allowed. Duplicate contributors are coalesced, zero contributions omitted, and a zero-weight contributor may have a missing ID. Explicit channel-level incidence overrides the effect-level mapping; group and incidence cannot both override the same channel. A channel-level group keeps the existing single-individual override. All nonzero contributor IDs must belong to the relationship domain. Compiled and packed genetic products use the existing sparse record design; packed preparation uses all records. Observation predictions use the same mixture, while level predictions remain per individual. Group-mean mixture mappings and weighted PEV requests are not exposed.

Sparse external natural relationship covariance

Complete example (models/multitrait-sparse-relationship.yml):

inputs:
  relationship_levels:
    type: table
    file: data/relationship_levels.csv
    columns: {individual: sample_id}
  relationship_covariance:
    type: sparse_covariance
    file: data/relationship_sparse.csv
    levels: relationship_levels
relationships:
  E: {type: sparse_covariance, input: relationship_covariance}

The CSV header is row,column,value; IDs use the declared levels-table order. Supply unique lower-triangle entries sorted by column index then row index, including every positive diagonal. Missing off-diagonal entries mean zero. The input is covariance, not precision. Symmetry follows the triangle contract; native sparse factorization certifies positive definiteness. The frontend streams entries to disk and native PCG retains a sparse covariance factor as its inverse-action operator. It does not form a dense covariance inverse. This supports extra compiled-model effects, not a new raw-genomic H route. Use threads.relationship and an explicit sparse backend.

Storage depends on factor fill, not only input nonzeros. Exact precision diagonals use selected inversion with CHOLMOD; the current PARDISO adapter uses bounded batches of unit solves. That fallback avoids quadratic storage but can be expensive for large relationships. Industrial setup time and arbitrary out-of-core relationship plugins remain unqualified.

Supplied nonseparable trait-by-visit residual covariance

Complete example (models/multitrait-nonseparable-residual.yml):

model:
  residual:
    type: supplied_blocks
    group: individual
    position: time
    covariance: residual
    coordinates:
      - {trait: trait1, position: 20}
      - {trait: trait2, position: 20}
      - {trait: trait1, position: 21}
      - {trait: trait2, position: 21}

The supplied SPD matrix spans these coordinates in exactly this order. It can have visit-specific trait covariance and any admissible serial cross-covariance. Each group receives one independent copy. Missing trait/visit observations use the observed marginal submatrix, with no imputation or conditioning on missing responses. Duplicate observed coordinates within a group are errors. Wide and long inputs are equivalent. Named covariance row/column labels, if used, are trait@position, with numeric positions formatted without unnecessary zeros. This compiled adapter is intended for moderate within-group coordinate grids; large dense grids and their parameterization storage can be costly. It does not add a parametric nonstationary serial kernel, packed residual-block route, or serial-parameter REML.

Focused qualification of the extended adapters

Local Windows validation on 2026-10-16 passed 19 focused Python reference and affected-route regression tests, the native compiled-MME optimizer check and all eight added complete YAML examples. Balanced univariate and multivariate REML recover analytic covariance estimates; whole fixed matrices remain unchanged. Other references assemble the observation design, covariance and mixed-model equations independently, including combined exact/APY H models, weighted packed products, missing nonseparable residual coordinates and both sparse backends. Nonconvergence, indefinite covariance and stale source inputs are rejected. Existing packed, direct/maternal, FA, serial and raw APY cases remain included. These checks do not constitute a new industrial benchmark, Linux qualification, broad package qualification or resumed PEV precision work.

The existing R wrapper was rebuilt against the extended native optimizer control structure, and canonical package runtime verification passed. A separate test-correlated-formula.R check has 29 passing and two failing expectations: it assumes a one-iteration nonconverged REML fit has successful status and score output. The committed native optimizer already reports iteration exhaustion as failure, and the unchanged R wrapper returns no score output for that failure. Those existing expectations and the R interface were left unchanged; the local R test suite is not claimed fully green.

The reviewed fixed-block REML path checks average-information curvature only on free transformed coordinates, including when the unrestricted covariance components would be confounded. The unrestricted model still rejects an unusable information matrix. The portable 50K validation recipe exports 14 complete YAMLs from the preserved fixture and checks seven equivalent model pairs, a bounded analytic REML reference and immutable input provenance. Full-fixture prediction and moderate-sized REML qualification remain separate.

The final review additionally passed all 12 extended model tests with only PARDISO selected, the two portable fixture-workflow tests (including deliberate output-tamper rejection), and the native optimizer regression check. Both starts on the actual fixture’s bounded REML submodel agree with its independent residual-variance reference, 1.9994769667549. The full 50K paired campaign is provided for GenomeDK execution; it has not been run locally as part of this bundle refresh. Existing PEV qualification limits remain in force.

FA and serial REML, coordinate fixing and precision allocation

Format 8 adds runnable FA REML (models/factor-analytic-reml.yml), multivariate serial REML (models/multivariate-serial-reml.yml), packed fixed groups (models/multitrait-packed-fixed-groups.yml), packed combined bases (models/multitrait-packed-combined-base.yml) and independent-pilot PEV (models/multitrait-single-step-pilot-uncertainty.yml) examples. These additions leave the R interfaces unchanged. Reprepare compiled inputs with the matching frontend and driver; frozen earlier bundles are unchanged.

Covariance parameters and fixed coordinates

analysis.estimate selects complete parameter blocks. For a serial residual, the block includes its covariance and its correlation/rate. For FA it contains the identified free loadings and positive specific variances. Omitted blocks stay fixed at their supplied natural starts. A selected block can additionally hold named optimizer coordinates at their starts:

analysis:
  task: reml
  estimate: [genetic, residual]
  constraints:
    genetic:
      fixed: ['loading[trait2,shared]', 'specific[trait1]']
    residual:
      fixed: [correlation]
  optimizer:
    maximum_iterations: 100
    maximum_step: 1.0
    score_tolerance: 1.0e-6

FA names use the declared channel and factor labels. Fixing a positive loading or specific variance fixes its log coordinate, hence the same natural value. An unstructured multivariate block names lower Cholesky entries in column order, such as cholesky[trait2,trait1]; diagonal entries are optimized on the log scale. This holds a Cholesky element, not a covariance matrix entry. Other covariance entries can change. A scalar covariance uses variance[trait1]. A scalar serial residual uses variance and correlation/rate; a multivariate serial block uses Cholesky names followed by correlation/rate. describe lists the exact names in analysis.parameter_coordinates. Unknown, duplicate, entirely fixed or inapplicable declarations are rejected. Arbitrary scalar interval bounds, cross-block equality tying and zero-specific/zero-variance boundary fitting are not exposed.

The optimizer applies the Jacobian, AI submatrix and stationarity test to the free coordinates. Fixed whole natural blocks are restored exactly. A fixed direction with a nonzero score does not prevent a valid free-submodel fit. parameters.csv separates natural coordinate names/values from optimizer coordinate names/fixed status. variance_components.csv reports matrix entries and their whole-block status. For FA REML, covariance_models.json distinguishes the starting model from fitted parameters and implied covariance.

Serial parameter estimation

For a scalar model, reference a named 1×1 covariance in variance:

covariances:
  residual:
    - [1.0]
model:
  residual: {type: ar1, group: individual, position: time, variance: residual, correlation: 0.4}

For multivariate data, use covariance: residual with a natural trait covariance matrix. AR1 uses integer positions and a correlation strictly within (-1,1); exponential uses numeric positions and a strictly positive rate. Native REML rebuilds the exact observed trait-by-time covariance and its analytic derivatives at each trial. Missing traits/visits are marginalized by the observed submatrix. It does not freeze the prediction whitening transform while estimating parameters. Trait covariance and serial parameters can be estimated jointly or selectively fixed through the same block controls.

This is deterministic sparse-direct REML for moderate compiled problems. It stores dense covariance/derivative blocks within each series, and the usual global sparse derivative arrays. The supplied-covariance packed prediction route retains its existing Markov preparation. Industrial packed REML and long-series REML scalability have not been qualified. Failed fits retain diagnostics and publish no qualified predictions; no convergence or likelihood gate is relaxed automatically.

For ordinary pedigrees with unstructured genetic covariance and independent residual blocks, the separate packed EM-REML route uses analysis.traces to select exact or stochastic PCG moment traces. It shares the pedigree factor across covariance updates and reports scores/MCSE rather than a restricted likelihood. Whole-block fixing is supported; FA, serial parameters, groups and metafounders are outside this initial packed REML scope.

Packed fixed and combined base models

The fixed-group syntax and ordinary inheritance definition above also work with formulation: ss_snp, either pedigree factorization and either supported packed preconditioner. Combined models keep the metafounder covariance pedigree and the original ordinary group-inheritance pedigree separate. The small references compare all coefficients, group contributions and total predictions with dense equations under both pedigree factorization routes.

For these additions, preparation computes F once in parent-first order and writes a fingerprinted NumPy array, with 8 × ordinary-individuals × groups bytes per group effect, plus metadata/array overhead. describe.group_proportion_bytes reports this storage. Fixed reference contrasts or random inherited means are lowered into record columns; random group coefficients retain the supplied independent group prior. Packed genotype products and the exact genetic model are unchanged. The original ordinary random-group route retains implicit pedigree sweeps. Models with many group levels can make the explicit preparation array expensive; no application RAM ceiling is imposed. Packed inherited random means currently require the common, unweighted individual incidence used in these references; their combination with weighted or separate channel mappings is rejected pending validation.

Independent-pilot PEV allocation

PEV precision work has resumed. Existing fixed-budget estimators remain available. For unstructured stochastic PEV with exact prior diagonals, optionally add:

output:
  uncertainty:
    method: stochastic
    effect: genetic
    individuals: genotyped
    probes: 28
    pilot_probes: 8
    deflation_rank: 4
    seed: 4401
    maximum_solves: 128
    maximum_reliability_mcse: 0.02
    prior: {method: exact, maximum_individuals: 0}

The independent, nondeflated pilot estimates raw posterior sampling variance. For each target its pilot reliability MCSE predicts the required final count using inverse sample-count variance scaling. The maximum across requested targets is rounded up, with probes as the minimum and maximum_solves - pilot_probes - deflation_rank as the maximum. The pilot can be conservative because the final estimator may use deflation.

Pilot samples are discarded. Final probes use a separate stream and a fixed count conditional on the pilot, avoiding reuse and repeated-stopping selection of the final diagonal/MCSE. Diagnostics separately count pilot queries, allocated and completed final probes, sketch queries and all posterior solves. Both equation residual gates still apply to every solve. The precision target can remain unmet; such finite candidates remain unqualified and are never clipped. MCSE is Monte Carlo sampling error, not a confidence-coverage or model-calibration certificate. More probes can improve precision at extra runtime. This does not establish a better estimator at the same total solve budget or industrial performance.

Focused Windows qualification passed 25 affected model/reference tests, four native compiled-MME gates and the five complete new examples. Additional CHOLMOD checks cover FA/scalar constraints and serial likelihoods. Final direct dense comparisons also check complete natural scores and AI matrices for FA and both serial families. The nonlinear serial covariance curvature term is included, and both strict 1e-8 and default-control serial references converge on PARDISO and CHOLMOD. Iteration-limited nonconvergence remains a rejected fit. The canonical package-only rebuild and R runtime smoke checks passed; this is not a full R CMD check or GenomeDK qualification. Receipts are under build/validation/model-parameter-* and build/validation/parameter-model-examples-20261016/. In the small independent pilot comparison, mean reliability MCSE changed from 0.25277 with four final probes to 0.05464 with 64 final probes plus four discarded pilot queries. The final values exactly match a separate fixed-64-probe run. The requested precision target was still unmet and remained unqualified; this demonstrates work allocation and error reduction, not production precision or coverage.

Optional prediction-truth evaluation

The optional top-level evaluation mapping references a named individual-level output.predictions request, a long truth CSV and its column mappings:

evaluation:
  predictions: genetic_values
  truth:
    file: true_breeding_values.csv
    columns: {individual: animal, trait: trait, value: true_breeding_value}
  subset: heldout
  individual: individual
  genotypes: genotypes

predictions, truth, subset and individual are required. genotypes is an optional PLINK input reference for stratified reporting. individual names a mapped model-data column. subset is heldout, phenotyped or all; cohort membership uses individuals with any observed response after wide/long input lowering. Evaluation does not mask records or change the model. Record/grid prediction requests are not evaluation targets. Unknown keys, missing or nonfinite truth values, duplicate labels and unknown target levels/traits are rejected. validate checks settings and targets; describe includes resolved paths and cohort counts. Only run performs scoring after numerical acceptance.

Outputs, metric definitions, qualification and the frontend version boundary are described in the standalone guide.