APY preparation from PLINK genotypes

The complete runnable YAML example (models/multitrait-ssgblup-apy-plink.yml) connects raw PLINK data to native streamed APY preparation and compiled ssGBLUP PCG. Input covariance matrices remain natural covariances. This is prediction with supplied parameters; PEV precision work remains paused.

inputs:
  genotypes: {type: plink, prefix: data/genotypes, id: IID}
  core:
    type: table
    file: data/apy_core_ids.csv
    columns: {individual: sample_id}
relationships:
  A: {type: pedigree, input: pedigree}
  G:
    type: apy
    genotypes: genotypes
    core: core
    centering: sample_mean
    normalization: sum_2p1mp
    scale: 0.95
    offset: 0.0
    independent_ridge: 0.05
    minimum_conditional_variance: 0.0
    row_tile: 256
    marker_tile: 256
  H: {type: single_step, pedigree: A, genomic: G}
execution:
  formulation: ss_gblup
  solver: pcg
  sparse_backend: pardiso
  preconditioner: block_jacobi
  relationship_metadata: prediction_only
  threads: {factor: 1, relationship: 1, decode: 4, genomic_products: 2}

Add the phenotype/pedigree bindings, trait covariances, model and outputs as shown in the complete example. The core CSV contains a header sample_id and one individual ID per row. Its nonempty, unique IDs determine core order; all other mapped FAM samples become noncore in physical FAM order. PLINK id_map supports explicit FID/IID-to-individual mapping, including repeated IIDs in distinct families. Every mapped individual must be an ordinary pedigree individual. Core selection is explicit; automatic selection and marker selection/weights are outside this spelling.

Genomic covariance

All physical BIM markers have unit weight. Let c be their observed mean BIM allele-1 dosages across all FAM samples, p=c/2 and Z the centered dosages. Missing calls are replaced by c. The covariance before APY approximation is

G = scale * Z Z' / sum(2*p*(1-p)) + offset * 11' + independent_ridge * I.

Only sample_mean and sum_2p1mp are accepted. Scale defaults to 1 and must be positive; offset defaults to 0 and must be nonnegative. Ridge must be explicit and nonnegative. It is an independent genomic ridge. No pedigree moment matching, A22 blending, automatic ridge or conditional-variance clipping is inferred. All-missing markers and zero total normalization fail; individual monomorphic markers contribute zero after centering. APY conditional variances must strictly exceed minimum_conditional_variance, which defaults to zero.

APY keeps the core covariance and core/noncore covariance, replacing noncore conditional off-diagonals with zero. Thus it changes the genomic covariance model; PCG convergence does not eliminate that approximation. Selecting all genotyped individuals as core gives exact G with quadratic core storage.

Storage, threads and reporting

Python’s centering pass reads one marker at a time. Native preparation uses bounded row/marker tiles, threads.decode for decoding and threads.genomic_products for dense products/core solves. The latter also controls streamed cross products during prediction. Both default to one; multithreaded dense products require native MKL support. The centering pass is serial. threads.relationship controls the sparse pedigree factor; threads.factor retains its compiled direct-factor meaning.

Row and marker tiles default to 256. They control work granularity; there is no hardcoded application RAM ceiling or timeout. Numeric genomic preparation RAM is O(core^2) plus bounded tiles and linear metadata. The core Cholesky factor is reused; no full G, genomic inverse, decoded panel, resident cross block or projection is constructed. Cross storage is O(noncorecore) on disk. Core factorization costs O(core^3), and conditional solves cost O(noncorecore^2). Compiled ssGBLUP prediction defaults to execution.relationship_metadata: prediction_only. It retains the exact H precision diagonal, determinant, pedigree factor and preconditioner, while omitting H covariance diagonals that supplied-parameter prediction never uses. This removes repeated per-non-genotyped-individual covariance actions. An omitted diagonal is absent, not approximated. relationship_metadata: complete restores exact covariance metadata for diagnostics and matched setup comparisons. Automatic covariance starts/reliability calculations that require that metadata must use a complete native relationship; this prediction-only YAML route does not request them. Native constructors keep complete metadata by default. Prepared compiled files use format 6 and must be regenerated from YAML with the refreshed driver.

Genomic preparation and prediction still have independent scaling limits, including core size, cross-file traffic, precision-diagonal preparation and PCG convergence; no industrial APY timing claim follows from this optimization.

Results include native/apy_N_cross.bin and native/apy_N_metadata.csv for each raw genomic preparation. Keep the cross files with the results. The CSV contains numeric pedigree IDs, covariance/precision diagonals and noncore conditional variances; model.json maps these back to source IDs and records core/noncore order. run.json contains phase timing, peak native RSS, thread requests, metadata policy and omitted-column count, core size, cross-file size, normalization and minimum conditional variance. Sources and prepared centers are fingerprinted. Failed preparation can leave partial files in a fresh rejected output directory.

Fixed-covariance compiled prediction also omits unused global REML first/second covariance-derivative matrices and RHS arrays. It retains the same numerical MME, residual precision caches and preconditioner; full estimation keeps its existing derivative assembly. This removes quadratic-in-parameter-count storage from prediction without imposing a memory ceiling. The refreshed native capability is prediction_value_only_assembly.