Bronchial Aspirate-Based Profiling Identifies MicroRNA Signatures Associated With COVID-19 and Fatal Disease in Critically Ill Patients

Background The pathophysiology of COVID-19-related critical illness is not completely understood. Here, we analyzed the microRNA (miRNA) profile of bronchial aspirate (BAS) samples from COVID-19 and non-COVID-19 patients admitted to the ICU to identify prognostic biomarkers of fatal outcomes and to define molecular pathways involved in the disease and adverse events. Methods Two patient populations were included (n = 89): (i) a study population composed of critically ill COVID-19 and non-COVID-19 patients; (ii) a prospective study cohort composed of COVID-19 survivors and non-survivors among patients assisted by invasive mechanical ventilation (IMV). BAS samples were obtained by bronchoaspiration during the ICU stay. The miRNA profile was analyzed using RT-qPCR. Detailed biomarker and bioinformatics analyses were performed. Results The deregulation in five miRNA ratios (miR-122-5p/miR-199a-5p, miR-125a-5p/miR-133a-3p, miR-155-5p/miR-486-5p, miR-214-3p/miR-222-3p, and miR-221-3p/miR-27a-3p) was observed when COVID-19 and non-COVID-19 patients were compared. In addition, five miRNA ratios segregated between ICU survivors and nonsurvivors (miR-1-3p/miR-124-3p, miR-125b-5p/miR-34a-5p, miR-126-3p/miR-16-5p, miR-199a-5p/miR-9-5p, and miR-221-3p/miR-491-5p). Through multivariable analysis, we constructed a miRNA ratio-based prediction model for ICU mortality that optimized the best combination of miRNA ratios (miR-125b-5p/miR-34a-5p, miR-199a-5p/miR-9-5p, and miR-221-3p/miR-491-5p). The model (AUC 0.85) and the miR-199a-5p/miR-9-5p ratio (AUC 0.80) showed an optimal discrimination value and outperformed the best clinical predictor for ICU mortality (days from first symptoms to IMV initiation, AUC 0.73). The survival analysis confirmed the usefulness of the miRNA ratio model and the individual ratio to identify patients at high risk of fatal outcomes following IMV initiation. Functional enrichment analyses identified pathological mechanisms implicated in fibrosis, coagulation, viral infections, immune responses and inflammation. Conclusions COVID-19 induces a specific miRNA signature in BAS from critically ill patients. In addition, specific miRNA ratios in BAS samples hold individual and collective potential to improve risk-based patient stratification following IMV initiation in COVID-19-related critical illness. The biological role of the host miRNA profiles may allow a better understanding of the different pathological axes of the disease.

optimal discrimination value and outperformed the best clinical predictor for ICU mortality (days from first symptoms to IMV initiation, AUC 0.73). The survival analysis confirmed the usefulness of the miRNA ratio model and the individual ratio to identify patients at high risk of fatal outcomes following IMV initiation. Functional enrichment analyses identified pathological mechanisms implicated in fibrosis, coagulation, viral infections, immune responses and inflammation.

INTRODUCTION
The characteristic feature of patients with severe COVID-19 is the development of acute respiratory distress syndrome (ARDS). Consequently, more than 15% of hospitalized patients are eventually transferred to the intensive care unit (ICU), where 70-90% of these patients receive invasive ventilation and 20-50% present fatal outcomes (1)(2)(3). Although mortality has been mainly concentrated in patients admitted to the ICU, especially in those who require invasive mechanical ventilation (IMV), there is a lack of tools for predicting the progression of the disease in critically ill patients (4). The development of prognostic assays reflecting the risk of clinical decompensation will be useful for medical decision-making. This approach will also provide valuable molecular information on the pathological mechanisms linked to the severe clinical courses of COVID-19. In this scenario, the use of relatively novel omic technologies, such as transcriptomics, emerges as an interesting approach for in-depth molecular phenotypic analyses.
MicroRNAs (miRNAs) are non-coding RNAs (ncRNAs) consisting on transcripts of 19-25 nucleotides that regulate gene expression posttranscriptionally (5). MiRNAs are finetuning regulators of diverse biological processes, including differentiation, proliferation and migration; and specially, the regulation of homeostasis and the response to stress (6). Experimental evidence has suggested that these small transcripts are critical components in viral infections and the host defense, particularly among RNA viruses (7). MiRNAs are also found in the extracellular environment. The release of miRNAs into the extracellular space can be passive, resulting from tissue damage, or active, which is due to secretion mechanisms (8). Actively secreted miRNAs participate in cell-to-cell communication as signaling molecules by regulating the gene expression of recipient cells (9). In addition, the extracellular forms of miRNAs have been recognized as biomarkers with great potential for patient management. Their concentration and expression profile vary in response to physiological stress or pathological insults (10,11). MiRNAs are highly stable and can be quantified in accessible specimens with techniques currently used in clinical laboratories, i.e., quantitative PCR (qPCR). The role as clinical indicators of miRNAs quantified in a wide array of bodily fluids has been reported (12)(13)(14). Although the incorporation of miRNAbased biomarkers to standard clinical care must address several methodological limitations (15), a great progress has been made in the past decade. Some cost-effective tests are currently available in the market to assist in the management of malignant diseases and osteoporosis (16,17).
While the ICU patient is assisted by IMV, tracheobronchial toilet is crucial to treat the airways obstructed by secretions causing lung atelectasis and the consequent gas exchange deterioration. The specimen extracted by flexible bronchoscopy is known as bronchial aspirate (BAS). BAS is considered a sample type of great interest at the research and clinical levels since it can reflect pathophysiological processes in the local lung environment (18). Furthermore, BAS sampling is a component of standard care, and its collection is less invasive than that of other respiratory specimens, such as bronchoalveolar lavage fluid (BALF).
We hypothesized that critical life-threatening infections induce characteristic molecular changes that can be detected in BAS samples. We first aimed to compare the miRNA profile of BAS samples from critically ill COVID-19 and non-COVID-19 patients to identify molecular pathways associated with the disease. Then, we analyzed the expression profile of miRNAs in BAS samples from survivors and non-survivors of ICU stays. In this second objective, we sought to identify miRNA signatures as prognostic biomarkers of a fatal outcome in COVID-19-related critical illness and to define molecular pathways involved in the development of adverse events. Eighteen COVID-19-positive patients confirmed by nasopharyngeal swab PCR and 14 non-COVID-19 patients admitted to the ICU from April 2020 to August 2020 were included. Both study groups were matched by age.

Study Design and Patients
Comprehensive clinical data were manually extracted from the electronic medical records by specialized clinical research assistants.
The samples and data from patients included were provided by Biobank IdISBa B.0000527 (www.idisba.es) and CIBERES Pulmonary Biobank Consortium B.0000471, a network currently formed by twelve tertiary Spanish hospitals (www.ciberes.org) integrated in the Spanish National Biobanks Network. Samples were processed following standard operating procedures with the appropriate approval of the Ethics and Scientific Committees (s007-BBCOV) and with the collaboration of the healthcare services of the Hospital Universitario Son Espases and Hospital Son Llatzer (Palma, Spain).

Study 2. BAS miRNA Profiling of Survivors and Non-survivors to COVID-19-Related Critical Illness
This was a prospective, observational, single-center cohort study. The cohort was composed of consecutive patients 18 years or older with a positive nasopharyngeal swab RT-qPCR test result for SARS-CoV-2 diagnosed with ARDS secondary to COVID-19 infection and under IMV by orotracheal intubation or tracheostomy hospitalized in the ICU from the Hospital Arnau de Vilanova (Lleida, Spain) from March 2020 to January 2021. Demographic, clinical and pharmacological data and information regarding hospital admission, ICU admission, initiation of IMV and outcomes were extracted from the medical records by specialized clinical research assistants. The clinical endpoint of the study was ICU mortality, days under IMV, duration of ICU stay and duration of hospital stay.
Samples were obtained with support by IRBLleida Biobank (B.0000682) and "Plataforma Biobancos PT17/0015/0027." The procedure to obtain the samples was approved by the Research Ethics Committee (CEIC 2273) and was performed in compliance with the ethical requirements according to the Declaration of Helsinki. The participating patients or their relatives were informed about the research and gave their written informed consent before the use of their biological samples and clinical information in the study.

Bronchial Aspirates Collection
BAS samples were obtained during the ICU stay as part of the clinical management with two objectives: i) to perform toilet bronchoscopy to aspirate retained secretions and to prevent or to revert lung atelectasis and ii) to rule out superinfections. Bronchoaspiration was performed by flexible bronchoscopy with an Ambu R aScope TM 4 Broncho Large 5.8/2.8 bronchoscope (Ambu, Ballerup, Denmark) connected to a vacuum. Patients undergoing this procedure were under continuous sedoanalgesia; in addition, extra boluses of short-acting sedatives (midazolam 5-10 mg) and analgesics (fentanyl 0.05-0.15 mcg) were required to ensure proper tolerance. Bronchial secretion samples (2-5 mL) were obtained. In some cases, simple bronchial washing with 0.9% saline solution (Braun, Melsungen, Germany) was required (19). The samples were immediately aliquoted, frozen and stored at −80 • C. For study population 1, the frozen aliquots were shipped on dry ice to the Lleida Institute for Biomedical Research (Lleida, Spain). For study population 2, the samples were stored at the biobank of the Lleida Institute for Biomedical Research.

RNA Isolation
RNA isolation and miRNA quantification were performed by experienced staff in a blinded fashion. All experiments were performed under standardized conditions in the same laboratory.
A total of 39 miRNAs were selected according to preset criteria by experienced researchers after a comprehensive literature search (Supplementary Table 1). MiRNAs have been previously described as biomarkers of and/or mediators in molecular pathways related to respiratory viral infection, lung damage/fibrosis, immune response, inflammation and coagulation.
BAS samples were thawed at 4 • C. RNA isolation was performed in aliquots that had undergone their first thaw. Total RNA was isolated from 100 µL of BAS using the miRNeasy Mini Kit (Qiagen, Hilden, Germany). Five volumes of QIAzol lysis reagent was mixed with one volume of BAS and incubated for 5 min at room temperature. The RNA Spike-In Kit (synthetic UniSp2, UniSp4, and UniSp5) (Qiagen) and Caenorhabditis elegans cel-miR-39-3p (1.6 × 10 8 copies/µL) (Qiagen) were added to all extractions to monitor RNA isolation efficiency. The mixture was supplemented with 1 µg of MS2 bacteriophage RNA (Roche, Merck, Darmstadt, Germany), an RNA carrier not containing miRNAs, to improve RNA yield. All reagents were spiked into samples during RNA isolation after incubation with the denaturing solution. Subsequently, one volume of chloroform was added, and after 3 min at room temperature, the mixture was centrifuged at 12,000 g and 4 • C for 15 min. The upper aqueous phase was transferred to a fresh reagent tube, and 1.5 volumes of 100% ethanol were added. Purification of RNA was performed with RNeasy MinElute spin columns according to the manufacturer's instructions. RNA was eluted in 30 µL RNase-free H 2 O and stored in a −80 • C freezer.

RT-qPCR
For isolated total RNA, retrotranscription (RT) was performed using the miRCURY R LNA R RT Kit (Qiagen) according to the manufacturer's instructions. One microliter of RNA was reverse transcribed in 10 µL reactions. Additional spike-in (synthetic UniSp6) (Qiagen) was added to the cDNA synthesis reaction to check for RT robustness and PCR inhibitors. RT reaction was performed with the following conditions: incubation for 60 min at 42 • C, heat inactivation for 5 min at 95 • C, and immediate cooling to 4 • C. Then, cDNA was stored at −20 • C. miRNAs were quantified by qPCR using the miRCURY LNA SYBR R Green PCR Kit (Qiagen) in 384-well miRCURY LNA miRNA Custom PCR Panels (Qiagen) with precoated primers, according to the manufacturer's instructions. This technology is sensitive and specific for miRNA quantification in samples with low RNA yield and allows to discriminate between the target and other closely related mature miRNAs or miRNA precursors. The system provides optimal reproducibility and analytical accuracy (20). cDNA was diluted 40x, and 4 µL were used in 10 µL qPCR reactions. Synthetic UniSp3 was included as an interplate calibrator and qPCR control. The qPCR was performed using a QuantStudio TM 7 Flex Real-Time PCR System under the following conditions: 2 min at 95 • C, 40 cycles of 10 s at 95 • C and 60 s at 56 • C and, finally, melting curve analysis.

Expression Analysis
QuantStudio Software v1.3 (Thermo Fisher Scientific, Massachusetts, USA) was used for both the determination of the quantification cycle (Cq) and the melting curve analysis. Cq was defined as the fractional cycle number at which the fluorescence exceeded a given threshold. To ensure the highest data quality, RNA and cDNA spike-in levels were analyzed during miRNA quantification. The synthetic spike-in UniSp3 was also used as an interplate calibrator. The miRNAs were considered to be expressed when Cq values < 35. Samples that were below the limit of detection were censored at the lowest level observed for each miRNA. The miRNA ratios were used to reduce sample-to-sample variation in miRNA expression levels. This approach minimizes technical variability caused by sample characteristics, collection and processing (21).

Pathway and Gene Ontology Enrichment Analyses
Pathway and Gene Ontology (GO) analyses were performed by using the web-based computational tool DIANA-miRPath v3.0 (22) (April 2021). DIANA-miRPath v3.0 utilizes experimentally validated miRNA interactions derived from TarBase v7.0 and combines the results with KEGG (Kyoto Encyclopedia of Genes and Genomes) and GO analyses (biological process). The level of significance was set at a false discovery rate (FDR)-adjusted p-value < 0.001.

Statistical Analysis
The statistical analysis was performed using R software, version 4.0.2 (www.r-project.org). Descriptive statistics were used to summarize the characteristics of the study populations. The normality of the distribution was analyzed using the Shapiro-Wilk test. Data are presented as the mean (standard deviation) or median [25th percentile; 75th percentile] for continuous variables, according to normality, and as the frequency (percentage) for categorical variables. Continuous variables and categorical variables were compared between groups using Student's t-test, the Mann-Whitney U test and Fisher's exact test, respectively. Spearman's rho coefficient was used to assess the correlation between continuous variables. The miRNA ratios were calculated using the expression levels: ratio A/B = 2 −CqA / 2 −CqB . All possible miRNA ratios were generated. The selection of ratios for subsequent statistical analysis was based on the following criteria: (1) higher area under the receiver operating characteristic (ROC) curve (AUC); (2) p-value < 0.05; and (3) the same miRNA could not be part of different ratios selected for the statistical analysis. Differences in miRNA ratios between groups were evaluated using linear models for arrays (23). The same models were used to adjust miRNA ratio levels for age and sex. The principal component analysis (PCA) and hierarchical clustering included the differentially expressed miRNA ratios.
A process of miRNA ratio selection for the prediction of mortality based on a random forest algorithm was performed (24, 25). Finally, a predictive model (miRNA ratio score) with the selected miRNA ratios was fitted using a logistic regression model. ROC curves and precision-recall curves were constructed to test the ability of miRNA ratios and clinical predictors to distinguish between survivors and non-survivors. The AUC was used as the global discrimination value measure. For selected features, a miRNA ratio cutoff point was established to fit mortality risk using a maximally selected log-rank statistic (26). Finally, independent models of mortality risk were constructed with a Cox regression model including dichotomized levels of the miRNA ratios. The same analysis was performed with the predictive model. Kaplan-Meier curves were used to illustrate differences among groups in the observed time-to-event outcome. The p-value threshold defining significance was set at <0.05.

COVID-19 Is Associated With a Specific MicroRNA Profile in Bronchial Aspirates From Critically Ill Patients
We first analyzed the impact of COVID-19 on the miRNA profile in the BAS samples from critically ill patients. The main demographic, clinical, pharmacological and biochemical data of the study groups are summarized in Table 1. The causes of ICU admission among non-COVID-19 patients are displayed in Supplementary Table 2. Compared to non-COVID-19 patients, COVID-19 patients showed a significantly lower PaO 2 at the time of ICU admission. The duration of the ICU stay, the requirement for non-invasive MV and prone positioning and the use of antibiotics and hydroxychloroquine were also greater in COVID-19 patients.

ICU admission
Oxygen signaling pathway" and "ErbB signaling pathway." One hundred and thirteen GO terms were identified. Processes included "viral process, " "viral life cycle, " "viral transcription, " "positive regulation of viral transcription, " "immune system process, " "innate immune response, " "leukocyte migration, " several signaling pathways linked to Toll-like receptors, "cell death" and several pathways related to apoptotic signaling, "blood coagulation, " "platelet activation" and "platelet degranulation." Next, we sought to determine whether the miRNA signature of BAS samples can predict patient outcome. Figure 2 displays the study flowchart. The patient demographics, clinical and pharmacological characteristics according to ICU survival are presented in Table 2. ICU mortality was 31.6% (n = 18/57). Non-survivors were typically older. At the time of ICU admission, non-survivors showed lower glucose concentrations. The days until IMV ventilation since the beginning of the symptoms and the hours under prone positioning were higher in those who eventually died. The causes of death were respiratory insufficiency (66.7%), multiorgan failure (27.8%) and cardiovascular accidents (5.5%). The comparison between those patients with BAS samples and all patients under IMV or patients under IMV without BAS samples is displayed in Supplementary Tables 6, 7. Microbiological and FIGURE 2 | Study flowchart. Among 226 patients admitted to the ICU, 197 patients were positive for SARS-CoV-2. Eighteen were negative for SARS-CoV-2, four were referred to other centers early, and seven did not provide written informed consent. One hundred and eight patients were under invasive mechanical ventilation by tracheostomy or orotracheal intubation during the ICU stay. Among these, 77 patients were subjected to bronchoscopy. Fifty-seven samples were available for microRNA quantification.
fungal examinations showed no differences between the study groups (Supplementary Table 8).
A predictor selection procedure based on random forest was used to select the optimal combination of ratios for the prediction of mortality during the ICU stay. The multivariable analysis selected a 3-miRNA ratio prediction model (miRNA ratio score) composed of miR-125b-5p/miR-34a-5p, miR-199a-5p/miR-9-5p and miR-221-3p/miR-491-5p ratios (Figure 4E). The discriminatory capacity of the score composed of the three ratios was higher than that observed for individual miRNA ratios or clinical predictors (AUC 0.85) (Figure 4D). The performance of the miRNA ratio score was the highest in the precisionrecall curves (Figure 4F). The Kaplan-Meier analysis showed that patients with high levels in the score based on the 3-miRNA ratio model and the individual miR-199a-5p/miR-9-5p ratio were at significantly higher risk of in-ICU death (HR 4.42 and HR 3.8, respectively) ( Figure 4G). An inverse correlation between the ratio miR-1-3p/miR-124-3p and the number of days under IMV was observed (Supplementary Figure 2).
In silico analysis was performed to identify the molecular pathways and biological processes most closely linked to the miRNA profile from non-survivors (Supplementary Figures 3A,B and Supplementary Tables 10,  11). Twenty-five KEGG pathways and 122 GO terms were enriched with the verified targets of miRNAs that comprised the ratios. Pathways and GO terms related to COVID-19 pathophysiology were identified, including the "TGF-beta signaling pathway, " "viral process, " "viral life cycle, " "viral transcription, " "blood coagulation, " "platelet activation, " "platelet degranulation, " "immune system process, " "positive regulation of type I interferon production, " "innate immune response" and "leukocyte migration, " "cell death" and several signaling pathways linked to Toll-like receptors and apoptosis.

DISCUSSION
Significant advancements in the outcomes and treatment of patients with COVID-19 have been made since the beginning of the pandemic. Nevertheless, there is considerable room to improve the risk stratification and therapeutic management of patients with COVID-19-related critical illness (4). As such, we analyzed the host miRNA profile of BAS samples from COVID-19 and non-COVID-19 patients admitted to the ICU. We report that (i) the miRNA signature of critically ill COVID-19 patients is essentially different from that observed in non-COVID-19 patients; (ii) non-survivors of an ICU stay show a characteristic miRNA profile when compared to that of survivors; (iii) miRNA profiling in respiratory specimens allows risk stratification in COVID-19 patients under IMV; and iv) the miRNA profile provides information about the factors that mediate the disease and fatal outcomes among critically ill COVID-19 patients.

Clinical Implications for the Management of Critically Ill Patients With COVID-19 Assisted by Invasive Mechanical Ventilation
Data regarding the host miRNA profile within the lung in COVID-19 patients are limited. Aberrant individual miRNA or miRNA profiles in the lung compartment have been observed in respiratory viral infections, such as influenza A virus (27). Since the beginning of the pandemic, most studies have explored  Continuous variables are expressed as the median [P25; P75] or mean (SD) and categorical variables are expressed as n (%). ALT, alanine aminotransferase; AST, aspartate aminotransferase; COPD, chronic obstructive pulmonary disease; FiO 2 , fraction of inspired oxygen; ICU, intensive care unit; IMV, invasive mechanical ventilation; LDH, lactate dehydrogenase; PaCO 2 , carbon dioxide partial pressure; PaO 2 , oxygen partial pressure; SaO 2 , oxygen saturation.
alterations in the miRNA profile of matrices such as whole blood, plasma and serum (28)(29)(30)(31). However, the identification of biomarkers for lung disease using blood samples is controversial, at least in pneumonia-related ARDS. Lung compartmentalization exists for specific pathological mechanisms, e.g., inflammation (32). From a practical standpoint, miRNA assessment in BAS emerges as an interesting tool to define phenotypes among ICU patients under IMV. Of note, BAS may reflect in situ host response against SARS-CoV-2 infection in the lung, is relatively easy to collect and can replace other procedures that are not well-tolerated by the patient, such as BALF collection (33). Our results demonstrated for the first time the biomarker potential of miRNAs in BAS samples. To our knowledge, this is also the first study that reports host miRNA profiles in respiratory specimens as predictors of clinical outcomes in COVID-19related critical illness. Critically ill patients with COVID-19 under IMV can be classified by multi miRNA signatures quantified in the respiratory tract. In particular, five miRNA ratios segregated between ICU survivors and non-survivors. Simultaneous assessment of miRNA profiles or miRNA ratios may hold promise to provide more comprehensive information about the clinical evolution of the patient. This hypothesis, which seems to be true for circulating miRNAs (34), can also be translated to respiratory specimens such as BAS. Using multivariable analysis for multimarker analysis, we constructed a miRNA ratio-based prediction model (miRNA ratio score) to optimize the best combination of ratios for risk prediction. The score based upon miRNA expression ratios (miR-125b-5p/miR-34a-5p, miR-199a-5p/miR-9-5p, and miR-221-3p/miR-491-5p) and the individual miR-199a-5p/miR-9-5p ratio showed excellent discrimination value; AUC 0.85 and 0.80, respectively. The survival analysis confirmed the great potential of the miRNA ratio score and the miR-199a-5p/miR-9-5p ratio for identifying patients at very high risk of fatal outcomes following IMV initiation.
To date, the prognosis of critically ill patients is mostly based on clinical characterization. Consequently, the estimated ratios were also compared with sociodemographic variables, clinical parameters or laboratory parameters, including a comprehensive set of indicators collected at the time of ICU admission or initiation of IMV. Some of these parameters, including the PaO 2 :FiO 2 ratio, leukocyte counts and D-dimer, C-reactive protein (CRP) and creatinine concentrations, among others, have been previously reported as mortality predictors in ventilated and/or ICU COVID-19 patients (2,(35)(36)(37). The miRNA ratio signature and the miR-199a-5p/miR-9-5p ratio outperformed the best clinical predictor for ICU mortality (days from first symptoms to IMV initiation, AUC 0.73).   Based on these findings, miRNA quantification in BAS samples may support patient stratification and constitute an objective tool for individualizing and guiding clinical care during the ICU stay, e.g., treatment replacement and the use of novel medications in addition to conventional medical treatment or more intensive monitoring and allocation of hospital resources.

Molecular Mechanisms Associated With COVID-19 and Fatal Courses of the Disease
Previous investigations have proposed crucial roles of ncRNAs in the pathophysiology of COVID-19 (29). Since these RNAs, such as miRNAs, can be measured in respiratory secretions, ncRNAs constitute an interesting tool to molecularly phenotype patients. The current results may improve the limited understanding of SARS-CoV-2 pathogenesis and the host biological response. This information is particularly relevant due to the atypical features linked to critically ill COVID-19 patients. Notably, miRNA targeting has been previously reported as an interesting approach for the treatment of viral infections (38) and chronic disease (39).
We report that SARS-CoV-2 infection induces a unique miRNA profile in BAS samples from COVID-19 patients. Five miRNA ratios distinguished between laboratory-confirmed infected and non-infected patients. Functional enrichment analyses using experimentally validated miRNA:gene interactions suggested the implication of these transcripts in different pathophysiological axes of the disease and are in line with previous findings that demonstrated a characteristic immunoinflammatory profile of COVID-19 patients (40). A number of pathways and GO terms were associated with mechanisms such as immune cell differentiation and activation, cytokine and chemokine synthesis and the regulation of the inflammatory process. Several GO terms were linked to the coagulation process, a hallmark of COVID-19. Notably, the ratio containing miR-155-5p, proposed as a key regulator of T-cell maturation and the innate immune response (41), showed significant differences between the study groups. Some of these miRNAs are also known as mediators of viral infections. Hepatitis C virus uses miR-122 to facilitate its replication (38). In addition, miR-155 expression is required for the growth of Epstein-Barr virus-infected B cells (42). Interestingly, infection by SARS-CoV in bronchoalveolar stem cells causes miR-214 overexpression, presented in the ratio that showed higher differences, which has been described as a mechanism for evading immune elimination until successful transmission takes place (43).
Additionally, miRNA ratios may be informative about the driving factors and triggers of fatal COVID-19 forms. miR-199-5p has been described as a miRNA with pro-viral functions. For instance, it enhances hepatitis C virus replication, and its downregulation suppresses virus replication in an in vitro knockdown model (44). Furthermore, miR-199-5p, together with miR-214-3p, has also been proposed as a major regulator of lung fibrosis. Independent investigations using in vitro and in vivo models reported this miRNA as a key effector of TGF-β-induced lung myofibroblast activation (45,46). Pulmonary miR-199-5p expression is increased in patients with idiopathic forms of fibrotic lung disease (46), and its extracellular levels are elevated in bodily fluids (e.g., sputum) in respiratory conditions (47). The possible role of miR-9-5p in SARS-CoV-2 infectivity should not be disregarded. Angiotensin-converting enzyme 2 (ACE2) mRNA presents high-probability miRNA binding sites for this transcript in the 3 ′ -UTR (48). Supporting the participation of miR-9-5p in the regulation of coronavirus infections, the human coronavirus OC43 nucleocapsid protein binds miR-9, a negative regulator of the multifunctional transcription factor NF-κB, which ultimately impacts the regulation of the inflammatory response to viral infection (49). In general, there is a preponderance of miRNAs related to inflammatory mechanisms. miR-16-5p, miR-124-3p, miR-125b-5p, miR-126-3p, and miR-221-3p are able to target molecules belonging to pathways with a key role in inflammatory responses, including receptors (TLRs), signaling molecules (MYD88), transcription factors (NF-κB) and cytokines/chemokines (TNF, IL-6 and IL-8) (50)(51)(52)(53). Some miRNAs, such as miR-126-5p and miR-221-3p, may be associated with the severe endothelial injury and coagulopathic features observed in lung samples from fatal cases (54). Both miRNAs, in addition to miR-16-5p, may also be linked to the vascular angiogenesis previously described in the lungs of patients who have succumbed to the disease (55,56). In the context of coagulation, recent findings demonstrated that miR-34a-5p is downregulated in postmortem lung biopsies of COVID-19 patients with severe respiratory injuries and thrombotic events (57). The miRNA target sites (MTSs) within SARS-CoV-2 have been identified for this miRNA (58). The same authors suggested that the virus genome may modulate the levels of miR-34a-5p, acting as a miRNA sponge to regulate unfolded protein response (UPR)-related apoptosis. The KEGG and GO analyses support the participation of these host miRNAs in the pathogenesis of the disease.
Overall, the findings add additional information to genomic, proteomic, metabolic and cellular data to identify new drug targets or therapeutic agents (59). This study provides useful hypothesis-generating data and constitutes a basis for future investigations in critically ill patients. The miRNAs described in the current study deserve additional analyses as therapeutic targets/agents. The role as mediators or consequence of the disease and the pathophysiological importance in fatal outcomes should be explored in detail using alternative approaches, including mechanistic in vitro and in vivo studies.

Strengths and Limitations of the Study
The strengths of the study include the rigorous quality control of clinical and miRNA data, the exploration of a real-world setting that included consecutive patients who underwent bronchoscopy and the evaluation of the miRNA ratios in conjunction with comprehensive clinical information. Nevertheless, the conclusions should be interpreted in the context of several limitations. First, mortality was evaluated in a single-center study that included a specific subpopulation of ICU patients, i.e., patients assisted by IMV and subjected to bronchoscopy. Although the highest mortality rates during the pandemic have been concentrated in the group of patients who require IMV secondary to ARDS, the predictive power of the miRNA signature remains to be validated in alternative patient subgroups. Second, the impact of the disease on the miRNA profile and the potential of these small transcripts as prognostic biomarkers need to be confirmed in larger populations. Third, the non-COVID-19 group was composed of patients presenting heterogeneous causes of ICU admission. The main causes of ICU admission among non-COVID-19 patients were pneumonia (35.7%) and brain edema (14.3%) (Supplementary Table 2). The inclusion of a group of patients with viral pneumonia constitutes an interesting approach to characterize with more detail the molecular impact of SARS-CoV-2 infection. Forth, we could not exclude the impact of complications and treatments during ICU stay that were not recorded. Fifth, the source for the miRNA profile in BAS is not fully defined and could reflect the intracellular profile of or release from immune, pulmonary, epithelial or vascular cells. It is worth noting that the contribution of tissue damage to the BAS miRNA profile should not be disregarded. MiRNAs are passively released to the extracellular space during cell death, with cardiomyocytespecific miRNAs; i.e., miR-208a-3p and miR-208b-3p, as familiar examples (60). Additional approaches are fundamental to infer any causal association. Sixth, miRNA profiling in BAS is challenging because of the extraordinary complexity and broad dynamic range of miRNA concentrations. BAS shows a high sample-to-sample variation, which is due to the collection of the sample during bronchoaspiration (e.g., bronchial washing with saline solution) that cannot be accounted for using standard normalization methods. The comparison of individual miRNAs between BAS samples is not possible. To overcome this limitation, we analyzed miRNA ratios, which precludes the need for normalization factors instead of individual miRNAs.

CONCLUSIONS
COVID-19 induces a specific host miRNA profile in BAS samples from critically ill patients. Furthermore, specific miRNA ratios in respiratory secretions, particularly BAS, show individual and collective potential for risk-based patient stratification following IMV initiation in COVID-19-related critical illness. The biological role of the miRNA profiles in BAS samples allows a better understanding of the different pathological axes of the disease. The role of miRNAs as mediators and biomarkers should be further explored.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Files, further inquiries can be directed to the corresponding author/s.

ETHICS STATEMENT
Samples were processed following standard operating procedures with the appropriate approval of the Ethics and Scientific Committees (s007-BBCOV) and with the collaboration of the healthcare services of the Hospital Universitario Son Espases and Hospital Son Llatzer (Palma, Spain). and DdG-C contributed to the data analysis and interpretation. DdG-C is the guarantor of the paper, had full access to all of the data in the study, and takes responsibility for the integrity of the data and the accuracy of the data analysis. All authors contributed to the manuscript draft, critically revised the manuscript for important intellectual content, and approved the final version.

FUNDING
This work was supported by Instituto de Salud Carlos III (COV20/00110), co-funded by European Regional Development Fund (ERDF)/A way to make Europe. CIBERES is an initiative of the Instituto de Salud Carlos III. Supported by: Programa de donaciones estar preparados; UNESPA (Madrid, Spain). DdG-C (Miguel Servet 2020: CP20/00041) and MM (PFIS: FI21/00187) have received financial support from the Instituto de Salud Carlos III, co-funded by the European Social Fund (ESF)/Investing in your future. MM is the recipient of a predoctoral fellowship (PERIS, PIF-Salut, SLT017/20/000142) from the Department de Salut (Generalitat de Catalunya). MG-H is the recipient of a predoctoral fellowship from University of Lleida.
participation of the Units of Intensive Care, Clinical Analysis and Pulmonology of Hospital Universitario Son Espases and Hospital Son Llatzer for their collaboration. This work was also supported by IRBLleida Biobank (B.0000682) and Plataforma Biobancos PT17/0015/0027/.