Multi-trait and multi-ancestry genetic analysis of comorbid lung diseases and traits improves genetic discovery and polygenic risk prediction

Largest genome-wide association study of lung function in the East Asian population identifies novel loci

Spirometry is a pulmonary function test that evaluates lung function. As continuous measures of lung function, they are useful traits for studying the heritable basis of respiratory diseases. However, the largest genetic studies of spirometry to date have been Eurocentric. To identify novel loci associated with lung function in understudied populations, we performed the largest GWAS to date of East Asian (EAS) ancestry individuals of FEV1, FVC, and FEV1/FVC, combining GWAS summary statistics from the Korean Cancer Prevention Study-II (KCPS2) and Taiwanese Biobank (TWB) (Supplementary Table 1, Methods)^23–25. We conducted fixed-effects inverse-variance weighted meta-analysis for the three continuous spirometry measures of lung function across 132,200 total individuals. We applied quality control filters (Methods) that resulted in 8 million unique single-nucleotide polymorphisms (SNPs) for meta-analysis.

The meta-analysis identified 44, 73, and 31 independent loci for FEV1, FVC, and FEV1/FVC, respectively (Figure 2, Supplementary Figure 1-3). Of these, 44 are novel and have not been previously associated with each trait (Supplementary Table 2), including 25 novel loci for FVC, 17 for FEV1, and 2 for FEV1/FVC. Among these novel loci, several were previously associated with height, including rs3782886 in BRAP (P=4.20×10^-8, β=0.0219 with FEV1), rs7290267 near FLJ27365 (P=1.49×10^-9, β=-0.0272 with FEV1), rs724016 in ZBTB38 (P=1.35×10^-10, β=-0.0178 with FVC), rs58744877 near GNAS (P=2.71×10^-14, β=-0.033 with FVC), rs3176466 near CDKN2C (P=4.59×10^-8, β=0.0236 with FVC), rs28839214 near RP11-361D14.2 (P=1.91×10^-13, β=0.0197 with FVC), rs149580940 in FIBIN (P=7.89×10^-10, β=0.0449 with FVC)^26–28. rs3782886 has also shown previous associations with other traits including alcohol use disorder, high-density lipoprotein, and type 2 diabetes^29–31. rs724016 has also been shown to be associated with Crohn’s disease³², and rs150971595 near GABBR1 (P=1.49×10^-8, β=-0.0307 with FVC) associated with total cholesterol level and triglyceride levels³³. Previous reported associations for each locus can be found in Supplementary Table 3.

• Download figure
• Open in new tab

Figure 2. Meta-analysis results of spirometry in East Asian cohorts.

Frequency and effect size of risk alleles of 44, 73, and 31 index variants for FVC (a), FEV1 (b), and FEV1/FVC (c). Novel associations are highlighted in red. Genes with large effect sizes are highlighted. d) Results of our combined EAS meta-analysis of FEV1/FVC with the largest multi-ancestry GWAS to date. The Manhattan plot of the previous multi-ancestry GWAS is shown in red on the bottom, with yellow dots representing previously identified loci associated with FEV1/FVC. The Manhattan plot of the meta-analyzed results of the previous study with our new EAS meta-analysis is shown in blue on the top, with new loci annotated and indicated with black stars. Gene names indicate the nearest gene to novel loci.

We combined the results of our EAS meta-analysis of lung function with the largest multi-ancestry GWAS of lung function to date³⁴. Incorporating results from KCPS2 and TWB resulted in 519, 526, and 476 total loci discovered for FEV1, FVC, and FEV1/FVC, respectively, of which 131, 163, and 149 were not originally reported to be associated with each trait (Figure 2, Supplementary Figure 4-5). Significantly associated loci are reported in Supplementary Table 4-6.

Pervasive pleiotropic effects across respiratory traits and diseases inform shared components of heritable risk

Measures of lung function are known clinical risk factors for respiratory diseases such as chronic obstructive pulmonary disease (COPD), asthma, and lung cancer. Using the largest multi-ancestry GWAS summary statistics to date, we found that lung function (FEV1, FVC, FEV1/FVC), smoking behavior (smoking status and cigarettes/day), and respiratory diseases (asthma, COPD, and lung cancer) were all genetically significantly correlated with each other (P<0.05) with the exception of FVC and lung cancer (Figure 3a).

• Download figure
• Open in new tab

Figure 3: Shared and distinct heritable components inform the molecular basis of trait differences.

a) Genetic correlations between traits (top) and heritability estimates (bottom) analyzed in this study based on S-LDSC correlations. Colors indicate magnitude of correlation and stars indicate degree of significance. Comparison of effect sizes of variants from GWAS on EUR samples for (b) COPD vs. Cigarettes smoked per day (Cig/day) and (c) Lung cancer vs. Cig/day. Effect sizes of variants on diseases (COPD and lung cancer) are on the x-axis and effect sizes of variants on Cig/day are on the y-axis. Three distinct groups of variant effects were identified by a Bayesian classifier in shared variants analysis, indicated by (colors) as follows: variants significantly associated with both traits (red), variants significantly associated with disease only (orange), and variants only significantly associated with Cig/day (navy). Solid dots are variants confidently assigned to the group indicated by the color of the dots with posterior probability >95%, while crosses are variants not confidently assigned to any of the groups. The colored ellipse range indicates the 95% probability regions of the fitted bivariate effect size distributions with each class. Nearest gene names are labeled on disease-specific variants with the highest posterior probability for each gene.

To examine the distinct and shared variants among respiratory diseases and auxiliary traits, we first evaluated the consistency of lead SNPs between the traits. Many variants demonstrated pleiotropic effects associated with multiple traits. In particular, we explored the shared and distinct genetic etiology among the three respiratory disease traits (asthma, COPD, and lung cancer), risk factors (smoking status, smoking heaviness), and a measure of lung function (FEV1/FVC). For each trait pair, we selected variants significantly associated with either or both of the traits. We computed the linear slopes of variant effect sizes within each category using an expectation-maximization (EM) algorithm, which typically revealed more than one linear trend and divergent slopes (Methods, Figure 3, Supplementary Figure 6-8). For example, genetic variants associated with the amount and frequency of smoking have mostly distinct effects from asthma. This suggests that there is no uniform explanation of the genetic mechanisms of the two traits across those variants.

We then implemented a Bayesian algorithm, linemodels³⁵, to classify the variant effect sizes into confident associations with trait 1, trait 2, or both traits. We used ancestry-specific (AFR, EAS, EUR) as well as meta-analyzed multi-ancestry GWAS results to identify variants having confident associations with a posterior probability above 0.95. We found 31 variants that were confidently associated with COPD and 108 with lung cancer, which are distinct from those associated with smoking heaviness, as measured by Cigarettes/Day, in EUR. No variants were confidently associated primarily with both traits or with smoking heaviness only (Supplementary Table 7-8). Two variants on chromosome 8 near the GULOP pseudogene, rs113623975 and rs56223946, were predominantly associated with both COPD and lung cancer.

When we further investigated variants with predominant effects on lung cancer across EAS, EUR, and the multi-ancestry meta-analysis GWAS results, we identified 7 variants from TERT appearing in all three groups and 6 variants significant only in the meta-analysis results (Supplementary Figure 9-10). These include rs34517439 from DNAJB4, rs11778371 from CHRNA2, and 4 other variants within or near TERT.

PRSxtra jointly models multi-ancestry and multi-trait effects to predict diseases and exacerbations

To develop PRSxtra, we combined the results of our new GWAS meta-analysis of lung function in EAS with the largest and most diverse GWAS of COPD, asthma, lung cancer, and smoking (Figure 1) to conduct ancestry-specific multi-trait GWAS (MTAG)¹³. We then used PRS-CSx to jointly model trait-specific MTAG results across ancestry groups to derive candidate scores^36,37 to include in the final PRSxtra.

First, MTAG resulted in a greater number of significant associations and loci discovered across traits (Supplementary Table 9). AMR and EUR ancestry-specific MTAG analysis resulted in the largest gains. For example, FVC in AMR and COPD in EUR cohorts had the greatest proportional increases in the number of lead SNPs discovered, from 2 to 34 (17-fold) and 27 to 441 (16-fold), respectively. In total, MTAG identified 609 more lead SNPs across all traits and ancestry groups. We assessed the gain in power for each run of MTAG by the increase in mean χ² statistic, where all but five GWAS (all of which were ancestry-specific spirometry GWAS) increased in power (Supplementary Table 9).

Next, we aggregated the summary statistics from MTAG to derive multi-trait and multi-ancestry polygenic risk scores. For each trait, we used PRS-CSx to jointly model trait-specific MTAG summary statistics across ancestry groups to derive a total of 39 candidate scores for individuals in the All of Us Research Program (AoU), a longitudinal cohort study continuously enrolling adults in the US with genotype, self-reported, and linked health record data. The program places a strong emphasis on including diverse populations that have traditionally been underrepresented in biomedical research. The individuals from AoU are independent from the cohorts in which the GWAS was obtained. For each trait, we randomly split participants who passed quality control into 70% for training and 30% for validation (Figure 1, Methods, Supplementary Table 10-13). We then used ridge regression to calculate linear combinations of the scores in the training set (Supplementary Table 14-16, Supplementary Figure 11-13).

We evaluated the performance of PRSxtra in the held-out validation cohort for predicting asthma (N_cases=9450, N_controls=64070), COPD (N_cases=4560, N_controls=66605), and lung cancer (N_cases=578, N_controls=70600) compared to a single ancestry- and trait-matched PRS (Figure 1, Supplementary Table 10). The validation cohort includes individuals of AFR, AMR, EAS, EUR, MID, and SAS genetic ancestries. In this multi-ancestry validation cohort, ancestry- and trait-matched PRS and PRSxtra were significantly correlated with each other (P<0.0001), with r=0.444, r=0.185, and r=0.120 for asthma, COPD, and lung cancer, respectively.

PRSxtra alone predicted asthma, COPD, and lung cancer more accurately than the trait- and ancestry-matched PRS alone (P<0.0001) in the multi-ancestry validation cohort (Supplementary Table 17-19). For asthma, the AUC improved from 0.543 (95% CI=[0.537, 0.549]) to 0.563 (95% CI=[0.557, 0.569]) using PRSxtra versus PRS. For COPD, AUC improved from 0.540 (95% CI=[0.532, 0.549]) to 0.589 (95% CI=[0.581, 0.597]). For lung cancer, AUC improved from 0.539 (95% CI=[0.516, 0.561]) to 0.592 (95% CI=[0.569, 0.616]) (Figure 3).

Among individuals with each disease, the median percentile of PRSxtra was much higher than PRS (Supplementary Table 20). We observed the largest difference between scores in lung cancer, where cases had a median PRSxtra in the 64th percentile compared to the 57th percentile for PRS. Controls had a median PRSxtra and PRS both in the 50th percentile. The AMR population had the largest improvement in disease prediction between PRS and PRSxtra. For asthma, AUC improved from 0.509 (95% CI=[0.493, 0.524]) to 0.630 (95% CI=[0.616, 0.645]). For COPD, AUC improved from 0.540 (95% CI=[0.532, 0.549]) to 0.589 (95% CI=[0.581, 0.597]). These differences could reflect the increased EAS GWAS power for spirometry and the relatively low genetic divergence between EAS and AMR groups. However, these new ancestry-specific GWAS do not explain the most phenotypic variation when predicting these traits. Specifically, we further assessed the robustness of PRSxtra by conducting leave-one-out analyses for predicting COPD and asthma in the AMR cohort, where there were the largest improvements in prediction. Removing any single score did not significantly change the prediction ability of PRSxtra. For example, removing the EUR lung cancer candidate PRS resulted in the largest decrease in AUC of -0.0175 for predicting COPD, and removing EUR FVC candidate PRS resulted in the largest decrease in AUC of -0.000421 for predicting asthma.

PRSxtra improved overall stratification in the multi-ancestry validation cohort compared to PRS. In the bottom and top deciles of PRSxtra, asthma prevalence was 8.4% and 16.8% (+8.4%), respectively, versus 9.0% and 14.1% (+5.1%) for PRS. For COPD, the prevalence in the bottom and top deciles of PRSxtra was 2.8% and 9.1% (+6.3%), respectively, compared to 4.6% and 7.6% for PRS (+3.0%). Lung cancer prevalence in the bottom and top deciles of PRSxtra was 0.5% and 1.4% (+0.9%), respectively, compared to 0.46% to 0.90% (+0.44%) for PRS (Supplementary Figure 14, Supplementary Table 21-23). PRSxtra was also better at identifying individuals with the highest risk of disease. For example, individuals in the top decile of PRSxtra had 3.5 times the odds of COPD compared to those in the first decile, while PRS had 1.7 times the odds. The improvement was most apparent in the AMR population for asthma and COPD. Individuals in the top decile of PRSxtra had 4-fold odds of asthma and COPD compared to the first decile, whereas the top decile of PRS had almost the same odds of asthma and COPD compared to the first decile (Supplementary Figure 14, Supplementary Table 24-26).

We compared the performance of PRS and PRSxtra alone to a joint model of the strongest clinical risk factors for each disease, i.e. family history of asthma, and smoking for COPD and lung cancer. Clinical risk factors generally provided the largest improvement in performance from sex and age. In the multivariable model, PRSxtra remained significantly better than PRS for predicting asthma in the multi-ancestry validation cohort and for predicting asthma and COPD in the AMR subgroup. For lung cancer, the risk scores do not add additional information to the multivariable model of sex, age, and smoking status (Figure 4, Supplementary Table 17-19).

• Download figure
• Open in new tab

Figure 4. PRSxtra significantly improves prediction for several respiratory diseases compared to PRS and clinical risk factors.

Diseases include: asthma (a), COPD (b), and lung cancer (c) in the full multi-ancestry held-out validation set and in ancestry-specific subgroups. Green dots indicate the performance of PRS or PRSxtra alone to predict diseases, blue dots represent the performance of the combined model with incremental variables added from sex and age, family history (for asthma) or smoking status (for COPD and lung cancer), and PRS or PRSxtra in gradually darker shades. The full validation set included individuals of AFR, AMR, EAS, EUR, MID, and SAS predicted genetic ancestry background.

To investigate performance in clinical subgroups, we tested the association of PRSxtra with COPD and lung cancer in subgroups of never and ever smokers, and of asthma in individuals with and without family history. We found similar performance between subgroups (Supplementary Table 17-19, Supplementary Figure 15-17). PRSxtra was positively associated with COPD and asthma exacerbations. PRSxtra predicted COPD exacerbation significantly better than PRS in the validation cohort, with AUC increasing from 0.572 (95% CI=[0.558-0.586]) to 0.600 (0.587-0.614), However, PRSxtra did not significantly improve prediction of asthma exacerbation compared to PRS (Supplementary Table 27).

Meta-analysis of lung function across East Asian populations

We performed a fixed-effects meta-analysis with inverse variance weighting, as implemented in METAL v2001-03-25 software, for FEV1, FVC, and FEV1/FVC across two East Asian cohorts using published summary statistics from Korean Cancer Prevention Study-II (KCPS2)^24,25, Taiwanese Biobank (TWB)^23,33. The details of their study design have been described previously. In summary, the Korean Cancer Prevention Study-II Biobank (KCPS2) is a prospective cohort study based in Korea of 153,950 subjects with genotype data and phenotype measurements between 2004 and 2013. The Taiwanese Biobank is a prospective cohort study of the Taiwanese population with 149,894 participants between the ages of 30-70 years old at recruitment (as of April 2021). FEV1 was defined as the total air blown between zero and 1 second (measured in Liters). FVC was defined as the vital capacity during forced expiration (measured in Liters). In KCPS2, FEV1 and FVC were measured using pulmonary/metabolic systems Vmax 20, Carefusion, USA. For each trait, samples with measurements that were more than 6 standard deviations away from the sample average were excluded.

Altogether, these cohorts had a total sample size of 132,200. In conducting our metanalysis, we excluded genetic variants with minor allele frequency < 0.01. We used FUMA v. 1.5.2 to annotate and functionally map variants in the meta-analysis. Genome-wide significance was defined using a threshold of P<5×10^-8. Genomic risk loci included variants correlated with the most significant variant at R²>0.6 using the 1000G Phase 3 EAS reference panel. Genome positions are reported in build hg37 for index variants. To designate a locus as previously known or potentially novel, the index variants, or the most significant variants in each locus, were at least 1 Mb in distance from a previously discovered genome-wide significant variant associated with the trait. Previously discovered variants were compiled from Shrine et al³⁴.

Comparison of pleiotropic effect analysis

We used S-LDSC to estimate heritability and genetic correlation of and between phenotypes using summary statistics of EUR populations. We then applied the linemodels package (https://github.com/mjpirinen/linemodels) to the GWAS summary statistics of three respiratory diseases (asthma, COPD, and lung cancer) and the three featured environmental and genetic risk factors (smoking status, smoking heaviness (Cigday), and FEV1/FVC) across three ancestry groups (AFR, EAS, and EUR) plus the corresponding meta-analysis results. We focused on comparisons for 12 pairs of traits selected from above. Three of these were between diseases phenotypes: asthma and COPD, asthma and lung cancer, COPD and lung cancer, and nine pairs were between disease phenotype and smoking or lung function. We considered variants significantly associated with either one of the traits being compared (Supplementary Table 28). We classified the variants into three classes (two when there is no variant associated with both traits) based on their association patterns: associated with trait 1 only, trait 2 only, and both. The classes were represented by linear trends whose slopes were estimated using an EM algorithm. Conditioning on these classes, we ran the linemodels package on the GWAS effect sizes and standard errors of overlapped variants of the two traits, where we set the scale parameters determining the magnitude of effect sizes to 0.2, the correlation parameters determining the allowed deviation from the lines to 0.99 as default, and the slope parameter to the estimates from the previous EM step. The membership probabilities in the three classes were computed separately for each variant by assuming that the classes were equally probable a priori. We assumed no overlapping samples between the two GWASs being compared and set the correlation of their effect estimators to 0. Confident associations are defined as having a posterior probability above 0.95.

Multi-trait analysis

We conducted multi-trait genome-wide association studies as implemented in MTAG v. 1.0.7 for each ancestry group (AFR, AMR, EAS, EUR) by combining the ancestry-specific GWAS summary statistics for GBMI asthma, GBMI COPD, GWMA lung cancer, spirometry meta-analysis, and GSCAN smoking behaviors. MTAG performs a joint analysis of GWAS results from related traits to improve the number of genetic loci identified and the predictive power of polygenic scores. For each population, we used the ancestry-specific LD reference panel from the gnomAD reference panels v2.1.1. For each ancestry-specific MTAG, we included traits with χ² > 1.02.

Study population

The All of Us Research Program is a longitudinal cohort study that has continuously enrolled US adults 18 years or older since May 2017. The program aims to engage in one million or more US participants and places a strong emphasis on including diverse populations that have traditionally been underrepresented in biomedical research. Details of the All of Us cohort have been previously described⁵⁴. In summary, participants of the program opt to provide self-reported data, linked health record data, and biospecimen data to be made available for research uses. The program’s primary objective is to build a resource to help researchers understand individual differences in biological, clinical, social, and environmental determinants of health and disease to advance precision health care.

Informed consent for all the participants in the All of Us Research Program are conducted in person or through an eConsent platform. The protocol was reviewed by the Institutional Review Board (IRB) of the All of Us Research Program. Data can be accessed through the All of US Research Workbench, a secure cloud-based analytic platform. Whole genome sequencing, genotyping array variant data, variant annotations, computed ancestry, and quality reports are accessible through the Controlled Tier of the AoU. This project is registered in the All of Us program under the workspace name “PRSxtra AoU”. In our analysis, we included individuals with whole genome data in the v7 Data Release, self-reported sex, and date of birth along with additional disease-dependent filtering criteria: For COPD and lung cancer, we excluded individuals who did not self-report report smoking status. For COPD, we additionally excluded individuals with homozygous polymorphism of SERPINA1, which encodes for a serine protease inhibitor alpha 1 antitrypsin, as this is a known risk allele associated with COPD⁵⁶. Participants were randomly split 70% for training and 30% for validation.

Phenotype ascertainment

We curated clinical phenotypes from All of Us using a combination of electronic health record data, and/or self-reported personal history data from the All of Us v7 Data Release. ICD codes for each phenotype and exacerbation are detailed in Supplementary Table 29-33.

We define smoking status and family history based on self-reported data. Previous smokers are individuals who smoked more than 100 pack years but do not currently smoke. Individuals who have never smoked more than 100 cigarettes are considered never smokers. Family history included mother, father, and siblings with the same record of disease as the participants. No family history included those who did not explicitly self-report a family history.

PRSxtra construction

We constructed PRSxtra in a three-layer process. Layer 1 consisted of performing multi-trait meta-analysis across related traits for each ancestry population as previously described. In layer 2, we used PRS-CSx, which leverages linkage disequilibrium across discovery samples to jointly model the genetic effects across populations via a shared continuous shrinkage prior. We used the default parameters on PRS-CSx on AFR, AMR, EAS, and EUR ancestry-specific meta-analysis results. Only HapMap3 variants—a set of 3 million variants compiled by the International HapMap Project which capture common patterns of variation in a variety of human populations—were included in calculating scores. In layer 3, we use ridge regression, as implemented by the “glmnet” R package⁵⁷, to jointly model the 39 standardized PRS (with mean 0 and variant 1) generated in layer 2 to construct an ancestry-specific PRSxtra for the three disease phenotypes: COPD, asthma, and lung cancer. In ridge regression, we used 10-fold cross-validation and minimum lambda value to estimate the weights of each PRS. PRSxtra was validated in the held-out multi-ancestry cohort from AoU.

PRS construction

As a baseline comparison to PRSxtra, we derived trait- and ancestry-matched PRS using PRS-CS, which uses a Bayesian regression framework to infer posterior effect sizes of SNPS. We ran PRS-CS with default parameters on the trait- and ancestry-specific GWAS summary statistic (Supplementary Table).

Statistical analysis

For COPD, asthma, and lung cancer, we placed individuals into bins by their PRS and PRSxtra deciles. In each decile, we calculated the prevalence of disease and disease exacerbation. We calculated the risk of disease and exacerbation for each decile compared to the lowest decile of PRS and PRSxtra using logistic regression models. We evaluated the performance of predicting diseases based on PRS and PRSxtra alone, as well as in a joint multi-variable model with covariates. The baseline model included age and sex. We then subsequently added clinical risk factors (smoking, family history), and PRS or PRSxtra. We evaluated the predictive performance of each model using the area under the receiving operating curve. In the full population, we used the trait- and EUR-specific PRS as the baseline for comparison. All statistical analyses were two-sided and performed with the use of R software, version 3.5 (R Project for Statistical Computing).