Immunoglobulin E-virus phenotypes of infant bronchiolitis and risk of childhood asthma

Background: Bronchiolitis is the leading cause of infant hospitalization in U.S. and is associated with increased risk for childhood asthma. Immunoglobulin E (IgE) not only plays major roles in antiviral immune responses and atopic predisposition, but also offers a potential therapeutic target.

Objective: We aimed to identify phenotypes of infant bronchiolitis by using total IgE (tIgE) and virus data, to determine their association with asthma development, and examine their biological characteristics.

Methods: In a multicenter prospective cohort study of 1,016 infants (age <1 year) hospitalized for bronchiolitis, we applied clustering approaches to identify phenotypes by integrating tIgE and virus (respiratory syncytial virus [RSV], rhinovirus [RV]) data at hospitalization. We examined their longitudinal association with the risk of developing asthma by age 6 years and investigated their biological characteristics by integrating the upper airway mRNA and microRNA data in a subset (n=182).

Results: In infants hospitalized for bronchiolitis, we identified 4 phenotypes: 1) tIgE^lowvirus^RSV-high, 2) tIgE^lowvirus^RSV-low/RV, 3) tIgE^highvirus^RSV-high, and 4) tIgE^highvirus^RSV-low/RV phenotypes. Compared to phenotype 1 infants (resembling “classic” bronchiolitis), phenotype 4 infants (tIgE^highvirus^RSV-low/RV) had a significantly higher risk for developing asthma (19% vs. 43%; adjOR, 2.93; 95% CI, 1.02–8.43; P=.046). Phenotypes 3 and 4 (tIgE^high) had depleted type I interferon and enriched antigen presentation pathways; phenotype 4 also had depleted airway epithelium structure pathways.

Conclusions: In this multicenter cohort, tIgE-virus clustering identified distinct phenotypes of infant bronchiolitis with differential risks of asthma development and unique biological characteristics.

Introduction

Bronchiolitis is the leading cause of infant hospitalization in the U.S., accounting for 110,000 hospitalizations annually (1). Bronchiolitis also puts infants at increased risk for subsequent morbidity—~30% of infants hospitalized for bronchiolitis develop asthma (2). While bronchiolitis has traditionally been considered a single disease with a similar mechanism (3), growing evidence supports disease heterogeneity (4–6). Indeed, recent studies have reported subgroups (phenotypes and endotypes) of infant bronchiolitis with different risks of asthma by using clinical (7, 8) and molecular (9–11) data. Yet, no subgroup-specific approach for asthma prevention has been developed in this high-risk population.

In childhood asthma, immunoglobulin E (IgE) plays an important role (12) and is a major therapeutic target (13, 14). Research has also shown that, among infants (age <12 months) with parental atopy (15) and hospitalized for bronchiolitis (16), IgE-mediated sensitization is associated with a higher risk of developing asthma. Additionally, considering that IgE-mediated sensitization is uncommon in early infancy (16), a recent study of infant bronchiolitis has demonstrated the relationship of a higher total IgE (tIgE) level with an increased risk of developing asthma (17). Furthermore, research has found that IgE impairs innate antiviral immune responses (e.g., type I interferon) through cross-linking of high-affinity IgE receptor (FcϵRI) (18–20), thereby potentially leading to severe viral infection and subsequent asthma development. Despite the public health and research significance, little is known about the relationship among IgE, viral infection, immune responses in infants with bronchiolitis, and their integrated contributions to the subsequent development of asthma.

To address this knowledge gap, we analyzed data from a multicenter cohort study of infants hospitalized for bronchiolitis to 1) identify tIgE-virus phenotypes, 2) examine their longitudinal associations with asthma development, and 3) determine their biological characteristics by integrating upper airway mRNA and microRNA (miRNA) data. A better understanding of tIgE-driven phenotypes may inform early-life prevention strategies (e.g., anti-IgE antibody) against the development of asthma.

Methods

Study design, setting, and participants

We analyzed data from the 35th Multicenter Airway Research Collaboration (MARC-35) study—a multicenter prospective cohort study. Details of the study design, setting, participants, data collection, testing, and statistical analysis may be found in the Supplementary Methods. Briefly, investigators enrolled infants (<12 months) hospitalized with attending physician-diagnosis of bronchiolitis at 17 sites across 14 U.S. states (Table E1) in 2011–2014. Bronchiolitis was defined by the American Academy of Pediatrics (AAP) guidelines—acute respiratory illness with some combination of rhinitis, cough, tachypnea, wheezing, crackles, and retractions, regardless of previous breathing problem episodes (3). We excluded infants with a known heart-lung disease, immunodeficiency, immunosuppression, or gestational age < 32 weeks. All patients were treated at the discretion of the treating physician. The institutional review board at each participating hospital approved the study with written informed consent obtained from the parent or guardian.

Data collection

Clinical data (demographic characteristics; medical, environmental, and family history; and details of the bronchiolitis course) were collected via structured interview and chart review using a standardized protocol (21). All data were reviewed at the EMNet Coordinating Center (Boston, MA, U.S.), and site investigators were queried about missing data and discrepancies identified by manual data checks. Serum, nasopharyngeal, and nasal specimens were collected using standard protocols (21) within 24 hours of hospitalization. Serum tIgE was measured by using ImmunoCAP Total IgE (Thermo Fisher Scientific, Waltham, MA, USA); the range of this assay is 2 to 5,000 kU/L. Upper airway (nasopharyngeal and nasal) specimens were tested for respiratory viruses (e.g., respiratory syncytial virus [RSV] and rhinovirus [RV]) by real-time polymerase chain reaction assays (16), mRNA profiling by RNA sequencing (RNA-seq), and miRNA profiling by small RNA-seq.

Nasopharyngeal mRNA profiling

The details of RNA extraction, RNA-seq, and quality control are described in Supplementary Methods and previous studies (9). Briefly, after applying total RNA extraction, DNase treatment, and rRNA to randomly-selected 244 specimens, we performed RNA-seq with a NovaSeq6000 (Illumina, San Diego, CA, U.S.) using an S4 100bp PE Flowcell (Illumina). After quality control, all RNA-seq samples had high sequence coverage (mean of 8,067,019 pair-end reads/sample). We estimated transcript abundances in Salmon (22) using the human genome (hg38) and the mapping-based mode.

Nasal miRNA profiling

The details of RNA extraction, small RNA-seq, and quality control are described in Supplementary Methods. Briefly, after total RNA extraction, DNase treatment, and rRNA reduction, we used 624 specimens with sufficient RNA quantity and quality to perform small RNA-seq with a NovaSeq6000 using an S2 50bp PE Flowcell (Illumina). All small RNA-seq samples had sufficient sequence depth (mean of 26,736,490 pair-end reads/sample) to obtain a high sequence coverage. From clean small RNA-seq reads, we estimated miRNA detection and abundance by sMETASeq (23). We mapped trimmed reads against human miRNA sequences from miRBase V22 (24).

Outcome measure

The outcome of interest was the development of asthma by age 6 years. Asthma was defined using a commonly used epidemiologic definition (25): physician diagnosis of asthma, with either asthma medication use (e.g., inhaled bronchodilators and inhaled corticosteroids) or asthma-related symptoms (e.g., wheezing and nocturnal cough) during the year before the evaluation at age 6 years.

Statistical analysis

The analytic workflow is summarized in Figure 1. The details of the statistical analysis can be found in the Supplementary Methods. Briefly, we first identified mutually exclusive clusters for each of the tIgE and virus datasets collected at the index hospitalization. For the tIgE data, we generated the low and high tIgE clusters by using the median value. In contrast to tIgE clustering with the use of median cut-off value, for the virus data (including the genomic load [i.e., cycle threshold value] of RSV and RV), we computed a Gower distance and derived the virus clusters by using a consensus clustering algorithm with partitioning around medoids (PAM) method. To choose an optimal number of the virus clusters, we used a combination of separations of the consensus matrix and relative change of the area under the cumulative distribution function (CDF) curve, in addition to the cluster size and clinical plausibility (16, 26). Second, we combined the tIgE and virus clusters to derive a fused matrix, computed a Gower distance, and identified 4 mutually exclusive phenotypes by using a consensus clustering algorithm with PAM method. Third, to interpret the clinical characteristics of the phenotypes, we developed chord diagrams. Fourth, we determined the longitudinal association of the phenotypes with the asthma risk. Of the 1,016 MARC-35 infants, we used 182 infants with mRNA, miRNA, and asthma outcome data as the analytic cohort. Then, we constructed unadjusted logistic regression and multivariable mixed-effects logistic regression models. In the multivariable mixed-effects model, we adjusted for site effects and patient-level potential confounders (i.e., age, sex, parental history of asthma, prematurity [<37 weeks], pre-hospitalization use of inhaled and/or systemic corticosteroids) based on clinical plausibility and a priori knowledge (27, 28). Fifth, we examined the difference in biological characteristics between the derived phenotypes by using the nasopharyngeal mRNA data. We conducted differential expression gene analysis between the phenotypes and Gene Set Enrichment Analysis (GSEA) (29) based on Biological Processes in Gene Ontology (false discovery rate [FDR] <0.10). Sixth, we also examined the difference in biological characteristics by integrating the mRNA and miRNA data. As with the mRNA data, we conducted differential expression miRNA analysis and miRNA set enrichment analysis based on Biological Processes in Gene Ontology derived from the miRNA-target databases, miRTarBase 8.0 (30). We identified pathways that are enriched in the mRNA data and depleted in the miRNA data or that are depleted in the mRNA data and enriched in the miRNA data (FDR <0.10).

FIGURE 1

Figure 1 Analytic workflow of phenotyping with total immunoglobulin E and virus data. 1. Clustering of individual variables. By using the multicenter prospective cohort (MARC-35) data of 1,016 infants (age <1 year) hospitalized for bronchiolitis, we first identified mutually exclusive clusters for each of the total immunoglobulin E (tIgE) and virus datasets. We generated the low and high tIgE clusters by the median value and identified the 7 virus clusters by computing a Gower distance, followed by applying consensus clustering algorithms to the distance matrices. 2. Clustering of the fused matrix. By integrating these derived clusters from the tIgE and virus datasets, we generated a fused matrix and computed a Gower distance. We identified 4 mutually exclusive tIgE-virus phenotypes by applying a consensus clustering algorithm to the distance matrix. To select an optimal number of profiles, we used a combination of consensus matrix and consensus cumulative distribution function, in addition to phenotype size and clinical plausibility. 3. Clinical interpretation of phenotypes. To interpret the clinical characteristics of the 4 phenotypes, we developed chord diagrams on major clinical and virus factors. 4. Clinical significance of phenotypes. Of MARC-35, 182 infants with nasopharyngeal mRNA and nasal microRNA (miRNA) data at enrollment with the asthma outcome data were selected as the analytic cohort. To examine the clinical significance of the 4 phenotypes, we determined the longitudinal relationship of the phenotypes with the risk of developing asthma by constructing logistic regression models. 5. Biological characteristics of phenotypes. To examine the difference in biological characteristics between the phenotypes, we first conducted a functional pathway analysis by using the nasopharyngeal mRNA data. Next, we conducted a functional pathway analysis by integrating the mRNA data with the nasal miRNA data. RSV, respiratory syncytial virus; RV, rhinovirus.

In the sensitivity analysis, we first computed E-values to determine the robustness of causal inference to potential unmeasured confounding. We next examined the phenotype-outcome associations after excluding infants with a previous history of breathing problems. We also examined the robustness of the phenotype-outcome associations under different definitions of exposure (31) by repeating the analysis using a different number of phenotypes. We conducted statistical analyses using R version 4.1.2 (R Foundation, Vienna, Austria). All P-values were two-tailed, with p < 0.05 considered statistically significant. We corrected for multiple hypothesis testing using the Benjamini-Hochberg FDR with FDR <0.10 considered statistically significant.

Results

Among 1,016 infants hospitalized for bronchiolitis, the median age was 3 (interquartile range 2–6) months, 60% were male, and 42% were non-Hispanic white. Overall, 81% had RSV, 21% had RV, and 16% underwent intensive treatment use (the use of positive pressure ventilation and/or admission to intensive care unit) during the index hospitalization (Table E2).

Integrated clustering of the tIgE and virus data identified distinct phenotypes among infants with bronchiolitis

First, the low and high tIgE clusters by the median value were generated. Second, by applying clustering approaches to the virus data, a 7-class model provided an optimal fit (Figure E1). Next, by integrating these cluster data using the consensus clustering approach, the combination of the consensus matrix and relative change of the area under the CDF curve in addition to phenotype size (n = 225–283) and clinical plausibility found that a 4-class model provided an optimal fit (Figure E2) with the 4 tIgE-virus phenotypes called 1, 2, 3, and 4. The 4 distinct phenotypes were 1) tIgE^lowvirus^RSV-high (23%), 2) tIgE^lowvirus^RSV-low/RV (27%), 3) tIgE^highvirus^RSV-high (22%); and 4) tIgE^highvirus^RSV-low/RV (28%; Table E2; Figure 1).

Descriptively, infants with a phenotype 1 or 2 were similarly characterized by young age at the index hospitalization with a low proportion of previous breathing problems, tIgE level, and proportion of IgE sensitization (Table E2; Figure 2A). Regarding the virus data, infants with a phenotype 1 were characterized by RSV infection (with a high genomic load), while those with a phenotype 2 were characterized by RSV (with a low genomic load) and RV infection. As the phenotype 1 clinically resembled “classic” bronchiolitis (3), this group served as the reference group for the following analyses. In contrast, infants with a phenotype 3 or 4 were similarly characterized by older age at the index hospitalization with a high tIgE level and proportion of IgE sensitization (Table E2; Figures 2B, C). Infants with a phenotype 3 were characterized by a low proportion of previous breathing problems and RSV infection (with a high genomic load), while those with a phenotype 4 were characterized by a high proportion of previous breathing problems and RSV (with a low genomic load) and RV infection.

FIGURE 2

Figure 2 Clinical and virus characteristics of infants hospitalized for bronchiolitis in MARC-35, according to phenotypes. To interpret clinical and virus characteristics of the 4 phenotypes, we constructed chord diagrams that represent the comparison between phenotypes (A) 1 and 2, (B) 1 and 3, and (C) 1 and 4. Ribbons connect each of the phenotypes (phenotypes 1–4) with major clinical and virus characteristics. The width of the ribbon represents the proportion of infants or the mean value within the phenotypes who have the corresponding clinical or virus characteristic, which was scaled to a total of 100% or 1. Of the clinical and virus characteristics, we selected 8 characteristics that may be related to “classic” bronchiolitis and coordinated the presence or increase of those in the same direction. Hence, the phenotype with the narrower bundle of ribbons resembles “classic” bronchiolitis (phenotype 1). CT, cycle threshold; IgE, immunoglobulin E; RSV, respiratory syncytial virus; RV, rhinovirus.

tIgE-virus phenotypes had differential risks for developing asthma

Of the 1,016 infants, the following analyses focused on 182 infants with the nasopharyngeal mRNA and nasal miRNA data at enrollment with the asthma outcome data at age 6 years. The analytic and non-analytic cohorts did not differ in the clinical characteristics (p ≥ 0.05; Table E3), except for the proportion of daycare use, oxygen saturation at the presentation, and proportion of RSV infection. Furthermore, the significantly different characteristics among the 4 phenotypes were almost consistent between MARC-35 and the analytic cohort (Table 1). The identified phenotypes had differential risks for developing asthma. Compared with phenotype 1, phenotype 4 (tIgE^highvirus^RSV-low/RV) had a significantly higher risk of developing asthma (19% vs. 43%; adjusted odds ratio [adjOR] 2.93; 95% confidence interval [CI], 1.02–8.43; p = 0.046; E-value = 2.82; Figure 3). In contrast, the risk of asthma was not significantly different in phenotype 2 or 3 (Figure 3).

TABLE 1

Table 1 Baseline patient characteristics and clinical course of infants hospitalized for bronchiolitis in the analytic cohort, according to 4 phenotypes.

FIGURE 3

Figure 3 Association of phenotypes of infant bronchiolitis with risk of developing childhood asthma. To examine the association of bronchiolitis phenotypes (phenotype 1 as the reference) with the risk of developing asthma, logistic regression models were constructed. ^†The E-value (with its lower 95% confidence interval [CI] bound) represents how strongly unmeasured confounder(s) represents how strongly a set of unmeasured confounders would be associated with the exposure and outcome to fully eliminate the observed association. ^‡Multivariable mixed-effects logistic regression model accounting for patient clustering by site and adjusted for potential confounders (i.e., age, sex, parental history of asthma, prematurity [<37 weeks], previous history of breathing problems, and pre-hospitalization use of inhaled and/or systemic corticosteroids).

tIgE-virus phenotypes had distinct biological characteristics

The tIgE-virus phenotypes of infant bronchiolitis had distinct biological characteristics. With the mRNA data, compared with phenotype 1, phenotype 2, 3, and 4 had different gene expression signatures (Figure E3) and significantly enriched or depleted pathways (FDR < 0.10; Figures 4, E4). Of the identified pathways, phenotypes 3 and 4 (i.e., tIgE^high phenotypes) shared the antigen presentation, type I interferon, and airway epithelium structure pathways (Figures 4A, B). Phenotype 3 (tIgE^highvirus^RSV-high) had enriched antigen presentation pathways (e.g., peptide antigen assembly with major histocompatibility complex [MHC] class II protein complex) and depleted type I interferon pathways (e.g., positive regulation of interferon-beta production; Figures 4A, C). In contrast, phenotype 4 (tIgE^highvirus^RSV-low/RV) had depleted airway epithelium structure pathways (e.g., cilium movement; Figures 4B, D). Likewise, with the miRNA data, phenotypes 2, 3, and 4 had different miRNA expression signatures (Figure E5). By integrating the significantly enriched or depleted miRNA pathways (FDR < 0.10) with the identified mRNA pathways, we determined the enriched or depleted pathways in both data (Figures 5, E6). Similar to the results of mRNA data, phenotypes 3 and 4 shared 4 common pathways (e.g., response to virus and type I interferon signaling pathway; Figures 5A, B). Likewise, phenotype 3 had depleted virus infection response (e.g., defense response to virus; Figure 5C) and type I interferon pathways (e.g., positive regulation of interferon-beta production), while phenotype 4 had depleted airway epithelium structure pathways (e.g., cell junction assembly; Figure 5D).

FIGURE 4

Figure 4 Between-phenotype differences (1 vs. 3 and 4) in nasopharyngeal mRNA pathways among infants hospitalized for bronchiolitis. To examine the difference in the biological characteristics between phenotypes (1 [the reference] vs. 3 and 4 [i.e., tIgE^high phenotypes]), we applied the gene set enrichment analysis based on Biological Processes in Gene Ontology to the nasopharyngeal mRNA data. We identified enriched (normalized enrichment score [NES] ≥ 0 and false discovery rate [FDR] <.10) or depleted (NES <0 and FDR <.10) pathways. Of the identified pathways shared between the phenotype 1 vs. 3 and 1 vs. 4 comparisons, we selected the 20 pathways with the lowest FDR for (A) phenotype 1 vs. 3 and (B) phenotype 1 vs. 4. Likewise, of the identified pathways unique to the phenotype (C) 1 vs. 3 and (D) 1 vs. 4 comparisons, we selected the 20 pathways with the lowest FDR. MHC, major histocompatibility complex; NF-kappaB, nuclear factor κ-light-chain-enhancer of activated B cells.

FIGURE 5

Figure 5 Between-phenotype differences (1 vs. 3 and 4) in nasopharyngeal mRNA pathways and nasal microRNA pathways among infants hospitalized for bronchiolitis. To examine the difference in the biological characteristics between phenotypes (1 [the reference] vs. 3 and 4 [i.e., tIgE^high phenotypes]), we applied the gene set enrichment analysis based on Biological Processes in Gene Ontology to the nasopharyngeal mRNA data and the nasal microRNA (miRNA) data. We identified the pathways “enriched in the mRNA data and depleted in the miRNA data” or “depleted in the mRNA data and enriched in the miRNA data” (false discovery rate [FDR] <.10), accounting for suppressive nature of miRNA to mRNA. The identified pathways shared between phenotype 1 vs. 3 and 1 vs. 4 comparisons are shown in (A, B). Of the identified pathways unique to (C) phenotype 1 vs. 3 and (D) phenotype 1 vs. 4, we selected the 15 pathways with the lowest FDR in the mRNA data.

Sensitivity analysis

In the analysis excluding infants with a previous history of breathing problems with limited statistical power (n=155), the phenotype-asthma associations did not show statistical significance in the unadjusted and adjusted models yet remained consistent. Compared with phenotype 1, phenotype 4 had a non-significantly higher risk of developing asthma (17% vs. 30%; adjOR 2.28; 95% CI, 0.68–7.62; p = 0.18; Figure E7). In phenotyping with 3- and 5-class models, the alluvial plot (Figure E8) demonstrates a consistency of the original phenotypes (phenotypes 1–4) across the different number-class chosen. For example, in a 5-class model, phenotype A, B, C, or D had 100% concordance with the original phenotype 1, 2, 3, or 4, respectively (Figure E9; Table E4). Compared to phenotype A, phenotype E—which is concordant with phenotypes 2 and 4 (i.e., virus^RSV-low/RV phenotypes)—had a significantly higher asthma risk (19% vs. 60%; adjOR, 2.21; 95% CI, 1.27–3.85; p = 0.005; E-value = 2.34; Figure E10).

Discussion

By integrating tIgE and virus data from a multicenter prospective cohort of infants hospitalized for bronchiolitis, we identified 4 clinically distinct phenotypes. For example, compared with the reference phenotype 1, infants with phenotypes 3 and 4 (tIgE^high) had distinct biological characteristics (e.g., depleted type I interferon pathways and enriched antigen presentation pathways) based on upper airway mRNA and miRNA data. Furthermore, infants with phenotype 4 (tIgE^highvirus^RSV-high/RV) had a significantly higher risk of developing asthma and depleted pathways related to airway epithelium structure. The sensitivity analysis revealed the robustness of our findings. To the best of our knowledge, this is the first investigation that has identified IgE-virus phenotypes of infant bronchiolitis and their longitudinal relationship with the risks of developing asthma.

Bronchiolitis has conventionally deemed a single disease with similar pathobiological mechanisms (3). However, recent research has indicated that bronchiolitis is a syndrome with high heterogeneity (4–6). Our phenotypes derived through tIgE and virus data are in agreement with studies that have evaluated the heterogeneity in terms of atopic predisposition and respiratory virus (7, 8, 16, 17). For example, although a previous single-center study has shown, in 206 infants with RSV bronchiolitis, no association between tIgE levels and developing asthma (32), our previous multicenter study has revealed, in 921 infants hospitalized for bronchiolitis, the association of a higher tIgE level with an increased risk of developing childhood asthma (17). Additionally, another analysis has found that among infants hospitalized for bronchiolitis, those with IgE-mediated sensitization (to either food or aeroallergens) and RV type C infections had increased risks of developing childhood asthma (16). Furthermore, research using three cohort data of infants with bronchiolitis in the U.S. and Finland has reported the association of profiles based on atopic predisposition (e.g., previous history of eczema and parental history of asthma) and RV infection with different risks of developing asthma (8). Likewise, recent studies using omics data have also found endotypes based on the proteome (10), metabolome (11), and microbiome (metatranscriptome) (9) data. For example, an analysis of the nasopharyngeal metabolome (11) and microbiome (9) data of infants with bronchiolitis has reported that an endotype with downregulated type I interferon pathway had a higher risk of developing asthma. The current study builds on these earlier findings and extends them by identifying distinct tIgE-virus phenotypes with differential risks of developing asthma and their biological characteristics by using upper airway mRNA and miRNA data. A better understanding of tIgE-driven phenotypes may inform potential prevention strategies (e.g., the use of monoclonal anti-IgE antibody) against asthma in young infants.

Although the exact mechanisms underlying the observed phenotypes warrant further investigation, the current study found that tIgE^high phenotypes 3 and 4 had depleted type I interferon pathways and enriched antigen presentation pathways. Consistently, earlier studies of children and adults with atopic predisposition have demonstrated an inverse association between tIgE levels and virus-stimulated interferon-α production in peripheral blood plasmacytoid dendritic cells (18–20). Research has also revealed that FcϵRI expression is inversely related to interferon-α production and further inhibited by IgE-mediated cross-linking of FcϵRI (19, 20). Furthermore, a previous study has demonstrated, in RV infection to bronchial epithelial cells of adults with asthma, decreased interferon-β production, impaired apoptosis, and increased virus replication (33). In addition to dysregulated type I interferon pathways, phenotypes 3 and 4 had enriched antigen presentation pathways (e.g., MHC class II). Human MHC class II is highly expressed on the surface of airway epithelial cells and antigen-presenting cells to regulate the immune system, especially in allergic diseases and viral infections (34). Consistent with these phenotypes, genome-wide studies have also reported the association of single-nucleotide polymorphisms related to human MHC class II with tIgE level (35), IgE-mediated sensitization to aeroallergens (36), and asthma (35, 37). Additionally, another study has reported increased MHC class II expression in bronchial epithelial cells of adults with asthma (38). Furthermore, previous studies have demonstrated increased MHC class II expression in pulmonary mononuclear cells of RSV-infected mice (39), not in RV-infected alveolar epithelial cell lines (40).

Of the tIgE^high phenotypes, phenotype 4 (tIgE^highvirus^RSV-low/RV), which is at a higher risk of asthma, had depleted airway epithelial structure (e.g., cilia and intercellular junction) pathways. Consistent with our observations, growing evidence also suggests that airway epithelial abnormality and dysfunction were observed in adults and children with asthma (41, 42). For example, ciliated epithelium in the airway serves as a physical barrier by facilitating the clearance of mucus (43). Previous research using bronchial epithelial cells of adults (44) and children (45) with asthma has shown ciliary ultrastructural abnormality (e.g., loss of ciliated cells and ciliary disorientation) and dysfunction (e.g., low beat frequency). Likewise, intercellular junctions (e.g., tight junctions) also serve as an airway epithelial barrier by regulating permeability (46). In adults with asthma, tight junctions are disrupted, and epithelial permeability in bronchial epithelial cells is increased (47). Additionally, in airway epithelial cells of children with asthma, RV infection reduces the expression of tight junctional protein and increases permeability (48). Regardless of the complexity of these potential mechanisms, the current study has revealed that tIgE- and virus-derived phenotypes have distinct biological characteristics and differential asthma risks. Our data should advance research into developing strategies for asthma prevention based on tIgE-virus phenotyping.

The current study has several potential limitations. First, the current study did not have non-bronchiolitis “controls”. Nevertheless, this study did not aim to elicit infant phenotypes related to the incidence of bronchiolitis but to examine bronchiolitis phenotypes and their risk of developing asthma. Second, while the analytic and non-analytic cohorts did not differ in most clinical characteristics, the proportion of RSV infection was significantly different. Therefore, the phenotypes in the analytic cohort may not completely represent those in MARC-35. Third, although bronchiolitis involves inflammation in both upper and lower airways, the current study relied on upper airway mRNA and miRNA data. Nonetheless, research has shown that upper airway data represent inflammatory profiles in the lower airway, which is more invasive and difficult to collect in young infants (49). Fourth, the mRNA and miRNA data were derived from different upper airway specimens. Additionally, these specimens comprised heterogeneous cell populations, limiting pathways in the upper airway among the phenotypes. Fifth, these mRNA and miRNA data were measured only at the timing of bronchiolitis hospitalization. Although longitudinal measurement would be helpful, the single measurement successfully characterized the biologically distinct phenotypes. Sixth, in the analysis using the mRNA and miRNA data, the threshold of statistical significance was set at FDR <0.10. However, in the mRNA data, most identified pathways met FDR of <0.05. Seventh, in the sensitivity analysis excluding infants with a previous history of breathing problems, the phenotype-outcome associations did not show a statistical significance. However, the associations were directionally consistent between the main and sensitivity analysis. The consistency could support the robustness of the main analysis. Eighth, the asthma diagnosis may have been misclassified, and some children are going to develop asthma later. To address the potential limitations, the cohort is currently being followed until age 9 years. Ninth, while our findings are biologically and clinically plausible, they warrant further validation in an independent cohort. Lastly, our inferences may not be generalizable to populations other than infants hospitalized for bronchiolitis. Nevertheless, our data are highly relevant to >110,000 infants hospitalized for bronchiolitis annually in the U.S.—a population with a substantial morbidity burden.

Conclusions

By applying an unsupervised clustering approach to the tIgE and virus data from a prospective multicenter cohort study of infants hospitalized for bronchiolitis, we identified 4 clinically meaningful phenotypes. Specifically, the phenotype characterized by a high tIgE level and RSV (with a low genomic load) and RV infection had the greatest risk of developing asthma. Additionally, by using upper airway mRNA and miRNA data, these phenotypes revealed biologically distinct airway responses (e.g., type I interferon, antigen presentation, and airway epithelium structure pathways) during bronchiolitis. These findings offer an evidence-base for the early identification of high-risk infants during a critical period of airway development. Furthermore, these findings should advance asthma prevention strategies (e.g., modulating immune response by biologics targeting IgE, such as omalizumab) in this large patient population with a high morbidity burden.

Data availability statement

The data that support the findings of this study are available on request from the corresponding author through controlled access to be compliant with the informed consent forms of the MARC studies. The data are not publicly available due to privacy and ethical restrictions.

Ethics statement

The studies involving human participants were reviewed and approved by Massachusetts General Hospital. Written informed consent to participate in this study was provided by the participants’ legal guardian/next of kin.

Author contributions

Dr. RS carried out the main statistical analysis, drafted the initial manuscript, and approved the final manuscript as submitted. Drs. ZZ and TO reviewed and revised the initial manuscript, and approved the final manuscript as submitted. Drs. RF, JM, and ST collected the study data, reviewed and revised the initial manuscript, and approved the final manuscript as submitted. Drs. IR-T and MP-L carried out the bioinformatic analyses of the genomic data, reviewed and revised the initial manuscript, and approved the final manuscript as submitted. Drs. KH and CC conceptualized the study, obtained funding, supervised the statistical analysis, reviewed and revised the initial manuscript, and approved the final manuscript as submitted. All authors contributed to the article and approved the submitted version.

Funding

This study was supported by grants (UG3/UH3 OD-023253, R01 AI-127507, R01 AI-134940, and R01 AI-137091) from the National Institutes of Health (Bethesda, MD, U.S.). The content of this manuscript is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health. The funding organization was not involved in the collection, management, or analysis of the data; preparation or approval of the manuscript; or decision to submit the manuscript for publication.

Conflict of interest

Dr. ST receives royalties from UpToDate as a Section Editor and honoraria for lecture from Medscape.

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

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2023.1187065/full#supplementary-material

Abbreviations

adjOR, adjusted odds ratio; CDF, cumulative distribution function; CI, confidence interval; FDR, false discovery rate; GSEA, Gene Set Enrichment Analysis; IgE, immunoglobulin E; MARC, multicenter airway research collaboration; MHC, major histocompatibility complex; miRNA, microRNA; PAM, partition around medoids; RNA-seq, RNA sequencing; RSV, respiratory syncytial virus; RV, rhinovirus; tIgE, total immunoglobulin E.

References

1. Fujiogi M, Goto T, Yasunaga H, Fujishiro J, Mansbach JM, Camargo CA Jr, et al. Trends in bronchiolitis hospitalizations in the united states: 2000-2016. Pediatrics (2019) 144(6):e20192614. doi: 10.1542/peds.2019-2614

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Hasegawa K, Dumas O, Hartert TV, Camargo CA Jr. Advancing our understanding of infant bronchiolitis through phenotyping and endotyping: clinical and molecular approaches. Expert Rev Respir Med (2016) 10:891–99. doi: 10.1080/17476348.2016.1190647

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Ralston SL, Lieberthal AS, Meissner HC, Alverson BK, Baley JE, Gadomski AM, et al. Clinical practice guideline: the diagnosis, management, and prevention of bronchiolitis. Pediatrics (2014) 134:e1474–502. doi: 10.1542/peds.2014-2742

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Dumas O, Mansbach JM, Jartti T, Hasegawa K, Sullivan AF, Piedra PA, et al. A clustering approach to identify severe bronchiolitis profiles in children. Thorax (2016) 71:712–18. doi: 10.1136/thoraxjnl-2016-208535

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Arroyo M, Salka K, Perez GF, Rodríguez-Martínez CE, Castro-Rodriguez JA, Gutierrez MJ, et al. Phenotypical Sub-setting of the first episode of severe viral respiratory infection based on clinical assessment and underlying airway disease: a pilot study. Front Pediatr (2020) 8:121. doi: 10.3389/fped.2020.00121

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Niu H, Chang AB, Oguoma VM, Wang Z, McCallum GB. Latent class analysis to identify clinical profiles among indigenous infants with bronchiolitis. Pediatr Pulmonol (2020) 55:3096–103. doi: 10.1002/ppul.25044

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Dumas O, Erkkola R, Bergroth E, Hasegawa K, Mansbach JM, Piedra PA, et al. Severe bronchiolitis profiles and risk of asthma development in Finnish children. J Allergy Clin Immunol (2022) 149:1281–5.e1. doi: 10.1016/j.jaci.2021.08.035

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Fujiogi M, Dumas O, Hasegawa K, Jartti T, Camargo CA. Identifying and predicting severe bronchiolitis profiles at high risk for developing asthma: analysis of three prospective cohorts. EClinicalMedicine (2022) 43:101257. doi: 10.1016/j.eclinm.2021.101257

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Raita Y, Pérez-Losada M, Freishtat RJ, Hahn A, Castro-Nallar E, Ramos-Tapia I, et al. Nasopharyngeal metatranscriptome profiles of infants with bronchiolitis and risk of childhood asthma: a multicentre prospective study. Eur Respir J (2022) 60(1):2102293. doi: 10.1183/13993003.02293-2021

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Ooka T, Raita Y, Fujiogi M, Freishtat RJ, Gerszten RE, Mansbach JM, et al. Proteomics endotyping of infants with severe bronchiolitis and risk of childhood asthma. Allergy (2022) 77:3350–61. doi: 10.1111/all.15390

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Zhu Z, Camargo CA Jr, Raita Y, Fujiogi M, Liang L, Rhee EP, et al. Metabolome subtyping of severe bronchiolitis in infancy and risk of childhood asthma. J Allergy Clin Immunol (2022) 149:102–12. doi: 10.1016/j.jaci.2021.05.036

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Gould HJ, Sutton BJ. IgE in allergy and asthma today. Nat Rev Immunol (2008) 8:205–17. doi: 10.1038/nri2273

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Teach SJ, Gill MA, Togias A, Sorkness CA, Arbes SJ Jr, Calatroni A, et al. Preseasonal treatment with either omalizumab or an inhaled corticosteroid boost to prevent fall asthma exacerbations. J Allergy Clin Immunol (2015) 136:1476–85. doi: 10.1016/j.jaci.2015.09.008

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Phipatanakul W, Mauger DT, Guilbert TW, Bacharier LB, Durrani S, Jackson DJ, et al. Preventing asthma in high risk kids (PARK) with omalizumab: design, rationale, methods, lessons learned and adaptation. Contemp Clin Trials (2021) 100:106228. doi: 10.1016/j.cct.2020.106228

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Brockow I, Zutavern A, Hoffmann U, Grübl A, von Berg A, Koletzko S, et al. Early allergic sensitizations and their relevance to atopic diseases in children aged 6 years: results of the GINI study. J Investig Allergol Clin Immunol (2009) 19:180–7.

PubMed Abstract | Google Scholar

16. Hasegawa K, Mansbach JM, Bochkov YA, Gern JE, Piedra PA, Bauer CS, et al. Association of rhinovirus c bronchiolitis and immunoglobulin e sensitization during infancy with development of recurrent wheeze. JAMA Pediatr (2019) 173:544–52. doi: 10.1001/jamapediatrics.2019.0384

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Shibata R, Fujiogi M, Nanishi M, Ooka T, Mansbach JM, Teach SJ, et al. Total immunoglobulin e in infant bronchiolitis and risk of developing asthma. J Allergy Clin Immunol Pract (2022) 10:2761–3.e2. doi: 10.1016/j.jaip.2022.06.043

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Tversky JR, Le TV, Bieneman AP, Chichester KL, Hamilton RG, Schroeder JT. Human blood dendritic cells from allergic subjects have impaired capacity to produce interferon-alpha via toll-like receptor 9. Clin Exp Allergy (2008) 38:781–8. doi: 10.1111/j.1365-2222.2008.02954.x

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Gill MA, Bajwa G, George TA, Dong CC, Dougherty II, Jiang N, et al. Counterregulation between the FcepsilonRI pathway and antiviral responses in human plasmacytoid dendritic cells. J Immunol (2010) 184:5999–6006. doi: 10.4049/jimmunol.0901194

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Durrani SR, Montville DJ, Pratt AS, Sahu S, DeVries MK, Rajamanickam V, et al. Innate immune responses to rhinovirus are reduced by the high-affinity IgE receptor in allergic asthmatic children. J Allergy Clin Immunol (2012) 130:489–95. doi: 10.1016/j.jaci.2012.05.023

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Hasegawa K, Mansbach JM, Ajami NJ, Espinola JA, Henke DM, Petrosino JF, et al. Association of nasopharyngeal microbiota profiles with bronchiolitis severity in infants hospitalised for bronchiolitis. Eur Respir J (2016) 48:1329–39. doi: 10.1183/13993003.00152-2016

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Patro R, Duggal G, Love MI, Irizarry RA, Kingsford C. Salmon provides fast and bias-aware quantification of transcript expression. Nat Methods (2017) 14:417–9. doi: 10.1038/nmeth.4197

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Mjelle R, Aass KR, Sjursen W, Hofsli E, Sætrom P. sMETASeq: combined profiling of microbiota and host small RNAs. iScience (2020) 23:101131. doi: 10.1016/j.isci.2020.101131

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Kozomara A, Birgaoanu M, Griffiths-Jones S. miRBase: from microRNA sequences to function. Nucleic Acids Res (2019) 47:D155–62. doi: 10.1093/nar/gky1141

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Camargo CA Jr, Ingham T, Wickens K, Thadhani R, Silvers KM, Epton MJ, et al. Cord-blood 25-hydroxyvitamin d levels and risk of respiratory infection, wheezing, and asthma. Pediatrics (2011) 127:e180–7. doi: 10.1542/peds.2010-0442

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Fujiogi M, Zhu Z, Raita Y, Ooka T, Celedon JC, Freishtat R, et al. Nasopharyngeal lipidomic endotypes of infants with bronchiolitis and risk of childhood asthma: a multicentre prospective study. Thorax (2022) 77:1059–69. doi: 10.1136/thorax-2022-219016

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Mansbach JM, Piedra PA, Teach SJ, Sullivan AF, Forgey T, Clark S, et al. MARC-30 investigators. prospective multicenter study of viral etiology and hospital length of stay in children with severe bronchiolitis. Arch Pediatr Adolesc Med (2012) 166:700–6. doi: 10.1001/archpediatrics.2011.1669

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Hasegawa K, Jartti T, Mansbach JM, Laham FR, Jewell AM, Espinola JA, et al. Respiratory syncytial virus genomic load and disease severity among children hospitalized with bronchiolitis: multicenter cohort studies in the united states and Finland. J Infect Dis (2015) 211:1550–9. doi: 10.1093/infdis/jiu658

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U.S.A. (2005) 102:15545–50. doi: 10.1073/pnas.0506580102

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Huang H-Y, Lin Y-C-D, Li J, Huang K-Y, Shrestha S, Hong H-C, et al. miRTarBase 2020: updates to the experimentally validated microRNA-target interaction database. Nucleic Acids Res (2020) 48:D148–54. doi: 10.1093/nar/gkz896

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Thabane L, Mbuagbaw L, Zhang S, Samaan Z, Marcucci M, Ye C, et al. A tutorial on sensitivity analyses in clinical trials: the what, why, when and how. BMC Med Res Methodol (2013) 13:92. doi: 10.1186/1471-2288-13-92

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Bacharier LB, Cohen R, Schweiger T, Yin-Declue H, Christie C, Zheng J, et al. Determinants of asthma after severe respiratory syncytial virus bronchiolitis. J Allergy Clin Immunol (2012) 130:91–100.e3. doi: 10.1016/j.jaci.2012.02.010

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Wark PAB, Johnston SL, Bucchieri F, Powell R, Puddicombe S, Laza-Stanca V, et al. Asthmatic bronchial epithelial cells have a deficient innate immune response to infection with rhinovirus. J Exp Med (2005) 201:937–47. doi: 10.1084/jem.20041901

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Wosen JE, Mukhopadhyay D, Macaubas C, Mellins ED. Epithelial MHC class II expression and its role in antigen presentation in the gastrointestinal and respiratory tracts. Front Immunol (2018) 9:2144. doi: 10.3389/fimmu.2018.02144

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Moffatt MF, Gut IG, Demenais F, Strachan DP, Bouzigon E, Heath S, et al. A large-scale, consortium-based genomewide association study of asthma. N Engl J Med (2010) 363:1211–21. doi: 10.1056/NEJMoa0906312

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Gheerbrant H, Guillien A, Vernet R, Lupinek C, Pison C, Pin I, et al. Associations between specific IgE sensitization to 26 respiratory allergen molecules and HLA class II alleles in the EGEA cohort. Allergy (2021) 76:2575–86. doi: 10.1111/all.14820

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Zhu Z, Lee PH, Chaffin MD, Chung W, Loh P-R, Lu Q, et al. A genome-wide cross-trait analysis from UK biobank highlights the shared genetic architecture of asthma and allergic diseases. Nat Genet (2018) 50:857–64. doi: 10.1038/s41588-018-0121-0

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Vignola AM, Campbell AM, Chanez P, Bousquet J, Paul-Lacoste P, Michel FB, et al. HLA-DR and ICAM-1 expression on bronchial epithelial cells in asthma and chronic bronchitis. Am Rev Respir Dis (1993) 148:689–94. doi: 10.1164/ajrccm/148.3.689

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Beyer M, Bartz H, Hörner K, Doths S, Koerner-Rettberg C, Schwarze J. Sustained increases in numbers of pulmonary dendritic cells after respiratory syncytial virus infection. J Allergy Clin Immunol (2004) 113:127–33. doi: 10.1016/j.jaci.2003.10.057

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Papi A, Stanciu LA, Papadopoulos NG, Teran LM, Holgate ST, Johnston SL. Rhinovirus infection induces major histocompatibility complex class I and costimulatory molecule upregulation on respiratory epithelial cells. J Infect Dis (2000) 181:1780–4. doi: 10.1086/315463

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Heijink IH, Kuchibhotla VNS, Roffel MP, Maes T, Knight DA, Sayers I, et al. Epithelial cell dysfunction, a major driver of asthma development. Allergy (2020) 75:1902–17. doi: 10.1111/all.14421

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Carlier FM, de Fays C, Pilette C. Epithelial barrier dysfunction in chronic respiratory diseases. Front Physiol (2021) 12:691227. doi: 10.3389/fphys.2021.691227

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Wanner A, Salathé M, O’Riordan TG. Mucociliary clearance in the airways. Am J Respir Crit Care Med (1996) 154:1868–902. doi: 10.1164/ajrccm.154.6.8970383

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Thomas B, Rutman A, Hirst RA, Haldar P, Wardlaw AJ, Bankart J, et al. Ciliary dysfunction and ultrastructural abnormalities are features of severe asthma. J Allergy Clin Immunol (2010) 126:722–9.e2. doi: 10.1016/j.jaci.2010.05.046

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Cokuğraş H, Akçakaya N, Seçkin, Camcioğlu Y, Sarimurat N, Aksoy F. Ultrastructural examination of bronchial biopsy specimens from children with moderate asthma. Thorax (2001) 56:25–9. doi: 10.1136/thorax.56.1.25

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Hartsock A, Nelson WJ. Adherens and tight junctions: structure, function and connections to the actin cytoskeleton. Biochim Biophys Acta (2008) 1778:660–9. doi: 10.1016/j.bbamem.2007.07.012

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Xiao C, Puddicombe SM, Field S, Haywood J, Broughton-Head V, Puxeddu I, et al. Defective epithelial barrier function in asthma. J Allergy Clin Immunol (2011) 128:549-56.e1-12. doi: 10.1164/ajrccm-conference.2011.183.1_MeetingAbstracts.A6367

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Looi K, Buckley AG, Rigby PJ, Garratt LW, Iosifidis T, Zosky GR, et al. Effects of human rhinovirus on epithelial barrier integrity and function in children with asthma. Clin Exp Allergy (2018) 48:513–24. doi: 10.1111/cea.13097

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Poole A, Urbanek C, Eng C, Schageman J, Jacobson S, O’Connor BP, et al. Dissecting childhood asthma with nasal transcriptomics distinguishes subphenotypes of disease. J Allergy Clin Immunol (2014) 133:670–8.e12. doi: 10.1016/j.jaci.2013.11.025

PubMed Abstract | CrossRef Full Text | Google Scholar