GWAS summary statistics (p values) were downloaded from the ImmunoBase website ([www.immunobase.org](http://www.immunobase.org)). Raw T1D genotype data are available from dbGaP (<http://www.ncbi.nlm.nih.gov/projects/gap/cgi-bin/study.cgi?study_id=phs000180.v2.p2>).

Introduction {#sec001}
============

Genome-wide association studies (GWAS) have enabled identification of genetic variants associated with a wide range of complex phenotypes, but in many cases these variants explain only a proportion of the known heritability \[[@pgen.1004926.ref004]\]. There is increasing evidence that this is due to the combined contribution of small effects arising from multiple distinct variants \[[@pgen.1004926.ref005]\]. The testing of a large number of potential variants in parallel, with a comparatively low number of samples, mandates a stringent threshold for significance in order to limit false positives (type 1 errors), meaning that discovery of variants responsible for small effects requires very large sample sizes. Detection of such variants by increasing numbers of samples in studies is time-consuming and expensive, particularly for rare phenotypes, but it may be possible to improve detection by re-analysis of existing data \[[@pgen.1004926.ref006]\]. One promising strategy is to co-analyse GWAS results from similar phenotypes to exploit potential similarities in genetic aetiology. This has been attempted using several different methods \[[@pgen.1004926.ref002], [@pgen.1004926.ref007], [@pgen.1004926.ref008]\].

The assumption that GWAS for similar diseases may yield overlapping sets of disease-associated variants is based on the phenomenon of pleiotropy, in which a genetic variant is associated with more than one trait or disease \[[@pgen.1004926.ref009]\]. Pleiotropy is common in human genes: even when exclusively considering single nucleotide polymorphisms (SNPs) with strong evidence of association, around 15% of those associated with at least one trait are associated with multiple traits \[[@pgen.1004926.ref010]\]. Elements of shared genetic aetiology may be suspected in diseases with similar symptomatology, such as bipolar disorder and schizophrenia \[[@pgen.1004926.ref011]\] or in diseases with common risk factors, such as type 2 diabetes and obesity \[[@pgen.1004926.ref012]\]. If two diseases are known or suspected to share associated genetic variants, a degree of association of a locus with one disease may increase the likelihood of association with the other. Use of external covariates in this way can alleviate some of the effect of multiple testing \[[@pgen.1004926.ref008]\], meaning that phenotypic similarity may lead to improved detection of disease-associated variants. Correspondingly, discovery and specification of shared genetic aetiology between two diseases may suggest some shared pathophysiology \[[@pgen.1004926.ref012]\].

A technique for improved discovery of disease variants using pleiotropy between pairs of diseases has been successfully developed and applied by Andreasson et al \[[@pgen.1004926.ref003], [@pgen.1004926.ref013], [@pgen.1004926.ref014]\]. The technique extends the empirical Bayesian false discovery rate \[[@pgen.1004926.ref015]\] to a two-phenotype scenario, in which association with one phenotype is tested conditional on varying degrees of association with another. We denote the phenotype for which association is being tested as the 'principal phenotype' and the other as the 'conditional phenotype'.

By successively restricting attention to SNPs with a given strength of association in the conditional phenotype, the number of parallel tests to perform for association with the principal phenotype is reduced. If the two phenotypes share common associated variants, this restriction will retain disease-associated SNPs at a higher rate than null SNPs, resulting in a higher proportion of disease-associated SNPs in the restricted group than in the whole. The 'conditional false discovery rate' (cFDR), defined as the probability that a SNP is not associated the principal phenotype given its p values for the principal and conditional phenotypes are below some thresholds, exploits this effect. By computing cFDRs for schizophrenia conditioned on bipolar disorder and vice versa, Andreasson et al \[[@pgen.1004926.ref003]\] identified multiple previously undiscovered loci for both. In a separate study computing cFDRs for hypertension conditioned on 12 related traits \[[@pgen.1004926.ref013]\], 42 new loci associated with hypertension were reported. These constituted considerable improvement on existing results using single GWAS, albeit using a rather relaxed threshold of estimated *cFDR* ≤ 0.01.

A major disadvantage of the algorithm developed and used by Andreasson et al is the requirement that control groups for the two GWAS be distinct, in order to ensure that observed effect sizes are uncorrelated at null SNPs. This requires splitting a pool of potential controls between studies, with the summary statistics for each GWAS computed from only the controls allocated to that study. This may be impractical as it requires access to raw genotype data. More importantly, accuracy of effect size estimates improves with larger control groups, and consequently splitting controls in this way weakens the effect size estimates for individual studies. For this reason, many researchers employ a study design in which controls are pooled into a large group; for example, the Wellcome Trust Case Control and ImmunoChip consortia \[[@pgen.1004926.ref016], [@pgen.1004926.ref017]\].

Here we extend the cFDR approach to studies with overlapping control groups, exploiting an approach developed by Zaykin et al, following Lin et al \[[@pgen.1004926.ref018], [@pgen.1004926.ref019]\] to adjust for the effect of shared controls. This allows the strongest available estimates of effect sizes to be used for calculation, and consequently strengthens the power of the technique. Our technique additionally allows cFDR rates to be computed from summary statistics alone, without the need to recalculate effect sizes after re-allocating controls. We demonstrate the improvement arising from sharing controls in a type 1 diabetes data set.

We also identify a previously undiscussed difficulty with the technique potentially leading to a falsely low estimate of false discovery rate amongst SNPs declared non-null. Multiple overlapping sets of SNPs may be defined each of which has *cFDR* ≤ *α*. However, the union of these sets does not necessarily have an expected false-discovery rate less than *α* and is generally higher. An implication of this is that if we declare non-null all SNPs for which estimated cFDR is less than *α*, the expected overall false-discovery rate amongst SNPs declared non-null is greater than *α*. We describe an upper bound on the false discovery rate amongst such SNPs based on areas of regions of the unit square.

We apply our method to summary SNP association statistics for ten phenotypically distinct autoimmune diseases: type 1 diabetes (T1D) \[[@pgen.1004926.ref020]\], autoimmune thyroid disease (ATD) \[[@pgen.1004926.ref021]\], coeliac disease (CEL) \[[@pgen.1004926.ref022]\], multiple sclerosis (MS) \[[@pgen.1004926.ref023]\], narcolepsy (NAR) \[[@pgen.1004926.ref024]\], primary biliary cirrhosis (PBC) \[[@pgen.1004926.ref025]\], psoriasis (PS) \[[@pgen.1004926.ref026]\], rheumatoid arthritis (RA) \[[@pgen.1004926.ref027]\], ulcerative colitis \[[@pgen.1004926.ref028]\], and Crohn's disease \[[@pgen.1004926.ref028]\]. All were genotyped using a common SNP array: the ImmunoChip, designed to provide dense genotype coverage of regions associated with autoimmune disease. Many autoimmune traits are known to have significant heritability, much of which remains unexplained \[[@pgen.1004926.ref029]\]. We hypothesised that our method can improve detection of disease-associated variants in these diseases without the need for distinct control groups.

Results {#sec002}
=======

Overview of method {#sec002a}
------------------

The unconditional false discovery rate for a set of SNPs with p values \< *p* ~*i*~ is defined as the probability that a random SNP from this set is null. We denote this as *uFDR*(*p* ~*i*~), and our estimate as $\hat{uFDR}(p_{i})$.

The conditional false discovery rate (*cFDR*) is defined \[[@pgen.1004926.ref003], [@pgen.1004926.ref014]\] as the probability that a random SNP is null for a phenotype *i* given that the observed p values at that SNP for phenotypes *i* and *j* are less than (*p* ~*i*~,*p* ~*j*~); that is, $Pr(H_{0}^{(i)} \mid P_{i} \leq p_{i},P_{j} \leq p_{j})$, where $H_{0}^{(i)}$ is the null hypothesis that the SNP is not associated with phenotype *i*. We denote this quantity as *cFDR*(*p* ~*i*~∣*p* ~*j*~), and call phenotype *i* the 'principal phenotype' and phenotype *j* the 'conditional phenotype'.

We first apply genomic control to allow the assumption that, globally, P values for null SNPs are uniformly distributed on \[0, 1\]. We compute an estimate of the *cFDR*, which we denote $\hat{cFDR}(p_{i} \mid p_{j})$, in a similar manner to that proposed by Andreasson et al, but incorporating expected non-uniformity in the distribution of *P* ~*i*~ due to the sharing of controls. As $\hat{uFDR}(p_{i})$ is monotonically related to *p* ~*i*~, we set a significance cutoff at the maximum value of $\hat{uFDR}(p_{i})$ with *p* ~*i*~ \< 5 × 10^−8^. Correspondingly, we set a significance cutoff for $\hat{cFDR}(p_{i},p_{j})$ at the maximum $\hat{cFDR}(p_{i},p_{j})$ with *p* ~*i*~ \< 5 × 10^−8^. Implementation of these steps in R is available from <https://github.com/jamesliley/cFDR-common-controls>.

Sharing of control subjects {#sec002b}
---------------------------

If no controls are shared between studies, it is reasonable to assume that observed effect sizes for the two phenotypes are independent under a null hypothesis for the principal phenotype. This implies that the expected quantile of a given SNP's p value for the principal phenotype is simply the p value itself regardless of its p value for the conditional phenotype. However, when control samples are shared, this assumption is invalid. Shared controls induce a positive correlation on estimated effect sizes for the principal and conditional phenotype \[[@pgen.1004926.ref018], [@pgen.1004926.ref019]\], meaning that when attention is restricted to SNPs with a given degree of association with the conditional phenotype, the p values for the principal phenotype will be falsely low; that is, the probability $Pr(P_{i} \leq p_{i} \mid P_{j} \leq p_{j},H_{0}^{(i)})$ will not in general be equal to *p* ~*i*~ ($= Pr(P_{i} \leq p_{i} \mid H_{0}^{(i)})$); in fact it will usually be higher.

When controls are shared, the distribution of p values for the principal phenotype given p values for the conditional phenotype depends on the underlying effect of each SNP on the conditional phenotype. For any given SNP, this underlying effect size, which we denote *η*, is not known. However, across all SNPs, *η* may be considered to be realisations of a random variable *H* whose distribution is mirrored by the distribution of observed effect sizes for the conditional phenotype. By integrating over this unknown true effect size for the conditional phenotype, allowance can be made for shared controls, and the 'expected quantile' of a p value for the principal phenotype, defined as $Pr(P_{i} \leq p_{i} \mid P_{j} \leq p_{j},H_{0}^{(i)})$, can be calculated, as detailed in the Methods section.

We assume that *H* has a mixture distribution defined by two parameters (*π* ~0~,*σ* ^2^), such that *H* = 0 with probability *π* ~0~ and *H* ∼ *N*(0,*σ* ~2~) with probability 1−*π* ~0~. The parameters (*π* ~0~,*σ*) are estimated from the observed distribution of effect sizes for the conditional phenotype. In order to show the effect of our p value adjustment, we simulated p values for 20,000 SNPs for a principal and conditional phenotype, with controls shared between simulated studies. All SNPs were null for the principal phenotype, and were variably null or non-null for the conditional phenotype with probability 0.9, 0.1 respectively. *Z* scores at non-null SNPs for the conditional phenotype were distributed as *N*(0,*σ* ^2^), as per our assumption. A value of *σ* = 3 was used, which was similar to the values of *σ* in real data estimated by our E-M algorithm.

We considered the set of simulated SNPs with p values for the conditional phenotype less than 0.05 ([Fig. 1](#pgen.1004926.g001){ref-type="fig"}). In the absence of shared controls, we expect the distribution of *p* ~*i*~ amongst this set to be uniform, and hence expect the black dots to lie along the x-y line. However, we see the principal p values are biased downward in this set (black dots, [Fig. 1](#pgen.1004926.g001){ref-type="fig"}). Our computed expected quantile (blue dots) agrees closely with the observed quantile. In a sense, this constitutes 'adjusting' the p values for the principal phenotype so that the expected distribution is uniform under the null hypothesis. Software to generate this simulation is available at <https://github.com/jamesliley/cFDR-common-controls>.

![Correction for shared controls.\
Simulation of GWAS summary statistics for 20000 SNPs, all null for phenotype *i* and variably null or non-null for phenotype *j*, with association tested using a shared control group. Black dots show p values for phenotype *i* at all SNPs with p value for phenotype *j* less than 0.05, with evident downward bias. Blue dots show our adjustment to expected quantile of p values. The red dots show the expected quantile we would compute if we were to assume incorrectly that all SNPs were null for the conditional phenotype *i*. We see that this quantity overestimates the true quantile.](pgen.1004926.g001){#pgen.1004926.g001}

Our formula can easily be adapted to arbitrary distributions of *H* at the cost of increased computational time, but the form of the distribution of *H* is not generally known. We show in [S1 Text](#pgen.1004926.s001){ref-type="supplementary-material"} (section A) that, even for distributions of *H* which differ markedly from our assumption of normality, the error in the estimate is not large, and generally translates to a negligible difference in the set of SNPs declared non-null using the $\hat{cFDR}$ method. While our assumption has the potential to be anti-conservative if *H* is bimodal, nonparametric estimates for distributions of effect sizes suggest they have a uni-modal distribution centred on zero\[[@pgen.1004926.ref030]\]. Reassuringly, our assumption is conservative if *H* has heavier tails than a normal.

Comparison to split control approach {#sec002c}
------------------------------------

We compared Andreasson's approach to SNP discovery which advocated splitting controls into non-overlapping subsets to our extended shared-control approach using a type 1 diabetes dataset with a total 12,175 cases and 15,171 controls. Controls and cases were each split into two sets (control sets had size 7,585 and 7,586, cases 6087 and 6088). 'Split' p values were computed using one set of controls and one set of cases and corresponding 'shared' p values were computed using the complete set of controls. As expected, more shared p values reached genome-wide significance than did split p values ([Fig. 2](#pgen.1004926.g002){ref-type="fig"}).

![Validation of the shared-control approach and the p-value adjustment due to shared controls.\
Panel A shows the effect of splitting controls on power to detect association. The number of SNPs with p values less than a given cutoff are shown for split-control and shared-control approaches. For all p-value cutoffs, fewer SNPs reach significance when using a split-control design. Panel B shows the number of SNPs with $\hat{cFDR}$ values less than a given cutoff using the existing method on a split-control design, our extended method on a shared control design with the adjustment for shared controls, or using the split-control approach naively on the shared-control design; that is, assuming incorrectly that $Pr(P_{i} \leq p_{i} \mid P_{j} \leq p_{j},H_{0}^{(i)}) = p_{i}$. The second figure shows that failing to correctly calculate $Pr(P_{i} \leq p_{i} \mid P_{j} \leq p_{j},H_{0}^{(i)})$ leads to a subtle increase in the number of SNPs declared non-null at all cutoffs, due to the incorrect underestimation of $\hat{cFDR}$.](pgen.1004926.g002){#pgen.1004926.g002}

We computed *cFDR* values by labelling one set of cases 'conditional' and the other 'principal' using the split-control p values using Andreasson's approach and using the shared control p values using our method. For reference, we compared these to a naive application of Andreasson's method on the shared-control p values ([Fig. 2B](#pgen.1004926.g002){ref-type="fig"}). More SNPs can be declared significant according to *cFDR* using the shared-control than split-control approach at all reasonable thresholds, and naive application of Andreasson's approach to shared-control p values again increases the number declared significant.

Because the quantity $Pr(P_{i} \leq p_{i} \mid P_{j} \leq p_{j},H_{0}^{(i)})$ is systematically underestimated when using this naive method (by assuming it is equal to *p* ~*i*~) as shown in [S1 Text](#pgen.1004926.s001){ref-type="supplementary-material"} (section C), it leads to a falsely low $\hat{cFDR}$. The increase in observed number of SNPs declared significant when using the naive method shows that it can indeed lead to false discoveries.

For principal phenotype p values in the range 5 × 10^−6^ - 5 × 10^−8^ - effectively the region from which 'new' SNPs may be discovered by cFDR rather than p value alone---the naive *cFDR* is frequently underestimated by 2--3 fold ([S1 Fig.](#pgen.1004926.s013){ref-type="supplementary-material"}, left panel). For lower p values, the naive *cFDR* may underestimate by hundreds- or thousand-fold, with the potential fold underestimation increasing with decreasing p value ([S1 Fig.](#pgen.1004926.s013){ref-type="supplementary-material"}, right panel). Because of the relatively high ratio of number of controls to number of cases, the correlation between effect sizes is lower in this constructed case (c. 0.22) than between most phenotypes in our study (c. 0.5). The underestimation of *cFDR* using the 'naive' method worsens with higher correlation, so we would expect that the fold-underestimate we see here is less severe than that which would be observed if applying this to other studies.

An upper bound on the false discovery rate of all declared SNPs {#sec002d}
---------------------------------------------------------------

An important property of our method is the control of the expected false discovery rate (FDR): the expected proportion of false positives among the SNPs found by our method. The p values at a SNP for the principal and conditional phenotype correspond to a point in the unit square. In this sense, we can define the expected FDR of a region *R* of the unit square as the ratio of the expected number of null SNPs whose p values are in *R* divided by the expected total number of SNPs whose p values are in *R*. From a result of Benjamini and Hochberg \[[@pgen.1004926.ref031]\], the expected FDR when *R* is a rectangle with vertices at (0,0),(*p* ~*i*~,0),(*p* ~*i*~,*p* ~*j*~),(0,*p* ~*j*~) is at most $\alpha = \hat{cFDR}(p_{i} \mid p_{j})$. If we denote by *L* the closed region defined by the set of p value pairs *p* ~*i*~,*p* ~*j*~ such that $\hat{cFDR}(p_{i} \mid p_{j}) \leq \alpha$, then *L* has the property that the FDR of any rectangle of this form contained within *L* is less than *α*.

However, the expected FDR over *L* is not necessarily bounded by *α* ([Fig. 3](#pgen.1004926.g003){ref-type="fig"}). This can be seen most easily in the extreme scenario in which all non-null SNPs are concentrated in the lower left corner of the unit square, and all null SNPs are also null for the conditional phenotype. In this case, the expected number of null SNPs in a rectangle is proportional to its area, so L is the union of all rectangles of the form above of a given area containing the lower-left corner of the unit square; that is, a hyperbola. Clearly the area of L is larger than the area of these constituent rectangles, yet it contains the same number of non-null SNPs, so it has a lower FDR.

![*L* is the locus of a set of points with $\hat{cFDR} = \alpha$.\
The FDR is the ratio of null SNPs to total SNPs in *L*. If all the non-null SNPs were concentrated in the lower left corner, then the number of non-null SNPs in *L* would be equal to that in any individual rectangle with vertices at the origin and on *L*, but the number of null SNPs would be greater, meaning that the expected FDR of all SNPs in *L* would be greater than *α*. *M* ^\*^ is the largest rectangle by area contained within *L*. The false discovery rate within *M* ^\*^ is less than *α* ^\*^, the value of $\hat{cFDR}$ at the upper right vertex, which is usually equal to *α*, as in this case. The FDR of *L* is bounded by *α* ^\*^ *v*(*L*)/*v*(*M*), where *v*(*L*) and *v*(*M* ^\*^) are the expected number of null SNPs in *L* and *M* ^\*^ respectively ([S1 Text](#pgen.1004926.s001){ref-type="supplementary-material"}, section B).](pgen.1004926.g003){#pgen.1004926.g003}

In the original method \[[@pgen.1004926.ref003]\], SNPs were declared significant if they were contained within any rectangular regions with a $\hat{cFDR}$ value of less than 0.01. Our reasoning demonstrates that the expected false-discovery rate amongst all such SNPs was higher than 0.01. We can derive an upper bound for the expected FDR of *L* by considering *M* ^\*^, the largest rectangle in *L*. We show in [S1 Text](#pgen.1004926.s001){ref-type="supplementary-material"} (section B) that the bound may be expressed simply as $\frac{v(L)}{v(M^{*})}\alpha^{*}$, where *v*() denotes the expected number of null SNPs contained within *L* or *M* ^\*^ (approximately the area of *L* and *M* ^\*^) and *α* ^\*^ is the cFDR at the upper-right vertex of *M* ^\*^ ([Fig. 3](#pgen.1004926.g003){ref-type="fig"}).

Application to ten immune mediated diseases {#sec002e}
-------------------------------------------

We obtained summary statistics in the form of p values for ten immune mediated diseases from ImmunoBase ([www.immunobase.org](http://www.immunobase.org), accessed 19/3/14). For each pair of diseases, the number of shared controls was estimated according to the description of the control samples in each paper. The numbers of cases, controls and our estimated numbers of shared controls for each study are shown in [Table 1](#pgen.1004926.t001){ref-type="table"}. Uniform quality control criteria were applied to all SNPs, and the MHC region, which exhibits both strong LD and strong association with immune mediated diseases was excluded. P values were corrected within each trait for genomic inflation using a standard algorithm \[[@pgen.1004926.ref032]\] applied to SNPs included on the ImmunoChip to replicate a GWAS study of reading and maths ability (Steve Eyre and Cathryn Lewis, personal communication), unlikely to be related to any immune mediated disease studied here.

10.1371/journal.pgen.1004926.t001

###### Study sizes.

![](pgen.1004926.t001){#pgen.1004926.t001g}

  Disease   Controls                                                          Cases   Estimated number of pairwise shared controls                                                                  
  --------- ----------------------------------------------------------------- ------- ---------------------------------------------- ---- ------ ------- ------ ------ ------ ------ ------ ------- -------
  T1D       \[[@pgen.1004926.ref020]\][\*](#t001fn001){ref-type="table-fn"}   12175   15171                                          \-   9364   12228   8430   4289   8514   4822   8430   4020    4020
  ATD       \[[@pgen.1004926.ref021]\]                                        9364    2733                                                \-     9364    8430   4289   8514   4822   8430   4020    4020
  CEL       \[[@pgen.1004926.ref022]\]                                        12228   12041                                                      \-      8430   4289   8514   4822   8430   4020    4020
  MS        \[[@pgen.1004926.ref023]\]                                        24091   14498                                                              \-     4289   8430   4822   8430   10102   10102
  NAR       \[[@pgen.1004926.ref024]\]                                        10421   1886                                                                      \-     4289   4289   4289   4020    4020
  PBC       \[[@pgen.1004926.ref025]\]                                        8514    2861                                                                             \-     4822   8430   4020    4020
  PSO       \[[@pgen.1004926.ref026]\]                                        22806   10588                                                                                   \-     4822   4020    4020
  RA        \[[@pgen.1004926.ref027]\]                                        15870   11475                                                                                          \-     4020    4020
  UC        \[[@pgen.1004926.ref028]\]                                        15977   10920                                                                                                 \-      15977
  CRO       \[[@pgen.1004926.ref028]\]                                        15977   14763                                                                                                         \-

Number of cases and controls for each study, and relevant references, together with our estimates of the number of controls shared between studies.

\*P values for T1D are from a meta-analysis of case-control and TDT data, with effective numbers of cases computed as shown in the methods section.

P values for each principal phenotype were adjusted to *p* ^′^ as described above in order to account for the effect of shared controls. For each ordered pair of phenotypes, a Q-Q plot was generated as per Andreasson et al \[[@pgen.1004926.ref003]\]. A Q-Q plot is a graph of the observed distribution of a random variable against the expected distribution. We overlaid Q-Q plots for *log* ~10~(*p* ^′^) values for the principal phenotype for subsets of SNPs exhibiting successively smaller p values for the conditional phenotype. [Fig. 4](#pgen.1004926.g004){ref-type="fig"} shows QQ plots for T1D conditional on RA and PSO; plots for all other pairwise comparisons may be found in [S4](#pgen.1004926.s016){ref-type="supplementary-material"}--[S13](#pgen.1004926.s025){ref-type="supplementary-material"} Figs. Notably, if lines shift further left with more stringent cutoffs on association with the conditional phenotype, then SNPs which are associated with the conditional phenotype are more likely to be associated with the principal phenotype, indicating pleiotropic effects of SNPs on the two phenotypes. In many cases, the Q-Q plots demonstrate considerable leftward shift with conditioning on association with a second disease, and we see strong evidence for pleiotropy for T1D conditioned on RA and little or no evidence for pleiotropy for T1D conditioned on PSO.

![Q-Q plots for T1D conditional on RA (Panel A) and PSO (Panel B).\
Y axes show $log_{10}(p_{T1D}^{\prime})$; X axes show log quantile (rank) of p values in various sets of SNPs. The degree of leftward shift of a black point from the diagonal is proportional to the unconditional FDR of that p value for the principal phenotype, and the degree of leftward shift of a coloured point is proportional to the conditional FDR of the p value for the principal phenotype and the p-cutoff corresponding to the colour for the conditional phenotype. As expected, a leftward shift is seen even for the unconditional Q-Q plots (black line) owing to the use of the ImmunoChip, which focuses on potential autoimmune-associated regions. Each colour corresponds to the Q-Q plot for *p* ~*T*1*D*~ amongst a subset of SNPs with *p* ~*RA*~ or *p* ~*PSO*~ less than the indicated cutoff. P values for T1D are adjusted for the effect of shared controls between studies. A leftward shift with decreasing *p* ~*RA*~ or *p* ~*PSO*~ cutoff indicates that SNPs which are associated with the conditional phenotype (RA or PSO) are more likely to be associated with the principal phenotype (T1D), presumably due to pleiotropic effects on phenotypes. Good enrichment is seen for T1D conditioning on RA (Panel A), and little or no enrichment conditioning on PSO (Panel B).](pgen.1004926.g004){#pgen.1004926.g004}

We estimated the unconditional and conditional false discovery rates, $\hat{uFDR}(p_{i})$ and $\hat{cFDR}(p_{i} \mid p_{j})$, at each SNP for each phenotype and each ordered pair of phenotypes respectively. [Fig. 5](#pgen.1004926.g005){ref-type="fig"} shows $\hat{cFDR}$ for T1D conditioned on RA. The advantage gained by $\hat{cFDR}$ can be seen in the left-shift of the region in which a SNP can be declared significant (blue dots), corresponding to a higher p-value cutoff for significance for T1D among SNPs with low p values for RA. Indeed, if only SNPs with a p value for RA less than some threshold *ζ* are considered, a p value cutoff for significance for T1D is given by the leftmost border of the blue dots on the line *P* ~*j*~ = *ζ*.

![$\hat{cFDR}$ (red-yellow) for T1D conditioned on RA.\
White dots signify the region for which $\hat{uFDR}$ is less than *α* corresponding to *p* \< 5 × 10^−8^. Blue dots signify the region for which cFDR is less than the same *α*. Note the leftward shift of blue points and the general leftward shift of colours corresponding to an increased p-value threshold for association with T1D for SNPs with low p values for RA. Black dots show a random sample of the observed p value pairs.](pgen.1004926.g005){#pgen.1004926.g005}

The degree of leftward shift in the Q-Q plots clearly contains information about the degree of pleiotropy between diseases. We defined a statistic summarizing some aspects of this evidence for pleiotropy and used it to visualise the set of pairwise relationships between diseases as a network ([Fig. 6](#pgen.1004926.g006){ref-type="fig"}). The network encouragingly reflects several pathophysiological associations: UC is linked to CRO, and T1D to ATD. Strong linkage is also seen both ways for MS and PBC, and between T1D and RA, findings which can also be seen in the Q-Q plots ([S4](#pgen.1004926.s016){ref-type="supplementary-material"}--[S13](#pgen.1004926.s025){ref-type="supplementary-material"} Figs.). One way relationships suggest the presence of a larger total number of associated SNPs for the disease at the start of the arrow than at the end.

![Network of degree of pleiotropy between phenotypes.\
An arrow runs from vertex *i* to vertex *j* if and only if by conditioning on *p* \< 5 × 10^−6^ for the conditional phenotype *j* we can increase the threshold for significance for the p value for the principal phenotype *i* from 5 × 10^−8^ to 4 × 10^−7^ or greater. Edges are thickened if the cutoff could be increased more than this. The threshold 4 × 10^−7^ was selected as the minimum value for which the network is weakly connected; that is, having an arrow to or from each edge.](pgen.1004926.g006){#pgen.1004926.g006}

Discovery of novel associations {#sec002f}
-------------------------------

The numbers of SNPs deemed significant for each phenotype by analysis using unconditional and conditional approaches are shown in [table 2](#pgen.1004926.t002){ref-type="table"}, with details in [S2](#pgen.1004926.s003){ref-type="supplementary-material"}--[S11](#pgen.1004926.s012){ref-type="supplementary-material"} Tables. $\hat{cFDR}$ allows certain SNPs with p values as high as 3 × 10^−6^ to be declared significant while controlling the false discovery rate at a relatively low value. Fifty-one of the 59 SNPs we identify uniquely through $\hat{cFDR}$ have previously been reported to be associated with the relevant disease through use of alternative significance thresholds, other genomic control procedures, other GWAS or additional samples not genotyped by ImmunoChip, a useful verification of our technique. Eight of the SNPs we discover uniquely through cFDR were in regions not previously known to be associated with the corresponding disease ([table 3](#pgen.1004926.t003){ref-type="table"}). These will require replication in independent samples to be declared truly associated, but they contain some potentially interesting signals, such as an association for RA at SNP rs72928038 near existing MS, ATD and T1D associations in *BACH2*, a transcriptional regulator involved in transcription repression and activation by *MAFK* \[[@pgen.1004926.ref033]\]

10.1371/journal.pgen.1004926.t002

###### Number of association signals found by unconditional and conditional methods.

![](pgen.1004926.t002){#pgen.1004926.t002g}

        Univariate   Conditional     Efficiency ratio                   
  ----- ------------ --------------- ------------------ --------------- ------
  T1D   44           3.79 × 10^−6^   4                  7.32 × 10^−6^   0.11
  ATD   3            1.04 × 10^−5^   4                  1.02 × 10^−5^   0.03
  CEL   46           3.97 × 10^−6^   7                  5.77 × 10^−6^   0.24
  MS    43           5.97 × 10^−6^   12                 1.19 × 10^−5^   0.15
  NAR   3            1.61 × 10^−4^   0                  3.05 × 10^−4^   0.83
  PBC   25           8.86 × 10^−6^   2                  2.09 × 10^−5^   0.12
  PSO   44           2.33 × 10^−5^   2                  5.45 × 10^−5^   0.06
  RA    16           2.04 × 10^−5^   10                 6.00 × 10^−5^   0.11
  UC    49           4.49 × 10^−6^   6                  6.73 × 10^−6^   0.14
  CRO   80           2.08 × 10^−6^   12                 2.73 × 10^−6^   0.28

Number of SNPs with *p* ≤ 5 × 10^−8^ after genomic found by analysis of the principal phenotype alone and the estimated equivalent FDR for this set of SNPs. Conditional analysis shows the number of additional SNPs found through conditional FDR analysis, and the upper bound for the FDR of all SNPs selected by this method, including adjustment for the multiple phenotypes conditioned upon. Finally, we summarise the performance of the cFDR approach by the FDR ratio---the ratio of the upper bound for the FDR for all SNPs selected by the cFDR approach, and the estimated FDR of all SNPs with *p* less than the maximum *p* value within this set if we had not conditioned. Note, we list only the most associated SNPs in each LD block by pruning according to LD, as described in Methods.

10.1371/journal.pgen.1004926.t003

###### Novel SNP-disease associations.

![](pgen.1004926.t003){#pgen.1004926.t003g}

  Chr        Pos         SNP                  Disease   p value         C. phen.   MAF     Nearby Genes
  ---------- ----------- -------------------- --------- --------------- ---------- ------- ------------------
  1p31.3     61791863    rs6691768: A \> G    CEL       8.56 × 10^−8^   CRO        0.373   *NFIA*
  1p13.1     117076399   rs1034920: T \> C    PBC       1.43 × 10^−6^   MS         0.100   *CD58*
  2p21       43359275    rs6705577: G \> C    CEL       3.61 × 10^−7^   MS         0.272   *THADA*
  2q32.1     185501065   rs79248157: T \> C   UC        1.53 × 10^−7^   CRO        0.031   *ZNF804A*
  6q15       90976768    rs72928038: G \> A   RA        5.89 × 10^−7^   T1D        0.176   *BACH2*
  16q12.1    51080214    rs12924003:C \> T    CRO       5.82 × 10^−8^   NAR        0.321   *NOD2*
  20q13.12   44596207    rs6032606:G \> C     CEL       1.31 × 10^−7^   MS         0.052   *ZNF355*, *MMP9*
  22q13      37633851    rs9610686:G \> A     CEL       1.90 × 10^−7^   ATD        0.387   *RAC2*

Eight SNP associations with the indicated disease discovered by cFDR but not previously published to our knowledge, together with the phenotype upon which they have been conditioned (C. phen.) and nearby genes. SNPs are shown by RSID: major\>minor alleles. The disease (prinicipal phenotype) p value has been corrected for genomic inflation. Note that a SNP reaching significance by $\hat{cFDR}$ for the principal phenotype does not constitute evidence of association with the conditional phenotype (C. phen.). Chr = chromosome, Pos = position (build GRCh37), SNPs are shown with major\>minor alleles, MAF = minor allele frequency in UK controls.

The SNP rs1034290 in region 1p13.1, which we found to be associated with PBC, is in intron three of *CD58*, which is a surface receptor involved in binding and activation of T-lymphocytes. The protective effect of the MS-associated allele is postulated to arise from upregulation of the transcription factor *FOXP3* \[[@pgen.1004926.ref034]\] and the patterns of association in the region suggest the two diseases may share a causal variant here (<http://www.immunobase.org>).

Discussion {#sec003}
==========

We have extended a technique for computing conditional Bayesian False Discovery Rates to GWAS for independent diseases with shared control groups. This technique enables improved detection of disease-associated SNPs compared to conventional methods. By enabling larger control groups for each study, our method uses data more efficiently than in corresponding study designs in which control groups are independent, and is applicable to a wider range of GWAS datasets for which only summary statistics are available.

Combination of GWAS by analysis of pleiotropy in this sense has several attractive advantages over single-phenotype analysis. The most obvious advantage is improved detection of disease-associated SNPs using GWAS without the need for additional samples. A secondary advantage arises from understanding of the pleiotropic structure between phenotypes: if a SNP is known to exhibit pleiotropy between two conditions, it may be causative for a shared risk factor or pre-disease state. Analysis of such SNPs has the potential to yield information on disease aetiology, with implications for preventative medicine and development of treatment.

A further potential use for this technique could be the genomic analysis of diseases with complex phenotypes. In many cases, distinction between two diseases may be difficult; for instance, Crohn's disease and Ulcerative Colitis \[[@pgen.1004926.ref035]\]. Additionally, many diseases, including narcolepsy (<http://www.uptodate.com/contents/clinical-features-and-diagnosis-of-narcolepsy>, accessed 20/6/14), are definitively diagnosed on clinical grounds. This implies that these diseases may constitute a range of biochemical and genetic states. Inclusion criteria based on objective biochemical grounds, such as that used for narcolepsy in the context of this paper \[[@pgen.1004926.ref024]\] are unlikely to characterise all patients with these diseases, and conclusions drawn from studies will not necessarily be medically applicable to the whole patient population. Given this, diseases defined phenotypically with potential genomic diversity may be better analysed by separate consideration of biochemically-defined subtypes, with a collective analysis performed by a method such as $\hat{cFDR}$, avoiding the assumption that the genomic bases of disease subtypes are identical.

We identify a counter intuitive property that the FDR in the union of all regions with $\hat{cFDR}$ less than a given *α* may be greater than *α*, and propose a method to overcome this problem. Our methods for adjusting cutoffs to control FDR and account for multiple testing demonstrate the geometrical elegance of the theory of these techniques, with the possibility for further improvements and understanding. They are complex to apply, but could be much simplified if interest was directed to SNPs with conditional p values less than some threshold *p* ~0~. Our method would ensure that the expected false discovery rate at SNPs with $\hat{cFDR}(p_{i} \mid p_{0}) \leq \alpha$ would indeed be controlled at *α*. Our more complicated method to control FDR is necessary if the variable *p* ~*j*~ is used in place of the constant *p* ~0~.

An important consideration in both our method and the original Andreasson method is that a $\hat{cFDR}(p_{i} \mid p_{j})$ value which reaches significance does not constitute genome-wide evidence of association with the conditional phenotype *j*; indeed, the probability of association with the conditional phenotype relates to $\hat{cFDR}(p_{j} \mid p_{i})$ and in general $\hat{cFDR}(p_{i} \mid p_{j}) \neq \hat{cFDR}(p_{j} \mid p_{i})$. In some cases, where the principal p value is very close to genome-wide significance, even conditioning on *p* ~*j*~ ≤ 0.5 can theoretically be enough to reach the relevant $\hat{cFDR}$ threshold. This is not a weakness of the $\hat{cFDR}$ method as such, but a consequence of using a discrete technique (a significance cutoff) on a variable which essentially continuous in two dimensions ($\hat{cFDR}$). Principal p values greater than 5 × 10^−8^ which can be declared significant conditioning on large conditional p value cutoffs correspond to an increase in the area of the region *L* (see [results](#sec002){ref-type="sec"} section), which is accounted for by our FDR-controlling method.

Our method enables improved detection of SNPs compared to analysis of unconditional FDR (principal p value alone). However, the improvement is smaller than that reported by Andreasson et al \[[@pgen.1004926.ref003], [@pgen.1004926.ref013], [@pgen.1004926.ref014]\], who detected almost twice as many SNPs using *cFDR* as they would have detected with *uFDR*. This is expected for two reasons. Firstly, the gain in power from *cFDR* essentially comes from an increase in the total number of controls and the effective number of cases. If controls are shared, the only information gain can come from increasing the number of effective cases. Consequently, the difference in power between *cFDR* and *uFDR* will not be as large when controls are shared, although both outperform their counterparts when controls are split. Secondly, we were careful to use stringent cutoffs for FDR which were chosen to mirror the established genomewide significance threshold of *p* ≤ 5 × 10^−8^, generally equivalent to a false discovery rate around 5 × 10^−6^ to 5 × 10^−5^, compared to Andreasson et al who declared non-null all SNPs with $\hat{cFDR} < 0.01$.

One alternative way to exploit pleiotropic relationships is by meta-analysing two related diseases together, as though the diseases were the same. Our method confers several advantages over this approach. The most important of these is that our method borrows strength from other SNPs according to the level of genome wide pleiotropy between diseases; that is, if the two GWAS suggest extensive pleiotropy (such as [Fig. 4](#pgen.1004926.g004){ref-type="fig"} for T1D --- RA), a low p value for a conditional phenotype will 'sway' our judgement of association with the principal phenotype more than the same p value for a conditional phenotype with poor pleiotropy (such as [Fig. 4](#pgen.1004926.g004){ref-type="fig"}, for T1D --- PSO). A meta-analysis would not distinguish these two scenarios. A secondary advantage of our technique is that SNP detection is not systematically weakened if the two diseases do not exhibit pleiotropy, as would be the case in meta-analysis; this arises because we are testing association with only one of the two phenotypes at a time.

Methods {#sec004}
=======

Ethics statement {#sec004a}
----------------

This paper re-analyses previously published datasets. All patient data were handled in accordance with the policies and procedures of the participating organisations.

Datasets {#sec004b}
--------

We obtained SNP summary statistics from ten studies on autoimmune diseases from ImmunoBase ([www.immunobase.org](http://www.immunobase.org)). Inclusion and exclusion criteria for the studies are described in detail in the original publications (\[[@pgen.1004926.ref020]--[@pgen.1004926.ref028], [@pgen.1004926.ref028]\]. Generally, some or all controls from different studies were obtained from common data sources, resulting in overlapping control groups. All studies used the ImmunoChip array \[[@pgen.1004926.ref017]\].

P values for type 1 diabetes were from a meta-analysis of a case-control study and familial study using the transmission disequilibrium test (TDT). In order to calculate the correlation between *p* values for different diseases, we needed to calculate effective numbers of cases and controls for the combined T1D study. For a case control study, under the assumptions of Hardy-Weinberg and the null hypothesis, the variance of the log odds ratio may be expressed as $$\frac{n_{0} + n_{1}}{n_{0}n_{1}}\frac{1}{f\left( 1 - f \right)}$$ where *n* ~0~ and *n* ~1~ are the numbers of cases and controls and *f* is the minor allele frequency in controls.

Given the standard error of a log OR for the TDT study, *σ̂*, and a minor allele frequency, we estimated *M* = *σ̂* ^2^ *f*(1−*f*) for all ImmunoChip SNPs which did not show deviation from the null hypothesis (*p* \> 0.5). The distribution of log(*M*) is shown in [S3 Fig](#pgen.1004926.s015){ref-type="supplementary-material"}. By equating the median of *M* with $\frac{n_{0} + n_{1}}{n_{0}n_{1}}$, and assuming that each TDT family contributed the equivalent information to one control in a case-control study, ie *n* ~0~ = 2943, we estimated an equivalent number of cases as 4126. This seemed reasonable, given that there are a total of 5505 (dependent) cases across those families.

SNPs were excluded on the basis of QC summaries calculated on 12,888 common controls: call rate less than 99%, minor allele frequency less than 0.02, or deviation from Hardy-Weinberg equilibrium (∣*Z*∣ \> 5). Given the strong association of immune mediated diseases with the MHC and the extended LD in the region, we were concerned that MHC SNPs might cause inaccurate estimation of pleiotropy. We therefore excluded SNPs in a wide band around the MHC region on chromosome 6 (co-ordinates 24500000: 34800000, build NCBI36). After quality control, genotype data was available for at least one phenotype at a total of 110677 SNPs.

Genomic control {#sec004c}
---------------

P values were corrected for genomic inflation using a genomic control algorithm \[[@pgen.1004926.ref032]\]. A set of SNPs known to be unassociated with autoimmune disease was obtained from the Wellcome Trust Case Control Consortium (WTCCC) study on reading and mathematics ability. These SNPs were pruned so that none were in LD with *r* ^2^ \> 0.2, and any SNPs within 500 kb of known autoimmune-associated regions were removed. The average degree of inflation was computed for each disease at the remaining 1761 SNPs, and all effect sizes and p values were adjusted accordingly.

Procedures for uFDR and cFDR {#sec004d}
----------------------------

We assume that the p-values for a phenotype *i* across all SNPs are instances of a random variable *P* ~*i*~. If *p* ~*i*~ is an instance of this random variable corresponding to a SNP of interest, the unconditional false discovery rate *uFDR*(*p* ~*i*~) is defined as $$\begin{matrix}
\left. uFDR(p_{i}) = Pr(H_{0}^{(i)} \middle| P_{i} \leq p_{i}) \right. \\
{= Pr\left( H_{0}^{(i)} \right)\frac{\left. Pr(P_{i} \leq p_{i} \middle| H_{0}^{(i)}) \right.}{Pr(P_{i} \leq p_{i})}} \\
{= Pr\left( H_{0}^{(i)} \right)\frac{p_{i}}{Pr(P_{i} \leq p_{i})}} \\
\end{matrix}$$ where $H_{0}^{(i)}$ is the null hypothesis that the SNP of interest is not associated with phenotype *i*. Given a set of observed p values $\{ p_{i}^{1},p_{i}^{2}\ldots p_{i}^{N}\}$ for a phenotype *i* at *N* different SNPs, and an observed p value *p* ~*i*~ for a SNP of interest, we estimate this quantity as $$\begin{array}{ccl}
{\hat{uFDR}\#\left( p_{i} \right)} & = & \frac{p_{i}}{\#\left( \text{p}\;\text{values}\; p_{i}^{k}\;\text{with}{p_{i}^{k} \leq p_{i}} \right)/N} \\
 & = & \frac{{Expected}\;{quantile}\;{of}\; p_{i}\;{under}\; H_{0}^{(i)}}{{Observed}\;{quantile}\;{of}\; p_{i}} \\
\end{array}$$ Because we make the approximation $Pr(H_{0}^{i}) = 1$, the estimate $\hat{uFDR}$ is a an upwards-biased estimate of *uFDR*; that is, its expected value is greater than the true *uFDR*, making it a conservative estimator.

We compute the quantity ([1](#pgen.1004926.e050){ref-type="disp-formula"}) for each SNP at each phenotype, declaring any SNP for which $\hat{uFDR}(p_{i}) \leq \alpha$ as non-null for phenotype *i*. Defining *V* as the number of SNPs falsely declared non-null, *R* as the total number of SNPs declared non-null, and *Q* = *V*/*R*, a theorem of Benjamini and Hochberg \[[@pgen.1004926.ref031]\] shows the expected false discovery rate *E*(*Q*) among SNPs with $\hat{uFDR} \leq \alpha$ is less than *α*.

The *cFDR* constitutes a natural extension of this idea. We assume that the p-values for two phenotypes *i* and *j* across all SNPs are instances of a pair of random variables *P* ~*i*~,*P* ~*j*~. If *p* ~*i*~ and *p* ~*j*~ are instances of these variables corresponding to a SNP of interest then the conditional false discovery rate *cFDR* is defined for the set of SNPs with p values for each phenotype less than or equal to those at this SNP (as per Andreasson et al \[[@pgen.1004926.ref003]\]) as $$\begin{array}{ccl}
{cFDR\left( p_{i} \middle| p_{j} \right)} & = & {Pr\left( H_{0}^{(i)} \middle| P_{i} \leq p_{i},P_{j} \leq p_{j} \right)} \\
 & = & {Pr\left( H_{0}^{(i)} \middle| P_{j} \leq p_{j} \right)\frac{Pr\left( P_{i} \leq p_{i} \middle| P_{j} \leq p_{j},H_{0}^{(i)} \right)}{Pr\left( P_{i} \leq p_{i} \middle| P_{j} \leq p_{j} \right)}} \\
\end{array}$$ The estimation of this quantity proceeds in a similar way to *uFDR*. Given a set of observed p value pairs $\{(p_{i}^{1},p_{j}^{1}),(p_{i}^{2},p_{j}^{2})\ldots(p_{i}^{N},p_{j}^{N})\}$ for two phenotypes *i* and *j* at *N* different SNPs, and an observed p value pair (*p* ~*i*~, *p* ~*j*~) for a SNP of interest, we define *N* ~1~ as the number of p value pairs with *P* ~*j*~ ≤ *p* ~*j*~, and estimate the *cFDR* as $$\begin{array}{ccl}
{\hat{cFDR}\left( p_{i} \middle| p_{j} \right)} & = & \frac{Pr\left( P_{i} \leq p_{i} \middle| P_{j} \leq p_{j},H_{0}^{(i)} \right)}{\#\left( \text{pairs;}{\left( p_{i}^{k},p_{j}^{k} \right) \in \left( P_{i},P_{j} \right)}\;\text{with}\;{p_{i}^{k} \leq p_{i}}\;\text{and}\;{p_{j}^{k} \leq p_{j}} \right)/N_{1}} \\
 & = & \frac{{Expected}\;{quantile}\;{of}\; p_{i}\;{under}\; H_{0}^{(i)}\;{amongst}\; p_{i}^{k}\;{with}\; k\;{satisfying}\; p_{j}^{k} \leq p_{j}}{{Observed}\;{quantile}\;{of}\; p_{i}\;{amongst}\; p_{i}^{k}\;{with}\; k\;{satisfying}\; p_{j}^{k} \leq p_{j}} \\
\end{array}$$ Again, this estimate is conservative, due to the approximation $Pr(H_{0}^{(i)} \mid P_{j} \leq p_{j}) = 1$.

We compute the quantity ([2](#pgen.1004926.e057){ref-type="disp-formula"}) for each SNP at each pair of phenotypes, declaring any SNP for which $\hat{cFDR}(p_{i} \mid p_{j}) \leq \alpha^{\prime}$ as non-null for phenotype *i*. However, as noted earlier, this does not guarantee that the expected false discovery rate amongst such SNPs is less than *α* ^′^. We show that the FDR is controlled at a higher level dependent on the region of the unit square defined by rectangles for which $\hat{cFDR}(x \mid y) \leq \alpha^{\prime}$.

Our method here diverges from the original method proposed by Andreasson et al, in the use of the expected quantile $Pr(P_{i} \leq p_{i} \mid P_{j} \leq p_{j},H_{0}^{(i)})$ in place of the p-value *p* ~*i*~. If studies share no controls, it can be reasonably assumed that, for a SNP which is null for phenotype *i*, the p values (*p* ~*i*~,*p* ~*j*~) are independent, so $p_{i}^{\prime} = Pr(P_{i} \leq p_{i} \mid P_{j} \leq p_{j},H_{0}^{(i)}) = p_{i}$. This is the approach taken by Andreasson et al \[[@pgen.1004926.ref003]\]. We propose a method for computing $p_{i}^{\prime}$ when controls are shared between studies, and the independence assumption above is not valid.

Our approach is to compute the related quantity $Pr(P_{i} \leq p_{i} \mid P_{j} \leq p_{j},H_{0}^{(i)},H_{\eta}^{(j)})$, where *η* is the (unobserved) effect size we would observe for a given SNP for phenotype *j* if the observed MAFs agreed exactly with the population MAFs for that SNP, and $H_{\eta}^{(j)}$ is the hypothesis that *Z* ~*j*~ ∼ *N*(*η*,1) for that SNP. This quantity can be thought of as the 'expected quantile' of *p* ~*i*~; that is, the proportion of p values we expect to be less than *p* ~*i*~.

Computation of expected quantile {#sec004e}
--------------------------------

From the first part of ([2](#pgen.1004926.e057){ref-type="disp-formula"}), we have: $$cFDR\left( p_{i} \middle| p_{j} \right) = \frac{Pr\left( H_{0}^{(i)} \middle| P_{j} \leq p_{j} \right)Pr\left( P_{i} \leq p_{i} \middle| P_{j} \leq p_{j},H_{0}^{(i)} \right)}{Pr\left( P_{i} \leq p_{i} \middle| P_{j} \leq p_{j} \right)}$$ As per Andreasson et al \[[@pgen.1004926.ref003]\], the quantity $Pr(H_{0}^{(i)} \mid P_{j} \leq p_{j})$ is set conservatively at 1, and the quantity *Pr*(*P* ~*i*~ ≤ *p* ~*i*~∣*P* ~*j*~ ≤ *p* ~*j*~) is estimated empirically as the proportion of pairs of observed p values $({p^{\prime}}_{i},\mspace{180mu}{p^{\prime}}_{j})$ with ${p^{\prime}}_{j}\mspace{180mu} \leq \mspace{180mu} p_{j}$ which also satisfy ${p^{\prime}}_{i}\mspace{180mu} \leq \mspace{180mu} p_{i}$.

For a given SNP, let *η* denote the standardised mean allele frequency (MAF) difference; that is, the *Z* value we would compute if the observed MAFs agreed exactly with the population MAFs. We consider *η* for a random SNP as being an instance of a random variable *H*, and that the observed *z* value for that SNP *Z*∣*H* = *η* is distributed as $$\left. Z \middle| H = \eta \sim N\left( \eta,1 \right) \right.$$ We further assume that *H* follows a mixture distribution taking the value 0 with probability $\pi_{0}^{(j)}$ and a normal *pdf* with probability $1 - \pi_{0}^{(j)}$: $$H \sim \begin{cases}
{0,} & {p = \pi_{0}^{(j)}} \\
{N\left( 0,\sigma^{2} \right),} & {p = 1 - \pi_{0}^{(j)}} \\
\end{cases}$$

This implies $$Z_{j} \sim \begin{cases}
{N\left( 0,1 \right),} & {p = \pi_{0}^{(j)}} \\
{N\left( 0,1 + \sigma^{2} \right),} & {p = 1 - \pi_{0}^{(j)}} \\
\end{cases}$$ Thus, given the observed distribution of *Z* ~*j*~, the parameters $\pi_{0}^{(j)}$ and *σ* ~*j*~ may be estimated by an expectation - maximisation algorithm (<https://gist.github.com/chr1swallace/11421212>).

We assume as per Zaykin \[[@pgen.1004926.ref018]\] that the distribution of pairs of observed z values (*Z* ~*i*~,*Z* ~*j*~) for a single given SNP is bivariate normal. Denote by $H_{\eta}^{(j)}$ the event that, for a given SNP, the values *Z* ~*j*~ are distributed as *N*(*η*,1), with *η* depending on the SNP.

Under our assumption of the null hypothesis $H_{0}^{(i)}$ for the principal phenotype and a population MAF difference corresponding to *η* for the conditional phenotype, we have $$\left( Z_{i},Z_{j} \middle| H_{0}^{(i)},H_{\eta}^{(j)} \right) \sim N\left( {\begin{pmatrix}
0 \\
\eta \\
\end{pmatrix},\begin{pmatrix}
1 & \rho \\
\rho & 1 \\
\end{pmatrix}} \right)$$

The correlation *ρ* arises from the shared controls between groups \[[@pgen.1004926.ref018], [@pgen.1004926.ref019]\] and is asymptotically equal to $$\rho = \frac{1}{\sqrt{\left( {1 + \frac{N_{0i}}{N_{0}}} \right)\left( {1 + \frac{N_{0j}}{N_{0}}} \right)\left( {1 + \frac{N_{0i}}{N_{i}} + \frac{N_{0}}{N_{i}}} \right)\left( {1 + \frac{N_{0j}}{N_{j}} + \frac{N_{0}}{N_{j}}} \right)}}$$ where *N* ~*i*~ and *N* ~*j*~ are the numbers of cases, *N* ~0*i*~ and *N* ~0*j*~ are the numbers of non-shared controls, and *N* ~0~ is the number of shared controls for the original GWAS for the principal and conditional phenotypes respectively. There is good agreement with the asymptotic correlation when group sizes are greater than 100 \[[@pgen.1004926.ref018]\].

Given equations ([5](#pgen.1004926.e074){ref-type="disp-formula"})--([8](#pgen.1004926.e080){ref-type="disp-formula"}), the joint distribution of *Z* ~*i*~ and *Z* ~*j*~ can be computed under only the assumption $H_{0}^{(i)}$. The value of the partial PDF of $(Z_{i},Z_{j} \mid H_{0}^{(i)})$ at (*x*,*y*) can be derived in a similar way to ([6](#pgen.1004926.e075){ref-type="disp-formula"}): $$\left. \left( Z_{i},Z_{j} \middle| H_{0}^{(i)} \right) \right.\sim\begin{cases}
{N\left( {(\begin{array}{l}
0 \\
0 \\
\end{array}),(\begin{array}{ll}
1 & \rho \\
\rho & 1 \\
\end{array})} \right),} & {p = \pi_{0}^{(j)}} \\
{N\left( {(\begin{array}{l}
0 \\
0 \\
\end{array}),(\begin{array}{ll}
1 & \rho \\
\rho & {1 + \sigma^{2}} \\
\end{array})} \right),} & {p = 1 - \pi_{0}^{(j)}} \\
\end{cases}$$

We now compute the final probability in equation ([3](#pgen.1004926.e066){ref-type="disp-formula"}). Define $$P_{\eta}\left( X \right) = Pr\left( X \middle| H_{0}^{(i)},H_{\eta}^{(j)} \right)$$ as the probability of observing events *X* for a particular SNP with true effect size *η* (which may be 0, corresponding to the general null). Then, $$Pr\left( P_{i} \leq p_{i} \middle| P_{j} \leq p_{j},H_{0}^{(i)} \right) = \frac{Pr\left( P_{i} \leq p_{i},P_{j} \leq p_{j} \middle| H_{0}^{(i)} \right)}{Pr\left( P_{j} \leq p_{j} \middle| H_{0}^{(i)} \right)}$$ $$= \frac{\pi_{0}^{(j)}P_{0}\left( P_{i} \leq p_{i},P_{j} \leq p_{j} \right) + \left( 1 - \pi_{0}^{(j)} \right)\int_{- \infty}^{\infty}P_{\eta}\left( P_{i} \leq p_{i},P_{j} \leq p_{j} \right)f\left( \eta \right)d\eta}{\pi_{0}^{(j)}P_{0}\left( P_{j} \leq p_{j} \right) + \left( 1 - \pi_{0}^{(j)} \right)\int_{- \infty}^{\infty}P_{\eta}\left( P_{j} \leq p_{j} \right)f\left( \eta \right)d\eta}.$$

If the distribution of *H* is estimable by other means, quantity ([11](#pgen.1004926.e086){ref-type="disp-formula"}) can be calculated numerically without the assumption that the non-null component of *H* be normally distributed, at the cost of higher computation time. Under our assumptions, equations ([6](#pgen.1004926.e075){ref-type="disp-formula"}) and ([9](#pgen.1004926.e083){ref-type="disp-formula"}) enable the fast computation of quantity ([11](#pgen.1004926.e086){ref-type="disp-formula"}) by normal CDFs; writing $$\begin{array}{l}
{\Lambda_{(\rho,\sigma^{2})}(z_{i},z_{j})\mspace{180mu} = \mspace{180mu}{\int_{|x| > |z_{i}|,|y| > |z_{j}|}N_{(\begin{array}{l}
0 \\
0 \\
\end{array}),(\begin{array}{ll}
1 & \rho \\
\rho & {1 + \sigma^{2}} \\
\end{array})}}(x,y)dxdy} \\
{\lambda_{\sigma^{2}}(z_{j})\mspace{180mu} = \mspace{180mu}{\int_{|y| > |z_{j}|}N_{(0,1 + \sigma^{2})}}(y)dy} \\
\end{array}$$ we have $$\begin{array}{r}
{Pr\left( P_{i} \leq p_{i} \middle| P_{j} \leq p_{j},H_{0}^{(i)} \right) = \frac{\pi_{0}^{(j)}\Lambda_{(\rho,0)}\left( z_{i},z_{j} \right) + \left( 1 - \pi_{0}^{(j)} \right)\Lambda_{(\rho,\sigma^{2})}\left( z_{i},z_{j} \right)}{\pi_{0}^{(j)}\lambda_{0}\left( z_{j} \right) + \left( 1 - \pi_{0}^{(j)} \right)\lambda_{\sigma^{2}}\left( z_{j} \right)}} \\
\end{array}$$

Point expected quantile {#sec004f}
-----------------------

Because the formula for $Pr(P_{i} \leq p_{i} \mid P_{j} \leq p_{j},H_{0}^{(i)})$ is differentiable on the unit square, an expression for the expected quantile of *p* ~*i*~ given an exact value for *p* ~*j*~ can be computed by taking the partial derivative with respect to *p* ~*j*~: $$Pr\left( P_{i} \leq p_{i} \middle| P_{j} = p_{j},H_{0}^{(i)} \right) = \frac{\partial}{\partial p_{j}}Pr\left( P_{i} \leq p_{i} \middle| P_{j} \leq p_{j},H_{0}^{(i)} \right) = \frac{A + B}{C}$$ where $$\begin{array}{l}
{A\mspace{180mu} = \pi_{0}^{(j)}N_{(0,1)}(z_{j}){\int_{|x| \geq z_{i}}N_{(\rho z_{j},1 - \rho^{2})}}(x)dx\mspace{180mu}} \\
{B\mspace{180mu} = (1 - \pi_{0}^{(j)})N_{(0,1 + \sigma^{2})}(z_{j}){\int_{|x| \geq z_{i}}N_{(\frac{\rho z_{j}}{1 + \sigma^{2}},\frac{1 - \rho^{2} + \sigma^{2}}{1 + \rho^{2}})}}(x)dx} \\
{C\mspace{180mu} = \pi_{0}^{(j)}N_{(0,1)}(z_{j}) + (1 - \pi_{0}^{(j)})N_{(0,1 + \sigma^{2})}(z_{j})} \\
\end{array}$$ where *N* ~*μ*,*σ*^2^~(*x*) denotes the value of the normal *pdf* with mean *μ* and variance *σ* ^2^ at *x*.

Significance thresholds {#sec004g}
-----------------------

Because $\hat{uFDR}$ values are monotonically related to p values, the widely accepted GWAS p value cutoff of 5 × 10^−8^ corresponds naturally to a cutoff for $\hat{uFDR}$. For each phenotype *i*, we set a significance threshold *β* ^*i*^ for $\hat{uFDR}(p_{i})$ as the lowest possible value of *γ* for which $\left. \hat{uFDR}(p_{i}) \leq \gamma\Leftrightarrow p_{i} \leq 5 \times 10^{- 8} \right.$.

We then applied an analagous approach to $\hat{cFDR}$. For each pair of phenotypes (*i*,*j*), we set a significance threshold $\alpha_{j}^{i}$ as the lowest possible value of *γ* ^′^ for which $\left. \hat{cFDR}(p_{i} \mid p_{j}) \leq \gamma^{\prime}\Leftrightarrow p_{i} \leq 5 \times 10^{- 8} \right.$. Given the distribution of *P* ~*j*~, it is possible that this could lead to declaring SNPs with *p* ~*i*~ \> 5 × 10^−8^,*p* ~*j*~ ≈ 1 as significant. To avoid this, if $\alpha_{j}^{i}$ was larger (less stringent) than *β* ^*i*^, we set $\alpha_{j}^{i} = \beta^{i}$.

For each ordered pair of phenotypes (*i*,*j*), we declared all SNPs with $\hat{cFDR}(p_{i} \mid p_{j}) \leq \alpha_{j}^{i}$ as non-null for phenotype *i*. This included all SNPs with $\hat{uFDR}(p_{i}) \leq \beta^{i}$. We then used a technique described in [S1 Text](#pgen.1004926.s001){ref-type="supplementary-material"} (section B) to compute upper bounds $c_{j}^{(i)}$ on the false discovery rate amongst SNPs for which $\hat{cFDR}(p_{i},p_{j}) \leq \alpha_{j}^{i}$. For each phenotype, this gave nine upper bounds, corresponding to each of the nine conditional phenotypes.

Network and heatmap representation of pleiotropy {#sec004h}
------------------------------------------------

We compared the degree of pleiotropy between diseases by considering how much the p-value threshold for significance for the principal phenotype changed when conditioning on a small p-value threshold for the conditional phenotype. We used the $\hat{cFDR}$ algorithm to compute the number $p_{i}^{j*}$ such that $P(H_{0}^{(i)} \mid P_{i} \leq p_{i}^{j*},p_{j} \leq 5 \times 10^{- 6}) = P(H_{0}^{(i)} \mid P_{i} \leq 5 \times 10^{- 8})$; that is, $\hat{cFDR}(p_{i}^{j*} \mid 5 \times 10^{- 6}) = \hat{uFDR}(5 \times 10^{- 8})$. We then considered the ratio $p_{i}^{j*}/5 \times 10^{- 8}$; that is, the fold increase in significance cutoff after conditioning.

We note that because of the fixed value of *p* ~*j*~ = 5 × 10^−6^, the expected false discovery rate amongst the set of SNPs which satisfy $P(H_{0}^{(i)} \mid P_{i} \leq p_{i}^{j*},p_{j} \leq 5 \times 10^{- 6}) \leq P(H_{0}^{(i)} \mid P_{i} \leq 5 \times 10^{- 8})$ is bounded above by $P(H_{0}^{(i)} \mid P_{i} \leq 5 \times 10^{- 8})$, by the Benjamini-Hochberg result. Thus the expected false discovery rate amongst SNPs with $P_{i} \leq p_{i}^{j*}$ and *P* ~*j*~ ≤ 5 × 10^−6^ is bounded by the same value as the expected false discovery rate amongst SNPs with *P* ~*i*~ ≤ 5 × 10^−8^.

We visualised the ratio $p_{i}^{j*}/5 \times 10^{- 8}$ as a heatmap ([S2 Fig.](#pgen.1004926.s014){ref-type="supplementary-material"}). We also produced a network ([Fig. 6](#pgen.1004926.g006){ref-type="fig"}), with an edge from vertex *i* to vertex *j* if and only if, by conditioning on *P* ~*j*~ ≤ 5 × 10^−6^, the cutoff for significance for *P* ~*i*~ could be increased from 5 × 10^−8^ to 4 × 10^−7^. This cutoff was chosen as the smallest value such that the network was weakly connected; that is, each vertex had an arrow either to it or from it.

Discovery of novel SNP associations {#sec004i}
-----------------------------------

Cutoffs *β* ^*i*^ were chosen such that $\left. \hat{uFDR}(p_{i}) \leq \beta^{i}\Leftrightarrow p_{i} \leq 5 \times 10^{- 8} \right.$. SNPs were deemed significant for a principal phenotype *i* if $\hat{cFDR}(p_{i} \mid p_{j}) \leq \alpha_{j}^{i}$ for any conditional phenotype *j* and $\alpha_{j}^{i} \leq \beta^{i}$.

The expected false discovery rate amongst SNPs for which $\hat{uFDR}(p_{i}) \leq \beta^{i}$ is less than *β* ^*i*^ due to a theorem of Benjamini and Hochberg. However, as discussed above and in [S1 Text](#pgen.1004926.s001){ref-type="supplementary-material"} (section B), the expected false discovery rate amongst SNPs for which $\hat{cFDR}(p_{i} \mid p_{j}) \leq \alpha_{j}^{i}$ is not necessarily lower than $\alpha_{j}^{i}$. For each ordered pair of phenotypes (*i*,*j*), an upper bound $c_{j}^{i}$ was computed for the expected false discovery rates amongst SNPs with $\hat{cFDR}(p_{i} \mid p_{j}) \leq \alpha_{j}^{i}$. The list of SNPs declared non-null for phenotype *i* was pruned to allow for linkage disequilibrium (LD) by listing all SNPs in increasing order of $\left. \mathit{\min}_{j}(\hat{cFDR}(p_{i} \middle| p_{j}) \right.$ and stepping through the list from left to right, at each stage removing all SNPs in LD with *r* ^2^ ≥ 0.1 to the right of the current SNP. This ideally leads to the inclusion of at most one SNP from each LD block.

Multiple testing {#sec004j}
----------------

A multiple testing problem arises from considering p values for one disease conditioned separately on nine others. Specifically, if the criterion for declaring a SNP non-null for phenotype *i* is that $\hat{cFDR}(p_{i} \mid p_{j}) \leq \alpha_{j}^{i}$ for at least one of the nine possible values of *j*, then the FDR for all SNPs declared non-null will be greater than the FDR among the smaller set of SNPs for which $\hat{cFDR}(p_{i} \mid p_{j}) \leq \alpha_{j}^{i}$ for only one value of *j*, due to multiple testing.

However, this excess FDR is not enough to warrant a Bonferroni (Sidak) correction; the $\hat{cFDR}(p_{i} \mid p_{j})$ values for a phenotype *i* are highly correlated, as all are in turn highly correlated with *p* ~*i*~. A Bonferroni correction tends to remove any advantage in SNP detection gained from $\hat{cFDR}$, though an advantage may still be seen when only considering one conditional phenotype *j*.

We opted to use a method proposed by Nyholt \[[@pgen.1004926.ref036]\] to correct for multiple testing in SNPs with high LD. We estimated a correlation matrix Ω for potentially non-null cFDR values using Spearman's rank correlation. The variance of the eigenvalues of Ω, *Var*(*λ* ~*obs*~), was computed and used to estimate the effective number of variables *M* ~*eff*~ according to the equation $$M_{eff} = 1 + 9\left( {1 - \frac{Var\left( \lambda_{obs}) \right.}{9}} \right)$$

Note that *Var*(*λ* ~*obs*~) is between 9 (completely correlated variables, effectively one test) and 0 (completely uncorrelated variables, essentially a Bonferroni correction).

Denoting by $n_{j}^{i}$ the number of SNPs with $\hat{cFDR}(p_{i} \mid p_{j}) \leq \alpha_{j}^{i}$, corresponding to an upper bound on the FDR of $c_{j}^{i}$, an upper bound for the FDR among all SNPs declared significant for phenotype *i* was then computed as $$c_{0}^{i} = M_{eff}\frac{\sum_{j = 1..10,j \neq i}c_{j}^{i}n_{j}^{i}}{\sum_{j = 1..10,j \neq i}n_{j}^{i}},$$ intuitively, multiplying the expected average number of false discoveries across conditional phenotypes ($c_{j}^{i}n_{j}^{i}$) by the effective number of tests. Values of *M* ~*eff*~ and $c_{0}^{i}$ are shown in [S1 Table](#pgen.1004926.s002){ref-type="supplementary-material"}.

Supporting Information {#sec005}
======================

###### Mathematical basis for estimation of cFDR and establishment of upper bound on expected FDR when using $\hat{cFDR}$ method.

(PDF)

###### 

Click here for additional data file.

###### Calculation of false discovery rates for SNPs reaching $\hat{cFDR}$ significance levels.

*M* ~*eff*~ gives the 'effective number of tests', relating to the multiple testing adjustment for the multiple phenotypes conditioned upon (see Methods section). Max p is the maximum principal p value at which a SNP was able to be declared significant using $\hat{cFDR}$. Eq FDR shows the false-discovery rate we would be forced to control at in order to detect all these SNPs the principal p value alone. The FDR bound (bold) is the false discovery rate at which the set of SNPs discovered by the $\hat{cFDR}$ method is controlled.

(CSV)

###### 

Click here for additional data file.

###### [S2](#pgen.1004926.s003){ref-type="supplementary-material"}--[S11](#pgen.1004926.s012){ref-type="supplementary-material"} Tables show associated SNPs for each phenotype, ordered by best cFDR.

P values shown are after adjustment for genomic inflation. Chromosome positions are from the NCBI36 assembly. The conditional phenotype shown is the phenotype for which the cFDR was most below the relevant cutoff. Column CP is the conditional phenotype for which corrected cFDR was lowest. SNPs with p value greater than 5 × 10^−8^ for the principal phenotype are asterisked, and SNP-disease associations not previously known are suffixed with a'+'. [S2 Table](#pgen.1004926.s003){ref-type="supplementary-material"} shows SNPs associated with T1D (type 1 diabetes).

(CSV)

###### 

Click here for additional data file.

###### SNPs associated with ATD (autoimmune thyroid disease).

See legend for [S2 Table](#pgen.1004926.s003){ref-type="supplementary-material"}.

(CSV)

###### 

Click here for additional data file.

###### SNPs associated with CEL (celiac disease).

See legend for [S2 Table](#pgen.1004926.s003){ref-type="supplementary-material"}.

(CSV)

###### 

Click here for additional data file.

###### SNPs associated with MS (multiple sclerosis).

See legend for [S2 Table](#pgen.1004926.s003){ref-type="supplementary-material"}.

(CSV)

###### 

Click here for additional data file.

###### SNPs associated with NAR (narcolepsy).

See legend for [S2 Table](#pgen.1004926.s003){ref-type="supplementary-material"}.

(CSV)

###### 

Click here for additional data file.

###### SNPs associated with PBC (primary biliary cirrhosis).

See legend for [S2 Table](#pgen.1004926.s003){ref-type="supplementary-material"}.

(CSV)

###### 

Click here for additional data file.

###### SNPs associated with PSO (psoriasis).

See legend for [S2 Table](#pgen.1004926.s003){ref-type="supplementary-material"}.

(CSV)

###### 

Click here for additional data file.

###### SNPs associated with RA (rheumatoid arthritis).

See legend for [S2 Table](#pgen.1004926.s003){ref-type="supplementary-material"}.

(CSV)

###### 

Click here for additional data file.

###### SNPs associated with UC (ulcerative colitis).

See legend for [S2 Table](#pgen.1004926.s003){ref-type="supplementary-material"}.

(CSV)

###### 

Click here for additional data file.

###### SNPs associated with CRO (Crohn's disease).

See legend for [S2 Table](#pgen.1004926.s003){ref-type="supplementary-material"}.

(CSV)

###### 

Click here for additional data file.

###### Effect of adjusting $\hat{cFDR}$ for shared controls.

Plots A and B show the ratio between the true $\hat{cFDR}$ (computed using our method) to the'naive' $\hat{cFDR}$ (computed by naively applying the existing split-control approach to shared-control data without adjustment) for a range of p values for the principal phenotype. The p values forming the x-coordinates were obtained from the shared-control design. The left-hand plot shows ratios of true to'naive' cFDR for p values greater than 1 × 10^?10^, demonstrating 2--3 fold underestimation. The right-hand plot shows log-ratios of trueto'naive' cFDR for smaller p values, demonstrating hundred- or thousand- fold underestimation.

(PDF)

###### 

Click here for additional data file.

###### Summary of pleiotropy between phenotypes.

The colour for phenotype *i* (horizontal) and phenotype *j* (vertical) corresponds to the p-value cutoff for significance for phenotype *i*, given that a p-value cutoff for phenotype *j* is less than 5 × 10^−6^.

(PDF)

###### 

Click here for additional data file.

###### Distribution of log(*M*) amongst null SNPs.

*M* is proportional to the variance of the log odds ratio from TDT data, defined as *σ̂* ^2^ *f*(1−*f*)), where *f* is the minor allele frequency amongst null SNPs, and *σ̂* is the standard error. Equating the median of *M* with a known expression for variance of the log odds ratio in a case-control study enables back-calculation of the effective number of cases and controls. This technique was used for computing the number of cases and controls in the T1D study, for which p values were obtained from a meta-analysis of case-control and TDT data.

(PDF)

###### 

Click here for additional data file.

###### [S4](#pgen.1004926.s016){ref-type="supplementary-material"}--[S13](#pgen.1004926.s025){ref-type="supplementary-material"} Figs. are Q-Q plots labeled "*i*∣*j*", where *i* is the principal phenotype and *j* the conditional phenotype.

Y axes show $log_{10}(p_{i}^{\prime})$; X axes show log quantile (rank) of p values in various sets of SNPs. Each colour corresponds to the Q-Q plot for *p* ~*i*~ amongst only SNPs such that *p* ~*j*~ is less than a certain cutoff, with the black line corresponding to the Q-Q plot for all SNPs. P values for the principal phenotype are adjusted for the effect of shared controls between studies. A leftward shift with decreasing *p* ~*j*~ cutoff indicates enrichment of SNP sets from conditioning on degrees of association with a conditional phenotype, probably due to pleiotropic effects between phenotypes. Because the studies used the ImmunoChip, which covers only potential autoimmune-associated regions, the black line also shows considerable enrichment compared to quantiles. [S4 Fig.](#pgen.1004926.s016){ref-type="supplementary-material"} shows Q-Q plots with T1D (type 1 diabetes) as the principal phenotype

(PDF)

###### 

Click here for additional data file.

###### Q-Q plots with ATD (autoimmune thyroid disease) as the principal phenotype.

See legend for [S4 Fig.](#pgen.1004926.s016){ref-type="supplementary-material"}

(PDF)

###### 

Click here for additional data file.

###### Q-Q plots with CEL (celiac disease) as the principal phenotype.

See legend for [S4 Fig.](#pgen.1004926.s016){ref-type="supplementary-material"}

(PDF)

###### 

Click here for additional data file.

###### Q-Q plots with MS (multiple sclerosis) as the principal phenotype.

See legend for [S4 Fig.](#pgen.1004926.s016){ref-type="supplementary-material"}

(PDF)

###### 

Click here for additional data file.

###### Q-Q plots with NAR (narcolepsy) as the principal phenotype.

See legend for [S4 Fig.](#pgen.1004926.s016){ref-type="supplementary-material"}

(PDF)

###### 

Click here for additional data file.

###### Q-Q plots with PBC (primary biliary cirrhosis) as the principal phenotype.

See legend for [S4 Fig.](#pgen.1004926.s016){ref-type="supplementary-material"}

(PDF)

###### 

Click here for additional data file.

###### Q-Q plots with PSO (psoriasis) as the principal phenotype.

See legend for [S4 Fig.](#pgen.1004926.s016){ref-type="supplementary-material"}

(PDF)

###### 

Click here for additional data file.

###### Q-Q plots with RA (rheumatoid arthritis) as the principal phenotype.

See legend for [S4 Fig.](#pgen.1004926.s016){ref-type="supplementary-material"}

(PDF)

###### 

Click here for additional data file.

###### Q-Q plots with UC (ulcerative colitis) as the principal phenotype.

See legend for [S4 Fig.](#pgen.1004926.s016){ref-type="supplementary-material"}

(PDF)

###### 

Click here for additional data file.

###### Q-Q plots with CRO (Crohn's disease) as the principal phenotype.

See legend for [S4 Fig.](#pgen.1004926.s016){ref-type="supplementary-material"}

(PDF)

###### 

Click here for additional data file.

We thank Professor John Todd for help in interpreting the results. Thanks to Professor Cathryn Lewis and Dr Steve Eyre for the list of SNPs used for genomic control. We gratefully acknowledge the following groups and individuals who provided biological samples or data for this study. The IIBDGC shared summary association data on Crohn's and ulcerative colitis. We would like to thank the UK Medical Research Council and Wellcome Trust for funding the collection of DNA for the British 1958 Birth Cohort (MRC grant G0000934, WT grant 068545/Z/02). We thank The Avon Longitudinal Study of Parents and Children laboratory in Bristol and the British 1958 Birth Cohort team, including S. Ring, R. Jones, M. Pembrey, W. McArdle, D. Strachan and P. Burton, for preparing and providing the control DNA samples. We acknowledge use of DNA from The UK Blood Services collection of Common Controls (UKBS collection). The collection was established as part of the Wellcome Trust Case-Control Consortium. We acknowledge use of DNA samples from the NIHR Cambridge BioResource. We thank volunteers for their support and participation in the Cambridge BioResource and members of the Cambridge BioResource SAB and Management Committee for their support of our study. Access to Cambridge BioResource volunteers and their data and samples is governed by the Cambridge BioResource SAB. Documents describing access arrangements and contact details are available at <http://www.cambridgebioresource.org.uk/>.

Acknowledgements for the individual studies whose summary data we accessed via ImmunoBase may be found at <http://www.immunobase.org/poster/immunochip-paper-acknowledgments/>.

[^1]: I have read the journal\'s policy and the authors of this manuscript have the following competing interests: The JDRF/Wellcome Trust Diabetes and Inflammation Laboratory receives funding from Hoffmann La Roche and Eli-Lilly and Company. ImmunoBase is funded in part by Eli-Lilly and Company.

[^2]: Conceived and designed the experiments: JL CW. Performed the experiments: JL. Analyzed the data: JL. Contributed reagents/materials/analysis tools: JL CW. Wrote the paper: JL CW.
