Computes, for each biallelic SNV on a single chromosome, the LD score \(\ell_j = \sum_k \tilde r^2_{jk}\) summed over biallelic SNVs \(k\) within \(\pm\)window bp of \(j\), from the cohort's own genotypes. This is the LD-score input to ldsc_regression.

compute_ld_scores(
  gds_file,
  variant_id = NULL,
  window = 1e+06,
  segment = 2000L,
  sample_n = 0L,
  seed = 1L,
  verbose = 1
)

Arguments

gds_file

Path to a GDS/aGDS file for one chromosome. Opened read-only and closed on exit.

variant_id

Optional integer vector of variant.ids to restrict the variant set (both targets and LD-score neighbours). Default NULL uses all biallelic SNVs in the file.

window

Integer LD-score window half-width in bp (default 1e6, i.e. \(\pm\)1 Mb).

segment

Integer number of core variants processed per BLAS block (default 2000L). A memory/speed knob only; results are invariant.

sample_n

Integer. If > 0, a fixed random subsample of that many individuals is used as the LD reference (default 0L = all samples).

seed

Integer seed for the sample_n subsample (default 1L); unused when sample_n = 0L.

verbose

Integer verbosity (0 = silent, 1 = a one-line summary).

Value

A data.frame with one row per biallelic SNV: chr, pos, ref, alt, variant_id, and ld (the LD score), ordered by pos.

Details

The unbiased \(r^2\) estimator (Bulik-Sullivan et al. 2015) is used, $$\tilde r^2_{jk} = r^2_{jk} - \frac{1 - r^2_{jk}}{n - 2},$$ so the LD score is \(\ell_j = \sum_{k:\,|\mathrm{pos}_k-\mathrm{pos}_j|\le w} \tilde r^2_{jk}\). Variants are processed in contiguous core segments of segment variants; for each segment the genotypes of \([\mathrm{segment} \pm \mathrm{window}]\) are read, mean-imputed, column-standardised, and a single BLAS cross-product gives the window-by-core correlation block, from which each core variant's LD score is a \(\pm\)window prefix-sum over \(\tilde r^2\). Peak memory is bounded by one segment plus its flanks, never a whole chromosome.

Single-chromosome contract. Positions repeat across chromosomes, so the \(\pm\)window mask is only meaningful within one chromosome; the function therefore requires gds_file to hold a single chromosome (and errors otherwise). For a genome-wide table, call it per chromosome and rbind (mirrors how marginal_scan is driven per-chr).

LD scores are a population property estimated from genotypes only; they do not depend on the phenotype or covariates, so one cached table serves every covariate config. sample_n > 0 uses a fixed (seeded) random reference subsample of individuals for speed; the unbiased estimator absorbs the resulting small-\(n\) bias.

Ancestry caveat. In an ancestry-mixed cohort, in-sample LD scores carry some cross-ancestry LD. The LDSC intercept (the confounding read-out that ldsc_regression returns) is robust to LD-score mis-scaling, which the slope absorbs.

References

Bulik-Sullivan B. et al. (2015). LD Score regression distinguishes confounding from polygenicity in GWAS. Nat Genet 47:291.

Examples

if (FALSE) { # \dontrun{
ld22 <- compute_ld_scores("chr22_essentialdb.gds")
# Genome-wide: loop per-chromosome GDS files and bind.
ld <- do.call(rbind, lapply(sprintf("chr%d.gds", 1:22), compute_ld_scores))
} # }