Longitudinal respiratory subphenotypes in patients with COVID-19-related acute respiratory distress syndrome: results from three observational cohorts.

BACKGROUND: Patients with COVID-19-related acute respiratory distress syndrome (ARDS) have been postulated to present with distinct respiratory subphenotypes. However, most phenotyping schema have been limited by sample size, disregard for temporal dynamics, and insufficient validation. We aimed to identify respiratory subphenotypes of COVID-19-related ARDS using unbiased data-driven approaches.
METHODS: PRoVENT-COVID was an investigator-initiated, national, multicentre, prospective, observational cohort study at 22 intensive care units (ICUs) in the Netherlands. Consecutive patients who had received invasive mechanical ventilation for COVID-19 (aged 18 years or older) served as the derivation cohort, and similar patients from two ICUs in the USA served as the replication cohorts. COVID-19 was confirmed by positive RT-PCR. We used latent class analysis to identify subphenotypes using clinically available respiratory data cross-sectionally at baseline, and longitudinally using 8-hourly data from the first 4 days of invasive ventilation. We used group-based trajectory modelling to evaluate trajectories of individual variables and to facilitate potential clinical translation. The PRoVENT-COVID study is registered with ClinicalTrials.gov, NCT04346342.
FINDINGS: Between March 1, 2020, and May 15, 2020, 1007 patients were admitted to participating ICUs in the Netherlands, and included in the derivation cohort. Data for 288 patients were included in replication cohort 1 and 326 in replication cohort 2. Cross-sectional latent class analysis did not identify any underlying subphenotypes. Longitudinal latent class analysis identified two distinct subphenotypes. Subphenotype 2 was characterised by higher mechanical power, minute ventilation, and ventilatory ratio over the first 4 days of invasive mechanical ventilation than subphenotype 1, but PaO2/FiO2, pH, and compliance of the respiratory system did not differ between the two subphenotypes. 185 (28%) of 671 patients with subphenotype 1 and 109 (32%) of 336 patients with subphenotype 2 had died at day 28 (p=0·10). However, patients with subphenotype 2 had fewer ventilator-free days at day 28 (median 0, IQR 0-15 vs 5, 0-17; p=0·016) and more frequent venous thrombotic events (109 [32%] of 336 patients vs 176 [26%] of 671 patients; p=0·048) compared with subphenotype 1. Group-based trajectory modelling revealed trajectories of ventilatory ratio and mechanical power with similar dynamics to those observed in latent class analysis-derived trajectory subphenotypes. The two trajectories were: a stable value for ventilatory ratio or mechanical power over the first 4 days of invasive mechanical ventilation (trajectory A) or an upward trajectory (trajectory B). However, upward trajectories were better independent prognosticators for 28-day mortality (OR 1·64, 95% CI 1·17-2·29 for ventilatory ratio; 1·82, 1·24-2·66 for mechanical power). The association between upward ventilatory ratio trajectories (trajectory B) and 28-day mortality was confirmed in the replication cohorts (OR 4·65, 95% CI 1·87-11·6 for ventilatory ratio in replication cohort 1; 1·89, 1·05-3·37 for ventilatory ratio in replication cohort 2).
INTERPRETATION: At baseline, COVID-19-related ARDS has no consistent respiratory subphenotype. Patients diverged from a fairly homogenous to a more heterogeneous population, with trajectories of ventilatory ratio and mechanical power being the most discriminatory. Modelling these parameters alone provided prognostic value for duration of mechanical ventilation and mortality. FUNDING: Amsterdam UMC.

Introduction

SARS-CoV-2 has infected 229 million individuals worldwide and caused more than 4·7 million deaths as of the end of September, 2021. Infection with SARS-CoV-2, referred to as COVID-19, frequently results in acute respiratory failure that might require intensive care unit (ICU) admission for respiratory support, including invasive ventilation, with high mortality.2, 3, 4

Since acute respiratory failure in COVID-19 is caused by a single pathogen, it can be postulated that critically ill patients with COVID-19 form a homogenous group with a single phenotype. However, multiple studies have identified various clinical subphenotypes of COVID-19.5, 6, 7, 8, 9, 10, 11, 12 Most studies have focused on variation in biochemical information, such as plasma creatinine, D-dimer, and IL-6. These studies typically used data obtained at a single timepoint and, to date, external validation has been rare.

Evidence before this study

We searched MEDLINE, Embase, CINAHL, and Web of Science on June 21, 2021, with the terms “phenotypes” AND (”coronavirus” OR ”COVID-19”), with no date or language restrictions. Studies including patients who were not receiving ventilation were excluded, as were those that reported on paediatric and single-centre populations. No studies identified by our search described subphenotypes of COVID-19-related acute respiratory distress syndrome (ARDS) in a multicentre setting. Several studies described subphenotypes of COVID-19 with varying severities, and most concluded that there are large differences between patients requiring different levels of respiratory support. One study reported on clinical subphenotypes in the specific population of critically ill patients with COVID-19 but was done in a single centre. In that study, a cross-sectional analysis was done at baseline using both laboratory data and respiratory variables. Two subphenotypes were identified: one with high D-dimer, IL-6, and creatinine and a second without apparent biochemical alterations. There are numerous letters and short reports speculating on time-related subphenotypes of COVID-19-related ARDS that did not use statistical methods to evaluate the existence of these subphenotypes and were limited to a sample size below 100 patients from a single centre. Most of these studies speculated that compliance of the respiratory system is an important time-varying variable.

Added value of this study

To our knowledge, this is the first study to systematically evaluate the existence of respiratory subphenotypes in patients with COVID-19-related ARDS. We found no statistical evidence for respiratory subphenotypes using static data at any time during the first 4 days of mechanical ventilation. Time-dependent analysis showed that two subphenotypes developed during the first days of mechanical ventilation. Patients with an upward trajectory of ventilatory ratio, a marker of dead space ventilation, had a higher risk of venous thrombotic events, more frequently developed acute kidney injury, required longer invasive mechanical ventilation, and had higher mortality. These findings were validated in two replication cohorts.

Implications of all the available evidence

We should reject the hypothesis that there are multiple static subphenotypes of COVID-19-related ARDS. Respiratory system compliance does not seem to be the primary driver of time-dependent heterogeneity. Rather, a subphenotype with increasing dead space ventilation is associated with important patient-related outcomes and might be targetable with therapeutic anticoagulation. The increasing mechanical power subphenotype might be indicative of inappropriate use of high positive end-expiratory pressure. Therefore, both subphenotypes could be studied as treatable traits in the future.

A classification of critically ill patients with COVID-19 based on respiratory system compliance was postulated to present a treatable trait and to be dynamic in nature, with patients progressing from normal compliance to decreased compliance,13, 14 although this trend was not observed in a meta-analysis. These classifications received criticism because they were based on retrospective inspection of small patient cohorts, mostly obtained in a single centre, without the use of statistical methods that are typically used to identify subclasses. Unbiased identification of respiratory subphenotypes of critically ill patients with COVID-19 infection requires analysis of large, multicentre cohorts with longitudinal data using appropriate methodological approaches.

In this study, using a derivation and two replication cohorts, we aimed to assess the existence of respiratory subphenotypes in critically ill patients with COVID-19. We hypothesised that using cross-sectional data, based on respiratory variables, no distinct subphenotypes of COVID-19 would be observed. Furthermore, we hypothesised that the richness of using longitudinal data for partitioning variables would be more informative in identifying respiratory phenotypes with distinct characteristics. Finally, we aimed to determine the extent to which trajectories of single variables can be used as prognostic indicators, because single variable trajectories are more clinically implementable than longitudinal multivariate subphenotyping approaches.

Methods

Derivation cohort

PRoVENT–COVID was an investigator-initiated, national, multicentre, prospective, observational cohort study at 22 ICUs in the Netherlands.4, 18 Consecutive patients aged 18 years or older were eligible for participation in the PRoVENT–COVID study if they were admitted to a participating ICU between March 1, 2020, and May 15, 2020, and had received invasive ventilation for respiratory failure related to COVID-19. COVID-19 was confirmed by positive RT-PCR. Patients who received non-invasive ventilation only, and patients who were transferred to a non-participating ICU within 1 h after intubation and start of invasive ventilation, were excluded from the current analysis.

A detailed description of study procedures in the derivation cohort has been reported elsewhere. Day 0 was defined as the first calendar day that a patient received invasive ventilation in a participating ICU, irrespective of hospital or ICU admission date. In the first hour of invasive ventilation and thereafter every 8 h at fixed timepoints up to day 3, ventilator settings and parameters and vital signs were entered into the database. ARDS was defined by the Berlin definition. Outcomes were duration of invasive mechanical ventilation, ICU mortality, hospital mortality, overall 28-day mortality, overall 90-day mortality, development of venous thrombotic events during ICU stay, acute kidney injury, need for renal replacement therapy during ICU stay, and length of ICU stay. All data were entered into a password-secured, internet-based, electronic case report form (Castor EDC; Amsterdam, Netherlands). Before analysis, the study coordinators screened all data for potentially erroneous or incomplete recordings and verified and corrected information as appropriate with the help of local doctors and data collectors. After cleaning, the database was closed for analysis.

The study protocol has been previously published. The institutional review boards of the participating centres approved the study protocol and the need for individual informed consent was waived. The ethics boards of the participating hospitals approved the collection of data for the study purposes. The study coordinators and trained data collectors assisted local doctors and monitored the study according to the International Conference on Harmonization Good Clinical Practice guidelines. Study coordinators ensured integrity and timely completion of data collection.

Replication cohorts

The two replication cohorts included patients 18 years or older who required invasive mechanical ventilation for COVID-19 pneumonia between March 1, 2020, and Dec 31, 2020, at Michigan Medicine (Ann Arbor, MI, USA) or at one of ten BJC HealthCare Hospitals (St Louis, MO, USA) in or around the St Louis metropolitan area. For both cohorts, patients transferred from outside hospitals were excluded if they were initiated on invasive ventilation before the date of transfer, or if they received less than 48 h of mechanical ventilation. COVID-19 was confirmed by positive RT-PCR. All data used in the replication cohorts were collected retrospectively and extracted from electronic health records. Analyses in the replication cohorts were done locally, such that no protected health information was exchanged between study sites.

Statistical analysis

Subphenotypes were studied cross-sectionally at the start of invasive ventilation in the intensive care unit and every 24 h after, and longitudinally using time-dependent analysis with 8-hourly data from the first 4 days of invasive ventilation. Subphenotypes were identified using latent class analysis and using group-based trajectory modelling. Figure 1A shows a visual representation of the differences in the data used in each model. We used group-based trajectory modelling to allow for external validation of our results as no other dataset with granular multivariate data could be identified, which is indicative for the fact that group-based trajectory modelling is more clinically implementable than longitudinal latent class analysis. A stepwise approach was used that included: data setup, model estimation, model evaluation, and interpretation of the optimal model.

Figure 1

Study statistical methods to identify subphenotypes

Summary of statistical analysis methods to model heterogeneity in respiratory variables. (A) The x-axis shows different timepoints. The y-axis shows different patients. The z-axis shows different variables. (B) Cross-sectional latent class analysis was done on each timepoint using multiple variables and did not yield any subphenotypes. (C) Longitudinal LCA was done on all timepoints combined using multiple variables and yielded two subphenotypes with differences in dynamics of mechanical power and ventilatory ratio. Subphenotype 2 had less ventilator-free days and more venous thrombotic events, but no difference in mortality. (D) GBTM was used to evaluate individual trajectories over all timepoints of a single variable. An upward trajectory of ventilatory ratio matched with the longitudinal LCA dynamics and was also associated with more venous thrombotic events, higher mortality, and fewer ventilator-free days. GBTM=group-based trajectory modelling. LCA=latent class analysis. MP=mechanical power. VR=ventilatory ratio. VTE=venous thrombotic events.

The following readily available clinical variables were used as input for the identification of latent classes: driving pressure, compliance of the respiratory system (tidal volume over driving pressure), minute volume ventilation, pH, ventilatory ratio (VEmeasured × PaCO2measured ÷ VEpredicted × PaCO2predicted)—where VE is the expired minute ventilation and PaCO2 is the arterial carbon dioxide tension—difference between arterial and end-tidal CO2, and mechanical power of ventilation (0·098 × VT [in L] × RR × [Pmax – Pdriving÷2])— where VT is the tidal volume, RR is the respiratory rate, Pmax is the maximum airway pressure on pressure-controlled ventilation, and Pdriving is the driving pressure. Variables with mathematical coupling were not considered (eg, not including respiratory rate, tidal volume, and minute volume ventilation).

Latent class analysis was done using data obtained directly after the start of mechanical ventilation to identify cross-sectional subphenotypes. This analysis was repeated for data collected every 24 h. Latent class analysis was repeated using 8-hourly data from the first 4 days of invasive ventilation to identify time-dependent subphenotypes that emerged during the first 96 h of ICU admission. The specific model parameters are given in the appendix (p 2).

The optimal number of latent classes was selected using the lowest Bayesian information criterion, integrated complete likelihood, and Akaike information criterion values. Entropy and probability of group membership were used to evaluate class separation and robustness of class membership. Entropy is a standardised index of model-based classification accuracy, with higher values indicating more precise assignment of individuals to latent profiles. Five models, comprising 1 to 5 classes, were fitted and if one or more models had similar goodness of fit, the model with the lowest number of classes was selected.

The prognostic value of the subphenotypes was evaluated by comparing differences in ventilator-free days, 28-day mortality, occurrence of venous thrombotic events, and use of renal replacement therapy for acute kidney injury. Venous thrombotic event diagnosis was based on venous Doppler ultrasound or CT angiography.

We used group-based trajectory monitoring to assess if the trajectory of a single variable could be used to identify trajectory subphenotypes with similar dynamics to those identified by time-dependent latent class analysis, since this method was straightforward enough to allow ease of analysis in the replication datasets. We applied group-based trajectory measurement on 8-hourly data from the first 4 days to identify the trajectory subphenotypes for PaO2/FiO2, ventilatory ratio, mechanical power, and respiratory system compliance. Group-based trajectory monitoring is a finite mixture model used to identify clusters of patients following similar trajectories of a variable of interest and is easier to apply to a smaller dataset than longitudinal latent class analysis. The group-based trajectory monitoring algorithm computes a unique equation of the variable of interest as a function of time for each of the subphenotypes. Patients are classified into the trajectory subphenotype whose function most closely matches their measurements. This approach has previously successfully been applied to temperature trajectories.

The fit of the group-based trajectory model from the derivation cohort was applied to the replication cohorts and the association with outcomes was evaluated.

To study the influence of treatment strategies on the prevalence of the identified subphenotype, we compared the following interventions before or at the day of initiation of invasive mechanical ventilation between the groups: high or low positive end-expiratory pressure strategy based on the selected positive end-expiratory pressure/FiO2 table as extensively described, prone positioning, and use of remdesivir, corticosteroids, and tociluzimab.

Patient data are presented as mean and SD for normally distributed continuous variables, median and IQR for variables that are not normally distributed, and numbers with percentages for categorical data. Differences between identified subphenotypes were tested using t-test, Mann-Whitney U test, and Fisher's exact test. The dynamics of each variable per subphenotype were fitted using a second-degree polynomial to show the trajectory over time between the groups. Differences between subphenotypes at baseline and differences in dynamic change over time were assessed by linear mixed effect model analysis using hours since intubation, subphenotype, and an interaction term between the two as fixed effects and a random intercept per patient. We used logistic and linear regression to assess the association between subphenotype classification and binary or continuous outcomes, respectively. Odds ratios (ORs) with 95% CIs for the assocation between subphenotype and 28-day mortality, adjusted for age, gender, and BMI were calculated. A sensitivity analysis was done on data collected from patients who remained on invasive mechanical ventilation for more than 96 h. Missing data were imputed for the latent class analysis using multivariate imputation by chained equations; the results from the first dataset are presented and the additonal datasets were evaluated for consistency.

Data analysis was done in R version 4.0.3 through the R studio interface, with the exception of group-based trajectory monitoring, which applied using the traj package command in Stata MP16.

The PRoVENT-COVID study is registered with is registered with ClinicalTrials.gov, NCT04346342.

Role of the funding source

The funder of the study had no role in study design, data collection, data analysis, data interpretation, or writing of the report.

Results

Between March 1, 2020, and May 15, 2020, 1122 patients were screened for study entry, of whom 115 were excluded because they previously received invasive ventilation, leaving 1007 patients included in the derivation cohort (appendix p 4). We recorded patient characteristics (table 1 ). Included patients were mostly men (814 [73%] of 1007 patients), with a mean age of 63·7 years (SD 10·8) and a mean BMI of 28·6 kg/m2 (SD 5·6). 294 (29%) of 1007 (29%) patients died by day 28 after ICU admission. 288 patients were included for analysis from replication cohort 1 and 326 from replication cohort 2 (appendix p 4). Patients in the replication cohorts had a higher mean BMI of 34·0 kg/m2 (SD 18·6) in replication cohort 1 and 32·6 kg/m2 (10·0) in replication cohort 2. At day 28, 78 (27%) patients in replication cohort 1 and 138 (42%) in replication cohort 2 had died.

Table 1

Patient characteristics at the start of mechanical ventilation (ie, day 0) and outcomes in the derivation cohort

			All (n=1007)	Time-dependent subphenotypes*		p value
				1 (n=671)	2 (n=336)
Age, years			63·7 (10·8)	63·9 (11·0)	63·2 (10·2)	0·35
Sex						<0·001
	Female		279 (28%)	226 (34%)	53 (16%)
	Male		728 (73%)	445 (66%)	283 (84%)
BMI, kg/m2			28·6 (5·6)	28·3 (4·5)	29·1 (7·4)	0·054
Comorbidities
	Arterial hypertension		342 (34%)	229 (34%)	113 (34%)	0·930
	Heart failure		42 (4%)	28 (4%)	14 (4%)	1·000
	Diabetes		225 (22%)	163 (24%)	62 (18%)	0·044
	Chronic kidney disease		45 (4%)	163 (24%)	62 (18%)	0·015
	Chronic obstructive pulmonary disease		83 (8%)	59 (9%)	23 (7%)	0·350
	Immunosuppression		24 (2%)	15 (2%)	9 (3%)	0·830
Day of admission†
	Tidal volume per kg predicted bodyweight, mL/kg		6·46 (1·55)	6·42 (1·28)	6·54 (2·00)	0·24
	Respiratory rate, breaths/min		21·00 (4·60)	20·95 (4·54)	21·96 (4·84)	0·001
	Minute ventilation, L/min		9·40 (2·70)	8·97 (2·40)	10·25 (3·12)	<0·001
	Positive end-expiratory pressure, cmH2O		12·6 (2·9)	12·4 (2·8)	13·2 (2·9)	<0·001
	High positive end-expiratory pressure strategy		269 (27%)	167 (25%)	102 (30%)	0·130
	Plateau pressure, cmH2O		27·3 (5·2)	26·9 (5·1)	28·0 (5·3)	0·003
	Driving pressure, cmH2O		14·9 (4·4)	14·8 (4·3)	15·2 (4·5)	0·200
	Compliance respiratory system, mL/cmH2O		32·3 (12·0)	31·8 (11·5)	33·3 (13·0)	0·082
	PaO2, mm Hg		93 (43)	91 (43)	98 (43)	0·026
	PaO2/FiO2, mm Hg		148 (75)	143 (77)	153 (71)	0·182
		Mild	262 (26%)	162 (24%)	100 (30%)	0·110
		Moderate	520 (52%)	350 (52%)	170 (51%)	..
		Severe	225 (22%)	159 (24%)	66 (20%)	..
	Mechanical power, J/min		18·4 (6·7)	17·3 (6·1)	20·6 (7·3)	<0·001
	PaCO2, mm Hg		45·2 (11·8)	43·5 (10·7)	48·8 (13·1)	<0·001
	Ventilatory ratio		1·64 (0·68)	1·52 (0·49)	1·89 (0·91)	<0·001
	Difference between arterial and end-tidal CO2		7·72 (11·14)	6·71 (10·40)	9·79 (12·31)	<0·001
	pH‡		7·36 (0·10)	7·36 (0·09)	7·34 (0·10)	<0·001
	Creatinine, μmol/L‡		74 (61–96)	73 (58–93)	77 (65–97)	0·015
	Urine output, mL/day		705 (370–1148)	670 (350–1125)	725 (415–1202)	0·071
	Fluid balance, mL/day		608 (29–1418)	539 (13–1268)	712 (71–1587)	0·048
	Prone position during first ICU day		300 (30%)	199 (30%)	101 (30%)	0·520
	Vasopressor in first ICU day		792 (79%)	520 (77%)	271 (81%)	0·340
Outcomes
	Acute kidney injury§		447 (44%)	270 (40%)	177 (53%)	0·001
	Renal replacement§		180 (18%)	109 (16%)	70 (21%)	0·088
	Venous thrombotic event§		285 (28%)	176 (26%)	109 (32%)	0·048
	ICU length of stay, days		15 (9–27)	15 (9–25)	16 (9–30)	0·36
	Duration of invasive mechanical ventilation, days		13 (8–23)	13 (8–23)	14 (8–24)	0·41
	Ventilator-free days and alive at day 28, days		2 (0–16)	5 (0–17]	0 (0–15)	0·016
	Day 7 mortality		112 (11%)	72 (11%)	40 (12%)	0·83
	Day 28 mortality		294 (29%)	185 (28%)	109 (32%)	0·10
	Day 90 mortality		349 (35%)	223 (33%)	123 (37%)	0·24

Data are mean (SD), n (%), or median (IQR). Mild, moderate, and severe PaO2/FiO2 were graded according to Berlin definition. ICU=intensive care unit.

Latent class analysis using 8-hourly data from the first 4 days of invasive ventilation was used to identify time-dependent subphenotypes that emerged during the first 96 h of ICU admission.

Date of ICU admission was the same as the day of start of mechanical ventilation in this cohort.

In blood samples.

During ICU stay.

Latent class analysis using cross-sectional data at initiation of invasive mechanical ventilation did not yield any latent classes in the derivation cohort at baseline (figure 1B). A one-class model had the lowest Akaike information criterion, Bayesian information criterion, integrated complete likelihood, and entropy values (appendix p 5). These analyses were repeated for data obtained 24 h, 48 h, 72 h, and 96 h after intubation and yielded similar results (appendix p 5). Therefore, this analysis did not yield any evidence for respiratory subphenotypes using cross-sectional data.

Using the longitudinal data from the first 4 days of invasive ventilation (figure 1C), a two-class model best fit the derivation cohort (appendix p 6). Entropy was 82% and the probability of class membership was good (appendix p 7). Sensitivity analysis including only the 826 (82%) patients who remained on invasive mechanical ventilation for more than 96 h showed similar results (appendix p 8), with only 22 (3%) patients changing class membership. We recorded standardised mean differences between classes over time (figure 2 ). Subphenotype 2 (336 [33%] of 1007 patients) was characterised by increasing minute ventilation, mechanical power, and ventilatory ratio over the first 4 days of invasive mechanical ventilation (appendix p 9).

Figure 2

Standardised mean differences between the two longitudinal respiratory subphenotypes

Venous thromboembolism was more common in subphenotype 2 (109 [32%] of 336 patients) than in subphenotype 1 (176 [26%] of 671 patients; p=0·048). The median numbers of ventilator-free days and patients alive at day 28 were lower in subphenotype 2 than in subphenotype 1 (5, IQR 0–17 for subphenotype 1 and 0, 0–15 for subphenotype 2; p=0·016; table 1; appendix p 10). Day 7, day 28, and day 90 mortality individually were not significantly different in subphenotype 2 compared with subphenotype 1 (table 1).

The trajectories of PaO2/FiO2 and respiratory compliance were distinct from those observed in the time-dependent latent class analysis-derived subphenotypes and were not considered further. The group-based trajectory model of ventilatory ratio and mechanical power showed overlapping trajectories with the time-dependent latent class analysis-derived subphenotypes (Figure 1, Figure 3 ). Two trajectories were observed, as follows: a stable value for ventilatory ratio or mechanical power over the first 4 days of invasive mechanical ventilation (trajectory A overlapping with subphenotype 1) or an upward trajectory (trajectory B overlapping with subphenotype 2). 28-day mortality and the rate of venous thrombotic events were higher in patients categorised to ventilatory ratio and mechanical power trajectory B than in their trajectory A counterparts (table 2 ). Trajectory B was independently associated with 28-day mortality for ventilatory ratio (OR 1·64, 95% CI 1·17–2·29) and mechanical power (1·82, 1·24–2·66) after adjusting for sex, age, and BMI. Thus, group-based trajectory model-derived trajectories for ventilatory ratio and mechanical power had prognostic potential in the derivation cohort and were evaluated in the replication cohorts.

Figure 3

Comparison of dynamic changes of time dependent latent class analysis subphenotypes and trajectory analysis

At each 8-hourly timepoint the median and IQR is plotted. The line shows second-degree polynomial regression.

Table 2

Outcomes for trajectories of ventilatory ratio and mechanical power in derivation cohort

		All (n=1007)	Ventilatory ratio trajectories		p value	Mechanical power trajectories		p value
			A (n=777)	B (n=230)		A (n=695)	B (n=312)
Subphenotype 2 (%)		336 (33%)	190 (24%)	146 (63%)	<0·001	118 (17%)	218 (70%)	<0·001
Age, years		63·7 (10·8)	63·63 (10·94)	63·89 (10·10)	0·750	64·12 (10·53)	62·72 (11·18)	0·057
Male sex		728 (72%)	573 (74%)	155 (67%)	0·071	464 (66·8)	264 (84·6)	<0·001
Body-mass index, kg/m2		28·6 (5·6)	28·2 (4·3)	29·9 (8·7)	<0·001	28·0 (4·3)	29·9 (7·6)	<0·001
High positive end-expiratory pressure strategy		269 (27%)	211 (27%)	58 (25%)	0·510	158 (23%)	111 (36%)	<0·001
Outcomes
	Acute kidney injury	447 (44%)	316 (41%)	131 (57%)	<0·001	262 (38%)	185 (59%)	<0·001
	Renal replacement	180 (18%)	120 (15%)	59 (26%)	0·001	92 (13%)	87 (28%)	<0·001
	Venous thrombotic event	285 (28%)	195 (25%)	90 (39%)	<0·001	180 (26%)	105 (34%)	0·014
	ICU length of stay, days	15 (9–27)	15 (9–25)	19 (10–32)	0·007	15 (9–24)	17 (9–31)	0·063
	Duration of invasive mechanical ventilation, days	13 (8–23)	13 (8–22)	16 (8–26)	0·008	13 (8–22)	15 (8–25)	0·039
	Ventilator-free days and alive at day 28, days	2 (0–16)	6 (0–18)	0 (0–9)	<0·001	7 (0–18)	0 (0–9)	<0·001
	Day 7 mortality	112 (11%)	74 (10%)	38 (17%)	0·011	60 (9%)	52 (17%)	0·001
	Day 28 mortality	294 (29%)	210 (27%)	84 (37%)	0·018	173 (25%)	121 (39%)	<0·001
	Day 90 mortality	349 (35%)	242 (31%)	104 (45%)	<0·001	212 (31%)	134 (43%)	<0·001

Data are n (%) or mean (SD). Trajectory A is a stable trajectory and trajectory B is an increasing trajectory.

After adjusting for sex, age, and BMI, ventilatory ratio trajectory B was independently associated with 28-day mortality in replication cohort 1 (OR 4·65, 95% CI 1·87–11·6) and replication cohort 2 (1·89, 1·05–3·37). Mechanical power trajectory B was independently associated with mortality in replication cohort 1 (OR 2·98, 95% CI 1·51–5·87), but not in replication cohort 2 (0·96, 0·54–1·72). Patients with mechanical power trajectory B had a longer duration of mechanical ventilation than patients in mechanical power trajectory A in the derivation and replication cohorts (Table 2, Table 3 ).

Table 3

Characteristics and outcomes in the replication cohorts

		All	VR trajectories		p value between VR trajectories	MP trajectories		p value between MP trajectories
			A	B		A	B
Replication cohort 1
Total		288	246	29	..	214	72	..
Age, years		60·2 (14·5)	60·6 (14·6)	56·2 (13·5)	0·13	61·5 (14·2)	56·2 (14·7)	0·007
Male sex		186 (65%)	161 (65%)	16 (55%)	0·38	126 (59%)	59 (82%)	0·001
BMI, kg/m2		34·0 (18·6)	33·6 (19·5)	39·8 (11·3)	0·097	32·7 (18·1)	35·9 (9·6)	0·158
Pharmacological treatments
	Remdesivir	94 (33%)	87 (35%)	7 (24%)	0·22	78 (36%)	18 (25%)	0·075
	Corticosteroids	163 (57%)	143 (58%)	20 (69%)	0·26	119 (56%)	49 (68%)	0·063
	Tociluzimab	83 (29%)	70 (28%)	13 (45%)	0·069	62 (29%)	25 (35%)	0·360
Outcomes
	ICU length of stay, days	15 (9–27)	15 (9–26)	18 (8–45)	0·41	14 (8·5–25)	18 (11–31·5)	0·015
	Duration of invasive mechanical ventilation, days	10·5 (6–21)	10 (6–20)	14 (7–39)	0·070	9 (6–20)	14 (7–27)	0·014
	Ventilator-free days and alive at day 28, days	7 (0–20)	10 (0–21)	0 (0–0)	<0·001	11 (0–21)	0 (0–15)	0·002
	Day 28 mortality	78 (27%)	63 (26%)	13 (45%)	0·049	53 (25%)	25 (35%)	0·137
Replication cohort 2
Total		326	255	71	..	86	186	..
Age, years		64·2 (14·4)	64·8 (14·7)	62·0 (12·9)	0·15	67·1 (14·9)	64·6 (13·3)	0·160
Male sex		192 (59%)	161 (63%)	31 (44%)	0·005	43 (50%)	116 (62%)	0·073
BMI, kg/m2		32·6 (10·0)	31·4 (9·2)	36·7 (11·7)	<0·001	28·4 (6·5)	32·4 (9·2)	<0·001
Pharmacological treatments
	Remdesivir	51 (16%)	33 (13%)	18 (25%)	0·011	6 (7%)	33 (18%)	0·018
	Corticosteroids	111 (34%)	77 (30%)	34 (48%)	0·005	20 (23%)	69 (37%)	0·023
	Tociluzimab	13 (4%)	7 (3%)	6 (8%)	0·002	1 (1%)	8 (4%)	0·180
Outcomes
	ICU length of stay, days	12 (6–22)	13 (6–24)	12 (7–20)	0·43	8 (3–17)	14 (7–24)	<0·001
	Duration of invasive mechanical ventilation, days	10 (5–16)	10 (5–16)	10 (7–15)	0·68	6·21 (3–11)	11 (7–17)	<0·001
	Ventilator-free days and alive at day 28, days	0 (0–16)	0 (0–16)	0 (0–11)	0·35	0 (0–20)	0 (0–14)	0·017
	Day 28 mortality	138 (42%)	102 (40%)	36 (51%)	0·134	37 (43%)	77 (41%)	0·900

Data are n, mean (SD), n (%), or median (IQR), unless otherwise indicated. Trajectory A was a stable trajectory and trajectory B was an increasing trajectory. BMI=body-mass index. ICU=intensive care unit. MP=mechanical power. VR=ventilatory ratio.

We found no difference in the use of a high positive end-expiratory pressure strategy between the time-dependent latent class analysis derived subphenotypes 1 and 2 (table 1), nor for the ventilatory ratio trajectories A and B (table 3). Patients with the upward trajectory (ie trajectory B) of mechanical power were more frequently exposed to a high positive end-expiratory pressure strategy (table 3). We found no difference in the use of remdesivir, corticosteroids, or tociluzimab between the A and B trajectories of mechanical power and ventilatory ratio in replication cohort 1, but there were differences in replication cohort 2 (table 3). However, use of these drugs did not explain the differences in prevalence of the trajectory subphenotypes between the cohorts.

Discussion

In these large observational cohort studies with granular clinical data of critically ill patients with acute respiratory failure due to severe COVID-19, we found no empirical evidence for the existence of respiratory subphenotypes at the start of invasive ventilation, nor at cross-sectional analysis in the succeeding 4 days. However, using time-dependent analysis, we identified two subphenotypes that developed during the first 4 days of invasive mechanical ventilation. Trajectories of ventilatory ratio and mechanical power were most discriminatory and modelling these parameters alone provided prognostic value for duration of mechanical ventilation and mortality.

Using respiratory and ventilatory data available directly after the start of invasive mechanical ventilation in patients with COVID-19-related ARDS, we found no suggestion for the existence of latent classes. In other words, the data were best explained by the presence of a single phenotype. This finding emphasises the importance of a data-driven approach to subphenotyping and contradicts earlier efforts that attempted to identify respiratory subphenotypes based on physiological reasoning and clinical observations.5, 6, 7, 8, 9 This conclusion does not necessarily mean that the physiological description of the cases that were reported were false, but rather that these fall at the extremes of a normal distribution and do not represent distinct subclasses. Insufficient evidence for latent classes based on respiratory and gas exchange variables alone is in line with a study in which classification of subphenotypes of COVID-19-related ARDS was driven by plasma biomarkers of organ failure, inflammation, and coagulation, whereas respiratory and gas exchange variables did not provide any discrimination.

Subphenotypes emerged in the derivation cohort using longitudinal data from the first 4 days of invasive ventilation. The importance of dynamic changes over time has been exemplified by the identification of temperature trajectory subphenotypes in patients with sepsis and COVID-19.23, 25 Patterns in temperature trajectories showed prognostic enrichment for mortality and were associated with differences in inflammatory markers in plasma. We used group-based trajectory modelling to identify subphenotypes based on the dynamic changes in PaO2/FiO2, compliance of the respiratory system, ventilatory ratio, and mechanical power. The major advantage of modelling one variable with group-based trajectory modelling compared with time-dependent latent class analysis is that the former is more suitable for application in a clinical setting and across multiple cohorts. Indeed, we were unable to identify any database with sufficiently granular data of all eight variables to externally validate the time-dependent latent class analysis approach. Furthermore, in our analyses, we found that both ventilatory ratio and mechanical power trajectories were predictive of duration of mechanical ventilation in two datasets, while only ventilatory ratio trajectories were indicative of an increased likelihood of 28-day mortality.

Trajectory analysis of the ventilatory ratio seems to a promising method for prognostication. The two trajectories had similar rates of baseline organ dysfunction but differentially developed complications such as venous thrombotic events and acute kidney injury requiring renal replacement therapy. Ventilatory ratio is easily calculated at the bedside when arterial blood gas analysis is available and is a good surrogate marker for dead space ventilation, although it can be influenced by other factors. The subphenotype with an upward ventilatory ratio trajectory showed a higher mortality independent of baseline risk factors. Patients with an upward ventilatory ratio trajectory also more frequently had venous thrombotic events, and pulmonary embolism is a likely contributor to ventilation-perfusion mismatch in this patient group, although we did not study this relationship directly in this study.27, 28 We speculate that patients with an upward ventilatory ratio trajectory show more pulmonary perfusion defects and might benefit from more intensive anticoagulatory treatment, even in the absence of pulmonary embolism. Perfusion defects in patients with COVID-19 can also be driven by immune-response related in-situ thrombosis, for which anticoagulation therapy might be less effective.29, 30, 31 However, immunomodulatory therapies such as dexamethasone and tociluzimab were not associated with a decreased prevalence of the upward trajectory subphenotype in this study, suggesting that these therapies do not eliminate the increase in ventilatory ratio over time.

Ventilatory ratio and mechanical power might be influenced by treatment strategies; therefore, we studied if variation in clinical practice could explain the existence of the identified subphenotypes. A difference in positive end-expiratory pressure strategy was observed between the mechanical power trajectory subphenotypes, with a higher positive end-expiratory pressure strategy more common in patients with an upward trajectory. Given that an increase in positive end-expiratory pressure without a decrease in driving pressure can result in an increase in mechanical power, this could be suggestive for inappropriate use of high levels of positive end-expiratory pressure in this subset of patients. This finding requires further evaluation in larger cohorts of critically ill patients with COVID-19-related ARDS. However, pharmacological treatment with remdesivir, corticosteroids, and tociluzimab did not consistently explain differences in the prevalence of the trajectory subphenotypes.

This study has several important strengths. Other studies have reassessed the prognostic accuracy of mechanical power, ventilatory ratio, and PaO2/FiO2 after 24–48 h26, 32, 33 but, to our knowledge, this study is the first to identify subphenotypes based on respiratory and gas exchange variables in a large cohort of patients, with external validation of the findings. We combined cross-sectional latent class analysis, time-dependent latent class analysis, and group-based trajectory modelling to move from a simplistic model of two cross-sectional subphenotypes to a classification based on the individual trajectories of ventilatory ratio and mechanical power during the first 96 h of mechanical ventilation. Use of a single variable is important for clinical applicability, which is shown by the fact that we were unable to identify a replication cohort with 8-hourly data for all variables using time-dependent latent class analysis. A major limitation of this study is the observational nature of the data, which prevented us from studying the potential for predictive enrichment. Furthermore, data were carefully curated in the derivation cohort, whereas there was an automated system to extract data from electronic health records in replication cohorts. Despite this difference in data curation, the relationship between trajectories and outcomes persisted. Furthermore, even though consecutive patients were included in all cohorts, the variation in mortality rates suggest that case mix varied. Both factors could explain the different incidence of the trajectory B of ventilatory ratio and mechanical power between cohorts.

In conclusion, patients with COVID-19-related ARDS consistently show a single respiratory phenotype at the start of invasive mechanical ventilation. The data suggest that there are at least two distinct trajectories during the first days of invasive mechanical ventilation, with one subphenotype showing increasing minute ventilation, mechanical power, and ventilatory ratio. Trajectories of mechanical power and ventilatory ratio were independently associated with outcome. COVID-19-related ARDS seems to diverge from a fairly homogenous respiratory physiology to a more heterogeneous population during the first 4 days of invasive mechanical ventilation. This finding reveals the importance of including time as a key variable in future efforts towards subphenotyping COVID-19. Anticoagulation and positive end-expiratory pressure selection should be considered as treatable traits in the identified subphenotypes.

Data sharing

Deidentified participant data with a data dictionary can be shared after approval of a proposal with a signed data access agreement and in collaboration with the study group.

This online publication has been corrected. The corrected version first appeared at thelancet.com/respiratory on November 30, 2021

Declaration of interests

LDJB reports grants from the Dutch Lung Foundation (Young Investigator grant), grants from the Dutch Lung Foundation and Health Holland (Public–Private Partnership grant), grants from the Dutch Lung Foundation (Dirkje Postma Award), grants from the IMI COVID19 iniative, and grants from Amsterdam UMC fellowship, outside the submitted work. All other authors declare no competing interests.