Mapping the human genetic architecture of COVID-19

Niemi, Mari E. K.; Karjalainen, Juha; Liao, Rachel G.; Neale, Benjamin M.; Daly, Mark; Ganna, Andrea; Pathak, Gita A.; Andrews, Shea J.; Kanai, Masahiro; Veerapen, Kumar; Fernandez-Cadenas, Israel; Schulte, Eva C.; Striano, Pasquale; Marttila, Minttu; Minica, Camelia; Marouli, Eirini; Karim, Mohd Anisul; Wendt, Frank R.; Savage, Jeanne; Sloofman, Laura; Butler-Laporte, Guillaume; Kim, Han-Na; Kanoni, Stavroula; Okada, Yukinori; Byun, Jinyoung; Han, Younghun; Uddin, Mohammed Jashim; Smith, George Davey; Willer, Cristen J.; Buxbaum, Joseph D.; Karjalainen, Juha; Mehtonen, Juha; Folkersen, Lasse; Moltke, Ida; Hinney, Anke; Wang, Chen; Ellinghaus, David; Brunak, Søren; Castro, Pedro; GuindoMartinez, Marta; Lee, Ji Yeon; Loos, Ruth J. F.; Zhao, Jing Hua; Sun, Yan V.; Wang, Bo; Parker, Robert; Garcia, Sara Mingo; Smith, Oliver; Davis, Christopher; COVID-19 Host Genetics Initiative; 23andMe COVID-19 Team; Norwegian SARS-CoV-2 Study Grp; Humanitas COVID-19 Task Force; Humanitas Gavazzeni COVID-19 Task; FHoGID; RegCOVID; P-PredictUs; SeroCOVID; CRiPSI; Genes & Hlth Res Team; UCLA Hlth ATLAS Data Mart Working

The COVID-19 pandemic, caused by infection with SARS-CoV-2, has resulted in an enormous health and economic burden worldwide. One of the most remarkable features of SARS-CoV-2 infection is the variation in consequences, which range from asymptomatic to life-threatening, viral pneumonia and acute respiratory distress syndrome 8 . Although established host factors correlate with disease severity (for example, increasing age, being a man and higher body-mass index 1 ), these risk factors alone do not explain all of the variability in disease severity observed across individuals.
Genetic factors contributing to COVID-19 susceptibility and severity may provide new biological insights into disease pathogenesis and identify mechanistic targets for therapeutic development or drug repurposing, as treating the disease remains a highly important goal despite the recent development of vaccines. Further supporting this line of inquiry, rare loss-of-function variants in genes involved in the type I interferon response may be involved in severe forms of COVID-19 [9][10][11] . At the same time, several genome-wide association studies that investigate the contribution of common genetic variation 12-15 to COVID-19 have provided robust support for the involvement of several genomic loci associated with COVID-19 severity and susceptibility, with the strongest and most robust finding for severity being at the 3p21.31 locus [12][13][14][15][16] . However, much remains unknown about the genetic basis of susceptibility to SARS-CoV-2 and severity of COVID-19.
The COVID-19 Host Genetics Initiative (COVID-19 HGI) (https://www. covid19hg.org/) 17 is an international, open-science collaboration to share scientific methods and resources with research groups across the world with the goal to robustly map the host genetic determinants of SARS-CoV-2 infection and the severity of the resulting COVID-19 disease. Here, we report the latest results of meta-analyses of 46 studies from 19 countries (Fig. 1) for COVID-19 host genetic effects. regardless of symptoms (Methods). Controls for all three analyses were selected as genetically ancestry-matched samples without known SARS-CoV-2 infection, if that information was available (Methods). The average age of the participants with COVID-19 across studies was 55 years (Supplementary Table 1). We report quantile-quantile plots in Supplementary Fig. 1 and ancestry principal component plots for contributing studies in Extended Data Fig. 2.
Across our three analyses, we reported a total of 13 independent genome-wide significant loci associated with COVID-19 (the threshold of P < 1.67 × 10 −8 is adjusted for multiple trait testing) (Supplementary Table 2), most of which were shared between two or more COVID-19 phenotypes. Two of these loci are in very close proximity within the 3p21.31 region, which was previously reported as a single locus associated with COVID-19 severity 12-16 (Extended Data Fig. 3). Overall, we find six genome-wide significant associations for critical illness due to COVID-19, using data from 6,179 cases and 1,483,780 controls from 16 studies (Extended Data Fig. 4). Nine genome-wide significant loci were detected for moderate to severe hospitalized COVID-19 (including five of the six critical illness loci) from an analysis of 13,641 cases of COVID-19 and 2,070,709 controls across 29 studies (Fig. 2a, top). Finally, seven loci reached genome-wide significance in the analysis using data for all available 49,562 reported cases of SARS-CoV-2 infection and 1,770,206 controls, using data from a total of 44 studies (Fig. 2a, bottom). The proportion of cases with non-European genetic ancestry for each of the three analyses was 23%, 29% and 22%. We report the results for the lead variants at the 13 loci in different ancestry-group meta-analyses in Supplementary Table 3. We note that two loci, tagged by lead variants rs1886814 and rs72711165, had higher allele frequencies in southeast Asian (rs1886814; 15%) and East Asian genetic ancestry (rs72711165; 8%) whereas the minor allele frequencies in European populations were less than 3%. This highlights the value of including data from diverse populations for genetic discovery. We discuss the replication of previous findings and the new discoveries from these three analyses in the Supplementary Note.

Variant effects on severity and susceptibility
We found no genome-wide significant sex-specific effects at the 13 loci. However, we did identify significant heterogeneous effects (P < 0.004) across studies for 3 out of the 13 loci (Methods), which probably reflects the differential ascertainment of cases (Supplementary Table 2). There was a small number of overlapping samples (n = 8,380 European ancestry; n = 745 East Asian ancestry) between controls from the genOMICC and the UK Biobank studies, but leave-one-out sensitivity analyses did not reveal any bias in the corresponding effect sizes or P values (Extended Data Fig. 5

and Supplementary Information).
We next wanted to better understand whether the 13 significant loci were acting through mechanisms that increased the susceptibility to infection or that affected the progression of symptoms towards more severe disease. For all 13 loci, we compared the lead variant (strongest association P value) odds ratios (ORs) for the risk-increasing allele across our different COVID-19 phenotype definitions.
Focusing on the two better powered analyses: all cases with a reported SARS-CoV-2 infection and all cases hospitalized due to COVID-19, we find that four of the loci have similar odds ratios between these two analyses (Methods and Supplementary Table 2). Such consistency suggests a stronger link to susceptibility to SARS-CoV-2 infection rather than to the development of severe COVID-19. The strongest susceptibility signal was the previously reported ABO locus (rs912805253) 12,13,15,16 . Notably, and in agreement with a previously reported study 15 , we also report a locus within the 3p21.31 region that was more strongly associated with susceptibility to SARS-CoV-2 than progression to more severe COVID-19 phenotypes. rs2271616 showed a stronger association with a reported SARS-CoV-2 infection (P = 1.79 × 10 −34 ; OR (95% confidence interval (CI)) = 1.15 (1.13-1.18)) than hospitalization (P = 1.05 × 10 −5 ; OR (95% CI) = 1.12 (1.06-1.19)). For this locus-which contains additional independent signals-the linkage-disequilibrium (LD) pattern is discordant with the P-value expectation (Extended Data Fig. 6 and Supplementary Note), pointing to a key missing causal variant or to a potentially undiscovered multi-allelic or structural variant in this locus.

Article
By contrast, 9 out of the 13 loci were associated with increased risk of severe symptoms with significantly larger odds ratios for hospitalized COVID-19 compared with the mildest phenotype of reported SARS-CoV-2 infection (eight loci were below the threshold of P < 0.004 (test for effect size difference) and, in addition, the lead variant rs10774671 had a clear increase in odds ratios despite not passing this threshold) (Supplementary Table 2). We further compared the odds ratios for these nine loci for critical illness due to COVID-19 versus hospitalized due to COVID-19, and found that these loci exhibited a general increase in effect risk for critical illness (Methods, Extended Data Fig. 7a and Supplementary Table 4), but the lower power for association analysis of critically ill COVID-19 means that these results should be considered as suggestive. Overall, these results indicated that these nine loci were more likely to be associated with progression of the disease and worse outcome from SARS-CoV-2 infection compared to being associated with susceptibility to SARS-CoV-2 infection.
For some of these analyses, the controls were simply existing population controls without knowledge of SARS-CoV-2 infection or COVID-19 status, which may bias effect size estimates as some of these individuals may have either become infected with SARS-CoV-2 or developed COVID-19. We perform several sensitivity analyses (Extended Data Fig. 7b, Supplementary Note and Supplementary Table 4) in which we show that using population controls can be a valid and powerful strategy for host genetic discovery of infectious disease, and particularly those that are widespread and with rare severe outcomes.
A lung-specific cis-expression quantitative trait loci (cis-eQTLs) from GTEx v.8 20 (n = 515) and the Lung eQTL Consortium 21 (n = 1,103) provided further support for a subset of loci (Supplementary  Table 7   rs2109069:G>A in DPP9 (chr. 19p13.3), which is positively associated with critical illness, was previously reported to be risk-increasing for interstitial lung disease (tag lead variant rs12610495:A>G (p.Leu8Pro); OR = 1.29, P = 2.0 × 10 −12 ) 5 . The COVID-19 lead variant rs1886814:A>C in the FOXP4 locus is correlated (r 2 = 0.64) with a lead variant of lung adenocarcinoma (tag variant is rs7741164; OR = 1.2, P = 6.0 × 10 −13 ) 6,22 and similarly with a lead variant reported for subclinical interstitial lung disease 23 . In severe COVID-19, lung cancer and interstitial lung disease, the minor, expression-increasing allele is associated with increased risk. We also found that intronic variants (chr. 1q22) and rs1819040:T>A in KANSL1 (chr. 17q21.31), associated with protection against hospitalization due to COVID-19, were previously reported for reduced lung function (for example, tag lead variant rs141942982:G>T; OR (95% CI) = 0.96 (0.95-0.97), P = 1.00 × 10 −20 ) 7 . Notably, the 17q21.31 locus is a well-known locus for structural variants containing a megabase inversion polymorphism (H1 and inverted H2 forms) and complex copy-number variations, in which the inverted H2 forms were shown to be positively selected in European individuals 24,25 . Lastly, there are two loci in the 3p21.31 region with varying genes prioritized by different methods for different independent signals. For the severity lead variant rs10490770:T>C, we prioritized CXCR6 with the Variant2Gene (V2G) algorithm 26 , although LZTFL1 is the closest gene. The CXCR6 has a role in chemokine signalling 27 and LZTFL1 has been implicated in lung cancer 28 . rs2271616:G>T, which is associated with susceptibility, tags a complex region including several independent signals (Supplementary Note) that are all located within the gene body of SLC6A20, which encodes a protein that is known to functionally interact with the SARS-CoV-2 receptor ACE2 29 . However, none of the lead variants in the 3p21.31 region has been previously associated with other traits or diseases in our PheWAS analysis. Although these results provide supporting in silico evidence for candidate causal gene prioritization, further functional characterization is needed. Detailed locus descriptions and LocusZoom plots are provided in Supplementary Fig. 2.

Polygenic architecture of COVID-19
To further investigate the genetic architecture of COVID-19, we used results from meta-analyses including samples from European ancestries (sample sizes are described in the Methods and Supplementary Table 1) to estimate the heritability explained by common single-nucleotide polymorphisms-that is, the proportion of variation in the two phenotypes that was attributable to common genetic variants-and to determine whether heritability of COVID-19 phenotypes was enriched in genes that were specifically expressed in certain tissues 30 from the GTEx dataset 31 . We detected low, but significant, heritability across all three analyses (<1% on observed scale, all P values were P < 0.0001) (Supplementary Table 8). The values are low compared to previously published studies 14 , but may be explained by differences in the reported estimate scale (observed versus liability), the specific method used, disease-prevalence estimates, phenotypic differences between patient cohorts or ascertainment of controls. Despite the low reported values, we found that heritability of a reported SARS-CoV-2 infection was significantly enriched in genes that were specifically expressed in the lung (P = 5.0 × 10 −4 ) (Supplementary Table 9). These findings, together with the genome-wide significant loci identified in the meta-analyses, suggest that there is a significant polygenic architecture that can be better leveraged with future, larger, sample sizes.
We found evidence (false-discovery rate (FDR) < 0.05) of significant genetic correlations between nine traits and hospitalized COVID-19 and reported SARS-CoV-2 infection (Fig. 3, Extended Data Fig. 9 and Supplementary Table 11). Notably, genetic liability to ischaemic as a medium square and P > 0.5 as a small square. Genetic correlations or causal estimates that are significantly different from zero at an FDR of 5% are marked with an asterisk. Two-sided P values were calculated using LDSC for genetic correlations and inverse-variance-weighted analysis for Mendelian randomization. ADHD, attention-deficit hyperactivity disorder; BMI, body mass index; CRP, C-reactive protein; eGFR, estimated glomerular filtration rate.
Article stroke was only significantly positively correlated with critical illness or hospitalization due to COVID-19, but not with a higher likelihood of reported SARS-CoV-2 infection (infection r g = 0.019 versus hospitalization r g = 0.41, z = 2.7, P = 0.006; infection r g = 0.019 versus critical illness r g = 0.40, z = 2.49, P = 0.013). We next used two-sample Mendelian randomization to infer potentially causal relationships between these traits. After correcting for multiple testing (FDR < 0.05), eight exposure-COVID-19 trait pairs showed suggestive evidence of a causal association (Fig. 3, Extended Data Fig. 10, Supplementary Table 12 and Supplementary Fig. 3). Five of these associations were robust to potential violations of the underlying assumptions of Mendelian randomization. Corroborating our genetic correlation results and evidence from epidemiological studies, genetically predicted higher body-mass index (OR (95% CI) = 1.4 (1.3-1.6), P = 8.5 × 10 −11 ) and smoking (OR (95% CI) = 1.9 (1.3-2.8), P = 0.0012) were associated with increased risk of COVID-19 hospitalization, with body-mass index also being associated with increased risk of SARS-CoV-2 infection (OR (95% CI) = 1.1 (1.1-1.2), P = 4.8 × 10 −7 ). Genetically predicted increased height (OR (95% CI) = 1.1 (1-1.1)), P = 8.9 × 10 −4 ) was associated with an increased risk of reported SARS-CoV-2 infection, whereas a genetically predicted higher red-blood-cell count (OR (95% CI) = 0.93 (0.89-0.96), P = 5.7 × 10 −5 ) was associated with a reduced risk of reported SARS-CoV-2 infection. Despite evidence of a genetic correlation between type II diabetes and COVID-19 outcomes, there was no evidence of a causal association in the Mendelian randomization analyses, which suggests that the observed genetic correlations are due to pleiotropic effects between body-mass index and type 2 diabetes. Further sensitivity analyses relating to sample overlap are discussed in the Supplementary Information.

Discussion
The COVID-19 HGI has brought together investigators from across the world to advance genetic discovery for SARS-CoV-2 infection and severe COVID-19 disease. We report 13 genome-wide significant loci associated with some aspect of SARS-CoV-2 infection or COVID-19. Many of these loci overlap with previously reported associations with lung-related phenotypes or autoimmune or inflammatory diseases, but some loci have no obvious candidate gene.
Four out of the thirteen genome-wide significant loci showed similar effects in the reported SARS-CoV-2 infection analysis (a proxy for disease susceptibility) and all-hospitalized COVID-19 (a proxy for disease severity). Of these, one locus was in close proximity to, yet independent of, the major genetic signal for COVID-19 severity at the 3p21.31 locus. Notably, this locus was associated with COVID-19 susceptibility rather than severity. The locus overlaps SLC6A20, which encodes an amino acid transporter that interacts with ACE2. Nonetheless, we caution that more data are needed to resolve the nature of the relationship between genetic variation and COVID-19 at this locus, particularly as the physical proximity, LD structure and patterns of association suggest that untagged genetic variation could drive the association signal in the region. Our findings support the notion that some genetic variants, most notably at the ABO and PPP1R15A loci, in addition to SLC6A20, can indeed affect susceptibility to infection rather than progression to severe COVID-19 once infected.
Several of the loci reported here-as noted in previous publications 12,14 -intersect with well-known genetic variants that have established genetic associations. Examples of these include variants at DPP9 and FOXP4, which show previous evidence of increasing risk for interstitial lung disease 5 , and missense variants within TYK2 that show a protective effect on several autoimmune-related diseases 32-35 . Together with the heritability enrichment observed in genes expressed in lung tissues, these results highlight the involvement of lung-related biological pathways in the development of severe COVID-19. Several other loci show no previously documented genome-wide significant associations, despite the high significance and attractive candidate genes for COVID-19 (for example, CXCR6, LZTFL1, IFNAR2 and OAS1/ OAS2/OAS3 loci). The previously reported associations for the strongest association for COVID-19 severity at the 3p21.31 locus and monocytes count are likely to be due to proximity and not a true co-localization.
Increasing the global representation in genetic studies enhances the ability to detect novel associations. Two of the loci that affect disease severity were only discovered by including the four studies of individuals with East Asian ancestry. One of these loci-close to FOXP4-is common particularly in East Asian participants (32%) as well as admixed American participants in the Americas (20%) and Middle Eastern participants (7%), but has a low frequency in most European ancestries (2-3%) in our data. Although we cannot be certain of the mechanism of action, the FOXP4 association is an attractive biological target, as it is expressed in the proximal and distal airway epithelium 36 and has been shown to have a role in controlling epithelial cell fate during lung development 37 . The COVID-19 HGI continues to pursue expansion of the datasets included in the analyses of the consortium to populations from underrepresented populations in upcoming data releases. We plan to release ancestry-specific results in full once the sample sizes allow for a well-powered meta-analysis.
Care should be taken when interpreting the results from a meta-analysis because of challenges with case and control ascertainment and collider bias (see Supplementary Note for a more detailed discussion on study limitations). Drawing a comprehensive and reproducible map of the host genetics factors associated with COVID-19 severity and SARS-CoV-2 requires a sustained international effort to include diverse ancestries and study designs. To accelerate downstream research and therapeutic discovery, the COVID-19 HGI regularly publishes meta-analysis results from periodic data freezes on the website https://www.covid19hg.org/ and provides an interactive explorer through which researchers can browse the results and the genomic loci in more detail. Future work will be required to better understand the biological and clinical value of these findings. Continued efforts to collect more samples and detailed phenotypic data should be endorsed globally, allowing for more thorough investigation of variable, heritable symptoms, particularly in light of the newly emerging strains of SARS-CoV-2, which may provoke different host responses that lead to disease.

Online content
Any methods, additional references, Nature Research reporting summaries, source data, extended data, supplementary information, acknowledgements, peer review information; details of author contributions and competing interests; and statements of data and code availability are available at https://doi.org/10.1038/s41586-021-03767-x.

Contributing studies
All of the participants were recruited following protocols approved by local Institutional Review Boards; this information is collected in Supplementary Table 1 for all 46 studies. All protocols followed local ethics recommendations and informed consent was obtained when required. Information about sample numbers, sex and age from for each contributing study is given in Supplementary Table 1. In total, 16 studies contributed data to the analysis of critical illness due to COVID-19, 29 studies contributed data to hospitalized COVID-19 analysis and 44 studies contributed to the analysis of all cases of COVID-19. Each individual study that contributed data to a particular analysis met a minimum threshold of 50 cases, as defined by the phenotypic criteria, for statistical robustness. The effective sample sizes for each ancestry group shown in Fig. 1 were calculated for display using the formula: ((4 × N case × N control )/(N case + N control )). Details of contributing research groups are provided in Supplementary Table 1.

Phenotype definitions COVID-19 disease status (critical illness and hospitalization status) was assessed following the Diagnosis and Treatment Protocol for Novel
Coronavirus Pneumonia 38 . The critically ill COVID-19 group included patients who were hospitalized owing to symptoms associated with laboratory-confirmed SARS-CoV-2 infection and who required respiratory support or whose cause of death was associated with COVID-19. The hospitalized COVID-19 group included patients who were hospitalized owing to symptoms associated with laboratory-confirmed SARS-CoV-2 infection. The reported SARS-CoV-2 infection group included individuals with laboratory-confirmed SARS-CoV-2 infection or electronic health record, ICD coding or clinically confirmed COVID-19, or self-reported COVID-19 (for example, by questionnaire), with or without symptoms of any severity. Genetic-ancestry-matched control individuals for the three case definitions were sourced from population-based cohorts, including individuals whose exposure status to SARS-CoV-2 was either unknown or infection-negative for questionnaire/electronic-health-record-based cohorts. Additional information regarding individual studies contributing to the consortium are described in Supplementary Table 1.

Genome-wide association studies and meta-analyses
Each contributing study genotyped the samples and performed quality controls, data imputation and analysis independently, but following the consortium recommendations (information is available at https:// www.covid19hg.org/). We recommended that genome-wide association study (GWAS) analyses were run using Scalable and Accurate Implementation of GEneralized mixed model (SAIGE) 39 on chromosomes 1-22 and X. The recommended analysis tool was SAIGE, but studies also used other software such as PLINK 40 . The suggested covariates were age, age 2 , sex, age × sex and the 20 first principal components. Any other study-specific covariates to account for known technical artefacts could be added. SAIGE automatically accounts for sample relatedness and case-control imbalances. Quality-control and analysis approaches for individual studies are reported in Supplementary Table 1.
Study-specific summary statistics were then processed for meta-analysis. Potential false positives, inflation and deflation were examined for each submitted GWAS. Allele frequency plots against gnomAD 3.0 genomes were manually inspected for each study. Standard error values as a function of the effective sample size were used to find studies that deviated from the expected trend. Summary statistics passing this manual quality control were included in the meta-analysis. Variants with an allele frequency of >0.1% and an imputation INFO score of >0.6 were carried forward from each study. Variants and alleles were lifted over to genome build GRCh38, if needed, and harmonized to gnomAD 3.0 genomes 41 by finding matching variants by strand flipping or switching the ordering of alleles. If multiple matching variants were included, the best match was chosen according to the minimum fold change in absolute allele frequency. Meta-analysis was performed using the inverse-variance-weighted (IVW) method on variants that were present in at least two-thirds of the studies contributing to the phenotype analysis. The method summarizes effect sizes across the multiple studies by computing the mean of the effect sizes weighted by the inverse variance in each individual study.
We report 13 meta-analysis variants that pass the genome-wide significance threshold after adjusting the threshold for multiple traits tested (P < 5 × 10 −8 /3). We report the unadjusted P values for each variant. We tested for heterogeneity between estimates from contributing studies using Cochran's Q-test 42,43 . This is calculated for each variant as the weighted sum of squared differences between the effects sizes and their meta-analysis effect, the weights being the inverse variance of the effect size. Q is distributed as a χ 2 statistic with k (number of studies) minus one degrees of freedom. Two loci reached genome-wide significance but were excluded from the significant results in Supplementary Table 2 due to heterogeneity between estimates from contributing studies and missingness between studies at chr. 6: 31057940-31380334 and chr. 7: 54671568-54759789; however, these regions are not excluded from the corresponding summary statistics in data release 5 (COVID-19 HGI (https://www.covid19hg. org/results/r5/) and GWAS Catalog (study code GCST011074)). For each of the lead variants reported in Supplementary Table 2, we aimed to find loci specific to susceptibility or severity by testing whether there was heterogeneity between the effect sizes associated with hospitalized COVID-19 (progression to severe disease) and reported SARS-CoV-2 infection. We used the Cochran's Q measure 42,43 , calculated for each variant as the weighted sum of squared differences between the two analysis effect sizes and their meta-analysis effect with the weights being the inverse variance of the effect size. A significant P value of P < 0.004 ((0.05/13 loci) for multiple tests) indicates that the effect sizes for a particular variant are significantly different in the two analyses (Supplementary Table 2). For the nine loci, in which the lead variant effect size was significantly higher for hospitalized COVID-19, we carried out the same test again but comparing effect sizes from hospitalized COVID-19 with critically ill COVID-19 (Supplementary Table 4). Furthermore, we carried out the same test comparing meta-analysed hospitalized COVID-19 (population as controls) and hospitalized COVID-19 (SARS-CoV-2-positive but non-hospitalized as controls) (Supplementary Table 4). For these pairs of phenotype comparisons, we generated new meta-analysis summary statistics to use; including only those studies that could contribute data to both phenotypes that were under comparison.

Principal component projection
To project every GWAS participant into the same principal component (PC) space, we used pre-computed PC loadings and reference allele frequencies. For reference, we used unrelated samples from the 1000 Genomes Project and the Human Genome Diversity Project and computed PC loadings and allele frequencies for the 117,221 single-nucleotide polymorphisms (SNPs) that (1) are available in every cohort; (2) have a minor allele frequency of >0.1% in the reference; and (3) are LD-pruned (r 2 < 0.8; 500-kb window). We then asked each cohort to project their samples using our automated script provided at https:// github.com/covid19-hg/. It internally uses the PLINK2 44 --score function with the variance-standardize option and reference allele frequencies (--read-freq); so that each cohort-specific genotype/dosage matrix is mean-centred and variance-standardized with respect to reference allele frequencies, but not cohort-specific allele frequencies. We further normalized the projected PC scores by dividing the values by a square root of the number of variants used for projection to account for a subtle difference due to missing variants.

Gene prioritization
To prioritize candidate causal genes reported in full in Supplementary Table 2, we used various gene prioritization approaches using both locus-based and similarity-based methods. Because we only describe the in silico gene prioritization results without characterizing the actual functional activity in vitro or in vivo, we aimed to provide a systematic approach to nominate potential causal genes in a locus using the following criteria.
(1) The closest gene: a gene that is closest to a lead variant by distance to the gene body.
(2) Genes in the LD region: genes that overlap with a genomic range containing any variants in LD (r 2 > 0.6) with a lead variant. For LD computation, we retrieved LD matrices provided by gnomAD v.2.1.1 41 for each population analysed in this study (except for admixed American, Middle Eastern and South Asian genetic ancestry populations, for whom data are not available). We then constructed a weighted-average LD matrix by per-population sample sizes in each meta-analysis, which we used as a LD reference.
(3) Genes with coding variants: genes with at least one loss-of-function or missense variant (annotated by VEP 45 v.95 with GENCODE v.29) that is in LD with a lead variant (r 2 > 0.6).
(4) eGenes: genes with at least one fine-mapped cis-eQTL variant (PIP > 0.1) that is in LD with a lead variant (r 2 > 0.6) (Supplementary Table 5). We retrieved fine-mapped variants from the GTEx v.8 20 (https:// www.finucanelab.org/) and eQTL catalogue 46 . In addition, we looked up significant associations in the Lung eQTL Consortium 21 (n = 1,103) to further support our findings in lung with a larger sample size (Supplementary Table 7). We note that, in contrast to the GTEx or eQTL catalogue, we only looked at associations and did not fine-map our data to the Lung eQTL Consortium data.
(5) V2G: a gene with the highest overall V2G score based on Open Targets Genetics (OTG) 26 . For each variant, the overall V2G score aggregates differentially weighted evidence of variant-gene associations from several data sources, including molecular cis-QTL data (for example, cis-protein QTLs from ref. 47 , cis-eQTLs from GTEx v.7 and so on), interaction-based datasets (for example, promoter capture Hi-C), genomic distance and variant effect predictions (VEP) from Ensembl. A detailed description of the evidence sources and weights used is provided in the OTG documentation (https://genetics-docs.opentargets. org/our-approach/data-pipeline) 26 .

Phenome-wide association study
To investigate the evidence of shared effects of 15 index variants for COVID-19 and previously reported phenotypes, we performed a phenome-wide association study. We considered phenotypes in OTG obtained from the GWAS catalogue (this included studies with and without full summary statistics, n = 300 and 14,013, respectively) 48 and from the UK Biobank. Summary statistics for UK Biobank traits were extracted from SAIGE 39 for binary outcomes (n = 1,283 traits) and Neale v.2 (n = 2,139 traits) for both binary and quantitative traits (http://www.nealelab.is/uk-biobank/) and FinnGen Freeze 4 cohort (https://www.finngen.fi/en/access_results). We report PheWAS results for phenotypes for which the lead variants were in high LD (r 2 > 0.8) with the 13 genome-wide significant lead variants from our main COVID-19 meta-analysis (Supplementary Table 6). This conservative approach allowed spurious signals primarily driven by proximity rather than actual colocalization to be removed (see Methods).
To remove plausible spurious associations, we retrieved phenotypes for GWAS lead variants that were in LD (r 2 > 0.8) with COVID-19 index variants.

Heritability
LD score regression v.1.0.1 49 was used to estimate the SNP heritability of the phenotypes from the meta-analysis summary statistic files. As this method depends on matching the LD structure of the analysis sample to a reference panel, the summary statistics of European ancestry only were used. Sample sizes were n = 5,101 critically ill cases of COVID-19 and n = 1,383,241 control participants, n = 9,986 hospitalized cases of COVID-19 and n = 1,877,672 control participants, and n = 38,984 cases and n = 1,644,784 control participants for the analysis of all casesall including the 23andMe cohort. Pre-calculated LD scores from the 1000 Genomes European reference population were obtained online (https://data.broadinstitute.org/alkesgroup/LDSCORE/). Analyses were conducted using the standard program settings for variant filtering (removal of non-HapMap3 SNPs, the HLA region on chromosome 6, non-autosomal, χ 2 > 30, minor allele frequency of <1%, or allele mismatch with reference). We additionally report SNP heritability estimates for the all-ancestries meta-analyses, calculated using European panel LD scores, in Supplementary Table 8.

Partitioned heritability
We used partitioned LD score regression 50 to partition COVID-19 SNP heritability in cell types in our summary statistics for European ancestry only. We ran the analysis using the baseline model LD scores calculated for European populations and regression weights that are available online (https://github.com/bulik/ldsc). We used the COVID-19 summary statistics for European ancestry only for the analysis.

Genome-wide association summary statistics
We obtained genome-wide association summary statistics for 43 complex-disease, neuropsychiatric, behavioural or biomarker phenotypes (Supplementary Table 10). These phenotypes were selected based on their putative relevance to COVID-19 susceptibility, severity or mortality, with 19 selected based on the Centers for Disease Control list of underlying medical conditions associated with COVID-19 severity 51 or traits reported to be associated with increased risk of COVID-19 mortality by OpenSafely 52 . Summary statistics generated from GWAS using individuals of European ancestry were preferentially selected if available. These summary statistics were used in subsequent genetic correlation and Mendelian randomization analyses.

Genetic correlation
LD score regression 50 was also used to estimate the genetic correlations between our COVID-19 meta-analysis phenotypes reported using samples of only European ancestry, and between these and the curated set of 38 summary statistics. Genetic correlations were estimated using the same LD score regression settings as for heritability calculations. Differences between the observed genetic correlations of SARS-CoV-2 infection and COVID-19 severity were compared using a z-score method 53 .

Mendelian randomization
Two-sample Mendelian randomization was used to evaluate the potential for causal association of the 38 traits on COVID-19 hospitalization, on COVID-19 severity and reported SARS-CoV-2 infection using samples of only European ancestry. Independent genome-wide significant SNPs robustly associated with the exposures of interest (P < 5 × 10 −8 ) were selected as genetic instruments by performing LD clumping using PLINK 40 . We used a strict r 2 threshold of 0.001, a 10-Mb clumping window, and the European reference panel from the 1000 Genomes Project 54 to discard SNPs in LD with another variant with a smaller P-value association. For genetic variants that were not present in the hospitalized COVID-19 analysis, PLINK was used to identify proxy variants that were in LD (r 2 > 0.8). Next, the exposure and outcome datasets were harmonized using the R package TwoSampleMR 55 . Namely, we ensured that the effect of a variant on the exposure and outcome corresponded to the same allele, we inferred positive-strand alleles and dropped palindromes with ambiguous allele frequencies, as well as incompatible alleles. Supplementary Table 10 includes the harmonized datasets used in the analyses.
The global test from Mendelian randomization pleiotropy residual sum and outlier (MR-PRESSO) 56 software was used to investigate overall horizontal pleiotropy. In brief, the standard IVW meta-analytic framework was used to calculate the average causal effect by excluding each genetic variant used to instrument the analysis. A global statistic was calculated by summing the observed residual sum of squares, that is, the difference between the effect predicted by the IVW slope excluding the SNP, and the observed effect of the SNP on the outcome. Overall horizontal pleiotropy was subsequently analysed by comparing the observed residual sum of squares, with the residual sum of squares expected under the null hypothesis of no pleiotropy. The MR-PRESSO global test was shown to perform well when the outcome and exposure GWASs are not disjoint (although the power to detect horizontal pleiotropy is slightly reduced by complete sample overlap). We also used the regression intercept in MR-Egger 57 to evaluate potential bias due to directional pleiotropic effects. This additional check was used in Mendelian randomization analyses with an I GX 2 index surpassing the recommended threshold (I > 90% GX 2 ) 58 . Contingent on the MR-PRESSO global test results we analysed the causal effect of each exposure on COVID-19 hospitalization by using a fixed-effect IVW meta-analysis as the primary analysis, or, if pleiotropy was present, the MR-PRESSO outlier-corrected test. The IVW approach estimates the causal effect by aggregating the single-SNP causal effects (obtained using the ratio of coefficients method-that is, the ratio of the effect of the SNP on the outcome over the effect of the SNP on the exposure) in a fixed-effects meta-analysis. The SNPs were assigned weights based on their inverse variance. The IVW method confers the greatest statistical power for estimating causal associations 59 , but assumes that all variants are valid instruments and can produce biased estimates if the average pleiotropic effect differs from zero. Alternatively, when horizontal pleiotropy was present, we used the MR-PRESSO outlier-corrected method to correct the IVW test by removing outlier SNPs. We conducted further sensitivity analyses using alternative Mendelian randomization methods that provide consistent estimates of the causal effect even when some instrumental variables are invalid, at the cost of reduced statistical power including: (1) Weighted median estimator (WME); (2) weighted mode-based estimator (WMBE); and (3) MR-Egger regression. Robust causal estimates were defined as those that were significant at an FDR of 5% and either (1) showed no evidence of heterogeneity (MR-PRESSO global test P > 0.05) or horizontal pleiotropy (Egger intercept P > 0.05); or (2) in the presence of heterogeneity or horizontal pleiotropy, the WME-, WMBE-, MR-Egger-or MR-PRESSO-corrected estimates were significant (P < 0.05). All statistical analyses were conducted using R v.4.0.3. Mendelian randomization analysis was performed using the 'TwoSampleMR' v.0.5.5 package 55 .

Website and data distribution
In anticipation of the need to coordinate many international partners around a single meta-analysis effort, we created the COVID-19 HGI website (https://covid19hg.org). We were able to centralize information, recruit partner studies, rapidly distribute summary statistics and present preliminary interpretations of the results to the public. Open meetings are held on a monthly basis to discuss future plans and new results; video recordings and supporting documents are shared (https://covid19hg.org/meeting-archive). This centralized resource provides a conceptual and technological framework for organizing global academic and industry groups around a shared goal. The website source code and additional technical details are available at https:// github.com/covid19-hg/covid19hg.
To recruit new international partner studies, we developed a workflow in which new studies are registered and verified by a curation team (https://covid19hg.org/register). Users can explore the registered studies using a customized interface to find and contact studies with similar goals or approaches (https://covid19hg.org/partners). This helps to promote organic assembly around focused projects that are adjacent to the centralized effort (https://covid19hg.org/projects). Visitors can query study information, including study design and research questions. Registered studies are visualized on a world map and are searchable by institutional affiliation, city and country.
To encourage data sharing and other forms of participation, we created a rolling acknowledgements page (https://covid19hg.org/ acknowledgements) and directions on how to contribute data to the central meta-analysis effort (https://covid19hg.org/data-sharing). Upon the completion of each data freeze, we post summary statistics, plots and sample size breakdowns for each phenotype and contributing cohort (https://covid19hg.org/results). The results can be explored using an interactive web browser (https://app.covid19hg.org). Several computational research groups carry out follow-up analyses, which are made available for download (https://covid19hg.org/in-silico). To enhance scientific communication to the public, preliminary results are described in blog posts by the scientific communications team and shared on Twitter. The first post was translated to 30 languages with the help of 85 volunteer translators. We compile publications and preprints submitted by participating groups and summarize genome-wide significant findings from these publications (https://covid19hg.org/ publications).

Reporting summary
Further information on research design is available in the Nature Research Reporting Summary linked to this paper.

Data availability
Summary statistics generated by the COVID-19 HGI are available at https://www.covid19hg.org/results/r5/ and are available in the GWAS Catalog (study code GCST011074). The analyses described here include the freeze-5 data. COVID-19 HGI continues to regularly release new data freezes. Summary statistics for non-European ancestry samples are not currently available due to the small individual sample sizes of these groups, but results for lead variants of 13 loci are reported in Supplementary Table 3. Individual level data can be requested directly from contributing studies, listed in Supplementary Table 1. We used publicly available data from GTEx (https://gtexportal.org/home/), the Neale lab (http://www.nealelab.is/uk-biobank/), Finucane lab (https:// www.finucanelab.org), the FinnGen Freeze 4 cohort (https://www. finngen.fi/en/access_results) and the eQTL catalogue release 3 (http:// www.ebi.ac.uk/eqtl/).

Code availability
The code for summary statistics lift-over, the projection PCA pipeline including precomputed loadings and meta-analyses are available on GitHub (https://github.com/covid19-hg/) and the code for the Mendelian randomization and genetic correlation pipeline is available on GitHub at https://github.com/marcoralab/MRcovid.

Corresponding author(s): Andrea Ganna
Last updated by author(s): Jun 1, 2021 Reporting Summary Nature Research wishes to improve the reproducibility of the work that we publish. This form provides structure for consistency and transparency in reporting. For further information on Nature Research policies, see our Editorial Policies and the Editorial Policy Checklist.

Statistics
For all statistical analyses, confirm that the following items are present in the figure legend, table legend, main text, or Methods section.

n/a Confirmed
The exact sample size (n) for each experimental group/condition, given as a discrete number and unit of measurement A statement on whether measurements were taken from distinct samples or whether the same sample was measured repeatedly The statistical test(s) used AND whether they are one-or two-sided Only common tests should be described solely by name; describe more complex techniques in the Methods section.
A description of all covariates tested A description of any assumptions or corrections, such as tests of normality and adjustment for multiple comparisons A full description of the statistical parameters including central tendency (e.g. means) or other basic estimates (e.g. regression coefficient) AND variation (e.g. standard deviation) or associated estimates of uncertainty (e.g. confidence intervals) For null hypothesis testing, the test statistic (e.g. F, t, r) with confidence intervals, effect sizes, degrees of freedom and P value noted

Software and code
Policy information about availability of computer code Data collection No code was used to collect data in the study.

Data analysis
Each individual study that contributed genetic-phenotype association summary statistics to the consortium carried out their association analyses independently of the consortium (study-specific information outlined in Supplementary Table 1). However, the consortium did release phenotyping and analysis guidelines as a recommendation (https://www.covid19hg.org/). For quality control of genotype data we recommended using the Ricopili pipeline (PMID: 31393554). For genotype phasing and imputation we recommended the TopMed Imputation Server (PMID: 27571263) or Michigan Imputation Server (PMID: 27571263). For genome-wide association study (GWAS), we recommended SAIGE (PMID: 30104761), but some studies used PLINK (PMID: 17701901). Each study then submitted their GWAS summary statistics to the consortium for meta-analysis. Code availability statement: The code for summary statistics liftover, projection PCA pipeline including precomputed loadings and metaanalysis are available at https://github.com/covid19-hg/ and the code for Mendelian randomization and genetic correlation pipeline at https://github.com/marcoralab/MRcovid. For manuscripts utilizing custom algorithms or software that are central to the research but not yet described in published literature, software must be made available to editors and reviewers. We strongly encourage code deposition in a community repository (e.g. GitHub). See the Nature Research guidelines for submitting code & software for further information.

April 2020
Data Policy information about availability of data All manuscripts must include a data availability statement. This statement should provide the following information, where applicable: -Accession codes, unique identifiers, or web links for publicly available datasets -A list of figures that have associated raw data -A description of any restrictions on data availability Data availability statement: Summary statistics generated by COVID-19 HGI are available at https://www.covid19hg.org/results/r5/ and are available on GWAS Catalog (study code GCST011074). The analyses described here utilize the freeze 5 data. COVID-19 HGI continues to regularly release new data freezes. Summary statistics for non-European ancestry samples are not currently available due to the small individual sample sizes of these groups, but results for 13 loci lead variants are reported in Supplementary Table 3. Individual level data can be requested directly from contributing studies, listed in Supplementary Table 1. We used publicly available data from GTEx (https://gtexportal.org/home/), the Neale lab (http://www.nealelab.is/uk-biobank/), Finucane lab (https://www.finucanelab.org), FinnGen Freeze 4 cohort (https://www.finngen.fi/en/access_results), and eQTL catalogue release 3 (http://www.ebi.ac.uk/eqtl/).

Field-specific reporting
Please select the one below that is the best fit for your research. If you are not sure, read the appropriate sections before making your selection.

Life sciences Behavioural & social sciences Ecological, evolutionary & environmental sciences
For a reference copy of the document with all sections, see nature.com/documents/nr-reporting-summary-flat.pdf

Life sciences study design
All studies must disclose on these points even when the disclosure is negative.

Sample size
The consortium meta-analysed genome-wide association study (GWAS) summary statistics from any individual study that had included a minimum of n=50 cases and n=50 controls in their analysis. The cutoff at n=50 cases and n=50 controls was aimed at reducing noise to the meta-analysis, but also to be inclusive of studies that had not yet accumulated large numbers of COVID-19 patient data. No statistical calculation for adequate sample size was performed, but the results identifying multiple genomic regions at genome-wide significance threshold indicates adequate power for genetic discovery.
Data exclusions Individual level phenotype and genotype data exclusions were performed by each individual study, following the consortium analysis plan recommendations (www.covid19hg.org). Possible reasons for sample exclusion included removing genetic ancestry outliers within a study (using principal components analysis), poor quality of genetic data or lack of phenotypic data for a sample.
The consortium manually examined GWAS summary statistics data submitted by each study (for each submitted analysis separately), including sample size used for analysis, allele frequency check against gnomad reference panel, and distribution of test statistics. After metaanalysis, the results were checked for heterogeneity variant effects between contributing studies, and Table 1 excludes two genome-wide significant loci that were deemed to have extremely heterogeneous effects, but these variants are reported in the released consortium summary statistics (with heterogeneity test values).

Replication
No replication was performed. The consortium meta-analysed GWAS summary statistics, bringing together as many studies as possible to achieve the largest possible sample size and statistical power for association. this meant that the consortium included most large studies of COVID-19 host genetics that have been performed to date, so it was not possible to perform replication analyses in external cohorts. Therefore we performed manual checks on each study contributing summary statistics before entering them into the meta-analysis. In addition, after meta-analysis, we performed a check for heterogeneity between variant association estimates across studies contributing data. This allowed us to better understand whether the variant effects differed much between individual studies.
Randomization No randomization was performed because there was no allocation of samples to experimental groups.

Blinding
Blinding was not relevant to the study. The case status and severity of symptoms was evaluated for each sample by investigators from each study respectively. The consortium recommended using covariates to control for confounding: age + age2 + sex + age*sex + 20 principal components (obtained using genetic data) + study specific covariates (if any). The consortium meta-analysed summary statistics from these case/control studies, not individual level data. Details of which variables each study used and how the calculated PCs for their analysis are available in Supplementary Table 1. Reporting for specific materials, systems and methods We require information from authors about some types of materials, experimental systems and methods used in many studies. Here, indicate whether each material, system or method listed is relevant to your study. If you are not sure if a list item applies to your research, read the appropriate section before selecting a response.