Multi-omic profiling to assess the effect of iron starvation in Streptococcus pneumoniae TIGR4

Bacterial strains and culture conditions

S. pneumoniae TIGR4 was grown without agitation at 37 °C in air with 5% CO₂ in Todd Hewitt Broth (THB) until mid-exponential phase (OD₆₀₀ = 0.3), and kept at −80 °C with 20% glycerol. Three different biological replicates were made for proteomics and metabolomics experiments, and four replicates for transcriptomics (20 mL, 100 mL, and 200 mL cultures to perform transcriptomics, metabolomics, and proteomics experiments, respectively). Each set of replicates was standardised by inoculum, using starter cultures from glycerol-kept vials, in order to prepare the standard inoculum that was further added to the different biological replicates for each “ome” extraction. Iron-depleted cultures were prepared by adding deferoxamine (DFO) mesylate salt (Sigma) dissolved in water, at 100 µM, according to the dose used for this chelator in other works for the pneumococcus (Trappetti et al., 2011).

Protein extracts

Total cellular proteins and secreted proteins were obtained as described (Mitsuwan et al., 2017). Briefly, to obtain cellular proteins, pellets were washed three times in sterile phosphate buffered saline (PBS) pH 7.4. Bacterial cell wall was digested at 37 °C with top-down agitation by adding 100 U mutanolysin (Sigma–Aldrich, St. Louis, MO, USA). Protoplasts were resuspended in 4 mL of rehydration buffer (7 M urea, 2 M thiourea, 4% CHAPS, 0.5% Triton X-100, 0.005% bromophenol blue, 0.5% Bio-lyte 3–10 ampholytes (Bio-Rad, Hercules, CA, USA)) and disrupted by sonication (6 cycles; 20-s pulses, 90% amplitude). Proteins were recovered by centrifugation (5,000× g, 7 min), dialyzed and concentrated using a centrifugal filter device (Amicon Ultra15, 10 kDa; Millipore, Burlington, MA, USA). To obtain secretome fractions, proteins were precipitated from the supernatants with 10% trichloroacetic acid, after removing cell debris through filter devices (0.22 µm, Millipore, Burlington, MA, USA). Protein pellets were washed twice with 1 mL of ice-cold absolute ethanol. Finally, proteins were air-dried and resuspended in 500 µL of rehydration buffer.

Two-dimensional polyacrylamide gel electrophoresis and image analysis

Protein samples were cleaned using the 2-D clean up kit (Life Sciences, Marlborough, MA, USA) according to manufacturer’s instructions. Proteins were resuspended in 200 µL of rehydration buffer and quantified by the Bradford method (Bradford, 1976). Five hundred µg of protein were subjected to isoelectric focusing (IEF) on 18 cm Immobiline DryStrips immobilized pH gradient (IPG) gel strips (4–7 pH linear gradient (Life Sciences, Marlborough, MA, USA)). The strips were loaded onto a Bio-Rad Protean IEF Cell system (Bio Rad, Hercules, CA, USA), and IEF was performed at 20 °C using the following conditions: 2 h of passive rehydration, 50 V for 10 h followed by a voltage-ramp (250 V for 15 min; 500 V for 30 min; 1,000 V for 1 h; 2,000 V for 1 h; 5,000 V for 1 h; 8,000 V for 2 h); finally, proteins were focused on 70,000 Vh. Before the second dimension, the IPG strips were first soaked for 15 min in equilibration solution (50 mM Tris- HCl buffer, pH 8.8, 6 M urea, 30% v/v glycerol, 2% SDS, and bromophenol blue traces) containing 2.5 mg/mL DTT, and subsequently soaked for 15 min in equilibration solution containing 45 mg/mL iodoacetamide. The second dimension was performed on 12% polyacrylamide gels, using the Protean plus Dodeca Cell system (Bio-Rad, Hercules, CA, USA). Gels were run at 90 V until the dye reached the bottom. Then, gels were stained with brilliant blue G-colloidal solution (Sigma-Aldrich, St. Louis, MO, USA) according to manufacturer’s instructions. Gels were scanned with a GS-800 densitometer (Bio-Rad, Hercules, CA, USA). Digitized images were analyzed with PD-Quest v8.1.0 (Bio-Rad, Hercules, CA, USA). Two analytical gels were made per sample (i.e., biological replicate and protein extraction). Consistent spots were considered as those whose presence remained constant at the three biological replicates. Spots showing a consistent change in the intensity value of at least twofold were included in the quantitative analysis.

Protein identification by MALDI-TOF/TOF MS/MS

Spots excision, protein digestion, peptide desalting and mass spectrometry analysis were performed as already described by our group (Mitsuwan et al., 2017) with slight modifications: after mass spectra acquisition using a MALDI-TOF/TOF (4800 Proteomics Analyzer, Applied Biosystems, Foster City, CA, USA) mass spectrometer in the m/z range 800 to 4,000, Mascot 2.0 search engine (Matrix Science Ltd., London, UK) was used for protein identification running on GPS Explorer™ software v3.5 (Applied Biosystems, Foster City, CA, USA) over the National Center for Biotechnology Information (NCBI) protein database (updated monthly). Search setting allowed one missed cleavage with the selected trypsin enzyme, Streptococcus pneumoniae for taxonomy restrictions, cysteine carbamidomethylation as a fixed modification, methionine oxidation as a variable modification, a MS/MS fragment tolerance of 0.2 Da, and a precursor mass tolerance of 10 ppm. Identifications with a Mascot score >70 (p-value <0.05) were considered as significant.

Inductively coupled plasma-mass spectrometry (ICP-MS) analysis

Iron concentration in THB medium was determined by ICP-MS. Five ml aliquots of samples were chemically digested with 2 mL concentrated nitric acid on a hot plate with a heating ramp of 20 °C until reaching 130 °C, maintaining this temperature for 2 hours. High purity deionized water was added to the digested samples to a final volume of 25 mL. Rh-solution was added as internal standard (final concentration of 10 µg/L). The isotope ⁵⁶Fe and the internal standard were analyzed with a NexION 350X instrument (PerkinElmer, Waltham, MA, USA), equipped with a PFA concentric microFlow nebulizer and a cyclonic PFA spray chamber, and operated at 1,600 W in He collision mode. Results were expressed as µg of Fe per liter of culture.

RNA isolation

Cells resuspended in 1 mL of Tri-Reagent (Sigma–Aldrich, St. Louis, MO, USA) were disrupted by vortexing (20 min) with 0.5 g of glass beads (Sigma–Aldrich, St. Louis, MO, USA). After recovering the supernatant by centrifugation (1 min, 12,000× g, 4 °C), 200 µL of chloroform were added. Samples were centrifuged again (12,000× g for 15 min; 4 °C) and 500 µL of ice-cold isopropanol were added. After 15 min incubation (4 °C), samples were centrifuged (30 min, 12,000× g, 4 °C) and washed with 500 µL of 70% ice-cold ethanol. RNA was air-dried and resuspended in 40 µL of distilled water previously treated with 1% of diethylpyrocarbonate (DEPC). Samples were treated with DNAse (Ambion, Austin, TX, USA) according to manufacturer’s instructions.

RNA amplification, labeling and hybridization to DNA microarrays

RNA quality was assessed using a TapeStation (Agilent Technologies, Santa Clara, CA, USA). The RNA integrity number (RIN) ranged between 7.0 and 9.2. Samples with RIN >7 were considered for analysis. RNA concentration and dye incorporation were measured using a UV–VIS spectrophotometer (Nanodrop 1000, Agilent Technologies, Santa Clara, CA, USA). Hybridization to custom 8 × 15 K Gene Expression Microarrays (ID 044371, Agilent Technologies, Santa Clara, CA, USA) containing the whole genome of S. pneumoniae TIGR4 was conducted following manufacturer’s protocol using a two-color (Cy3 and Cy5) Microarray-Based Gene Expression Analysis (v. 6.5, Agilent Technologies, Santa Clara, CA, USA). Microarrays were then washed and scanned using a DNA Microarray Scanner (Model G2505C).

Gene expression analysis

Microarray hybridization data were obtained with the Feature Extraction Software v. 10.7 (Agilent Technologies, Santa Clara, CA, USA), using the default variables. Data analysis was performed using the R Limma Bioconductor package (Ritchie et al., 2015), according to a direct two color design. A total of eight microarrays were done, corresponding to four biological replicates for each condition using swapped and random mixtures to cope with batch effects. Functional annotation of the differentially expressed genes was done using embedded BlastX included into the Blast2GO program (Conesa et al., 2005) using the public NCBI nr database. Raw feature intensities were corrected using the Normexp background correction algorithm. An initial within-array normalization was done using spatial and intensity-dependent Loess method, followed by a between-array Aquantile normalization. Normalized data are shown in Figs. S1 and S2. Differential expression was ordered according to their adjusted p-values, and the expression of each gene is reported as the log₂ ratio of the value obtained for each condition compared to control condition. A gene was considered differentially expressed if it displayed an adjusted p-value less than 0.05 by the Student t-test. Finally, over- and under- expressed genes were analyzed in terms of gene ontology by using a hypergeometric analysis (GOStats package). Prediction of operons was obtained from the DOOR database (Dam et al., 2007; Mao et al., 2009).

Preparation of intracellular metabolite samples

Cell pellets were washed twice with PBS and resuspended in lysis buffer (PBS and 30% sucrose) containing 100 U mutanolysin. Samples were incubated at 37 °C overnight. Metabolic quenching was achieved by adding ice-cold 50% methanol. Cells were disrupted by sonication (6 cycles: 20 s, 90% amplitude). After centrifugation (5,000× g, 7 min), cell debris was separated through 0.22 µm membrane filters (Millipore, Burlington, MA, USA). Then, supernatants were collected and ultracentrifuged (100,000× g; 1.5 h, 4 °C). Finally, metabolite samples were kept at −80 °C before analysis.

Analysis of intracellular metabolites

Intracellular metabolite samples were analyzed by LC–QTOF MS/MS using an Agilent 1200 Series LC system coupled to an Agilent 6540 UHD Accurate-Mass QTOF hybrid mass spectrometer equipped with dual electrospray (ESI) source as described (Mitsuwan et al., 2017) without modifications. Briefly, chromatography was performed using a C18 reverse-phase analytical column (Mediterranean, 50 mm × 0.46 mm i.d., 3 µm particle size; Teknokroma, Barcelona, Spain), thermostated at 25 °C. The mobile phases were 5% ACN (phase A) and 95% ACN (phase B) both with 0.1% formic acid as ionization agent. The LC pump was programmed with a flow rate of 0.8 mL/min with the following elution gradient: 3% phase B was kept as initial mobile phase constant from min 0 to 1; from 3 to 100% of phase B from min 1 to 13. A post-time of 5 min was set to equilibrate the initial conditions for the next analysis. The injection volume was 3 µL and the injector needle was washed for 10 times between injections with 80% methanol. The parameters of the electrospray ionization source, operating in negative and positive ionization mode, were as follows: the capillary and fragmentor voltage were set at ±3.5 kV and 175 V, respectively; N₂ in the nebulizer was flowed at 40 psi; the flow rate and temperature of the N₂ as drying gas were 8 L/min and 350 °C, respectively. MS and MS/MS data were collected in both polarities using the centroid mode at a rate of 2.6 spectra per second in the extended dynamic range mode (2 GHz). Accurate mass spectra in auto MS/MS mode were acquired in bot MS and MS/MS m/z ranges of 60–1,100 Da. To assure the desired mass accuracy of recorded ions, continuous internal calibration was performed during analyses by using the signals at m/z 121.0509 (protonated purine) and m/z 922.0098 (protonated hexakis (1H, 1H, 3H-tetrafluoropropoxy) phosphazine or HP-921) in positive ion mode; while in negative ion mode ions with m/z 119.0362 (proton abstracted purine) and m/z 966.0007 (formate adduct of HP-921) were used. The auto MS/MS mode was configured with two maximum precursors per cycle and an exclusion window of 0.25 min after two consecutive selections of the same precursor. The collision energy selected was 20 V.

Metabolomics data processing and identification

MassHunter Workstation software (version 5.00 Qualitative Analysis, Agilent Technologies, Santa Clara, CA, USA) was used to process all data obtained by LC–QTOF in auto MS/MS mode. Treatment of raw data files was initiated by extraction of potential molecular features (MFs) with the suited algorithm included in the software. For this purpose, the extraction algorithm considered all ions exceeding 1,000 counts with a single charge state. Additionally, the isotopic distribution to consider a molecular feature as valid should be defined by two or more ions (with a peak spacing tolerance of 0.0025 m/ z, plus 10.0 ppm in mass accuracy). Adduct formation in the positive (+Na) and negative (+HCOO) modes, as well as neutral loss by dehydration were included to identify features corresponding to the same potential metabolite. Then, MFs, characterized by their retention time, intensity in the apex of chromatographic peak and accurate mass, were exported in compound exchange format files (.cef files) and imported into Mass Profiler Professional (MPP) software package (version 12.1, Agilent Technologies, Santa Clara, CA, USA) for alignment and further processing. MPP allowed statistical analysis by Volcano plot (a combination of analysis of variance and fold change analysis). The MS/MS METLIN Personal Compound and Database Library (PCDL) was used to identify significant compounds using both MS and MS/MS information to assure metabolite identification. To identify compounds with no MS/MS information in the METLIN Database, MetaCyc database, and the MassBank database were used.

Statistical analysis

All the quantitative analyses (transcriptomics, proteomics, and metabolomics) were performed from three or four independent biological replicates, and the results are expressed as the mean ± standard deviation. Paired data were analyzed by univariate analysis using the Student’s t-test. Principal component analysis (PCA) was done with the web-based software NIA array analysis tool (http://lgsun.grc.nia.nih.gov/anova/index.html) (Sharov, Dudekula & Ko, 2005). p-Values lower than 0.05 were considered statistically significant.

Effect of iron starvation on pneumococcal growth in vitro

S. pneumoniae TIGR4 was grown in THB, in which the iron concentration available for the bacteria, measured by ICP-MS, was 734 µg/L. In order to perform the further transcriptomic, proteomic and metabolomic analyses, we firstly monitored the growth of this strain in the presence of 100 µM DFO. As observed for other microorganisms, iron depletion caused a slight lag phase in the growth of pneumococcus (control growth rate: 1.93 h⁻¹; DFO-treated culture growth rate: 1.39 h⁻¹), and instead lead to a decrease in the OD at the plateau phase (Fig. 1). We chose the mid-exponential phase for sampling RNA, proteins and metabolites in curve points such that almost no delays in the time points of collecting cultures took place. These corresponded to OD = 0.3.

Transcriptomics

Limma analysis of the signals obtained after hybridizing the RNAs with the arrays revealed a number of differentially expressed genes when comparing the DFO treatment with its control (Supplemental Information 1): 338 genes changed significantly after DFO exposure, of which 118 genes increased their expression and 220 were down-regulated. Figure S3 shows the differentially expressed genes as a Volcano plot.

The log₂ fold change values, resulting from averaging four independent biological replicates for each condition, ranged between 1.78 and −2.13. We noticed a consistent differential expression in which all of the genes pertaining to a same operon (according to the DOOR annotation database) changed in the same trend and even experienced almost identical log₂ fold change values. We grouped the genes into operons according to the DOOR database (Supplemental Information 1). In total, 29 operons were clearly up-regulated, including those coding for the iron-compound ABC transporter system (genes SP_1869, SP_1870 and SP_1871), the closely neighbour gene SP_1872 coding for the iron-compound ABC transporter, the manganese ABC transporter system (genes SP_1648, SP_1649 and SP_1650), the ATP synthase complex (genes SP_1509 to SP_1513), and others including several transporter systems. Also, an operon containing genes for the riboflavin biosynthetic system (genes SP_0175 to SP_0178) was clearly up-regulated. Forty-four operons were clearly down-regulated, including the branched-chain amino acid transporter system (genes SP_0750 to SP_0753) and the large operon of 13 genes (SP_0419 to SP_0431) involved in fatty acid biosynthesis, as well as two operons involved in competence (genes SP_0042 and SP_0043, and genes SP_2235 and SP_2236). In all the cases, the log₂ FC values were almost identical for all the genes of the same operon. Table 1 shows the main operons that were differentially regulated.

The differentially expressed genes were annotated with Gene Ontology terms, using Blast2GO (Fig. 2). According to the molecular function level, the most enriched term among the over-expressed genes corresponded to “ion binding” (GO:0043167), representing around 25% of the genes. At the biological process level, the most enriched term was “cell metabolic process” (GO:0044237). For down-expressed genes, the most prevalent term at the biological process level was the same as for over-expressed genes, i.e., “cell metabolic process”, and at the molecular function level, three different terms, including also that of “ion binding” (GO:0043167), represented around 13% each.

We validated the results obtained in the microarrays with RT-q-PCR on a subset of 16 differentially expressed genes of the TIGR4 strain (Table S1). For all the genes, there was the same trend (either over- or down-expression) and a high correlation (r² = 0.93) in the fold changes measured both at the microarray and the RT-q-PCR, thus confirming that the microarray data were reliable.

Proteomics

We compared the proteomes of two different protein fractions, cell extracts and secreted proteins, by 2-DE under iron deprivation conditions (Fig. S4). We then analyzed both the qualitative (i.e., absolute appearance or disappearance of a protein in a condition compared to the other one) and quantitative changes (i.e., changes in protein abundance on spots present in both conditions) when comparing for each protein fraction the iron starvation condition with its non-iron deprived control (Tables 2 and 3).

For the two analyzed protein fractions, most of the changes corresponded to predicted cytoplasmic proteins, mainly affecting enzymes of the primary metabolism. Glyceraldehyde-3-phosphate dehydrogenase was less abundant under iron deprivation. This enzyme was strongly reduced in total extracts after DFO treatment (fold-change decrease > 4), as well as not detected in the secretome fraction of DFO-exposed culture (but detected in the control secretome). Another enzyme of the glycolysis pathway that also decreased its abundance was the phosphoglycerate kinase. Other two Mg²⁺-dependent enzymes of the glycolysis pathway, enolase and pyruvate kinase, also decreased in secretomes of DFO-exposed cells. We also found that one enzyme of the fatty acid biosynthesis pathway was more abundant under iron deprivation: 3-oxoacyl-[acyl-carrier protein] reductase in total extracts (fold change around 7). Also, two major and very abundant extracellular proteins were detected to decrease in DFO-treated cultures: the choline-binding protein A (SP_2190) and the PcsB protein (SP_2216).

We also found a decrease in the abundance levels of an enzyme taking part in the recycling pathway of amino sugar compounds: the divalent cation-dependent N-acetylglucosamine-6-phosphate deacetylase, SP_2056 (NagA), which was found to be almost ninefold less abundant in DFO-exposed secretomes.

Other cation-containing or dependent enzymes decreased their abundances under iron deprivation conditions. Thus, lower levels of zinc-containing alcohol dehydrogenase and manganese-dependent inorganic pyrophosphatase were detected, as well as of hypoxanthine-guanine phosphoribosyltransferase, which requires Mg ²⁺ as a cofactor.

Metabolomics

Finally, we studied the changes in the metabolic profile of intracellular metabolite fractions from both pneumococcal strains, by using LC–MS/MS, which allowed detection of 719 and 826 different chromatographic peaks, in negative and positive ionization mode, respectively. Statistical analysis by Volcano plot revealed that 64 entities presented a p-value below 0.05 and a fold change, in terms of relative concentrations, higher than 2 for discrimination between treated and non-treated TIGR4 samples (data not shown).

A multivariate statistical analysis was carried out with significant entities to evaluate whether the iron deprivation treatment had an effect on the metabolite profile. Figure S5 shows the principal component analysis (PCA) of identified metabolites when comparing iron-depleted and non-depleted cultures, revealing that the first principal component (X-axis) clearly grouped the three control biological replicates, which were clearly separated from two out of the three DFO-treated samples. However, there was dispersion in the three biological replicates of the DFO treatment, as the first principal component did not group them.

The tentative identification of significantly changing entities led to a panel of 17 compounds (Table 4). Although the number of changing metabolites positively identified was low, we clearly found an increase in the concentration of intermediate amino sugar metabolites involved in the cell wall peptidoglycan metabolism: there was an increase in uridine-5′-diphosphoglucuronic acid (2.2-fold), UDP-N-acetylmuraminate (>20-fold), N-acetylglucosamine (3.8-fold), and UDP-N-acetylglucosamine (2.4-fold).