Role and Diagnostic Performance of Host Epigenome in Respiratory Morbidity after RSV Infection: The EPIRESVi Study

Background: Respiratory syncytial virus (RSV) infection has been associated with the subsequent development of recurrent wheezing and asthma, although the mechanisms involved are still unknown. We investigate the role of epigenetics in the respiratory morbidity after infection by comparing methylation patterns from children who develop recurrent wheezing (RW-RSV), subsequent asthma (AS-RVS), and those experiencing complete recovery (CR-RSV).

Methods: Prospective, observational study of infants aged < 2 years with RSV respiratory infection admitted to hospital and followed-up after discharge for at least three years. According to their clinical course, patients were categorized into subgroups: RW-RSV (n = 36), AS-RSV (n = 9), and CR-RSV (n = 32). The DNA genome-wide methylation pattern was analyzed in whole blood samples, collected during the acute phase of the infection, using the Illumina Infinium Methylation EPIC BeadChip (850K CpG sites). Differences in methylation were determined through a linear regression model adjusted for age, gender and cell composition.

Results: Patients who developed respiratory sequelae showed a statistically significant higher proportion of NK and CD8T cells (inferred through a deconvolution approach) than those with complete recovery. We identified 5,097 significant differentially methylated positions (DMPs) when comparing RW-RSV and AS-RVS together against CR-RSV. Methylation profiles affect several genes involved in airway inflammation processes. The most significant DMPs were found to be hypomethylated in cases and therefore generally leading to overexpression of affected genes. The lead CpG position (cg24509398) falls at the gene body of EYA3 (P-value = 2.77×10^-10), a tyrosine phosphatase connected with pulmonary vascular remodeling, a key process in the asthma pathology. Logistic regression analysis resulted in a diagnostic epigenetic signature of 3-DMPs (involving genes ZNF2698, LOC102723354 and RPL15/NKIRAS1) that allows to efficiently differentiate sequelae cases from CR-RSV patients (AUC = 1.00). Enrichment pathway analysis reveals the role of the cell cycle checkpoint (FDR P-value = 4.71×10^-2), DNA damage (FDP-value = 2.53×10^-2), and DNA integrity checkpoint (FDR P-value = 2.56×10^-2) in differentiating sequelae from CR-RSV patients.

Conclusions: Epigenetic mechanisms might play a fundamental role in the long-term sequelae after RSV infection, contributing to explain the different phenotypes observed.

Introduction

Respiratory syncytial virus (RSV) is a common pathogen that infects children by two years of age and is the leading global cause of hospitalization of infants (1). It is the principal cause of acute lower respiratory infections (ALRI) in young children, and it is associated with morbidity and mortality in childhood (2). Approximately 34 million new ALRI episodes in children worldwide were attributable to RSV (1), a huge number that resulted in 3.2 million hospital admissions and almost 60,000 global childhood deaths each year. No effective vaccines are available to treat RSV yet, although some good candidates are being tested in large-state clinical trials (3, 4). For the time being, palivizumab is the only monoclonal antibody approved to prevent severe RSV in infants and children at high risk for severe disease (5, 6).

Besides the acute burden of RSV, a growing body of evidence from epidemiology data supports that RSV infection in the first three years of life can be directly correlated with long-term respiratory morbidities, such as recurrent wheezing and asthma (7, 8). It has been observed that RSV, like other respiratory viruses, causes a “hit and run” phenomenon, characterized by the increased risk of developing recurrent wheezing and asthma in childhood after the infection, as a permanent phenotype that persists long after the virus clearance (9). Wheezing is the typical high-pitched, whistling sound made during breathing. Asthma, on the other hand, is characterized by abnormalities in lung function that include variable airway obstruction and increased bronchial reactivity (10). Recognizing asthma is usually obvious and most of the time asthmatic patients also report wheezing episodes; however, it is very challenging to predict and distinguish which children will present only early-life symptoms, from those whose symptoms persist, and those who may develop definitive wheeze or asthma (11).

The risk of wheezing and/or asthma incidence has been increasingly related to a combination of genetic and environmental factors as well as the severity of the respiratory infection (12–14). Several efforts have been made to establish the ultimate causes underlying the apparition of respiratory sequelae as a consequence of RSV infection, but further complementary studies are still necessary to understand the underlying genetic and molecular mechanisms. It has been hypothesized that prenatal interaction between the maternal and the child immune system as well as the high cytokines production after the infection can affect the lung structure and function inducing changes in the regulation of immune response (15).

The immune system plays a pivotal role in the sequelae reported in patients suffering from respiratory infections. The immune response shows many functional differences in neonates and adults, and these differences might be associated with epigenetic modifications of genes that control inflammation and immune response (16). Thus, the modulation of epigenetic mechanisms governs the immune cell phenotype and function allowing the external environment to influence the immune response outcome (17).

Epigenetic regulation behaves as a dynamic interface between genome and environment; in the case of viral infection, this regulation in host defense cells is directly related to disease development (18). Pieces of evidence support the involvement of epigenetic mechanisms in the modulation of the interaction between host and pathogen (19), as it also occurs as a consequence of vaccination [e.g. rotavirus infection (20)]. For instance, it has been shown that some viral infections (e.g. by herpesvirus, KSHV, EBV) can lead to the modification of host epigenetic marks, and this probably contributes to the establishment of latency and some pathogenic roles (21, 22).

DNA methylation involves the addition of methyl groups to the DNA molecule and is the major epigenetic factor influencing gene activities. It is generally associated with a decrease of gene expression when it occurs in promoters; however, according to recent evidence (23, 24) the effect of DNA methylation occurring in low CpG density regions, such as the gene body, or intergenic regions, might have also great impact on the regulation of gene expression, with repressive or permissive effects. Recently, Elgislouli and colleagues (25) observed an alteration of the methylation pattern of a particular enhancer region in the perforin 1 gene (PRF1), an essential cytotoxic protein for the control of viral infection, occurring after severe RSV infection. Infants hospitalized with exacerbated RSV-induced bronchiolitis exhibited decreased methylation of the perforin-1 enhancer after four years of follow-up, suggesting that immune response changes due to RSV could persist even years after infection (26). Besides, Wang et al. (18) showed that cultured bronchial human epithelial cells (BECs) infected with RSV presented a high expression level of the gene NODAL (a member of the TGF-beta superfamily), whose promoter was found to be hypermethylated in normal BECs, resulting in increased Th2 and Th17 skewing of T cells. In murine models, it was observed an overexpression of the demethylase genes Kdm5b and H3K4, in dendritic cells after RSV infection (27); expression of these genes has the potential to repress transcription of type I IFN and other innate cytokines, causing a decrease of pro-inflammatory cytokines and expression of a Th2 phenotype. Methylation of H3K4 in regulatory T cells by the histone methyltransferase SMYD3 is necessary to control inflammation in the lungs after RSV infection (28).

The aim of the present study is to analyze the DNA methylation changes related to the development of wheezing and/or asthma induced after RSV infection. To the best of our knowledge, this is the first time that an epigenomic approach is applied to the investigation of epigenetic modifications due to RSV infection and the subsequent long-term respiratory sequelae.

Material and Methods

Study Subjects and Design

EPIRSVI is a prospective, transversal, observational study of children admitted to the Hospital Clínico Universitario de Santiago de Compostela, and the Complejo Hospitalario Universitario de Ourense, part of the GENDRES consortium (www.gendres.org), for a respiratory infection due to RSV (Supplementary Table 1). The recruitment of patients was performed from 2010 to 2015, during their acute phase of RSV infection (within seven days from the beginning of symptoms). Patients, independently from their age, were all followed up from their recruitment for a time-lapse of at least three years after discharge to specifically monitor the onset of wheezing or asthma. Patients were categorized according to their clinical course into three different subgroups (n = 77): (a) recurrent wheezing RSV cases (RW-RSV; n = 36), (b) asthma RSV cases (AS-RSV; n = 9); and (c) not-wheezing/asthma RSV cases with complete recovery (CR-RSV; n = 32). The three groups matched perfectly for age, gender, and disease severity at the time of sample collection.

Wheezing was determined here as a respiratory episode occurring with wheezing lasting more than one day. The interval between two episodes was defined as a period of at least seven days without respiratory symptoms. Instead, recurrent wheeze was defined as three or more episodes of wheezing during the first year of life. The diagnosis of asthma in preschool children was carried out by a pediatric pulmonologist based on the reiterated presence of compatible symptoms (several episodes of bronchial obstruction, generally witnessed by the physician), adequate response to regular treatment, and the exclusion of other alternative diagnoses, following the international practice guidelines. From the clinical perspective, the diagnosis of asthma in these children involves a series of peculiarities inherent to patients of this age, such as limitations for performing pulmonary function tests, greater attention to differential diagnoses, limited response to common asthma treatments, and a high probability of symptoms remission during childhood.

DNA Isolation and Methylation Profiling

Blood samples were collected in EDTA tubes, within 24/48 hours from the hospital admission of each patient. Genomic DNA was isolated using the Wizard^® Genomic DNA Purification Kit (PROMEGA) QIAamp DNA mini kit (Qiagen) according to the manufacturer’s instructions. After DNA extraction, DNA was purified, and concentration was evaluated using a Nanodrop and Picogreen assay.

Genomic DNA from whole blood samples was bisulfite-treated using the EZ-96 DNA Methylation kit (Zymo Research Corp) following the manufacturer’s recommendations for Infinium assays. The technique consists of the conversion of unmethylated cytosine residues to uracil through deamination, while leaving the methylated cytosine (5-mC) unaffected.

Following bisulfite conversion, the samples were hybridized to Illumina Infinium MethylationEPIC BeadChip. This microarray allows the examination of >850K methylation sites quantitatively across the genome at single-nucleotide resolution. The microarray coverage includes 99% of the RefSeq genes, 95% of CpG islands, non-CpG islands, and differentially methylated sites, high coverage of transcription factors binding sites, miRNA promoter regions, and enhancer sequences. We compared the methylome of RW-RSV, AS-RSV, with CR-RSV to evaluate whether the divergent epigenetic marks could be specific for the development of recurrent wheezing and/or asthma.

Data Processing, Cell Deconvolution and Statistical Analysis

DNA methylation quality control, processing, normalization, and statistical analyses were performed in the statistical R software (29), using different Bioconductor packages, and following the workflow of EPIC methylation analysis described recently by Fortin et al. (30). The raw intensity files (IDAT), were imported into R and were preprocessed and transformed into β and M-values using the minfi (v1.30.0) package (31, 32). The β-values are defined as the intensity ratio of the methylated signals over the total (methylated and unmethylated) signals for each site. Their values range between 0, illustrating that the position is not methylated at all in all the cells and 1, indicating that all cells in the sample have a given position methylated (33). The M-values were calculated as the log₂ ratio of the intensities of methylated probes vs. unmethylated probes; their values range from –1 to 1, where values close to 0 define a similar intensity between the methylated and unmethylated probes, positive M-values indicate that there are more methylated than unmethylated molecules, while negative M-values mean the opposite (34).

Counting blood cells is a necessary step to adjust for individual differences in cellular heterogeneity in the blood sample from which genomic DNA was extracted. Because of the extreme cell-type specificity of DNA methylation (35), variations in cell-type composition between phenotypes can confound analyses and bring to a wrong interpretation of the results. Therefore, before checking for differentially methylated positions (DMPs), we used the package FlowSorted.Blood.EPIC and its modified function estimateCellCounts2 for the estimation of blood cell type composition. At the base of the function, there is a modification of the Houseman algorithm (36), that allows estimating the relative proportions of white blood cell subtypes, namely, CD4+ T-lymphocytes, CD8+ T-lymphocytes, natural killer (NK) cells, B-lymphocytes, monocytes, and neutrophils.

Infinium 850K has two kinds of probes, Infinium I and Infinium II, that are not directly comparable, therefore, a correction must be applied to control for design bias. The normalization procedures were performed using the function preprocessQuantile available in the minfi package (37). Moreover, probes underwent several filtering processes, by removing probes with P-value > 0.01 and probes located in the sex chromosomes. Additionally, sites containing SNPs or with a minor allele frequency (MAF) < 0.05 were also excluded from the data analysis because probe binding might be affected by genetic variation in the binding area. Finally, we removed probes known to have cross-reaction.

DMPs were identified with limma (38) assuming a linear model, where M-values of each probe were used as quantitative dependent variables in all analyses, and including cell composition, age, and gender as covariates in all the models. After running the linear model, we applied the statistical analysis using an empirical Bayes method to moderate standard errors. DMPs were filtered using a significance threshold of adjusted P-value < 0.01 [Benjamini–Hochberg method (39)].

We used a Principal Component Analysis (PCA) to compare DNA methylation of children with complete recovery (CR-RVS) and those with respiratory sequelae after an RSV infection (RW-RVS + AS-RVS). The threshold Delta β of 0.1 and adjusted P-value < 0.01 were used to identify candidate DMPs among the comparison groups. The diagnostic efficiency of the most significant candidate DMPs was evaluated using receiver operating characteristic (ROC) curve analyses. Furthermore, the candidate DMPs were included in the logistic regression analysis using the Parallel Regularised Regression Model Search (PReMS) (40) method that allows identifying the minimum epigenetic signature with the highest diagnostic performance to differentiate between RW/AS-RSV and CR-RSV groups. This logistic regression model balances small positions number with accurate discrimination, minimizing the number of biomarkers selected in the signatures. The optimal model size (with the lowest out-of-sample log-likelihood) was determined using 20 cross-validation folds. In addition to the best model, results for model sizes with predictive log-likelihood within one-standard error of the best model is also shown. We randomly split the cohort into 75% training and 25% test sets, ensuring equal proportions of RW/AS-RSV and CR-RSV in each set. The epigenetic signature was identified in the training set and afterwards validated in the test set. Wilcoxon test was used to assess statistical significance between patient groups.

To deeply analyze the functional interpretation of our reported differential probes, we performed pathways enrichment analysis on the candidate DMPs applying the methylglm and methylRRA functions of the R methylGSA package (41). To adjust for the number of CpGs, the first method employs a logistic regression model incorporating the number of CpGs as a covariate, while the second approach makes use of a Robust Rank aggregation (RRA) procedure. Indeed, the number of CpGs is variable for genes of similar length and DNA methylation studies give rise to multiple CpG association P-values per gene.

To find co-methylated modules related to respiratory sequelae, we built a signed weighted correlation network using the WGCNA package (42) with the CpGs showing the most different β-values between samples (top 25% with a higher variance; n = 123,617). We selected a soft-thresholding power of 12 based on the criterion of scale-free topology after testing a set of candidate powers (Supplementary Figure 1A). As module detection parameters, we chose a minimum module size of 30, a medium sensitivity for cluster splitting and a 0.25 as dendrogram cut heigh threshold for module merging. We calculated gene significance (GS) to detect significant associations between modules/genes and phenotype. Module membership (MM), as a measurement of intramodular connectivity, was also calculated by correlating the methylation profile with the eigengene of a given module. We explored the correlation between GS and MM and calculated the average absolute gene significance for all genes within a module to find the most important modules. The top hub genes within the most relevant module were selected using both MM > 0.8 and GS ≤ –0.7 as thresholds. We studied the biological significance of the most important module related to the trait by performing an over-representation analysis through the Clusterprofiler R package (43), and using Gene Ontology (GO) database as a reference.

Results

Patient Characteristics

Table 1 summarizes the clinical characteristics of the subjects studied and maintained for the downstream statistical analysis after filtering raw data (n = 68), and it provides information regarding the most important risk factors related to the development of asthma in early childhood. For the clinical characteristics of the full cohort see Supplementary Table 1. More males than females represented sequelae and complete recovery groups (60.7% and 70.0% for CR-RSV and RW/AS-RSV, respectively). Concerning disease severity, all children were hospitalized because of RSV infection and their admission lasted > 5 days. Most subjects were diagnosed with bronchiolitis and required oxygen during their hospital admission. Most children were under six months of age; however, a significant difference (P < 0.001) was observed in the mean age of both groups, being children with respiratory sequelae slightly younger than patients with complete recovery after RSV. No comorbidities were present in the study cohort apart from prematurity (see Table 1). Thus, admissions prior to RSV were due to prematurity, bronchiolitis (mainly), rotavirus gastroenteritis in one case and hypernatremia dehydration in another case. Most of the children admitted below the age of one year were normally more common to have oxygen requirement, suffer apneas or need of intravenous fluids due to lack of appetite, therefore the number of admissions is generally higher compared to older ages. According to the previous medical history of these patients, the only clinical variable found to be statistically significant was the personal history of atopic dermatitis (Fisher exact test; P-value = 0.017), while alimentary and stational allergies showed no significant difference between the two groups. The family history of asthma and respiratory problems were not found to be significantly different between groups. Finally, we observed a significant difference (Fisher exact test; P-value = 0.005) in the number of bronchitis episode prior to hospitalization in patients with respiratory sequelae, and a higher number of suspected bacterial infection (Fisher exact test; P-value <0.001) in children with normal recovery.

TABLE 1

Table 1 Clinical characteristics of patients left for downstream analysis, classified as recurrent wheezing RSV (RW-RSV), asthma RSV (AS-RSV), and complete recovery RSV cases (CR-RSV).

DMPs Between Complete Recovered and Sequelae Cases

We interrogated the genome-wide DNA methylation profiles of 77 children during their acute phase of RSV infection that experienced complete recovered or developed respiratory sequelae.

During the data preprocessing nine samples were discarded because they did not pass the initial quality control; therefore, 68 out of 77 samples were included in the statistical analysis. After removing poorly performing probes, namely, those overlapping with SNPs or hybridized to multiple locations in the genome and those located on the X and Y chromosomes, a total of 811,035 probes were preserved for downstream analysis.

Cell deconvolution analysis indicates a significant difference in the number of NK and CD8+T cells between sequelae and recovered cases (P-value = 0.03 and P-value = 0.04, respectively); Figure 1A. When we compared RW/AS-RSV against CR-RSV patients, adjusting the linear model for cell-type composition, age, and gender, we found 5,097 significantly DMPs (FDR P-value < 0.01), corresponding to 3,278 unique genes according to the Illumina Human Methylation EPIC manifest annotation file (44) (Supplementary Table 2). Among these, 1,155 DMPs were hypomethylated (22.7%) and 3,942 hypermethylated (77.3%) in CR-RSV. PCA of these significant DMPs showed a clear separation of the two groups along with the first principal component (PC1; accounting for 43.2% of the variance) (Figure 1B).

FIGURE 1

Figure 1 (A) Boxplot showing the proportion of leukocyte cell type in RW/AS-RSV and CR-RSV groups. Red rectangles highlight the two types of cells that show statistically significant differences between groups. (B) PCA of the significant DMPs (FDR P-value < 0.01) between RW/AS-RSV and CR-RSV cases.

There are regional differences between the distribution of statistical different hypomethylated and hypermethylated CpGs when examined in a genomic context. Thus, hypermethylated positions are more concentrated in the gene body (59.2%), while the hypomethylated positions are distributed more homogeneously along the gene body (29.8%), promoter regions (45.6%), and transcriptions start sites [TSS1500 (22.1%), TSS200 (23.5%)], with a lower amount falling at ‘UTR (15.4%) and 1^st Exon (7.2%) (Figure 2A). Island’s regions are overrepresented by hypomethylated sites (38.8%), while the OpenSea regions are more characterized by hypermethylated DMPs (61%); the distribution of DMPs at shelf and shore was quite similar between the two comparison groups (Figure 2B). Moreover, all chromosomes were more represented by hypermethylated DMPs than the hypomethylated ones (Figure 2C).

FIGURE 2

Figure 2 (A) Pie chart showing the percentage of the significant DMPs in the genomic context, and (B) according to their distribution in the Island context. (C) Barplot showing the distribution of hypomethylated and hypermethylated DMPs in the chromosomes.

To narrow down the analysis and reduce the large number of significant methylated positions observed, we focused on finding the most consistent methylation changes by selecting only significant CpGs with an absolute Delta β > 0.10 when comparing CR-RSV against RW/AS-RSV (Figure 3A). This analysis found 28 CpGs (Table 2) associated with 23 unique genes; five of them were moderately hypermethylated, the remaining showed a lower methylation level for sequelae cases.

FIGURE 3

Figure 3 (A) Volcano plots for differential DNA methylation status. The x-axis shows the mean DNA methylation (β -value) difference, whereas the y-axis indicates the –log₁₀ of the adjusted P-value for each CpG site. The most significant DMPs (n = 28) in children with respiratory sequelae when compared with children complete recovered (threshold: Delta β > 0.10, FDR P-value < 0.01) are plotted in red (hypomethylated sites) and blue (hypermethylated CpGs). (B) Receiver operating characteristic (ROC) curves indicate that the accuracy of the test based on the reported CpGs is very high (AUC > 90%) when comparing RW/AS-RSV and CR-RSV groups. (C) Boxplot of the nine most significant DMPs with a Delta β > 0.10 in the comparison RW/AS-RSV vs. CR-RSV cases.

TABLE 2

Table 2 The 28 most significant DMPs in the comparison RSV sequelae patients vs. CR-RSV cases (FDR P-value < 0.01, absolute average β-value > 0.10).

Methylation Signature for RSV Sequelae

A total of 9 out of the 28 significant positions with absolute Delta β > 0.10 showed an FDR P-value < 1×10^-4, and were selected to evaluate their ability, in terms of methylation levels, to discriminate RW/AS-RSV from CR-RSV groups (Figure 3B). All the selected positions exhibit hypomethylation pattern in RW/AS-RSV patients when compared to the CR-RSV group (Figure 3C). Two of the nine positions, namely, the cg00044440 (Chromosome 5; OpenSea; FDR P-value = 6.77×10^-6; Delta β = -0.12), and the cg06579481 (Chromosome 7; North-Shelf; FDR P-value = 8.57×10^-6; Delta β = -0.12), were not annotated to a gene. The remaining seven DMPs, instead, were associated to eight different genes: cg24509398 (FDR P-value = 2.77×10^-10; Delta β = –0.14) falls in the body of EYA3; cg15619333 (FDR P-value = 4.39×10^-9; Delta β = –0.116×10^-1) locates in the body of ZNF268; cg12216772 (FDR P-value = 3.33×10^-8; Delta β = –0.12) is in the 5’UTR of ANUBL1; cg15822108 (FDR P-value = 4.28×10^-7; Delta β = –0.11) settles in the body of NAA15; cg02737727 (FDR P-value = 1.84×10^-9; Delta β = –0.11) falls in the body of LOC102723354; cg14819618 (FDR P-value = 2.26×10^-6; Delta β = –0.10) resides within the body of PRAGMIN; and finally, cg23621438 (FDR P-value = 1.84×10^-9; Delta β = –0.10), locates within the TSS1500 of MIR6817 and the body of CRYBB2P1.

A logistic regression analysis performed with the 28-CpG lead positions found a minimum epigenetic signature composed of three methylation markers (cg15619333, cg02737727, and cg06259441) as the best diagnostic model to differentiate RW/AS-RSV from CR-RSV in the training set (n = 52) (Figure 4A). The cg15619333 and the cg02737727 were previously described among the candidate CpGs with the lowest FDR P-value, while the North Shore cg06259441 resides within two overlapping genes RPL15 and NKIRAS1 in chromosome 3. Interestingly, this last position showed an opposite pattern of methylation, being highly methylated in RW/AS-RSV cases.

FIGURE 4

Figure 4 (A) Optimal model size according to logistic regression analysis. The x-axis represents the training set predictive log-likelihood, while the y-axis shows the number of genes of the signature. Solid grey bars indicate one standard error (SE) of the predictive log-likelihood. Vertical dashed lines show the model with the best predictive log-likelihood and the model within one-SE of the best model. (B) Violin and boxplots of the predicted values from the posterior mean of the optimal model in training and test cohorts. (C) ROC curves of the 3-CpGs signature from both training and test cohorts show the area under the curve (AUC) and 95% confidence intervals (CIs). AUC values for the whole cohort are also displayed.

This 3-position signature showed high diagnostic performance with the training set (P-value = 1.00×10^-14; Figure 4B), with AUC value of 100% (95% CI: 1-1) (Figure 4C). When the 3-signature model was applied to the test set (n = 16), we obtained similar prediction results (P-value = 2.50×10^-4) with an AUC value of 100% (95% CI: 1–1) as in the training set (Figure 4C). The diagnosis accuracy of the model was not affected when all samples were analysed together (Figure 4C).

Enrichment Pathways Analysis

To address the potential biological significance of DMPs, a gene set enrichment and pathways analysis was carried out. When comparing the sequelae group against CR-RSV cases, many significantly enriched pathways (FDR P-value < 0.05) in the GO database was detected (Figure 5A). Among the top overrepresented categories, there were those related to several cell cycle processes, such as the cell cycle checkpoint (FDR P-value = 4.71×10^-2), DNA damage (FDR P-value = 2.53×10^-2), and DNA integrity checkpoint (FDR P-value = 2.56×10^-2). The enrichment pathways analysis performed exclusively on positions overlapping promoters (TSS200 and TSS1500), revealed significant enrichment in the GO ubiquitin ligase binding (FDR P-value = 1.21×10^-2), and response to transforming growth factor beta (FDR p.value =1.21x10^-2), among other pathways (Figure 5B).

FIGURE 5

Figure 5 (A) Dot plot of the top GO pathways (FDR P-value < 0.05) obtained between RW/AS-RSV and CR-RSV (methylglm approach). (B) Bar plot of the enrichment pathways analysis performed considering only CpGs within the promoter regions (methylRRA approach). Size along the x-axis indicates the number of genes involved in each pathway. Colors correspond to the different FDR P-values associated with the pathways.

Co-Methylation Modules

We detected 99 co-methylation modules (Supplementary Figure 1B) when comparing methylation patterns of RW/AS-RSV against CR-RSV patients. The top one (showing the highest correlation with the RW/AS-RSV phenotype; P-value = 1.67×10^-9; Supplementary Table 3; Supplementary Figure 1C) includes 3,147 CpGs (Supplementary Table 3). This module shows a negative correlation with the trait (–0.65), pointing to a hypomethylation state in RW/AS-RSV concerning the CR-RSV group. This pattern can also be inferred from the heatmap and samples eigengenes plots, with overall higher methylation values in the CR-RSV cohort (Supplementary Figures 2A, B). A highly significant correlation between gene significance (GS) and module membership (MM) was observed for this module, meaning that CpGs strongly associated with RW/AS-RSV were also core elements of the module (Supplementary Figure 2C). Interestingly, the vast majority of the main hub CpGs within this module corresponds with the previously identified 5.097 DMPs (~92%). Three of them are among the 28-CpGs previously described as most remarkable DMPs between sequelae and non-sequelae groups (cg24509398 [EYA3], cg15822108 [NAA15], cg15619333 [ZNF268]; Supplementary Table 3). Genes linked to the most significant module positions were mainly involved in pathways related to GTPase activity, epigenetic modifications, and cell cycle; and molecular functions related to actin, cadherin, and calmodulin binding and also, with focal adhesion components (Supplementary Table 3).

Identification of DMPs Associated With Sequelae Sub-Phenotypes

We next investigated RW-RSV and AS-RSV separately using the same linear model and thresholds described above for the merged sequelae group. When comparing RW-RSV vs. CR-RSV we detected 3,432 significant DMPs, of which 2,663 were hypermethylated, and 769 were hypomethylated in CR-RSV cases. Likewise, the comparison of AS-RSV patients vs. CR-RSV allowed us to identify 1,358 DMPs, with the vast majority being hypermethylated (n = 1,059), and 299 being hypomethylated in CR-RSV. A total of 512 significant DMPs were shared by the two contrasts, and also overlap with the most significant DMPs found when comparing respiratory sequelae cases vs. CR-RSV.

As shown above, a PCA considering DMPs between sequelae cases vs. CR-RSV patients allows also to separate the two sequelae phenotypic groups in its second component (Figure 2: see PC2; accounting for ~9% of the variation). By running a linear model with the contrast AS-RSV vs. RW-RSV and using a less stringent threshold (FDR P-value <0.05), a total of 47 significantly DMPs associated with 41 genes could be detected. A PCA was carried out exclusively with these 41 DMPs and the two sequelae phenotypes clearly separate the two groups in its PC1 (49% of the full variance) (Figure 6). Of the significant CpGs, 30 positions were found to be hypermethylated in RW-RSV when compared to AS-RSV, while the remaining 17 exhibit the opposite pattern (Table 3). The position with the highest difference in β value (FDR P-value = 1.8×10^-2; Delta β >0.14) was the cg18873878 located within the TP73 gene, while the most significant CpGs (FDR P-value = 0.11×10^-1) was the cg05838113 located within the body of ADAM8 gene.

FIGURE 6

Figure 6 PCA of the most different DMPs (FDR P-value < 0.05) shows an almost clear separation between RW-RSV and AS-RSV group.

TABLE 3

Table 3 PDMPs between AS-RSV and RW-RSV groups are ordered by Delta β-values (FDR P-value < 0.05).

Discussion

It is known that host factors might impact subsequent respiratory morbidity in RSV infection, but the underlying biological mechanisms through which recurrent wheezing and asthma are prone to emerge after RSV infection are still unknown.

The analysis of the DNA methylome in RSV infected children allowed us to find >5.000 DMPs mainly located in inflammatory genes when comparing patterns in individuals with sequelae from those completely recovered (Figure 7). Overall, the findings suggest that epigenetic mechanisms might play a fundamental role in the long-term sequelae after RSV infection. The association between RSV infection and the subsequent development of wheezing and/or asthma in infants has been widely reported (7, 8, 45, 46). This correlation seems to be very complex owed to the several host environmental factors that influence the expression of respiratory sequelae developed after the infection, such as genetics, age, prematurity, RSV bronchiolitis, atopic dermatitis, and maternal asthma. The patients analyzed in the present study had a very similar clinical course, with the particularity that those developing respiratory sequelae were two months younger than those completely recovered and show more incidence of atopic dermatitis; similar findings were already described to be related with respiratory morbidity after RSV infection (47).

FIGURE 7

Figure 7 Study design and overview of main results.Supplementary Material Legend.

Estimation of the relative proportion of cell types in peripheral blood samples of patients revealed a statistically significantly higher proportion of NK and CD8T cells in children with respiratory sequelae than those with complete recovery (Figure 7). This is in line with evidences demonstrating that RSV infection can induce severe acute lung immune injury promoting the accumulation of lung NK cells at the early stage of infection in mice, as a consequence of the increased production of IFN-γ (48). Similar characteristics were described for CD8+ T cells, which are activated after RSV infection to produce inflammatory cytokines (49), and seem to enhance airway inflammation and airway dysfunction in mice (50). However, there is growing body of evidence suggesting that NK cells, as well as CD8+T cells are involved in both the promotion and inhibition of allergic lung inflammation and airway diseases (51).

After correcting the linear model for cell composition, age, and gender, the outcomes of the preliminary analysis revealed 5,097 DMPs when comparing CR-RSV vs. RW/AS-RSV. The most significant CpG (FDR P-value = 2.77×10^-10; Delta β = -0.14) was cg24509398 within the EYA3 gene, a tyrosine phosphatase involved in DNA repair and in distinguishing apoptotic and repair responses to genotoxicity. Recently, Wang et al. (52) suggested that the over-expression for this gene in patients with pulmonary arterial hypertension (PAH) might stimulate the survival of pulmonary vascular cells in the presence of DNA damage causing vascular remodeling, a typical feature of PAH. Although PAH and asthma are considered two different and unrelated clinical phenomena, there are pieces of evidence demonstrating that they share some pathological features such as inflammation, smooth muscle constriction, and proliferation (53). Although DNA methylation in the CpG selected might not be exerting strong changes in the expression of EYA3, the putative involvement of this gene in the vascular remodeling (a key process in the asthma pathology), might indicate a possible connection between the altered methylation status of the gene and its uncontrolled gene expression.

To test the diagnostic potential of epigenetic markers in the context of RSV sequelae, we carried out a logistic regression analysis using the top 28 DMPs. We found a 3- CpGs position model (related to ZNF2698, LOC102723354 and the overlapping RPL15/NKIRAS1 genes; AUC = 1.00) as the best epigenetic signature to distinguish RW/AS-RSV phenotypes from CR-RSV in both training and test sets as well as in the whole cohort (Figure 7). While LOC102723354 is related to an uncharacterized long non-coding RNA (lncRNA), whose function is still completely unknown, ZNF268 and NKIRAS1 are known to regulate, at a different level and in a different way, the nuclear factor ‘kappa-light-chain-enhancer’ of activated B-cells (NF-κB) (54, 55), an important transcription factor that plays a critical role in the production of many inflammatory cytokines. NF-κB is involved in activating and coordinating both innate and adaptive immune responses; it is found to be associated with allergic airway diseases, and it is activated in bronchial asthmatic patient biopsies and airway epithelium from mice (56).

Among the pathways found to be significantly enriched from the functional analysis, those involved in the cell cycle checkpoint are the most significantly enriched. These finding might suggest that cell proliferation may play an important role in the pathogenesis of wheezing and asthma after RSV infection. Besides, significant enrichment in pathways involved in DNA damage and integrity checkpoint was also found, probably suggesting that the cellular DNA damage machinery is activated and exploited by viruses which have acquired the ability to manipulate the key regulators of these pathways to promote their own replication. Considering the positions residing within promoters, we identified a significant enrichment for the nuclear ubiquitin ligase complex pathway. It has been shown that allergic airway inflammation associated with rhinovirus infection leads to the upregulation of the E3 ubiquitin ligase midline 1 (MID1) in mouse bronchial epithelium, thus, the inhibition of MID1 attenuates rhinovirus-induced airway inflammation and asthma exacerbations (57). In addition, the transforming growth factor-beta (TGF-β) signaling pathway was also found to be significant (FDR P-value = 1.21x10^-2) in promoters’ enrichment pathways. This cytokine’s family seems to play a critical role in the development of airway inflammation and remodeling in asthma (58, 59). In asthmatic airways, TGF-β can induce an antiapoptotic effect in the airway epithelial cells through the SMAD signaling pathway and can trigger apoptosis activating MAPK signaling pathways (60–62). Significant DMPs fall within genes involved in the TFG-β related pathways (e.g. BMP2, FLCN, ING2, PMEPA1, PRDM16, TGFB1I1 and TGFBR3), in line with previous studies that reported an association between these genes and the increased asthma risk (63–68). All over, these results would suggest that the regulation of the TGF-β could be involved in the induction of a moderate or severe course of asthma, a common respiratory sequela observed in children after RSV infection.

In addition, we found a set of co-methylated positions that were significant negatively correlated with the sequelae phenotype, pointing to a hypomethylation state of this CpGs block in RW/AS-RSV patients. Genes associated with these positions pointed out to GTPase signaling as the most significant route in sequelae after RSV infection. Other related significant pathways were those involved in Ca²⁺ metabolism, actin binding and cell junctions. The Rho family of GTPases has been proposed as a promising therapeutic target for asthma, and it is known to play an important role in the pathophysiology of asthma, including airway smooth muscle contraction, airway hyper‐responsiveness and bronchial epithelial barrier dysfunction and recently, in mesenchymal stem cell differentiation and migration for airway remodeling and repairing (69).

A moderate number of statistically significant DMPs emerged when analyzing wheezing and asthma phenotypes separately. The highest difference in Delta β (0.148) corresponds to the cg18873878 observed within the TSS200 region of the TP73 gene being hypermethylated in RW-RSV children in comparison with those with asthma. TP73 gene encodes a member of the p53 family of transcription factors involved in cellular responses to stress and development; its role in airway epithelium is unknown, even if its homolog, TP63, is found to be essential for tracheobronchial epithelium development and homeostasis (70).Another interesting position with a lower methylation pattern in children with asthma was the cg05838113 located within the body of the ADAM8 gene, an ADAM Metallopeptidase Domain which is linked to asthma. It seems that mice with allergic airway inflammation (AAI) show higher levels of expression of ADAM8 in airway epithelium and airway inflammatory cells (71). In addition, increased expression of ADAM8 was observed in the sputum and endobronchial biopsies of patients with moderate and severe asthma (72).

There are a few limitations in the present study. First, the main limitation is the limited number of patients analyzed, and the relative unbalanced in the available samples for the different RSV phenotypes studied. This fact led us to consider a less conservative threshold for the contrast analysis RW-RSV vs. AS-RSV which, in consequence, might increase the probability for false positives. In turns, our study is better powered for the recurrent wheezing outcome. Second, the clinical history of patients is incomplete for some co-variates of interest. Third, the three years of follow up could be not sufficient to discriminate between recurrent wheezing and asthma, being patients with recurrent wheezing at risk to develop asthma in the following years; however, a three-year follow-up should highly capture the diagnosis of asthma in those children with predisposition to develop the disease. Finally, in our cohort, there is higher bacterial superinfection rates in CR-RSV (71.9%) vs. the RSV sequelae cohort (28.9%). The presence of bacteria could have impacted the epigenetics results if the changes observed were due to bacterial presence rather than differences in response to the viral infection. However, it is important to note that bacterial superinfection, as defined in our study, was not always based on a bacterial isolation in a sterile site, but rather frequently based on suggestive clinical symptoms, radiological findings and/or analytical values such as elevated biomarkers (i.e procalcitonin).To account for the mentioned limitations, statistical analyses were carried out considering all possible confounding factors. Even so, further studies will be needed to corroborate the main findings of the present study, as well as other complementary epigenetic studies that takes into account other epigenomic modifications e.g. in the chromatin.

Conclusion

In conclusion, there is suggestive evidence in the present study indicating that epigenetic factors might contribute to the susceptibility to develop recurrent wheezing and asthma after RSV infection. Many DMPs associated with respiratory sequelae developed after RSV infection were detected, with remarkable patterns of methylation profiles observed for genes involved in airway inflammation processes. Functional analysis considering all significant CpGs between sequelae and complete recovery patients revealed significant enrichment for pathways involved in cell cycle checkpoint, DNA damage and integrity checkpoint. In addition, we reported a 3-CpGs epigenetic signature that might be of interest as a diagnostic tool for RSV sequelae. DNA methylation might play a fundamental role in the development of asthma and/or wheezing after RSV infection and could explain the different post-infection sequelae observed. Further investigation using larger cohorts, as well as transcriptomic studies, would be needed to further disentangle the role of epigenomics in asthma and other clinical manifestations (e.g. wheezing), as well as to validate epigenetic biomarkers that allow predicting with precision the different respiratory sequelae that can emerge after infection.

Data Availability Statement

The data presented in the study are deposited on Gene Expression Omnibus (GEO) repository, Accession number: GSE199334.

Ethics Statement

The studies involving human participants were reviewed and approved by the Ethics Committee of Clinical Investigation of Galicia (CEIC 2010/015; updated version 2017/07/26, reg 2017/398). Written informed consent to participate in this study was provided by the participants’ legal guardian/next of kin.

Author Contributions

FM-T, AS, and JG conceived and designed the experiments. IR-C, MC-L and MTG-L were involved in sample recruitment and collection of clinical data. SP, AC, AD-L performed the experiments. SP, and AG-C analyzed the data. SP, AG-C, and AS wrote the first draft of the manuscript, and was contributed by FM-T. All authors read and approved the final manuscript.

Funding

This study received support from Instituto de Salud Carlos III (ISCIII): GePEM (PI16/01478/Cofinanciado FEDER; AS), DIAVIR (DTS19/00049/Cofinanciado FEDER, AS), Resvi-Omics (PI19/01039/Cofinanciado FEDER, AS), Agencia Gallega de Innovación (GAIN): Grupos con Potencial de Crecimiento (IN607B 2020/08, AS); Agencia Gallega para la Gestión del Conocimiento en Salud (ACIS): BI-BACVIR (PRIS-3, AS), and CovidPhy (SA 304 C, AS); ReSVinext (PI16/01569/Cofinanciado FEDER, FM-T), Enterogen (PI19/01090/Cofinanciado FEDER, FM-T) and consorcio Centro de Investigación Biomédica en Red de Enfermedades Respiratorias (CB21/06/00103; FM-T; GEN-COVID (IN845D 2020/23, FM-T) and Grupos de Referencia Competitiva (IIN607A2021/05, FM-T). AD-L is funded by a contract “Juan Rodés” from ISCIII (JR17/00016). ABC is a Miguel Servet researcher (ISCIII; CP17/0008). The funders were not involved in the study design, collection, analysis, interpretation of data, the writing of this article or the decision to submit it for publication.

Conflict of Interest

IR-C has received honoraria from GSK, Pfizer, Sanofi Pasteur and MSD for taking part in advisory boards and expert meetings and for acting as a speaker in congresses outside the scope of the submitted work. IR-C has also acted as subinvestigator in randomized controlled trials of Ablynx, Abbot, Seqirus, Sanofi Pasteur MSD, Sanofi Pasteur, Cubist, Wyeth, Merck, Pfizer, Roche, Regeneron, Jansen, Medimmune, Novavax, Novartis and GSK. FM-T has received honoraria from GSK group of companies, Pfizer Inc, Sanofi Pasteur, MSD, Seqirus, Biofabri and Janssen for taking part in advisory boards and expert meetings and for acting as a speaker in congresses outside the scope of the submitted work. FM-T has also acted as principal investigator in randomized controlled trials of the above-mentioned companies as well as Ablynx, Gilead, Regeneron, Roche, Abbott, Novavax, and MedImmune, with honoraria paid to his institution.

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.

Acknowledgments

We greatly thank all the patients that kindly contributed to this study. The authors thank the study investigators of GENDRES network (www.gendres.org) (Annex/Appendix), and the nursery and laboratory service from the Hospital Clínico Universitario de Santiago de Compostela for their effort.

Supplementary Material

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

Supplementary Figure 1 | (A) Selection of the soft-thresholding power. Plots showed a correlation between soft-thresholding powers and both the scale-free fit index (upper) and the mean connectivity (lower). (B) Clustering dendrogram of CpGs and co-methylation modules detected represented by different colors. (C) Average significance of all CpGs in each module.

Supplementary Figure 2 | (A) β-values heatmap of the CpGs from the pink module. Samples eigengenes values are also represented (B) Differences in samples eigengenes values from pink module between RW/AS-RSV and CR-RSV groups. (C) Correlation between gene significance (GS) for sequelae phenotype and module membership (MM) in the pink module. (D) Over-representation analysis of gene ontology terms using genes in the pink module (BP, biological processes; MF, molecular functions; CC, cellular components).

References

1. Shi T, McAllister DA, O'Brien KL, Simoes EAF, Madhi SA, Gessner BD, et al. Global, Regional, and National Disease Burden Estimates of Acute Lower Respiratory Infections Due to Respiratory Syncytial Virus in Young Children in 2015: A Systematic Review and Modelling Study. Lancet (2017) 390(10098):946–58. doi: 10.1016/S0140-6736(17)30938-8

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Liu L, Oza S, Hogan D, Chu Y, Perin J, Zhu J, et al. Global, Regional, and National Causes of Under-5 Mortality in 2000–15: An Updated Systematic Analysis With Implications for the Sustainable Development Goals. Lancet (2016) 388(10063):3027–35. doi: 10.1016/S0140-6736(16)31593-8

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Mazur NI, Martinón-Torres F, Baraldi E, Fauroux B, Greenough A, Heikkinen T, et al. Lower Respiratory Tract Infection Caused by Respiratory Syncytial Virus: Current Management and New Therapeutics. Lancet Resp Med (2015) 3(11):888–900. doi: 10.1016/S2213-2600(15)00255-6

CrossRef Full Text | Google Scholar

4. Hammitt LL, Dagan R, Yuan Y, Baca Cots M, Bosheva M, Madhi SA, et al. Nirsevimab for Prevention of RSV in Healthy Late-Preterm and Term Infants. (2022) N Engl J Med 386(9):837–46. doi: 10.1056/NEJMoa2110275

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Groothuis JR, Nishida H. Prevention of Respiratory Syncytial Virus Infections in High-Risk Infantsby Monoclonal Antibody (Palivizumab). Pediatr Int (2002) 44(3):235–41. doi: 10.1046/j.1442-200X.2002.01558.x

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Resch B. Product Review on the Monoclonal Antibody Palivizumab for Prevention of Respiratory Syncytial Virus Infection. JHv Immunotherapeutics (2017) 13(9):2138–49. doi: 10.1080/21645515.2017.1337614

CrossRef Full Text | Google Scholar

7. Sigurs N, Gustafsson PM, Bjarnason R, Lundberg F, Schmidt S, Sigurbergsson F, et al. Severe Respiratory Syncytial Virus Bronchiolitis in Infancy and Asthma and Allergy at Age 13. Am J Respir Crit Care Med (2005) 171(2):137–41. doi: 10.1164/rccm.200406-730OC

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Sigurs N, Bjarnason R, Sigurbergsson F, Kjellman B. Respiratory Syncytial Virus Bronchiolitis in Infancy Is an Important Risk Factor for Asthma and Allergy at Age 7. Am J Respir Crit Care Med (2000) 161(5):1501–7. doi: 10.1164/ajrccm.161.5.9906076

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Holtzman MJ, Shornick LP, Grayson MH, Kim EY, Tyner JW, Patel AC, et al. “Hit-And-Run” Effects of Paramyxoviruses as a Basis for Chronic Respiratory Disease. Pediatr Infect Dis J (2004) 23(11):S235–S45. doi: 10.1097/01.inf.0000144674.24802.c1

PubMed Abstract | CrossRef Full Text | Google Scholar

10. James AL, Wenzel S. Clinical Relevance of Airway Remodelling in Airway Diseases. Eur Respir J (2007) 30(1):134–55. doi: 10.1183/09031936.00146905

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Clough JB, Keeping KA, Edwards LC, Freeman WM, Warner JA, Warner JO. Can We Predict Which Wheezy Infants Will Continue to Wheeze? Am J Respir Crit Care Med (1999) 160(5):1473–80. doi: 10.1164/ajrccm.160.5.9807019

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Everard ML. The Relationship Between Respiratory Syncytial Virus Infections and the Development of Wheezing and Asthma in Children. Curr Opin Allergy Clin Immunol (2006) 6(1):56–61. doi: 10.1097/01.all.0000200506.62048.06

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Martinez FD. Development of Wheezing Disorders and Asthma in Preschool Children. Pediatrics (2002) 109(2 Suppl):362–7. doi: 10.1542/peds.109.SE1.362

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Shi T, Ooi Y, Zaw EM, Utjesanovic N, Campbell H, Cunningham S, et al. Association Between Respiratory Syncytial Virus-Associated Acute Lower Respiratory Infection in Early Life and Recurrent Wheeze and Asthma in Later Childhood. J Infect Dis (2019):S628–33. doi: 10.1093/infdis/jiz311

CrossRef Full Text | Google Scholar

15. Busse WW, Banks-Schlegel SP, Larsen GL. Effects of Growth and Development on Lung Function. Models for Study of Childhood Asthma. Am J Respir Crit Care Med (1997) 156(1):314–9. doi: 10.1164/ajrccm.156.1.9612121

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Fonseca W, Lukacs NW, Ptaschinski C. Factors Affecting the Immunity to Respiratory Syncytial Virus: From Epigenetics to Microbiome. Front Immunol (2018) 9:226. doi: 10.3389/fimmu.2018.00226

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Lee GR, Kim ST, Spilianakis CG, Fields PE, Flavell RA. T Helper Cell Differentiation: Regulation by Cis Elements and Epigenetics. Immun (2006) 24(4):369–79. doi: 10.1016/j.immuni.2006.03.007

CrossRef Full Text | Google Scholar

18. Wang L, Wu G, Qin X, Ma Q, Zhou Y, Liu S, et al. Expression of Nodal on Bronchial Epithelial Cells Influenced by Lung Microbes Through DNA Methylation Modulates the Differentiation of T-Helper Cells. Cell Physiol Biochem (2015) 37(5):2012–22. doi: 10.1159/000438561

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Syn G, Blackwell JM, Jamieson SE. Epigenetics in Infectious Diseases. Epigenet Biomarkers Diagnostics (2016) 19:377–400. doi: 10.1016/B978-0-12-801899-6.00019-X

CrossRef Full Text | Google Scholar

20. Gómez-Carballa A, Barral-Arca R, Cebey-López M, Currás-Tuala MJ, Pischedda S, Gómez-Rial J, et al. Host Transcriptomic Response Following Administration of Rotavirus Vaccine in Infants’ Mimics Wild Type Infection. Front Immunol (2021) 11:3517. doi: 10.3389/fimmu.2020.580219

CrossRef Full Text | Google Scholar

21. Chen HS, Lu F, Lieberman PM. Epigenetic Regulation of EBV and KSHV Latency. Curr Opin Virol (2013) 3(3):251–9. doi: 10.1016/j.coviro.2013.03.004

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Paschos K, Allday MJ. Epigenetic Reprogramming of Host Genes in Viral and Microbial Pathogenesis. Trends Microbiol (2010) 18(10):439–47. doi: 10.1016/j.tim.2010.07.003

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Jjingo D, Conley AB, Soojin VY, Lunyak VV, Jordan IKJO. On the Presence and Role of Human Gene-Body DNA Methylation. Oncotarget (2012) 3(4):462. doi: 10.18632/oncotarget.497

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Yang X, Han H, De Carvalho DD, Lay FD, Jones PA, Liang G. Gene Body Methylation Can Alter Gene Expression and Is a Therapeutic Target in Cancer. Cancer Cell (2014) 26(4):577–90. doi: 10.1016/j.ccr.2014.07.028

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Elgizouli M, Logan C, Grychtol R, Rothenbacher D, Nieters A, Heinzmann A. Reduced PRF1 Enhancer Methylation in Children With a History of Severe RSV Bronchiolitis in Infancy: An Association Study. BMC Pediatr (2017) 17(1):65. doi: 10.1186/s12887-017-0817-9

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Elgizouli M, Logan C, Nieters A, Brenner H, Rothenbacher D. Cord Blood PRF1 Methylation Patterns and Risk of Lower Respiratory Tract Infections in Infants: Findings From the Ulm Birth Cohort. Med (Baltimore) (2015) 94(1):e332. doi: 10.1097/MD.0000000000000332

CrossRef Full Text | Google Scholar

27. Ptaschinski C, Mukherjee S, Moore ML, Albert M, Helin K, Kunkel SL, et al. RSV-Induced H3K4 Demethylase KDM5B Leads to Regulation of Dendritic Cell-Derived Innate Cytokines and Exacerbates Pathogenesis In Vivo. PloS Pathog (2015) 11(6):e1004978. doi: 10.1371/journal.ppat.1004978

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Nagata DE, Ting HA, Cavassani KA, Schaller MA, Mukherjee S, Ptaschinski C, et al. Epigenetic Control of Foxp3 by SMYD3 H3K4 Histone Methyltransferase Controls Itreg Development and Regulates Pathogenic T-Cell Responses During Pulmonary Viral Infection. Mucosal Immunol (2015) 8(5):1131–43. doi: 10.1038/mi.2015.4

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Team RC. R: A Language and Environment for Statistical Computing. Austria: Vienna (2013).

Google Scholar

30. Fortin JP, Triche TJ Jr., Hansen KD. Preprocessing, Normalization and Integration of the Illumina HumanMethylationEPIC Array With Minfi. Bioinf (2017) 33(4):558–60. doi: 10.1093/bioinformatics/btw691

CrossRef Full Text | Google Scholar

31. Aryee MJ, Jaffe AE, Corrada-Bravo H, Ladd-Acosta C, Feinberg AP, Hansen KD, et al. Minfi: A Flexible and Comprehensive Bioconductor Package for the Analysis of Infinium DNA Methylation Microarrays. Bioinf (2014) 30(10):1363–9. doi: 10.1093/bioinformatics/btu049

CrossRef Full Text | Google Scholar

32. Bibikova M, Le J, Barnes B, Saedinia-Melnyk S, Zhou L, Shen R, et al. Genome-Wide DNA Methylation Profiling Using Infinium(R) Assay. Epigenomics (2009) 1(1):177–200. doi: 10.2217/epi.09.14

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Bibikova M, Lin Z, Zhou L, Chudin E, Garcia EW, Wu B, et al. High-Throughput DNA Methylation Profiling Using Universal Bead Arrays. Genome Res (2006) 16(3):383–93. doi: 10.1101/gr.4410706

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Du P, Zhang X, Huang CC, Jafari N, Kibbe WA, Hou L, et al. Comparison of Beta-Value and M-Value Methods for Quantifying Methylation Levels by Microarray Analysis. BMC Bioinf (2010) 11:587. doi: 10.1186/1471-2105-11-587

CrossRef Full Text | Google Scholar

35. Ziller MJ, Gu H, Muller F, Donaghey J, Tsai LT, Kohlbacher O, et al. Charting a Dynamic DNA Methylation Landscape of the Human Genome. Nat (2013) 500(7463):477–81. doi: 10.1038/nature12433

CrossRef Full Text | Google Scholar

36. Houseman EA, Accomando WP, Koestler DC, Christensen BC, Marsit CJ, Nelson HH, et al. DNA Methylation Arrays as Surrogate Measures of Cell Mixture Distribution. BMC Bioinf (2012) 13:86. doi: 10.1186/1471-2105-13-86

CrossRef Full Text | Google Scholar

37. Touleimat N, Tost J. Complete Pipeline for Infinium® Human Methylation 450k BeadChip Data Processing Using Subset Quantile Normalization for Accurate DNA Methylation Estimation. Epigenomics (2012) 4(3):325–41. doi: 10.2217/epi.12.21

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, et al. Limma Powers Differential Expression Analyses for RNA-Sequencing and Microarray Studies. Nucleic Acids Res (2015) 43(7):e47. doi: 10.1093/nar/gkv007

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Benjamini Y, Hochberg Y. Controlling the False Discovery Rate: A Practical and Powerful Approach to Multiple Testing. J R Stat Soc: Ser B (Methodol) (1995) 57(1):289–300. doi: 10.1111/j.2517-6161.1995.tb02031.x

CrossRef Full Text | Google Scholar

40. Hoggart CJJB. PRReMS: Parallel Regularised Regression Model Search for Bio-Signature Discovery. bioRxiv (2018) 355479. doi: 10.1101/355479

CrossRef Full Text | Google Scholar

41. Ren X, Kuan PF. methylGSA: A Bioconductor Package and Shiny App for DNA Methylation Data Length Bias Adjustment in Gene Set Testing. Bioinf (2019) 35(11):1958–9. doi: 10.1093/bioinformatics/bty892

CrossRef Full Text | Google Scholar

42. Langfelder P, Horvath S. WGCNA: An R Package for Weighted Correlation Network Analysis. BMC Bioinform (2008) 9(1):1–13. doi: 10.1186/1471-2105-9-559

CrossRef Full Text | Google Scholar

43. Wu T, Hu E, Xu S, Chen M, Guo P, Dai Z, et al. Clusterprofiler 4.0: A Universal Enrichment Tool for Interpreting Omics Data. The Innovation (2021) 2(3):100141. doi: 10.1016/j.xinn.2021.100141

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Hansen K, Aryee M. IlluminaHumanMethylationEPICmanifest: Manifest for Illumina's EPIC Methylation Arrays. R package version 03 0. (2016).

Google Scholar

45. Coutts J, Fullarton J, Morris C, Grubb E, Buchan S, Rodgers-Gray B, et al. Association Between Respiratory Syncytial Virus Hospitalization in Infancy and Childhood Asthma. Pediatr Pulmonol (2020) 55(5):1104–10. doi: 10.1002/ppul.24676

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Mejias A, Wu B, Tandon N, Chow W, Varma R, Franco E, et al. Risk of Childhood Wheeze and Asthma After Respiratory Syncytial Virus Infection in Full-Term Infants. Pediatr Allergy Immunol (2020) 31(1):47–56. doi: 10.1111/pai.13131

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Koponen P, Helminen M, Paassilta M, Luukkaala T, Korppi MJERJ. Preschool Asthma After Bronchiolitis in Infancy. Eur Respir J (2012) 39(1):76–80. doi: 10.1183/09031936.00040211

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Li F, Zhu H, Sun R, Wei H, Tian Z. Natural Killer Cells Are Involved in Acute Lung Immune Injury Caused by Respiratory Syncytial Virus Infection. J Virol (2012) 86(4):2251–8. doi: 10.1128/JVI.06209-11

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Schmidt ME, Varga SM. The CD8 T cell response to respiratory virus infections. Front Immunol (2018) 9:678. doi: 10.3389/fimmu.2018.00678

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Hamelmann E, Oshiba A, Paluh J, Bradley K, Loader J, Potter TA, et al. Requirement for CD8+ T Cells in the Development of Airway Hyperresponsiveness in a Marine Model of Airway Sensitization. J Exp Med (1996) 183(4):1719–29. doi: 10.1084/jem.183.4.1719

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Kaiko GE, Phipps S, Angkasekwinai P, Dong C, Foster PS. NK Cell Deficiency Predisposes to Viral-Induced Th2-Type Allergic Inflammation via Epithelial-Derived IL-25. J Immunol (2010) 185(8):4681–90. doi: 10.4049/jimmunol.1001758

PubMed Abstract | CrossRef Full Text | Google Scholar

52. Wang Y, Pandey RN, York AJ, Mallela J, Nichols WC, Hu YC, et al. The EYA3 Tyrosine Phosphatase Activity Promotes Pulmonary Vascular Remodeling in Pulmonary Arterial Hypertension. Nat Commun (2019) 10(1):4143. doi: 10.1038/s41467-019-12226-1

PubMed Abstract | CrossRef Full Text | Google Scholar

53. Said SI, Hamidi SA, Gonzalez Bosc L. Asthma and Pulmonary Arterial Hypertension: Do They Share a Key Mechanism of Pathogenesis? Eur Respir J (2010) 35(4):730–4. doi: 10.1183/09031936.00097109

PubMed Abstract | CrossRef Full Text | Google Scholar

54. Wang W, Guo M, Hu L, Cai J, Zeng Y, Luo J, et al. The Zinc Finger Protein ZNF268 Is Overexpressed in Human Cervical Cancer and Contributes to Tumorigenesis via Enhancing NF-κb Signaling. J Biol Chem (2012) 287(51):42856–66. doi: 10.1074/jbc.M112.399923

PubMed Abstract | CrossRef Full Text | Google Scholar

55. Sarais F, Rebl H, Verleih M, Ostermann S, Krasnov A, Köllner B, et al. Characterisation of the Teleostean κb-Ras Family: The Two Members NKIRAS1 and NKIRAS2 From Rainbow Trout Influence the Activity of NF-κb in Opposite Ways. Fish Shellfish Immunol (2020) 106:1004–13. doi: 10.1016/j.fsi.2020.08.052

PubMed Abstract | CrossRef Full Text | Google Scholar

56. Rico-Rosillo G. Vega-Robledo GBJRAdM. The Involvement of NF-κb Transcription Factor in Asthma. Rev Alerg Mex (2011) 58(2):107–11.

PubMed Abstract | Google Scholar

57. Collison A, Hatchwell L, Verrills N, Wark PA, De Siqueira AP, Tooze M, et al. The E3 Ubiquitin Ligase Midline 1 Promotes Allergen and Rhinovirus-Induced Asthma by Inhibiting Protein Phosphatase 2A Activity. Nat Med (2013) 19(2):232–7. doi: 10.1038/nm.3049

PubMed Abstract | CrossRef Full Text | Google Scholar

58. Ohno I, Nitta Y, Yamauchi K, Hoshi H, Honma M, Woolley K, et al. Transforming Growth Factor Beta 1 (TGF Beta 1) Gene Expression by Eosinophils in Asthmatic Airway Inflammation. Am J Respir Cell Mol Biol (1996) 15(3):404–9. doi: 10.1165/ajrcmb.15.3.8810646

PubMed Abstract | CrossRef Full Text | Google Scholar

59. Gagliardo R, Chanez P, Gjomarkaj M, La Grutta S, Bonanno A, Montalbano AM, et al. The Role of Transforming Growth Factor-Beta1 in Airway Inflammation of Childhood Asthma. Int J Immunopathol Pharmacol (2013) 26(3):725–38. doi: 10.1177/039463201302600316

PubMed Abstract | CrossRef Full Text | Google Scholar

60. Halwani R, Al-Muhsen S, Al-Jahdali H, Hamid Q. Role of Transforming Growth Factor-Beta in Airway Remodeling in Asthma. Am J Respir Cell Mol Biol (2011) 44(2):127–33. doi: 10.1165/rcmb.2010-0027TR

PubMed Abstract | CrossRef Full Text | Google Scholar

61. Tirado-Rodriguez B, Ortega E, Segura-Medina P, Huerta-Yepez S. TGF- Beta: An Important Mediator of Allergic Disease and a Molecule With Dual Activity in Cancer Development. J Immunol Res (2014) 2014:318481. doi: 10.1155/2014/318481

PubMed Abstract | CrossRef Full Text | Google Scholar

62. Aschner Y, Downey GP. Transforming Growth Factor-Beta: Master Regulator of the Respiratory System in Health and Disease. Am J Respir Cell Mol Biol (2016) 54(5):647–55. doi: 10.1165/rcmb.2015-0391TR

PubMed Abstract | CrossRef Full Text | Google Scholar

63. DeVries A, Wlasiuk G, Miller SJ, Bosco A, Stern DA, Lohman IC, et al. Epigenome-Wide Analysis Links SMAD3 Methylation at Birth to Asthma in Children of Asthmatic Mothers. J Allergy Clin Immunol (2017) 140(2):534–42. doi: 10.1016/j.jaci.2016.10.041

PubMed Abstract | CrossRef Full Text | Google Scholar

64. Lund RJ, Osmala M, Malonzo M, Lukkarinen M, Leino A, Salmi J, et al. Atopic Asthma After Rhinovirus-Induced Wheezing Is Associated With DNA Methylation Change in the SMAD3 Gene Promoter. Allergy (2018) 73(8):1735–40. doi: 10.1111/all.13473

PubMed Abstract | CrossRef Full Text | Google Scholar

65. Nakao A, Sagara H, Setoguchi Y, Okada T, Okumura K, Ogawa H, et al. Expression of Smad7 in Bronchial Epithelial Cells Is Inversely Correlated to Basement Membrane Thickness and Airway Hyperresponsiveness in Patients With Asthma. J Allergy Clin Immunol (2002) 110(6):873–8. doi: 10.1067/mai.2002.129236

PubMed Abstract | CrossRef Full Text | Google Scholar

66. Kim H-K, Jang T-W, Jung M-H, Park H-W, Lee J-E, Shin E-S, et al. Association Between Genetic Variations of the Transforming Growth Factor β Receptor Type III and Asthma in a Korean Population. Exp Mol Med (2010) 42(6):420–7. doi: 10.3858/emm.2010.42.6.043

PubMed Abstract | CrossRef Full Text | Google Scholar

67. Kariyawasam HH, Xanthou G, Barkans J, Aizen M, Kay AB, Robinson DS, et al. Basal Expression of Bone Morphogenetic Protein Receptor Is Reduced in Mild Asthma. Am J Respir Crit Care Med (2008) 177(10):1074–81. doi: 10.1164/rccm.200709-1376OC

PubMed Abstract | CrossRef Full Text | Google Scholar

68. Trinh HKT, Pham DL, Choi Y, Kim HM, Kim SH, Park HSJC, et al. Epithelial Folliculin Enhances Airway Inflammation in Aspirin-Exacerbated Respiratory Disease. Clin Exp Allergy (2018) 48(11):1464–73. doi: 10.1111/cea.13253

PubMed Abstract | CrossRef Full Text | Google Scholar

69. Zhang Y, Saradna A, Ratan R, Ke X, Tu W, Do DC, et al. RhoA/Rho-Kinases in Asthma: From Pathogenesis to Therapeutic Targets. Clin Transl Immunol (2020) 9(5):e1134. doi: 10.1002/cti2.1134

CrossRef Full Text | Google Scholar

70. Hogan BL, Barkauskas CE, Chapman HA, Epstein JA, Jain R, Hsia CC, et al. Repair and Regeneration of the Respiratory System: Complexity, Plasticity, and Mechanisms of Lung Stem Cell Function. Cell Stem Cell (2014) 15(2):123–38. doi: 10.1016/j.stem.2014.07.012

PubMed Abstract | CrossRef Full Text | Google Scholar

71. Knolle MD, Owen CA. ADAM8: A New Therapeutic Target for Asthma. Expert Opinion on Therapeutic Targets. Expert Opin Ther Targets (2009) 13(5):523–40. doi: 10.1517/14728220902889788

PubMed Abstract | CrossRef Full Text | Google Scholar

72. Foley SC, Mogas AK, Olivenstein R, Fiset PO, Chakir J, Bourbeau J, et al. Increased Expression of ADAM33 and ADAM8 With Disease Progression in Asthma. J Allergy Clin Immunol (2007) 119(4):863–71. doi: 10.1016/j.jaci.2006.12.665

PubMed Abstract | CrossRef Full Text | Google Scholar