Comparing GREML, Haseman–Elston, and LDSC Estimates

Why Compare the Estimators?

GREML, Haseman–Elston regression, and LD score regression can all produce estimates called SNP heritability, but they do not use the same data or fit the same statistical model. Comparing their numerical outputs is informative only after checking whether the estimates refer to comparable SNP sets, populations, phenotypes, and scales. Differences can reflect sampling error, different estimands, or violations of method-specific assumptions; a discrepancy by itself does not identify which explanation applies.

The common target is usually the proportion of phenotypic variance attributable to additive effects tagged by a specified set of SNPs, in a specified population and on a specified phenotype scale. “Tagged” matters: a SNP set may capture effects through linkage disequilibrium even when it does not contain every causal variant. The target is therefore not automatically the proportion explained by all genetic differences.

The methods are most directly comparable when they analyze the same phenotype in the same target population, use compatible SNP coverage and allele-frequency handling, and express heritability on the same scale. Even then, they estimate the target in different ways, so their estimates need not coincide in a finite sample.

What Each Method Estimates

GREML uses individual-level genotypes and phenotypes. It models the phenotype as fixed covariate effects plus a random additive genetic effect and residual error. A genomic relationship matrix, constructed from the analyzed SNPs, describes the covariance expected between individuals under the genetic random-effects model. In the simplest model,

Here \(\mathbf{X}\) contains modeled covariates, \(\mathbf{K}\) is the relationship matrix for the analyzed SNPs, and the genetic and residual effects are assumed independent. Restricted maximum likelihood (REML) estimates the variance components while accounting for the fitted fixed effects. Under this model, the corresponding heritability estimate is \(\hat{h}^2_{\mathrm{SNP}}=\hat{\sigma}_g^2/(\hat{\sigma}_g^2+\hat{\sigma}_e^2)\). It relies on the relationship matrix representing the relevant tagged genetic covariance and on adequate control of structure and other sources of covariance.

Haseman–Elston (HE) regression is a moment-based approach using individual-level genotypes and phenotypes. It relates products of phenotypes from pairs of individuals to their genetic similarity. In the simple model, for distinct individuals \(i\) and \(j\), the expected phenotype product changes with their genetic relatedness \(K_{ij}\) at a rate determined by the genetic variance. With a standardized phenotype, an estimated slope can be interpreted as an estimate of genetic variance, and hence of heritability if total phenotypic variance is one. HE does not require fitting the same likelihood as GREML, but it still depends on the genetic similarity measure, phenotype model, and appropriate handling of pair dependence and confounding.

LD score regression (LDSC) uses GWAS summary statistics and LD scores, typically calculated from a reference panel. It regresses association-statistic inflation across SNPs on their LD scores. Under its model, SNPs tagging more variants tend to have larger expected association statistics when those variants contribute to polygenic signal. An intercept is included to accommodate inflation that is approximately independent of LD score; interpreting it as a diagnostic requires the assumptions of the regression model.

These methods therefore differ in the information they use. GREML and HE use the covariance among individual phenotypes and genotypes. LDSC uses the pattern of GWAS association statistics across SNPs and an LD reference. LDSC is not simply GREML applied to summary data: it identifies signal from a relationship between association strength and LD score.

Worked Estimates

Why the LDSC Slope Tracks Tagged Heritability

A short derivation explains the role of LD score. For SNP \(j\), let its LD score be \(l_j=\sum_k r_{jk}^2\), summing squared correlations with the SNPs in the relevant set. Under a simplified polygenic model, suppose the per-SNP effect variances are equal, with \(\operatorname{Var}(\beta_k)=h^2_{\mathrm{SNP}}/M\). The noncentral part of the expected association statistic at SNP \(j\) is approximately \(N\sum_k r_{jk}^2\operatorname{Var}(\beta_k)\). Substitution gives

Adding the null contribution and a term for approximately LD-score-independent inflation yields the basic regression form

where \(a\) is the intercept under this parameterization. The derivation is an approximation: its interpretation depends on the genetic architecture, the LD scores, the GWAS statistics, and the treatment of confounding. The intercept is not a universal correction for any bias that happens to affect the GWAS.

Uncertainty and Comparability

All three methods are subject to sampling uncertainty, but their standard errors arise from different calculations. GREML commonly obtains uncertainty from the fitted likelihood or an approximation to it; variance-component estimates can be imprecise near the parameter boundary. HE uncertainty must reflect the dependence among pairs, rather than count pairs as independent. LDSC uncertainty is based on the regression and commonly uses genomic blocks to account for correlation among nearby SNP statistics. These standard errors do not necessarily include uncertainty from a mismatched reference panel or an incorrect model.

A numerical comparison should begin with the target, not the point estimates. Check the phenotype scale, analyzed SNPs and their coverage, ancestry and sample composition, sample size, and covariates. For binary traits, confirm that the estimates are on the same scale; a liability-scale conversion is not interchangeable with an observed-scale estimate. Check method-specific assumptions too: a valid relationship matrix and covariance model for individual-level methods, and suitable summary statistics and LD scores for LDSC.

Diagnosing Disagreement

A useful diagnosis separates differences in the target from differences in the estimator. Start by asking whether the methods use the same trait definition and scale, the same population, and comparable SNP coverage. A different SNP set or LD reference can change what genetic variance is tagged. Differences in sample ascertainment, phenotype measurement, and covariate adjustment can also change the effective target or the covariance pattern being estimated.

Next consider assumptions and data quality. For GREML and HE, unmodeled structure or other sources of phenotype covariance related to genetic similarity can affect estimates. For LDSC, inaccurate or ancestry-mismatched LD scores, problematic GWAS calibration, and a poor fit of the assumed relationship between association statistics and LD score can affect the slope. These possibilities do not imply that one method must always be higher or lower than another.

Finally, compare uncertainty, not just point estimates. A difference between two point estimates may be small relative to their uncertainty, and the estimates may be correlated. If results disagree materially, report the methods, SNP sets, scales, uncertainty procedures, and relevant diagnostics together. A single estimate should not be selected merely because it is larger, more precise-looking, or closer to a preferred value.