Genetic Modifiers of Cystic Fibrosis Lung Disease Severity: Whole-Genome Analysis of 7,840 Patients

Abstract

Rationale: Lung disease is the major cause of morbidity and mortality in persons with cystic fibrosis (pwCF). Variability in CF lung disease has substantial non-CFTR (CF transmembrane conductance regulator) genetic influence. Identification of genetic modifiers has prognostic and therapeutic importance. Objectives: Identify genetic modifier loci and genes/pathways associated with pulmonary disease severity. Methods: Whole-genome sequencing data on 4,248 unique pwCF with pancreatic insufficiency and lung function measures were combined with imputed genotypes from an additional 3,592 patients with pancreatic insufficiency from the United States, Canada, and France. This report describes association of approximately 15.9 million SNPs using the quantitative Kulich normal residual mortality-adjusted (KNoRMA) lung disease phenotype in 7,840 pwCF using premodulator lung function data. Measurements and Main Results: Testing included common and rare SNPs, transcriptome-wide association, gene-level, and pathway analyses. Pathway analyses identified novel associations with genes that have key roles in organ development, and we hypothesize that these genes may relate to dysanapsis and/or variability in lung repair. Results confirmed and extended previous genome-wide association study findings. These whole-genome sequencing data provide finely mapped genetic information to support mechanistic studies. No novel primary associations with common single variants or rare variants were found. Multilocus effects at chr5p13 (SLC9A3/CEP72) and chr11p13 (EHF/APIP) were identified. Variant effect size estimates at associated loci were consistently ordered across the cohorts, indicating possible age or birth cohort effects. Conclusions: This premodulator genomic, transcriptomic, and pathway association study of 7,840 pwCF will facilitate mechanistic and postmodulator genetic studies and the development of novel therapeutics for CF lung disease.

Keywords: GWAS/TWAS; cystic fibrosis; lung disease severity; pathway analyses; whole-genome sequencing.

Figure 1.

Distributions of the Kulich normal residual mortality-adjusted (KNoRMA) age-adjusted lung function phenotype by site cohort. The line inside each box is the median KNoRMA, and the box represents the interquartile range, or distance between the first and third quartiles (the 25th and 75th percentiles). Violin plots show the overall population distribution. Sample sizes are shown, with symbol areas proportional to sample size. CGS = Canadian CF Gene Modifier Study (population-based); FrGMS = French CF Gene Modifier Consortium (population-based); JHU = Johns Hopkins University (twin-siblings design); UNC = University of North Carolina (extremes of phenotype); UW = University of Washington (longitudinal study for effect of Pseudomonas aeruginosa acquisition on lung disease).

Figure 2.

Genetic loci significantly associated with Kulich normal residual mortality-adjusted (KNoRMA) lung phenotype. Genome-wide Manhattan plot of associations with KNoRMA in all 7,840 patients. Red line shows genome-wide significance of P < 5 × 10⁻⁸. Blue line shows a suggestive significance of P < 5 × 10⁻⁷.

Figure 3.

Forest plots for SNP association effect size by cohort at significant loci. β (coefficient) refers to the average change in Kulich normal residual mortality-adjusted (KNoRMA) phenotype for each copy of the protective allele. Square sizes are proportional to the sample size (n) of each cohort, and the line segments are 95% confidence intervals of each β. The most significant SNP from each locus was chosen. For each SNP, the protective allele is listed on the left, and frequencies of the protective alleles are shown in parentheses. Cohorts are arranged by increasing mean age (at KNoRMA). CFGP = Cystic Fibrosis Genome Project; CGS = Canadian CF Gene Modifier Study (population-based); FF = patients with CF who are F508del homozygous in the CFTR gene; FrGMS = French CF Gene Modifier Consortium (population-based); JHU = Johns Hopkins University (twin-siblings design); non-FF = patients with CF who are not homozygous for F508del in the CFTR gene; UNC = University of North Carolina (extremes of phenotype); UW = University of Washington (longitudinal study for effect of Pseudomonas aeruginosa acquisition on lung disease).

Figure 4.

Significant genes based on transcriptome-wide association (TWAS) evidence for expression versus lung function (Kulich normal residual mortality-adjusted phenotype). Genome-wide Manhattan plot of TWAS associations. The red line corresponds to transcriptome-wide false-discovery q lower than 0.10, with significant genes labeled. Red-colored text corresponds to increased expression associated with improved lung function, and blue-colored text corresponds to increased expression associated with decreased lung function. Regions of significant genome-wide phenotype–genotype association are marked with black arrows on the x-axis.

Figure 5.

Genes that drive core enrichment-significant results for this branching morphogenesis pathway (Gene Ontology 0048754). This VEGAS2 analysis Gene Set Enrichment Analysis plot includes 32 genes that are in three key signaling pathways (Shh [25]; TGFb [26]; and Wnt [27]) for lung development (including branching morphogenesis) and/or interact with genes in those three signaling pathways and/or have other roles in lung development (shown in bold). The 18 genes that are associated with lung repair and/or play a role in molecular pathogenic aspects of lung disorders (e.g., chronic obstructive pulmonary disease, asthma, lung fibrosis, cellular morphogenesis) are shown with asterisks. The remaining 11 genes are reported to have a role in development and/or morphogenesis in other tissues.