Analysis of Differentially Expressed MicroRNAs in Serum and Lung Tissues from Individuals with Severe Asthma Treated with Oral Glucocorticoids

Abstract

Nowadays, microRNAs (miRNAs) are increasingly used as biomarkers due to their potential contribution to the diagnosis and targeted treatment of a range of diseases. The aim of the study was to analyze the miRNA expression profiles in serum and lung tissue from patients with severe asthma treated with oral corticosteroids (OCS) and those without OCS treatment. For this purpose, serum and lung tissue miRNAs of OCS and non-OCS asthmatic individuals were evaluated by miRNAs-Seq, and subsequently miRNA validation was performed using RT-qPCR. Additionally, pathway enrichment analysis of deregulated miRNAs was conducted. We observed altered expression by the next-generation sequencing (NGS) of 11 miRNAs in serum, of which five (hsa-miR-148b-3p, hsa-miR-221-5p, hsa-miR-618, hsa-miR-941, and hsa-miR-769-5p) were validated by RT-qPCR, and three miRNAs in lung tissue (hsa-miR-144-3p, hsa-miR-144-5p, and hsa-miR-451a). The best multivariate logistic regression model to differentiate individuals with severe asthma, treated and untreated with OCS, was to combine the serum miRNAs hsa-miR-221-5p and hsa-miR-769-5p. Expression of hsa-miR-148b-3p and hsa-miR-221-5p correlated with FEV₁/FVC (%) and these altered miRNAs act in key signaling pathways for asthma disease and the regulated expression of some genes (FOXO3, PTEN, and MAPK3) involved in these pathways. In conclusion, there are miRNA profiles differentially expressed in OCS-treated individuals with asthma and could be used as biomarkers of OCS treatment.

1. Introduction

Asthma is a heterogeneous disease of the lower airways that affects approximately 300 million people worldwide [1]. Chronic airway inflammation as well as airway damage and remodeling are hallmark features of the disease and induce bronchial hyper-reactivity and variable expiratory airflow limitation [2]. Individuals with asthma exhibit a broad spectrum of disease severity. Severe asthma, defined as that which requires treatment with high-dose inhaled corticosteroids plus a second controller and/or systemic corticosteroids [3], affects about 5–10% of the asthmatic population, and is associated with increased mortality and morbidity, diminished quality of life, and greater health expenditure [4].

The pathophysiology of asthma is highly complex, comprising distinct phenotypes and endotypes [5]. In addition, a broad range of triggers including environmental risk factors, genetic polymorphisms, and epigenetic changes contribute to the onset of asthma, variations in phenotypes, and response to steroid treatment [6]. Severe asthma comprises two predominant inflammatory endotypes: Type 2 (T2)-high and T2-low disease [7]. T2-high asthma is more severe and difficult to treat, and many patients with this endotype require high-dose treatments and/or biological therapy for better control of the symptoms and exacerbations, and to avoid adverse reactions caused by the administration of oral corticosteroids (OCS) [8,9].

To help guide the selection of treatment in patients with severe asthma, currently available biomarkers such as peripheral blood and induced sputum eosinophil count, fraction of exhaled nitric oxide (FeNO) level, and immunoglobulin-E (IgE) levels have been used [10]. However, interest in using microRNA (miRNA) profiles as biomarkers for diseases is currently on the rise [11]. These small, single-stranded, non-coding RNA molecules [12] participate in the regulation of gene expression by inhibiting protein translation [13] and are useful as biomarkers as they are stable and disease-related particles [14]. Previous studies have demonstrated that miRNAs hold great potential as asthma mediators and biomarkers and can further be of use in asthma endotyping and phenotyping, thereby making it possible to offer personalized therapies for patients [15,16]. Among them, Mirra et al. recently identified a new expression profile of circulating miRNAs including a member of the miR-181 family, related to lung inflammation, which could be used as a clinical marker of bronchial asthma and response to pharmacological treatment [17].

The aim of this study was to search for a profile of differentially expressed miRNAs between individuals with severe asthma receiving OCS treatment versus those without OCS therapy using next-generation sequencing (NGS). The miRNAs found could serve as lung-tissue or systemic biomarkers of oral corticosteroid treatment and facilitate the selection of the most appropriate management approach.

2. Results

2.1. Characteristics of the Study Population

Descriptive data reflecting the demographic, inflammatory, functional, and clinical characteristics of all individuals studied are summarized in Table 1 and Table 2. The subjects were distributed into two different analyses: 20 OCS-treated versus (vs.) 26 non-OCS-treated severe asthmatics (serum samples analysis) and six OCS-treated vs. seven non-OCS-treated asthmatic subjects (lung biopsy samples analysis).

As shown in Table 1, regarding inflammatory parameters, the OCS-treated and non-OCS-treated severe asthmatics from which serum samples were used had significant differences in IgE levels (366.0 and 96.8 IU, p < 0.05). Additionally, as expected, we observed significant differences in treatment with systemic corticosteroids (p < 0.0001). Other characteristics (i.e., demographic, inflammatory, functional, clinical) did not present statistically significant differences.

Moreover, comparing the treated and non-treated with OCS groups in the analysis of lung biopsy samples (Table 2), the only significant differences were detected in the percentage of monocytes (10.0% vs. 5.9%, p < 0.05) and, as planned, in the oral corticosteroid treatment (p < 0.01).

2.2. MiRNA Profile in Serum Samples

Analysis of the data obtained by NGS (miRNAs-seq) showed 11 differentially expressed serum miRNAs (hsa-miR-148b-3p, hsa-miR-221-5p, hsa-miR-618, hsa-miR-200b-3p, hsa-miR-548I, hsa-miR-200a-3p, hsa-miR-941, hsa-miR-181d-5p, hsa-miR-769-5p, hsa-miR-133a-3p, and hsa-miR-3614-5p) between the OCS-treated and non-OCS-treated patients with severe asthma (adjusted p < 0.05) (Figure 1a). All of these 11 miRNAs were overexpressed in the OCS-treated subjects.

When PCA was performed using the miRNA expression values, the results revealed that both groups were clearly differentiated (Figure 1b).

After the miRNAs-Seq analysis, RT-qPCR validation of miRNAs was performed in a larger cohort of serum samples. Of the 11 differentially expressed miRNAs obtained by miRNAs-Seq between the OCS-treated and non-OCS-treated individuals with severe asthma, we found significantly increased levels in five of them (hsa-miR-148b-3p, hsa-miR-221-5p, hsa-miR-618, hsa-miR-941, and hsa-miR-769-5p) in the OCS-treated patients (Figure 1c). Hsa-miR-3614-5p was not detected by RT-qPCR, and the remaining miRNAs evaluated did not reach statistical significance (Figure 1c).

These results confirm that serum hsa-miR-148b-3p, hsa-miR-221-5p, hsa-miR-618, hsa-miR-941, and hsa-miR-769-5p were differentially expressed miRNAs between the OCS-treated and non-OCS-treated subjects with severe asthma, suggesting that these could serve as markers of treatment.

2.3. MiRNAs Correlation and Signaling Pathways

In order to establish some relation between the expression levels of these serum miRNAs (ΔCt) and different clinical characteristics, we performed correlation analyses. We observed a significant direct correlation between FEV₁/FVC and hsa-miR-148b-3p in the two groups (correlation coefficient (r) = 0.63 for the OCS-treated patients with severe asthma or correlation coefficient (r) = 0.59 for non-OCS-treated subjects with severe asthma) and hsa-miR-221-5p among the OCS-treated individuals with severe asthma (correlation coefficient (r) = 0.56) (Figure 2a). Since there is a correlation between FEV₁/FVC and hsa-miR-221-5p in the OCS-treated severe asthmatics, it can explain the difference in treatment between the two groups, whereas, since FEV₁/FVC correlated with hsa-miR-148b-3p in both groups, we can say that, in this case, there is no explanation for the difference in treatment, but there was a link with asthma pathology. Furthermore, despite not being a result of statistical significance, it is worth noting the opposite trend between the groups of the OCS-treated and non-OCS-treated asthmatic individuals, in the correlation of hsa-miR-148b-3p with the FeNO levels, perhaps because OCS treatment modified the FeNO/miR-148b-3p relation.

ROC curves were generated, and we calculated the areas under the curve (AUC). Multivariate logistic regression models (continuous predictors) were created for better discrimination of the OCS-treated and non-OCS-treated patients with severe asthma (Supplementary Table S1). Hsa-miR-221-5p and hsa-miR-769-5p showed AUC values of 0.75 and 0.72, respectively, indicating that both are acceptable as univariate predictors. Furthermore, when both serum miRNAs were combined, the AUC was even better (AUC = 0.77), as this was the best multivariate logistic regression model for differentiating OCS-treated and non-OCS-treated severe asthmatics (Figure 2b).

When an in silico analysis was carried out to determine their involvement in the biological processes with the differently expressed miRNAs, we determined that these miRNAs are involved in the regulation of several crucial asthma pathways, which are extracellular matrix (ECM)–receptor interaction, fatty acid biosynthesis, steroid biosynthesis, the Hippo signaling pathway, and adherens junction (Figure 3a). As TGFRB1, FOXO3, and PTEN are target genes of hsa-miR-148b-3p, and MAPK3 is the objective gene of hsa-miR-221-5p, all are involved in the signaling pathways related to asthma disease. The RT-qPCR assay was conducted in lung tissue samples from asthmatic patients and healthy subjects, revealing that FOXO3, PTEN, and MAPK3 were significantly increased in asthmatic individuals compared to the controls (Figure 3b).

Furthermore, to verify the association between miRNAs and the expression levels of their respective target genes (ΔCt), correlation analyses were carried out. We found a significant positive correlation for the total population between hsa-miR-221-5p and their target gene, MAPK3 (Spearman r = 0.35) (Figure 3c).

Finally, we examined the relation between the miRNA expression levels and eosinophil counts. We showed a statistically significant direct correlation between hsa-miR-221-5p with the percentage of eosinophils in the different population subsets (i.e., total subjects and asthmatic individuals) (Figure 3d).

2.4. MiRNA Profile in Lung Samples

Therefore, validated serum miRNAs were analyzed in the target tissue of asthma disease. We found no statistically significant differences in the expression of these miRNAs in the lung biopsies. For this reason, we decided to conduct miRNA-Seq in these samples from asthmatic patients in order to understand the miRNA profile in lung tissue. Here, three miRNAs were differentially expressed between the OCS-treated and non–OCS-treated asthmatics (q ≥ 8). Hsa-miR-144-3p, hsa-miR-144-5p, and hsa-miR-451a were overexpressed in the OCS-treated patients (Figure 4).

We then performed RT-qPCR on asthmatic lung biopsies from a larger cohort to validate the miRNAs obtained by sequencing analysis. None of the three miRNAs were validated, however, we observed a similar expression tendency in OCS-treated asthmatic patients (2^−ΔΔCt for each hsa-miR: 144-3p = 6.5, 144-5p = 5.2 and 451a = 2.0), as obtained in miRNA-Seq.

3. Discussion

This article draws on data from a comprehensive study of the miRNAs differently expressed by NGS in the serum (of which 5 were validated: hsa-miR-148b-3p, hsa-miR-221-5p, hsa-miR-618, hsa-miR-941, and hsa-miR-769-5p) and in lung tissue (hsa-miR-144-3p, hsa-miR-144-5p and hsa-miR-451a) between the OCS-treated and non-OCS-treated patients. These miRNAs could be used as epigenetic biomarkers to distinguish individuals with severe asthma under OCS treatment, although we do not know whether this miRNAs alteration is due to treatment or a characteristic of patients that need OCS in their asthma control.

Asthma is a heterogeneous disease comprising different phenotypes and endotypes, and as a result, patient stratification by inflammatory endotype is a central component of the algorithm used to assess and manage severe disease [18,19]. Identification and application of biomarkers used to identify phenotypes and endotypes of severe asthma and guide therapeutic strategy is of increasing interest to clinicians, as individual approaches to disease management and personalized medicine are imperative [20]. The severe asthma endotype is currently defined by biomarkers of varying accuracy including blood and sputum eosinophils, periostin, FeNO, and IgE for T2 asthma [16,21,22]. As reported in the literature, we found that serum samples of OCS-treated asthmatic patients had elevated IgE levels compared to the non-treated asthmatics [23]. Here, it is worth mentioning that the expression of biomarkers may be modified or even suppressed by pharmacological treatment (i.e., oral corticosteroids for blood eosinophils), making it advisable to examine biomarkers in combination. For example, we observed that the combination of two serum miRNAs used as biomarkers (hsa-miR-221-5p and hsa-miR-769-5p) resulted in a better AUC than when these were analyzed individually, suggesting that the profile of these two miRNAs was greater than using them separately.

Biomarkers are promising tools for the recognition and management of non-adherence to inhaled and oral corticosteroid therapy in asthma and can help to identify potential responders to therapeutic options, thus enabling their use as predictors of treatment response [24]. In recent years, miRNAs have been proposed as diagnostic and prognostic biomarkers for diseases [25]. NGS is a valuable technique in these efforts and can be employed for whole genome analysis including mRNA and small RNA expression such as miRNAs [26]. Our group previously used NGS to describe 15 miRNAs differentially expressed between eosinophilic and non-eosinophilic asthmatic patients, two of which were validated and could serve as instruments to classify patients into different phenotypes/endotypes [15]. In addition, hsa-miR-144-3p was previously described by our group to present an increase in miRNA levels in corticosteroid treated subjects by RT-qPCR [27], which validates the result obtained in miRNAs-Seq from the lung biopsy samples of asthmatics.

In this study, we found different profiles of miRNAs in the serum and lung tissue. Dysregulation of hsa-miR-148b-3p, hsa-miR-221-5p, hsa-miR-618, hsa-miR-941, hsa-miR-769-5p, hsa-miR-144-3p, hsa-miR-144-5p, and hsa-miR-451a has been previously associated with different aspects of asthma pathogenesis and/or other respiratory diseases [28,29,30,31,32]. Nevertheless, these profiles have never been described as indicators of treatment. These miRNAs could be relevant molecules with implications for the phenotypic distinction of patients based on whether they are treated or untreated and may be indicated for use in asthma pathogenesis, since they were significantly correlated with clinical parameters. Interestingly, we found different miRNA profiles in asthmatics regarding the sample origin, which underscore the roles of these molecules as epigenetic regulators, as found in previous studies of miRNA tissue specificity [33].

Furthermore, we observed an enrichment in the ECM–receptor interaction pathway, which is widely associated with asthma due to the central role it plays in the development of airway inflammation and remodeling associated with asthma [2]. Moreover, miRNAs have been shown to affect the expression of multiple genes directly or indirectly, and for this purpose, we studied target genes related to differential miRNA expression. Here, FOXO3, PTEN, and MAPK3 exhibited increased expression in asthmatic samples, which correlated with the miRNA expression studied, which emphasizes the implication of these inflammatory genes in the disease, as previously described by others [34,35,36]. In addition, in this study, we observed a positive correlation between hsa-miR-221-5p and MAPK3, its target gene, instead of detecting, as usual, an inverse miRNA-target gene correlation. Given that miRNAs generally behave as silencers of their target mRNAs, a possible explanation would be that miRNA–mRNA relationships are complex, and the lung tissue samples in which they have been studied differ from each other (influence of genetic and/or environmental factors), while there may be other miRNAs involved in the modulation of a particular mRNA (MAPK3 mRNA), and there could also be interaction of other target genes of this miRNA (hsa-miR-221-5p) that could be inhibitors of MAPK3, in this case, remarking the intricacy of miRNA–mRNA post-transcriptional regulation in whole tissue.

One of the main limitations of this study is the lack of the same-subject pre-treatment samples, which made it impossible to determine whether the difference in the miRNA profiles between the groups of patients treated with OCS therapy was a consequence of the treatment (reflects treatment response) or if, in contrast, the existence of differentially expressed miRNAs between the groups was a baseline characteristic of the disease and selected the treatment. Another limitation is the low number of lung biopsy samples available from asthmatic patients and, in fact, some of these samples were derived from individuals with other severe lung pathology, which could probably explain the reason why no significant differences were observed. Due to the limited number of analyzed lung biopsy samples available from the asthmatic patients, it is a preliminary level study in this case. Finally, it would be very interesting to evaluate the levels of protein expression by ELISA of the genes regulated by two of the altered miRNAs (hsa-miR-148b-3p and hsa-miR-221-5p) in this study.

4. Materials and Methods

4.1. Study Subjects and Sample Collection

Serum samples were obtained from 46 subjects with a severe asthma diagnosis who were recruited from allergy and pulmonology units of a series of hospitals in Spain; from these, we selected samples from 20 patients undergoing treatment with OCS and 26 not receiving OCS. MiRNA-sequencing (miRNA-Seq) was performed on serum samples from six of the OCS-treated and seven non-OCS-treated patients, respectively. Of the remaining 33 samples, 14 of those obtained from OCS-treated and 19 from non-OCS-treated subjects were included for miRNA validation by means of quantitative real-time polymerase chain reaction (qPCR). Descriptive data representing the demographic, inflammatory, functional, and clinical characteristics of the study subjects were collected.

All patients were participants in the MEGA project, which studies a cohort of asthmatic individuals with varying grades of severity [37]. The inclusion criteria were as follows: (i) acceptance to participate by providing signed informed consent; (ii) asthma diagnosis following the 2021 GINA criteria [38]; and (iii) age between 18 and 75 years. The definition of OCS-treated and non-OCS-treated severe asthmatic was established based on their treatment.

The study was conducted in accordance with the tenets of the Declaration of Helsinki, and the protocol was approved by the participating hospital ethics committees.

In addition, we analyzed the lung biopsy samples and relevant descriptive data corresponding to 13 asthmatic patient donors; these were also subdivided depending on whether the individuals were treated with OCS (OCS-treated) or not (non-OCS-treated). Furthermore, we had 20 lung biopsy samples from the control subjects, required for the gene expression analysis in lung tissue. Lung biopsy samples from these asthmatic patients and control individuals were provided by the CIBERES Pulmonary Biobank Consortium (PT13/0010/0030), a network currently formed by 12 tertiary Spanish hospitals (www.ciberes.org, accessed on 10 January 2023) detailed in the Acknowledgements section, and integrated in the Spanish National Biobanks Network. Lung biopsies were processed in accordance with the standard operating procedures, and subsequent approval was granted for processing by the local ethics and scientific committees.

Serum samples were obtained by blood clotting in anticoagulant-free tubes and subsequent centrifugation at 3000 rpm for 10 min at 4 °C and stored at −80 °C until use. Lung tissue samples were preserved in the RNAlater stabilization solution at −80 °C until use.

4.2. Library Preparation, MiRNAs-Seq, Bioinformatic Analysis

Serum RNA (including miRNAs) was extracted from 200 μL of serum using the miRNeasy Serum/Plasma Advanced Kit (Qiagen, Hilden, Germany), as indicated by the manufacturer. Lung tissue RNA was purified using the QIAzol Lysis Reagent (Qiagen, Hilden, Germany), followed by the application of the acid guanidinium thiocyanate-phenol-chloroform extraction method [39]. In all cases, the RNA enriched in miRNAs was eluted by adding 20 μL of RNase-free water.

Then, small RNA (miRNA-enriched RNA) was converted to Illumina sequencing libraries using the NEXTFLEX^® Small RNA-Seq Kit v3 (Bioo Scientific Corporation, Austin, TX, USA), strictly according to the manufacturer’s instructions. The size profile of the individual libraries was quantified using High Sensitivity D1000 Screentape on a 4200 TapeStation System (both Agilent, Santa Clara, CA, USA). Quantified libraries were sequenced on an Illumina MiniSeq 550 platform (Illumina, San Diego, CA, USA) using a MiniSeq 500/550 75-Cycle High Output Kit.

Bioinformatic analyses including quality control, pre-processing, and statistical analysis of small RNA-Seq data were carried out by the Bioinformatics Unit of IIS-Fundación Jiménez Díaz. Adaptor removal and trimming of raw reads were performed using Cutadapt [40] by following the instructions for the NEXTflex small RNA-Seq Kit. Adaptor-trimmed reads between 17 and 25 nt were retained and aligned to the reference genome (GRCh38 assembly) using Bowtie2 [41]. Mapping of reads to known miRNAs was performed with HTSeq-count2 [42] using mature miRNA annotation retrieved from the miRBase database (miRBase v22). Raw miRNA counts across samples were normalized by sequencing depth and RNA composition using the trimmed mean of M-values (TMM) function from the NOISeq Bioconductor R package [43]. Subsequent principal component analysis (PCA) on normalized and scaled values was applied using the prcomp function from the R stats package [44]. Comparison of the normalized expression levels across groups was performed following two alternative methods for testing differential expression in the sequencing data: NOISeq [43] and DESeq2 [45]. Fold-change and p-values adjusted by false discovery rate (FDR) were calculated and used to identify significant differentially expressed miRNAs. MiRNAs were considered biologically relevant if they were differentially expressed (adjusted p < 0.05) and presented a Log2 fold change ≥1.5 between two groups, or if q ≥ 0.8 in the case of NOISeq.

4.3. RT-qPCR and Pathway Enrichment Analyses

For validation of the miRNA-Seq results or the expression analysis of the validated miRNAs, 4 µL of serum miRNA samples or 30 ng of lung tissue miRNA-enriched RNA were reverse transcribed to cDNA using the miRCURY LNA RT Kit (Qiagen), following the manufacturer’s protocol. The synthetic miRNAs SP6 and cel-miR-39-3p were used to control for correct reverse transcription. The cDNA obtained was stored at −20 °C until use. Then, miRNA expression was evaluated by qPCR using the miRCURY LNA SYBR Green PCR Kit (Qiagen), as indicated in the instructions. For this purpose, we used 3 μL of cDNA from the serum or lung tissue miRNAs diluted 1:60 or 1:30, respectively, in RNase-free water. Based on the results of the miRNA-Seq, the probes (Qiagen) used for the validation analysis of miRNAs in the serum and lung tissue were the following: hsa-miR-148b-3p, hsa-miR-221-5p, hsa-miR-618, hsa-miR-200b-3p, hsa-miR-548I, hsa-miR-200a-3p, hsa-miR-941, hsa-miR-181d-5p, hsa-miR-769-5p, hsa-miR-133a-3p, hsa-miR-451a, hsa-miR-144-3p, and hsa-miR-144-5p. Additionally, hsa-miR-103a-3p, hsa-miR-191-5p, SP6, cel-miR-39-3p, and U6 (Qiagen) were selected as housekeeping miRNAs. All samples were run in triplicate, and reactions were performed in a Light Cycler^® 96 thermocycler (Roche, Basel, Switzerland). Cycle threshold (Ct) values were analyzed with LightCycler^® 96 SW 1.1 (Roche) software and miRNA expression was calculated by normalizing to the endogenous miRNA controls by applying the 2^−ΔCt method [46], where ΔCt = Ct_miRNA − Ct_{housekeeping miRNAs}. The relative quantification of differences in expression (RQ = 2^−ΔΔCt; where ΔΔCt = ΔCt_OCS-treated − ΔCt_{non-OCS-treated}) was carried out by the ΔΔCt method [46].

For gene expression analysis in lung tissue, 500 ng of RNA quantified by a Nanodrop ND-1000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA) was reverse transcribed using the High-Capacity cDNA Reverse Transcription Kit (Applied Biosystems, Foster City, CA, USA), followed by qPCR according to the manufacturer’s guidelines on a StepOne Real-Time PCR System (Applied Biosystems). TaqMan^TM gene expression probes were purchased for TGFBR1, FOXO3, PTEN, MAPK3, and GAPDH using TaqMan^TM Gene Expression MasterMix (Applied Biosystems). Gene expression was calculated by normalizing to the endogenous gene GAPDH control by applying the 2^−ΔCt method as previously reported (RQ = 2^−ΔΔCt; where ΔΔCt = ΔCt_asthmatics − ΔCt_{healthy subjects}) [46].

In order to identify target genes linked to asthmatic pathology of hsa-miR-148b-3p, hsa-miR-221-5p, hsa-miR-618, hsa-miR-941, and hsa-miR-769-5p, differentially expressed miRNAs in the serum between the OCS-treated and non-OCS-treated patients with severe asthma, pathway enrichment analysis of dysregulated miRNAs was performed using the DIANA-miRPath v3.0 bioinformatic online resource [47]. Relevant pathways for asthma disease were represented when a p-value < 0.05, and the related genes were analyzed by RT-PCR, as previously mentioned.

4.4. Statistical Analysis

Statistical analyses and graphs were created with GraphPad Prism^® v6-8.0 (GraphPad Software Inc., San Diego, CA, USA) and R software^® v4.1.0 [44] (R Foundation for Statistical Computing, Vienna, Austria).

Results are shown as the mean (standard deviation, SD) or median (interquartile range, IQR) values. For all statistical analyses, differences showing p < 0.05 were considered significant. Normality was analyzed by means of the Shapiro–Wilk test. For continuous variables, comparisons of normally distributed data between non-paired groups were performed via an unpaired t test (equal SD) and t test with Welch’s correction (different SD), and non-normally distributed data and non-paired groups were compared by the Mann–Whitney test.

Additionally, correlations between the miRNA expression levels (ΔCt) and some clinical parameters (quantitative variables) were estimated by Spearman (non-normally distributed data) or Pearson (normally distributed data) correlation; the Fisher exact test was applied to a 2 × 2 contingency table to test the null hypothesis of the independence of groups and some clinical characteristics (qualitative variables). Finally, the expression profile (ΔCt) of each differentially expressed miRNA was used to create receiver operator characteristic (ROC) curves, and logistic regression models were developed to evaluate the performance of miRNA as biomarkers; an area under the curve (AUC) of 0.7 indicated an acceptable biomarker.

5. Conclusions

In summary, we describe the significant differences in the expression of eight miRNAs, hsa-miR-148b-3p, hsa-miR-221-5p, hsa-miR-618, hsa-miR-941, hsa-miR-769-5p, hsa-miR-144-3p, hsa-miR-144-5p, and hsa-miR-451a (the first five in serum and the last three in lung tissue). These miRNAs could be used as biomarkers of oral corticosteroid treatment, allowing for the differentiation of patients treated with OCS from the patients not treated with OCS.

Author Contributions

V.d.P. conceived the manuscript and designed the study. M.G.-M., C.L.-S., J.M.R.-M. and M.Á.L. performed the experiments and data analysis. M.G.-M. and C.L.-S. edited and wrote the manuscript. G.N.-M., L.d.l.F., P.M. and I.M.-F. performed the bioinformatics and statistical data analysis. J.S., M.V.-M., S.Q., M.L.C., F.J.G.-B., E.A., I.B., A.V., X.M., M.J.C., C.M.-R., V.P. and J.M.O. recruited the patients and obtained the samples and clinical data. All authors contributed to the article and approved the submitted version. All authors have read and agreed to the published version of the manuscript.

Institutional Review Board Statement

The study protocol was approved by the Ethics Committee of each different hospital, and the study was conducted in accordance with the principles set forth in the Declaration of Helsinki.

Informed Consent Statement

Signed informed consent was obtained from all patients. The study was conducted in accordance with the tenets of the Declaration of Helsinki, and the protocol was approved by the participating hospital ethics committees.

Data Availability Statement

The data that support the findings of this study are available from the corresponding author, V.d.P., upon reasonable request.

Acknowledgments

We wish to particularly acknowledge the patients and the CIBERES Pulmonary Biobank Consortium (PT13/0010/0030) on behalf of Hospital Universitario 12 de Octubre, Hospital del Mar-IMIM, Hospital Universitario de Getafe, Consorcio Universitario Hospital General de Valencia, Fundación Jiménez-Díaz, Hospital Clínic de Barcelona, Hospital Universitario de Bellvitge, Hospital Universitario Joan XXIII, Hospital Universitario Germans Trias i Pujol, Hospital Universitario Virgen del Rocío, Hospital Vall d’Hebron, and Hospital Universitario Son Espases, integrated in the Spanish National Biobanks Network, for their collaboration. The authors also acknowledge Oliver Shaw for the editorial support.

Conflicts of Interest

Authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as potential conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. miRNA expression profile in the serum samples from asthmatic patients. (a) Heatmap of differentially expressed serum miRNAs between the OCS-treated and non-OCS-treated patients with severe asthma (adjusted p < 0.05; Log2 fold change ≥1.5). (b) Principal component analysis (PCA) of severe asthmatic serum samples showing the two different treatment groups (OCS-treated vs. non-OCS-treated severe asthmatics). (c) RT-qPCR validation of the differentially expressed serum miRNAs between OCS-treated and non-OCS-treated individuals with severe asthma. 2−∆∆Ct values are represented in the graph (**, p < 0.01; *, p < 0.05).

Figure 2. Correlation analyses and the ROC curves of differentially expressed serum miRNAs. (a) Table showing the correlation coefficient (r) between the distinct clinical parameters and serum miRNA expression levels. Red color intensity indicates a higher negative correlation, whereas green color intensity shows higher positive correlations. Bold and underlined values were statistically significant (p < 0.05). (b) ROC curves and multivariate logistic regression models (continuous predictors) of serum hsa-miR-221-5p and hsa-miR-769-5p. Bold values indicate an AUC over 0.70.

Figure 3. In silico pathway analysis, gene expression, and correlation data of serum miRNA results. (a) Graph showing enrichment (as −log p value) of the signaling pathways obtained using the DIANA-miRPath v3.0 bioinformatic tool. (b) RT-qPCR in the lung biopsies of target genes of the differentially expressed serum miRNAs (hsa-miR-148b-3p and hsa-miR-221-5p). GAPDH was used as a housekeeping gene (*, p < 0.05). (c) The relative miRNA expression of hsa-miR-221-5p (ΔCt) correlated directly with MAPK3 in the lung tissue samples (p < 0.05). (d) Correlation graphs of hsa-miR-221-5p in the lung biopsy samples with blood eosinophil percentage (%) in different population groups (asthmatic individuals and total subjects) (p < 0.05). Dots (red and green) are individual values, black dashed line indicates error bars and red/green solid line shows linear regression.

Figure 4. miRNA expression profile in the lung tissue biopsies from asthmatic individuals. Graph showing the mean NOIseq values obtained after miRNAs-Seq bioinformatic analyses (*, q ≥ 0.8).

Table 1

Demographic, inflammatory, functional, and clinical characteristics of the study subjects (serum sample analysis).

		OCS-Treated(n = 20)	Non-OCS-Treated(n = 26)	p-Value
Age (years) †		54.5 (44.3–59.5)	53.0 (42.8–61.3)	N.S.
Sex (%)	Female	17 (85.0%)	23 (88.5%)	N.S.
BMI †		26.7 (23.1–30.4)	28.7 (24.2–33.4)	N.S.
Smoking habit (%)	Smokers	1 (5.0%)	4/25 (16.0%)	N.S.
	Passive	2 (10.0%)	2/25 (8.0%)	N.S.
	Ex-smokers	7 (35.0%)	9/25 (36.0%)	N.S.
	Non-smokers	10 (50.0%)	10/25 (40.0%)	N.S.
Blood eosinophils (cells/µL) †		300.0 (100.0–500.0)	300.0 (200.0–525.0)	N.S.
Sputum eosinophils (%) †		0.0% (0.0–35.0)	3.1% (2.0–22.1)	N.S.
Atopy (%)		12 (60.0%)	17 (65.4%)	N.S.
IgE (IU) †		366.0 (110.0–690.0)	96.8 (36.6–285.8)	*
FEV1/FVC (%) #		67.4% (±15.9)	73.3% (±23.5)	N.S.
FeNO (ppb) †		23.0 (13.5–78.3)	32.0 (13.8–62.0)	N.S.
Exacerbations (%)		17 (85.0%)	16 (61.5%)	N.S.
Severity (%)	Severe	20 (100.0%)	26 (100.0%)	N.S.
	Moderate	0 (0.0%)	0 (0.0%)	N.S.
	Mild	0 (0.0%)	0 (0.0%)	N.S.
	Intermittent	0 (0.0%)	0 (0.0%)	N.S.
ACT #		16.8 (±5.9)	17.2 (±5.5)	N.S.
ICS and LABA (%)		20 (100.0%)	26 (100.0%)	N.S.
Systemic corticosteroid (%)		20 (100.0%)	0 (0.0%)	****

Results are expressed as ^# mean (± SD) or ^† median (IQR); N.S., non-significant; ****, p < 0.0001; *, p < 0.05; BMI, body mass index; FEV₁, forced expiratory volume measured during the first second; FVC, forced vital capacity; FeNO, fractional exhaled nitric oxide; ppb, parts per billion; ACT, asthma control test; ICS and LABA, inhaled corticosteroids and long-acting β2-agonists; OCS, oral corticosteroids.

Table 2

Demographic, inflammatory, and clinical characteristics of the study subjects (lung biopsy sample analysis).

		OCS-Treated(n = 6)	Non-OCS-Treated(n = 7)	p-Value
Age (years) †		63.0 (56.8–66.3)	34.0 (27.0–74.0)	N.S.
Sex (%)	Female	3 (50.0%)	1 (14.3%)	N.S.
Smoking habit (%)	Smokers	1 (16.7%)	1/6 (16.7%)	N.S.
	Ex-smokers	3 (50.0%)	4/6 (66.7%)	N.S.
	Non-smokers	2 (33.3%)	1/6 (16.7%)	N.S.
Neutrophils (%) †		61.9% (50.3–72.5)	83.6% (70.6–89.5)	N.S.
Lymphocytes (%) †		24.1% (16.5–35.5)	8.7% (5.9–16.5)	N.S.
Monocytes (%) †		10.0% (8.9–12.2)	5.9% (4.4–9.0)	*
Eosinophils (%) †		2.1% (0.3–3.0)	0.5% (0.1–3.1)	N.S.
Basophils (%) †		0.4% (0.3–0.7)	0.2% (0.1–0.5)	N.S.
Atopy (%)		1 (16.7%)	1 (14.3%)	N.S.
OCS (%)		6 (100.0%)	0 (0.0%)	**

Results are expressed as ^† median (IQR); N.S., non-significant; **, p < 0.01; *, p < 0.05; OCS, oral corticosteroids.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/ijms24021611/s1.

References

1. James, S.L.; Abate, D.; Abate, K.H.; Abay, S.M.; Abbafati, C.; Abbasi, N.; Abbastabar, H.; Abd-Allah, F.; Abdela, J.; Abdelalim, A. et al. Global, regional, and national incidence, prevalence, and years lived with disability for 354 diseases and injuries for 195 countries and territories, 1990–2017: A systematic analysis for the Global Burden of Disease Study 2017. Lancet; 2018; 392, pp. 1789-1858. [DOI: https://dx.doi.org/10.1016/S0140-6736(18)32279-7] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30496104]

2. Hough, K.P.; Curtiss, M.L.; Blain, T.J.; Liu, R.M.; Trevor, J.; Deshane, J.S.; Thannickal, V.J. Airway Remodeling in Asthma. Front. Med.; 2020; 7, 191. [DOI: https://dx.doi.org/10.3389/fmed.2020.00191] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32509793]

3. Chung, K.F.; Wenzel, S.E.; Brozek, J.L.; Bush, A.; Castro, M.; Sterk, P.J.; Adcock, I.M.; Bateman, E.D.; Bel, E.H.; Bleecker, E.R. et al. International ERS/ATS guidelines on definition, evaluation and treatment of severe asthma task force report ERS/ATS guidelines on severe asthma Executive Summary. Eur. Respir. J.; 2014; 43, pp. 343-373. [DOI: https://dx.doi.org/10.1183/09031936.00202013] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/24337046]

4. Côté, A.; Godbout, K.; Boulet, L.-P. The management of severe asthma in 2020. Biochem. Pharmacol.; 2020; 179, 114112. [DOI: https://dx.doi.org/10.1016/j.bcp.2020.114112] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32598948]

5. Agache, I.; Akdis, C.A. Precision medicine and phenotypes, endotypes, genotypes, regiotypes, and theratypes of allergic diseases. J. Clin. Investig.; 2019; 129, pp. 1493-1503. [DOI: https://dx.doi.org/10.1172/JCI124611] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30855278]

6. Sheikhpour, M.; Maleki, M.; Ebrahimi Vargoorani, M.; Amiri, V. A review of epigenetic changes in asthma: Methylation and acetylation. Clin. Epigenet.; 2021; 13, 65. [DOI: https://dx.doi.org/10.1186/s13148-021-01049-x] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/33781317]

7. Kuruvilla, M.E.; Lee, F.E.H.; Lee, G.B. Understanding Asthma Phenotypes, Endotypes, and Mechanisms of Disease. Clin. Rev. Allergy Immunol.; 2019; 56, pp. 219-233. [DOI: https://dx.doi.org/10.1007/s12016-018-8712-1] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30206782]

8. Busse, W.W. Biological treatments for severe asthma: A major advance in asthma care. Allergol. Int.; 2019; 68, pp. 158-166. [DOI: https://dx.doi.org/10.1016/j.alit.2019.01.004] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30792118]

9. Suehs, C.M.; Menzies-Gow, A.; Price, D.; Bleecker, E.R.; Canonica, G.W.; Gurnell, M.; Bourdin, A. Expert consensus on the tapering of oral corticosteroids for the treatment of asthma: A delphi study. Am. J. Respir. Crit. Care Med.; 2021; 203, pp. 871-881. [DOI: https://dx.doi.org/10.1164/rccm.202007-2721OC]

10. Li, H.; Zhang, Q.; Wang, J.; Gao, S.; Li, C.; Wang, J.; Zhang, S.; Lin, J. Variability of Type 2 inflammatory markers guiding biologic therapy of severe asthma: A 5-year retrospective study from a single tertiary hospital. World Allergy Organ. J.; 2021; 14, 100547. [DOI: https://dx.doi.org/10.1016/j.waojou.2021.100547]

11. Weidner, J.; Bartel, S.; Kılıç, A.; Zissler, U.M.; Renz, H.; Schwarze, J.; Schmidt-Weber, C.B.; Maes, T.; Rebane, A.; Krauss-Etschmann, S. et al. Spotlight on microRNAs in allergy and asthma. Allergy Eur. J. Allergy Clin. Immunol.; 2021; 76, pp. 1661-1678. [DOI: https://dx.doi.org/10.1111/all.14646] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/33128813]

12. Ying, S.Y.; Chang, D.C.; Lin, S.L. The MicroRNA (miRNA): Overview of the RNA Genes that Modulate Gene Function. Mol. Biotechnol.; 2008; 38, pp. 257-268. [DOI: https://dx.doi.org/10.1007/s12033-007-9013-8] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/17999201]

13. Dexheimer, P.J.; Cochella, L. MicroRNAs: From Mechanism to Organism. Front. Cell Dev. Biol.; 2020; 8, 409. [DOI: https://dx.doi.org/10.3389/fcell.2020.00409] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32582699]

14. Mori, M.A.; Ludwig, R.G.; Garcia-Martin, R.; Brandão, B.B.; Kahn, C.R. Extracellular miRNAs: From Biomarkers to Mediators of Physiology and Disease. Cell Metab.; 2019; 30, pp. 656-673. [DOI: https://dx.doi.org/10.1016/j.cmet.2019.07.011] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/31447320]

15. Gil-Martínez, M.; Rodrigo-Muñoz, J.; Sastre, B.; Cañas, J.; García-Latorre, R.; Redondo, N.; de la Fuente, L.; Mínguez, P.; Mahíllo-Fernández, I.; Sastre, J. et al. Serum microRNAs catalog asthmatic patients by phenotype. J. Investig. Allergol. Clin. Immunol.; 2021; 32, pp. 471-478. [DOI: https://dx.doi.org/10.18176/jiaci.0753] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34489228]

16. Popović-Grle, S.; Štajduhar, A.; Lampalo, M.; Rnjak, D. Biomarkers in Different Asthma Phenotypes. Genes; 2021; 12, 801. [DOI: https://dx.doi.org/10.3390/genes12060801] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34070316]

17. Mirra, D.; Cione, E.; Spaziano, G.; Esposito, R.; Sorgenti, M.; Granato, E.; Cerqua, I.; Muraca, L.; Iovino, P.; Gallelli, L. et al. Circulating MicroRNAs Expression Profile in Lung Inflammation: A Preliminary Study. J. Clin. Med.; 2022; 11, 5446. [DOI: https://dx.doi.org/10.3390/jcm11185446]

18. Carr, T.F.; Bleecker, E. Asthma heterogeneity and severity. World Allergy Organ. J.; 2016; 9, 41. [DOI: https://dx.doi.org/10.1186/s40413-016-0131-2]

19. Licari, A.; Castagnoli, R.; Brambilla, I.; Marseglia, A.; Tosca, M.A.; Marseglia, G.L.; Ciprandi, G. Asthma endotyping and biomarkers in childhood asthma. Pediatr. Allergy Immunol. Pulmonol.; 2018; 31, pp. 44-55. [DOI: https://dx.doi.org/10.1089/ped.2018.0886]

20. Xtara, X.; Carr, F.; Xmd, X.; Kraft, X. Use of biomarkers to identify phenotypes and endotypes of severe asthma. Ann. Allergy Asthma Immunol.; 2018; 121, pp. 414-420. [DOI: https://dx.doi.org/10.1016/j.anai.2018.07.029]

21. Tiotiu, A. Biomarkers in asthma: State of the art. Asthma Res. Pract.; 2018; 4, 10. [DOI: https://dx.doi.org/10.1186/s40733-018-0047-4] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30598830]

22. Canonica, G.W.; Blasi, F.; Crimi, N.; Paggiaro, P.; Papi, A.; Fanelli, F.; Stassaldi, A.; Furneri, G. Defining type 2 asthma and patients eligible for dupilumab in Italy: A biomarker-based analysis. Clin. Mol. Allergy; 2021; 19, 5. [DOI: https://dx.doi.org/10.1186/s12948-021-00146-9] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34020658]

23. Zieg, G.; Lack, G.; Harbeck, R.J.; Gelfand, E.W.; Leung, D.Y. In vivo effects of glucocorticoids on IgE production. J. Allergy Clin. Immunol.; 1994; 94, pp. 222-230. [DOI: https://dx.doi.org/10.1053/ai.1994.v94.a54936] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/8064074]

24. FitzGerald, J.M.; Lemiere, C.; Lougheed, M.D.; Ducharme, F.M.; Dell, S.D.; Ramsey, C.; Yang, M.C.L.; Côté, A.; Watson, W.; Olivenstein, R. et al. Recognition and management of severe asthma: A Canadian Thoracic Society position statement. Can. J. Respir. Crit. Care Sleep Med.; 2017; 1, pp. 199-221. [DOI: https://dx.doi.org/10.1080/24745332.2017.1395250]

25. Wu, Y.; Li, Q.; Zhang, R.; Dai, X.; Chen, W.; Xing, D. Circulating microRNAs: Biomarkers of disease. Clin. Chim. Acta; 2021; 516, pp. 46-54. [DOI: https://dx.doi.org/10.1016/j.cca.2021.01.008]

26. Potla, P.; Ali, S.A.; Kapoor, M. A bioinformatics approach to microRNA-sequencing analysis. Osteoarthr. Cartil. Open; 2021; 3, 100131. [DOI: https://dx.doi.org/10.1016/j.ocarto.2020.100131]

27. Rodrigo-Muñoz, J.M.; Gil-Martínez, M.; Lorente-Sorolla, C.; García-Latorre, R.; Valverde-Monge, M.; Quirce, S.; Sastre, J.; del Pozo, V. miR-144-3p Is a Biomarker Related to Severe Corticosteroid-Dependent Asthma. Front. Immunol.; 2022; 13, 858722. [DOI: https://dx.doi.org/10.3389/fimmu.2022.858722]

28. Ong, J.; Faiz, A.; Timens, W.; van den Berge, M.; Terpstra, M.M.; Kok, K.; van den Berg, A.; Kluiver, J.; Brandsma, C.A. Marked TGF-β-regulated miRNA expression changes in both COPD and control lung fibroblasts. Sci. Rep.; 2019; 9, 18214. [DOI: https://dx.doi.org/10.1038/s41598-019-54728-4]

29. Guan, Y.; Ma, Y.; Tang, Y.; Liu, X.; Zhao, Y.; An, L. MiRNA-221-5p suppressed the Th17/Treg ratio in asthma via RORγt/Foxp3 by targeting SOCS1. Allergy Asthma Clin. Immunol.; 2021; 17, 123. [DOI: https://dx.doi.org/10.1186/s13223-021-00620-8]

30. Li, C.; Hu, X.; Li, L.; Li, J. hui Differential microRNA expression in the peripheral blood from human patients with COVID-19. J. Clin. Lab. Anal.; 2020; 34, e23590. [DOI: https://dx.doi.org/10.1002/jcla.23590] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32960473]

31. Liu, P.F.; Yan, P.; Zhao, D.H.; Shi, W.F.; Meng, S.; Liu, Y.; Liu, B.; Li, G.F.; Xie, L.X. The effect of environmental factors on the differential expression of miRNAs in patients with chronic obstructive pulmonary disease: A pilot clinical study. Int. J. Chron. Obstruct. Pulmon. Dis.; 2018; 13, pp. 741-751. [DOI: https://dx.doi.org/10.2147/COPD.S156865] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/29535514]

32. Rodrigo-Muñoz, J.M.; Cañas, J.A.; Sastre, B.; Rego, N.; Greif, G.; Rial, M.; Mínguez, P.; Mahíllo-Fernández, I.; Fernández-Nieto, M.; Mora, I. et al. Asthma diagnosis using integrated analysis of eosinophil microRNAs. Allergy Eur. J. Allergy Clin. Immunol.; 2019; 74, pp. 507-517. [DOI: https://dx.doi.org/10.1111/all.13570] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30040124]

33. Jovičić, A.; Roshan, R.; Moisoi, N.; Pradervand, S.; Moser, R.; Pillai, B.; Luthi-Carter, R. Comprehensive Expression Analyses of Neural Cell-Type-Specific miRNAs Identify New Determinants of the Specification and Maintenance of Neuronal Phenotypes. J. Neurosci.; 2013; 33, pp. 5127-5137. [DOI: https://dx.doi.org/10.1523/JNEUROSCI.0600-12.2013] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/23516279]

34. Imraish, A.; Abu-Thiab, T.; Zihlif, M. IL-13 and FOXO3 genes polymorphisms regulate IgE levels in asthmatic patients. Biomed. Rep.; 2021; 14, 55. [DOI: https://dx.doi.org/10.3892/br.2021.1431] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/33884198]

35. Wen, X.; Yan, J.; Han, X.R.; Zheng, G.H.; Tang, R.; Liu, L.F.; Wu, D.M.; Lu, J.; Zheng, Y.L. PTEN gene silencing contributes to airway remodeling and induces airway smooth muscle cell proliferation in mice with allergic asthma. J. Thorac. Dis.; 2018; 10, pp. 202-211. [DOI: https://dx.doi.org/10.21037/jtd.2017.12.104]

36. Alam, R.; Gorska, M.M. MAPK Signaling and ERK1/2 bistability in Asthma. Clin. Exp. Allergy; 2011; 41, pp. 149-159. [DOI: https://dx.doi.org/10.1111/j.1365-2222.2010.03658.x]

37. Muñoz, X.; Álvarez-Puebla, M.J.; Arismendi, E.; Arochena, L.; del Pilar Ausín, M.; Barranco, P.; Bobolea, I.; Cañas, J.A.; Cardaba, B.; Crespo, A. et al. Estudio de los mecanismos implicados en la génesis y evolución del asma (proyecto MEGA): Creación y seguimiento a largo plazo de una cohorte de pacientes asmáticos. Arch. Bronconeumol.; 2018; 54, pp. 378-385. [DOI: https://dx.doi.org/10.1016/j.arbres.2017.12.012]

38. Global Initiative for Asthma-Global Initiative for Asthma-GINA. Available online: https://ginasthma.org/ (accessed on 10 May 2021).

39. Chomczynski, P.; Sacchi, N. Single-step method of RNA isolation by acid guanidinium thiocyanate-phenol-chloroform extraction. Anal. Biochem.; 1987; 162, pp. 156-159. [DOI: https://dx.doi.org/10.1016/0003-2697(87)90021-2]

40. Martin, M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet. J.; 2011; 17, pp. 10-12. [DOI: https://dx.doi.org/10.14806/ej.17.1.200]

41. Langmead, B.; Salzberg, S.L. Fast gapped-read alignment with Bowtie 2. Nat. Methods; 2012; 9, pp. 357-359. [DOI: https://dx.doi.org/10.1038/nmeth.1923]

42. Anders, S.; Pyl, P.T.; Huber, W. HTSeq—A Python framework to work with high-throughput sequencing data. Bioinformatics; 2015; 31, pp. 166-169. [DOI: https://dx.doi.org/10.1093/bioinformatics/btu638] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/25260700]

43. Tarazona, S.; García-Alcalde, F.; Dopazo, J.; Ferrer, A.; Conesa, A. Differential expression in RNA-seq: A matter of depth. Genome Res.; 2011; 21, pp. 2213-2223. [DOI: https://dx.doi.org/10.1101/gr.124321.111] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/21903743]

44. R Core Team. R: A Language and Environment for Statistical Computing; R Foundation for Statistical Computing: Vienna, Austria, 2013; Available online: http://www.R-project.org/ (accessed on 10 January 2023).

45. Love, M.I.; Huber, W.; Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol.; 2014; 15, 550. [DOI: https://dx.doi.org/10.1186/s13059-014-0550-8]

46. Livak, K.J.; Schmittgen, T.D. Analysis of relative gene expression data using real-time quantitative PCR and the 2^−ΔΔCT method. Methods; 2001; 25, pp. 402-408. [DOI: https://dx.doi.org/10.1006/meth.2001.1262] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/11846609]

47. Vlachos, I.S.; Zagganas, K.; Paraskevopoulou, M.D.; Georgakilas, G.; Karagkouni, D.; Vergoulis, T.; Dalamagas, T.; Hatzigeorgiou, A.G. DIANA-miRPath v3.0: Deciphering microRNA function with experimental support. Nucleic Acids Res.; 2015; 43, pp. W460-W466. [DOI: https://dx.doi.org/10.1093/nar/gkv403] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/25977294]