Estimate Combined Percent Familial Risk Explained by GWAS Loci

Most current genome-wide association studies (GWAS) include a calculation of the percent familial risk the discovered loci explain.  This calculation indicates how much of the familial risk can be accounted for by the known loci and is usually used as evidence there are additional yet undetected loci that remain to be discovered.  Looking through references, it can be a bit difficult to find exactly how this calculation is performed.  One reference I found that includes a formula for the calculation is by Cox et al. 2007 (PMID:17293864), but I am sure there are plenty others that also include a formula.   While the Cox et al. formula for calculating familial relative risk due to each locus is arranged differently than the ones below, the two are equivalent.  I just find this arrangement less cumbersome to use.  The overall equation is to compare the cumulative risk of the known loci (sum of log lambda k) to the estimated risk of a first degree relative (log lambda 0).

where:

      p is the risk allele frequency for locus k

      r is the per allele odds ratio for locus k.

To make calculations easy, I made a simple R script that does all the calculations automatically.  The input for the script is a file with 3 columns:

(1) Annotation for the SNP - this can be anything, for example: RS number, chromosomal coordinates, etc.
(2) Risk allele frequency - this is the frequency of the risk allele (range: 0-1) equal to p in the above equation.
(3) Per allele odds ratio - odds ratio for every one unit increase in the number of risk alleles.

Note, the risk allele frequency is the frequency of the risk allele and not the minor allele frequency.  The program also needs an estimate of the familial relative risk (lambda 0).  This can usually be done by looking for previous familial studies for the disease.

Here is the R script:

	# Import command line arguments
	args <- commandArgs(trailingOnly=T)
	snp_filename <- args[1]
	fd_fam_risk <- as.numeric(args[2])
	# Import SNP data
	snps <- read.table(snp_filename, as.is=T, header=F)
	names(snps) <- c("id", "p", "r")
	attach(snps)
	# Observed familial risk to first degree relatives
	print("Familial Risk Estimate (OR):")
	fd_fam_risk
	# Calculate familial risk due to each locus
	snps$lambda <- (p*r**2+(1-p))/((p*r+(1-p))**2)
	snps
	# Contribution of known SNPs to familial risk
	print("Familial Risk Explained (%):")
	sum(log(snps$lambda))/log(fd_fam_risk)*100

view raw familial_risk_snps.R hosted with ❤ by GitHub

It can be run from the command line by the example command:

Rscript familial_risk_snps.R snp_lst.txt 4

where:
      familial_risk_snps.R is the name of the script.
      snp_lst.txt is the input file with three columns described above.
      4 is the estimate of the familial relative risk of the disease.