key: cord-0814845-k5569i0l authors: Murray, Liisa M; Thillaiyampalam, Gayathri; Xi, Yang; Cristino, Alexandre S; Upham, John W title: Whole transcriptome analysis of high and low IFN‐α producers reveals differential response patterns following rhinovirus stimulation date: 2021-11-17 journal: Clin Transl Immunology DOI: 10.1002/cti2.1356 sha: 9e7e21bad9f2b17a6286a9b06f94d7d9094431ce doc_id: 814845 cord_uid: k5569i0l OBJECTIVES: Viral respiratory infections cause considerable morbidity and economic loss. While rhinoviruses (RV) typically cause little more than the common cold, they can produce severe infections and disease exacerbations in susceptible individuals, such as those with asthma. Variations in the regulation of key antiviral cytokines, particularly type I interferon (IFN‐α and IFN‐β), may contribute to RV susceptibility. To understand this variability, we compared the transcriptomes of high and low type I IFN producers. METHODS: Blood mononuclear cells from 238 individuals with or without asthma were cultured in the presence or absence of RV. Those samples demonstrating high or low RV‐stimulated IFN‐α production (N = 75) underwent RNA‐sequencing. RESULTS: Gene expression patterns were similar in samples from healthy participants and those with asthma. At baseline, the high IFN‐α producer group showed higher expression of genes associated with plasmacytoid dendritic cells, the innate immune response and vitamin D activation, but lower expression of oxidative stress pathways than the low IFN‐α producer group. After RV stimulation, the high IFN‐α producer group showed higher expression of genes found in immune response biological pathways and lower expression of genes linked to developmental and catabolic processes when compared to the low IFN‐α producer group. CONCLUSIONS: These differences suggest that the high IFN‐α group has a higher level of immune system readiness, resulting in a more intense and perhaps more focussed pathogen‐specific immune response. These results contribute to a better understanding of the variability in type I IFN production between individuals. Rhinoviruses (RV) remain the most common respiratory virus causing common colds, even during the SARS-CoV-2 pandemic. 1 While RV infection is seemingly innocuous in most healthy people, it can cause complications in vulnerable groups with chronic respiratory diseases such as cystic fibrosis, idiopathic pulmonary fibrosis, chronic obstructive pulmonary disease and asthma. 2 Indeed, RV causes more asthma exacerbations than the influenza virus. 3 RV infection worsens airway inflammation and may predispose to secondary bacterial infections, though whether these are due to abnormal host immune defence, to RV-induced immunopathology or to both remains unclear. 4 Rhinoviruses-stimulated airway structural cells and circulating leukocytes from people with asthma may produce insufficient amounts of type I interferons (IFN-a and IFN-b), as reviewed recently, 5 although not all investigators have been able to confirm these findings. 6, 7 IFN production appears to vary with asthma severity 8 and between different inflammatory phenotypes, 9 suggesting that IFN insufficiency may be characteristic of specific asthma subtypes, rather than being common to all with asthma. Recent study has linked both high and low IFN production to acute wheezing in children with respiratory illness, further supporting variable IFN responses to viruses in asthma patients. 10 Antiviral immune responses depend on pathogen pattern receptors, such as melanoma differentiation-associated gene 5 (MDA-5), retinoic acid-inducible protein I (RIG-1) and various Toll-like receptors (TLRs) that recognise specific viral molecular patterns and induce type I IFN production. While all nucleated cells produce IFN-a and IFN-b, plasmacytoid dendritic cells (pDC) have the greatest capacity to produce these cytokines and are central to antiviral immune responses. 11 Optimal IFN-a and IFN-b production is critical for the antiviral immune response, as too little or too much IFN has adverse health consequences. 12 Hence, type I IFN production needs to be tightly regulated, involving multiple checkpoints and complex patterns of gene expressions. 13 At a molecular level, the interferon regulatory factor (IRF)7 is a major hub regulating type I IFN production, 14 while at a cellular level, pDC are recognised for their capacity to rapidly produce large amounts of IFN-a, with recent publications highlighting key mechanisms. 15, 16 Although such studies have shed important light on type I IFN biology, it is also important to understand how IFN-a production varies between individuals. Investigating differences between people exhibiting efficient and weak antiviral responses may provide important insights into infection susceptibility. Differences occurring before and during infection may contribute to efficient/weak antiviral responses. While a variety of researchers have made important observations concerning antiviral immune response variability, [17] [18] [19] [20] [21] [22] [23] these studies have not examined RV, which is known to induce distinct host responses relative to other viruses. 24 To address this important knowledge gap, the primary study aims were to define the folllowing: 1. In the absence of RV stimulation which transcriptional variations differentiate high IFN-a and low IFN-a producers, and 2. During RV stimulation which transcriptional variations differentiate high IFN-a and low IFN-a producers. The secondary study aim was to determine whether transcriptional variations differed between those with asthma and healthy participants in unstimulated and RV-stimulated conditions. We made use of samples obtained during our recent study of the immunological and clinical variables associated with cold frequency in 301 individuals. 25 RV-stimulated peripheral blood mononuclear cells (PBMC) showed a very broad range of IFN-a production in these individuals. From these samples, we selected a subsample of 75 individuals with contrasting high IFN-a and low IFN-a production in order to perform differential gene expression analysis in these 'extreme phenotypes'. Because the study groups differed in mean body mass index (BMI), statistical methods were used to assess whether BMI might be affecting gene expression. 1690.6). In order to provide more detailed analysis of gene expression patterns, we focussed on the extreme phenotypes of IFN-a production: 75 participants with the highest IFN-a and lowest IFN-a production were selected from asthma and healthy groups to contrast differences in wholegenome transcriptomes ( Figure 2 ). The analysis was restricted to participants with European ancestry, in order to minimise confounding by ancestral background. Four samples did not pass quality control during RNA sequencing or principal component analysis (PCA) including one control asthma IFN-a-low sample, one control healthy IFN-a-high sample, one RV healthy IFN-ahigh sample and one RV healthy IFN-a-low sample. These four samples were excluded from RNA sequencing, thereby reducing the total number of samples used in the complete analysis from 150 to 146 (Supplementary table 1) . Further details regarding study design are shown in Figure 2 . The only criteria used to select samples for transcriptome analysis were high or low IFN-a production ( Figure 2 ) and the presence or absence of asthma. It is noteworthy that the IFN-a-low groups had higher BMI (24.0 and 27.0 kg m À2 vs 23.0 and 24.0 kg m À2 ; P = 0.008) than the IFN-ahigh group (see Supplementary table 1) . There were age differences between the four groups (P = 0.042): the asthma IFN-a-low group was oldest (median age 47 years), whereas the asthma IFN-a-high group was youngest (median age 27 years) with the two healthy groups in between (healthy IFN-a-high median 32 years; healthy IFNa-low median 34 years). Minor differences in sex ratios were observed between groups, but these were not statistically significant. There were no significant differences in asthma control, asthma severity or self-reported respiratory infection frequency between IFN-a-high and IFN-a-low groups (Supplementary table 2) . Gene expression signatures of IFN-a producer groups are not determined by BMI As the IFN-a-low groups had significantly higher BMI than the high IFN-a producer groups, this raised a concern that BMI differences could confound the analysis. We used statistical two-way ANOVA Relative frequency 414 2200 Figure 1 . Variation in RV-stimulated IFN-a production. Frequency distribution of RV-stimulated IFN-a production in 238 PBMC samples tested in two biological replicates. RV-stimulated samples with low IFN-a production (lowest 15% of IFN-a production) produced less than 414 pg mL À1 IFN-a (orange-shaded), and samples with high IFN-a production (highest 15% of IFN-a production) produced more than 2200 pg mL À1 IFN-a (green-shaded testing to identify which differentially expressed genes (DEGs) between high and low IFN-a producer groups were likely confounded by BMI. The participants were divided into two groups with BMI below or at/above 25 kg m À2 , the threshold for normal weight and overweight condition. The resulting interaction between IFN-a producer groups and BMI groups identified the ability of BMI to confound the signatures of IFN-a producer groups. We observed that only a small fraction of genes show a significant interaction between BMI and IFN-a producer groups (baseline -2%, RV treatment -5%) and are likely to behave differently in obese and normal weight samples; however, majority of genes are used reliably to differentiate gene expression between high and low IFN producer groups (Supplementary figure 1) . Further details on the DEGs that associate with BMI regardless of IFN-a producing ability are listed in Supplementary figure 2. High IFN-a producers have more circulating plasmacytoid dendritic cells and produce more TNF As pDC are the most potent IFN-a-producing cells in the circulation, variations in blood pDC numbers are likely to impact on IFN-a production in response to RV stimulation. Hence, we compared pDC quantity and cytokine production by RV-stimulated PBMC across the 75 control samples stratified by IFN-a production and presence or absence of asthma. We have previously shown that the gene expression of Ctype lectin domain family 4 member C (CLEC4C) can be used to quantify circulating pDC in whole blood. 26 As expected, CLEC4C gene expression in whole blood was significantly lower in the IFN-alow groups signifying lower quantities of the IFNa-producing pDC ( Figure 3 ; P-value < 0.001). We also assessed production of other cytokines by RV-stimulated PBMC. The production of the TLR8-stimulated pro-inflammatory cytokine TNF was also lower in the IFN-a-low groups (Figure 3 ; P-value = 0.006), whereas the production of another TLR8-activated pro-inflammatory cytokine IL-12 was similar in all four groups. The whole-genome transcriptomes were compared across the 75 participants stratified by IFN-a production and the presence or absence of asthma. PCA transcriptomes indicated a marked difference between unstimulated and RVstimulated samples (not shown), and a modest grouping of samples by IFN-a producer groups across PC1 [ Figure 4 unstimulated (a), RVstimulated (b)]. Because PCA revealed no apparent distinction between asthmatic and non-asthmatic samples, it seemed reasonable to undertake a pooled analysis of the IFN-a-high and IFN-a-low groups, regardless of asthma status, to address factors associated with variations in IFN-a productionthe primary objective of the study. We next performed a separate analysis of the transcriptomes using unstimulated and RVstimulated PBMC from the IFN-a-high and IFN-alow producer groups. With a cut-off of false discovery rate (FDR) q-value < 0.05 and log foldchange > 1 of gene expression, unpaired differential gene expression analysis in the unstimulated samples returned 39 genes with high expression (Supplementary table 3 ) and 10 genes with low expression (Supplementary table 4) in the IFN-a-high group compared to the IFN-a-low group. The heatmap visualisation in Figure 5 reveals a major subcluster of genes with high expression that have consistently high expression in the IFN-a-high group and low expression in the IFN-a-low group. Sample clustering with Ward hierarchical cluster analysis positions asthma samples randomly, which supports our results from the PCA ( Figure 4 ) that asthma status does not have a major impact on differential gene expression in our samples. We searched the manually annotated records in the protein knowledgebase UniProt 27 and PubMed (NIH, USA) for any additional literature supporting the understanding of the functions of the DEG products. Common functional groups that arose from the literature are annotated for each gene in Supplementary tables 3 and 4. In comparison with the IFN-a-low group, the genes with high expression in the IFN-a-high group were most often related to immune function, whereas the genes with low expression were most often related to oxygen transport. The immune-related genes included those associated with the complement system (C1QA, C1QB, C1QC, C3, CD93), antigen presentation (CLEC4C, CPVL, RNASE1, ASGR2), B-cell chemokine (CXCL13), immune regulation (IDO1, CHI3L1, UBD, VCAN), macrophage and pDC surface markers (CLEC4C, MS4A4A, MS4A6A). Three potassium channel encoding genes (KCNJ10, KCNJ15 and KCNMA1) and a transcription factor gene ETV3L associated with vitamin D activation also had high expression in the IFN-a-high samples. Several genes with high expression in IFN-a-high have been associated with allergy [CD93 28 ; Figure 3 . CLEC4C gene expression in whole blood (a) and cytokine production in PBMC by study group (b, c). Individual datapoints are overlaid with median (bold line) and interquartile range markers (thin lines) for each group. Statistically significant Kruskal-Wallis ranks sum tests between the four groups and significant post hoc Mann-Whitney U-tests are indicated as ***P < 0.001; **P < 0.01; *P < 0.05. A group of oxygen transport genes that have low expression encodes haemoglobin subunits (HBA1, HBA2, HBB, HBG2) and an enzyme involved in haem biosynthesis (ALAS2). The gene group forms a subcluster with low expression in most IFN-a-high samples. Four other genes with low expression have a role in wound healing (MET, DPYSL3), stress response (HSPA1B) and as inflammatory markers (HSPA1B and NR4A2). We next assessed differential gene expression in RV-stimulated PBMC. With a cut-off of FDR qvalue < 0.05 and log fold-change > 1, unpaired differential gene expression analysis returned 55 genes with high expression (Supplementary table 5) , and 73 genes with low expression (Supplementary table 6) in the IFN-a-high group. Like the clustering analysis of the unstimulated samples ( Figure 5 ), the RV-stimulated PBMC samples show expression profiles that cluster primarily by IFN-a-producer group, rather than by the presence or absence of asthma as presented in Figure 6 . This clustering analysis reveals three distinct clusters of genes with similar expression patterns ( Figure 6 ). The gene cluster in the middle contains the upregulated genes, including IFN genes. The expression of those genes is uniformly high in the IFN-a-high group and low in the IFN-a-low group. In contrast, the genes with low expression separate into two main clusters at the top and bottom of the heatmap. The expression of those genes is consistently low in the IFN-a-high group samples but is more varied in the IFN-a-low group samples. The top gene cluster contains several subclusters of genes with high expression in distinct IFN-a-low sample clusters. In contrast, the entire bottom cluster of genes has high expression in only one IFN-a-low sample cluster consisting of mostly non-asthma samples. Some clustering with asthma and healthy samples is evident, and these sample subclusters have a specific gene subcluster expression pattern. We searched the manually annotated records in the protein knowledge database UniProt and any additional literature available for the function of the DEG products. The most common relevant functions are annotated for each gene in Supplementary tables 5 and 6. As expected, IFNB1, IFNL1 and IFNW1 had high expression in the IFN-a-high samples (Supplementary table 4 Conversely, under RV stimulation, cells from the IFN-a-high group showed low expression of gene products related to structural or extracellular matrix (ECM) function. Three matrix metalloproteinases (MMP9, MMP7, MMP8) involved in ECM breakdown were prominent, while the matrix metalloproteinase inactivator TIMP3 showed only low expression. SERPINE1 and COL23A1 are also associated with ECM organisation. The gene products of COL23A1, GREM1, STAB1, CEACAM8, CEACAM6, DSC1, FLRT2, ITGA11 and TIMP3 have a role in cell-cell or cell-ECM adhesion. Similarly, SEMA4C, PDGFC, MMP9, MMP7, MMP8 and CEACAM6 upregulate cell migration, whereas PODN, CYP1B1, SERPINE1 and DPYSL3 downregulate cell migration. The gene products of PDGFC, CEACAM6 and GREM1 upregulate cell proliferation, whereas PODN and CD9 downregulate cell proliferation. Also, GPC4, MET and SDC2 may be involved in cell proliferation. Three genes had high expression in both unstimulated and RV-stimulated samples, while seven genes showed only low expression. UBD, which was highly expressed, encodes an ubiquitinlike protein that tags proteins for degradation and regulates TNF-mediated NFKB signalling and dendritic cell (DC) activation/maturation. CXCL13 and PPFIA4 genes were highly expressed; they encode a B-cell chemokine and focal adhesion regulator, respectively. IFN-a-high samples showed low expression of the haemoglobin subunits (HBA1, HBA2, HBB, HBG2) and haem biosynthesis enzyme (ALAS2), cytoskeleton remodelling (DPYSL3) genes and MET, a gene encoding a transmembrane receptor involved in cellular proliferation and fibrosis. Biological pathway analysis comparing differentially expressed genes in the IFN-a high and low producers Gene set enrichment analysis (GSEA) was used to search for functional enrichment within the DEG sets. Figure 7a presents the gene ontology (GO) pathways associated with DEGs from unstimulated samples for the IFN-a-high group compared to the IFN-a-low group. The pathways with high expression predominantly reflect upregulation of multiple immune functions, whereas the pathways with low expression are much more varied in their function. A network representation visualises shared genes in selected pathways (Figure 7c ), demonstrating considerable gene overlap within both up-and downregulated pathways. Several downregulated pathways were enriched in the list of DEGs. Prominent among these important pathways were oxidative stress response genes, with high statistical significance, driven by a relatively small number of genes, including haemoglobin subunits (HBA1, HBA2, HBB) and inflammatory markers (HSPA1B and NR4A2). Together, the downregulated gene functions and the enriched GO pathways indicate a differential state of cellular stress between the IFN-a producer groups, such that the cells from high IFN-a producing individuals showed efficient downregulation of oxidative stress gene expression, whereas a subset of samples from low IFN-a producing individuals exhibited high expression of oxidative stress genes. The oxidative stress pathway was not enriched in the RVstimulated samples. Figure 6 . RV-stimulated differential gene expression. Gene expression of the DEGs in 73 RV-stimulated PBMC samples IFN-a-high group vs IFN-alow group. Log 2 -transformed gene expression is presented as a difference from the median red denoting high expression and blue low expression. Samples are hierarchically clustered with the Ward method and genes with the weighted pair group method with arithmetic mean. Samples are colour-coded for asthma status and IFN-a producer group. As expected, the significantly upregulated GO pathways during RV stimulation relate to activation of antiviral immune responses (Figure 7b) , whereas the genes with low expression are involved in cell proliferation, cell-cell communication and cell migration, which is reflected in the enriched pathways. A selection of a representative up-and downregulated pathways is presented in a network shown as Figure 7d . The RV-stimulated pathways show high interconnectivity and several shared genes between the upregulated pathways, whereas the downregulated pathways share less genes and are less connected. Validation of the DEGs using microarray results from a separate cohort To validate the results from the differential gene expression analysis, we accessed a similar cohort of RV-stimulated PBMC samples from a separate cohort of 17 asthma samples and 17 controls collected in 2014, as previously reported. 37 The comparison between the validation dataset and the current dataset is shown in Supplementary tables 1 and 2. We used the microarray results from the validation cohort to confirm the association between the expression of the DEGs and IFN-a production. Linear regression was used to test associations between gene expression and IFN-a production as IFN-a production in the validation dataset was lower than in the current study, and there were insufficient numbers to restrict the analysis to high and low IFN-a producers. Because of lower numbers of genes analysed by the microarray experiment, only 11 out of the 45 unstimulated DEGs from the main study were detected in the validation microarray (Figure 8a ), and 36 out of the 129 RV-stimulated DEGs were detected in the validation microarray (Figure 8b ). After correcting for multiple testing, the unstimulated DEG CLEC4C gene expression was validated as significantly associated with IFN-a production in the validation cohort and RVstimulated DEGs IFNW1 and GCOM1 gene expression was validated as significantly associated with IFN-a production (Supplementary table 7) . The primary aim of this study was to examine transcriptional variations associated with RVinduced IFN-a production. Having shown the broad range of IFN-a production in a large group of study participants, we focussed on the extreme phenotypes, contrasting whole-genome transcriptomes in the 15% with the highest or the 15% with the lowest IFN-a production. While others have comprehensively profiled variability in antimicrobial immune responses, 17 a strength of the current study was the large cohort in which weak and strong IFN-a producers were identified prior to transcriptomic analysis. Previous studies have not comprehensively examined variations in the RV-induced immune response, which is particularly important, given the frequency with which RV induces exacerbations of asthma and other respiratory illnesses and its predominance relative to other common respiratory viruses. Contrary to our expectation, we saw only modest differences in gene transcription patterns between those with asthma and healthy participants. A subgroup of people with both asthma and low IFN-a production showed high expression of two gene clusters ( Figure 6 ). This supports the notion that only a subgroup of people with asthma have specific alterations in IFN-a production and antiviral immunity, rather than this being a general characteristic of all people with asthma. 8, 9 The unstimulated samples were used as a baseline measurement with the potential to reveal biological factors predicting the strength of IFN-a response. The high and low IFN-a-producing groups demonstrated multiple DEGs, enriched in immune-related genes and GO pathways. The DEGs expressed by specific immune cells may point towards differences in the function of those cells or merely quantities. Conventional DC are essential for detecting invading pathogens and activating and directing the adaptive immune response, while pDC have a specialised role as potent IFN-a producers. 11, 37 The pDC-associated gene CLEC4C regulates antigen presentation and type I IFN production, while its consistent expression profile enables its use for pDC quantification. 26 CLEC4C had higher expression in unstimulated samples in the IFN-a-high group than in the IFN-a-low group in whole blood samples and in the PBMC samples used for whole transcriptome analysis (Figures 3 and 5) . This was confirmed in the validation study in which CLEC4C was the only DEG in unstimulated samples significantly associated with IFN-a production ( Figure 8a ). Thus, CLEC4C expression appears as a key variable predicting the magnitude of subsequent RV-induced IFN-a production. Other antigen-presenting cell (APC)-associated genes, including CD93, CPVL, IDO1, MS4A4A and VCAN, also had higher expression in the high IFN-a producer group in unstimulated samples supporting the close relationship between APC and an efficient IFN-a response. While MS4A4A is found on the surface of alternatively activated M2 macrophages, 38 MS4A4A also forms part of a type I IFN signature in early rheumatoid arthritis. 39 Thus, our finding of high MS4A4A expression in the high IFN-a producer group may be a consequence of systemic inflammatory signals, rather than fitting into a simple M1/M2 paradigm. The endosome protein encoded by ASGR2 is important for APC function for sampling the extracellular environment by mediating the endocytosis and lysosomal degradation of glycoproteins. 40 Moreover, the complement system is increasingly recognised to play an important role in antiviral immunity, 41 and many of the genes with high expression in the IFN-a-high group were members of the complement system whose genes are enriched in the innate immune response GOpathway. Notably, others have reported that this GO pathway contains multiple genes whose expression shows a high degree of inter-individual variability. 17 The RV-stimulated samples were used to examine RV infection in vitro. A triggered antiviral response results in a substantial increase in the expression of genes involved in defence against the pathogen. The larger amount of type I IFN produced in response to RV and the IFN genes with higher expression in the IFN-a-high group compared to IFN-a-low group reflect that. The higher expression of multiple Th1 and antiviral immune response genes in the IFN-a-high group than the IFN-a-low group further indicates a more extensive antiviral response. Notably, the IFN-ahigh group also showed higher TLR8-induced TNF production. Given that TLR8-induced cytokine production is dominated by TNF and IL-12 production, rather than IFN-a production, this suggests the IFN-a-high group has a greater capacity to respond to viruses and viral nucleic acids that extends beyond IFN-a production. The difference in developmental and metabolic process pathways between the two groups also supports the notion the IFN-a-high group devotes their transcriptome to a robust antiviral immune response and suppresses non-critical pathways, whereas the IFN-a-low group transcriptome retains the transcription of those other functions. Of the 56 genes with higher expression in the IFN-a-high group, ten IFN and three non-IFN genes correspond to significant cytokine expression quantitative trait loci (eQTL) genes. 20 One of these genes with high expression was IL18, a TLR8-induced cytokine that activates type 1 cytokine production in natural killer cells. 42 The association of IL18 with immune variation has previously been documented: an IL18 eQTL modulates influenza-induced IFNb production by DC 20 ; IL-18 binding protein (IL-18BP) restricts IL-18 availability and is inversely associated with cytokine production. 17 During experimental RV infection, IL-18 appears protective, with IL-18 concentrations in nasal and bronchial lining fluid inversely proportional to infection severity, both in healthy people and those with asthma. 43 Our own findings described herein support the conclusion that IL-18 production is an important component of effective host defence against severe rhinovirus infections and offers an attractive target for future research. The high IFN-a producer group showed lower RV-associated expression of multiple antibacterial genes such as the antibacterial peptides bactericidal permeability increasing protein (BPI), hepcidin antimicrobial peptide (HAMP) and lactotransferrin (LTF). The IFN-a-low group upregulated genes for lactotransferrin and several haem biosynthesis and haemoglobin components, which are important in binding iron. Bacteria require iron for survival, and binding iron to those molecules is a powerful antibacterial mechanism. 44 We speculate that this might reflect an immune response that is unspecific to pathogen type, whereas high IFN-a producers respond with a more highly focussed antiviral response. There are clearly complex interactions between viruses and bacteria in the lung. 45 Rhinoviruses are known to induce degradation of antimicrobial peptides and impair macrophage antibacterial responses. 46, 47 Outside the immediate immunological functions, cells from high IFN-a producers also showed efficient downregulation of the oxidative stress pathway at baseline, whereas a subset of low IFN-a producers exhibited upregulation of the oxidative stress pathway. Oxidative stress typically results from external factors such as cigarette smoke or intrinsic factors including reactive oxygen species (ROS) produced in the context of inflammation. Exogenous oxidants can reduce IFNa production by pDC, while mitochondrial reactive oxygen species can inhibit TLR7 function 48, 49 ; hence, oxidative stress might be a factor constraining virus-induced IFN production. Oxidative stress is known to be present in asthma and correlates with clinical severity, 50,51 which might contribute to the weak RV-induced IFN-a production reported in some studies. Addressing oxidative stress in asthma might be beneficial by promoting a stronger antiviral immune response. The groups showed differential expression of genes related to nutrient availability and function. The high IFN-a producer group showed higher baseline expression of ETV3L, a vitamin Dassociated transcription factor 52 and higher RVstimulated expression of CYP27B1, and the protein encoded by this gene converts 25-OH vitamin D3 to its active form. The expression of these genes indicates that the vitamin D availability may be necessary for a strong antiviral response. In contrast, TCN1, a vitamin B12 binding protein, was downregulated in IFN-a-high RVstimulated samples. By downregulating the B12 binding protein, high IFN-a producers may be making more B12 available. 53 Interestingly, Khoo and colleagues identified molecular phenotypes with some similarities to those identified herein in their study of upper airway specimens collected from children with acute asthma exacerbations. 10 In their cluster analysis, they identified two distinct molecular phenotypes, one characterised by high IRF7 and high IFN expression, and a second characterised by low IRF7 expression, growth factor signalling and downregulation of IFN. 10 Clinical features differed between these two phenotypes. There are a number of differences in study design, and Khoo and colleagues measured gene expression in upper airway cells collected from acutely unwell children without in vitro stimulation, whereas our study assessed circulating immune cells collected when participants were relatively well; gene expression was assessed in both stimulated and unstimulated cells. Nonetheless, we think that the two studies provide important complementary findings. Association tests between gene expression and IFN-a production in the validation dataset confirmed CLEC4C as an important baseline gene and confirmed that RV-stimulated IFNW1 and GCOM1 are significantly associated with IFN-a production. GCOM1 codes for GRINL1A complex locus 1 that interacts with the N-methyl Daspartate (NMDA) receptor in the nervous system. 54 There is some evidence that NMDA-type glutamate receptors are expressed on lymphocytes and neutrophils and that its activation has functional consequences. 55 In the central nervous system, interactions have been described between type I IFN function and NMDA receptors, 56 though whether this is relevant in relation to RV requires further investigation. The GCOM1 gene has been reported to be involved with transcription elongation, and because of its nature as a complex transcription unit, several transcriptional variants are produced that could potentially have a wide range of functions. 57 The study groups were not matched for BMI, age and gender, which may have influenced the findings in this study. However, our analysis indicated that only a handful of genes showed a significant interaction between IFN-a producer group and BMI. Piasecka et al. 58 showed that age influences influenza-induced IFN-a production, but we found little evidence in the current study that age was affecting RV-induced IFN-a production, although we acknowledge that the participants did not cover a wide age range. Sex hormones influence antiviral immunity 59,60 and we recently observed that age is associated with respiratory infection frequency in women but not men, implying a role for sex hormones in antiviral immunity. 25 Since the ratio of men and women in each of the groups was similar in that study and the current study, we think it is unlikely that gender imbalance had a major influence on our key findings. For the differential gene expression analysis, we chose to combine the asthma and healthy groups because IFN-a production was similar in both groups and asthma was not a significant factor in our principal component analysis. Despite this, our gene expression findings are in line with previous studies of healthy individuals. 17, 20 Finally, to validate our results, we accessed a similar cohort of asthma and control samples gene expression data performed with microarray technology. However, gene matches to DEGs were limited because of the restricted sensitivity of the microarray technology compared to RNAseq. Because of the study design, it would have been impractical to collect a second cohort with high and low producers, and as the validation samples are closer to the IFN-a-low group in their IFN-a production than IFN-a-high, the number of false negatives may have increased in the validation results. Collectively, these results show that at baseline, high IFN-a responders express essential components of the innate immunity (especially pDC-related CLEC4C) at a higher level than the low IFN-a responders. During the antiviral immune response, high IFN-a producers devote more of their transcriptome to the production of antiviral cytokines and effector proteins while downregulating expression of genes encoding antibacterial proteins and cellular processes, such as development and metabolism pathways. Differences between the high and low IFN-a antiviral response can be partly attributed to host factors regulating oxidative stress and availability of vitamin D and B12. This suggests that the high and low IFN-a groups have varying levels of immune system readiness and capacity to translate this into a pathogen-specific immune response. Our findings described herein contribute to a better understanding of the inter-individual variability in type I IFN production in the context of respiratory virus infections, vaccination, asthma and autoimmune disease. An asthma case-control study was designed to recruit 300 participants, but actually recruited 301. Details of the study cohort and participant characteristics were described recently. 25 The experiments described herein relate to PBMC samples from the subset of 238 participants with European ancestry (Figure 1 ), selected in order to minimise genetic variations related to ancestry. RV-induced IFN-a production was measured and two contrasting IFN-a producing groups were selected for more detailed analysis of gene expression, comprising the samples from the highest 15% or lowest 15% of IFN-a production in the European subset. Samples in the high IFN-a producer group were obtained from 18 participants with asthma and 19 healthy participants, while samples in the low IFN-a producer group were obtained from 19 participants with asthma and 19 healthy participants, comprising a total of 75 participants. The high and low IFN-a producer groups will henceforth be referred to as IFN-a-high and IFN-a-low groups. To validate the results, we used samples from a separate cohort recruited in 2015. 37 We had access to samples from 18 healthy participants and 17 participants with asthma from this latter study. One healthy sample that was part of both cohorts was excluded for the validation test reducing the number of healthy samples to 17 Peripheral blood mononuclear cells were isolated and freshly stimulated with RV16 or the TLR8 agonist VTX-2337 (1 lM; Sapphire Bioscience, Waterloo Australia) for 24 h as described in detail elsewhere. 25, 37 Control cells were cultured in media alone with no added stimuli. IFN-a, TNF and IL-12 concentrations in culture supernatants were measured with enzyme-linked immunosorbent assay (ELISA; pan-specific IFN-a, Mabtech Ab, Sweden; IL-12 (p70), BD OptEIA, BD Biosciences, USA; and TNF BD Biosciences, USA) as described recently. 25 For the current study, RNA in PBMC samples was preserved in RNAprotect TM (Qiagen, Hilden, Germany) and stored at À80°C. Samples were extracted for total RNA and sequenced at Macrogen Inc. (Seoul, South Korea) with TruSeq mRNA kit using the NovaSeq6000 platform, with a minimum of 20 million, 100 bp paired-end reads per sample. Two PBMC RNA samples were sequenced per individual: one unstimulated (baseline) and one RVstimulated. The data are available from Gene Expression Omnibus repository (accession number GSE166292). The Genome Informatics Group at QIMR Berghofer, Queensland, Australia, aligned the sequence reads. Sequence reads were trimmed for adapter sequences using Cutadapt [version 1.11 61 65 ]. Gene expression of the pDC marker CLEC4C 26 was quantified in whole blood for the main cohort using relative quantitative RT-PCR as described. 25 RNA from the validation samples was extracted, and transcriptional profiling was performed using Illumina Human HT-12 microarrays (SanDiego, CA, USA) at the microarray facility at the Diamantina Institute, University of Queensland as described. 37 Both RV and unstimulated PBMC samples were profiled. Illumina BeadStudio summary probe and summary control probe profiles were read into R 66 using the read.ilmn() function available in the limma package. 67 Background correction and normalisation were performed using the neqc() function 68 available in limma, which uses the normal exponential convolution model for background correction followed by quantile normalisation. Distribution of probes for each array was assessed using boxplots, before and after normalisation, and no outlying arrays were identified. Probes that are not expressed were filtered out prior to analysis using a threshold of probes that are expressed in at least three arrays according to the detection of P-values of 5%. 44 samples and 10 907 probes were included for further association analysis. EdgeR package 69 was utilised for the discovery of DEGs in the IFN-a-high group compared with the IFN-a-low group. Library size was corrected using counts per million (CPM), which involves dividing each sample gene count by the total number of mapped reads. Trimmed mean of M-values was used to normalise differences in RNA composition between samples with the function calcNormFactors() from the edgeR package. Function glmQLFit() fits a quasilikelihood negative binomial generalised log-linear model to count data to identify DEGs. Significant genes were filtered with multiple testing corrected FDR < 0.05 and log fold-change > 1. RV-stimulated and unstimulated samples were tested separately. Identification of significantly enriched biological pathways was performed using GSEA 70 focussing on DEGs between IFN-a-high and IFN-a-low groups in unstimulated and RVstimulated samples. Fold change in gene expression was used for ranked list of input to GSEA, and GO pathways with FDR < 0.25 were considered as biological pathways significantly related to the DEGs as recommended by GSEA. 70 All statistical analyses were performed with R (version 3.4.4 66 ). Variables were tested for normality and consequently treated as nonparametric variables. CreateTableOne() 71 function in the R package tableone was used to create the demographics table of the study group and perform statistical tests. The difference in sample distribution was tested with a nonparametric Mann-Whitney U-test. Linear regression was used to test association. A P-value < 0.05 was considered statistically significant after multiple testing adjustment where necessary. To determine the genes most likely associated with differences in BMI, two-way ANOVA was first used to exclude interactions between IFN-a producer groups and groups with BMI below/above 25 kg m À2 . Then, least absolute shrinkage and selection operator (LASSO) regression analysis was used to test BMI association with gene expression. Boxplots and before-after plots were created with GraphPad Prism software (version 7, La Jolla, USA). GOpathway networks were generated in Cytoscape interface [San Diego, USA 72 ], and other graphs were created with R using packages ggfortify, 73 ggplot2, 74 Pheatmap 75 and dendsort. 76 JW Upham received an unrestricted grant from AstraZeneca that partially funded this study. Within the last 5 years, he has received personal fees from GSK, Novartis, Sanofi and Boehringer Ingelheim that are unrelated to the submitted work. To the best of our knowledge, the other authors have no conflict of interest, financial or otherwise. Coinfection in SARS-CoV-2 infected patients: Where are influenza virus and rhinovirus/enterovirus? Chapter 6 -Exacerbations of chronic respiratory diseases Distinction between rhinovirusinduced acute asthma and asthma-augmented influenza infection How rhinovirus infections cause exacerbations of asthma Viral infections in allergy and immunology: How allergic inflammation influences viral infections and illness Interferon response and respiratory virus control are preserved in bronchial epithelial cells in asthma Rhinovirusinduced interferon production is not deficient in well controlled asthma Impaired innate interferon induction in severe therapy resistant atopic asthmatic children Reduced antiviral interferon production in poorly controlled asthma is associated with neutrophilic inflammation and highdose inhaled corticosteroids Upper airway cell transcriptomics identify a major new immunological phenotype with strong clinical correlates in young children with acute wheezing Rhinovirus stimulated IFN-a production: how important are plasmacytoid DCs, monocytes and endosomal pH? The goldilocks zone of type I IFNs: lessons from human genetics Induction and function of type I and III interferon in response to viral infection Interferon regulatory factor 7 is a major hub connecting interferon-mediated responses in virusinduced asthma exacerbations in vivo Single-cell analysis reveals that stochasticity and paracrine signaling control interferon-alpha production by plasmacytoid dendritic cells Self-priming determines high type I IFN production by plasmacytoid dendritic cells Integration of multi-omics data and deep phenotyping enables prediction of cytokine responses Variation in the human immune system is largely driven by non-heritable influences The cellular composition of the human immune system is shaped by age and cohabitation Common genetic variants modulate pathogen-sensing responses in human dendritic cells A functional genomics approach to understand variation in cytokine production in humans Genetic variants regulating immune cell levels in health and disease Host and environmental factors influencing individual human cytokine responses Respiratory viral infections and host responses; insights from genomics Risks for cold frequency vary by sex: role of asthma, age, TLR7 and leukocyte subsets CLEC4C gene expression can be used to quantify circulating plasmacytoid dendritic cells UniProt: a worldwide hub of protein knowledge Soluble CD93 in serum as a marker of allergic inflammation Genetic variation along the histamine pathway in children with allergic versus nonallergic asthma Transcriptional heterogeneity of mast cells and basophils upon activation Effect of variation in CHI3L1 on serum YKL-40 level, risk of asthma, and lung function Sputum hyaluronan and versican in severe eosinophilic asthma Rare SOX2 + airway progenitor cells generate KRT5 + cells that repopulate damaged alveolar parenchyma following influenza virus infection Female mice lacking Pald1 exhibit endothelial cell apoptosis and emphysema Fatty acid elongase 7 catalyzes lipidome remodeling essential for human cytomegalovirus replication Cell-based screen identifies human interferon-stimulated regulators of listeria monocytogenes infection Critical role of plasmacytoid dendritic cells in regulating gene expression and innate immune responses to human rhinovirus-16 MS4A4A: a novel cell surface marker for M2 macrophages and plasma cells Novel biomarkers of a peripheral blood interferon signature associated with drug-naive early arthritis patients distinguish persistent from self-limiting disease course Pattern-specific transcriptomics identifies ASGR2 as a predictor of hematogenous recurrence of gastric cancer Novel mechanisms and functions of complement Distinct indirect pathways govern human NK-cell activation by TLR-7 and TLR-8 agonists Interleukin-18 is associated with protection against rhinovirus-induced colds and asthma exacerbations The role of iron and siderophores in infection, and the development of siderophore antibiotics Viral-bacterial interactions in the respiratory tract Rhinovirus infection induces degradation of antimicrobial peptides and secondary bacterial infection in chronic obstructive pulmonary disease Rhinovirus exposure impairs immune responses to bacterial products in human alveolar macrophages Modulatory effects of low-dose hydrogen peroxide on the function of human plasmacytoid dendritic cells Regulation of type I interferon responses by mitochondria-derived reactive oxygen species in plasmacytoid dendritic cells Oxidative stress in asthma: Part of the puzzle Induced sputum 8-isoprostane concentrations in inflammatory airway diseases Mapping variation in cellular and transcriptional response to 1,25-dihydroxyvitamin D3 in peripheral blood mononuclear cells Vitamin B12 transport from food to the body's cells-a sophisticated, multistep pathway GRINL1A colocalizes with N-methyl Daspartate receptor NR1 subunit and reduces N-methyl D-aspartate toxicity NMDA receptors in immune competent cells Interferon-b1a modulates glutamate neurotransmission in the CNS through CaMKII and GluN2A-containing NMDA receptors Control of the RNA polymerase II phosphorylation state in promoter regions by CTD interaction domain-containing proteins RPRD1A and RPRD1B Distinctive roles of age, sex, and genetics in shaping transcriptional variation of human immune responses to microbial challenges The TLRmediated response of plasmacytoid dendritic cells is positively regulated by estradiol in vivo through cellintrinsic estrogen receptor a signaling Estrogen modulation of endosome-associated toll-like receptor 8: an IFNaindependent mechanism of sex-bias in systemic lupus erythematosus Cutadapt removes adapter sequences from high-throughput sequencing reads STAR: ultrafast universal RNA-seq aligner RNA-SeQC: RNA-seq metrics for quality control and process optimization RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome R: A language and environment for statistical computing limma powers differential expression analyses for RNA-sequencing and microarray studies Optimizing the noise versus bias trade-off for Illumina whole genome expression BeadChips edgeR: a Bioconductor package for differential expression analysis of digital gene expression data Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles tableone: Create 'Table 1' to Describe Baseline Characteristics Cytoscape: a software environment for integrated models of biomolecular interaction networks ggfortify: unified interface to visualize statistical results of popular R packages Elegant Graphics for Data Analysis pheatmap: Pretty Heatmaps dendsort: Modular Leaf Ordering Methods for Dendrogram Nodes Additional supporting information may be found online in the Supporting Information section at the end of the article.This is an open access article under the terms of the Creative Commons Attribution-NonCommercial-NoDerivs License, which permits use and distribution in any medium, provided the original work is properly cited, the use is non-commercial and no modifications or adaptations are made.