Summary
Asthma is a complex disease that affects millions of people and varies in prevalence by an order of magnitude across geographic regions and populations. However, the extent to which genetic variation contributes to these disparities is unclear, as studies probing the genetics of asthma have been primarily limited to populations of European (EUR) descent. As part of the Global Biobank Meta-analysis Initiative (GBMI), we conducted the largest genome-wide association study of asthma to date (153,763 cases and 1,647,022 controls) via meta-analysis across 18 biobanks spanning multiple countries and ancestries. Altogether, we discovered 180 genome-wide significant loci (p < 5×10−8) associated with asthma, 49 of which are not previously reported. We replicate well-known associations such as IL1RL1 and STAT6, and find that overall the novel associations have smaller effects than previously-discovered loci, highlighting our substantial increase in statistical power. Despite the considerable range in prevalence among biobanks, from 3% to 24%, the genetic effects of associated loci are largely consistent across the biobanks and ancestries. To further investigate the polygenic architecture of asthma, we construct polygenic risk scores (PRS) using a multi-ancestry approach, which yields higher predictive power for asthma in non-EUR populations compared to PRS derived from previous asthma meta-analyses and using other methods. Additionally, we find considerable genetic overlap between asthma and chronic obstructive pulmonary disease (COPD) across ancestries but minimal overlap in enriched biological pathways. Our work underscores the multifactorial nature of asthma development and offers insight into the shared genetic architecture of asthma that may be differentially perturbed by environmental factors and contribute to variation in prevalence.
Introduction
Asthma is a complex and multifactorial disease that affects millions of people worldwide, yet much of its genetic architecture has eluded discovery. Genetic factors contribute substantially to asthma risk, with heritability estimates ranging from 35% to 95% (Ober and Yao, 2011; 2015). The most recent meta-analysis of asthma discovered 212 asthma-associated loci across the genome, confirming the polygenic nature of asthma (Han et al., 2020). However, these risk loci only account for a small proportion of the total heritability for asthma. Furthermore, the discovery GWAS, like the majority of previous asthma GWAS, were primarily conducted in populations of European ancestry. A major exception is the Trans-National Asthma Genetic Consortium (TAGC) which included modest sample sizes from populations of African, Japanese and Latino ancestries (Demenais et al., 2018). Surveys of asthma worldwide have found that prevalence can vary by as much as 21-fold among countries (Asher et al., 2020; Sembajwe et al., 2010). Prevalence varies greatly within countries as well (Akinbami et al., 2012) and this variation cannot be attributed to any single known risk factor such as air pollution. Rather, genetic and environmental contributions to the wide variation in asthma prevalence are complex. Therefore, assessing the genetic architecture of asthma in diverse cohorts is particularly important.
This heterogeneity in prevalence is mirrored by, and may be a consequence of, the heterogeneity of the disease itself. Asthma is commonly viewed as not a single disease, but rather a syndrome composed of overlapping phenotypes and endotypes ((Borish and Culp, 2008; Dharmage et al., 2019)). These are driven by a myriad of risk factors, both genetic and non-genetic (Borish and Culp, 2008; Dharmage et al., 2019). Asthma also shares genetic components with various comorbid diseases, including other respiratory diseases like chronic obstructive pulmonary disease (COPD), allergic diseases, obesity, and neuropsychiatric disorders (Maselli and Hanania, 2018; Postma and Rabe, 2015; Zhu et al., 2018, 2019, 2020, 2021). The complex etiology and clinical presentations of asthma in turn complicate standards for defining the phenotype; for example, nearly 60 different definitions of “childhood asthma” were used across more than 100 studies in the literature (Van Wonderen et al., 2010). This likely also contributes to the observed variation in asthma prevalence across populations.
High quality clinical models of asthma are necessary to comprehensively and quantitatively investigate risk factors within and among populations. However, with relatively few large and diverse datasets of asthma, genetic risk predictors for asthma have been rare in the literature. For a highly heterogeneous disease like asthma, polygenic risk scores (PRS) may be particularly useful tools for predicting subtypes, disease severity, and the development of comorbidities in the clinical setting, and for investigating potential gene-environment interactions in the research setting. Existing PRS for asthma were generated primarily from UK Biobank (UKBB), a large population-based cohort study, and BioBank Japan (BBJ), a similarly large-scale hospital-based study (Amariuta et al., 2020), or cohorts of smaller sample sizes (Song et al., 2020; Sordillo et al., 2021). The PRS had limited predictive ability in these studies and lower performance compared to PRS for other common complex diseases. This underscores the genetic complexity of asthma and highlights the need for more large-scale, genomic studies of asthma.
To more deeply interrogate the genetic architecture of asthma across different populations, we analyzed paired phenotypic and genetic data from the Global Biobank Meta-analysis (GBMI). Participating biobanks shared summary statistics for the meta-analyses of 14 phenotypes, including asthma (Zhou et al., 2021). Compared to previous asthma resources and studies, this collaborative effort increased both the sample size and diversity of asthma cases by many folds, covered more variants with higher imputation quality, and harmonized phenotypes using consistent electronic health record definitions for asthma across datasets. Harnessing this resource, we identified 49 loci not previously associated with asthma. Despite prevalence differences of nearly an order of magnitude, we also demonstrated remarkable consistency of genetic effects across the biobanks and ancestries in GBMI. Further, we showed that the increased diversity of data from GBMI improves genetic risk prediction accuracy in multiple populations. Finally, we provided additional evidence for shared genetic architectures between asthma and known coexisting diseases such as COPD and hay fever. Our findings highlight the need for future investigations into how genetic effects shared with different diseases contribute to the heterogeneity of asthma.
Results
Meta-analysis for asthma across 18 biobanks in GBMI
To identify novel loci associated with asthma, we performed fixed-effects inverse-variance weighted meta-analysis using the harmonized GWAS summary statistics for asthma from 18 biobanks participating in GBMI. The combined sample size from all studies was 153,763 cases and 1,647,022 controls, spanning individuals of European (EUR), African (AFR), Admixed American (AMR), East Asian (EAS), Middle Eastern (MID), and Central and South Asian (CSA) ancestry (Fig. 1). Despite the standardized phenotype definitions used by each biobank, which included the asthma PheCode and/or self-reported data (Supplementary Table 3), the prevalence of asthma varies widely across these biobanks, ranging from 3% in the Taiwan Biobank to 24% in the Mass General Brigham Biobank. We applied pre- and post-GWAS quality control filters that resulted in 70.8 million SNPs for meta-analysis; for downstream analyses we analyzed SNPs present in at least 2 biobanks (Zhou et al., 2021). The meta-analysis identified 180 loci of genome-wide significance (p < 5×10−8), 49 of which have not been previously reported to be associated with asthma (Fig. 2A). The potentially novel associations had smaller effect sizes on average compared to the previously reported loci, across the allele frequency spectrum (Fig. 2B). This illustrates that with the increased power and effective sample size of GBMI, we were able to uncover SNPs with more modest effects on asthma.
Because the GBMI meta-analysis includes data from UKBB, we compared our results to the TAGC meta-analysis results that did not include the UKBB GWAS to facilitate analyses that require independent samples (Demenais et al., 2018). Of the 180 lead variants in GBMI, 122 were in the TAGC meta-analysis or had a tagging variant in high LD (r2 > 0.8) in the TAGC study; 76 of these had p < 0.05 in the TAGC results. We compared the effect sizes of these 76 SNPs in the GBMI and the TAGC meta-analyses using a previously proposed Deming regression method that considers standard errors in both effect size estimates (Deming, 1943). We found that all 76 SNPs were directionally consistent and aligned across the studies (Supplementary Table 4, Supplementary Fig. 1).
Among the 49 novel loci, six were missense variants or in high LD (r2 > 0.8) with a missense variant (Supplementary Table 2). One of these SNPs, chr10:94279840:G:C (pmeta-analysis = 2.5×10−9), resides in PLCE1, an autosomal recessive nephrotic syndrome gene (Jefferson and Shankland, 2007); high prevalence of atopic disorders, like asthma, among children with nephrotic syndrome has long been observed in the clinic, suggesting potential shared pathways underlying asthma and nephrotic syndrome (Riar et al., 2019). The asthma risk allele has also been previously linked to lower blood pressure (UK Biobank — Neale lab). The lead SNPs chr14:100883117:G:T (pmeta-analysis = 2.6×10−8) and chr19:56222056:C:A (pmeta-analysis = 2.4×10−8) also implicate novel genes, RTL1 and ZSCAN5A respectively. RTL1 has been found to play a role in muscle regeneration (Loo et al., 2019), and ZSCAN5A has been linked to monocyte count (Chen et al., 2020). Additionally, three of the novel lead SNPs co-localized with a fine-mapped cis-eQTL (Supplementary Table 2). One of these, chr19:49513502:C:T (pmeta-analysis = 7.98×10−9), implicates FCGRT, which regulates IgG recycling and is a potential drug target for autoimmune neurological disease therapies (Dalakas and Spaeth, 2021). The other previously-reported missense variants replicated previous findings; among these, chr4:102267552:C:T (p.Ala391Thr, p = 2.5×10−12) is a highly pleiotropic variant in SLC39A8 that has been associated with many psychiatric, neurologic, inflammatory and metabolic diseases (Huang et al., 2017; International Consortium for Blood Pressure Genome-Wide Association Studies et al., 2011; Li et al., 2016; Nebert and Liu, 2019; Pickrell et al., 2016; Schizophrenia Working Group of the Psychiatric Genomics Consortium, 2014; Speliotes et al., 2010) and has been demonstrated to disrupt manganese homeostasis (Nakata et al., 2020). Variants implicating well-known asthma-associated genes with large effects, like IL1RL1, IL2RA, STAT6, IL33, GSDMB, and TSLP, were replicated in the meta-analysis as well.
Shared genetic architecture of asthma across biobanks and ancestry groups
Given that sample size, disease prevalence, ancestry, and sampling approaches differed across the 18 biobanks, we next investigated the consistency of the asthma-associated loci across the biobanks and their attributes. We first implemented an approach that estimates the correlation (r b) between the effects of the 180 lead variants in each biobank GWAS and the corresponding meta-analysis excluding that biobank (Qi et al., 2018). Most of the r b estimates were highly correlated (i.e. did not differ significantly from 1) (Supplementary Table 5). To further interrogate the consistency of the lead variants in all biobanks, we computed the ratio of the effect size of these SNPs in the biobank-specific GWAS over that in the corresponding leave-that-biobank-out meta-analysis. We found that the average per-biobank ratios were almost evenly split between those greater than and less than 1 (Supplementary Fig. 2). This indicates that any significant difference in effects likely does not reflect technical artifacts in the meta-analysis or GWAS procedures. We also applied Deming regression to assess the alignment of the lead SNP effects in each biobank-specific GWAS with the effects in the corresponding leave-that-biobank-out meta-analysis (Deming, 1943) and observed that the effect sizes were comparable across the biobanks (Fig. 3). Furthermore, the genome-wide genetic correlations between the biobanks with non-zero heritability estimates and the respective leave-that-biobank-out meta-analyses were all close to 1 (Zhou et al., 2021). Taken together, these analyses indicate that the genetic architecture of asthma is largely shared across diverse cohorts, even when cohort characteristics like sample size and disease prevalence differ.
We also found little evidence of heterogeneity in the ancestry-specific effect sizes for the lead SNPs. One SNP, chr10:9010779:G:A, was significantly heterogeneous (p-value for Cochran’s Q test < 0.0003, the Bonferroni-corrected p-value threshold) across the ancestry-specific meta-analyses of AFR, AMR, CSA, EAS, and EUR individuals (Fig. 4A, Supplementary Table 6). One known SNP that nearly reached the Bonferroni-corrected p-value threshold for heterogeneity, chr16:27344041:G:A, displayed different effects in the EUR and EAS cohorts. This SNP lies within an intron of IL4R (Fig. 4B), which has known associations with asthma (El-Husseini et al., 2020; Massoud et al., 2016). Previous studies have investigated the association of IL4R with asthma in different populations, with inconsistent results, so future studies on the population-specific effects of this gene will be needed (Battle et al., 2007; Kousha et al., 2020; Nie et al., 2013). Our findings demonstrate that despite broad consistency of effect sizes across ancestries among the lead variants, the increased power and diversity of GBMI enabled the detection of SNPs with significantly different effects across ancestries.
Polygenic risk score analyses
We next explored the impact of the increased sample sizes and diversity in GBMI on genome-wide risk prediction of asthma. To establish a baseline understanding of PRS performance for the 14 phenotypes analyzed in GBMI, we implemented PRS-CS (Ge et al., 2019) to construct PRS using the leave-one-biobank-out meta-analyses as discovery data (Wang et al., 2021). The PRS derived from the GBMI leave-one-biobank-out meta-analyses of asthma had higher predictive accuracy, as measured by R2 on the liability scale , compared to the PRS constructed from the TAGC meta-analysis (Demenais et al., 2018) across all biobanks tested (Fig. 5).
To expand on these analyses, we tested a recently-developed extension of PRS-CS, PRS-CSx (Ruan et al., 2021). This method jointly models multiple summary statistics from different ancestries to enable more accurate effect size estimation for prediction. For input to PRS-CSx, we used the AFR, AMR, EAS, CSA, and EUR ancestry-specific meta-analyses from GBMI; the discovery meta-analysis that matched the ancestry of the target cohort excluded the target cohort (Supplementary Fig. 5). With the posterior SNP effect size estimates from PRS-CSx, we tested the multi-ancestry PRS in the following target populations: AFR ancestry individuals in UKBB, CSA ancestry individuals in UKBB, a holdout set of EAS ancestry individuals in BBJ, and a holdout set of EUR ancestry individuals in UKBB. The final prediction models tested in these target populations were the optimal linear combinations of the population-specific PRS. In both the EAS and EUR target cohorts, the approached the SNP-based heritability , estimated to be 0.085 for asthma using the all-biobank meta-analysis (Wang et al., 2021) (Fig. 5). However, we acknowledge that estimates may differ across biobanks and ancestries given differences in disease prevalence, environmental exposures, phenotype definitions, and other factors; these differences may contribute to the PRS in EAS individuals performing similarly to PRS in EUR individuals in our analyses, despite the smaller sample size of the EAS discovery cohort. The across the target populations for the PRS-CSx scores were roughly the same as the of the PRS derived from the PRS-CS analyses. It is important to note that the discovery data used in the PRS-CS analyses differed slightly in sample size and composition, since the leave-one-biobank-out approach was used for PRS-CS, but the target cohorts in which the PRS were evaluated were the same (Supplementary Table 14). To investigate why improvement in performance using PRS-CSx was only incremental in most of the target cohorts, we examined the performances of each population-specific PRS. We found that across all target cohorts, PRS derived from either the EUR or EAS set of posterior effec size estimates outperformed the linear combination, and the of these PRS were also higher compared to that of the PRS-CS scores (Supplementary Fig. 6). This suggests that the addition of more discovery GWAS to PRS-CSx can improve the accuracy of PRS based on a single set of posterior effect size estimates, but the linear combination of PRS from multiple GWAS does not necessarily yield higher accuracy. This may be due to the considerably smaller sample sizes of some of the input discovery meta-analyses in our analyses and thus varying signal to noise ratios.
Genetic correlation between asthma and significantly heritable phenotypes across biobanks
Previous studies have shown that asthma is genetically correlated with a wide range of diseases and traits (Maselli and Hanania, 2018; Postma and Rabe, 2015; Zhu et al., 2018, 2019, 2020, 2021). We aimed to do a comprehensive genetic correlation analysis and identify all phenotypes that are significantly correlated with asthma from the GBMI biobanks. Among the 14 phenotypes analyzed in GBMI, COPD, a late-onset respiratory disease, had the highest genetic correlation with asthma (rg (se) = 0.67 (0.021), p = 1.55×10−226). This genetic correlation estimate is higher than estimates from previous studies, which ranged from 0.38-0.42 (Hobbs et al., 2017; Sakornsakolpat et al., 2019). This may be a result of greater power in the GBMI meta-analysis of COPD and different phenotype definitions used. Of note, only 6 of the 180 asthma-associated index variants (3%) had a genome-wide significant p-value in the COPD meta-analysis. Conversely, 12 of the 46 COPD-associated index variants (26%) had a genome-wide significant p-value in the asthma meta-analysis (Supplementary Table 7); 8 of these 12 variants were previously reported COPD associations. This shows that despite the strong genetic overlap between the COPD and asthma meta-analyses, the COPD meta-analysis largely identified variants with COPD-specific effects.
We next expanded these genetic correlation analyses beyond GBMI to measure correlations between asthma and the full breadth of phenotypes in UKBB. Of the 7,142 phenotypes for which GWAS were conducted in the EUR ancestry cohort in UKBB, 1,008 were significantly heritable (heritability Z score > 4) (Pan UKBB). Applying linkage-disequilibrium score correlation (LDSC) to these GWAS and the GBMI leave-UKBB-out meta-analysis of asthma, we obtained pairwise genetic correlation estimates between the heritable UKBB phenotypes and asthma, and observed strong correlations (|rg| > 0.4) with 95 of these phenotypes, which spanned prescriptions, PheCodes, and other categories (Supplementary Table 8). Many of these replicated previously-found correlations with asthma. Digestive system disorders, including gastritis and gastroesophageal reflux disease (GERD), emerged as a disease category with significant and strong genetic correlations with asthma. Although the association between asthma and digestive disorders has not been as well studied, this does reinforce a previous finding of shared genetics between asthma and diseases of the digestive system (Demenais et al., 2018), indicating that the commonly-observed co-presentation of asthma and gastroesophageal disease in the clinic may be partially due to pleiotropic genetic effects. Our results also showed moderate and significant correlations (rg ranging from 0.2-0.3) between asthma and neuropsychiatric diseases, like anxiety and depression, and obesity-related traits, like body mass index, which is similarly consistent with previous findings (Zhu et al., 2019, 2020).
To assess the consistency of the correlations in another population, we computed genetic correlation estimates between the GBMI EAS meta-analysis of asthma and other phenotypes in BBJ (Supplementary Table 9, Supplementary Fig. 4). Of the 48 available diseases in BBJ, 8 were significantly heritable in both BBJ and UKBB. Among these phenotypes, COPD showed the strongest and most significant correlation with asthma in BBJ (rg = 0.29, p = 6.41×10−6), although this was notably lower than across the full GBMI meta-analysis; differences in phenotype definition and curation may potentially contribute to variation in correlation estimates. Among all significantly heritable phenotypes in BBJ, pollinosis, also known as hay fever, showed moderate correlation with asthma as well (rg = 0.28, p = 0.0004). These were directionally consistent with the correlation results from UKBB, which showed strong and significant correlation with COPD (rg = 0.71, p = 3.88×10−57) and moderate correlation with pollinosis (rg = 0.39, p = 4.60×10−3).
Gene and pathway enrichment analyses for asthma
To further evaluate the extent of genetic overlap between asthma and COPD, we applied a gene prioritization method, MAGMA, to the GBMI EUR, AFR, EAS, and CSA meta-analyses of asthma as well as the GBMI EUR, AFR, and EAS meta-analyses of COPD (de Leeuw et al., 2015). After Bonferroni correction, we found that 442, 149, and 6 genes were significantly associated with asthma in the EUR (p < 2.50×10−6), EAS (p < 2.50×10−6), and CSA (p < 2.52×10−6) populations, respectively, with no significantly associated genes in the AFR cohort (all p > 2.51×10−6) (Supplementary Table 10). The majority of the genes associated with asthma identified in the EAS meta-analysis overlapped with the genes from the EUR meta-analysis (126 out of 149 genes), and all 6 genes associated with asthma as identified in the CSA meta-analysis were also significantly associated in the EUR and EAS meta-analyses. We identified 46 and 33 genes significantly associated with COPD in the EUR (p < 2.50×10−6) and EAS (p < 2.50×10−6) cohorts, respectively, and similarly to asthma, no significantly associated genes from the AFR meta-analysis (all p > 2.51×10−6) (Supplementary Table 11). Of the 75 genes associated with COPD across the EUR and EAS meta-analyses, 24 overlapped with the asthma-associated genes.
Utilizing these sets of genes, we also adopted MAGMA for gene-set enrichment based on the curated and ontology gene sets from the Molecular Signatures Database (MSigDB) (Subramanian et al., 2005). We found hundreds of gene sets that were significantly enriched (FDR < 0.05) by the asthma-associated genes discovered in the EUR and EAS meta-analyses (Supplementary Table 12). In contrast, only a handful of gene sets were significantly enriched among COPD-associated genes discovered in the AFR meta-analysis, likely reflecting the smaller overall sample size of COPD (Supplementary Table 13). The top-ranked asthma pathways from the EUR meta-analysis included cytokine and interleukin signaling and T-cell activation. Consistently biologically, the EAS meta-analysis identified autoimmune thyroid disease and graft vs. host disease pathways. The top-ranked COPD pathways from the EUR meta-analysis, although not significant, included several pathways related to nicotine receptor activity. These results reinforce that despite the substantial genetic overlap, asthma and COPD are governed by distinct biological processes as well. Future investigations will be required to fully parse out the etiology and comorbidities of asthma, like COPD, that develop later on in adulthood.
Discussion
Assembling the largest and most diverse collection of asthma cohorts to date, we conducted a GWAS meta-analysis of 18 biobanks around the world and identified 49 novel associations. Despite the substantial sample sizes of previous meta-analyses of asthma (Demenais et al., 2018), our results indicated that the heterogeneity and complexity of asthma necessitate even larger sample sizes for genomic discovery. We interrogated the overall consistency of genetic effects across the cohorts and found that in spite of variability in recruitment continent, sampling strategy, health system design, and disease prevalence, the effects of the loci discovered in the meta-analysis were mostly concordant across the biobanks. Additionally, genetic correlation estimates across ancestries, which ranged from 0.65 to 0.99 for the well-powered ancestry groups, as well as specific genes prioritized using MAGMA, strongly supported the finding that the genetic architecture of asthma is largely shared across the ancestry groups studied. This study provided further evidence for substantial genetic overlap between asthma and well-known, immune-related comorbidities like COPD and allergic diseases. We also identified genetic correlations between asthma and less well-studied comorbidities like digestive system disorders, while highlighting additional complexity in the etiology and comorbidities of asthma. For example, gene set enrichment analyses using MAGMA did not yield many shared pathways for asthma and COPD despite the strong genetic correlation.
We also showed that applying novel Bayesian PRS construction methods like PRS-CSx and PRS-CS to association data from larger and more diverse cohorts can improve prediction in asthma, particularly for populations of non-European ancestry. However, we found that imbalances in the sample sizes of the discovery cohorts may need to be taken into careful consideration when using these methods. Previous studies have shown that imbalanced sample sizes across ancestries contribute somewhat unpredictably to varying prediction performances, with a low signal-to-noise ratio in ancestry-matched target populations reducing prediction performance (Majara et al., 2021). Therefore, further investigations are needed to fully understand the interplay between sample size and ancestry in the context of polygenic prediction. Ultimately, these analyses highlight the pressing need for more well-powered and ancestrally-diverse resources that will help reduce these imbalances.
We have highlighted the harmonization of the phenotype definitions across biobanks, but it is important to acknowledge that the criteria used, which allowed for both self-reported and PheCode information, are vulnerable to imprecision and variability in the data collected. Self-reported data for asthma is particularly susceptible to imprecision, since it relies on personal recollection of asthma diagnoses that are often given in childhood. On the other hand, PheCodes, which are based on ICD codes, may fail to capture diagnoses made earlier in the lifetime of individuals in hospital-based cohorts. Therefore, including both self-reported and PheCode data -- an approach adopted by some but not all biobanks -- may be optimal for association analyses for asthma. In the UKBB we found that the genetic correlation between the EUR GWAS based on the asthma PheCode and the EUR GWAS based on either PheCode or self-reported data was nearly 1 (rg = 0.97). So, while we have demonstrated that the variation of phenotype definition used does not significantly influence the genetic association results in this case, we cannot confirm the same pattern for all biobanks in GBMI and especially for other diseases. However, given the relative alignment of genetic effects across the biobanks, we would expect that minor differences in phenotype definition would not substantially change the association results.
Additionally, we acknowledge that since the definitions used here for asthma and COPD do not exclude individuals with concurrent diagnoses, we are not able to fully distinguish the distinct biological pathways affecting asthma and COPD. Comorbidity rates of asthma and COPD reported in the literature range across studies but population-based estimates generally are low, around 2-3% (Hosseini et al., 2019; Kumbhare et al., 2016), while hospital-based prevalence estimates tend to be higher, around 13% (Akmatov et al., 2020). Among biobanks participating in GBMI, for example, 15.5% of all individuals with asthma in UKBB have a concurrent COPD diagnosis, 21% in BioVU, and 7.4% in BBJ. A previous study found that using stricter definitions of asthma, such as excluding subjects with COPD, resulted in stronger association signals for some of the asthma-associated loci (Han et al., 2020). However, it is important to note that this case exclusion would introduce an additional source of ascertainment bias. We also note that estimates of genetic correlation by LDSC are not biased by sample overlap (Bulik-Sullivan et al., 2015). In fact, this has been explored in the context of asthma and allergic diseases, where rg estimates from LDSC were shown to be robust to overlapping cases and/or controls (Zhu et al., 2018).
We also recognize the importance of analyzing environmental factors in conjunction with genetic factors for a disease that is heavily influenced by the environment. Our genetic analyses offer insight into the potential shared biological pathways that may be differentially affected by non-genetic factors, but we were not able to explicitly investigate environmental effects given the lack of available environmental exposure data among the biobanks. The high degree of alignment among genetic associations, coupled with the large variability in asthma prevalence, points to a particularly important role of the environment for asthma risk across populations. Gaining a greater understanding of the specific non-genetic factors that contribute to asthma development in different environments may help guide more accurate disease prediction across populations.
This study, and importantly the data sharing across biobanks facilitated by this initiative, have laid the groundwork for deeper dives into the shared and distinct genetic signatures of asthma subtypes. Previous studies have categorized the UKBB individuals with asthma into childhood-versus adult-onset subtypes based on their ages at first diagnosis, discovering partly distinct genetic architectures (Ferreira et al., 2019; Zhu et al., 2020). Data from other biobanks in GBMI make it possible to perform similar stratifications and enable multiple downstream analyses. For example, future studies can evaluate the genetic overlap between subtypes, further validate previously reported subtype-specific variants in different populations, and test the power of PRS to discern different subtypes, empowered by the meta-analysis conducted here. Additionally, with access to biobanks with a wide range of phenotypes beyond the 14 analyzed in this initiative, we can begin to tease out the underlying biological relationships between asthma subtypes and other correlated phenotypes, particularly immune-related pulmonary diseases. Harnessing these cross-trait correlations for prediction may also be a fruitful approach to improving the accuracy of polygenic prediction models for asthma.
Methods
Asthma phenotype definitions for association analyses
The phenotype definition guidelines that were developed by GBMI and shared with all participating biobanks can be found in Zhou et al. (Zhou et al., 2021). Disease endpoints, including asthma, were defined following the PheCode maps, which maps ICD-9 or ICD-10 codes to PheCodes (Denny et al., 2013). Asthma cases were all study participants with the asthma PheCode, and controls were all study participants without the asthma PheCode (or asthma-related PheCodes). Biobanks that did not have ICD codes primarily used self-reported data (Supplementary Table 3).
Meta-analysis for asthma in GBMI
We performed fixed-effects meta-analysis with inverse variance weighting for 18 biobanks in GBMI: China Kadoorie Biobank (CKB), Generation Scotland (GS), Lifelines, QSKIN, East London Genes & Health (GNH), HUNT, UCLA Precision Health Biobank (UCLA), Colorado Center for Personalized Medicine (CCPM), Mass General Brigham (MGB), BioVU, BioMe, Michigan Genomics Initiative (MGI), BioBank Japan (BBJ), Estonian Biobank (ESTBB), deCODE Genetics (DECODE), FinnGen, Taiwan Biobank (TWB), and UK Biobank (UKBB). Basic information on the biobanks are described in Zhou et al., as well as details on the genotyping, imputation, GWAS, post-GWAS quality control, and meta-analysis procedures (Zhou et al., 2021). Genetic variants with minor allele count (MAC) < 20 and imputation score < 0.3 were excluded from the analyses. Altogether, these cohorts had a total sample size of 153,763 cases and 1,647,022 controls (Supplementary Table 1). GWAS meta-analyses were first conducted within continental ancestry groups to control for population stratification. 5,051 cases and 27,607 controls were of African (AFR) ancestry; 4,069 cases and 14,104 controls were of Admixed American (AMR) ancestry; 18,549 cases and 322,655 controls were of East Asian (EAS) ancestry; 121,940 cases and 1,254,131 were of European (EUR) ancestry; 139 cases and 1,434 controls were of Middle Eastern (MID) ancestry; and 4,015 cases and 27,091 controls were of Central and South Asian (CSA) ancestry.
Lead SNP and locus definitions
We used a threshold of p < 5e-8 to identify SNPs with a genome-wide significant effect. To identify lead variants, we used a window size of 500 kb upstream and downstream of the SNPs with the strongest evidence of association in the meta-analysis, and merged overlapping regions until no genome-wide significant variants were detected within the ± 500 kb region. To designate lead SNPs as previously discovered or potentially novel, we compiled a list of known asthma-associated SNPs (p < 5×10−8) from the associations collected by El-Husseini et al. (El-Husseini et al., 2020) and listed in the GWAS catalog (as of 11/14/2021) (Buniello et al., 2019). We extended 500 kb upstream and downstream of each of these variants to define a locus, and intersected these regions with the loci defined by the lead SNPs in our meta-analysis to identify any overlaps. We annotated genetic variants with the nearest genes using ANNOVAR (Wang et al., 2010) and putative loss-of-function using VEP (McLaren et al., 2016) with the LOFTEE plug (Karczewski et al., 2021) as implemented in Hail (Zhou et al., 2021). We also annotated whether the lead or tagging variants (r2 > 0.8) of asthma were fine-mapped in any of the cis-eQTL fine-mapping resources. We retrieved cis-eQTL fine-mapped variants with posterior inclusion probability (PIP) > 0.9 in any tissues and cell types from the GTEx v8 (Aguet et al., 2020) and eQTL catalogue release 4 (Kerimov et al., 2021). Fine-mapping was conducted using SuSiE (Wang et al., 2020) with summary statistics and covariate-adjusted in-sample LD matrix as described previously (Kanai et al., 2021; Ulirsch and Kanai) (Supplementary Table 2).
Lead SNP effects across biobanks
To estimate the correlation of SNP effects for the 180 top-associated SNPs between one specific biobank and the leave-that-biobank-out meta-analysis, we used the method proposed by Qi et al. using GWAS summary statistics (Qi et al., 2018) (Supplementary Table 5).
Specifically, the method directly calculates SNP effect correlation as: where and denote the estimated SNP effects from GWAS conducted in one specific biobank and from GWAS performed in the leave-that-biobank-out meta-analysis, respectively. The is calculated as the sampling covariance between and . The and are the estimated variances of and , separately. The and are the estimated variance of the estimation errors of and , which are approximated as the mean of the squared standard errors of estimated SNP effect ( and ) across all the top-associated SNPs, respectively. The standard error of is obtained through the jackknife approach by leaving one SNP out each time. SNPs with large standard errors in CKB and HUNT (chr12:123241280:T:C and chr17:7878812:T:C, respectively) were excluded from these analyses.
Then, for the lead SNPs present in each biobank, we computed: for the biobank and leave-that-biobank-out pair. We took the average ratio across the lead SNPs for each biobank and leave-that-biobank-out pair. We then used the regression method introduced in Deming et al., which considers the errors in both the X- and Y-variables, to compare the effect sizes of these SNPs in each biobank GWAS with their effects in the leave-that-biobank-out meta-analysis (Deming, 1943). We set the intercept equal to 0 for these analyses.
Ancestry-specific heterogeneity analyses
To assess heterogeneity of per-SNP effect sizes for the 180 lead SNPs across ancestries in GBMI, we conducted ancestry-specific meta-analyses of the five most well-powered ancestry groups in GBMI (EUR, AFR, AMR, EAS, and CSA). We applied the Cochran’s Q test to the SNP effects in the ancestry-specific meta-analyses and identified SNPs with significant heterogeneity based on a Bonferroni-corrected p-value cut-off of 0.05/169 = 0.0003, accounting for the number of SNPs present in all studies (Supplementary Table 6). Regions displaying heterogeneity in effects across ancestry groups were visualized using the LocalZoom tool (Boughton et al., 2021).
Polygenic risk score analyses
A description of the PRS analyses conducted using PRS-CS (Ge et al., 2019), as well as the leave-one-biobank-out meta-analysis strategy applied, is provided in Wang et al. (Wang et al., 2021).
We used PRS-CSx, which jointly models GWAS summary statistics from populations of different ancestries and returns posterior SNP effect size estimates for each input population (Ruan et al., 2021). We applied this method to the AMR, AFR, CSA, EAS, and EUR ancestry-specific meta-analyses, which served as the discovery data for PRS construction. For the ancestry-specific meta-analysis that matched the ancestry of the target cohort, we excluded the target cohort. We evaluated the predictive performance of the PRS in 4 target cohorts: 1) AFR ancestry individuals in UKBB (849 cases, 5190 controls), 2) CSA ancestry individuals in UKBB (1232 cases, 6744 controls), 3) EAS ancestry individuals in BBJ that were part of a randomly-selected 1k holdout set (500 cases, 500 controls), and 4) EUR ancestry individuals in UKBB that were part of a randomly-selected 10k holdout set (1164 cases, 7577 controls). As an example, for the AFR ancestry individuals, the full set of discovery data for PRS construction consisted of the AMR, CSA, EAS, and EUR ancestry-specific meta-analyses, as well as the AFR ancestry-specific meta-analysis excluding the AFR ancestry individuals in UKBB. The same strategy was applied to the other 3 target cohorts (Supplementary Fig. 5). We used ancestry-matched LD reference panels from UKBB data and the default PRS-CSx settings for all input parameters. We evenly and randomly split cases and controls in the target cohorts into validation and testing subsets. Using the posterior SNP effect size estimates from PRS-CSx, we computed one PRS per discovery population for the validation subsets to learn the optimal linear combination of the ancestry-specific PRS using PLINK v1.9 (Chang et al., 2015; Purcell and Chang). Then, with these weights, we evaluated the prediction accuracy of this linear combination of PRS in the testing subset. We reported the average prediction accuracy, measured by variance explained on the liability scale , estimated using the prevalence of asthma in the BBJ biobank for the EAS target cohort and in the UKBB biobank for the other target cohorts, across 100 random splits.
Genetic correlation analyses in UKBB and BBJ
We used linkage-disequilibrium score correlation (LDSC) to compute pairwise genetic correlations (rg) (Bulik-Sullivan et al., 2015). We estimated rg between all EUR-ancestry UKBB phenotypes with heritability Z-score > 4 and (1) the GBMI leave-UKBB-out meta-analysis for asthma and (2) the UKBB EUR-ancestry GWAS of asthma (PheCode ID 495 in UKBB). The heritability Z-scores were obtained from the stratified-LDSC (S-LDSC) computations of heritability reported by the Pan-UK Biobank team (Finucane et al., 2015, 2018; Pan UKBB). Summary statistics from the UKBB EUR GWAS were obtained from the Pan-UK Biobank team as well (Pan UKBB).
We also used LDSC (Bulik-Sullivan et al., 2015) to compute rg between 48 phenotypes in BioBank Japan (BBJ) and (1) the GBMI leave-BBJ-out meta-analysis for asthma and (2) the BBJ GWAS of asthma. We used publicly available GWAS summary statistics for all traits (Akiyama et al., 2017; Ishigaki et al., 2020; Kanai et al., 2018). Genetic correlation results were visualized using the R corrplot package (Wei and Simko, 2021).
Gene- and pathway-based enrichment analyses for asthma and COPD
Fixed-effects meta-analysis with inverse variance weighting was also performed for 16 biobanks in GBMI with COPD data: BBJ, BioMe, BioVU, CCPM, CKB, ESTBB, FinnGen, GNH, GS, HUNT, Lifelines, MGB, MGI, TWB, UCLA, and UKBB. The same processing and methods were applied here as for the asthma meta-analysis. These cohorts had a total sample size of 81,568 cases and 1,310,798 controls. COPD cases were defined based on the COPD PheCode, and controls were all study participants without the COPD or COPD-related PheCodes. Biobanks that did not have ICD codes available used spirometry data (Lifelines) or self-reported data (TWB). Details can be found in Zhou et al. (Zhou et al., 2021). Meta-analyses were also conducted within continental ancestry groups: 19,044 cases and 310,689 controls of EAS ancestry, 1,978 cases and 27,704 controls of AFR ancestry, and 58,559 cases and 937,358 controls of EUR ancestry.
We used MAGMA v1.09b for gene prioritization and gene-set enrichment analyses, applying this method to the GBMI asthma EUR, AFR, EAS, and CSA ancestry-specific meta-analyses and the GBMI COPD EUR, AFR, and EAS ancestry-specific meta-analyses (de Leeuw et al., 2015). For the gene-level analyses in MAGMA, we first mapped the SNPs to the provided list of genes using a window size of 20kb, and then performed gene analysis using the ancestry-matched 1000G LD reference panels to account for LD structure. Gene-set enrichment was performed using the default settings to correct for gene length, gene density, and the inverse mean minor allele count. The gene sets used were the c2, “curated gene sets,” and c5, “ontology gene sets,” obtained from the Molecular Signatures Database v7.4 (Subramanian et al., 2005).
Data Availability
The all-biobank GWAS summary statistics are publicly available for downloading at https://www.globalbiobankmeta.org/resources and can be browsed at the PheWeb Browser (http://results.globalbiobankmeta.org).
Data and Code Availability
The all-biobank GWAS summary statistics are publicly available for downloading at https://www.globalbiobankmeta.org/resources and can be browsed at the PheWeb Browser (http://results.globalbiobankmeta.org). Custom scripts used for quality control, meta-analysis, and loci definition are available at https://github.com/globalbiobankmeta. Other analyses utilized publicly available tools: the R deming package for Deming regression (Therneau, 2018), PRS-CSx for polygenic prediction (https://github.com/getian107/PRScsx), LDSC for genetic correlation (https://github.com/bulik/ldsc), and MAGMA v1.09b for gene-set enrichment (https://ctg.cncr.nl/software/magma).
Author Contributions
Study design: K.T., A.R.M., M.J.D., B.M.N.
Data collection/contribution: W.Z., Y.O., T.M.
Data analysis: K.T., W.Z., Y.W., M.K., S.N., R.G., L.M., L.N.
Writing: K.T., A.R.M., S.N., R.G.
Revision: K.T., A.R.M., Y.W., R.G., M.J.D.
Declaration of Interests
M.J.D. is a founder of Maze Therapeutics. B.M.N. is a member of the scientific advisory board at Deep Genomics and consultant for Camp4 Therapeutics, Takeda Pharmaceutical, and Biogen.
Supplementary Information
Supplementary Figures
Supplementary Tables
STable 1: Description of 18 biobanks in GBMI that contributed summary statistics for asthma meta-analysis with sample size, ancestry, and recruitment strategy information.
STable 2: 180 lead variants discovered by GBMI asthma meta-analysis with annotations for nearby genes, missense variants, and fine-mapped cis-eQTLs.
STable 3: Description of asthma definition used by each biobank.
STable 4: 122 of the 180 lead variants present or with tagging variant in the TAGC study with effect sizes and p-values from GBMI and TAGC meta-analyses.
STable 5: Correlations of SNP effects of the 180 lead variants between each biobank and the corresponding leave-that-biobank-out meta-analysis.
STable 6: Heterogeneity p-values, computed from Cochran’s Q statistic, for lead variants using effect sizes from the AFR, AMR, EAS, EUR, and CSA meta-analyses.
STable 7: 46 lead variants discovered by GBMI COPD meta-analysis with corresponding effect sizes, standard errors, and p-values from asthma meta-analysis.
STable 8: Genetic correlations estimated by LDSC between GBMI leave-UKBB-out meta-analysis and UKBB EUR GWAS of heritable phenotypes with phenotype descriptions.
STable 9: Genetic correlations estimated by LDSC between GBMI asthma meta-analyses and UKBB and BBJ GWAS of several diseases.
STable 10: Results from MAGMA gene analysis using GBMI EAS, CSA, and EUR asthma meta-analyses. Genes with p-values < Bonferroni-corrected p-value thresholds are reported.
STable 11: Results from MAGMA gene analysis using GBMI EAS and EUR COPD meta-analysis. Genes with p-value < Bonferroni-corrected p-value thresholds are reported.
STable 12: Results from MAGMA gene-set enrichment analyses for asthma. Gene sets with FDR < 0.05 are reported.
STable 13: Results from MAGMA gene-set enrichment analyses for COPD. Gene sets with FDR < 0.05 are reported.
STable 14: Description of discovery data used for PRS-CSx and PRS-CS analyses.
Acknowledgements
We acknowledge helpful comments from Cristen Willer, Hailiang Huang, Yunfeng Ruan, and Tian Ge. A.R.M is funded by the K99/R00MH117229. S.N. is supported by Takeda Science Foundation. Y.O. is supported by JSPS KAKENHI (19H01021, 20K21834), and AMED (JP21km0405211, JP21ek0109413, JP21ek0410075, JP21gm4010006, and JP21km0405217), JST Moonshot R&D (JPMJMS2021, JPMJMS2024), Takeda Science Foundation, and Bioinformatics Initiative of Osaka University Graduate School of Medicine, Osaka University. R.G. is supported by the T32AG000222.
Footnotes
Lead Contact: Kristin Tsuo (ktsuo{at}broadinstitute.org)