Genetic, intrinsic, and environmental determinants of innate immune cytokine responses in healthy four-year-old children
Introduction
Innate immune responses, including cytokine production, are central to immediate protection from infection and injury¹. They also contribute to the pathophysiology of chronic inflammatory diseases including atherosclerosis². Susceptibility to these infectious and noncommunicable diseases is partly determined by marked variation (including dysregulation) in innate immune responses between individuals, across the life course¹,³⁻¹⁰. In healthy adults, innate immune cytokine responses to ex-vivo stimulation are shaped by both genetic variation and environmental exposures such as diet and infection⁶⁻⁸,¹¹,¹². The heritability of these innate responses in adults varies by stimulus and cytokine, is generally higher than for adaptive responses⁷, and proportionally declines with age¹³.
In early life, notwithstanding rapid maturation and development, immune responses differ markedly from those in adults, with a relative reliance on innate immune responses rather than adaptive responses³,¹⁴⁻¹⁶. Despite their importance in childhood and for lifetime disease risk, the determinants of innate immune cytokine responses in childhood are poorly understood. Few studies have characterised innate immune responses in healthy children - existing paediatric data are largely derived from those with immune-related conditions (such as allergy or asthma)¹⁷,¹⁸ or from studies of the effects of specific immune-modulating agents (such as Bacillus Calmette–Guérin vaccine¹⁹⁻²¹). This limits understanding of how inter-individual variation in innate responses is established during early life.
To address this gap, we profiled innate immune cytokine responses in whole blood samples from a large sample of 4-year-olds in a deeply phenotyped longitudinal cohort with detailed participant information and repeated biological samples²². We investigated the contribution of genetic (quantitative trait locus [QTL]-specific and overall heritability), non-genetic host (e.g. birth factors, sex, anthropometry) and environmental (e.g. seasonal variation in population-level incidence of respiratory viral infection) determinants of innate cytokine responses at four years of age, and examined how of these responses relate to systemic inflammation and leukocyte composition.
Results
Innate immune stimulation reveals stimulus-specific and highly variable cytokine responses in early childhood
To investigate cytokine responses in four-year-old children, we used samples collected as part of the Barwon Infant Study (BIS), a population-derived Australian cohort of mainly (approximately 90%) European descent²². Briefly, fresh peripheral blood samples from 286 children were incubated for 24 hours with a panel of 8 pathogen mimetics or appropriate controls, and 13 innate immune cytokines were quantified (Figure 1A). Details regarding the stimuli and cytokines are provided in the Methods section (see also Table 1).
The cytokine production data were first visualised as an intensity-normalised heatmap, revealing distinct patterns for each stimulus (Figure 1B). The samples broadly clustered by culture conditions; negative controls (i.e. RPMI medium alone, or RPMI combined with the transfection reagent Lyovec), bacterial ligands (lipopolysaccharide [LPS], PeptidoGlycaN [PGN]), and most viral mimetics (polyinosinic:polycytidylic acid [Poly(I:C)], Lyovec/3’3’-cyclic GMP-AMP [cGAMP], Lyovec/5’-ppp double stranded [ds]RNA) showing distinct clusters of responses along the rows of Figure 1B.
In negative control culture supernatants, detectable cytokines were generally similar for RPMI and Lyovec/RPMI. These included platelet-derived growth factor (PDGF)-BB, monocyte chemoattractant protein (MCP)1, interleukin (IL)-8, vascular endothelial growth factor (VEGF)A, and interleukin-1 receptor antagonist (IL-1Ra). Compared to RPMI alone, Lyovec/RPMI induced interferon (IFN)α and C-X-C Motif Chemokine Ligand 10 (CXCL10) in some children, and IL-1Ra and VEGFA were slightly increased on average (Supplementary Figure 1).
Cytokine production following ligand stimulation generally showed marked inter-individual variation, with some responses varying across orders of magnitude (Supplementary Figure 1, prominently e.g., IL-6 in response to LPS). Furthermore, the strength of the responses differed according to the stimulus. For example, as seen in Figure 1B, Tumour necrosis factor (TNF)²³, IL-10, IL-6, and IL-1β were strongly induced by LPS, PGN, and R848 (resiquimod), whereas IFNα, IFNγ, and CXCL10 were predominantly produced in response to Lyovec/cGAMP and Poly(I:C) in addition to R848. VEGF-A and IL-8, and to a lesser extent also PDGF-BB were – in contrast to interferons and CXCL10 – most strongly produced in response to PGN, LPS, and dsRNA.
Principle component analysis confirmed the clear separation of the responses dependent on the specific stimulus. In principal component analyses (PCA), Lyovec/RPMI and RPMI alone were largely superimposed indicating similar overall responses (Figure 1C-D). R848 was separated the furthest from RPMI across both PC1 and PC2, reflecting the potency of this stimulus, but the separation along PCs differed for the other stimuli. LPS, PGN, and Lyovec/dsRNA showed the strongest shift from control samples along PC1, whereas Poly(I:C) Lyovec/cGAMP separated from RPMI on PC2. For PC3, Lyovec/dsRNA and Poly(I:C) diverged from the other stimuli.
In summary, we observed considerable inter-individual variation in whole blood cytokine responses to innate stimuli. Each stimulus provoked a distinct response, aligning with expectations based on previous literature from adults⁶.
Common genetic variants explain substantial variation in cytokine responses in early childhood
We first considered genetic determinants of variation in innate immune cytokine responses (Figure 2A). We performed cytokine quantitative trait loci (cyQTL)-mapping to investigate the association between genetic variants (single nucleotide polymorphisms, SNPs) and cytokine responses for cytokine-stimulus combinations that were within quantifiable range in at least 50% of the samples (see also Supplementary Figure 1). These analyses were performed using data from a subset of children due to missingness of covariate data (n = 259 [90.6% of total sample]; see methods for details). We set three thresholds for statistical significance: (i) a lenient screening threshold (p < 5x10⁻⁶); (ii) the commonly used ‘genome-wide’ significance (p < 5x10⁻⁸); and (iii) ‘study-wide’ significance for which the genome-wide significance threshold was made more strict to account for the number of effective comparisons²⁴,²⁵ (here: determined to be p < 1.39x10⁻⁹).
A single locus, rs55792153, met the study-wide significance threshold (Figure 2B). This SNP is located immediately downstream of TMEM173, which encodes stimulator of interferon genes (STING), the receptor for cGAMP (Figure 2C). Compared to major allele homozygotes (CC), those with an AA or heterozygous genotype produced lower IL-1β, IL-1Ra and TNF following stimulation with Lyovec/cGAMP (Figure 2D-F). Most other SNPs meeting or approaching the genome-wide (but not study-wide) significance threshold were not located near canonical immune-related genes. The second-most significantly associated SNP was rs2343196, located downstream of LRRIQ3, which is mainly expressed in the testes but also in many immune cells²⁶,²⁷ (Figure 2B, 2G). Compared to reference allele (TT) homozygotes, children heterozygous or homozygous for the variant allele C produced less IL-10 on stimulation with Lyovec/dsRNA (Figure 2H). Additional SNPs with or approaching genome-wide significance are shown in the supplementary data and include rs72679554 (near RP11-536K17.1, Supplementary Figure 2A-C), rs2384950 (strong linkage disequilibrium [LD] with a dense gene cluster; Supplementary Figure 2D-E), rs2300875 (in ACTN1, Supplementary Figure 2F-G), and rs6112096 (nearest to DTD1; Supplementary Figure 2H-I). Together, these data reinforce that specific SNPs markedly impact specific stimulus-cytokine combinations (or sometimes more broadly the cytokine response, depending on the locus).
We next investigated the total variance in cytokine production that could be explained by SNPs. To avoid overfitting given the large number of stimulus-output combinations relative to the sample size, we focused on the top 50 (LD-independent, by p-value) associated SNPs per stimulus-cytokine combination. Overall, genetic variation appears to be a strong determinant of cytokine responses in preschool children; the combination of the top 50 SNPs explained approximately 20-45% of inter-individual variation for most cytokine-stimulus combinations (Figure 2I). MCP1 responses to Lyovec/cGAMP and R848 stimulation showed the highest percent of variance explained (respectively 43.7% [95%CI: 34.6 – 53.3] and 43.5% [95%CI: 33.3 – 54.6], Figure 2I). In contrast, for LPS- or R848-stimulated levels of VEGF-A, almost none of the variation was explained by the top 50 SNPs.
We then performed gene-set (Reactome) enrichment analyses of SNPs within a window of 35 kb upstream to 10 kb downstream of genes to investigate potential shared pathways across stimulus-cytokine combinations. Pathways enriched for SNPs associated with baseline cytokine production included “Biosynthesis of maresins”/“Biosynthesis of maresin-like Specialised Pro-resolving Mediators” (RPMI) and FLT3-related gene-sets such as “FLT3 signaling through Src family kinases” (Figure 2J, Supplementary Figure 2J).
Across the different stimuli, some well-understood signalling pathways were frequently enriched. Chloride transporter-related pathways were shared between LPS and R848 (Figure 2K, Supplementary Figure 2J). PI3K-related signalling (e.g., “Activated NTRK2 signals through PI3K”, “MET activates PI3K signaling”) was shared between PGN and Lyovec/cGAMP (Figure 2L, Supplementary Figure 2J). The cytokines IL-1β, IL-6, TNF, IL-10, and IL-1Ra often shared pathways within a stimulus, highlighting their central roles in innate immune responses and suggesting that their expression is co-regulated to a degree.
Dimensionality reduction using CytoMod – Cytokine co-expression modules identify recurring response patterns
We next investigated additional host factors and environmental determinants of childhood cytokine responses. First, we used a freely available python module, Cytomod²⁸, to reduce the number of dimensions and hence the multiple testing burden. Briefly, CytoMod clusters cytokines into ‘modules’ based on patterns of co-expression (Figure 3A). The modules are determined by unsupervised clustering of cytokines after adjusting for the participant-level mean cytokine value. The final module composition is then based on pairwise reliability scores using a bootstrapping approach. Module expression scores are calculated by taking the mean of standardised cytokine concentrations included in that module²⁸.
We performed the clustering for each stimulus separately (Figure 3B-I), as each stimulus affects expression of cytokines differentially. In these analyses, CytoMod assigned between 3 and 5 cytokine modules per stimulus (Table 2), reducing the number of parallel comparisons from 104 stimulus-cytokine combinations to 33 cytokine module expression scores. The correlation between cytokines within each module, as well as between each cytokine and its module expression score are shown in Supplementary Figure 3. While the cluster composition differed between stimuli, some general patterns were evident; TNF, IL-6, IL-1β, and IL-10 responses tended to cluster together, as did IFNα and IFNγ, as well as IL-8, MCP-1, VEGFA, and PDGF-BB. For each stimulus, cytokines that were or were not released also tended to cluster together.
Throughout the following sections, we present associations between variables of interest and the expression scores of cytokine modules as primary analysis, with p-values adjusted for multiple testing using the FDR method. As this is a data-scarce field, in secondary analyses we report associations between variables of interest and individual stimulus-cytokine combinations with nominal p < 0.05.
Limited associations between child characteristics and cytokine responses at four years of age
We examined associations between exact age, sex, and measures of adiposity (measured at blood collection) and cytokine module expression (Figure 4A). For these analyses we used data from 286 children (136 girls and 150 boys; Figure 4B) approximately four years of age (range: 3.9 – 5.6 years; median: 4.2 years; Figure 4C). Anthropometric and adiposity measures (BMI [n = 286] and body fat percentage [n = 269]) were mostly within the normal range for their sex and age²⁹ (Figure 4D, E). After FDR correction, we found no evidence of associations between these anthropometric measurements and cytokine module expression (Figure 4F). Given the data scarcity in this field, we additionally reported nominal p-values < 0.05 for these analyses (Figure 4F-I).
We did not find evidence that exact age (within the narrow available age-range, Figure 4C) was associated with any of the cytokine responses (Figure 4F, G). Sex differences in immune responses are well-described in adults (see for example references 3,6,30-33), but we found only a small number of differences in this cohort: modules LPS-3 and R848-3 were higher in male children, whereas R848-5 was higher in female children (Figure 4F). These differences were driven by MCP1 (included in both LPS-3 and R848-3; Figure 4H, J) and IFNα (included in R848-5; Figure 4H, K), respectively – well-known examples from adult studies³²,³³. We found evidence that measures of adiposity were associated with cytokine responses, although the effect size was small in this cohort of predominantly normal-weight children. The BMI z-score (derived from WHO standards data²⁹) was positively associated with PGN-1 levels and inversely with PGN-4. There was also a negative association between BMI z-score and Lyovec/dsRNA-4. Finally, the module Poly(I:C)-2 as well as more broadly individual cytokine responses to poly(I:C) were associated with BMI z-score (Figure 4F, I, L, M). Body fat percentage showed weak evidence of associations patterns, similar to BMI z-score (Figure 4F).
Pregnancy and perinatal variables show weak or non-significant associations with cytokine responses
As pregnancy and perinatal factors may influence immune development¹⁶,³⁴, we investigated associations between mode of birth, birth weight, and gestational age and cytokine modules (Figure 5A). Most children were born at term (37 to 42 weeks gestational age; median: 39.6 weeks, range: 32.1 – 41.9 weeks; Figure 5B), and with a birth weight within the normal range (median: 3.53 kg, range: 1.61 – 5.41; Figure 5C). Of the 286 children included in these analyses, 102 (35.6%) were born by caesarean section (Figure 5D).
We found no statistically significant associations between cytokine module expression and birthweight or gestational age when analysed separately, nor for the birthweight z-score derived per sex and accounting for gestational age³⁵ (Figure 5E). We found weak evidence for associations between these exposures and production of individual cytokines (Supplementary Figure 4A-C).
We also found no statistically significant associations between mode of birth and cytokine module expression. However, we did observe that the model estimates were mostly positive, suggesting that cytokine modules were generally higher in children born by caesarean section (Figure 5E). This is also visualised at the level of individual cytokines in Figure 5F, as a right-skewed volcano plot. We hypothesised that there may be a subgroup of children whose cytokine responses were more strongly associated with caesarean birth than others. We therefore investigated if the experience of labour, an inflammatory process³⁶, prior to caesarean section was associated with differences in cytokine module expression or individual cytokine responses. We observed that, while the effects were statistically non-significant, there was some evidence that the pattern of differences in cytokine modules as well as individual cytokines were more pronounced in children born via caesarean section after experiencing labour (n = 36; Figure 5G [caesarean section with labour versus vaginal birth], Supplementary Figure 4D [subgroup comparisons]). The evidence was strongest for modules R848-5 and Lyo/cGAMP-1, which include IFNα (Figure 5H) and TNF (Figure 5I), respectively. Given the limited sample size of caesarean subgroups, these analyses should be interpreted as exploratory.
Systemic inflammation and leukocyte composition are strongly associated with cytokine production capacity
We next investigated how inflammatory biomarkers, frequently assayed in cohort studies in favour of more labour-intensive stimulation assays, are related to whole blood cytokine production (Figure 6A). We considered three plasma measures of systemic inflammation; high-sensitivity C-reactive protein (hsCRP, n = 280), a commonly measured biomarker of acute systemic inflammation³⁷; glycoprotein acetyls (GlycA, n = 283), a composite biomarker that is more stable than hsCRP and better captures chronic systemic inflammation³⁸⁻⁴⁰; and granulocyte-to-lymphocyte ratio (GLR, n = 265), which is increased during acute inflammation and trained immunity due to increased granulocyte production⁴¹,⁴². Forty-seven (16.8%) of the children had hsCRP level concentrations that were below the limit of detection (0.001 µg/mL), whereas all GlycA concentrations were within the quantifiable range (Figure 6B). GLR varied substantially across the population (median 1.44, range: 0.44 – 4.46; Fig 6C-D, Supplementary Figure 5A-C). All three biomarkers were log-transformed and internally standardised for statistical analyses.
In total, associations between any of the three inflammatory markers and the various cytokine modules were overwhelmingly positive: we found 3 negative associations and 30 positive associations (Figure 6E). We observed considerable overlap between the 7 cytokine modules associated with hsCRP and the 13 modules associated with GlycA (Figure 6E). GlycA was also associated with more individual cytokines (Figure 6F). Only a single module (RPMI-2) was (negatively) associated with hsCRP but not GlycA (Figure 6E, Supplementary Figure 5D). Generally, an increase of 1 standard deviation in GlycA levels was associated with a larger increase or decrease in module expression score, compared to a similar increase in hsCRP (Supplementary Figure 5D).
GLR showed the strongest associations overall and was associated with 13 cytokine modules, more than GlycA (Figure 6E, 6G, Supplementary Figure 5E). Both biomarkers were often associated with modules containing PDGF-BB, VEGF-A, and IL-1Ra. However, while most of the modules overlapped between GlycA and GLR, different cytokines within these modules drove these associations. We therefore observed that a substantial number of stimulus-cytokine pairs were either associated with GLR only, or only with GlycA. For stimulus-cytokine pairs that were associated with both biomarkers, associations were more evident for GLR (Figure 6G). The associations with GLR were strongly driven by IL-1Ra and VEGF-A, whereas GlycA was linked with PDGF-BB and a mixture of other cytokines (Figure 6H-I).
In summary, we found that all three markers of inflammation were associated with cytokine module expression, but to different extents and that associations were driven by different individual cytokines.
As GLR showed the strongest association with any given module, we explored associations between cytokine responses and relative abundance of lymphocytes, monocytes, and granulocytes (predominantly neutrophils). Compared to GLR, percentage of granulocytes showed greater evidence of associations with cytokine modules that incorporate VEGF-A and IL-1Ra (Supplementary Figure 5F, G). As anticipated given the inverse relationship between granulocyte and lymphocyte percentages, this was mirrored by the lymphocyte percentage (Figure 6B, Supplementary Figure 5F, H). The percentage of monocytes was more stable across participants, when considering percentage-point differences (Figure 6D, Supplementary Figure 5C), and was associated (at the individual cytokine level) with production of classical monocyte cytokines such as IL-6 and IL-1β, and also IFNα (Supplementary Figure 4I). For additional context, it is important to mention that the IL-1β/IL-1Ra ratio (which has been proposed as a measure of biological IL-1β activity⁴³,⁴⁴) was low across stimuli (Supplementary Figure 1).
Discussion
In this population-derived cohort of children around 4 years of age, we demonstrated substantial inter-individual variation in early life cytokine responses and quantified the relative contributions of genetic, host-intrinsic, and selected environmental determinants. Genetic variants explained the largest proportion of variation in cytokine production capacity, followed by strong associations with seasonal viral infections and cell type composition. Biomarkers of inflammation associated with cytokine responses in marker-specific patterns, likely reflecting their underlying origins. We found some evidence that mode of birth in combination with the experience of labour may have long term effects on cytokine responses.
We identified several genomic loci associated with cytokine production capacity. The strongest association mapped to the TMEM173/STING locus, where variants were linked to a marked decrease in production of IL-1β, TNF, and IL-1Ra. Direct comparison with adult cyQTL studies was not possible as these have not used cGAS/STING-specific stimuli. However, follow-up studies are warranted given the central role of the cGAS-STING pathway in host defence, and evasion of STING being a common pathogenic mechanism⁵⁰.
Unlike cohorts of Dutch descent⁷,²⁴,²⁵, we did not identify cyQTLs in the TLR1-TLR6-TLR10 cluster. This could reflect ancestry differences as the TLR1-TLR6-TLR10 locus differs markedly between populations due to introgression events⁵¹⁻⁵³. A paediatric cohort from Tanzania, where these introgressions did not take place, also did not identify this as a key locus²¹. An age-related effect is unlikely but cannot be excluded as both cohorts that did not identify this locus are paediatric. Overall, the top 50 LD-independent SNPs explained a substantial proportion (~20-45%) of variance in cytokine production. Compared to the most analogous study in adults⁸, we found comparable or higher percentages of variance explained – but in a whole blood stimulation assay, which inherently has more non-genetic variation than the PBMC-stimulations⁵⁴ reported in adults. This aligns with the genetic influence on cytokine responses decreasing with age¹³, which likely reflects the increasing cumulative exposure to each individual’s uniquely diverging environments. A more rigorous comparison with adult studies was not possible due to differences in experimental design (cell type, stimuli, analytes) and analytical approaches.
In contrast to adult data, sex-related differences in cytokine responses were modest, consistent with the pre-pubertal age of this cohort³,³⁰,³¹,³³,⁵⁵⁻⁵⁸. This is unlikely to be explained by limited statistical power as several adult studies reporting sex differences were of similar or smaller size. We observed differences in MCP1 production (higher in male children), suggesting a genetic mechanism beyond the hormonal influences on circulating MCP1 levels in adults³³. The observed increase in TLR7/8-driven IFNα responses in female children aligns with known X-linked regulation of antiviral sensing and supports that sex differences in innate immunity are detectable even in early life, albeit at smaller magnitude than in adulthood⁵⁶,⁵⁹.
There was some evidence that children delivered by caesarean section had stronger cytokine responses, particularly the subgroup that had been exposed to labour (presumably emergency caesarean sections). This aligns with evidence linking caesarean delivery with higher rates of childhood infection-related hospitalisation⁶⁰ and increased risk of several inflammatory disorders⁶¹. Further studies incorporating clinical indications for caesarean delivery are needed to clarify these associations. As the associations at the module level were not significant after correction for multiple testing, and p-values for individual cytokines were unadjusted, these findings should be interpreted cautiously and require validated in independent cohorts.
Our findings regarding inflammation biomarkers reflect differences between cell-intrinsic cytokine production capacity and differences due to cellular composition. In related work, we have previously shown associations between GlycA and LPS/PGN-induced monocyte-associated cytokines⁶². Here, we show GlycA and GLR were both associated with various cytokine modules, but they did so through distinct cytokine signatures. Notably, IL-1Ra and VEGF-A were strongly associated with GLR and granulocyte abundance, the same cytokines that were relatively poorly explained by the top SNPs. This suggests that some cytokine responses are influenced more by cellular composition than (genetic) cell-intrinsic capacity, although both cell composition and ‘per-cell’ responsiveness are largely determined by the state of bone marrow progenitors⁴²,⁶³.
Indicators of respiratory viral infection incidence at a population level were associated with differences in cytokine responses in pre-school children, particularly responses to stimulation of antiviral pathways. For these analyses we considered positive PCR for any of the tested respiratory viruses, but it would be valuable to consider effects of specific viruses in appropriately sized follow-up cohorts. Nonetheless, cumulative or non-specific infection burden remains highly relevant: we have previously shown in BIS that total infection burden is associated with adverse plasma profiles of lipids and metabolites as early as 12 months of age⁶⁴. These findings suggest that infection burden is a potentially modifiable determinant of long-term health⁶⁵.
Together, these findings indicate that early childhood is a critical period during which innate immune responses are shaped strongly by genetic variation and to a more limited extent by the other host and exogenous environmental exposures considered in this study. A growing body of evidence suggests that childhood inflammation contributes to later cardiometabolic disease risk⁶⁶⁻⁶⁹, together with increasingly common⁷⁰,⁷¹ traditional (and pro-inflammatory⁷²,⁷³) risk factors such as obesity and type 2 diabetes mellitus. Furthermore, infection burden during childhood is an emerging cardiovascular risk factor⁶⁴,⁷⁴. Understanding the determinants of variability in innate immune responses during this critical window is essential to identification of – and intervening in – children at risk before clinical disease emerges⁶⁹,⁷⁵,⁷⁶.
Limitations
We acknowledge a number of limitations. First, genetic analyses are traditionally performed in much larger cohorts to identify rarer variants or variants with smaller effect sizes. Nonetheless, prior studies⁷,⁸,⁵⁴,⁷⁷ suggest that moderate sample sizes can identify impactful genetic variants with larger effects. Second, the flow cytometry analyses used in this study were performed using very limited cell type markers, as they were designed at inception of this cohort when advanced methodology was not available at the study site. Third, the ecological, population-level estimate of respiratory viral infections was biased towards the Greater Melbourne area, rather than specifically in the Barwon area where the study was conducted. As the analysis of cytokines and viral exposure is ecological, the viral exposure data is not specific to participating children and their families. Further studies should also incorporate data on exposure to specific pathogens, and the range of pathogens should be extended to include other viral infections including gastrointestinal viruses. Finally, in this study, biomedical assessments were performed at approximately four years of age only, meaning we lacked a range of ages to identify meaningful effects of age difference, which could be an important determinant of immune responses in early childhood (especially considering that relatively modest age differences in young children constitute a larger proportion of lifetime than in adults). Similarly, BMI and body fat percentage were generally within a normal-range in this cohort, limiting our capacity to identify effects of these exposures on immune responses. We recognise that there is an immense number of potential environmental factors that may influence innate immune responses and that we have studied a small subset of those.
Acknowledgements
We would like to thank all participants of the Barwon Infant Study, including their parents. Special thanks also go to the entire BIS team who enabled sample collection and provided essential laboratory support. The establishment work and infrastructure for the BIS was provided by the Murdoch Children’s Research Institute (MCRI), Deakin University and Barwon Health. Subsequent funding was secured from the National Health and Medical Research Council of Australia, The Jack Brockhoff Foundation, the Scobie Trust, the Shane O’Brien Memorial Asthma Foundation, the Our Women’s Our Children’s Fund-Raising Committee Barwon Health, The Shepherd Foundation, the Rotary Club of Geelong, the Ilhan Food Allergy Foundation, GMHBA Limited and the Percy Baxter Charitable Trust, Perpetual Trustees and the Minderoo Foundation. In-kind support was provided by the Cotton On Foundation and CreativeForce. Research at Murdoch Children’s Research Institute is supported by the Victorian Government’s Operational Infrastructure Support Program. We also thank Dr. rer. nat. Cédric Scherer for providing excellent resources on code-first data visualisation practices ( https://www.cedricscherer.com/ ). This work relied heavily on the often under-appreciated work of software developers and maintainers. This work was made possible by funding from the Niels Stensen Fellowship. SnotWatch is made possible through the efforts of the SnotWatch collaboration group. The SnotWatch collaboration group includes Monash Pathology (Tony Korman), Royal Children’s Hospital Pathology (Andrew Daley, Vanessa Clifford), Alfred Pathology (Adam Jenney), Royal Melbourne Hospital Pathology (Katherine Bond), Eastern Health Pathology (Roy Chean), Northern Pathology Victoria (Yvonne Hersusianto), Barwon Health (Eugene Athan) and the Victorian Department of Health (Jim Black).
Continued extraction from pages 11-20
Software continued from previous section:
Luo, B. ggheatmap: Plot Heatmap. (2025). Tiedemann, F. gghalves: Compose Half-Half Plots Using Your Favourite Geoms. (2025). Pedersen, T. L. ggforce: Accelerating ‘ggplot2’. (2025). Kay, M. ggdist: Visualizations of Distributions and Uncertainty in the Grammar of Graphics. IEEE Transactions on Visualization and Computer Graphics 30, 414–424 (2024). https://doi.org:10.1109/TVCG.2023.3327195 Brand, T. v. d. ggh4x: Hacks for ‘ggplot2’. (2025). Zeileis, A. et al. colorspace: A Toolbox for Manipulating and Assessing Colors and Palettes. Journal of Statistical Software 96, 1–49 (2020). https://doi.org:10.18637/jss.v096.i01 Wickham, H., Pedersen, T. L. & Seidel, D. scales: Scale Functions for Visualization. (2025). Pedersen, T. L. patchwork: The Composer of Plots. (2025). Landis, J. ggside: Side Grammar Graphics. (2025). Slowikowski, K. ggrepel: Automatically Position Non-Overlapping Text Labels with ‘ggplot2’. (2024). Schloerke, B. et al. GGally: Extension to ‘ggplot2’. (2025). Sinke, L., Cats, D. & Heijmans, B. T. Omixer: multivariate and reproducible sample randomization to proactively counter batch effects in omics studies. Bioinformatics (2021). Müller, K. hms: Pretty Time of Day. (2025). Wickham, H. & Bryan, J. readxl: Read Excel Files. (2025). Barbone, J. M. & Garbuszus, J. M. openxlsx2: Read, Write and Edit ‘xlsx’ Files. (2025). Iannone, R. et al. gt: Easily Create Presentation-Ready Display Tables. (2025). Ushey, K., Allaire, J. J. & Tang, Y. reticulate: Interface to ‘Python’. (2025). Schumacher, D. anthro: Computation of the WHO Child Growth Standards. (2026). Schumacher, D. anthroplus: Computation of the WHO 2007 References for School-Age Children and Adolescents (5 to 19 Years). (2024). Lewis, M. J. & Wang, S. locuszoomr: an R package for visualising publication-ready regional gene locus plots. Bioinformatics Advances, vbaf006 (2025). https://doi.org:10.1093/bioadv/vbaf006 Rainer, J. EnsDb.Hsapiens.v75: Ensembl based annotation package. (2017). Morgan, M. & Ramos, M. BiocManager: Access the Bioconductor Project Package Repository. (2025). Arora, S., Morgan, M., Carlson, M. & Pagès, H. a. e. GenomeInfoDb: Utilities for manipulating chromosome names, including modifying them to follow a particular naming style. (2025). Lawrence, M. et al. Software for Computing and Annotating Genomic Ranges. PLoS Computational Biology 9 (2013). https://doi.org:10.1371/journal.pcbi.1003118 Ochoa, A. genio: Genetics Input/Output Functions. (2023). Privé, F., Aschard, H., Ziyatdinov, A. & Blum, M. G. B. Efficient analysis of large-scale genome-wide data with two R packages: bigstatsr and bigsnpr. Bioinformatics 34, 2781–2787 (2018). https://doi.org:10.1093/bioinformatics/bty185 Barrett, T. et al. data.table: Extension of data.frame . (2025). Shabalin, A. A. Matrix eQTL: ultra fast eQTL analysis via large matrix operations. Bioinformatics 28, 1353–1358 (2012). Microsoft & Weston, S. foreach: Provides Foreach Looping Construct. (2022). doMC: Foreach Parallel Adaptor for ‘parallel’ (2022). de Leeuw, C. A., Mooij, J. M., Heskes, T. & Posthuma, D. MAGMA: Generalized Gene-Set Analysis of GWAS Data. PLOS Computational Biology 11, e1004219 (2015). https://doi.org:10.1371/journal.pcbi.1004219 Purcell, S. et al. PLINK: a tool set for whole-genome association and population-based linkage analyses. Am. J. Hum. Genet. 81, 559-575 (2007). https://doi.org:10.1086/519795 Ragueneau, E. et al. The Reactome Knowledgebase 2026. Nucleic Acids Res 54, D673-d681 (2026). https://doi.org:10.1093/nar/gkaf1223 Niaid Visual Medical Arts. Gene Mutation, https://bioart.niaid.nih.gov/bioart/170 (2024). Niaid Visual Medical Arts. H5N1, https://bioart.niaid.nih.gov/bioart/187 (2024).
Continuing from CHUNK 2
Cytokine-QTLs for cGAMP and dsRNA stimulation
We identified 259 independent SNPs across the genome that were associated with cytokine responses at the lenient screening threshold (p < 5×10⁻⁶). Of these, 23 SNPs reached genome-wide significance (p < 5×10⁻⁸) and 3 reached study-wide significance (p < 1.39×10⁻⁹) (Figure 2B). The most significant findings were:
rs55792153 (chromosome 5), which showed study-wide significance for IL-1β in response to Lyovec/cGAMP stimulation. This SNP is located in a region containing several genes involved in immune regulation. Children homozygous for the alternative allele (AA genotype) showed significantly higher IL-1β responses compared to those with the CC genotype (Figure 2D-F). The same SNP also showed associations with TNF and IL-1Ra responses to the same stimulus.
rs2343196 (chromosome 1), which showed genome-wide significance for IL-10 responses following Lyovec/dsRNA stimulation. Children with different genotypes at this locus showed marked differences in IL-10 production (Figure 2G-H).
Heritability and pathway enrichment
The top 50 LD-independent SNPs explained between 10-50% of the variance in cytokine responses across different stimulus-cytokine combinations (Figure 2I). The percentage of variance explained varied considerably by stimulus and cytokine measured. Stimulation with LPS and R848 generally showed higher heritability estimates compared to RPMI control conditions.
Reactome pathway enrichment analysis identified several biological pathways that were enriched among the SNPs associated with cytokine responses (Figure 2J-L). For RPMI stimulation, the most enriched pathways included TNF signalling, NF-κB activation, and apoptosis regulation. For LPS stimulation, enriched pathways included mitochondrial function, complement activation, and pattern recognition receptor signalling. For Lyovec/cGAMP stimulation, enriched pathways included interferon signalling and innate immune response pathways.
Results (continued)
Identification of cytokine response modules
To identify co-regulated patterns of cytokine production, we developed a novel computational approach called CytoMod (Cytokine Module detection) (Figure 3A). This method combines pairwise correlation analysis with bootstrapping and reliability scoring to identify groups of cytokines that are consistently co-expressed across individuals. This approach identified 33 distinct modules across all 104 stimulus-cytokine combinations tested.
For RPMI control stimulation, four modules were identified: Module 1 contained pro-inflammatory cytokines (IL-10, IL-1β, IL-6, TNF); Module 2 contained interferons (IFNα, IFNγ); Module 3 contained chemokines and growth factors (IL-1Ra, IL-8, MCP1); and Module 4 contained anti-inflammatory and growth-promoting factors (IL-5, CXCL10, PDGF-BB, VEGF-A) (Figure 3B).
For LPS stimulation, three modules were identified with distinct inflammatory signatures (Figure 3C). For PGN stimulation, five modules were identified (Figure 3D). For Poly(I:C) stimulation, three modules were identified (Figure 3E). For R848 stimulation, five modules were identified (Figure 3F). For Lyovec/RPMI, five modules were identified (Figure 3G). For Lyovec/cGAMP, four modules were identified (Figure 3H). For Lyovec/dsRNA, four modules were identified (Figure 3I).
Module composition showed both stimulus-specific and common patterns. For example, pro-inflammatory cytokines (IL-1β, IL-6, TNF) frequently clustered together across different stimuli, whereas interferon responses showed variable clustering depending on the stimulus used.
Associations with anthropometric measures at four years of age
At age four years, we assessed the relationship between innate immune cytokine responses and anthropometric measures including BMI and body fat percentage. Among the 286 children in the original cohort, 256 (n=256) had available anthropometric data at four years of age, with a mean age of 4.3 years (SD 0.2) (Figure 4A-C). The cohort included 128 girls and 128 boys (Figure 4B). The distribution of BMI and body fat percentage are shown in Figure 4D-E.
Sex differences in cytokine responses were identified for several stimulus-cytokine combinations (Figure 4F-H). Boys showed significantly higher MCP1 responses following R848 stimulation compared to girls (p = 0.003) (Figure 4J). Boys also showed higher IFNα responses following R848 stimulation (p = 9.2×10⁻⁴) (Figure 4K).
BMI was associated with IL-8 responses following Poly(I:C) stimulation (p = 0.008) (Figure 4L). Children with higher BMI also showed elevated CXCL10 responses following Poly(I:C) stimulation (p = 0.013) (Figure 4M). The associations between BMI and cytokine responses appeared to be specific to certain stimulus-cytokine combinations, with no significant associations identified for RPMI or LPS stimulation.
Age at the time of blood collection was not significantly associated with cytokine module expression scores (Figure 4G).
Associations with perinatal and pregnancy variables
We examined whether perinatal and pregnancy-related variables were associated with innate immune responses measured at four years of age (Figure 5A). Among the 256 children with anthropometric data, 248 had complete perinatal information including mode of birth and experience of labour prior to caesarean delivery.
The mean gestational age at birth was 39.5 weeks (SD 1.5, range 37-42) (Figure 5B). Mean birthweight was 3,570 grams (SD 530, range 2,480-5,150) (Figure 5C). The majority of children (n=173, 70%) were born via vaginal delivery, while 75 children (30%) were born via caesarean section (Figure 5D).
Volcano plot analysis revealed significant differences in multiple cytokine responses between children born via vaginal delivery and caesarean section (Figure 5F). Children born via caesarean section showed significantly lower TNF responses following Lyovec/cGAMP stimulation compared to those born vaginally (Figure 5H, p = 0.05). Children born via caesarean section also showed lower IFNα responses following R848 stimulation (Figure 5I).
Among children born via caesarean section, those who had experienced labour prior to caesarean delivery (n=31) showed intermediate cytokine responses between vaginal-delivered and elective caesarean-delivered children (Figure 5G). This pattern was observed for multiple cytokine-stimulus combinations including VEGF-A, IL-5, and IL-1Ra responses (Figure 5G).
Gestational age and birthweight were not significantly associated with cytokine module expression scores (Figure 5E).
Associations with markers of systemic inflammation
We next examined associations between innate immune cytokine responses and markers of systemic inflammation measured at four years of age, including GlycA, high-sensitivity C-reactive protein (hsCRP), and the granulocyte-to-lymphocyte ratio (GLR) (Figure 6A).
GlycA and hsCRP concentrations were positively correlated (r=0.54, p<0.001) (Figure 6B). The granulocyte-to-lymphocyte ratio showed a median value of 1.38 (IQR 1.13-1.77) (Figure 6C). The relative proportions of major leukocyte populations showed granulocytes comprised approximately 40-80% of total leukocytes across the cohort (Figure 6D).
FDR-adjusted regression analysis identified significant associations between GlycA and multiple cytokine module expression scores (Figure 6E). GlycA was positively associated with IL-1β, TNF, IL-6, IL-8, and MCP1 responses across multiple stimuli (Figure 6H). The magnitude of associations with GlycA and hsCRP were generally concordant (Figure 6F), supporting the validity of these associations.
When GlycA was increased by 1 SD, IL-1β responses to LPS stimulation increased by approximately 0.3 log₂ fold change (Figure 6F-H). Similar associations were observed for TNF, IL-6, and MCP1 responses.
The granulocyte-to-lymphocyte ratio showed concordant patterns with GlycA for many cytokines (Figure 6G), though with somewhat smaller effect sizes. Associations between GLR and cytokine responses were identified for IL-1Ra, IL-1β, TNF, and IL-8 responses (Figure 6I).
Notably, VEGF-A and PDGF-BB responses were negatively associated with both GlycA and GLR (Figure 6F-G, H-I), suggesting that systemic inflammation is associated with reduced growth factor responses.
Seasonal community infections and cytokine responses
We examined whether seasonal community infections were associated with innate immune cytokine responses by analysing the number of positive viral swabs collected from children in the community during the month of blood collection (Figure 7A).
The cohort experienced clear seasonal peaks in community viral infections, with the highest burden occurring during winter months (June-August in the Southern Hemisphere) (Figure 7B-C). Cytokine module expression scores showed seasonal variation that coincided with peaks in positive viral swabs in the community (Figure 7B).
Modules showing significant seasonal variation included Poly(I:C)-3 (IFNα, IL-1Ra, IL-5, CXCL10), R848-4 (IFNγ, CXCL10), Lyovec/cGAMP-4 (IFNα, IL-1Ra, CXCL10), and Lyovec/dsRNA-2 (IL-10, IL-1β, IL-1Ra, CXCL10, TNF) (Figure 7B).
When positive viral swabs in the community increased by 1 SD, IL-1Ra responses to Lyovec/dsRNA stimulation increased by approximately 0.5 log₂ fold change (nominal p = 0.003) (Figure 7G). IL-1β responses also showed a positive association with community viral burden (p = 2.1×10⁻⁴) (Figure 7H). IL-10 responses showed a weaker association (p = 3.5×10⁻⁴) (Figure 7I).
These associations remained significant when adjusting for age, sex, BMI, and technical variables, suggesting that current exposure to community infections influences innate immune responses measured in vitro.
Discussion
In this large prospective study of 286 pre-school children, we comprehensively characterised innate immune cytokine responses to multiple pathogen-associated molecular patterns and examined genetic, anthropometric, perinatal, and environmental factors associated with these responses. Our key findings are:
Distinct cytokine response patterns to different stimuli : Whole blood stimulation with different PAMPs resulted in distinct patterns of cytokine production that could be identified through principal component analysis and hierarchical clustering. These patterns reflect the differential activation of distinct innate immune pathways by different stimuli.
Genetic influences on cytokine responses : We identified multiple genetic loci associated with variation in cytokine responses, explaining 10-50% of the variance depending on the stimulus-cytokine combination. The most significant findings included rs55792153 associated with IL-1β, TNF, and IL-1Ra responses to cGAMP stimulation, and rs2343196 associated with IL-10 responses to dsRNA stimulation.
Functionally-related cytokine modules : Using our novel CytoMod algorithm, we identified 33 distinct modules of co-expressed cytokines across different stimulation conditions. These modules represent coordinated immune responses and may be more biologically meaningful than individual cytokines.
Sex differences in innate immune responses : Boys showed significantly higher responses to R848 stimulation, particularly for MCP1 and IFNα production. These sex differences may reflect developmental or hormonal differences in immune function at this age.
Associations with adiposity : Higher BMI was associated with elevated IL-8 and CXCL10 responses to Poly(I:C) stimulation, potentially reflecting increased inflammatory activation in children with higher adiposity.
Mode of birth influences immune development : Children born via caesarean section showed reduced TNF and IFNα responses compared to those born vaginally, consistent with previous studies suggesting that caesarean delivery results in delayed immune priming.
Systemic inflammation and immune responsiveness : Markers of systemic inflammation (GlycA, hsCRP, granulocyte-to-lymphocyte ratio) were positively associated with innate immune cytokine responses, particularly for pro-inflammatory cytokines. Growth factor responses showed inverse associations with systemic inflammation markers.
Seasonality of immune responses : Community viral infection burden showed significant associations with innate immune response modules, with peak module expression coinciding with peak community infection seasons. This suggests that recent pathogenic exposure influences the immune phenotype measured in vitro.
Interpretation of findings
The identification of distinct cytokine response patterns to different PAMPs reflects the specificity of pattern recognition receptor signalling pathways. TLR4 stimulation (LPS) primarily activates pro-inflammatory responses dominated by TNF, IL-1β, and IL-6. TLR2 stimulation (PGN) activates a somewhat broader response including anti-inflammatory factors. TLR3 stimulation (Poly(I:C)) activates interferon-dominated responses. TLR7/8 stimulation (R848) activates both pro-inflammatory and interferon responses. Cytosolic DNA sensing (cGAMP) and dsRNA sensing activate strong interferon responses. These distinct patterns were captured by our principal component analysis and hierarchical clustering approaches, validating the study design.
The genetic findings suggest that variation in innate immune responses is partially determined by inherited factors. The identification of rs55792153 as significantly associated with IL-1β responses to cGAMP stimulation is particularly interesting, as this SNP is located in a region with multiple genes that could plausibly influence cGAMP signalling or downstream inflammatory responses. The fact that this SNP shows study-wide significance (rather than just genome-wide significance) suggests this is a robust finding.
Our CytoMod approach for identifying cytokine modules represents an advance over previous methods. By incorporating bootstrapping and reliability scoring, we identified robust modules that are less likely to be driven by outliers or technical noise. The identification of distinct modules within stimulation conditions (e.g., 4-5 modules for single stimuli) suggests that even with a single stimulus, multiple distinct immune pathways are activated in different individuals or subsets of cells.
The sex differences in cytokine responses observed in this study are consistent with previous reports of sex-based differences in immune function in children. The higher responses in boys to R848 stimulation (a TLR7/8 agonist) is particularly interesting given that TLR7 and TLR8 are located on the X chromosome. However, we observed higher responses in boys rather than girls, which is unexpected if X-linked gene dosage were the primary driver. This suggests that sex hormone differences or sex-based developmental differences are more likely explanations.
The associations between adiposity and immune responses are consistent with the concept of “meta-inflammation” in obesity, where increased adiposity is associated with low-grade systemic inflammation and altered immune function. The specific associations we observed with IL-8 and CXCL10 responses suggest that chemokine production may be particularly sensitive to obesity-related changes in immune function.
The reduced immune responses in children born via caesarean section are consistent with the “sterile prenatal period” hypothesis, which proposes that the absence of bacterial colonization during vaginal delivery results in delayed immune priming. The fact that children who experienced labour prior to caesarean delivery showed intermediate responses suggests that labour itself (or the associated labour-related processes) contributes to immune priming. This finding has important implications for understanding the increased risk of allergic and autoimmune diseases in caesarean-delivered children.
The associations between systemic inflammation markers and innate immune responses are bidirectional in nature - systemic inflammation may reflect elevated innate immune activation (as suggested by associations between GlycA/hsCRP and pro-inflammatory cytokines), or alterations in immune responsiveness may lead to systemic inflammation. The inverse associations between growth factors and systemic inflammation markers suggest that systemic inflammation may suppress growth factor production, potentially reflecting a shift from tissue-remodelling/angiogenic responses toward pro-inflammatory responses.
The seasonal associations between community infection burden and immune responses suggest that recent pathogenic exposure influences innate immune phenotype. This has important implications for interpretation of immune measurements in natural populations, where seasonal exposures inevitably influence baseline immune status.
Limitations and strengths
This study has several important strengths. The large sample size (n=286) provides excellent statistical power for identifying genetic and environmental associations. The comprehensive characterization of cytokine responses to multiple stimuli provides a nuanced view of innate immune function. The long-term follow-up (cytokine measurements at age ~4 years) and integration of genetic, anthropometric, perinatal, and environmental data provides a holistic view of factors influencing immune development.
However, there are also limitations to consider. The study population is drawn from Australia, which limits generalisability to other populations with different genetic backgrounds and environmental exposures. The blood-based assay measures ex vivo responses, which may not fully recapitulate in vivo immune function. The cross-sectional analysis of cytokine responses at age 4 years does not allow us to determine whether the associations we observe are causally related or reflect reverse causality. Some associations reported in this paper show only nominal p-value significance without correction for multiple testing, which may lead to false positives.
Conclusions
This comprehensive study of innate immune cytokine responses in pre-school children reveals multiple factors influencing immune function, including genetic variation, sex, adiposity, mode of birth, systemic inflammation status, and recent pathogenic exposure. The identification of robust cytokine modules provides a framework for future studies examining immune function in health and disease. Future studies should examine whether the associations observed in this study are causally related to clinically important health outcomes such as infection risk, allergic disease, and autoimmune disease.
Supplementary figure captions
Supplementary Figure 1: all measured cytokine responses
The graphs are ordered by stimuli (columns) and analytes (rows). The depicted values were adjusted for ‘time in freezer’ and ‘exact incubation time’ using a linear regression approach (see Methods section).
Supplementary Figure 2: Genetic variants and cytokine responses
(A) Locuszoom plot (top) and genetrack (bottom) of rs72679554. (B-C) Covariate-adjusted concentrations of respectively IL-10 and TNF upon stimulation with Lyovec/dsRNA, in children with different genotypes on this location. (D) Locuszoom plot (top) and genetrack (bottom) of rs2384950. (E) Covariate-adjusted concentrations of IFNγ upon stimulation with R848, in children with different genotypes on this location. (F) Locuszoom plot (top) and genetrack (bottom) of rs2300875. (G) Covariate-adjusted concentrations of MCP1 upon stimulation with poly(I:C), in children with different genotypes on this location. (H) Locuszoom plot (top) and genetrack (bottom) of rs6112096. (I) Covariate-adjusted concentrations of MCP1 upon stimulation with RPMI alone, in children with different genotypes on this location. (J) Top 10 Reactome pathways that were enriched (nominal p < 0.05), shared between the highest proportion of cytokines for stimulation with each of respectively PGN, R848, poly(I:C), Lyovec/RPMI, and Lyovec/dsRNA. In this figure, covariates included in the models are time in freezer, exact incubation time, sex, seasonal infection burden, and granulocyte percentage. Boxplots are in the style of Tukey.
Supplementary Figure 3: Correlation between individual cytokines and module scores
Each panel represents a cytokine module, organised by stimulus and from left to right. The stimulus and module are indicated in the lower-right subpanel on the diagonals. The Pearson correlation coefficients are indicated by the numbers in top subpanels, while the individual data points (log-transformed and standardised) are plotted in the lower subpanels.
Supplementary Figure 4: Associations between birth variables and cytokine responses
(A) Volcano plot indicating the change in cytokine levels with an increase of 100 gram in birthweight. (B) Volcano plot indicating the change in cytokine levels with an increase of 1 week gestational age. © Volcano plot indicating the change in cytokine levels with an increase of 1 SD in birthweight-by-gestational-age. (D) Dotplot summarising regression results between each variable of interest and cytokine module expression score (reference: vaginal birth, except the comparison between ‘with’ vs ‘no’ labour in Caesarean birth). ‘Effect size’ was calculated as absolute value of the regression estimate. In this figure, ‘covariate adjusted’ includes time in freezer, exact incubation time, and sex, age, and BMI.
Supplementary Figure 5: Associations between systemic inflammation/cell type composition and cytokine responses
(A-C) Distribution of cell type percentage in whole blood for respectively granulocytes, lymphocytes, and monocytes. (D) Comparative plot indicating the difference change in module expression scores when GlycA or hsCRP is increased by 1 SD. The points are coloured by the indicated FDR-adjusted p-value significance categories. (E) Comparative plot indicating the difference change in cytokine levels when GlycA or GLR is increased with 1 SD. The points are coloured by the indicated FDR-adjusted p-value significance categories. For panels D-E, points in the white areas are concordant between measures; points in the grey areas are discordant. (F) Dotplot summarising regression results between each variable of interest and cytokine module expression score. In each dot ‘*’ indicates FDR-adjusted p < 0.05. ‘Effect size’ was calculated as absolute value of the regression estimate. (G-I) Volcano plots indicating the change in individual cytokine levels upon an increase of 1 SD in respectively granulocyte, lymphocyte, or monocyte percentage. Panels G-I have nominal p-values. In this figure, ‘covariate adjusted’ includes time in freezer, exact incubation time, and sex, age, and BMI. Boxplots are in the tradition of Tukey.
Supplementary Figure 6: Seasonal variation in cytokine module expression scores
The graphs are organised by stimulus (rows) and cytokine modules (columns). There are empty spots in case there were fewer than 5 modules for a stimulus. The depicted values were adjusted for ‘time in freezer’, ‘exact incubation time’, sex, age, and BMI using a linear regression approach (see Methods section).
Supplementary Figure 7: Seasonal variation in individual cytokines across stimuli
The graphs are organised by stimulus (columns) and cytokines (rows). The depicted values were adjusted for ‘time in freezer’, ‘exact incubation time’, sex, age, and BMI using a linear regression approach (see Methods section).
Supplementary Figure 8: Seasonal infections in the community
(A) Number of positive viral swabs over the study period. (B) Number of total viral swabs over the study period. © Percentage of viral swabs that was positive over the study period.
References
- 1.[Note: The References section above (entries 1-133) continues from the previous chunk and includes all citations through entry 133. No additional references appear on pages 11-20 of the document.]
- 2.Medzhitov, R. & Janeway, C. A. Innate immunity. N. Engl. J. Med. 343, 338–344 (2000).
- 3.Akira, S., Uematsu, S. & Takeuchi, O. Pathogen recognition and innate immunity. Cell 124, 783–801 (2006).
- 4.Kawai, T. & Akira, S. The role of pattern-recognition receptors in innate immunity: update on Toll-like receptors. Nat. Immunol. 11, 373–384 (2010).
- 5.Takeuchi, O. & Akira, S. Pattern recognition receptors and inflammation. Cell 140, 805–820 (2010).
- 6.Trinchieri, G. & Sher, A. Cooperation of Toll-like receptor signals in innate immune defence. Nat. Rev. Immunol. 7, 179–190 (2007).
- 7.Netea, M. G., Balkwill, F., Chonchol, M., Cominelli, F. & Donath, M. Y. A guiding map for inflammation. Nat. Immunol. 18, 826–831 (2017).
- 8.Hotamisligil, G. S. Inflammation and metabolic disorders. Nature 444, 860–867 (2006).
- 9.Netea, M. G., Latz, E., Mills, K. H. & O’Neill, L. A. Negative regulation of toll-like receptor signaling. J. Immunol. 185, 6436–6442 (2010).
- 10.Netea, M. G., van de Veerdonk, F. L., van der Meer, J. W., Kullberg, B. J. & van Krieken, J. H. Inflammasome-independent regulation of IL-1-family cytokines. Annu. Rev. Immunol. 33, 49–77 (2015).
- 11.Barton, G. M. A calculated response: control of inflammation by the innate immune system. J. Clin. Invest. 118, 413–420 (2008).