GWAS-Assisted Best Linear Unbiased Prediction (GWABLUP)

GWAS-Assisted Best Linear Unbiased Prediction (GWABLUP)

Description

Genomic prediction using GWAS-weighted relationship matrix. Implements Meuwissen et al. (2024) GWABLUP: weights all markers by their posterior probability of having non-zero effects, derived from EMMAX GWAS results.

Usage

gwablup(
  y,
  markers,
  gwas_result,
  X = NULL,
  method = "auto",
  solver = "auto",
  max_iter = 100L,
  tol = 1e-06,
  n_threads = NULL,
  min_weight = 1e-04,
  ref_markers = NULL,
  trait = "continuous",
  link = "logit"
)

Arguments

y numeric named vector of phenotypes (length n)
markers

list of raw marker inputs:

  • snp_add: matrix (n x m), raw genotype 0/1/2

  • mh_add: list of chr data.frames

gwas_result object of class “gwas_result” from run_gwas(). Must be run with the same markers and individuals.
X fixed effects design matrix (n x c). If NULL, intercept only
method character, REML algorithm (default "auto"). Passed to masreml()
solver character, EBV solver (default "auto"). Passed to masreml()
max_iter integer, maximum REML iterations (default 100)
tol numeric, convergence tolerance (default 1e-6)
n_threads integer, number of threads (default: all cores)
ref_markers list of raw marker data for training individuals only, same format as markers. If provided, allele frequencies for building G_wa are computed from ref_markers instead of markers. Use in train/test context to avoid data leakage. If NULL (default), allele frequencies computed from markers.
trait character, trait distribution: “continuous” (default) or “binary”. Passed to masreml().
link character, link function for binary trait: “logit” (default) or “probit”. Passed to masreml().

Details

User workflow:

  1. Fit standard GBLUP: fit <- masreml(y, markers)

  2. Run GWAS: gwas <- run_gwas(markers, y, fit)

  3. Inspect GWAS results (optional): summary(gwas)

  4. Run GWABLUP: fit_wa <- gwablup(y, markers, gwas)

Value

Object of class c(“gwablup”, “masreml”) — identical structure to masreml() output, with additional $gwas slot containing the gwas_result. All S3 methods from masreml (summary, print, plot) work directly.

References

Meuwissen et al. (2024) GWABLUP: genome-wide association assisted BLUP. Genet Sel Evol. 56:17.

See Also

masreml, run_gwas, compute_accuracy

Examples

Code
library("masreml")

d       <- load_data("small")
y       <- d$pheno$y_cont_qtl_snp; names(y) <- d$pheno$id
W_train <- d$snp[d$train_idx, ]
W_test  <- d$snp[d$test_idx, ]
y_train <- y[d$train_idx]

# ── SNP GWABLUP train/test ──────────────
fit_tr  <- masreml(y_train, markers = list(snp_add = W_train))
gwas_tr <- run_gwas(
  markers     = list(snp_add = W_train),
  y           = y_train,
  masreml_fit = fit_tr,
  ref_markers = list(snp_add = W_train)
)
fit_wa_tr <- gwablup(
  y           = y_train,
  markers     = list(snp_add = W_train),
  gwas_result = gwas_tr,
  ref_markers = list(snp_add = W_train)
)

summary(fit_wa_tr)

G_full <- build_G_snp(d$snp, ref_W = W_train)
pred   <- predict(fit_wa_tr,
                  G_full    = list(snp_add = G_full),
                  train_ids = rownames(W_train),
                  test_ids  = rownames(W_test))
head(pred$GEBV)