INTRODUCTION
An estimated 14.3 million cases of typhoid and paratyphoid fever, collectively referred to as enteric fever, occur each year (
1). Systemic infection is caused by bacterial pathogens
Salmonella enterica serovars Typhi and Paratyphi A-C, acquired from contaminated food and water. While a protein-polysaccharide conjugate vaccine has been recently demonstrated to be highly effective against typhoid fever in children, no vaccine against paratyphoid fever has yet been licensed (
2,
3).
S. Typhi is thought to invade the intestinal epithelium and disseminate systemically within the first 24 hours after infection (
4). Infection is asymptomatic for up to 2 weeks (
5), after which a secondary bacteremia occurs, often accompanied by fever and fatigue (
6).
It is unclear whether vaccination and natural immunity exert their protective effect on bacterial invasion of the intestinal epithelium at a later stage in infection or at multiple points during infection. Furthermore, mechanisms of protection, which may differ between natural immunity and vaccines, are unclear. The enteric fever human challenge model, where volunteers are orally exposed to typhoidal
Salmonella, offers a unique opportunity to investigate host-pathogen interactions during the incubation period. Previous studies using the human challenge model have shown that live
S. Typhi can be detected in stool in the first days after oral challenge, both for participants who go on to develop disease and those who stay healthy (
7), suggesting that
S. Typhi reaches the lower gastrointestinal tract in both groups. Likewise,
S. Typhi DNA has been detected in the blood in the first few days after the challenge, both in those who later develop typhoid and those who do not (
8).
In this study, we apply transcriptomics as an unbiased approach to identifying early host responses to enteric fever, both in unvaccinated and vaccinated participants. Transcriptomics has previously been used to characterize pre-symptomatic host responses to challenge with influenza (
9), rhinovirus (
10), and malaria (
11), providing a broad snapshot of the pathways activated in the first days after exposure. Controlling for the confounding effect of circadian rhythms, we find enrichment of gene sets relating to innate immunity in both those who did and did not go on to develop enteric fever. Furthermore, estimating cell proportions from gene expression data suggested a greater number of monocytes may be in an activated state 12 hours post-challenge relative to baseline. We also find that plasma TNF-α is significantly raised in challenged participants but not unchallenged controls.
DISCUSSION
This is the first study to profile changes in the human blood transcriptome 12 hours after exposure to S. Typhi and S. Paratyphi A. 199 genes were transcriptionally upregulated 12 hours after challenge, suggesting a host response to exposure even at this early stage of infection. There was also an increase in plasma TNF-α and estimated activated monocytes 12 hours post-challenge. No differences in gene expression were identified between those who did and those who did not go on to develop enteric fever, indicating that this early host response occurs even in those who remain healthy despite exposure.
Specific genes upregulated 12 hours post-challenge included caspase
CASP4, which plays a key role in inflammasome-mediated immunity to
Salmonella infections (
22,
23);
IFITM2, which can inhibit intracellular bacterial growth in monocytes (
24); and
RBCK1, a component of a complex regulating NF-kB pathway activation (
25). The kinetics of these changes, returning to baseline by 24 hours post-challenge, are consistent with data from historical challenge studies: participants treated with chloramphenicol 24 hours after challenge still developed blood-culture positive enteric fever, suggesting that by this time
S. Typhi has reached an intracellular environment (
26).
As transcriptional profiling was carried out on whole blood, changes in gene expression were driven not only by transcriptional regulation but also by changes in the proportions of different immune cells. Previously, single-cell RNA-sequencing has been used to characterize cell-specific responses to
ex vivo infection with the related serovar
S. Typhimurium (
21). Exposed monocytes developed a distinct gene expression profile with 90% containing intracellular bacteria. We applied the resultant deconvolution algorithm to our data and found a small rise in monocytes estimated to be in this
Salmonella-exposed activated state. However, as the rise was very small, and cell expression profiles were determined using a different serovar and
ex vivo, this observation needs to be corroborated by future experimental data such as flow cytometry or single cell-RNA-sequencing. Nonetheless, it does support murine data suggesting that infected monocytes are responsible for both systemic
S. Typhi dissemination and early innate responses to infection (
6,
27). Furthermore, human data suggest
S. Typhi interacts with and is found within monocytes at later stages of infection (
14,
28). These results, together with previous research describing a rise in the ability of monocytes to bind
S. Typhi 24 hours after challenge (
29), could imply that monocytes are also involved in early human-
S. Typhi interactions. The rise in the predicted number of activated monocytes was only marked for unvaccinated participants who remained healthy. If confirmed, this could suggest the early response may be protective, and that an alternative antibody-mediated mechanism of protection may take place in vaccinated participants; for example, antibody-dependent phagocytosis and killing (
30) in the lamina propria.
In addition, we observed a rise in TNF-α 12 hours post-challenge in samples from a typhoid toxin human challenge study. A similar but non-significant rise has previously been observed in the typhoid dose-finding study [1.11-fold in the typhoid dose finding study (
P = 0.07) (
14) compared with 1.17-fold here in the typhoid toxin study (
P = 0.01)]. We found that following
in vitro stimulation with
S. Typhi, TNF-α was predominantly produced by monocytes. However, further
ex vivo experiments will be required to ascertain whether TNF-α is primarily expressed by monocytes or other cell types 12 hours post-challenge
in vivo. One possibility is secretion of cytokines from the intestinal mucosa into the blood. However, organotypic models have found TNF-α to be primarily released by monocyte/macrophages rather than lymphocytes, fibroblasts, endothelial cells, or epithelial cells following
S. Typhi or
S. Paratyphi A stimulation (
31,
32). Furthermore,
ex vivo stimulation of intestinal biopsies with
S. Typhi only elicited a non-significant increase in TNF-α (
33). Nonetheless,
in vivo, a number of cell types could be responsible for secretion, including CD8 T cells (
34) and mucosal-associated invariant T cells (
35).
Those who did not go on to develop enteric fever had broadly similar transcriptional responses to those who did in the hours after challenge. This implies that at a dose giving a 60%–75% attack rate, much of the protection against enteric fever may be mediated at a later stage of infection, although it is unclear if the change in the peripheral blood transcriptome reflects a mucosal or systemic host-pathogen interaction. While difficult to study in humans, mouse enteric fever models using
S. Typhimurium as a challenge agent have previously characterized variation in host responses. One possibility is that differences in susceptibility are due to innate genetic variation. In 32 genetically distinct mouse strains orally challenged with
S. Typhimurium, bacteria reached the spleen and liver of all strains but with considerable variation in bacterial load and survival after 7 days (
36). Another study found that whereas three mouse strains had similar bacterial loads in the liver/spleen at day 1 after intravenous infection, there was significant variation by day 3 (
37). Variation in pre-existing immunity may also exert an effect later in infection. For example, vaccinated and naïve mice had similar liver/spleen bacterial loads at day 1 after oral
S. Typhimurium challenge but diverged by day 5 (
38). Conversely, in another study where challenge was intravenously administered at a much lower dose, differences in liver/spleen load were already evident between vaccinated and naïve mice within 6 hours (
39).
Here, we carried out transcriptional profiling of the bulk peripheral blood transcriptome. This allowed us to include changes in granulocyte gene expression and reduced the effect of blood processing on gene expression as samples were immediately lysed by a stabilizing reagent. However, this approach prevented us from finding how early enteric fever challenge affects individual cell populations. Detailed longitudinal and single-cell transcriptional profiling of the incubation period, with participants challenged at varying times of day, would be hugely valuable in elucidating mechanisms of early immunity and disease pathogenesis. This study was also limited to peripheral blood due to the invasiveness of gut, spleen, liver, or bone marrow biopsies. It may be that analyzing host responses in these tissues would be necessary to identify the timepoint at which those who develop enteric fever diverge from those who do not.
Together, our results indicate an early innate response to enteric fever 12 hours following challenge. Monocytes may play a key role in early innate responses to enteric fever. However, early responses were broadly similar between those who remained healthy and those who developed enteric fever, as well as between vaccination groups. A possible explanation is that protection is mediated at a later stage of infection, following systemic dissemination. Further investigation of the mechanism of protection, and its implications for effective vaccine design, is warranted.
MATERIALS AND METHODS
Human challenge studies
Samples were collected from six enteric fever human challenge studies (
Table 1). Whereas, the typhoid- (
15) and paratyphoid-dose finding (
16) studies aimed to find the inoculum dose that would give an attack rate of 60%–75%, the typhoid oral vaccine trial (
7), typhoid parenteral vaccine trial (
17), typhoid toxin study (
18), and re-challenge study (
19) aimed to find the effect of oral vaccination, parenteral vaccination, typhoid toxin knockout, and prior challenge on attack rate, respectively. All studies were approved by the South Central-Oxford A Research Ethics Committee (10 /H0604/53, 11/SC/0302, 14/SC/0004, 14/SC/1427, 16/SC/0358, and 14/SC/1204), and all participants provided written informed consent. In the typhoid oral and parenteral vaccine trials, participants received a typhoid vaccine or placebo/control vaccine 28 days prior to challenge. Participants were challenged with
S. Typhi Quailes strain, toxin-deficient
S. Typhi SB6000 strain, or
S. Paratyphi A NVGH308 strain suspended in sodium bicarbonate. Unchallenged control participants were given sodium bicarbonate alone. Enteric fever was defined as either an oral temperature ≥38°C sustained for ≥12 hours, or a positive blood culture. Participants were treated with ciprofloxacin or azithromycin following diagnosis or at 14 days if they did not develop enteric fever.
Transcriptional profiling
Three milliliters of peripheral blood was collected in Tempus Blood RNA tubes before challenge (approximately 6–9 am) and approximately 12 and 24 hours following challenge. Total RNA was extracted using the Tempus Spin RNA Isolation Kit (Life Technologies).
RNA samples from individuals in the typhoid dose-finding study and typhoid oral vaccine trial were hybridized onto Illumina HT-12v4 bead-arrays at the Wellcome Trust Sanger Institute and Wellcome Trust Centre for Human Genetics, and fluorescent probe intensities captured with GenomeStudio software. Samples were hybridized in three batches. These two studies were grouped together for all following analyses. RNA samples from individuals in the paratyphoid dose-finding study and typhoid parenteral vaccine trial were poly-A selected (but not globin-depleted) and sequenced using a HiSeq V4 at the Wellcome Trust Sanger Institute. Fastq files from the same sample were concatenated, and prealignment quality control carried out using FASTQC. As all files had consistently high phred quality scores (>25) across their length, all were aligned to the human genome (GRCh38 Gencode version 26) using STAR-2.6.1c (
40). Total reads per sample ranged from 16 to 44 million. Reads per gene were counted using the STAR GeneCounts mode.
Publicly available circadian data sets
Two publicly available whole blood transcriptome data sets were downloaded to identify differential gene expression attributable to circadian rhythms: GSE48113 (
12), profiled using an Agilent Whole Human Genome 4 × 44K custom oligonucleotide microarray, and GSM3122578, profiled by RNA-sequencing (
13). Samples from GSM3122578 were aligned to the human genome (GRCh38 Gencode version 26) using STAR-2.6.1c and counted using STAR GeneCounts mode as above. As baseline and 24 hours post-challenge samples were taken at 6–9 am in the typhoid study and 12 hours post-challenge samples at 6–9 pm, data sets were subset to in-phase samples collected at these time points. Two samples from GSM3122578 were excluded on the basis of no mapped reads. Twenty morning samples and thirty-two evening samples were analyzed from GSE48113, and ten morning samples and nine evening samples from GSM3122578.
Data processing
Microarray data sets, including GSE48113, were background corrected and quantile normalized, then the package lumi used to identify outliers for exclusion. On this basis, two samples from one of the typhoid dose-finding/typhoid oral vaccine trial batches were excluded from the raw data, which were then re-normalized. Probes from each array were aligned to the human genome assembly (GRCh38 GenCode version 26) using STAR-2.6.1c (
40). Probes which did not align to genes were excluded, as were probes aligning to multiple genes. Where multiple probes corresponded to one gene symbol, the collapseRows function from the WGCNA R package was used to select the probe with the highest average mean as being representative, on the basis of this being the most reproducible approach (
41). Genes corresponding to probes which were not significantly expressed over the negative controls (
P > 0.05) in more than 37 samples (the number of baseline samples) in the typhoid dose-finding study/typhoid oral vaccine trial data set were excluded.
For the paratyphoid dose-finding study and typhoid parenteral vaccine trial RNA-sequencing data sets, principal component analysis was used for outlier detection, with no samples excluded on this basis. Non-protein coding and hemoglobin subunit genes were excluded. Count tables were filtered to exclude genes with <1 cpm in ≥31 samples (the number of baseline samples in control participants challenged with S. Typhi) and normalized using weighted trimmed mean of M-values scaling (edgeR).
Differential gene expression analysis
For microarray data, a linear regression model was fitted using the limma package. Challenge dose and batch were both incorporated into the model as categorical covariates, where those in the typhoid oral vaccine trial and escalated dose group in the typhoid dose-finding study [10–50 × 103 colony forming units (CFU)] were considered “high dose,” and those in the lower dose group (1–5 × 103 CFU) in the typhoid dose-finding study were considered “low dose.” Analysis was paired by using participant ID as a blocking variable. For RNA-sequencing data, the count matrix was first transformed using limma voom. A linear regression model was fitted with vaccination status, sequence pool, and dose as categorical covariates, where the deescalated dose group (0.5–1 × 103 CFU) in the paratyphoid challenge study was considered “low dose” and the higher dose group (1–5 × 103) as “high dose.” Analysis was paired using participant ID as a blocking variable.
Over-representation of genes relating to circadian rhythms
The results of differential expression analysis using limma were used to classify genes into differentially expressed after challenge (adjusted P value < 0.05) and not differentially expressed after challenge (adjusted P value > 0.05). Genes were also classified as potentially circadian (unadjusted P < 0.05 or absolute log2(fold change) >0.1 in the two circadian data sets) or not circadian. A contingency table, limited to genes present in both the enteric challenge and circadian data sets, was created. A fisher test was then used to test for over-representation.
Transcriptional changes 12 hours post-challenge
To stringently exclude any genes that may be differentially expressed due to delayed processing or circadian rhythms, genes with an unadjusted P value < 0.05 or log2(fold change) >0.1 in either GSE48113 or GSM3122578 were excluded from analysis. The package MetaVolcanoR was then used to combine the differential expression results from the typhoid dose-finding/typhoid oral vaccine trial, typhoid parenteral vaccine trial, and paratyphoid dose-finding study using a random effects model.
The gene list generated by the differential gene expression meta-analysis was ranked by fold change, such that the most upregulated genes were at the top of the list and the most downregulated at the bottom. The gene lists were then input into a gene set enrichment analysis algorithm (
42), in order to determine whether members of different blood transcriptional modules [groups of genes related by biological function or cell-specific expression (
43)], were unevenly distributed in the gene list.
Differential expression between subgroups
Genes that may be differentially expressed due to circadian rhythms were not excluded, on the basis that samples from different groups were subject to the same confounding factors. For different outcome groups, differences at baseline (Remained healthy baseline samples—Developed enteric fever baseline samples) and in response to challenge [(Remained healthy 12 hours post-challenge samples—Remained healthy baseline samples) − (Developed enteric fever 12 hours post-challenge samples – Developed Enteric Fever baseline samples)] were assessed by differential expression in limma. As above, the package MetaVolcanoR was then used to combine the differential expression results from the typhoid dose-finding/typhoid oral vaccine trial, typhoid parenteral vaccine trial, and paratyphoid dose-finding study using a random effects model.
For different dose groups, differences in response to challenge [(Higher dose 12 hours post-challenge samples – Higher dose baseline samples) – (Lower dose 12 hours post-challenge samples – Lower dose baseline samples)] were assessed by differential expression in limma for the typhoid dose-finding/typhoid oral vaccine trial and paratyphoid dose-finding study. Analysis of differences between vaccination groups was restricted to the typhoid parenteral vaccine trial, assessing differences at baseline, in the 12 hours response, and between outcome groups for vaccinated participants only.
Deconvolution of the transcriptional response
RNA-sequencing data sets were normalized by library size to give counts per million, then log-transformed. CIBERSORT was used to estimate changes in immune cell populations
in silico. Estimated cell fractions were compared with differential white blood cell counts taken at baseline in the typhoid dose-finding study and typhoid oral vaccine trial and absolute blood counts in the paratyphoid dose-finding study. Processed gene expression matrices were imported into MATLAB, and the dynamic deconvolution algorithm from
https://github.com/noabossel/Dynamic-deconvolution-algorithm applied to estimate changes in monocytes in an activated infection-induced state. A paired Wilcoxon test was used to assess whether monocytes in an infection-induced state were significantly changed from baseline. Comparing those who developed enteric fever with those who did not at baseline, samples were unpaired with systematic differences in levels of predicted infection-induced state monocytes between studies. Therefore, a logistic regression model with study as a covariate was used. A Mann-Whitney test was used to assess whether there was a difference in the log
2(fold change) in cells 12 hours post-challenge relative to baseline among different subgroups.
Plasma cytokine quantification
Heparinized blood samples were collected and centrifuged at 3,000 rpm for 10 minutes at 4°C. The supernatant was transferred to a fresh tube before centrifuging a second time. The supernatant was then aliquoted and protease inhibitor added in a 1:40 dilution before storage at −80°C. Plasma aliquots at baseline and 12 hours post-challenge or mock challenge were assayed in duplicate using a MILLIPLEX MAP Kit with nine analytes (EGF, TGF-α, CXCL1, CD40L, IL1RA, IL2, IL6, IL8, and TNF-α) following manufacturer instructions. For values with fluorescence intensities falling below that of the lowest standard (3.2 pg/mL), the concentration was set to half that of the lowest standard (1.6 pg/mL). Duplicates were averaged.
In vitro stimulation with S. typhi
To characterize cytokine secretion by different leukocyte subsets, 15 mL blood was collected from each of three healthy donors in a BD ethylenediaminetetraacetic acid vacutainer and layered onto an equal volume of Polymorphprep (Alere Technologies) then centrifuged at 500× g for 35 minutes at room temperature with the brake off. The granulocyte and peripheral blood mononuclear cells (PBMC) bands were transferred to separate Falcon tubes and topped up to 50 mL with phosphate-buffered saline (PBS). Cells were centrifuged at 250× g for 10 minutes then resuspended in 5 mL red blood cell lysis solution (eBioscience) for 5 minutes. The tube was once more topped up to 50 mL with PBS, centrifuged at 250 × g for 10 minutes, and the pellets resuspended in 4 mL PBS for cell counting. CD14 microbeads (Miltenyi Biotec) were diluted one in five in autoMACs running buffer and PBMCs suspended in 300 µL diluted beads. After 15 minutes at 4°C, 6 mL autoMACS running buffer was added, and the PBMCs centrifuged for 10 minutes at 250 G. PBMCs were then resuspended in 500 µL running buffer, and the CD14+ fraction (monocytes) and CD14– fraction (lymphocytes) separated using an autoMACS Pro Separator. Viability was assessed using trypan blue and was high. Monocytes, lymphocytes, and granulocytes were counted using a hemocytometer, suspended in Roswell Park Memorial Institute (RPMI) 1,640 Medium and seeded in duplicate (2 wells per donor per cell type) in a 96 wells plate at a density of 1,00,000 cells per well. Frozen stocks of 107 CFU/ml Quailes strain S. Typhi were thawed and washed twice with RPMI. Cells were inoculated with S. Typhi at a multiplicity of infection of 10 or RPMI as a negative control. After 50 minutes, gentamycin was added to a final concentration of 200 µg/mL. Four hours post-inoculation, the plate was centrifuged at 400× g for 4 minutes, and the supernatants filtered using 0.22 µm cellulose acetate Spin-X centrifuge tube filters. Samples were assayed in duplicate using a MILLIPLEX MAP kit with nine analytes (EGF, TGF-α, CXCL1, CD40L, IL1RA, IL2, IL6, IL8, and TNF-α). For values with fluorescence intensities falling below that of the lowest standard (3.2 pg/mL), the concentration was set to half that of the lowest standard (1.6 pg/mL). Duplicates were averaged.
ACKNOWLEDGMENTS
The authors would like to acknowledge all the volunteers for participating in the studies, the Data Safety and Monitoring Committee for providing safety oversight of the studies, the University of Maryland for providing S. Typhi Quailes challenge strain, and GSK Vaccines Institute for Global Health for providing the S. Paratyphi A NVGH308 strain. We would also like to thank Christina Dold for her work supporting the typhoid challenge studies as Biological Safety Officer and Henderson Zhu for his role in depositing the RNA-sequencing data.
This research was funded in whole, or in part, by the Wellcome Trust [092661/Z/10/Z and 090532/Z/09/Z]. For the purpose of Open Access, the author has applied a CC BY public copyright license to any Author Accepted Manuscript version arising from this submission.
This work was supported by the Wellcome Trust (grant numbers 092661/Z/10/Z awarded to A.J.P., 090532/Z/09/Z), The Bill & Melinda Gates Foundation (OPP1084259; Global Health Vaccine Accelerator Platform grant OPP1113682), European Commission FP7 grant “Advanced Immunization Technologies” (ADITEC), and the NIHR Oxford Biomedical Research Centre.