Genome-wide association and HLA region fine-mapping studies identify susceptibility loci for multiple common infections

Infectious diseases have a profound impact on our health and many studies suggest that host genetics play a major role in the pathogenesis of most of them. We perform 23 genome-wide association studies for common infections and infection-associated procedures, including chickenpox, shingles, cold sores, mononucleosis, mumps, hepatitis B, plantar warts, positive tuberculosis test results, strep throat, scarlet fever, pneumonia, bacterial meningitis, yeast infections, urinary tract infections, tonsillectomy, childhood ear infections, myringotomy, measles, hepatitis A, rheumatic fever, common colds, rubella and chronic sinus infection, in over 200,000 individuals of European ancestry. We detect 59 genome-wide significant (P < 5 × 10−8) associations in genes with key roles in immunity and embryonic development. We apply fine-mapping analysis to dissect associations in the human leukocyte antigen region, which suggests important roles of specific amino acid polymorphisms in the antigen-binding clefts. Our findings provide an important step toward dissecting the host genetic architecture of response to common infections.

I nfectious diseases, the second leading cause of death worldwide, represent persistent challenges to human health due to increasing resistance to established treatments, lack of life-saving vaccines and medications in developing countries, and the ongoing emergence of new pathogens 1 . Studies have linked susceptibility to infectious agents to cancers, autoimmune diseases and drug hypersensitivity. Human papillomaviruses have been associated with multiple cancers 2 ; rubella and mumps infection with type 1 diabetes (T1D) in children 3 ; and Varicella zoster (shingles) with multiple sclerosis 4 . Reactivation of chronic persistent human herpesviruses has been linked to drug-induced hypersensitivity 5 . Thus, infectious diseases have a profound impact on our health, both directly as well as through connections with other diseases. Nevertheless, genetic studies of common infectious diseases have lagged behind those of other major complex diseases. In particular, only a few genome-wide association studies (GWASes) have been undertaken for common infectious diseases of lower mortality or for infectious diseases for which vaccines are available (Supplementary Table 1). Previous work has linked variants in human leukocyte antigen (HLA) region to shingles defined from electronic medical record data 6 , to the risk for chronic hepatitis B (CHB) virus infection in Asian populations 7 and to the Epstein-Barr virus infection in Mexican American families 8 . Outside of HLA region, a variant near HORMAD2 at 22q12.2 (rs2412971, P = 1.48 × 10 −9 ) was recently reported to be associated with tonsillectomy 9 . Studies in a west African population have reported associations between two variants (rs4331426 and rs2057178) and susceptibility to tuberculosis (TB) 10,11 , and a later study found that the association with rs4331426 was not significant in an Asian population 12 . A prior study established associations between meningitis and variants in complement factor H (CFH; rs1065489) and in CFHrelated protein 3 (CFHR; rs426736) 13 . Neisseria meningitidis is known to evade complement-mediated killing by the binding of host CFH to the meningococcal factor H-b-bind protein (fHbp).
We conducted 23 imputed GWASes for 23 common infections and infection-associated procedures (Table 1) with data drawn from more than 200,000 23andMe research participants of European ancestry who completed a standardized questionnaire about their history of infections (Supplementary Notes 1 and 2). In this study, we detected 59 genome-wide significant (GWS) associations and in some instances our findings replicated previous GWASes. We further imputed and tested the HLA classical alleles and amino acid polymorphisms to dissect independent HLA signals for multiple common infections. Our fine-mapping studies identify novel HLA associations in European population and suggest important roles of specific amino acid polymorphisms in the antigen-binding clefts. By studying multiple infections together, we find that virus-induced infectious diseases are mainly associated with class I molecules and the bacteria-induced infectious diseases are mainly associated with class II molecules, while with notable exceptions, as discussed below.

Results
Study summary. In total, 59 GWS regions (P < 5 × 10 −8 ) were discovered from the 23 imputed GWASes for 17 common infections we studied ( Table 1). The diseases phenotypes were described in Supplementary Notes 1 and 2. Summary information for the index single-nucleotide polymorphism (SNP with lowest P value) is presented in Tables 2 and 3, and GWAS results are displayed as Manhattan plots in Supplementary Fig. 1. All association results were adjusted for the corresponding genomic control inflation factors, and the Q-Q plots do not exhibit significant inflation in the P values obtained ( Supplementary  Fig. 1). Regional association plots illustrate patterns of association for the index SNP and other variants within each associated genomic region ( Supplementary Fig. 2). The strength of association in the region, previous reports of associations from the National Human Genome Research Institute GWAS Catalog entries, nearby expression quantitative trait loci (eQTL) and coding SNPs are discussed below and are also summarized in Supplementary Tables 3 and 5. The HLA region is significantly associated (P < 5 × 10 −8 ) with 13 common infections, and we dissect independent signals in each association ( Fig. 1 and We estimated the genetic heritability and correlations of these common infections using the GWAS data. The proportion of variance attributable to genome-wide SNPs (SNP heritability) was estimated to be 6% on average for the common infections. Multiple infections showed significant positive genetic correlations (P < 2 × 10 −4 , Supplementary Fig. 4). The results are presented in Supplementary Tables 6 and 7 and discussed in detail in Supplementary Notes 3 and 4.
Chickenpox and shingles. We identified independent associations of adults with a history of chickenpox with HLA-A and HLA-B in the class I region. On conditional analysis within HLA-A, the amino acid polymorphism at position 107 of HLA-A protein (HLA-A Gly107, P = 3.77 × 10 −10 , odds ratio (OR) = 1.09) accounted for the HLA allele association (HLA-A*02:01, P = 1.08 × 10 −9 , OR = 0.92, conditional P value is 0.90) in this interval. Although HLA-A Gly107 is not located in the peptidebinding cleft, it is in high linkage disequilibrium (LD) with HLA-A Asp74, HLA-A Gly62 and HLA-A V95 (r 2 > 0.95), which are all in the peptide-binding cleft. In HLA-B, only rs9266089 (P = 1.00 × 10 −10 , OR = 1.12) met the threshold for GWS.
We identified a variant upstream of IFNA21 (rs7047299, P = 1.67 × 10 −8 , OR = 1.07) as a GWS association with shingles. None of the variants within 500 kb and in moderate LD (r 2 > 0.6) with rs7047299 were coding, nor were they reported as eQTL. IFNA21 encodes type I interferon and is mainly involved in innate immune response against viral infection. It has been shown to be involved in the pathogenesis of rubella 14 and may also influence susceptibility to asthma and atopy 15 .
Cold sores. We identified independent associations between two SNPs in the class I region and a history of cold sores. One is 2 kb upstream of POU5F1 indexed by rs885950 (P = 7.47 × 10 −13 , OR = 1.08) and the other is upstream of HCP5 indexed by rs4360170 (P = 3.41 × 10 −9 , OR = 1.13). We also found a GWS association with HLA-B-ThrGly45 in the peptide-binding cleft of the HLA-B protein (P = 4.91 × 10 −12 , OR = 0.91). Upon conditioning on HLA-B ThrGly45, the two SNPs had significant residual effects (conditional P values are 6.28 × 10 −9 and 1.32 × 10 −5 ). However, the effect of HLA-B ThrGly45 was largely removed (conditional P value is 0.001) after conditioning on the two SNPs.
Mononucleosis. We identified a variant upstream of HCP5 in the HLA class I region as a GWS association (indexed by rs2596465, P = 2.48 × 10 −9 , OR = 1.08). No GWS HLA classical allele or amino acid variant was identified.
FUT2 (rs516316, P = 9.63 × 10 −72 , OR = 1.25) and ST3GAL4 (rs3862630, P = 1.21 × 10 −8 , OR= 1.13) also have highly significant associations with mumps. Both are components of the glycosphingolipid (GSL) biosynthesis pathway, which was the most significant Meta-Analysis Gene-set Enrichment of variaNT Associations (MAGENTA) -analyzed pathway ( Table 6: P = 1 × 10 −4 and false discovery rate (FDR) = 0.003) identified using mumps GWAS data. This pathway also includes the ABO gene (rs8176643, P = 5.81 × 10 −5 , OR = 1.07). GSLs are found in cell membranes of organisms ranging from bacteria to humans and play important biological roles in membrane structure, host-pathogen interactions, cell-cell recognition and modulation of membrane protein function. Genetic variation in human FUT2 determines whether ABO blood group antigens are secreted into body fluids 16 . Detailed functional annotations of the genes  18 ('non-secretor' allele is linked to lower rate of childhood ear infections in our GWAS, Supplementary Fig. 5A), and increase risk to Crohn's disease 74 and T1D 19 . These results are interesting in light of studies showing that the number of T1D cases increase significantly after a mumps epidemic, suggesting that mumps virus might promote autoimmunity and trigger T1D 3 . The surface glycoprotein of mumps virus, hemagglutinin-neuraminidase, attaches to sialic acid receptors and promotes fusion and viral entry into host cells 20 . ABO antigens can modulate the interaction between pathogens and cell surface sialic acid receptors 21 ; this presents a plausible mechanism by which secretor status could modify susceptibility to mumps infections. ST3GAL4 encodes a member of the glycosyltransferase 29 family, a group of enzymes involved in protein glycosylation. We also identified a significant association between mumps and a variant on chromosome 14 (rs11160318, P = 4.56 × 10 −12 , OR = 1.1) that is located in the intergenic region upstream of BDKRB2. The role of BDKRB2 in mumps susceptibility is not clear, but the BDKRB2 product has high affinity for intact kinins, which mediate a wide spectrum of biological effects including inflammation, pain, vasodilation, and smooth muscle contraction and relaxation 22 .
Hepatitis A. Our GWAS showed a suggestive association with IFNL4 (rs66531907, P = 5.7 × 10 −8 , OR = 1.23, 1 kb upstream of IFNL4). IFNL4 is located upstream of IFNL3 (also known as IL28B). The high-risk allele rs66531907 (C) in our hepatitis A GWAS is in almost complete LD with rs8099917 (T) (r 2 = 0.992), which has been associated with progression to chronic hepatitis C infection (HCV) and poor response to HCV therapy in multiple studies 23 .
Hepatitis B. Our GWAS of hepatitis B virus (HBV) infection in Europeans identified a GWS association with rs9268652 (P = 3.14 × 10 −9 , OR = 1.32, in HLA-DRA), consistent with the association of the HLA class II region with CHB susceptibility. Fine-mapping analysis failed to identify GWS associations with HLA classical variants, but DQB1*06:02 was moderately associated with hepatitis B (P = 5.47 × 10 −6 , OR = 0.75).
We also identified a GWS association with rs6692209 (P = 5.25 × 10 −9 , OR = 1.08) near LCE3E in the Epidermal Differentiation Complexn 1q21. The EDC comprises a remarkable density of gene families that determine differentiation of the human epidermis 24 . The high-risk allele rs6692209 (T) is in high LD (r 2 = 0.86) with rs2105117 (A) (P = 2.98 × 10 −8 , OR = 1.07, in our GWAS), a missense variant in LCE5A. Variants in late cornified envelope (LCE) genes have been implicated in GWASes of atopic dermatitis (rs3126085-A, P = 6 × 10 −12 , r 2 = 0.047 with rs6692209) 25 and psoriasis (rs4085613-A, P = 7 × 10 −30 , r 2 = 0.002 with rs6692209) 26 , although the index SNPs identified in those GWASes were in low LD with our identified region. While studies have linked the EDC to various skin diseases, the role of the LCE proteins in epidermal biology has not been extensively studied. One study suggested that LCE proteins of groups 1, 2, 5 and 6 are involved in normal skin barrier function, while LCE3 genes encode proteins that are involved in barrier repair after injury or inflammation 27 .
We also found a GWS association with a small insertion, rs35395352, in 1p36.23 (P = 3.90 × 10 −8 , effect = 0.03). None of the variants within 500 kb and in moderate LD (r 2 > 0.6) with this SNP were coding. However, two variants have been reported as cis-eQTLs: rs7548511 (P = 7.84 × 10 −8 , effect = 0.03, risk allele A, r 2 = 0.96 with rs35395352) is an eQTL for ENO1 (eQTL P value is 8 × 10 −5 ) in lymphoblastoid 28 and rs11121129 (P = 9.62 × 10 −8 , effect = 0.03, risk allele G, r 2 = 0.96 with rs35395352) is an eQTL for TNFRSF9 (eQTL P value is 4.4 × 10 −4 ) in lymphoblastoid 29 . The product of ENO1, alpha enolase, has been identified as an autoantigen in Hashimoto's encephalopathy 30 and a putative autoantigen in severe asthma 31 and Behcet's disease 32 . The TNFRSF9-encoded protein is a member of the tumor necrosis factor (TNF) -receptor superfamily, which contributes to the clonal expansion, survival and development of T cells. The same 1p36.23 region has been implicated in GWASes of psoriasis (rs11121129, risk allele A, P = 1.7 × 10 −8 , r 2 = 0.96 with rs35395352) 33 . Interestingly, the strep throat risk allele rs11121129-A is protective against psoriasis. Studies found that exacerbation of chronic psoriasis can be associated with streptococcal throat infections, and the T cells generated by palatine tonsils can recognize skin keratin determinants in patients' blood 34 .
Pneumonia. We found a GWAS association between pneumonia and a variant in the HLA class I region (rs3131623, P = 1.99 × 10 −15 , OR = 1.1). Conditioning on rs3131623, the effects of HLA-B SerTrpAsn97 (P = 2.37 × 10 −12 , OR = 0.94, conditional P value is 9.33 × 10 −4 ), in the peptide-binding cleft, was significantly reduced. HLA-DRB1 LS11 in class II region also showed a significant independent effect (P = 4.10 × 10 −11 , OR = 0.94). Conditioning on rs3131623, HLA-DRB1 LS11 remained moderately significant (conditional P value is 3.89 × 10 −6 ). Bacterial meningitis. We identified a GWS association with CA10 (carbonic anhydrase X) (rs1392935, P = 1.19 × 10 −8 , OR = 1.78, intronic of CA10). None of the variants within 500 kb and in moderate LD (r 2 > 0.6) with rs1392935 were coding, nor were they reported as an eQTL. However, rs1392935 falls within an enhancer that is defined in the fetal brain. The protein encoded by CA10 is an acatalytic member of the alpha-carbonic anhydrase subgroup and it is thought to play a role in the central nervous system, especially in brain development 35 .
Vaginal yeast infection. We found GWS associations between vulvovaginal candidiasis and variants in DSG1 (rs200520431, P = 1.87 × 10 −16 , effect = 0.11). The index SNP rs200520431 is intronic and in high LD (r 2 > 0.8) with multiple missense mutations (rs8091003, rs8091117, rs16961689, rs61730306, rs34302455) in DSG1. The DSG1 gene product is a calciumbinding transmembrane glycoprotein component of desmosomes in vertebrate epithelial cells. It connects the cell surface to the keratin cytoskeleton and plays a key role in maintaining epidermal integrity and barrier function 36 . This glycoprotein has been identified as the autoantigen of the autoimmune skin blistering disease pemphigus foliaceus 37 and homozygous mutations in DSG1 have been showed to result in severe dermatitis, multiple allergies, and metabolic wasting syndrome 36 . A variant downstream of PRKCH (rs2251260, P = 3.46 × 10 −10 , effect = 0.05) was also significantly associated with yeast infections. None of the variants within 500 kb and in moderate LD (r 2 > 0.6) with the index SNP rs2251260 were coding, nor were they reported as an eQTL. However, rs2251260 falls within strong enhancers defined in hepatocellular carcinoma. PRKCH is a member of a family of serine-and threonine-specific protein kinase and is predominantly expressed in epithelial tissues. The PRKCH protein kinase can regulate keratinocyte differentiation 38 .
We also found a significant association with a variant in the 14q32.2 gene desert (rs7161578-T, P = 4.04 × 10 −8 , effect = 0.04). The index SNP is in high LD with many enhancer sequences defined in epidermal keratinocytes. A variant in the same region (rs7152623-A, P = 3 × 10 −15 , r 2 = (+) 0.98 with rs7161578-T) was implicated in aortic stiffness 39 .
Urinary tract infection. We found an association between urinary tract infections (UTIs) and variants in PSCA (rs2976388, P = 3.27 × 10 −10 , effect = 0.04). The index SNP rs2976388 falls within strong enhancers defined in epidermal keratinocytes. It is also in high LD (r2 = 0.95) with rs2294008 (P = 1.01 × 10 −9 ), which is in the 5'-untranslated region promoter sequence of PSCA and can cause changes in transcriptional repressor CCCTC-binding factor (CTCF) motif. The PSCA gene variants conferred risk of UTI only in females when we further tested the effect in female and male cohorts separately (P < 1e-11 in female, and P > 0.1 in male). Prior GWASes have found associations between the PSCA gene and duodenal ulcer (rs2294008-C, P = 2 × 10 −33 ) 75 and bladder cancer (rs2294008-T, P = 4 × 10 −11 ) 43 . PSCA encodes a glycosylphosphatidylinositol-anchored cell membrane glycoprotein of unknown function. This glycoprotein was initially identified as a prostate-specific cell-surface marker 42 and is overexpressed in a large proportion of prostate cancers, but also detected in bladder and pancreatic cancers 43 .  We also identified a variant intronic of FRMD5 (rs146906133, P = 2.02 × 10 −8 , effect = 0.38) that is in LD with rs138763871 (r 2 = 0.8, P = 1.75 × 10 −6 , effect = 0.33), a missense variant in STRC (R1521W). The role of FRMD5 or STRC in UTIs is unclear.
Tonsillectomy. Our tonsillectomy GWAS revealed 35 independent GWS associations. We observed strong HLA associations at HLA-B in the class I region. The amino acid polymorphism HLA-B Val97 (P = 1.44 × 10 −21 , OR = 1.21) in the peptide-binding cleft of the HLA-B protein, accounted for most of the SNP effect (rs41543314, P = 5.36 × 10 −21 , OR = 1.23, conditional P value is 0.06) and HLA allele effect (HLA-B*57:01, P = 6.29 × 10 −22 , OR = 1.22, conditional P value is 0.05) in this interval. We also noted additional independent signals in the class II and class III regions. Within HLA-DRB1 in the class II region, the effect of rs140177540 (P = 5.02 × 10 −15 , OR = 1.06, conditional P value is 0.001) was largely removed after conditioning on HLA-DRB1 LeuValGly11 (P = 1.02 × 10 −12 , OR = 1.06) in the peptide-binding cleft of the HLA-DRB1 protein.
We detected 34 additional GWS associations (Table 3 and  Supplementary Table 4 IgA is an antibody that plays a critical role in mucosal immunity. Tonsils belong to nasopharyngeal-associated lymphoid tissues and the generation of B cells is considered one of the major tonsillar functions; secretory dimeric IgA produced by the B cells is capable of preventing absorption and penetration of bacteria and/or viruses into the upper respiratory tract mucosa 44 . Tonsillectomy has been shown to significantly decrease levels of serum IgA and salivary secretory IgA levels 45 . Other identified significant signaling pathways include 'TACI and BCMA stimulation of B-cell immune system' (P = 9 × 10 −4 and FDR = 0.008), Although tonsillectomy has been performed for over 100 years, possible immunological effects of this procedure remain controversial. Many reports have demonstrated that in certain patients tonsillectomy is an effective therapy for psoriasis 34 and rheumatoid arthritis (RA) 46 ; the rationale for this effect is unknown. The 'intestinal immune network for IgA production' pathway was also identified as the most significant pathway in a recent GWAS of IgAN 47 . In that study, they also linked the intestinal mucosal inflammatory disorders and inflammatory bowel disease (IBD) with risk of IgAN. Our tonsillectomy GWAS showed significant overlap of the identified loci with those autoimmune and inflammatory disorders 40 . We found five risk loci that are also associated with RA (HLA, TNFIP3, CD40, SPRED2 and SH2B3), four loci associated with IBD (HLA, MTMR3, TNFAIP3 and CD40), four loci associated with ulcerative colitis (HLA, GNA12, MAPK3 and NFKB1), three loci associated with psoriasis (HLA, SLC12A8 and 1p36.23), and two loci that are shared with IgAN (HLA and MTMR3). We observed both concordant and opposing effects compared to these immune-mediated diseases. Our results may help to elucidate the connection between tonsillectomy and these diseases and may provide insight into clinical markers that could be used as indicators of tonsillectomy as a therapy for these diseases.
We found 13 additional GWS signals (Table 2 and Supplementary Table 5). Although MAGENTA pathway analysis did not identify significant pathways (Table 8) of interest, variants in FUT2 and ABO, which are involved in the GSL biosynthesis pathway and were implicated in our mumps GWAS, were also significantly associated with childhood ear infection. The risk allele rs681343(C) (P = 3.51 × 10 −30 , OR = 1.11) is a synonymous mutation in FUT2 and is in almost complete LD with rs601338 (G) (r 2 = 0.9993), the 'secretor' (se) allele that is associated with higher risk of childhood ear infection in our data. This is consistent with a previous report in which secretors were overrepresented among patients with upper respiratory infections 18 . Our second most significant association was with a variant in TBX1 (rs1978060, P = 1.17 × 10 −19 , OR = 1.09). The low-risk allele rs1978060 (A) is in LD (r 2 = 0.45) with rs72646967 (C), a missense (N397H) mutation in TBX1 that was also implicated in our tonsillectomy GWAS. T-box genes encode transcription factors that play important roles in tissue and organ formation during embryonic development. In mice, TBX1 haploinsufficiency in the DiGeorge syndrome region has been showed to disrupt the development of the pharyngeal arch arteries 41 and the middle and outer ear 59 . We also identified other genes involved in developmental processes, including MKX (rs2808290, P = 5.09 × 10 −16 , OR = 1.07), FGF3 (rs72931768, P = 2.63 × 10 −9 , OR = 1.09) and AUTS2 (rs35213789, P = 3.75 × 10 −9 , OR = 1.06).
The MKX product is an IRX family-related homeobox transcription factor and has been shown to play a critical role in tendon differentiation by regulating type I collagen production in tendon cells 49 . FGF3 is a member of the fibroblast growth factor family, which is involved in embryonic development, cell growth, morphogenesis, tissue repair and invasion. A study of similar genes in mouse and chicken suggested a role in inner ear formation 50 . AUTS2 has been implicated in neurodevelopment and is a candidate gene for autism spectrum disorders and developmental delay 76 . Finally, we identified missense mutations in CDHR3 (rs114947103, P = 5.40 × 10 −9 , OR = 1.07) and PLG (rs73015965, P = 3.78 × 10 −8 , OR= 1.43). The index SNP rs114947103 is in almost complete LD with rs6967330 (r2 = 0.97), a missense variant in CDHR3, that was recently identified in a GWAS as a susceptibility locus for asthma 52 . The biological role of CDHR3 is not clear, but it is highly expressed in airway epithelium and belongs to the cadherin family that is involved in cell adhesion and epithelial polarity 52 . Mutations in the PLG gene could cause congenital plasminogen deficiency, which results in inflamed growths on the mucous membranes. Studies in mice have showed that PLG plays an essential role in protecting against the spontaneous development of chronic otitis media (middle ear infection) and have also suggested the possibility of using PLG for clinical therapy of certain types of otitis media 53 . The same variant in TBX1 that was associated with tonsillectomy and childhood ear infection (see earlier) was also found to be significantly associated with myringotomy (rs1978060, P = 3.34 × 10 −10 , OR = 1.17).

Discussion
This is the first description of GWASes for numerous common infectious diseases in a single study. Our GWASes replicated the previously reported associations with shingles and tonsillectomy at much higher significance (Supplemental Table 1). In other instances our GWASes did not replicate prior findings and instead identified novel associations. Our GWASes of HBV infection, mononucleosis and positive TB test identified significant HLA associations; however, they are not in strong LD with any of the previously reported loci (Supplemental Table 1). For the associations that did not replicate, several factors may have reduced our power. First, there was genetic heterogeneity in different ethnic groups due to differences in allele frequency and LD structure, especially in the HLA region. Many complex diseases have been shown to be influenced by different genetic determinants in different ethnic groups 54 . Second, there were differences in how the disease phenotypes were determined. Our hepatitis B phenotype was defined by participants' self-declared history of HBV infection and thus combines cases of subclinical infection that were detected after the fact, acute infections that spontaneously cleared, successfully treated infections, and chronic persistent infection cases ('Survey and phenotype scoring logic' in Supplementary Notes 1). This definition differs from studies that focus solely on chronic HBV carriers. Similarly, our GWAS of patients with a history of a positive TB tests tested a different phenotype than published studies of patients with active TB.
Insufficient power due to small sample size may have in some instances contributed to our results differing from previous studies. Our GWASes of susceptibility to meningitis (based on a sample of 842 cases) failed to replicate the findings at CFH, but identified a significant association at CA10. Moreover, our case definition of self-reported 'meningitis' may be confounded by the fact that bacterial meningitis is frequently considered in patients with presenting to the hospital with severe headache and fever, a definitive diagnosis is made in a minority of cases 55 . In other words, patients initially suspected to have bacterial meningitis may ultimately be diagnosed with other syndromes such as viral meningitis, stroke, subdural empyema and cerebral abscess 55 .
Notwithstanding the limitations, our study identified numerous significant associations and also suggested crucial roles of innate immune responses in protective immunity to multiple common infections. For example, two genes involved in innate immunity of skin, LCE3E and DSG1, were significantly associated with plantar warts and yeast infections, respectively. We identified the greatest number of independent associations for tonsillectomy (n = 35) and childhood ear infection (n = 14), two relatively nonspecific phenotypes. The large number of associations may reflect, in part, the large sample sizes for these two traits (n = 173,412 and 121,810, respectively). There were many associations between both tonsillitis and middle ear infections with genes playing roles in the innate immune response. Lymphoid hyper-reaction throughout the mucosa-associated lymphoid tissue system has been suggested as a molecular mechanism underlying the genetic association with tonsillectomy 9 . In addition, we discovered associations between both syndromes and genes with roles in embryonic development, underscoring the important role of host anatomy in both syndromes. Patients undergo tonsillectomy for multiple reasons, including recurrent tonsillitis, obstructive sleep apnea and nasal airway obstruction 56 . Similarly, childhood ear infections are strongly influenced by the shape of the Eustachian tube 57 . The missense mutation (rs72646967-C, N397H) in TBX1, which is essential for inner ear development 58,59 , is associated with a lower risk of both tonsillectomy and inner ear infections ( Supplementary  Fig. 5B).
In the HLA region, we found that viral diseases-e.g., chickenpox, shingles, cold sores, mononucleosis (all caused by members of the human herpesvirus) as well as mumps (caused by mumps virus)-were mainly associated with variation in class I molecules (Fig. 1). The bacterial diseases-specifically having a positive TB test (caused by Mycobacterium tuberculosis), scarlet fever (caused by Streptococcus pyogenes), and childhood ear infection (primarily caused by bacteria)-were mainly associated with variation in class II molecules. Tonsillectomy and pneumonia, caused by either bacteria or viruses, were associated with both class I and class II molecules. These observations are consistent with previous knowledge about antigen presentation; viruses mostly replicate within nucleated cells in the cytosol and produce endogenous antigens that are presented by class I major histocompatibility complex (MHC) molecules, while-with notable exceptions-the bacteria responsible for the diseases in question grow primarily extracellularly and are taken up by endosomal compartments where they are processed for presentation by class II MHC molecules 60 . The two intracellular pathways of protein processing may not be completely separate. Activation of CD8 + , HLA class I-restricted T cells by exogenous antigens has been reported 61 and HLA class II-restricted cytotoxic T cells that recognize endogenously synthesized antigen, such as HBV envelope antigens, have also been described 62 .
Some diseases, however, do not follow the common principles of antigen presentation. HPV infections, causing plantar warts, are exclusively intraepithelial and there is no viremia or cytolysis in the infection cycle. The primary response to HPV antigens is more likely to be initiated by the antigen-presenting cells of squamous epithelia, the Langerhans cells (LCs). LCs capture antigens by macropinocytosis and receptor-mediated endocytosis, and then initiate class II processing of exogenous antigens 48 , which may explain why our plantar warts GWAS mainly identified associations with MHC class II molecules. In other cases, our results suggest even more complex interactions. Strep throat (caused by S. pyogenes) has both MHC class I and class II associations. Some of our reported strep throat cases may represent secondary bacterial infections following a viral upper respiratory infection, meaning that the discovered MHC class I associations actually reflect the initiating viral infection 63 . Some cases may have been misdiagnosed or misreported as strep throat, given that many sore throats are actually caused by viruses 63 . A final possibility is that HLA class I molecules can bind bacterial peptides derived from exogenous proteins that are internalized by endocytosis or phagocytosis, a process called crosspresentation 64 .
Although additional studies will be required to validate our associations, our findings are an important step toward dissecting the host genetic contribution to variation in response to infections. Research insights into infectious diseases will help to drive new diagnostic approaches and perhaps new therapies and preventions. Moreover, our results may also yield insights into autoimmune disorders that are associated with infectious triggers. A postulated mechanism for the relationship between infection and immunological disorders is that certain pathogens elicit immune responses to cross-reactive self epitopes, thereby initiating a cycle of damage and further immune activation 65 . Indeed, our study identified at least 10 infectious disease associations with genes that have been previously associated with autoimmune diseases.

Methods
Subjects. Participants used in our GWASes were drawn from the customer base of 23andMe, Inc, a personal genetics company. All individuals included in the analyses provided informed consent and answered surveys online according to our human subjects protocol, which was reviewed and approved by Ethical and Independent Review Services, a private institutional review board http://www. eandireview.com). Over 200,000 genotyped participants were included in analyses based on selection for having >97% European ancestry as determined through an analysis of local ancestry 51 and completing surveys about their history of common infections (Table 1, Supplementary Table 2, 'Survey and phenotype scoring logic' in Supplementary Notes 1). Except for the history of chickenpox infections, participants were also asked whether they have ever received the chickenpox vaccine, so that vaccinated individuals were excluded from the controls in the chickenpox study. The vaccination status for many other infections such as mumps, measles and rubella was not surveyed and the controls for them included vaccinated individuals ('Surveys and phenotype scoring logic' in Supplementary Note 1). For each study, a maximum set of unrelated individuals was chosen using a segmental identity-by-descent estimation algorithm 66 . Individuals were defined as related if they shared more than 700 cM of identity-by-descent, including regions where the two individuals share either one or both genomic segments identity-by-descent. This level of relatedness (involving~20% of the human genome) corresponds approximately to the minimal expected sharing between first cousins in an outbred population.
Genotyping and SNP imputation. DNA extraction and genotyping were performed on saliva samples by CLIA-certified and CAP-accredited clinical laboratories of Laboratory Corporation of America. Samples were genotyped on one of four genotyping platforms. The V1 and V2 platforms were variants of the Illumina HumanHap550 + BeadChip and contained a total of about 560,000 SNPs, including about 25,000 custom SNPs selected by 23andMe. The V3 platform was based on the Illumina OmniExpress + BeadChip and contained a total of about 950,000 SNPs and custom content to improve the overlap with our V2 array. The V4 platform in current use is a fully custom array of about 950,000 SNPs and includes a lower redundancy subset of V2 and V3 SNPs with additional coverage of lower-frequency coding variation. Samples that failed to reach 98.5% call rate were re-analyzed. Individuals whose analyses failed repeatedly were re-contacted by 23andMe customer service to provide additional samples, as is done for all 23andMe customers.
Participant genotype data were imputed using the March 2012 'v3' release of the 1000 Genomes Project reference haplotypes 67 . We phased and imputed data for each genotyping platform separately. First, we used BEAGLE 68 (version 3.3.1) to phase batches of 8000-9000 individuals across chromosomal segments of no more than 10,000 genotyped SNPs, with overlaps of 200 SNPs. We excluded SNPs with minor allele frequency < 0.001, Hardy-Weinberg equilibrium P < 1 × 10 −20 , call rate <95%, or large allele frequency discrepancies compared to the European 1000 Genomes Project reference data. Frequency discrepancies were identified by computing a 2 × 2 table of allele counts for European 1000 Genomes samples and 2000 randomly sampled 23andMe customers with European ancestry, and identifying SNPs with a χ 2 P < 10 −15 . We imputed each phased segment against all-ethnicity 1000 Genomes haplotypes (excluding monomorphic and singleton sites) using Minimac2 69 with five rounds and 200 states for parameter estimation.
For the non-pseudoautosomal region of the X chromosome, males and females were phased together in segments, treating the males as already phased; the pseudoautosomal regions were phased separately. We then imputed males and females together using Minimac as with the autosomes, treating males as homozygous pseudo-diploids for the non-pseudoautosomal region.
GWAS analysis. For case control comparisons, we tested for association using logistic regression, assuming additive allelic effects. For quantitative traits, association tests were performed using linear regression. For tests using imputed data, we use the imputed dosages rather than best-guess genotypes. We included covariates for age, gender, and the top five principal components to account for residual population structure. The association test P value was computed using a likelihood ratio test, which in our experience is better behaved than a Wald test on the regression coefficient. Results for the X chromosome were computed similarly, with men coded as if they were homozygous diploid for the observed allele.
SNP function annotation. To explore whether any of the significant SNPs identified might link to functional mutations or have potential regulatory functions, we used the online tools HaploReg (http://www.broadinstitute.org/ mammals/haploreg/haploreg.php) to confirm the location of each SNP in relation to annotated protein-coding genes and/or non-coding regulatory elements. We queried only the variants that were within 500 kb of and in moderate LD (r 2 > 0.6) with the index SNP.
Imputation of HLA classical alleles and respective amino acid variations. HLA imputation was performed with HIBAG 70 , an attribute bagging based statistical method that comes as a freely available R package and includes a pre-fit classifier. This classifier was trained from a large database of individuals with known HLA alleles and SNP variation within the HLA region. Over 98% of the tagging SNPs used in HIBAG were genotyped and passed quality control on 23andMe's platform. We imputed allelic dosage for HLA-A, B, C, DPB1, DQA1, DQB1 and DRB1 loci at four-digit resolution. We used the default settings of HIBAG and the run time for 100,000 samples was about 10 h on our cluster.
Using an approach suggested in 71 , we downloaded the files that map HLA alleles to amino acid sequences from https://www.broadinstitute.org/mpg/snp2hla/ and mapped our imputed HLA alleles at four-digit resolution to the corresponding amino acid sequences; in this way we translated the imputed HLA allelic dosages directly to amino acid dosages. We encoded all amino acid variants in the 23andMe European samples as biallelic markers, which facilitated downstream analysis. For example, position 45 of HLA-B protein had five different alleles (E: Glu, M: Met, T: Thr, K: Lys, G: Gly), we first encoded the position using five binary markers, each corresponding to the presence or absence of each allele (e.g., HLA-B Gly45 tags Gly at position 45 of HLA-B protein). For positions having three or more alleles, we also created markers that tag multiple alleles, each corresponding to the presence or absence of the multiple alleles (e.g., HLA-B ThrGly45 tags Thr and Gly at position 45 of HLA-B protein). Thus, we created binary indicators for all possible combinations of amino acid variants. We use this naming convention for amino acid polymorphisms throughout this paper.
We imputed 292 classical HLA alleles at four-digit resolution and 2395 biallelic amino acid polymorphisms. Similar to the SNP imputation, we measured imputation quality using r 2 , which is the ratio of the empirically observed variance of the allele dosage to the expected variance assuming Hardy-Weinberg equilibrium. The imputation quality (r 2 ) of the top associated HLA alleles and amino acids are in Supplementary Table 8.
HLA region fine mapping. To test associations between imputed HLA allele/ amino acid dosages and phenotypes, we performed logistic or linear regression using the same set of covariates used in the SNP-based GWAS for that phenotype. We applied a forward stepwise strategy, within each type of variant, to establish statistically independent signals in the HLA region. Within each variant type (e.g., SNP, HLA allele and HLA amino acid), we first identified the most strongly associated signals (lowest P value) for each disease and performed forward iterative conditional regression to identify other independent signals if the conditional P value was < 5 × 10 −8 . All analyses controlled for sex and five principal components of genetic ancestry. The P values were calculated using a likelihood ratio test. The iterative conditional regression dissected HLA signal into independent HLA associations. Within each identified HLA locus (HLA-A, B, C, DPB1, DQA1, DQB1 and DRB1), we further carried out reciprocal analyses, which are the conditional analyses across variants types, to see if the association can be attributed to the amino acid polymorphism within each HLA locus.
Pathway analysis. To better understand how multiple genes in the same pathway may contribute to certain infections, we performed pathway analysis using MAGENTA 72 , which tests for enrichment of genetic associations in predefined biological processes or sets of functionally related genes, using GWAS results as input. We used gene sets of 1320 canonical pathways from the Molecular Signatures Database (MsigDB) compiled by domain experts 73 and default settings of MAGENTA tool. We reported the most significant pathway (gene sets) with a gene set enrichment FDR < 0.01.
Data availability. The GWAS summary statistics for the top 8000 SNPs for each phenotype are available to download in Supplementary Data 1. The association test statistics for the HLA alleles and HLA amino acid with P value < 0.05 for each phenotype are available to download in Supplementary Data 2. Researchers can request to access the full set of GWAS statistics by applying to 23andMe research collaboration program.