Abstract

   The human nasal epithelium is the primary site of exposure to influenza
   virus, the initiator of host responses to influenza and the resultant
   pathologies. Influenza virus may cause serious respiratory infection
   resulting in major complications, as well as severe impairment of the
   airways. Here, we elucidated the global transcriptomic changes during
   H3N2 infection of human nasal epithelial cells from multiple
   individuals. Using RNA sequencing, we characterized the
   differentially-expressed genes and pathways associated with changes
   occurring at the nasal epithelium following infection. We used in vitro
   differentiated human nasal epithelial cell culture model derived from
   seven different donors who had no concurrent history of viral
   infections. Statistical analysis highlighted strong transcriptomic
   signatures significantly associated with 24 and 48 h after infection,
   but not at the earlier 8-h time point. In particular, we found that the
   influenza infection induced in the nasal epithelium early and altered
   responses in interferon gamma signaling, B-cell signaling, apoptosis,
   necrosis, smooth muscle proliferation, and metabolic alterations. These
   molecular events initiated at the infected nasal epithelium may
   potentially adversely impact the airway, and thus the genes we
   identified could serve as potential diagnostic biomarkers or
   therapeutic targets for influenza infection and associated disease
   management.

   Keywords: influenza, H3N2 virus, human nasal epithelial cells (hNECs),
   RNA sequencing, transcriptomics, epithelial responses, pathogenesis

1. Introduction

   The global burden of inter-pandemic influenza is high. It is estimated
   to affect 1 billion people annually, with 3–5 million severe cases
   requiring hospitalization or intensive care treatment, resulting in
   approximately 0.5 million deaths [[52]1]. Worryingly, drug-resistant
   influenza strains are emerging at a rapid rate that will severely
   hamper the ability of our healthcare systems to contain influenza
   outbreaks [[53]2]. Therefore, alternative strategies are needed against
   severe influenza infections during both seasonal and pandemic influenza
   outbreaks.

   The normal human airway epithelium is a pseudo-stratified layer of
   ciliated and non-ciliated columnar cells, goblet cells, club cells, and
   basal cells [[54]3]. The airway epithelium protects against airway
   infection via efficient mucociliary clearance (MCC), the production of
   inflammatory mediators and chemokines against viruses, and the
   recruitment of immune cells [[55]4]. When the influenza virus breaches
   the defense of the human airway epithelium, it causes a myriad of
   innate responses by the infected host in response to viral invasion
   [[56]5,[57]6]. Among these changes are critical factors that can
   determine disease severity, and which may lead to the development of
   diagnostic, prognostic prediction markers, or anti-influenza therapies
   [[58]7,[59]8,[60]9,[61]10]. However, few studies have hitherto been
   performed in relevant models [[62]11] and human models of influenza are
   not feasible due to potential severity of the infection. Therefore, the
   mechanistic study of viral-induced airway changes using relevant models
   can lead to better understanding of the development of severe
   complications. Additionally, we need greater clarity on the different
   immune responses in view of the rising prevalence of chronic diseases
   such as diabetes mellitus and asthma. Patients with these disorders are
   especially susceptible to severe influenza complications compared to
   healthy subjects [[63]12]. Thus, the establishment of a baseline
   response against influenza infection of “healthy” tissue is beneficial
   to facilitate future comparative studies to better manage influenza in
   patients with co-morbidities.

   Although the study of host responses in influenza infection is not new,
   current in vitro cell lines cannot accurately represent human airway
   infection due to the lack of key mucociliary features [[64]13]. Hence,
   we have previously developed an air–liquid interface (ALI) human airway
   epithelial cell culture that is able to sustain influenza infection
   [[65]5,[66]6]. We have also further compared the transcriptomic
   responses of our infected human nasal epithelial cells (hNECs) with 15
   other in vitro and in vivo influenza infection transcriptomic studies
   [[67]6]. The comparison revealed that at their peak responses against
   influenza, the differential transcriptome signature in hNECs was highly
   similar to the signatures from other influenza infection models
   [[68]6]. Interestingly, compared to the homogenous cell lines tested,
   our heterogenous hNEC model exhibited a more comparable response to the
   clinical influenza studies, indicating that most responses were
   initiated at the nasal epithelium [[69]6]. Therefore, in this study, we
   aim to further utilize the hNEC model as a physiologically relevant in
   vitro model to clarify the nasal epithelial responses against influenza
   H3N2 infection, which would then facilitate the identification of the
   key host factors that are significant for future studies.

   To establish host factors that are significantly altered in the nasal
   epithelium as a reference of early innate responses against influenza,
   the dynamic expression of the genes needs to be clearly elucidated.
   While there are many studies that utilize microarray analysis to
   identify the host responses against influenza, the limitation of the
   microarray is its inability to determine the full extent of gene
   changes due to its hybridization-based protocol [[70]14]. The aim of
   this study was to utilize RNA sequencing (RNAseq) technology to not
   only reveal the hNEC responses (from multiple individuals) against
   influenza infection, but also to identify those genes with high
   magnitude changes to serve as potential reference markers of the innate
   responses of influenza infection. Given that RNAseq functions by
   reading virtually all the RNAs present in the samples tested, we can
   also discern the magnitude of each RNA change and mark them as the
   canonical responses. In addition, as RNAseq is not constrained by probe
   usage as in microarrays, they are therefore more reliable in detecting
   novel interactions during influenza infections of hNECs. Hence, RNAseq
   analysis will further augment the transcriptomic data established
   previously by microarray analysis. The augmented baseline can then be
   applied to future clinical studies and practice against influenza
   infection, especially for comparison against patients with other
   underlying co-morbidities that may be affected by more severe disease.

2. Materials and Methods

2.1. Derivation of Human Nasal Epithelial Stem/Progenitor Cells (hNESPCs) and
In Vitro Differentiation of hNECs

   Approval to conduct this study was obtained from the National
   Healthcare Group Domain-Specific Board of Singapore (DSRB Ref:
   D/11/228) and the institutional review board of the National University
   of Singapore (IRB Ref: 13-509). Written consent was obtained from
   donors prior to the collection of the tissue biopsies. At the time of
   collection, all subjects were free of symptoms of URTI. The medical
   backgrounds of the subjects are summarized in[71]Table S1. The hNESPCs
   were isolated and enriched from the tissue biopsies according to a
   previously standardized protocol [[72]5,[73]15], which normalized the
   hNESPCs to a baseline state that differentiates into hNECs resembling
   healthy tissues if they pass the quality control checks for their
   differentiation [[74]5]. Following enrichment, the hNESPCs were
   expanded further and subjected to ALI culture in Transwells for in
   vitro differentiation according to previous protocol as well
   [[75]5,[76]15]. Briefly, primary cells were subjected to isolation for
   selection of hNESPCs, which were enriched and expanded with Dulbecco’s
   Modified Eagle Medium: Nutrient Mixture F-12 (DMEM/F12)
   (Gibco-Invitrogen, Carlsbad, CA, USA) containing 10 ng/mL of human
   epithelial growth factor (EGF, Gibco-Invitrogen, Carlsbad, CA, USA), 5
   μg/mL of insulin (Sigma, St. Louis, MO, USA), 0.1 nM of cholera toxin
   (Sigma, St. Louis, MO, USA), 0.5 μg/mL of hydrocortisone (Sigma, St.
   Louis, MO, USA), 2 nM/mL of 3,3′,5-triiodo-l-thyronine (T3) (Sigma, St.
   Louis, MO, USA), 10 μL/mL of an N-2 supplement (Gibco-Invitrogen) and
   100 IU/mL of antibiotic-Antimycotic (Gibco-Invitrogen, Carlsbad, CA,
   USA). The expanded hNESPCs were then transferred onto 12-well 0.4 μm
   Transwell inserts (Corning, Corning, NY, USA). Once confluent, growth
   medium was discarded and 700 μL of PneumaCult™-ALI Medium with inducer
   supplements (STEMCELL Technologies Inc., Vancouver, Canada) was added
   to the basal chamber to establish ALI conditions. The cells were
   cultured in ALI culture for 4 weeks, with media change every 2–3 days.
   After 3–4 weeks of differentiation, hNECs from a total of seven donors
   were then subjected to influenza H3N2 virus infection.

2.2. Inoculation of Human Influenza H3N2 Virus in Fully Differentiated hNECs
and Viral Quantification

   The Influenza A strain used in this study is of the H3N2 subtype
   (A/Aichi/2/1968) (ATCC, Manassas, VA, USA). The virus was propagated
   using embryonated egg culture and used for all the infection in the
   hNECs. Prior to infection, fully differentiated hNECs were washed with
   1× dPBS and infected with the H3N2 influenza virus at a multiplicity of
   infection (MOI) of 0.1 and incubated for 1 h at 35 °C. After the 1 h
   incubation, the viral inoculum was removed and the hNECs were incubated
   back in 35 °C. The control hNECs were harvested for apical wash and RNA
   prior to the infection at 0 h post-infection (hpi). The infected hNECs
   were then harvested for the apical wash and RNA following 8, 24, and 48
   hpi incubation at 35 °C.

2.3. Viral Plaque Assay

   At each infection time point, 150 µL of 1x dPBS was added and incubated
   in the apical chamber for 10 min at 35 °C to recover progeny viruses as
   the apical wash. The plaque assay for viral quantification was
   performed using overnight MDCK cultures (ATCC, Manassas, VA, USA) at
   85–95% confluence in 24-well plates. The MDCK cells were incubated with
   100 µL of serial dilutions (from 10^−1 to 10^−6) of virus from apical
   washes at 35 °C for 1 h, where plates were rocked every 15 min to
   ensure equal viral distribution. After incubation, the inocula were
   removed and replaced with 1 mL of Avicel (FMC Biopolymer, Philadelphia,
   PA, USA) overlay, and incubated at 35 °C for 65–72 h. After incubation,
   Avicel overlay were removed, and cells were fixed with 4% formaldehyde
   in 1× PBS for 1 h. Formaldehyde was then removed, and cells were washed
   with 1× PBS prior to staining with 1% crystal violet for 15 min before
   washing the stain away. The plaque-forming units (PFU) were calculated
   as follows: Number of plaques × dilution factor = number of PFU per 100
   µL.

2.4. Total RNA Extraction and Real-Time Quantitative PCR

   At each time point after the collection of apical wash, the hNECs were
   lysed using RNA lysis buffer. Total RNA was then extracted from the
   lysate using mirVana miRNA isolation kit (Life Technologies, Grand
   Island, NY, USA). The extracted total RNA was first subjected to
   nanodrop analysis to first ensure the RNA quality, before being
   submitted for RNAseq analysis. Then, 500 ng from the remaining RNA was
   subjected to cDNA synthesis using Maxima first-strand cDNA synthesis
   kit (ThermoScientific, Pittsburgh PA, USA). After this, qPCR analysis
   was performed to evaluate the transcriptional levels of host response
   genes selected based on previous microarray analysis using pre-designed
   primers (Sigma Aldrich). Each qPCR reaction was performed in duplicate
   using GoTaq-qPCR Master Mix kit (Promega, San Luis Obispo, CA, USA),
   and relative gene expression was calculated using the comparative
   method of 2-ΔΔCt normalized to the housekeeping gene PGK1. Relative
   gene expression levels were presented as median values and
   interquartile ranges, while statistical significance was determined
   using the Wilcoxon signed-rank test.

2.5. Library Preparation for RNAseq

   All human RNAs were analyzed on the Agilent Bioanalyzer (Agilent, Santa
   Clara, CA, USA) or the Perkin Elmer Labchip GX system (Perkin Elmer,
   Waltham, MA, USA) for quality assessment with RNA Integrity Number
   (RIN) or RNA Quality score range from 6.8–9.7 and median of 9.0. cDNA
   libraries were prepared using 2 ng of total RNA and 1 µL of a 1:50,000
   dilution of ERCC RNA Spike in Controls (Ambion^® Thermo Fisher
   Scientific, Waltham, MA, USA) using SMARTSeq v2 protocol [[77]16],
   except for the following modifications: (1) Use of 20 µM TSO, and (2)
   use of 250 pg of cDNA with 1/5 reaction of Illumina Nextera XT kit
   (Illumina, San Diego, CA, USA). The length distribution of the cDNA
   libraries was monitored using DNA High Sensitivity Reagent Kit on the
   Perkin Elmer Labchip GX system (Perkin Elmer, Waltham, MA, USA).
   Sixteen samples were subjected to an indexed PE sequencing run of 2 ×
   51 cycles on an Illumina HiSeq 2000 (16 samples/lane), and 65 samples
   to an indexed PE sequencing run of 2 ×151 cycles on an Illumina HiSeq
   4000 (30 samples/lane).

2.6. RNAseq Analysis

   FASTQ files were mapped to the human genome build hg38 using STAR. Gene
   counts were computed using featureCounts (part of the Subread package)
   using annotations from GENCODE version 26. Differential gene expression
   analysis was performed using edgeR in a paired fashion under R version
   3.3.3. Multiple testing correction was done using the method of
   Benjamini and Hochberg and p-values (False Discovery Rate; FDR) less
   than 0.05 was deemed to be significant.

2.7. Geneset Enrichment Analysis

   Geneset enrichment analysis using data from Gene Ontology (GO) was
   performed using the Bioconductor package topGO, while the analysis
   using Reactome Pathway was performed using the Vioconductor package
   ReactomePA. Both analyses were run in R version 3.3.3 using multiple
   testing-corrected significant differentially-expressed genes.

2.8. Luminex Assay for Cytokines and Chemokines

   Luminex assay was performed on apical supernatant collected from both
   uninfected control and infected hNECs at 8, 24, and 48 hpi using the
   Multiplexing Immunoassay Human Cytokine/Chemokine Panel 1, 41-plex
   (Merck, cat# HCYTOMAG-60K-41, Kenilworth, NJ, USA). The panel included
   the detection of: Human sCD40L, EGF, FGF-2, Flt-3 ligand, Fractalkine,
   G-CSF, GM-CSF, GRO, IFN-α2, IFN-γ, IL-1α, IL-1β, IL-1ra, IL-2, IL-3,
   IL-4, IL-5, IL-6, IL-7, IL-8, IL-9, IL-10, IL-12 (p40), IL-12 (p70),
   IL-13, IL-15, IL-17A, IP-10, MCP-1, MCP-3, MDC (CCL22), MIP-1α, MIP-1β,
   PDGF-AB/BB, RANTES, TGF-α, TNF-α, TNF-β, VEGF, Eotaxin/CCL11, and
   PDGF-AA. Samples and standards were incubated with fluorescent-coded
   magnetic beads which had been pre-coated with the respective capture
   antibodies. After an overnight incubation at 4 °C, plates were washed
   twice. Biotinylated detection antibodies were incubated with the
   complex for 1 h, and Streptavidin-PE was then added and incubated for
   another 30 min. Plates were washed twice again, then beads were
   re-suspended with sheath fluid before acquiring on the FLEXMAP^® 3D
   (Luminex) using xPONENT^® 4.0 (Luminex) acquisition software. Data
   analysis was done on Bio-Plex Manager™ 6.1.1 (Bio-Rad). Standard curves
   were generated with a 5-PL (5-parameter logistic) algorithm, reporting
   values for both MFI and concentration data. Results were then expressed
   as mean fold change compared with uninfected control, and p-values
   (FDR) of less than 0.05 were considered significant.

3. Results

3.1. Influenza Virus Infection Induces Similar Responses in hNECs from
Multiple Donors

   Prior to analysis, the responses of all seven hNECs donors following
   influenza infection were plotted on a principal component analysis
   (PCA) plot. The PCA plot indicated a degree of variability in the
   responses between donors and time points ([78]Figure 1). Nonetheless,
   the responses were clustered tightly enough following infection to
   signify their consistency of infection for further transcriptomic
   analysis—similar to those observed in our previous study [[79]6].

Figure 1.

   [80]Figure 1
   [81]Open in a new tab

   Principal component analysis (PCA) of hNECs from seven donors infected
   with influenza H3N2 virus. PCA plot of transcriptomic changes of
   different donor hNECs infected with influenza virus across four time
   points (0, 8, 24, and 48 h post-infection (hpi)). Each colored dot
   represents the sequencing ID of infected hNECs from a specific donor at
   a specific time point of infection. The majority of dots clustered
   closely following infection, thus signifying the consistency of
   infection for transcriptomic analysis.

3.2. Transcriptomic Changes Are Detected at 8 hpi, and Peak at 48 hpi

   Significant gene expression changes (FDR < 0.05) of infected hNECs were
   detected as early as 8 hpi, and further increased at 24 and 48 hpi
   ([82]Table 1; [83]Figure 2A). Also, the number of genes decreased in a
   linear fashion as the fold change in expression increased, as seen in
   the 10x fold change genes indicated in [84]Figure 2B, where about 10%
   of the significantly altered genes remained. At 8 hpi, there were 31
   upregulated genes and 13 downregulated genes. The major upregulated
   genes were the antiviral sensors and early response genes such as IFNs,
   IFITs, and IFIs. Interestingly enough, interferon lambda (IFNλ) gene
   IFNLs was the earliest response interferon of infected hNECs, as
   opposed to interferons alpha or beta, at 8 hpi. At later time points,
   the number of gene expression changes increased substantially, with
   upregulation of 704 and 1080 genes, and downregulation of 217 and 758
   genes at 24 and 48 hpi, respectively. There was augmented expression of
   antiviral effectors and inflammatory genes at both time points. IFNL
   remained the interferon gene with highest expression at both time
   points, while a marked elevation of cytokines such as CXCL10 and CXCL11
   was also observed. Considering downregulated genes, proliferative and
   transcriptomic functions appeared to be suppressed, with diminished
   expression of genes such as FMO2, KLK12, and FOSB. Genes associated
   with metabolism, cell cycle, and DNA repair were further suppressed
   following infection at 24 and 48 hpi. [85]Tables S2–S4 list the
   complete set of significant gene expression changes, arranged according
   to their fold change (Log[2]FC). In addition, we have also verified
   that the genes showing major expression changes by RNAseq generally
   concurred with RT-qPCR analyses. Of the 10 genes tested by qPCR at 48
   hpi, all of them showed the same directional changes in expression as
   observed by RNAseq. Hence, seven of these genes showed a p-value of
   <0.05 (IL4I1, IFNL1(IL29), CXCL10, TNFSF10, IFI6, CCL24, and CYP26A1),
   one gene had a p-value of <0.1 (CTGF), while only two genes were not
   statistically significant (TGFA and ANO5) ([86]Figure S1).

Table 1.

   Summary of RNAseq analysis across time points.
   Time, Hours Post Infection (hpi) 8 24 48
   Number of significant genes 44 (31↑; 13↓) 921 (704↑; 217↓) 1838 (1080↑;
   758↓)
   Number of significant genes over 10-fold (>Log[2]FC of 3.4) 3 (2↑; 1↓)
   99 (91↑; 8↓) 167 (112↑; 55↓)
   Highest upregulated gene IFNL2 IFNL1 IFNL1
   Highest downregulated gene KLK12 NAP1L3 CYP26A1
   [87]Open in a new tab

Figure 2.

   [88]Figure 2
   [89]Open in a new tab

   RNAseq analysis of significant gene expression changes following
   influenza infection. (A) Significant gene expression changes over time
   (FDR < 0.05). Red bars indicate the number of upregulated genes, while
   the blue bars indicate the number of downregulated genes. The number of
   significant gene expression changes peaked at 48 hpi. (B) Volcano plot
   showing individual points of gene expression changes over time (8, 24,
   and 48 hpi). The number of genes decreased in a linear fashion as the
   fold change increased. Genes with >10-fold difference in expression are
   depicted in dark red, and shown to include ~10% of the total number of
   significant gene changes.

3.3. Transcriptomic Change Alterations Correlate with Viral Titer, with
Intersecting Genes Remaining Consistent in Directional Changes of Expression

   We then further compared the transcriptomic alterations in the hNECs
   over time, following influenza H3N2 infection. The number of gene
   expression changes mirrored the viral titer changes, which peaked at 48
   hpi, and were consistent between donors ([90]Figure 3A). Approximately
   two thirds of genes at 8 and 24 hpi overlapped with other time points,
   while about one third of genes at 48 hpi overlapped ([91]Figure 3B).
   The overlapping genes displayed similar directional consistency at the
   significant time points. In addition, congruent with the consistent
   viral titer with most gene expression changes at 48 hpi, we also noted
   the most consistent alterations in expression of genes across donors.
   This is highlighted in [92]Figure 3C, which portrays the heatmaps of
   the top 10 genes with the smallest p-value, together with their
   direction and magnitude of change. Based on these analyses, we proposed
   that 48 hpi represents the optimal time point for the subsequent
   pathway analysis to ascertain influenza-specific pathway changes.

Figure 3.

   [93]Figure 3
   [94]Open in a new tab

   Viral titer and the transcriptomic change pattern post influenza H3N2
   infection. (A) Viral titer retrieved from apical wash of infected hNECs
   8, 24, and 48 hpi of influenza inoculation. The y-axis is presented in
   a Log[10] scale and viral titer is presented in PFU/100 µL. Each
   individual point at each time point represents viral titer retrieved
   from in vitro-differentiated culture of a hNECs donor. (B) Venn diagram
   of transcriptional change pattern over time post influenza infection.
   Venn diagram summarizes the transcriptional pattern changes over time.
   For 8 and 24 hpi, about two thirds of genes overlapped with other time
   points, while for 48 hpi, one third of genes overlapped. The direction
   of the gene change in overlapping regions were consistent across the
   time points it is significant in. A greater number of genes were
   upregulated than downregulated at all time points. (C) Heatmaps of top
   10 genes with smallest p-value at time points 8, 24 and 48 hpi (blue
   bar) compared to 0 hpi (red bar). The changes of the genes became more
   consistent overtime and showed same directional changes and clustering
   at 48 hpi.

3.4. Pathway Enrichment of Significant Gene Changes Revealed Potential
Epithelium-Initiated Pathways Contributing to Influenza Pathology and
Pathogenesis

   We then further subjected the significant gene changes to gene set
   enrichment using both GO and REACTOME databases. At time points 8, 24,
   and 48 hpi, there were 3, 41, and 30 significant (adjusted p-value <
   0.05) GO biological processes ([95]Table S5) and 3, 93, and 74
   significant (adjusted p-value < 0.05) REACTOME pathways ([96]Table S6),
   respectively. At the early time of 8 hpi, interferon-mediated antiviral
   responses were elevated as expected. At 48 hpi, the pathways appeared
   to be more stabilized and consistent for both GO and REACTOME analyses,
   despite displaying more gene expression changes at this time point.

   Responses to influenza virus skewing towards type I immunity were
   predominant in the GO analysis. The expected interferon-mediated
   functions by the epithelium validated the authenticity of our model,
   where we found enriched type I interferon (GO) and RIG-I (REACTOME)
   pathways with upregulation of nearly all significant gene members (data
   not shown). Besides the interferon and antiviral pathways, we
   identified several functions of interest initiated by the nasal
   epithelium that may contribute to the pathology and pathogenesis of
   influenza. At 48 hpi, GO pathway enrichment analysis revealed that the
   nasal epithelium was actively involved in initial IFNγ signaling
   (GO:0060333), despite not directly producing IFNγ. We also observed
   enriched function in apoptosis and necroptosis (GO:0008637 and
   GO:0060544), immune evasion (GO:0045824), and other pathways that may
   lead to complication events such as smooth muscle proliferation
   (GO:0048661) and response to fatty acid (GO:0070542) ([97]Table 2).

Table 2.

   48 hpi GO pathways potentially associated with influenza pathogenesis.
   GOid Description Annotated Significant Expected p-Value_adj Genes
   GO:0060333 Interferon-gamma-mediated signaling pathway 72 36 10.35
   9.184 × 10^−11 TRIM22 ↑; IRF9 ↑; TRIM38 ↑; GBP1 ↑; GBP2 ↑; HCK ↓; HLA-A
   ↑; HLA-B ↑; HLA-C ↑; HLA-E ↑; HLA-F ↑; HLA-G ↑; ICAM1 ↑; IRF1 ↑; IRF2
   ↑; IRF6 ↑; IRF7 ↑; JAK2 ↑; MT2A ↑; OAS1 ↑; OAS2 ↑; OAS3 ↑; PML ↑; B2M
   ↑; PTAFR ↑; SP100 ↑; TRIM21 ↑; STAT1 ↑; TRIM25 ↑; TRIM26 ↑; CAMK2D ↑;
   NLRC5 ↑; TRIM5 ↑; OASL ↑; NMI ↑; CD44 ↑
   GO:0008637 Apoptotic mitochondrial changes 81 22 11.64 0.015 CLU ↓;
   ERBB4 ↓; EYA2 ↓; IFI6 ↑; GCLC ↓; GCLM ↓; SFN ↑; SLC25A4 ↓; HK2 ↑; IFIT2
   ↑; MCL1 ↑; MMP9 ↑; PLAUR ↑; PMAIP1 ↑; BAK1 ↑; SOD2 ↑; BNIP3L ↓; PLA2G6
   ↓; CASP8 ↑; AIFM2 ↑; TNFSF10 ↑; NOL3 ↓
   GO:0048661 Positive regulation of smooth muscle cell proliferation 58
   19 8.34 0.036 NAMPT ↑; ABCC4 ↓; PPARGC1A ↓; HBEGF ↑; EREG ↑; FGFR2 ↓;
   ID2 ↓; IRAK1 ↑; JAK2 ↑; MMP9 ↑; NOTCH3 ↓; IRAK4 ↑; SERPINF2 ↓; PTAFR ↑;
   PTGS2 ↑; CCL5 ↑; STAT1 ↑; THBS1 ↑; CAMK2D ↑
   GO:0045824 Negative regulation of innate immune response 28 13 4.02
   0.015 ADAR ↑; TRAFD1 ↑; A2M ↓; HLA-E ↑; IFI16 ↑; LGALS9 ↑; SERPINB9 ↑;
   CEACAM1 ↑; SERPING1 ↑; TNFAIP3 ↑; DHX58 ↑; NLRC5 ↑; NMI ↑
   GO:0070542 Response to fatty acid 32 13 4.6 0.015 PPARGC1A ↓; CREB1 ↑;
   CTGF ↓; ALAD ↓; GIPR ↓; PDK4 ↓; PTAFR ↑; PTGS2 ↑; TLR2 ↑; UCP2 ↓;
   ZC3H12A ↑; CAT ↓; CD36 ↓
   GO:0060544 Regulation of necroptotic process 11 7 1.58 0.043 RBCK1 ↑;
   RIPK3 ↑; SLC25A4 ↓; BIRC2 ↑; BIRC3 ↑; CFLAR ↑; ARHGEF2 ↑
   [98]Open in a new tab

   For the REACTOME pathway analysis, we selected pathways that were
   enriched with more than 30 significant genes present in the enriched
   pathway, and these were generally in agreement with the GO analysis
   ([99]Table 3). In addition to IFNγ signaling (877300) and apoptosis
   (109581), it also revealed changes in epithelial-initiated B cell
   receptor signaling (983705 and 1168372) and amino acid metabolism
   (71291) following influenza infection. It is noteworthy that these
   pathways were initiated at the epithelial level without the
   participation of immune cells, thus highlighting the relevant genes of
   interest for future studies.

Table 3.

   48 hpi REACTOME pathways potentially associated with influenza
   pathogenesis.
   ID Description Gene Ratio Bg Ratio p-value p.adjust q-value geneID
   877300 Interferon gamma signaling 33/776 65/4637 2.081 × 10^−10 3.83 ×
   10^−8 3.32 × 10^−8 OAS2 ↑; OAS3 ↑;IRF7 ↑; OAS1 ↑; STAT1 ↑; TRIM21 ↑;
   TRIM5 ↑; PML ↑; HLA-F ↑; HLA-E ↑; SP100 ↑; OASL ↑; TRIM22 ↑; TRIM25 ↑;
   HLA-B ↑; GBP1 ↑; HLA-C ↑; B2M ↑; TRIM38 ↑; IRF6 ↑; HLA-A ↑; IRF9 ↑;
   JAK2 ↑; TRIM26 ↑; IRF1 ↑; GBP2 ↑; HLA-G ↑; IRF2 ↑; CAMK2D ↑; PTAFR ↑;
   CD44 ↑; MT2A ↑; ICAM1 ↑
   983705 Signaling by the B Cell Receptor (BCR) 38/776 119/4637 2.72153 ×
   10^−5 0.001440597 0.001249039 PSME2 ↑; BTC ↑; PSMB9 ↑; DAPP1 ↑; PSMB8
   ↑; RASGRP3 ↑; RICTOR ↑; EREG ↑; PSMA6 ↑; PSMA4 ↑; MOV10 ↑; HBEGF ↑;
   ITPR3 ↑; PSMF1 ↑; FGF2 ↑; CARD11 ↑; PSME1 ↑; CREB1 ↑; LYN ↑; PSMB7 ↑;
   FGFR3 ↓; MALT1 ↑; CDKN1A ↑; PIK3AP1 ↑; AGO1 ↓; FGF1 ↓; FGFR2 ↓; FGF5 ↑;
   ERBB4 ↓; ITPR1 ↓; UBC ↑; NCK1 ↑; ERBB2 ↓; SH3KBP1 ↑; PSMA3 ↑; PSMA5 ↑;
   TRIB3 ↑; GAB1 ↑
   109581 Apoptosis 30/776 88/4637 4.88082 × 10^−5 0.00227963 0.001976505
   PSME2 ↑; PSMB9 ↑; PMAIP1 ↑; TNFSF10 ↑; DSG3 ↑; PSMB8 ↑; BAK1 ↑; TICAM1
   ↑; H1F0 ↑; CASP7 ↑; TLR3 ↑; CFLAR ↑; FAS ↑; PSMA6 ↑; PSMA4 ↑; PSMF1 ↑;
   OCLN ↑; TRAF2 ↑; XIAP ↑; PSME1 ↑; CASP8 ↑; MAGED1 ↓; SFN ↑; PSMB7 ↑;
   TJP1 ↑; BIRC2 ↑; UBC ↑; PLEC ↑; PSMA3 ↑; PSMA5 ↑
   71291 Metabolism of amino acids and derivatives 40/776 134/4637 9.37168
   × 10^−5 0.004025804 0.003490489 IL4I1 ↑; IDO1 ↑; PSME2 ↑; SQRDL ↑;
   PSMB9 ↑; PSMB8 ↑; GAMT ↓; PSMA6 ↑; SMS ↑; PSMA4 ↑; KYNU ↑; SLC5A5 ↓;
   PSMF1 ↑; SHMT1 ↓; ALDH7A1 ↓; AZIN2 ↑; PSME1 ↑; AFMID ↓; BBOX1 ↓; GCLC
   ↓; IYD ↓; ALDH4A1 ↓; FAH ↓; PSMB7 ↑; HIBCH ↓; MCCC1 ↓; ALDH18A1 ↓; TST
   ↓; GCLM ↓; BCAT1 ↑; PSAT1 ↑; GPT ↓; ASNS ↑; PSMA3 ↑; ADI1 ↓; PSMA5 ↑;
   SLC25A10 ↓; ALDH6A1 ↓
   1168372 Downstream signaling events of B Cell Receptor (BCR) 31/776
   98/4637 0.000180538 0.004865597 0.004218614 PSME2 ↑; BTC ↑; PSMB9 ↑;
   PSMB8 ↑; RASGRP3 ↑; RICTOR ↑; EREG ↑; PSMA6 ↑; PSMA4 ↑; MOV10 ↑; HBEGF
   ↑; PSMF1 ↑; FGF2 ↑; CARD11 ↑; PSME1 ↑; CREB1 ↑; PSMB7 ↑; FGFR3 ↓; MALT1
   ↑; CDKN1A ↑; AGO1 ↓; FGF1 ↓; FGFR2 ↓; FGF5 ↑; ERBB4 ↓; UBC ↑; ERBB2 ↓;
   PSMA3 ↑; PSMA5 ↑; TRIB3 ↑; GAB1 ↑
   [100]Open in a new tab

3.5. Comparison with Established Influenza Gene Expression Signatures and
Transcriptomes Revealed Consistent Upregulation of Influenza Response Genes

   Given that RNAseq analysis facilitates more accurate expression changes
   following infection compared to hybridization technology such as
   microarray, we conducted further analysis on the levels of gene
   expression changes to enable more stringent and accurate transcriptomic
   analyses for future studies. By comparing these results to a previous
   study that identified 11 influenza-specific signatures, we verified
   that these 11 genes were all expressed in infected nasal epithelium
   later at 24 hpi, but not at 8 hpi. Furthermore, at both 24 and 48 hpi,
   all but one of the 11 gene signatures exhibited elevated expression of
   >2.5-fold change (>1.33 Log[2]FC) compared to uninfected control hNECs
   ([101]Table 4). When we applied the higher fold change cutoff, the
   number of significant genes decreased by approximately 50% ([102]Figure
   4A), which was also congruent with the linear association observed
   earlier. Therefore, future studies on early transcriptional alterations
   could consider adopting the 2.5-fold change in expression as a more
   stringent threshold, which may be more feasible, especially for large
   transcriptomic studies that yield large numbers of data points.

Table 4.

   Summary of alterations in expression of influenza-specific gene
   signatures in infected hNECs.
   Gene_Name logCPM FDR 8 h LogFC 8 h FDR 24 h LogFC 24 h FDR 48 h LogFC
   48 h
   ZBP1 5.467096 N.S. 2.71 × 10^−6 6.419275 2.70 × 10^−8 6.97086
   IFI6 7.740724 N.S. 2.77 × 10^−14 4.780717 1.80 × 10^−23 6.30199
   CD38 5.090161 N.S. 1.91 × 10^−8 5.232278 2.10 × 10^−10 5.699259
   MX1 10.7197 N.S. 2.36 × 10^−13 4.87848 3.48 × 10^−17 5.339009
   HERC5 6.743117 N.S. 1.05 × 10^−9 5.093646 6.89 × 10^−12 5.28321
   IFIH1 7.876294 N.S. 2.71 × 10^−4 3.566556 7.51 × 10^−6 3.972719
   HERC6 8.376077 N.S. 1.66 × 10^−11 3.339444 7.15 × 10^−6 3.843258
   RTP4 4.998082 N.S. 4.97 × 10^−11 2.846709 7.28 × 10^−15 3.25471
   PARP12 5.511677 N.S. 3.12 × 10^−2 2.203849 1.99 × 10^−2 2.11728
   LY6E 8.624757 N.S. 2.10 × 10^−6 1.451684 8.74 × 10^−12 1.989543
   LGALS3BP 11.0223 N.S. 1.30 × 10^−3 0.855231 1.57 × 10^−8 1.336304
   [103]Open in a new tab

Figure 4.

   [104]Figure 4
   [105]Open in a new tab

   Significant gene expression changes (>2.5-fold change) following
   influenza infection. (A) Significant gene expression changes over time
   (FDR < 0.05) at >2.5-fold change. Red bars indicate the number of
   upregulated genes, while blue bars indicate the number of downregulated
   genes. The number of significant gene expression changes peaked at 48
   hpi, with approximately half of the genes being retained following the
   application of the new threshold. (B) Venn diagram of significant genes
   found in RNAseq (blue) and microarray (yellow). A high degree of
   overlap was observed where 1/3 of RNAseq genes overlapped with half of
   the microarray genes.

   In addition, when compared to the previous microarray study on a
   similar hNEC model [[106]6], both RNAseq and microarray shared a high
   degree of overlap, with about one third and half of total genes from
   RNAseq and microarray overlapping, respectively ([107]Figure 4B). The
   overlap was generally observed in genes with highly altered expression,
   such as CXCL10, CXCL11, and RSAD2, which were changed to a similar
   magnitude in both RNAseq and microarray ([108]Table S7). When we
   compared the 11 influenza signature genes, RNAseq revealed a more
   consistent increase in magnitude, i.e., at 48 hpi, the magnitude of the
   gene change was generally higher than that of the microarray
   ([109]Table 5). In addition, RNAseq was also able to detect novel genes
   with expression changes of high magnitude that were generally higher
   than those found by microarray only (193 genes versus 22 genes with
   elevated expression greater than 2.5-fold). Genes such as HEATR9,
   PDCD1, IL4I1, ART3, and KCNH7 were altered to a higher magnitude than
   the 2.5-fold threshold. Hence, RNAseq-based transcriptomic analysis may
   augment transcriptomic findings to identify novel gene responses
   against influenza in the future.

Table 5.

   Expression of influenza signature genes in infected hNECs at 48 hpi
   detected by H3N2 microarray and RNAseq.
   Gene name Microarray LogFC 48 h Microarray FDR 48 h RNAseq LogFC 48 h
   RNAseq FDR 48 h
   ZBP1 1.30628 1.04 × 10^−6 6.97086 2.70 × 10^−8
   IFI6 3.78685 6.67 × 10^−6 6.30199 1.80 × 10^−23
   CD38 3.61175 5.88 × 10^−7 5.699259 2.10 × 10^−10
   MX1 3.17584 1.57 × 10^−8 5.339009 3.48 × 10^−17
   HERC5 3.53433 1.65 × 10^−8 5.28321 6.89 × 10^−12
   IFIH1 4.08252 5.01 × 10^−8 3.972719 7.51 × 10^−6
   HERC6 3.10932 9.16 × 10^−8 3.843258 7.15 × 10^−16
   RTP4 2.41668 6.25 × 10^−7 3.25471 7.28 × 10^−15
   PARP12 2.02756 4.75 × 10^−8 2.11728 1.99 × 10^−2
   LY6E 1.11967 6.02 × 10^−5 1.989543 8.74 × 10^−12
   LGALS3BP 0.863641 7.59 × 10^−4 1.336304 1.57 × 10^−8
   [110]Open in a new tab

3.6. Secretory Cytokine Profiles During Influenza Virus Infection of hNECs
Reveal Early Suppression of Cytokines in Apical Secretion

   After deriving the transcriptomes by RNAseq, we then further
   investigated whether the changes in expression of genes resulted in
   alterations in secretory cytokines and chemokines early in the
   infection of hNECs. Initially, we detected significant reductions in
   multiple cytokines at 8 hpi, with the exception of IL-15 which was
   increased ([111]Figure S2). This may reflect the initial immune
   suppression during influenza infection. However, at 24 and 48 hpi, less
   significant changes were observed, i.e., only increase in TNF-a and
   decrease in MDC and PDGF-AA were noted at 24 hpi. This was followed by
   increase in IP-10 (CXCL10) and TGF-a and decrease in PDGF-AA seen at 48
   hpi. This analysis highlights changes in IP-10, TGF-a, and PDGF-AA to
   be significant early responses in secretory cytokines/chemokines
   following influenza infection.

4. Discussion

   Our study has identified epithelium-initiated host responses which are
   found to be involved in both innate and adaptive responses. The finding
   is significant as we can now focus on the primary point of contact of
   influenza—the nasal epithelium in the study of early host responses for
   identifying host factors that can be utilized for diagnostic and
   therapeutic purposes [[112]7]. In addition, our study also showed that
   it is important for reference databases to use relevant human models
   like the hNECs model, which contains the mucociliary component of the
   airways, in order to provide closely representative host responses.
   While there exists a high number of microarray studies that showed the
   host responses using similar hNEC models, there are only a small number
   of equivalent RNAseq studies. Compared to microarrays, RNAseq analysis
   can provide a more comprehensive picture of the transcriptomic
   landscape, and is not limited by the hybrid library variant and
   concentrations [[113]14]. Hence, in order to derive accurate magnitude
   of gene expression changes, we performed an RNAseq analysis of H3N2
   infection using the hNECs model. H3N2 influenza virus was selected,
   given that it is a major circulating subtype over long periods of time.
   In addition, relatively lower efficacy of vaccines against this subtype
   prompted us to study its interactions with the primary host target to
   elucidate the immune responses and association with adaptive immunity
   [[114]17,[115]18]. This model has been previously evaluated to be a
   highly clinically-relevant model that can facilitate controlled
   infection of nasal cells from multiple individuals. In addition, we
   have also previously shown—by microarray analysis—that the nasal
   epithelium is responsible for the initiation of host responses
   following influenza infection [[116]6]. This renders the hNECs to serve
   as a valuable tool to analyze transcriptomics from different
   individuals infected under the same conditions to ensure consistent and
   relevant responses in humans.

   Once the magnitude of gene expression changes was considered, several
   interesting findings emerged. Firstly, the infected hNECs were observed
   with strong activation of antiviral genes and early inflammatory genes
   leading to type I immune responses. A large number of gene expression
   changes were of magnitude of over 100-fold difference (Log[2]7 to
   Log[2]9 fold change). Most of the genes with high-magnitude expression
   changes were verified by qPCR, with statistical significance congruent
   with the RNAseq analysis. Secondly, despite the absence of immune
   cells, the infected hNECs were able to generate strong type I responses
   that may likely aid the recruitment of cytotoxic cells to clear the
   infected cells. Thirdly, in early responses of the hNECs, IFNλ genes,
   which represent type III interferons, were more strongly induced than
   the more frequently observed type I interferons (IFNα and IFNβ), while
   type II interferons were not produced by hNECs, in agreement with
   previous studies [[117]5,[118]19]. The induction of type III
   interferons may reflect an important event within the hNECs where IFNλ,
   the initial responders against the infection, may be more beneficial in
   the antiviral response [[119]20,[120]21]. Moreover, we also observed
   notable suppression of expression of certain genes following influenza
   infection, including suppression of proliferation and DNA repair genes,
   which may contribute to the pathology and pathogenesis of influenza
   [[121]22]. Finally, RNAseq also unraveled expression changes of certain
   newly-discovered genes in response to influenza infection of the upper
   airway cells. Genes such as HEATR9 [[122]23], IL4I1 [[123]24], TNFSF13B
   (BAFF) [[124]25], and PDCD1 (PD-1) [[125]26] are recently implicated in
   influenza pathogenesis and mucosal defense, thereby signifying the role
   of the nasal epithelium against influenza infection. Furthermore,
   RNAseq identified altered expression of ART3 and KCNH7 genes that were
   not previously detected in influenza transcriptomes. These findings
   hence further reiterate the value of RNAseq in enhancing data on
   influenza transcriptomes for reference in future studies.

   Via pathway enrichment analysis, we have identified known antiviral
   pathways to validate the hNECs responses against influenza. In
   addition, we have also documented the potential pathways initiated by
   the nasal epithelium that may contribute to influenza pathogenesis as
   represented by the gene expression changes listed in [126]Table 2 and
   [127]Table 3. By analyses using literature-inferred GO and REACTOME
   databases, we have demonstrated that the nasal epithelium can play a
   role in the main antiviral signaling, i.e., IFNγ responses despite not
   being a direct producer of IFNγ. The pathway enrichment indicated that
   hNECs may serve as important regulators of type II interferons. Even
   though the effects of IFNγ are vital to the robust clearance of
   influenza viruses [[128]27], there are reports of unregulated IFNγ
   being a contributor to inflammatory damage [[129]28,[130]29].
   Therefore, the over-production of IFNγ response factors such as ICAM1
   and CD44 may contribute to inflammatory damage of the epithelium.
   Hence, production of factors such as STAT1 [[131]28] by the hNECs is
   also crucial in ensuring appropriate regulation of IFNγ-mediated
   expression of influenza response genes to modulate inflammation and to
   minimize damage. The primary contact of influenza virus with the nasal
   epithelium may subsequently lead to damage to the airway epithelium as
   well. This is apparent with the clear enrichment of the pathways of
   apoptosis, mitochondrial apoptotic processes, and necroptosis that
   contribute to cell death and mechanical barrier loss during infection
   [[132]5,[133]30]. Genes such as IFI6, BAK1, CAPS8, TNFSF10, and FAS
   suggest active apoptotic cell death that not only destroys cells in the
   epithelial barrier, but may also serve to propagate the virus and to
   perpetuate the damage [[134]31,[135]32,[136]33]. Furthermore, during
   virus infection, aberrant regulation of apoptosis may also lead to
   further injury to the epithelium and surrounding tissues [[137]34]. On
   the other hand, necroptosis pathways have also been observed to be
   enriched in influenza-infected hNECs. Compared to apoptosis, the study
   of necroptosis in influenza infection is relatively new with
   contradicting findings [[138]34]. RIPK3/necroptosis studies appear to
   generate contradictory results as to whether necroptosis protects
   against or is detrimental during influenza infection [[139]35,[140]36].
   Hence, its increased expression during infection of hNECs warrants
   further investigation on its role in influenza-induced damage. In
   addition, we also noted enrichment of B-cell signaling pathways by the
   infected hNECs which may be vital for B-cell responses during the
   adaptive immune response [[141]37]. We noted that most genes enriched
   in the B-cell pathways were related to antigen recognition such as
   proteasome subunits (PSME2, PSMB9, PSMA6, etc.) and B-cell
   receptor-associated genes such as DAPP1 and CARD11 [[142]38]. However,
   changes in expression of certain growth factors (including EREG and
   FGFs) following influenza infection may lead to complications involving
   airway remodeling and recruitment [[143]39,[144]40]. Further, the
   effects of the growth factors were further confirmed by the enrichment
   of pathways related to the proliferation of smooth muscle cells also
   induced by the infection. Changes to airway smooth muscle cells are
   usually implicated in airway remodeling [[145]41,[146]42,[147]43], and
   may also contribute to post-influenza complications. Hence, the genes
   found in this study may be crucial for elucidating the nasal-initiated
   responses that may contribute to the pathology and pathogenesis of
   influenza infection of the airways. Finally, another interesting
   pathway that may contribute to epithelial damage is the negative
   regulation of innate immune responses. These genes may serve as
   proviral factors and aid in immune evasion. For example, ADAR is a
   proviral factor that works in synergy with influenza NS1 to enhance
   viral replication [[148]44]. TRAFD1 is a negative regulator of
   toll-like receptor signaling which is upregulated in influenza-infected
   hNECs [[149]45]. DHX58 is a negative regulator of RIG-I/MDA5 signaling
   pathway [[150]46]. CEACAM1 is involved in regulation of liver
   inflammation [[151]47] and its expression appears to exert antiviral
   effects on influenza virus [[152]48]. NMI binds to influenza virus NS1
   and inhibits IRF7-mediated interferon signaling [[153]49,[154]50].
   Therefore, aberrant expression of genes in this signaling pathway may
   directly contribute to immune evasion of influenza, culminating in
   viral propagation and increased epithelial damage.

   We summarized the identified pathways (listed in [155]Table 2 and
   [156]Table 3) that alluded to immune evasion (negative regulation of
   innate immune responses), antigen processing (metabolism of amino acids
   and derivatives), and immunomodulation (interferon gamma signaling, B
   cell receptor signaling, and response to fatty acid) that may
   contribute to severity of influenza. There was evidence of direct
   pathway enrichment of potential influenza evasion strategies and/or
   immunomodulation with accompanying transcriptomic changes. The genes in
   the pathway may be analyzed for their immunomodulatory activity and
   whether their expression is beneficial to the virus (immune evasion) or
   the host (preventing cytokine storm). In addition, the infected hNECs
   also revealed modified responses associated with fatty acid, with many
   lipid signaling molecules such as leukotrienes that mediate antiviral
   responses and subsequent inflammation of the airway [[157]51,[158]52].
   Such modified responses may also determine the afforded in the airway
   and the severity of airway inflammation and damage. In addition, the
   modification may also affect the lower airway responses to inflammatory
   mediators; hence, the changes in these pathways may also suggest a
   potential mechanistic link to the pathogenesis of viral-induced
   exacerbation of chronic inflammatory diseases. Lastly, we also noted
   enrichment of pathways related to amino acid metabolism, which is
   important in antigen processing and proteasomal degradation of foreign
   protein. The changes in these genes at the hNECs, the target site of
   influenza infection, may determine the effectiveness of antiviral
   responses mounted and may therefore influenza disease severity.

   In addition, we also compared our RNAseq analysis against previously
   reported influenza-specific signatures in order to improve future
   transcriptomic analysis [[159]53]. In vitro transcriptomic analysis
   yields a large number of differentially-expressed genes that would
   require additional criteria to identify functionally significant genes.
   By means of this comparison, we discovered that almost all
   influenza-specific signatures exhibited differences in expression of
   above 2.5-fold. Hence, we propose applying fold change of >2.5 as a
   threshold for future in vitro transcriptomic systems analyses, in order
   to increase the stringency in detecting functionally significant gene
   changes.

   Finally, we observed that, unlike the transcriptome, there were notably
   fewer cytokines that were readily secreted during the acute phase of
   infection. Expression of cytokines was reduced at 8 hpi, except for
   IL-15, which interestingly is implicated in influenza-induced acute
   lung injury [[160]54]. Even fewer cytokines showed altered expression
   at later time points. Among them, only TGF-a, IP-10 (CXCL10), and
   PDGF-AA were significantly altered at 24 and 48 hpi. These may be
   significant markers that can be detected in the secretion of
   influenza-infected mucosal surface that may influence the severity of
   influenza. IP-10 is a well-established IFNγ response gene, and serves
   as a useful marker for response against influenza [[161]5,[162]6].
   TGF-a represents an important factor involved in the secretion of IL-8
   in response to influenza, and may determine the early appropriate
   innate responses to prevent severe disease [[163]55]. On the other
   hand, it is also involved in pulmonary fibrosis as a ligand of
   epidermal growth factor receptor (EGFR) and may contribute to
   complications in the lower airway [[164]56]. PDGF-AA was found to be
   elevated in the cerebrospinal fluid of influenza-associated
   encephalopathy [[165]57], but was consistently reduced in hNEC
   secretory fluid, thus warranting further investigation into its role in
   the infected nasal mucosa.

   The establishment of a reference transcriptome based on early responses
   of the human nasal epithelium model serves a key role in research on
   critical host factors involved in influenza. As the primary host
   contact with the virus, not only are immune responses against influenza
   important, but also the alterations in non-immune functions such as
   metabolism, cell content, and cell cycle, which may contribute to
   disease severity. In terms of translational potential, the model system
   identified gene expression changes of significant magnitude and
   pathways that impact responses against influenza and its severity.
   These genes may represent novel targets for future diagnostic and
   therapeutic development. Under controlled conditions, the hNECs
   clinically establish the baseline for “normal” innate immune responses
   of the nasal epithelium against influenza viral infection. Such a
   baseline can be particularly crucial when studying the changes in
   innate immune responses against influenza, especially in patients with
   underlying chronic diseases who may have aberrant airway responses
   against influenza. Their antiviral responses may differ from “normal”
   subjects, and this study thus provides the basis for comparing the
   differential responses that culminate in more severe infections in
   patients with co-morbidities such as diabetes and chronic airway
   inflammatory diseases. Such comparative clinical studies can
   potentially enhance the management of influenza viral infection in
   patients with chronic diseases.

5. Conclusions

   In conclusion, RNAseq technology allowed us to accurately quantify the
   magnitude of gene expression changes, as well as the relevant enriched
   pathways during H3N2 influenza virus infection of hNECs, which can
   serve as a baseline for future clinical studies. The establishment of
   this baseline under controlled condition elucidated the antiviral
   innate response by the infected nasal epithelium, and highlighted the
   molecular factors and abnormalities in the upper airway that may
   contribute to influenza severity. Furthermore, this study also
   culminated in the identification of novel gene signatures and host
   factors that may be harnessed for future research to develop influenza
   diagnostic markers and therapeutic targets.

Acknowledgments