Abstract
Rationale Idiopathic pulmonary fibrosis (IPF) is a devastating disease characterized by limited treatment options and high mortality. Novel therapies and prognostic biomarkers are needed.
Objective To identify and validate molecular determinants of IPF survival.
Methods A staged genome-wide association study (GWAS) was performed using paired genomic and survival data. Stage I cases were drawn from centers across the US and Europe and stage II cases from Vanderbilt University. Cox proportional hazards regression was used to identify gene variants associated with differential transplant-free survival (TFS). Stage I variants with nominal significance (p<5×10−5) were advanced for stage II testing and meta-analyzed to identify those reaching genome-wide significance (p<5×10−8). Downstream analyses were performed for genes and proteins associated with variants reaching genome-wide significance.
Main Results After quality controls, 1481 stage I cases and 397 stage II cases were included in the analysis. After filtering, 9,075,629 variants were tested in stage I, with 158 meeting advancement criteria. Four variants associated with TFS with consistent effect direction were identified in stage II, including one in an intron of proprotein convertase subtilisin/kexin type 6 (PCSK6) reaching genome-wide significance (HR 4.11; 95%CI 2.54-6.67; p=9.45×10−9). PCSK6 protein was highly expressed in IPF lung parenchyma and negatively correlated with survival. Peripheral blood PCSK6 gene expression and plasma concentration were associated with reduced transplant-free survival.
Conclusions We identified four novel variants associated with IPF survival, including one in PCSK6 that reached genome-wide significance. Downstream analyses suggested that PCSK6 protein may serve as prognostic biomarker in IPF and potential therapeutic target.
Introduction
Idiopathic Pulmonary Fibrosis (IPF) is a devastating disease characterized by progressive lung scarring and poor survival.1,2 Two anti-fibrotic therapies have been approved for the treatment of IPF after randomized controlled trials demonstrated efficacy in slowing lung function decline.3,4 Despite this advance, outcomes remain poor and anti-fibrotic therapy appears to provide only modest survival benefit.5 To improve IPF outcomes, novel therapeutic targets are needed.
We and others have identified molecular IPF risk factors through unbiased investigation of the genome, transcriptome, and proteome.6-15 Among the strongest molecular determinants of IPF is a common variant in the promoter region of MUC5B, which increases the odds of developing IPF by nearly 5-fold per risk allele.6-9 Despite this strong association with IPF onset, the MUC5B promoter was paradoxically associated with improved survival.16 Few other susceptibility-associated gene variants have been shown to reliably predict differential IPF survival, suggesting that molecular determinants of IPF susceptibility and progression may have limited overlap.
To better understand molecular drivers of IPF progression and identify new therapeutic targets, we conducted a two-stage, multi-center, international genome-wide association study (GWAS) of IPF survival, followed by downstream analysis of genes and proteins associated with top survival-associated variant to determine whether these circulating biomarkers also predicted differential survival.
Methods
Cohorts and case selection
All patients provided informed consent for research blood draw in accordance with protocols approved by the institutional review board at each participating institution. GWAS stage I cases consisted of unrelated IPF patients of European ancestry from three previously described case-control GWAS datasets from the United States (US),9 United Kingdom (UK),7 and a combined cohort from the US, UK and Spain (UUS)6. Available outcome data was gathered for all cases meeting international consensus criteria for IPF17 and survival plotted for individual cohorts within each dataset. Patients without available outcome data were excluded, as were clinical trial cohorts due to short follow-up (Supplemental Methods). Stage II cases consisted of previously described, unrelated IPF patients of European ancestry from Vanderbilt University.18
Genotyping and quality control
Genotypes were generated for stage I cases using SNP (single nucleotide polymorphism) genotyping arrays (Thermo Fisher Scientific, Waltham, MA) according to previously described methods.6,7,9 Imputation for stage I cases was performed using the Michigan Imputation Server using the Haplotype Reference Consortium panel (v1.1 2016). Genotypes for the stage II cases were determined by whole genome sequencing, as previously described.18 Stringent quality control measures were applied with a two-tier variant filtering scheme as follows: those ranging 0.5≤MAF<1% were retained when imputation R2≥0.8; and those with MAF≥1% when imputation R2 ≥ 0.5. In addition, variants deviating from Hardy-Weinberg equilibrium (p<1.0×10−6) were removed.
Genome-wide survival analysis
The primary endpoint assessed was transplant-free survival (TFS), defined as the time in months from site-determined date of IPF diagnosis to event (death or lung transplant) or censoring date. Variants associated with differential TFS were identified using a multivariable Cox proportional hazards regression model adjusted for age, sex, center, and first ten genetic principal components. Variant genotypes were treated as a continuous variable with each patient having an imputed genotype dosage between zero and two risk alleles. Variants present in at least two datasets (US, UK and UUS) were meta-analyzed using METAL (v2011-03-25) to generate stage I results.
Variants nominally associated with TFS in stage I were defined as those with Wald p<0.05 in at least two datasets with the same direction of effect and p<5.0×10−5 in stage I meta-analysis. Conditional analysis of these SNPs to deduce their independence was performed with GCTA-COJO v1.26. The proportional hazards assumption was then assessed for each independent variant meeting advancement criteria by testing whether Schoenfeld residual rank varied by genotype strata. Variants that satisfied the proportional hazards assumption were advanced for stage II testing. Stage I and II cases were then meta-analyzed using METAL with the genome-wide significance threshold set at p<5.0×10−8. In silico assessments were used to infer the biological effect of variants associated with TFS after stage II testing.
PCSK6 tissue expression
Formalin-fixed paraffin-embedded human lung tissue sections obtained from patients with IPF undergoing surgical lung biopsy were compared to control subjects undergoing lung resection for malignancy, with sections distal to areas of malignancy utilized.Immunohistochemistry was performed using standard methods (Supplemental Methods) and mean staining intensity of PCSK6 protein was compared between IPF cases and non-IPF controls using a Mann-Whitney U-test.
PCSK6 clinical outcome association
The association between circulating PCSK6 gene expression and TFS was assessed using three previously published microarray datasets from the COMET trial, Imperial College and the University of Chicago, which were analyzed separately with results meta-analyzed and presented as a forest plot.19 Circulating plasma PCSK6 protein concentration was then determined in patients with IPF from UC-Davis and UChicago (Supplemental Methods), log2 transformed, and tested for TFS association using Cox proportional hazards regression.20 The proportional hazards assumption was satisfied for all downstream survival analyses.
Results
Case selection for stage I
Patients comprising the US cohort included those from the University of Chicago (n=118) and University of Pittsburgh (n=200) (Figure E1, Figure E2). Those comprising the UK cohort included patients from the University of Edinburgh (n=119), Trent Lung Fibrosis Study (n=210), a subset of those participating in the prospective, multi-center PROFILE study (NCT 01134822) (n=175), and aggregated patients from smaller UK centers (Hull and Papworth) (n=61) (Figure E1, Figure E3). Patients comprising the UUS cohort included those from the University of Chicago (independent of those from the US cohort, n=187), PROFILE study (independent of those in the UK dataset, n=299), University of California (Davis and San Francisco) (n=84) and aggregated patients from centers in Spain (n=28) (Figure E1, Figure E4).
Baseline characteristics and outcomes
Following phenotypic exclusions, 1481 patients comprising stage I were included in the analysis. These included 318 patients from the US dataset, 565 from the UK dataset, and 598 from the UUS dataset. Baseline characteristics for each dataset are shown in Table 1. The mean age ranged from 67 to 72 years and males comprised 71-75% of each dataset. The mean percent predicted FVC and diffusion capacity of the lung for carbon monoxide (DLCO) was lowest in patients comprising the US dataset and highest in the UK dataset. A majority of patients in each dataset were classified as gender, age, physiology21 (GAP) stage I or II. Death events were highest in the UK cohort and lowest in the UUS cohort (Figure E5). Median survival was 39.3 months in the US dataset, 53.2 months in the UK dataset and 40.6 months in the UUS dataset. Median survival was 48 months in the Vanderbilt University validation cohort (Figure E6).
Genome-wide survival analysis
After filtering, 7,873,835 variants in the US dataset, 8,591,398 variants in the UK dataset, and 8,620,496 variants in the UUS dataset were tested for TFS association. Quintile-quintile plots for each stage I cohort suggested acceptable inflation (Figure E7). After stratifying stage I cohorts by minor allele frequency (MAF), inflation was higher for rare variants compared to low and high frequency variants, but within an acceptable range for each group (λ<1.1) (Figure E8). For meta-analysis, 9,075,629 variants were tested for TFS association in the aggregated stage I cohorts. One hundred and sixty-one independent SNPs had Wald p<0.05 in at least two datasets with the same direction of effect and p<5.0×10−5 in stage I meta-analysis (Figure 1). Of those, 158 satisfied the proportional hazards assumption and advanced for stage II testing (Table E1).
Manhattan plot of stage I gene variants associated with IPF survival. SNPs reaching nominal significance (p<5×10−5) fall above the blue line and those reaching genome-wide significance (p<5.0×10−8) fall above the red line.
Genotype data was available in the Vanderbilt University cohort for 154 of the 158 variants advanced from stage I. Six variants were associated with TFS in the Vanderbilt University cohort at p<0.05, including four with consistent effect direction that strengthened in TFS association after meta-analysis (Table 2; Table E2). These four were rs184498750 near Succinate-CoA Ligase GDP/ADP-Forming Subunit Alpha (SUCLG1), rs60514164 near ubiquitin-conjugating enzyme E2Q family member 2 (UBE2Q2), rs35647788 in an intron of Proprotein Convertase Subtilisin/Kexin Type 6 (PCSK6), and rs3893252 in an intron of Deleted In Azoospermia-Associated Protein 1 (DAZAP1) (Table 2). Of these, rs35647788 (PCSK6) showed strong TFS association in stage I (HR 4.76; 95% CI 2.62-8.64; p=2.96×10−7) and stage II (HR 3.12; 95% CI 1.37-7.11; p=6.70×10−3) cohorts and crossed the genome-wide significance threshold in meta-analysis (HR 4.11; 95% CI 2.54-6.67; p=9.45×10−9).With the exception of rs60514164 (UBE2Q2) (MAF=8%), these SNPs were low frequency, with MAF of ∼1% in the study population. Regional association plots for each of the four variants are shown in Figure E9.
In multivariable analysis, each variant except rs3893252 (DAZAP1) maintained survival association after adjustment for relevant confounders of IPF survival (Table E3). Among patients with the rs35647788 (PCSK6) variant, all were heterozygotes (Table E4) and were evenly distributed across centers comprising the UK and UUS cohorts. No rs35647788 (PCSK6) variants were observed in the US cohort despite good imputation quality (r2=0.74). In sensitivity analysis of the PCSK6 variant, results were consistent when censoring transplants (Table E5). In silico testing revealed functional effects for each of the four variants associated with TFS after stage II meta-analysis (Table E6). None of the variants identified had known association with fibrotic lung disease.
PCSK6 tissue expression
Morphologic assessment of histological sections from lung tissue in patients with IPF were compared with control subjects without fibrotic lung disease. In IPF lung, cytoplasmic PCSK6 expression localized to ciliated epithelial cells and alveolar epithelial cells and was markedly higher than PCSK6 expression in non-IPF control sections (Figure 2). Western Blot confirmed the presence of only a single PCSK6 band (Figure E10). Relative staining intensity was two-fold higher in IPF lung samples (n=86) compared with non-IPF controls (n=9) (p<0.001) (Figure 2).
PCSK6 immunohistochemistry showed increased cytoplasmic PCSK6 expression in ciliated epithelial cells (A) and alveolar epithelial cells (B-C) compared with normal lung control sections (D-E). Parallel IPF section sections (F-G) confirm increased PCSK6 staining (F) when compared to control section (G). Human kidney positive control (H) is provided for reference. Comparison of relative PCSK6 staining intensity (I) demonstrated significantly higher mean intensity in IPF lungs when compared to non-IPF control lungs (p<0.001). Magnification is 20um magnification for panels A-E and 50um for panels F-G.
PCSK6 clinical outcome association
Increased PCSK6 protein staining score was associated with reduced TFS in those with available survival data (n=71), with staining scores above the median associated with greater than 2-fold increased risk of death or lung transplant (HR 2.41; 95% CI 1.12-5.16; p=0.024) (Figure 3a). When assessing PCSK6 gene expression in the COMET (n=75), Imperial College (n=55) and University of Chicago (n=45) cohorts, increasing PCSK6 expression was associated with increased mortality risk in each cohort, with each one-unit increased associated with greater than three-fold increased risk of death or lung transplant in meta-analysis (HR 3.43; 95% CI 1.62-7.25; p=0.0012) (Figure 3b). When assessing PCSK6 plasma concentration in patients with IPF from UC-Davis (n=138) and UChicago (n=181), increasing plasma PCSK6 concentration was associated with reduced TFS, with each one-unit change in log-transformed plasma concentration associated with a nearly 50% increase in outcome risk (HR 1.47; 95% CI 1.14, 1.89; p=0.0031). These results were consistent across UC-Davis (HR 1.47, 95% CI 1.14-1.89; p=0.08) and UChicago (HR 1.34, 95% CI 0.92-1.94; p=0.12) cohorts. After stratification of the combined cohort by tertiles, those with PCSK6 concentration in the highest tertile displayed significantly worse survival than those in the second and third tertiles (p=0.0018) (Figure 3c).
Relationship between PCSK6 and clinically relevant IPF endpoints. Higher PCSK6 staining intensity (a), peripheral blood gene expression (b) and circulating plasma protein concentration (c) are associated with reduced transplant-free survival.
Overlap between IPF risk and transplant free survival
Variants previously associated with IPF risk6-10 were investigated for outcome association (Supplementary Methods). None of the 15 genetic variants with previously associated with IPF risk6-10 were associated with TFS after Bonferroni correction (p=0.0033) (Table E7). As previously reported,16,22 individuals with the MUC5B promoter polymorphism (rs3570590) displayed better overall survival, though this did not reach significance after adjustment for multiple testing. None of the four validated survival variants were associated with differential IPF risk (Table E8). When combing the effect of thousands of IPF risk variants in a polygenic risk score, this risk score was not significantly associated with TFS for any significance threshold used (Figure E11), again suggesting that variants that affect disease risk may have little impact on survival times after diagnosis.
Discussion
In this investigation, we conducted the first GWAS of IPF survival, identifying a variant intronic to PCSK6 that associated with differential TFS at genome-wide significance in two independent IPF cohorts totaling nearly two thousand patients. We subsequently found that PCSK6 protein was highly expressed in IPF lung tissue, localizing to alveolar epithelial cells, which play a key role in IPF onset and progression.23 Finally, we found that PCSK6 lung staining, peripheral blood gene expression and circulating plasma concentration negatively correlated with TFS across independent IPF cohorts. To our knowledge, this study is the first to systematically identify gene variants associated with IPF survival and the first to identify PCSK6 as a potential therapeutic target in IPF.
PCSK6, also called PACE4, encodes a widely expressed calcium-dependent serine endoprotease with known expression in alveoli and alveolar macrophages.24 PCSK6 is a critical mediator of TGF-β processing and is crucial for reproduction, embryological development and blood pressure regulation.25-29 A PCSK6 gene variant has been implicated in the development of hypertension29 and dysregulated PCSK6 gene expression has been linked to vascular disease,30-32 and cardiac remodeling following myocardial ischemia.33 These cardiovascular remodeling effects make PCSK6 of particular relevance to IPF, as PCSK6 overexpression can lead to increased collagen I and III deposition, TGF-β activation and extracellular matrix formation,33,34 which are cardinal features of IPF pathogenesis.2 Additionally, PCSK6 may bind tissue inhibitors of metalloproteinases (TIMPs),35 potentially counteracting pro-fibrotic metalloproteinases.36
PCSK6 dysregulation has also been implicated in the development of cancer of the lung,37 breast,38 ovary39 thyroid40 and prostate.41 PCSK6 has been shown to regulate apoptosis in prostate cancer42 and pancreatic cancer43 and also linked to increased cancer cell invasiveness by enhancing bioactivity of matrix metalloproteinases and cytokines.44 Accordingly, PCSK6 has been proposed as an anti-tumor therapeutic target45,46 and a bioavailable formulation of an anti-PCSK6 molecule is currently under investigation.47 In vitro PCSK6 inhibition has already been shown to reduce fibroblast proliferation, migration and invasion in rheumatoid arthritis-associated synovitis,48 suggesting it as a candidate target for the treatment of pulmonary fibrosis.
The PCSK6 variant identified in this study is rare and unlikely to singularly explain subsequent gene expression and protein findings. In silico analyses identified several nearby regulatory elements which may include common functional variants, with smaller effects, in linkage disequilibrium with this PCSK6 variant. While beyond the scope of this investigation, further research is needed to explore other rare variants in PCSK6 that may have an additive effect with the one identified, along with regulatory elements, expression quantitative trait loci and more complex structural variants that may be contributing to our findings.
This study also highlights important differences between genomic determinants of IPF susceptibility and survival, suggesting that genes involved in disease onset may be independent of those driving disease progression. No survival-associated variant identified in this analysis was associated with IPF risk. While replicating prior IPF risk association for the MUC5B promoter polymorphism, we found only weak association with improved survival, an observation that may be influenced by index event bias.22 As none of the survival-associated variants showed an association with disease risk, it is unlikely these survival results are affected by index event bias. These findings have potential implications for drug development, as genes associated with IPF survival may represent more effective therapeutic targets than those associated with IPF onset.
This study has several limitations. First, we acknowledge the relatively small sample size used to conduct this analysis. Despite being the largest genomic IPF outcome analysis reported to date, the modest size of this cohort limited our ability to identify higher frequency SNPs with modest effect sizes. Given sample size constraints for this rare disease, we pursued a two-stage approach with meta-analysis of candidate variants rather than a discovery/replication approach, which would have required substantially higher sample sizes in each cohort. The consistent effect association across cohorts and genome-wide significance for the PCSK6 variant after meta-analysis increases confidence that this represents a true association, as does the downstream clinical outcome analysis showing PCSK6 gene expression and protein concentration to be associated with differential TFS. Next, there were likely differences between cohorts used to conduct this study. Survival was best among patients in the UUS cohort, which may reflect selection bias, as these patients had less follow-up compared with those in the US and UK datasets. Additionally, a large proportion of patients comprising the UUS dataset were recruited after the US approval of pirfenidone and nintedanib, which might explain the better overall survival and patients recruited before 2012 may have been exposed to potentially harmful immunosuppression.49 The timing of study recruitment relative to IPF diagnosis was also unknown, which could have influenced results. Finally, we only assessed individuals of European ancestry. Validation of results is required in patients of other ancestral backgrounds.
Conclusion
Here we present results from the first GWAS of IPF survival conducted to date. This study sheds important light on the genetics of IPF progression and identified novel variants which may contribute to this process, including rs35647788 in an intron of PCSK6.Downstream analysis demonstrated PCSK6 protein lung staining, peripheral blood gene expression and circulating plasma concentration to be associated with reduced IPF survival, suggesting PCSK6 may serve as a potential therapeutic target in patients with IPF.
Data Availability
GWAS Summary statistics for this study can be accessed from https://github.com/genomicsITER/PFgenetics and will be made available at GWAS Catalog
Declaration of interests
JMO reports grants from the National Heart, Lung, and Blood Institute (NHLBI), American Lung Association and American Thoracic Society related to the submitted work and personal fees from Genentech, Boehringer Ingelheim, United Therapeutics, Lupin pharmaceuticals and AmMax Bio unrelated to the submitted work. PLM has received industry-academic funding from AstraZeneca and GSK and has received speaker and consultancy fees from Boehringer Ingelheim and Hoffman-La Roche outside the submitted work. JSK reports grants from the NHLBI and Pulmonary Fibrosis Foundation. AGN reports consultancy fees for Galapagos, Medical Quantitative Image Analysis and Boehringer Ingelheim, and fees for educational activities from Boehringer Ingelheim and UptoDate. AA reports grants from the National Heart, Lung, and Blood Institute (NHLBI), American College of Chest Physicians, and the Pulmonary Fibrosis Foundation, and speaking and advisory board fees from Genentech and Boehringer Ingelheim. MES has received grants from Boehringer Ingelheim and Galapagos and personal fees from Boehringer Ingelheim and Fibrogen. HP reports grants from Roche; personal fees from Boehringer Ingelheim, Roche and Pliant Therapeutics and is a Trustee for the charity Action for Pulmonary Fibrosis. IPH holds an NIHR Senior Investigator Award. AS and BLY are employees of Roche/Genentech with stock and stock options in Roche. TMM has received industry-academic funding from GlaxoSmithKline (GSK) R&D and UCB; and has received consultancy or speakers fees from Apellis, AstraZeneca, Bayer, Biogen Idec, Boehringer Ingelheim, Cipla, GSK R&D, InterMune, ProMetic, Roche, Sanofi-Aventis, Sanumed, and UCB. IN reports personal fees and consultancy, advisory board, speaker, and study contracts for Boehringer Ingelheim and Genentech from Boehringer Ingelheim and Genentech; personal fees from Immuneworks, Sanofi Aventis, Global blood therapeutics; and has a patent for Toll-interacting protein (TOLLIP) and pharmacogenetics pending. NK served as a consultant to Biogen Idec, Boehringer Ingelheim, Third Rock, Pliant, Samumed, NuMedii, Theravance, LifeMax, Three Lake Partners, Optikira, Astra Zeneca, Veracyte, Augmanity and CSL Behring, over the last 3 years, reports Equity in Pliant and a grant from Veracyte, Boehringer Ingelheim, BMS and non-financial support from MiRagen and Astra Zeneca. NK has IP on novel biomarkers and therapeutics in IPF licensed to Biotech. RGJ reports grants from GSK, UK MRC, Biogen, AstraZeneca, Galecto and NIHR; personal fees from Boehringer Ingleheim, Bristol Myers Squibb, Chiesi, Galapagos, GSK, Heptares, MedImmune, Pharmakea, Roche/Promedior, Pliant Therapeutics, and Veracyte; and is a Trustee for the charity Action for Pulmonary Fibrosis. LVW has received funding from GlaxoSmithKline. All other authors declare no competing interests.
Data Sharing Statement
GWAS Summary statistics for this study can be accessed from https://github.com/genomicsITER/PFgenetics and will be made available at GWAS Catalog
Acknowledgements
We thank all the patients that contributed biologic samples and clinical data for this research.
Footnotes
Funding NHLBI - K23HL138190 (JMO), R56HL158935 (JMO), K23HL150301 (JSK), K23HL146942 (AA), T32HL007605 (CTL), RO1HL130796 (IN), UG3HL145266 (IN)
Ministry of Science and Innovation (RTC-2017-6471-1; AEI/FEDER UE (CF)
Instituto de Salud Carlos III co-financed by the ERDF “A way of making Europe” from the EU (PI20/00876, CF), Instituto Tecnológico y de Energías Renovables (OA17/008 (CF) Consejería de Educación-Gobierno de Canarias and Cabildo Insular de Tenerife (BOC n.º 173, 24/08/2017 (JMLS)
UK National Institute for Health Research (NIHR) - RGJ holds an NIHR Research Professorship (RP-2017-08-ST2-014) and has MRC grants (G0901226 and MR/V00235X/1). TMM is supported by an NIHR Clinician Scientist Fellowship (NIHR reference CS-2013–13– 017). MDT is supported by an NIHR Senior Investigator Award (NIHR201371). The views expressed are our own and not necessarily those of the NHS, the NIHR, or the UK Department of Health.
British Lung Foundation Chair in Respiratory Research (C17-3) (TMM), GlaxoSmithKline/British Lung Foundation Chair in Respiratory Research (C17-1) (LWW) Wellcome Trust - WT202849/Z/16/Z (MDT), 221680/Z/20/Z (BGG). For the purpose of open access, the author has applied a CC BY public copyright license to any Author Accepted Manuscript version arising from this submission.
BREATHE — The Health Data Research Hub for Respiratory Health (MC_PC_19004) (LWW) Action for Pulmonary Fibrosis Research Fellows (RJA and PLM).
This research used the ALICE and SPECTRE High Performance Computing Facilities at the University of Leicester and TeideHPC (https://teidehpc.iter.es/en/home) at Instituto Tecnológico y de Energías Renovables.
Role of the Funding Source This study was funded by the agencies listed above, which had no role in the design or conduct of the study. All authors had full access to all data, and the corresponding author had final responsibility for the decision to submit for publication.
