Development and evaluation of a predictive algorithm for unsatisfactory response among patients with pulmonary arterial hypertension using health insurance claims data

Abstract

Objective

This study aimed to develop and validate a predictive algorithm for unsatisfactory response to initial pulmonary arterial hypertension (PAH) therapy using health insurance claims.

Methods

Adult patients with PAH initiated on a first PAH therapy (index date) were identified from Optum’s de-identified Clinformatics Data Mart Database (1/1/2010–12/31/2019). A random survival forest algorithm was developed using patient-month data and predicted the “survival function” (i.e. risk of not having unsatisfactory response) over time. For each patient-month observation, risk factors were assessed in the 12 months prior. Unsatisfactory response was defined as the first instance of (1) new PAH therapy, (2) PAH-related hospitalization or emergency room visit, (3) lung transplant or atrial septostomy, (4) PAH-related death or (5) chronic oxygen therapy initiation. To facilitate use in clinical practice, a simplified risk score was also developed based on a linear combination of the most important risk factors identified in the algorithm.

Results

In total, 4781 patients were included (median age = 69.0 years; 58.6% female). Over a median follow-up of 14.0 months, 3169 (66.3%) had an unsatisfactory response. The most important risk factors included in the algorithm were healthcare resource use (i.e. PAH-related outpatient visits, pulmonologist visits, cardiologist visits, all-cause hospitalizations), time since first PAH diagnosis, time since index date, Charlson Comorbidity Index, dyspnea, and age. Predictive accuracy was good for the full algorithm (C-statistic: 0.732) but was slightly lower for the simplified risk score (C-statistic: 0.668).

Conclusion

The present claims-based algorithm performed well in predicting time to unsatisfactory response following initial PAH therapy.

Keywords:

• Pulmonary arterial hypertension
• unsatisfactory response
• risk assessment
• machine learning
• predictive algorithm
• health-insurance claims
• combination therapy

Introduction

Pulmonary arterial hypertension (PAH) is a rare, progressive subtype of pulmonary hypertension (PH) characterized by the presence of pre-capillary PH¹. PAH has been traditionally defined as a mean pulmonary arterial pressure ≥25 mm Hg (although a threshold of >20 mm Hg has been recently proposed^2,3) pulmonary arterial wedge pressure ≤15 mm Hg, and pulmonary vascular resistance ≥3 wood units¹. It affects 15–50 persons per million within the United States (US), with a higher prevalence among women and a mean age at diagnosis of about 50 years^2,4–6. Symptoms of PAH are typically activity-related and include dyspnea, chest pain, fatigue, and dizziness^7,8. Due to the non-specific nature of these symptoms, patients often remain undiagnosed until the more advanced stages of the disease^8–10. As PAH progresses, symptoms may occur at rest, and patients are at increased risk of fatal right-heart failure^11,12. Despite recent advancements in treatment, the prognosis for patients diagnosed with PAH remains poor, with a 5-year survival rate of only about 60% in the US^4,13.

Treatment for PAH can be broadly categorized into general measures (e.g. dietary and lifestyle modifications), supportive therapies (e.g. diuretics and oral anticoagulants), calcium channel blockers, and PAH-specific drug therapies (i.e. targeting components of key PAH-relevant molecular pathways) ¹. The past two decades have seen a proliferation of pharmacological agents for the treatment of PAH. In particular, phosphodiesterase type 5 inhibitors (PDE5is) and endothelin-receptor antagonists (ERAs) are frequently used as first-line oral therapies for PAH¹. Prostacyclin pathway agents (PPAs) are generally used following disease progression, although a large proportion of patients with PAH are never treated with a PPA at any time^1,14,15. Additionally, combination therapy is recommended by current guidelines for the treatment of patients with PAH¹. Meta-analyses have shown that, compared with monotherapy, combination therapy significantly reduces the risk of progression by approximately 35% among patients with PAH^16,17.

Given that PAH is a progressive and life-threatening disease, delayed identification of patients at risk might lead to inadequate treatment and poor outcomes. A predictive algorithm that could identify patients at higher risk of an unsatisfactory response to initial PAH therapy might enable healthcare providers to target those who would benefit from early combination therapy⁹. Indeed, frequent and regular risk assessment is an essential aspect of PAH management per the guidelines of the European Society of Cardiology (ESC) and the European Respiratory Society (ERS), 6^th World Symposium on Pulmonary Hypertension (WSPH), and American College of Chest Physicians (CHEST) ^1,18,19. In current risk assessments, patients are typically stratified according to the risk of 1-year mortality using a comprehensive approach that takes multiple prognostic factors into account^1,20–24. For instance, the Registry to Evaluate Early and Long-Term PAH Disease Management (REVEAL) risk calculator was developed using data collected from patients enrolled in the REVEAL observational prospective registry^20,25,26. Risk factors for mortality at 1-year follow-up were identified by Cox regression models; these prognosis equations were then used to categorize patients into specific risk groups²⁵. To expedite risk assessment in the clinic, the REVEAL Lite 2 was developed, which uses only six variables including functional class (FC), vital signs (systolic BP [SBP] and heart rate), 6-min walk distance (6MWD), brain natriuretic peptide (BNP)/N-terminal prohormone of brain natriuretic peptide (NT-proBNP), and renal insufficiency (by estimated glomerular filtration rate [eGFR]) ²⁰.

The 2015 ESC/ERS guidelines recommend therapeutic approaches based on risk stratification, with the aim of achieving or maintaining a lower risk status among patients^1,14. However, current risk assessment tools, such as the REVEAL risk calculator, often require access to clinical values that are not widely available among patients (e.g. FC, 6MWD, NT-proBNP/BNP levels, measures of hemodynamics, and pulmonary function), thereby restricting their routine use in daily practice^1,20–24,27. To address this limitation, the present study aimed to develop and assess the performance of a predictive algorithm for unsatisfactory response to initial PAH therapy among patients diagnosed with PAH using health insurance claims data.

Methods

Data source

Patient-level data were retrieved from the Optum Clinformaticsⁱ Data Mart Database (1/1/2010–12/31/2019), which contains de-identified insurance claims from patients enrolled in commercial and Medicare Advantage health plans in the US. The information available in the database includes demographic characteristics, type of healthcare provider, diagnoses for all types of encounters within the network, detailed area of care during hospitalization, inpatient hospital procedures, inpatient medications, and pharmacy claims. The Optum database is de-identified and fully compliant with the patient confidentiality requirements of the Health Insurance Portability and Accountability Act (HIPAA). Therefore, no Institutional Review Board approval was required for the present study.

Study design

A retrospective study design was used to develop a predictive algorithm for unsatisfactory treatment response among patients with PAH who had been newly initiated on a first PAH therapy. The index date for a given patient was defined as the date of the first documented prescription claim for a PAH-specific PDE5i or ERA therapy. The baseline period was defined as the 12 months preceding the index date and served as a washout period to increase the likelihood of capturing the first PAH therapy as of the index date. The subsequent at-risk period spanned from the index date until the earliest of the end of data availability, the end of continuous healthcare plan enrollment, or death.

The data used to develop the algorithm were organized at the patient-month level, where an observation was included for each month between the index date and the end of the at-risk period or the first indicator of an unsatisfactory response to initial PAH therapy (see Study Outcome and Risk Factors). For each observation, the time to unsatisfactory response was assessed and risk factors were evaluated in the prior 12 months to allow for variation over time. This design allowed the algorithm to predict the time to unsatisfactory response at any given time point following the index date. Figure 1 depicts the study scheme.

Figure 1. Study design. Abbreviations. ERAs, endothelin-receptor antagonists; PAH, pulmonary arterial hypertension; PDE5is, phosphodiesterase type 5 inhibitors.

Study population

Eligible patients had ≥1 prescription fill for a PAH-specific PDE5i or ERA as part of the index regimen, which was defined as the treatment initiated on the index date along with any other treatment indicated for PAH initiated within a 60-day window after the index date. A 60-day window after the index date was used to account for any potential delays in initiation of certain agents from a combination therapy, for example, delays prior to receiving reimbursement authorization from health insurer or enrollment in a Risk Evaluation and Mitigation Strategy [REMS] program. Patients were also required to have ≥1 documented diagnosis of PAH (ICD-9-CM: 416.0, 416.8; ICD-10-CM: I27.0, I27.20, I27.21, I27.89) in an inpatient or emergency room setting or ≥2 documented diagnoses of PAH in an outpatient setting recorded during the baseline period or on the index date. Additionally, patients were required to be continuously enrolled in their healthcare plan for ≥12 months prior to the index date. Patients were excluded from the study if they were <18 years of age at the index date, had ≥1 prescription fill for a medication indicated for PAH (i.e. PDE5is, ERAs, PPAs, or riociguat) prior to the index date, received PPAs as part of their index regimen, or had ≥1 documented diagnosis of chronic thromboembolic pulmonary hypertension (CTEPH; ICD-9-CM: not available; ICD-10-CM: I27.24 available as of October 2017) at any time. Of note, PPAs were not considered for the index date, since they are typically used for more advanced stages of the disease, while the use of riociguat in monotherapy was not considered as it is not PAH-specific and could lead to the inclusion of patients with CTEPH in the study sample. All eligible observations were included in this retrospective analysis; therefore, sample size and power calculations were not performed.

Study outcome and risk factors

Patients were considered to have an unsatisfactory response if they had any of the following during the at-risk period: (1) initiation of a new treatment indicated for PAH not part of the index regimen (i.e. switch or add-on), (2) PAH-related or heart failure-related hospitalization or emergency room visit, (3) lung transplant or atrial septostomy, (4) PAH-related death, or (5) chronic oxygen therapy initiation (additional details are provided in Table 1). The time to unsatisfactory response was defined as the number of days between the index date and the date of the first occurrence of any of the components listed above. For each patient-month observation, the algorithm predicted the “survival function”, that is, the risk of not having the outcome at each subsequent time point. The predicted time to the outcome was then defined as the time at which the predicted survival function reached <50%. Candidate risk factors were selected based on a targeted literature review and included diagnosis codes, procedures codes, and pharmacological treatments received (Supplementary Table 1).

Table 1. Unsatisfactory clinical response during the at-risk period.

Components^a	Study sample
	N = 4781
(1) Initiation of a treatment indicated for PAH that was not part of the treatment regimen initiated on the index date, n (%)^b	932 (19.5%)
Agent(s) newly initiated, n (%)
ERA	340 (7.1%)
Ambrisentan	181 (3.8%)
Macitentan	115 (2.4%)
Bosentan	44 (0.9%)
PDE5i	329 (6.9%)
Tadalafil	197 (4.1%)
Sildenafil	60 (1.3%)
PPA	229 (4.8%)
Treprostinil	137 (2.9%)
Selexipag	61 (1.3%)
Iloprost	22 (0.5%)
Epoprostenol	11 (0.2%)
Riociguat	37 (0.8%)
(2) ≥1 PAH-related or heart failure-related hospitalization or ER visit, n (%)^c	2611 (54.6%)
(3) PAH-related or heart failure-related death, n (%)^d	933 (19.5%)
(4) Lung transplant or atrial septostomy, n (%)^e	24 (0.5%)
(5) Initiation of chronic oxygen therapy, n (%)^f	602 (12.6%)
Composite outcome
Any component	3169 (66.3%)
Duration of at-risk period (months), mean ± SD [median]	22.3 ± 23.7 [14.0]

Abbreviations. ER, emergency room; ERAs, endothelin-receptor antagonists; PAH, pulmonary arterial hypertension; PDE5is, phosphodiesterase type 5 inhibitors; PPA, prostacyclin pathway agent.

^aComponents were measured during the at-risk period. The at-risk period was defined as the period spanning from the index date until the earliest of the end of data availability, end of continuous healthcare plan enrollment, or death.

^bThe list of medications indicated for PAH included ERAs, PDE5is, PPA: prostacyclin pathway agents, and riociguat. The treatment regimen initiated on the index date was defined based on the agents received within the first 60 days of the index date.

^cPAH-related or heart failure-related hospitalization or ER visit was defined based on a diagnosis for PAH (ICD-9-CM: 416.0, 416.8; ICD-10-CM: I27.0, I27.20, I27.21, I27.89) or right heart failure (ICD-9-CM: 428.0; ICD-10-CM: I50.81) recorded in a medical claim in an inpatient or ER setting.

^dPAH-related or heart failure-related death was defined as having a diagnosis for PAH (ICD-9-CM: 416.0, 416.8; ICD-10-CM: I27.0, I27.20, I27.21, I27.89) or right heart failure (ICD-9-CM: 428.0; ICD-10-CM: I50.81) within 30 days of the date of death.

^eLung transplant was defined based on a lung or lung lobe transfer procedure (ICD-9-PCS: 33.5x, 33.6; ICD-10-CM: 0BYC, 0BYD, 0BYF, 0BYG, 0BYH, OBYJ, 0BYK, 0BYL, 0BYM; CPT: 32851, 32852, 32853, 32854, 33935). Atrial septostomy was defined based on an atrial septostomy procedure (ICD-9-PCS: 35.4x; ICD-10-PCS: 02QA4ZZ, 02QA3ZZ, 02QA0ZZ, 02B50ZZ, 02B53ZZ, 02B54ZZ; CPT: 92992, 92993, 33735, 33736, 33737).

^fInitiation of chronic oxygen therapy was defined based on medical claims for oxygen therapy (HCPCS: E0424, E0425, E0430, E0431, E0433, E0434, E0435, E0439, E0440, E0441, E0442, E0443, E0444, E0445, E0447, E0455, E1356, E1357, E1358, E1390, E1391, E1392, K0738, K0740, K0741, K0742, K0671, S8120, S8121), excluding claims recorded in inpatient or ER settings. It was further required that patients do not have a medical claim for oxygen therapy during the baseline period to be considered as initiating chronic oxygen therapy.

Statistical analyses

Descriptive statistics and univariate associations

Frequencies of the study outcome and its underlying components during the at-risk period were reported for patients meeting the inclusion and exclusion criteria. Univariate associations between the outcome and each candidate risk factor evaluated during the 12-month baseline period were assessed using Chi-Square test for categorical variables and Wilcoxon rank-sum test for continuous variables.

Algorithm development

The algorithm was developed using random survival forests, a tree-based machine learning strategy that can be used to build predictive models based on right-censored survival data²⁸. One advantage of random survival forests compared to traditional regression-based approaches for analyzing right-censored data is that it does not need to rely on strict underlying assumptions^28–30. For example, Cox regression models assume that the hazard ratio associated with a given risk factor should be constant across time (i.e. proportional hazards assumption), and that predictor variables have a linear relationship to outcomes. By contrast, random survival forests are nonparametric, and thus independent of model assumptions. Furthermore, common limitations of regression-based alternatives (e.g. Cox models), such as difficulty accommodating interactions between predictorsⁱⁱ, overfitting and inflated standard errors, are more naturally addressed by random survival forests. Consequently, random survival forests are well-suited for the analysis of high-dimensional right-censored survival data^28–30.

In the current random survival forests algorithm, an ensemble of 500 survival trees predicted a survival function for each patient-month observation, namely the probability that a patient does not have an unsatisfactory response at a given time point after that monthⁱⁱⁱ. The predictions of the 500 survival trees were averaged to produce a prediction for the random survival forest. Trees were developed using a bagging approach (i.e. bootstrap aggregating), where each survival tree was developed based on a random subsample of data (i.e. in-bag sample) and model performance was evaluated based on the observations that were left out (i.e. out-of-bag sample). Of note, each of the in-bag samples included no more than one observation for a given patient, and the remaining observations for patients included in a given in-bag sample were excluded from the corresponding out-of-bag sample.

Variable importance was determined based on the average contribution of each candidate risk factor to the algorithm’s performance. This was measured by permuting a single variable, refitting the model on the altered data, and calculating the decrease in model performance as measured by the C-statistic. This process was then repeated for all variables. Risk factors with zero or negative importance, which represented either variable with no prediction power or noise that reduced algorithm performance (non-predictive risk factors), were removed from the list of candidate risk factors¹¹. The algorithm was then retrained on the remaining subset of candidate risk factors until all risk factors had positive importance. The risk factors that were ultimately retained in the algorithm were then rank-ordered based on their importance for predicting time to unsatisfactory response.

The predictive accuracy of the final algorithm in out-of-bag samples was measured in terms of the concordance statistic (C-statistic), which indicates the degree of agreement between the observed and the predicted time to unsatisfactory response³¹. Based on prior literature, the C-statistic typically ranges between 0.6 and 0.75 for biomedical applications¹³. In addition, the sensitivity, specificity, negative predictive value, and positive predictive value of the algorithm were assessed at key time points based on the predicted time to the outcome.

Given that the random survival forest algorithm requires the averaging of predictions over 500 trees, this approach may not be easily applicable in real-world clinical practice. Therefore, a simplified linear risk score was also developed as a tool to assess the risk level for a given patient based on a numerical score (i.e. a higher score indicates a higher risk). It should be noted that this simplified risk score does not provide a prediction of the risk of unsatisfactory response at various time points (i.e. “survival curve”) for a given patient, as per the full algorithm. The simplified risk score was developed using an 80% random subsample (training set) and validated using the remaining 20% (validation set) ³². Cross-entropy optimization was used to maximize the predictive accuracy in the training sample with the following restrictions: (1) only the top 20 risk factors were considered based on variable importance derived from the random survival forest algorithm, (2) the score was a linear combination of risk factors, (3) coefficients were limited to integers between −10 and 10 and risk factors with a coefficient of 0 were dropped from the risk score, (4) risk factors with high collinearity with other risk factors included in the simplified risk score were further excluded (i.e. risk factors with a Pearson correlation >0.3; among collinear risk factors, those with the highest predictive accuracy were kept). For each patient month, observations were classified into risk groups based on the risk score terciles (i.e. each group comprised a third of the training sample). Given that risk factors could vary over time, a given patient could be classified into different risk groups at different time points. Kaplan-Meier analyses of time to unsatisfactory response over the duration of the at-risk period were conducted by the risk group.

The analyses in the present study were performed using R version 3.6.2; the random survival forest was developed using the Ranger package version 0.12 and the simplified risk score was developed using the CEoptim package version 1.2.

Results

A total of 17,026 patients with documented use of a PAH-specific PDE5i or ERA at any time were identified in the claims database. After applying the eligibility criteria, a total of 4781 adult patients with PAH were included in the study sample (Figure 2).

Figure 2. Sample selection flowchart. Abbreviations. CTEPH, chronic thromboembolic pulmonary hypertension; ERAs, endothelin-receptor antagonists; PAH, pulmonary arterial hypertension; PDE5is, phosphodiesterase type 5 inhibitors; PPA, prostacyclin pathway agent. Notes: 1. The index regimen was defined as the treatment initiated on the index date, along with any other treatment indicated for PAH initiated within a 60-day window after the index date.

Patient characteristics

Among patients in the study sample, the median duration of the at-risk period was 14.0 months. During the at-risk period, 66.3% of patients (n = 3169) experienced an unsatisfactory response (Table 1); the most common unsatisfactory response among these patients was ≥ 1 PAH-related or heart failure-related hospitalization or ER visit (54.6%; n = 2611).

Patient demographic and clinical characteristics during the baseline period are presented in Table 2. In the overall study sample, the median age was 69.0 years and 58.6% of patients were female. Patients in the overall study sample had a median Quan-CCI score of 4.0 and the most common comorbidities included valvular disease (60.4%), chronic obstructive pulmonary disorder (54.3%), and cardiomegaly (43.8%).

Table 2. Baseline characteristics.

	Study sample	Patients without outcome	Patients with outcome
Patient characteristics^a	N = 4781	N = 1612	N = 3169	P-value
Demographics
Age, mean ± SD [median]	67.8 ± 12.7 [69.0]	67.2 ± 12.8 [68.0]	68.1 ± 12.6 [70.0]	<.05
Age group, n (%)				<.05
18–29 years old	53 (1.1%)	20 (1.2%)	33 (1.0%)
30–39 years old	93 (1.9%)	30 (1.9%)	63 (2.0%)
40–49 years old	261 (5.5%)	102 (6.3%)	159 (5.0%)
50–59 years old	718 (15.0%)	247 (15.3%)	471 (14.9%)
60–69 years old	1275 (26.7%)	460 (28.5%)	815 (25.7%)
70–79 years old	1464 (30.6%)	468 (29.0%)	996 (31.4%)
80–89 years old	917 (19.2%)	285 (17.7%)	632 (19.9%)
Female	2800 (58.6%)	777 (48.2%)	2023 (63.8%)	<.001
Insurance type, n (%)				<.001
Exclusive provider organization	190 (4.0%)	77 (4.8%)	113 (3.6%)
Health maintenance organization	1100 (23.0%)	242 (15.0%)	858 (27.1%)
Indemnity	156 (3.3%)	39 (2.4%)	117 (3.7%)
Point-of-service	1037 (21.7%)	435 (27.0%)	602 (19.0%)
Preferred provider organization	318 (6.7%)	93 (5.8%)	225 (7.1%)
Other	1980 (41.4%)	726 (45.0%)	1254 (39.6%)
Geographic region, n (%)				<.001
Midwest	949 (19.9%)	272 (16.9%)	677 (21.4%)
Northeast	454 (9.5%)	169 (10.5%)	285 (9.0%)
South	2162 (45.3%)	775 (48.1%)	1387 (43.8%)
West	1210 (25.3%)	394 (24.5%)	816 (25.8%)
Time between first observed PAH diagnosis and the index date (months), mean ± SD [median]	5.2 ± 3.9 [4.2]	5.8 ± 3.9 [5.3]	4.9 ± 3.9 [3.7]	<.001
Treatment regimen^e as of the index date, n (%)
Monotherapy	4318 (90.3%)	1497 (92.9%)	2821 (89.0%)	<.001
PDE5i	3483 (72.9%)	1309 (81.2%)	2174 (68.6%)	<.001
Sildenafil	2656 (55.6%)	1012 (62.8%)	1644 (51.9%)	<.001
Tadalafil	827 (17.3%)	297 (18.4%)	530 (16.7%)	.153
ERA	835 (17.5%)	188 (11.7%)	647 (20.4%)	<.001
Ambrisentan	374 (7.8%)	88 (5.5%)	286 (9.0%)	<.001
Bosentan	274 (5.7%)	38 (2.4%)	236 (7.4%)	<.001
Macitentan	187 (3.9%)	62 (3.8%)	125 (3.9%)	.931
Combination therapy	463 (9.7%)	115 (7.1%)	348 (11.0%)	<.001
Quan-Charlson comorbidity index^b, mean ± SD [median]	3.9 ± 2.2 [4.0]	3.6 ± 2.2 [3.0]	4.1 ± 2.2 [4.0]	<.001
Medical encounters, mean ± SD [median]
All-cause hospital admissions	1.5 ± 2.0 [1.0]	1.2 ± 1.7 [1.0]	1.6 ± 2.0 [1.0]	<.001
PAH-related outpatient visits^c	5.6 ± 7.8 [4.0]	5.0 ± 8.2 [3.0]	5.9 ± 7.5 [4.0]	<.001
Cardiologist visits	9.7 ± 13.6 [6.0]	8.8 ± 13.6 [5.0]	10.1 ± 13.6 [6.0]	<.01
Pulmonologist visits	4.8 ± 7.1 [3.0]	3.7 ± 5.6 [2.0]	5.3 ± 7.7 [3.0]	<.001
PAH-related procedures, n (%)
Echocardiography	2111 (44.2%)	585 (36.3%)	1526 (48.2%)	<.001
RHC	1258 (26.3%)	309 (19.2%)	949 (29.9%)	<.001
Comorbidities, n (%)
Cardiomegaly	2092 (43.8%)	632 (39.2%)	1460 (46.1%)	<.001
Chronic obstructive pulmonary disease	2597 (54.3%)	810 (50.2%)	1787 (56.4%)	<.001
Connective tissue disorder	640 (13.4%)	184 (11.4%)	456 (14.4%)	<.01
Diabetes mellitus	2128 (44.5%)	675 (41.9%)	1453 (45.9%)	<.01
Pulmonary fibrosis	976 (20.4%)	241 (15.0%)	735 (23.2%)	<.001
Renal disease	1837 (38.4%)	560 (34.7%)	1277 (40.3%)	<.001
Thyroid disorders	1350 (28.2%)	412 (25.6%)	938 (29.6%)	<.01
Valvular disease	2886 (60.4%)	960 (59.6%)	1926 (60.8%)	.432
Symptoms, n (%)
Dyspnea	4241 (88.7%)	1349 (83.7%)	2892 (91.3%)	<.001
Chest pain	1942 (40.6%)	535 (33.2%)	1407 (44.4%)	<.001
Lower limb edema	1684 (35.2%)	467 (29.0%)	1217 (38.4%)	<.001

Abbreviations. AIDS, Acquired immunodeficiency syndrome; ER, emergency room; HIV, Human immunodeficiency virus; ERAs, endothelin-receptor antagonists; PAH, pulmonary arterial hypertension; PDE5is, phosphodiesterase type 5 inhibitors; PH, pulmonary hypertension; RHC, right heart catheterization.

^aUnless specified otherwise, characteristics were measured during the baseline period (i.e. 12-month period prior to the index date).

^bQuan et al. [41].

^cPAH-related is defined as having at least one diagnosis code for PAH (ICD-9-CM: 416.0, 416.8; ICD-10-CM: I27.0, I27.20, I27.21, I27.89) or right heart failure-related (ICD-9-CM: 428.0; ICD-10-CM: I50.81) recorded on a medical claim in an inpatient or ER setting.

Patients with an unsatisfactory response during the at-risk period had different demographic and clinical characteristics at baseline compared to those who did not have the study outcome (Table 2). In particular, patients with an unsatisfactory response were more likely to be older (median age: 70.0 vs. 68.0; p < .05) and female (63.8% vs. 48.2%; p < .001). During the baseline period, patients with an unsatisfactory response also had a greater comorbidity burden, as evidenced by a higher median Quan-CCI score (4.0 vs. 3.0; p < .001) and a greater proportion of patients with cardiomegaly (46.1% vs. 39.2%; p < .001), chronic obstructive pulmonary disorder (56.4% vs. 50.2%; p < .001), and pulmonary fibrosis (23.2% vs. 15.0%; p < .001). Patients with an unsatisfactory response also had a greater frequency of common PAH-related symptoms such as dyspnea (91.3% vs. 83.7%; p < .001), chest pain (44.4% vs. 33.2%; p < .001), and lower limb edema (38.4% vs. 29.0%; p < .001). Moreover, an unsatisfactory response to PAH therapy was associated with more frequent medical encounters. Compared to patients without the study outcome, patients with an unsatisfactory response had a greater mean number of all-cause hospital admissions (1.6 vs. 1.2; p < .001), PAH-related outpatient visits (5.9 vs. 5.0; p < .001), cardiologist visits (10.1 vs. 8.8; p < .001), and pulmonologist visits (5.3 vs. 3.7; p < .001), as well as more PAH-related procedures such as echocardiography (48.2% vs. 36.3%) and right heart catheterization (RHC; 29.9% vs. 19.2%; all p < .001).

Random survival forest algorithm

Predominant risk factors for unsatisfactory response identified by the algorithm are presented in Figure 3. Risk factors with the greatest importance during the prior 12-month period included HRU (i.e. PAH-related outpatient visits, pulmonologist visits, cardiologist visits, and all-cause hospital admissions), time since first PAH diagnosis, time since index date, Quan-CCI score, dyspnea, and age. Notable predictors of unsatisfactory response also included PAH-specific procedures (i.e. echocardiography and RHC), gender, comorbidities (i.e. cardiomegaly, pulmonary fibrosis, valvular disease, and chronic obstructive pulmonary disease), and other PAH-specific symptoms (i.e. chest pain, lower limb edema). The algorithm demonstrated a good ability to predict unsatisfactory response to initial PAH therapy (C-statistic: 0.732); the algorithm had 82% sensitivity and 55% specificity at 12 months (Table 3).

Figure 3. Top 20 most important risk factors based on random survival model algorithm. Abbreviations. PAH, pulmonary arterial hypertension; RHC, right heart catheterization.

Table 3. Algorithm performance at key time points.

Concordance of the model is 0.732
Key time point	Sensitivity	Specificity	PPV	NPV
3 months	0.41	0.89	0.45	0.88
6 months	0.65	0.75	0.51	0.84
9 months	0.75	0.64	0.57	0.80
12 months	0.82	0.55	0.62	0.78
15 months	0.87	0.48	0.67	0.77
18 months	0.90	0.43	0.71	0.75
21 months	0.93	0.40	0.74	0.74
24 months	0.94	0.37	0.77	0.73
3 years	0.96	0.29	0.86	0.62
4 years	0.97	0.23	0.92	0.47
5 years	0.98	0.16	0.95	0.29

Abbreviations. NPV, negative predictive value; PPV, positive predictive value.

Simplified risk score

The simplified risk score was developed based on the 20 top risk factors presented in Figure 3. A final set of seven risk factors were incorporated into the simplified risk score, with scores ranging from −10 to 27 (higher score indicates a higher risk of the outcome; Table 4). Specifically, the factors in the simplified risk score included: having ≥4 PAH-related outpatient visits (+8 points), ≥1 pulmonologist visit (+6 points), ≥1 hospitalization for any cause (+5 points), or a diagnosis of dyspnea (+5 points), pulmonary fibrosis (+2 points), or chest pain (+1 point). Additionally, the simplified risk score included ≥3 months since treatment initiation, which was associated with a reduced risk of unsatisfactory response (−10 points). The simplified risk score achieved good predictive accuracy in the training set (C-statistic: 0.645) and the validation set (C-statistic: 0.668; Table 4).

Table 4. Simplified risk score.

Risk factors in the prior 12 months	Points
Healthcare resource utilization
≥4 PAH-related outpatient visits	8
≥1 pulmonologist visits	6
≥1 hospitalization for any cause	5
Diagnosis
Dyspnea	5
Pulmonary fibrosis	2
Chest pain	1
≥3 months since treatment initiation	−10
Risk score range^a	−10–27
C-statistic (training set)	0.645
C-statistic (validation set)	0.668

^aRisk score is calculated as the sum of the points for each risk factor for given patient.

Survival curves by risk group

Based on risk score terciles, each patient-month observation was assigned to one of three risk groups as follows: (1) low-risk group (score ≤ 3), (2) medium-risk group (4 ≤ score ≤ 11), and (3) high-risk group (score ≥ 12). Compared to patients in the low-risk group, those in the high-risk group were more likely to have ≥4 PAH-related outpatient visits in the past year (79.6% vs. 7.4%; p < .001), ≥1 pulmonologist visit (95.1% vs. 65.4%; p < .001), and ≥1 hospitalization for any cause (71.6% vs. 9.3%; p < .001). They were also more likely to have dyspnea (96.2% vs. 44.8%; p < .001), pulmonary fibrosis (26.8% vs. 10.0%; p < .001), and chest pain (46.2% vs. 9.0%; p < .001; data not shown).

The Kaplan-Meier analysis shows that the three risk groups have a clear separation in terms of their probability of not having an unsatisfactory response at a given time point (Figure 4). In particular, the median time to the outcome was 0.7 years (8.7 months) for the high-risk group, 1.3 years (15.1 months) for the medium-risk group, and 3.1 years (37.5 months) for the low-risk group.

Figure 4. Kaplan-Meier curves for unsatisfactory response by risk groups.

Discussion

The present study used machine learning to develop and validate a predictive algorithm for unsatisfactory response to initial PAH therapy using health insurance claims data. An algorithm was developed to capture the impact of risk factors evaluated in the prior 12 months on unsatisfactory response during the subsequent at-risk period. Risk factors with the greatest importance included healthcare resource use (e.g. PAH-related outpatient visits, all-cause hospitalizations), time since first PAH diagnosis, time since index date, Quan-CCI score, dyspnea, and age. Based on the risk factors with a higher importance in the full algorithm, a simplified risk score was developed for ease of clinical application. The simplified risk score consisted of seven factors including HRU, diagnosis of symptoms/comorbidities, and time since treatment initiation. Despite the simplified risk score having a lower predictive accuracy than the full algorithm (C-statistic: 0.732 vs. 0.668), it provided a clear separation of risk among the three different risk categories (i.e. low, medium, and high risk). In light of these findings, the current risk assessments could potentially serve as part of a comprehensive strategy for the management of patients with PAH, alongside the use of clinical judgment and evidence-based guidelines. We caution that this model should not be used by payers to limit access to care for patients based on low-risk scores. Rather, it could serve as a complementary tool to identify patients at high risk that are not already flagged by physicians. In any case, a physician’s judgement should always prevail when determining the appropriate management of a given patient.

In the present study, increased HRU in the last 12 months was strongly associated with a higher risk of unsatisfactory response. Consistent with this, all-cause hospitalization was associated with an increased risk of 1-year mortality based on the REVEAL 2.0 and REVEAL Lite 2 risk scores^20,25. Similarly, another analysis of data from the REVEAL registry found that all-cause hospitalization was a strong marker for poor 1-year survival in both previously diagnosed and newly diagnosed PAH patients³³. Notably, in an analysis of data from the Phase III GRIPHON trial³⁴, the primary endpoint component of hospitalization for worsening of PAH was found to be a significant risk factor for subsequent mortality over a 1-year follow-up³⁵. Treatment guidelines recommend that healthcare providers perform a comprehensive assessment of patients with PAH at diagnosis/treatment initiation, at least every 6–12 months routinely in the event of clinical worsening, and 3–6 months after a change in therapy^1,14. Overall, evidence suggests that multi-parametric assessments of patients are being performed by physicians in the context of routine consultations²⁷. Although regular PAH-related outpatient visits are therefore to be expected as part of patient management, the present study found that ≥4 visits in the prior year is highly predictive of an increased risk of unsatisfactory response to treatment. Unfortunately, the information available in the database does not allow us to determine the reasons why certain patients might have had ≥4 PAH-related outpatient visits in the prior year. However, there are plausible explanations based on the prior literature. First, the common problem of delayed diagnosis could have resulted in the diagnosis of patients at a more advanced stage of disease³⁶, which would, in turn, necessitate more frequent PAH-related assessments. Second, patients diagnosed with PAH could have had an event of clinical worsening that was not captured in the current database, but which required more frequent assessments^1,14,27. Finally, suboptimal therapeutic management of PAH could have led to poor symptom management and thus an increased number of outpatient visits.

Diagnoses of certain PAH-related symptoms and comorbidities were also important predictors of unsatisfactory response in the present study. In particular, recorded diagnoses of dyspnea or chest pain in the last year emerged as salient risk factors in both the full algorithm and the simplified risk score. Since PAH often goes undetected until the later stages due to the non-specificity of its symptoms, it is possible that patients seeking medical care for such symptoms had greater disease severity than those who did not. Prior evidence suggests that a rapid rate of symptomatic progression is associated with a worse prognosis among patients with PAH^1,37. Therefore, it is critical to identify patients at higher risk of an unsatisfactory response at an early stage, when they have less severe disease and may respond more favorably to the appropriate therapy³⁷. Of the comorbidities observed in the present study, pulmonary fibrosis was associated with an increased risk of unsatisfactory response, which corroborates previous reports of an association between pulmonary diseases and poor clinical outcomes in patients with PAH^1,38,39.

The predictive accuracy of the present random survival forest algorithm is comparable to other risk scores (e.g. C-statistics of 0.76 for REVEAL 2.0 and 0.73 for REVEAL Lite 2) ²⁰, while the simplified algorithm had slightly lower predictive accuracy than other risk scores. However, it was not possible to directly compare the performance of these various risk tools in our sample due to the lack of availability of risk factors such as the FC, 6MWD, and NT-proBNP/BNP levels in our claims database. Future studies are warranted to compare the performance of our claims-based algorithms with that of other risk scores in the same population. Given machine learning algorithms allow for non-linearity and interactions between risk factors, they may have better performance than regression-based risk assessments in certain circumstances, such as when the regression model’s assumptions (e.g. linearity, proportionality) are not met or when there are multiple interactions between variables²⁴. However, variations in the predictive accuracy between risk assessments may also reflect additional factors, such as differences with respect to outcomes, study populations, and the availability of certain patient-level information. For instance, previously validated risk assessments mainly predicted survival^20,25,26, whereas the present study aims to predict unsatisfactory response to treatment. While prior risk assessments were designed to predict outcomes among incident patients^22,23, treatment naïve patients²¹, or any patients with PAH^20,24–26, the present study aimed to predict outcomes specifically among patients with PAH initiated on the first treatment with PAH-specific PDE5is or ERAs. Finally, many risk factors commonly included in other risk assessments (e.g. FC, 6MWD, and NT-proBNP/BNP levels,) are not available in claims data, as they require physicians’ assessments and laboratory tests. An advantage of the risk assessments developed in this study is that they rely solely on information available in claims data, which could facilitate their application on a large scale by payers or quality-of-care organizations.

Potential applications

The full algorithm and simplified risk score developed in this study could have interesting applications in clinical practice, especially in light of the potential benefits associated with earlier combination therapy among higher-risk patients. In particular, the 2015 ESC/ERS guidelines and the 6th WSPH proceedings recommend initial or upfront PPA-based combination therapy for treatment-naïve patients deemed at high risk of 1-year mortality^1,19. For patients with an intermediate risk status who are already being treated with a combination of PAH-specific PDE5i and ERA, these guidelines advise escalation to triple therapy by the addition of an oral or parenteral PPA^1,19.

The full algorithm and simplified risk score could be adapted to the specific needs and capabilities of different end-users in the healthcare system. For instance, the full algorithm could be used by organizations on a large scale to track temporal changes in risk levels in the patient population. A case in point would be the integration of the full algorithm into an automated warning system, which healthcare payers could use to detect when a patient’s predicted time to unsatisfactory response falls below a certain value. Alternatively, the simplified risk score could enable clinicians to readily discern who may be a plausible candidate for further monitoring and/or treatment escalation. In particular, patients deemed at higher risk could benefit from more frequent monitoring of clinical outcomes that are not routinely available, such as FC, 6MWD, NT-proBNP/BNP levels, and measures of hemodynamics and pulmonary function, which are in turn useful for determining the severity of PAH¹. Table 5 shows examples of risk scores for hypothetical patients. For example, a patient who has been treated for 5 months with a limited number of PAH-related outpatient visits and no hospitalizations in the past 12 months would be deemed at low risk of an unsatisfactory response, despite having dyspnea (Case 1). In contrast, a patient who has been treated for only one month and has required frequent PAH-related outpatient care and hospitalizations in the past 12 months would be identified as having a high risk of unsatisfactory response (Case 3).

Table 5. Predicted risk of an unsatisfactory response for hypothetical patient profiles.

Risk factors in the prior 12 months	Case 1	Score	Case 2	Score	Case 3	Score
Healthcare resource utilization
PAH-related outpatient visits	2 visits	0	3 visits	0	6 visits	8
Pulmonologist visits	0 visit	0	1 visit	6	2 visits	6
Hospitalization for any cause	0 visit	0	1 visit	5	2 visits	5
Diagnosis
Dyspnea	Yes	5	Yes	5	No	0
Pulmonary fibrosis	No	0	Yes	2	No	0
Chest pain	No	0	Yes	1	No	0
Time since treatment initiation	5 months	−10	12 months	−10	1 month	0
Total score	−5		9		19
Risk group	Low		Medium		High

Strengths and limitations

The present study has several methodological strengths. Notably, the use of machine learning, and specifically random survival forests, provides a number of advantages over regression-based approaches for risk prediction. First, the random survival forest algorithm is independent of specific functional form assumptions³⁰. Second, in the case of high dimensional data, limitations of regression approaches such as overfitting, inflated standard errors, or convergence problems are less likely to arise when using random survival forests²⁹. Third, random survival forests are naturally robust to outliers because predictors enter the model only ordinally; for a given observation, once the observation has the highest value in the sample of a particular predictor, further increases in the value will have no effect on the model prediction⁴⁰.

The present study should also be viewed in the context of certain limitations. The study outcome and risk factors were identified based on information derived from health insurance claims data, which may have been subject to misclassification errors or omissions. For instance, an ICD-9/ICD-10 diagnosis code on a medical claim is not an attestation that the patient had the diagnosis, as it may have been recorded incorrectly or represent a ruled-out diagnosis. Moreover, information regarding certain potentially important predictors of unsatisfactory response (e.g. FC, 6MWD, NT-proBNP/BNP plasma level) are not available in claims data, which may have limited the predictive accuracy of the algorithm developed. In addition, although the algorithm is predictive of unsatisfactory response, many of the risk factors identified suggest that patients with more severe comorbidity profiles are at greater risk; future studies focused on patients presenting certain comorbidities would be helpful in identifying any additional risk factors within specific subpopulations. Another limitation of this study was that the index regimen was defined based on PAH treatment initiated within a 60-day window after the index date. While this specific time window was selected based on the distribution of the data from our study sample (with most patients receiving a new PAH-specific regimen within 60 days), it is possible that early unsatisfactory response to certain medications might not have been fully captured. Furthermore, the algorithm predictions at longer duration post-index should be interpreted with caution given that a limited number of patients had an at-risk period of such duration. Finally, the present study population was commercially insured and, therefore, it may not be representative of the overall population of patients with PAH; external validation in different databases is warranted to confirm the generalizability of the current risk score to other populations. Finally, findings from this study may be interpreted as measures of association; however, no causal inference can be drawn from this analysis.

Conclusions

The present study used machine learning to develop and validate a health insurance claims-based algorithm for predicting the time to unsatisfactory response among patients following initial PAH therapy with a PDE5i or ERA. The algorithm reliably discriminated between patients based on their risk of unsatisfactory response following the initiation of PAH-specific therapy. In clinical practice, healthcare providers could utilize a simplified risk score derived from this algorithm to identify higher-risk patients with PAH who may require further confirmation with other risk score assessments or other monitoring tools. The timely identification of these patients could enable early interventions, particularly with combination therapy, which may help to mitigate the risk of an unsatisfactory response, reduce hospitalizations, and improve prognosis.

Transparency

Declaration of funding

This work was supported by Janssen Scientific Affairs, LLC. (Grant number: N/A).

Declaration of financial/other relationships

Janssen Scientific Affairs, LLC sponsored this study. Yuen Tsang, Peter Agron, Karimah S. Bell Lynum, and Sumeet Panjabi are employees of Janssen Scientific Affairs, LLC and may own stocks. Marjolaine Gauthier-Loiselle, Patrick Lefebvre, and Jimmy Royer are employees of Analysis Group, Inc., and Lucas Bennett is a former employee of Analysis Group, Inc., which has received consultancy fees from Janssen Scientific Affairs, LLC for the conduct of this study. The study sponsor was involved in all aspects of the research, including the collection of data, its analysis and interpretation, and approval of the final manuscript for publication. Peer reviewers on this manuscript have no relevant financial or other relationships to disclose.

Author contributions

All authors participated in the design of the study, data interpretation, and critically revised the intellectual content of this manuscript. Marjolaine Gauthier-Loiselle, Patrick Lefebvre, Jimmy Royer, and Lucas Bennett participated in the data analysis. All authors approved the final manuscript.

Acknowledgements

The authors would like to thank Emma Russon, who is a former employee of Analysis Group, Inc., for her contribution to the development of the algorithm. The authors would also like to thank Tom Cornwall, who is currently employed at Analysis Group, Inc., for his contributions to data analysis and interpretation. Medical writing assistance was provided by Mona Lisa Chanda, PhD, an employee of Analysis Group, Inc. Funding for this assistance was provided by Janssen Scientific Affairs, LLC.

Data availability statement

The data that support the findings of this study are available from Optum, but restrictions apply to the availability of these data, which were used under license for the current study, and so are not publicly available. Any researchers interested in obtaining the data used in this study can access database through Optum, under a license agreement, including the payment of appropriate license fee.

Notes

i Optum Clinformatics is a registered trademark of OptumInsight Life Sciences, Inc., Eden Prairie, MN, USA.

ii While interaction terms can be added to a Cox regression model, the difficulty arises from high dimensionality. Accommodating all two-way interactions requires k² parameters, while all three way-interactions requires approximately k³, and so on. The number of parameters very quickly exceeds the number of observations, such that the model cannot be properly estimated.

iii For the development of such algorithms, 500 trees are often used due to the marginal improvement in performance observed with an increasing number of trees and for computational burden considerations. Based on the literature, 128 trees are presumed to be sufficient (e.g. Oshiro, Thais Mayumi, Pedro Santoro Perez, and José Augusto Baranauskas. "How many trees in a random forest?" International workshop on machine learning and data mining in pattern recognition. Springer, Berlin, Heidelberg, 2012).