Identification of factors directly linked to incident chronic obstructive pulmonary disease: A causal graph modeling study

Abstract

Background: Beyond exposure to cigarette smoking and aging, the factors that influence lung function decline to incident chronic obstructive pulmonary disease (COPD) remain unclear. Advancements have been made in categorizing COPD into emphysema and airway predominant disease subtypes; however, predicting which healthy individuals will progress to COPD is difficult because they can exhibit profoundly different disease trajectories despite similar initial risk factors. This study aimed to identify clinical, genetic, and radiological features that are directly linked-and subsequently predict-abnormal lung function.

Methods and findings: We employed graph modeling on 2,643 COPDGene participants (aged 45 to 80 years, 51.25% female, 35.1% African Americans; enrollment 11/2007-4/2011) with smoking history but normal spirometry at study enrollment to identify variables that are directly linked to future lung function abnormalities. We developed logistic regression and random forest predictive models for distinguishing individuals who maintain lung function from those who decline. Of the 131 variables analyzed, 6 were identified as informative to future lung function abnormalities, namely forced expiratory flow in the middle range (FEF25-75%), average lung wall thickness in a 10 mm radius (Pi10), severe emphysema, age, sex, and height. We investigated whether these features predict individuals leaving GOLD 0 status (normal spirometry according to Global Initiative for Obstructive Lung Disease (GOLD) criteria). Linear models, trained with these features, were quite predictive (area under receiver operator characteristic curve or AUROC = 0.75). Random forest predictors performed similarly to logistic regression (AUROC = 0.7), indicating that no significant nonlinear effects were present. The results were externally validated on 150 participants from Specialized Center for Clinically Oriented Research (SCCOR) cohort (aged 45 to 80 years, 52.7% female, 4.7% African Americans; enrollment: 7/2007-12/2012) (AUROC = 0.89). The main limitation of longitudinal studies with 5- and 10-year follow-up is the introduction of mortality bias that disproportionately affects the more severe cases. However, our study focused on spirometrically normal individuals, who have a lower mortality rate. Another limitation is the use of strict criteria to define spirometrically normal individuals, which was unavoidable when studying factors associated with changes in normalized forced expiratory volume in 1 s (FEV1%predicted) or the ratio of FEV1/FVC (forced vital capacity).

Conclusions: This study took an agnostic approach to identify which baseline measurements differentiate and predict the early stages of lung function decline in individuals with previous smoking history. Our analysis suggests that emphysema affects obstruction onset, while airway predominant pathology may play a more important role in future FEV1 (%predicted) decline without obstruction, and FEF25-75% may affect both.

Fig 1. Overview of the methods used to predict GOLD 0 stage progression.

(A) Shows the COPD disease axes used to categorize individuals into different GOLD stages. (B) Depicts the transition in GOLD stage between the first and second COPDGene visit (a 5-year difference on average); 79.6% of individuals maintain GOLD 0 status and the remaining progress to more severe disease. (C) Details the workflow for data processing. Combined COPDGene data sets include demographic, spirometric, radiological, genetic, and survey information. Baseline variables only include measurements from the first visit. Imputation was performed using K-nearest neighbors. Categorical variables with more than 3 levels were merged. Correlation was determined using coefficient of determination for continuous variables and Cramer’s V-score for categorical variables. (D) Lays out the process for parameter tuning. Visit 1 data is split into a testing and training set, the latter is further divided for k-fold cross validation. For each fold and choice of ⍺ (graph sparsity), FCImax is used to determine variables causally linked to COPD progression (i.e., the Markov blanket). These variables are then used as predictors in a classification model where ⍺ is chosen by maximizing AUROC. AUROC, area under receiver operator characteristic curve; FEV₁, forced expiratory volume in 1 s; FVC, forced vital capacity; GOLD, Global Initiative for Obstructive Lung Disease.

Fig 2. Discovered Markov blanket that maximizes AUROC.

The Markov blanket encompasses every variable within the data set that can be used to infer information about our target variable: change in GOLD0 status between the first and second visit (ΔGOLD0). Each node in the graph corresponds to a measured variable in COPDGene and each edge represents a possible causal relationship that satisfies every conditional independence test performed. Arrows between variables show a direct causal link and/or an unmeasured latent confounder that causes both variables. (A) Depicts the optimal Markov blanket for the data set with variables related to FEV₁ and FVC removed. (B) Shows the same model, but all spirometric variables are removed. CT, computed tomography; FEF_25-75%, forced expiratory flow in the middle range; GOLD, Global Initiative for Obstructive Lung Disease; Perc15, 15th percentile cut-off for CT lung density in Hounsfield units; Pi10, average lung wall thickness in 10 mm radius.

Fig 3. Markov blanket composition across different graph sparsities.

The horizontal axis shows the log-transformed ⍺ value that controls graph sparsity (larger values lead to denser graphs), and the vertical axis shows variables that appeared in a Markov blanket. The color gradient counts the number of cross validating folds a variable appeared in the Markov blanket at that given graph sparsity. A lighter color indicates a variable appeared more frequently in the Markov blanket. (A) Shows the Markov blanket composition for the model with some spirometry measurements included. (B) Contains variables from the Markov blanket with no spirometry measurements included (showing the top 20, including ties). BMI, body mass index; CT, computed tomography; FEF_25-75%, forced expiratory flow in the middle range; Perc15: 15th percentile cut-off for CT lung density in Hounsfield units; Pi10: average lung wall thickness in 10 mm radius; SGRQ, St George Respiratory Questionnaire.

Fig 4. Classifier performance for predicting change in GOLD 0 status.

Receiver operator characteristic curves (ROC) illustrate model performance for FCImax + Logistic Regression and random forest. The shaded areas in the training data set designate ±1 standard deviation based on 10× cross-validations. Area under the curve measurements are displayed in the legend (no significant difference). (A–C) Include “limited spirometry” models. (D–F) Include “no spirometry” models. AUROC, area under the receiver operator characteristic curve; GOLD, Global Initiative for Obstructive Lung Disease.

Fig 5. Classifier performance for predicting change in GOLD 0 status in the SCCOR cohort.

Plotted are the ROC curves for COPDGene test data set and the external validation SCCOR cohort. The final model predicts change in GOLD 0 status without Pi10 as a predictor in the limited spirometry model. Each legend entry provides the area under the receiver operator curve (AUROC). AUROC, area under the receiver operator characteristic curve; SCCOR, Specialized Center for Clinically Oriented Research.

Fig 6. Exploring variables importance to the random forest model (“limited spirometry” model).

(A) Variables are ordered by importance using mean absolute Shapley values. (B–D) The distribution of Shapley values across measured demographics (age, biological sex, and race). Positive Shapley values on vertical axes indicate the random forest model was more likely to predict that individual to leave the GOLD 0 status (and vice versa). In (B) we do not display 8 (out of 2,114) individuals that were <45 years old at baseline, since they did not match the inclusion criteria of COPDGene. (E) SNP contributions to the random forest model prediction. Colors differentiate individuals with and without a given SNP. CT, computed tomography; FEF_25-75%, forced expiratory flow in the middle range; GOLD, Global Initiative for Obstructive Lung Disease; Pi10, average lung wall thickness in 10 mm radius; SGRQ, St George Respiratory Questionnaire; SNP, single-nucleotide polymorphism.