Skip to contents

Bioconductor-style infrastructure for breeding and quantitative genomics: genotypes, phenotypes, pedigree and marker coordinates in one object, with the relationship matrices that quantitative genetics runs on.

Bioconductor has excellent containers for expression, variants and single-cell data. It has nothing for the shape of data a breeding programme actually produces: a marker matrix, a phenotype table, a pedigree that reaches back further than the genotyping did, and a marker map. This package provides that container, and the matrices computed from it.

Why a container

Keeping genotypes, phenotypes and pedigree in separate objects means every filter or reordering is a chance to desynchronise them. BreedingExperiment extends RangedSummarizedExperiment, so:

  • markers are rows, individuals are columns, and subsetting keeps everything aligned;
  • marker positions are a GRanges, so a data set can be cut to a region, overlapped with annotation, or handled with any Bioconductor tool that speaks genomic ranges;
  • the pedigree travels with the data, and may include ancestors that were never genotyped.

Installation

From R-universe:

install.packages("BreedingExperiment",
                 repos = c("https://mqfarooqi1.r-universe.dev",
                           "https://bioconductor.org/packages/release/bioc",
                           "https://cloud.r-project.org"))

Or from GitHub:

# install.packages("BiocManager")
BiocManager::install("SummarizedExperiment")
remotes::install_github("mqfarooqi1/BreedingExperiment")

The package is under review for Bioconductor; once accepted it will install with BiocManager::install("BreedingExperiment").

Use

library(BreedingExperiment)

data(demoBreeding)
demoBreeding
#> class: BreedingExperiment
#> markers: 900  individuals genotyped: 120
#> assays(1): genotype
#> phenotypes(6): generation sex yield stature trueBV_yield trueBV_stature
#> sequences(10): chr1 chr2 chr3 chr4 chr5 chr6
#> pedigree: 180 individuals (120 genotyped, 60 not; 30 founders)

demoBreeding is a simulated five-generation population shaped like a real programme: the two earliest generations were born before genotyping began, so 60 of the 180 individuals have no markers; one trait is recorded on females only; a few sires dominate; and about one per cent of calls are missing. Use simulateBreeding() to generate your own.

Your own data goes in the same way:

be <- BreedingExperiment(
  genotypes  = dosage_matrix,        # markers x individuals, coded 0/1/2
  phenotypes = phenotype_table,      # one row per individual
  pedigree   = pedigree_table,       # id, sire, dam
  rowRanges  = marker_granges)       # marker positions

Relationship matrices

Function Matrix Source
Amatrix() Pedigree numerator relationships Henderson (1976)
Gmatrix() Genomic relationships from markers VanRaden (2008)
Dmatrix() Dominance relationships Vitezica et al. (2013)
Hmatrix() Single step: pedigree and markers combined Legarra et al. (2009); Christensen & Lund (2010)
A <- Amatrix(be)     # every individual in the pedigree
G <- Gmatrix(be)     # the genotyped ones
H <- Hmatrix(be)     # all of them, using both sources of information

The single-step matrix

Real programmes genotype only part of the population. Using only G throws away every ungenotyped relative; using only A throws away the markers. Hmatrix() keeps both, so genotyped and ungenotyped individuals can be analysed in one model — the relationship structure behind single-step genomic prediction (ssGBLUP).

G is first tuned to the scale of A22 and then blended with it, which also guarantees the result is invertible when markers are fewer than individuals.

Quality control

be <- filterMarkers(be, min_maf = 0.05, min_call_rate = 0.9)
be <- imputeMarkers(be)
alleleFrequency(be)      # freq, maf, call rate, heterozygosity per marker
inbreeding(be)           # from the diagonal of A

Correctness

The relationship matrices are checked against values that follow from theory, not against a previous run: parent–offspring 0.5, grandparent 0.25, full sibs 0.5, the offspring of a full-sib mating inbred at F = 0.25, H positive definite, and Hmatrix(blend = 1) reproducing A22 exactly.

Honest limitations

  • No mixed-model solver. This package supplies the relationship structure; fit the model with whichever solver you prefer. That boundary is deliberate — solvers are a crowded field, containers are not.
  • The A matrix is dense. The tabular method is O(n²) in memory, which is fine for tens of thousands of individuals but not for a national evaluation with millions. A sparse inverse is the obvious next step.
  • Imputation is mean imputation, adequate only for the low missingness left after quality control. Impute properly upstream for anything more.
  • Diploid, biallelic markers only.

References


MIT licensed. Muhammad Farooqi.