Single Nucleotide Polymorphism Spectral Decomposition (SNPSpD)

NOTE:
SNPSpD is not designed to analyse >100 SNPs - use SNPSpDlite or preferably SNPSpDsuperlite!
If you wish to estimate the effective number of independent markers on larger datasets then I recommend
i) download the SNPSpD scripts from the link provided below to perform SNPSpD analysis on your local machine,
ii) my matSpDlite approach (which only requires users to upload a correlation matrix),
iii) download the matSpDlite.R R script to perform matSpDlite analysis on your local machine, or
iiv) use the Genetic Type I error calculator (GEC).

Users may also be interested in the leaner version of SNPSpD, (SNPSpDlite) which ONLY calculates Meff and MeffLi (i.e., does not perform time-consuming varimax/promax rotations) thus allowing users to obtain Meff/MeffLi values for large numbers (1000s) of SNPs.

Users may also be interested in the related site (matSpD) which only requires users to upload a correlation matrix.

Please use the following reference when reporting results based on SNPSpD:
Nyholt DR (2004) A simple correction for multiple testing for SNPs in linkage disequilibrium with each other. Am J Hum Genet 74(4):765-769.
 

To run SNPSpD using all fully genotyped family members:

Notes:

• May 7, 2013: LD (r) values are now rounded to 3 decimal places because LD smaller than 3 decimal places is considered negligible.

• March 15, 2007: SNPSpD was installed on this FASTER server.

• Please note: a significant (or nearly significant) result using the SNPSpD Meff/MeffLi adjusted alpha which is not already significant using a standard Bonferroni (or Sidak) correction (for M tests), would indicate substantial non-independence (i.e., strong intermarker LD) among the tested SNPs and one may therefore benfit from performing permutations. Permutations are still considered the gold standard for multiple test correction.

• For further discussion on the utility of determining an approximate estimate of the effective number of independent SNPs please see Salyakina et al (2005), Nyholt  (2005) and Li and Ji (2005).

• SNPSpD has been updated (August 10, 2005) to also estimate the reportedly more accurate estimate of Meff [MeffLi] proposed by Li and Ji (2005). I recommend using MeffLi only if it is smaller than Meff (view a problematic example). I've found MeffLi to be highly comparable to permutations for datasets containing many SNPs.

• SNPSpD has been updated (July 1, 2005) to only use individuals genotyped at all loci and automatically remove redundant SNPs (see following point).

• As noted in the above paper, the SNPSpD approach can be conservative in the presence of higher order intermarker LD (i.e., very strong LD across >2 SNPs).  Upon closer inspection, it can be seen that by removing redundant SNPs [i.e., SNPs in complete LD (r=1) with another SNP] will greatly minimise the conservativeness of the approach when you have a combination of perfectly (r=1) and imperfectly (r?1) correlated SNPs [i.e., this 'problem' is essentially due to collinearity, where variables are so highly correlated that it is impossible to come up with reliable estimates of their individual coefficients].  For example, in the Keavney dataset there are 2 completely redundant sets {A240T, T93C} and {G2215A, I/D, G2350A}.  Clearly, these represent only 2 effective SNPs.  Therefore one should remove the redundant SNPs (e.g., T93C, I/D and G2350A) and re-run the data through SNPSpD.  To demonstrate: using the following "trimmed" pre and map files, produced the following "corrected" results for the Keavney data.  The resulting experiment-wide significance threshold (required to keep the Type I error rate at 5%) of 0.014 compares favourably with the permuted threshold of 0.015.
• A recent paper by Horne and Camp (2004), recommended an oblique rotation (rather than an orthogonal rotation) to remove inter-factor correlation. In their paper they say "...genetic theory suggests that groups of SNPs within a gene will have some degree of inter-relatedness due to origin on a common ancestral haplotype (or on a limited set of haplotypes). Thus, since oblique rotation can account for such relatedness, it likely is a more appropriate method for genetic analysis".  Hence, I have updated SNPSpD to also include results based on oblique rotation utilising the PROMAX function of R (however, these two rotation methods don't seem to indicate different htSNPs).
• To assist users (and myself) identify whether there is a problem with input files, SpD and/or rotation, I have separated the results into 4 sections [Meff, varimax rotation, promax rotation, and LDMAX output].

These files are run through a slightly altered version of Gon?alo Abecasis' ldmax program - part of the GOLD (Command Line Tools) package [gold-1.1.0.tar.gz].  Please refer to the ldmax documentation for specific details on the format of these input files.

The following links show example input files for a subset of the Keavney et al. (1998) data:
keavney.pre
keavney.map
Results from post-July 1 2005 SNPSpD analysis of Keavney data.
Results from pre-July 1 2005 SNPSpD analysis of Keavney data.
View original paper.

The following links show example input files for the Moffatt et al. (2000) data:
moffatt.pre
moffatt.map
Results from post-July 1 2005 SNPSpD analysis of Moffatt data.
Results from pre-July 1 2005 SNPSpD analysis of Moffatt data.
View original paper.

As described in the program's documentation, ldmax uses the expectation-maximization (EM) based approach of Slatkin and Excoffier (1995) Mol Biol Evol 12:921-7, to estimate haplotype frequencies.  Using these haplotype frequencies, ldmax outputs a table summarising pairwise linkage disequilibrium (LD) statistics to a file named LD.XT.  From this output, a perl script is used to create a matrix of the disequilibrium measure delta () , where

, and

the notation for estimated haplotype and marker allele frequencies in a 22 table are

Risch and Devlin (1995) Genomics 29:311-22, provide a summary of common LD measures and note  is the correlation coefficient for a 2 x 2 table [Hill and Robertson (1968) Theor Appl Genet 38:226-231].

Eigenvalues (s) for this LD correlation matrix may be calculated by principal components (factor) analysis, or more generally by spectral decomposition (SpD).  SNPSpD obtains s by spectral decomposition using the EIGEN function of R [R is a language and environment for statistical computing (R Development Core Team 2003).

The  for a given factor measures the variance in all the SNP-SNP correlations which is accounted for by that factor.

Higher correlation among the SNPs leads to higher s.

For example, if all SNPs are maximally correlated, the first  is equal to the number of SNPs represented in the matrix while the rest of the s are equal to zero.  In this situation the variance of the s is at its maximum and equal to the number of SNPs (M) in the matrix.

Conversely, when there is no correlation among SNPs, all of the s are equal to 1 and the set ofs has no variance.  Therefore, the variance of thes will range between 0, when all the SNPs are independent, and M.

Hence, the ratio of observed  variance to its maximum (M) gives the proportional reduction in the number of independent SNPs in a set, and the effective number of SNPs (Meff) may be calculated as follows:

Selecting a subset of SNPs:

The SpD of the LD correlation matrix also allows investigation of which SNPs best represent a group of SNPs in high LD with eachother.  That is, by examining the factor "loadings" for each , one may determine which SNP best captures the information contained in each group.

Although during the development of this web interface, Meng et al. (2003) Am J Hum Genet 73:115-130, independently came to the same conclusion, I have nonetheless extended SNPSpD to report results after orthogonal rotation.  The results include the s, principal component coefficients, and factor "loadings" after orthogonal rotation utilising the VARIMAX function of R.

I have maximised interpretability of these results by flagging the SNP(s) contributing the MOST to each rotated factor (group of SNPs).  These flagged SNPs may be viewed as "haplotye tagging SNPs".  Indeed, even in data with strong LD the rotated factors correspond well with haplotypes obtained via traditional methods.  For example, the 7 haplotypes reported in the Keavney et al. (1998) paper correspond to the 7 factors produced by SNPSpD after varimax rotation.

Finally, in order to select the minimum subset of SNPs which maximise information, I propose using the effective number of SNPs (Meff) as a guide to choosing the total number of haplotyes to tag [i.e., choose the tagging SNPs which best represent the factors with the Meff highest rotateds].  For example, the 23 SNPs in the Moffatt data have low intermarker LD, producing an Meff of 22.53, thus indicating the SNP loading highest on the factor with the smallest  [i.e., factor 23; =0.0483; SNP7 (hADV14S1A)] provides the least amount of independent information [i.e., it is highly correlated with SNP 9 (hADV14S1C)].  Analogously, the 10 SNPs in the Keavney data have very high intermarker LD, producing an Meff of 4.59, indicating the 5 SNPs loading highest on the factors with the 5 largest s [i.e., either one of SNP 7 (g2215a) / SNP 8 (i/d) / SNP9 (g2350a); SNP 1 (t-5991c); SNP 3 (t-3892c); SNP 6 (t-1237c); and SNP10 (4656ct)] provide MOST of the independent information of the 10 SNPs

Alternatively and perhaps more statistically appealing, you can choose the tagging SNPs which represent the factors explaining a selected proportion of variance.  For example, only the last factor explains ?1% of the variance in the Moffatt data.  Analogously, for the Keavney data; the first three factors explain 95% (94.65%) of the variance, while the first four factors explain 99% (98.74%) of the variance, etc.