Introduction

This example estimates heritability and genetic correlation for type 2 diabetes (GWAS1) and coronary artery disease (GWAS2) using gcorr and a gact database. LDSC, GNOVA and SumHer reuse stored LD scores. HESS uses the BED reference to prepare regional projections once, then reuses them for fitting.

Prepare the data

Use the database and filtered EUR reference from the sparse-LD example. gact supplies the harmonized statistics and sample sizes; Glist already contains the LD scores. The short interfaces below require gcorr 0.1.8 or later.

library(gact)
library(gcorr)

GAlist <- readRDS("path/to/GAlist_hsa.0.0.1.rds")
Glist <- readRDS("path/to/Glist_1000G_eur_filtered.rds")

studyIDs <- c("GWAS1", "GWAS2", "GWAS6")
stat <- getMarkerStat(GAlist=GAlist, studyID=studyIDs)

# Keep the common three-study marker selection; fit T2D and CAD.
stat <- list(z=stat$z[, 1:2, drop=FALSE],
             n=stat$n[, 1:2, drop=FALSE])

Genome-wide heritability and genetic correlation

fit_ldsc <- gcorr(Glist, stat, method="ldsc")
fit_ldsc$estimates

fit_gnova <- gcorr(Glist, stat, method="gnova")
fit_gnova$estimates

fit_sumher <- gcorr(Glist, stat, method="sumher")
fit_sumher$estimates

The wrapper matches marker names, reuses stored scores and supplies delete blocks. LDSC estimates intercepts and excludes squared Z above 80. GNOVA uses the LDSC cross-trait intercept by default; SumHer uses an equal-marker baseline with multiplicative scaling. The extracted stat$n is unchanged.

For LDSC, scores can also be retrieved directly from gact:

ldscores <- getLDscoresDB(GAlist=GAlist, ancestry="EUR", version="1000G")
fit <- gcorr::ldsc(z=stat$z, n=stat$n, ldscores=ldscores, what="rg")
fit$rg

The direct interface leaves squared-Z filtering optional. Recorded example results, with the cutoff of 80, are estimate ± standard error:

Method T2D h² CAD h² Genetic correlation
LDSC 0.470 ± 0.025 0.311 ± 0.019 0.384 ± 0.028
GNOVA 0.498 ± 0.028 0.111 ± 0.019 0.604 ± 0.059
SumHer baseline 0.488 ± 0.029 0.362 ± 0.023 0.387 ± 0.027

The SumHer baseline is illustrative. GNOVA uses a unit marginal noise baseline; LDSC and this SumHer fit estimate intercepts or scalings.

Chromosomes as marker sets

Define each chromosome as a set of marker names. These disjoint sets can reuse the ordinary LD scores, assuming negligible LD between chromosomes.

sets <- getMarkerSets(GAlist=GAlist, feature="Chromosomes")
fit_sets <- gcorr(Glist, stat, sets=sets, method="ldsc",
  set_model="grouped", overlap_intercept=0)
head(fit_sets$estimates)

Here the cross-trait intercept is fixed to zero, assuming no sample overlap or cross-trait confounding. Heritabilities and covariances sum across sets; correlations do not. Other annotations need their corresponding tagging scores.

Regional HESS and rho-HESS

Use Glist with independent LD regions already prepared by gprep(task="eigen", sets=sets), and Z oriented to its reference alleles. HESS uses the saved regions and samples, projecting both traits together. The defaults assume no sample overlap.

fit <- gcorr(Glist, stat, method="hess", nthreads=4)
fit$estimates

For these data, raw CAD HESS heritability is poorly determined and the genetic correlation exceeds one; that ratio is not interpretable as a biological correlation. The estimates use the database’s sample-size convention and are not liability-scale heritabilities. Genomic-control adjustments remain sensitivity analyses. See the gcorr guide for these options and assumptions.

Time and memory

Recorded on an Intel Core i7-1365U with 16 GB RAM. Memory is the whole R process’s resident peak with the required analysis inputs loaded; data preparation, loading and saving are outside the timed step. These timings measure the numerical fits; a GNOVA workflow call also estimates its LDSC intercept unless one is supplied. HESS preparation uses four threads.

Step Time Peak resident memory
LDSC 9.3 s 3.02 GiB
GNOVA 6.4 s 3.09 GiB
SumHer baseline 12.9 s 3.53 GiB
Chromosome sets and fit 29.9 s 3.27 GiB
HESS preparation 18 min 49 s 3.18 GiB
HESS pair fit from prepared moments 3.3–3.6 s 0.43 GiB