Single-cell analysis of Mycobacterium tuberculosis with diverse drug resistance using surface-enhanced Raman spectroscopy (SERS)

Microbial preparation

For this study, we used mycobacteria from the library of Mtb clinical strains collected in the Saint-Petersburg State Research Institute of Phthisiopulmonology. Clinical strains were obtained by culturing sputum samples from patients earlier hospitalized at the Institute (see Table 1). All patients had signed a written consent to store anonymized processed clinical material. This material’s usage for the study has been approved by the Independent Ethics Committee at the Saint-Petersburg State Research Institute of Phthisiopulmonology. No human participants were involved in this research, which operated with the bacterial strains kept in the Institute’s library.

The samples were stored in the library in the frozen state at −80 °C as a dense suspension in a physiological solution with 15% glycerol. The details of procedures of sample collecting and preparing for storage can be found in Zyubin et al. (2019) as well as procedures applied to assign them a drug resistance status.

The strains designated for subsequent culturing, deactivation, and examination via Raman spectroscopy methods were retrieved from storage and made available as of September 15, 2021. No additional human participants were involved in the present study. Among the bacterial samples extracted from the library (anonymous respectively to the patients’ names), we used material of six strains belonging to the following types: (A) strains isolated from the respiratory material of three patients suffering from pulmonary tuberculosis and (B) strains isolated from the site of the bone destruction after the surgical treatment of patients suffering from the bone tuberculosis. The drug sensitivity to particular kinds of medications of the chosen samples was assigned to them in the library’s description based on the standard method of the growth quantification on the Löwenstein–Jensen medium and/or by the BACTEC MGIIT 960. Table 1 reports this property of the explored strains; their denotation as multidrug-resistant (MDR) and extensive drug-resistant (XDR) follow the WOS classification.

Following the protocol described in detail earlier (Zyubin et al., 2019), the samples were unfrozen suspended in $250\mu \mathrm{l}$ of distilled water, killed by heating at $80$ °C during $20\mathrm{m}\mathrm{i}\mathrm{n}$ and re-suspended in 100–150 $\mu \mathrm{l}$ of distilled water after centrifugation during 20 min at 2,000 rpm. A drop containing diluted bacterial suspension was placed on the SERS substrate’s surface using an automated pipette in an amount of 5 $\mu \mathrm{l}$. The bacterial content was adjusted in such a way that after drying one can clearly identify one separate bacterial cell in the field of the 100× objective of an optical microscope.

SERS substrates

Gold SERS substrates were bought from Silmeco (Denmark). They consist of gold-coated freestanding silicon nanopillars and adopted for the $\lambda =785$ nm excitation. Under this condition, these substrates allow the enhancement of the Raman signal up to ${10}^4$ times.

SERS instrumentation

SERS-based records were performed with Raman Senterra (Bruker, Billerica, MA, USA) spectrometer under the following operational conditions: DPSS laser excitation corresponds to 785 nm and generates the power of 50 mW on the sample; the range of the spot size varies ranged from $1\times 3$ to $1\times 5\mu \mathrm{m}$. Such a spot provides an opportunity for accurate positioning in different points of a single bacterial cell’s surface, see Figs. 1B and 2 using a confocal microscope with a 100× (NA 0.9) objective; the diffraction grating corresponds to $1\text{,}200\mathrm{g}\mathrm{r}/\mathrm{m}{\mathrm{m}}^{-1}$; Rayleigh scattering was eliminated by the notch filters.

The optical scheme was supplied with a CCD camera IDUS 416 (Andor, London, UK) with digital picture resolutions of $1\text{,}024\times 256$ pixels. The preliminary calibration was carried out using a silicon substrate at a static spectrum centred at $520.1\mathrm{c}{\mathrm{m}}^{-1}$ for 1 s, and the sample control measurements without SERS were carried out for bacterial spectra using chemically purified quartz glass as the substrate. Due to the enhancement of the SERS structure, the surface plasmon field-enhanced Raman spectrum of the bacterial cell was registered from points of a single cell using the Raman Senterra spectrometer.

For each biological replicate per strain, the spectrum was recorded with 70 s acquisition time; the final spectra (one for each position on the cell’s surface) are represented by averaging over three registered measurements. Figure 1 illustrates such raw (without extracting auto-fluorescence-based baselines) spectral curves in correspondence to the points of a bacterium, where they were recorded. Such raw spectra were preprocessed before analysis, see below. The total number of measurements contains 172 spectra. Among this set, 60 of them exhibit the best signal-to-noise ratio; they were chosen for the subsequent analysis of the six different strains of Mtb.

Spectral preprocessing

SERS spectra were recorded within the range of 500–1,800 $\mathrm{c}{\mathrm{m}}^{-1}$, which belongs to the so-called “fingerprint” region. All spectra were normalized on the maximum intensity. The chosen raw spectra were preprocessed by following the standard route with a homemade MATLAB code. The first step consisted of the removal of spike-like artefacts that emerged usually due to cosmic rays. For this goal, the algorithm proposed in the work by Whitaker & Hayes (2018) was applied. It reveals the range of sampled wavenumbers, where the indicator, proposed in the cited article exceeds some prescribed limit. Here we used the level, i.e., the tenfold overcoming the median value over the whole sample. In contrast to the approach reported by Whitaker & Hayes (2018), we replaced the raw values within such intervals, not by those obtained with a weighted moving averaging but weighted moving median filter. This replacement is substantiated by the character of noise in the studied biological samples. The length of the median filter was chosen as a minimal odd length, which does not lead to local zero weights. The sample’s values, which were not classified as artificial spikes, were kept in the output sample the same as in the input to prevent double data smoothing. The second step consisted of this smoothing carried out with the cubic Savitzky-Golay filter with the window length equal to 21 (when another value is not specified explicitly). This relatively large width was determined as preferred for two reasons: (i) reach high-frequency content of the noise and (ii) the goal to reduce the overall number of peaks for analysis to a relatively few numbers of principal ones, which are repeatable over all spatial locations and have biophysical meaning (see the discussion below). Finally, the baseline correction using the MATLAB function made available by Al-Rumaithi (2023), which is principally based on the algorithm proposed in the article developed in the work by Schulze et al. (2012).

Individual cells

During the experiment, we recorded spectra from sensitive (SENS), multidrug-resistant (MDR), and extensive drug-resistant (XDR) bacterial cells at 100× magnification for pulmonary and extrapulmonary antibacterial strains. The applied experimental configuration allowed recording at several localized points (from three to six locations) of an individual bacterium as shown in Fig. 2. There are the following cases of objects with different drug-resistant status there: SENS (A, B); MDR (C, D); XDR (E, F). Figure 3 depicts the spectra obtained from locations denoted in Fig. 2 by arrows.

We would like to highlight two features of this applied approach, which make available the analysis of spectral features of interest: (i) the intensity of optical signal allows for catching the response not only from the cell wall but also from intracellular metabolites, e.g., the nucleic acid content, and so on; (ii) collecting signals from closely spaced locations along the bacterial cell enables the identification of points where enhanced spectra are most clearly expressed.

For this purpose, the total range of wavenumbers 400–1,800 $\mathrm{c}{\mathrm{m}}^{-1}$ corresponding to the expected fingerprint region was subdivided for the subsequent analysis into parts within three spectral intervals: 400–1,000 $\mathrm{c}{\mathrm{m}}^{-1}$, 1,000–1,400 $\mathrm{c}{\mathrm{m}}^{-1}$, and 1,400–1,800 $\mathrm{c}{\mathrm{m}}^{-1}$. Figure 3 denotes this subdivision.

Based on our previous articles (Zyubin et al., 2019, 2021), which operated with macroscopic, respectively to a cell’s size and scale, we are focused on the regions denoted as 1 and 2 in Fig. 3 and do not consider the third region, since the latter contains mainly spectral bands associated with identifiers of proteins and amide groups, not evidently involved in the biochemistry of drug resistance.

Among the regions of interest, the most important is the spectral range (720–800) $\mathrm{c}{\mathrm{m}}^{-1}$, which contains bands indicating DNA methylation signatures, as shown in Fig. 4A. It has been observed earlier, e.g., see the arguments provided by Papaleo et al. (2022), that DNA methylation plays a role in the formation of bacterial antibiotic resistance. Specifically, the intensity of the spectral bands at 747 and $775\mathrm{c}{\mathrm{m}}^{-1}$ can be corresponded to thymine base bands in probes of methylated DNA (Kim et al., 2017). In comparison with the results reported by Kim et al. (2017), the second thymine peak ( $747\mathrm{c}{\mathrm{m}}^{-1}$) is shifted to $775\mathrm{c}{\mathrm{m}}^{-1}$ for XDR strains. Regarding the peaks associated with the methylation process, it is worth noting that a peak at $722\mathrm{c}{\mathrm{m}}^{-1}$ corresponds to the intermolecular vibrations of DNA/RNA, and its intensity is higher in XDR compared to MDR and SENS strains (Fig. 4A), which aligns with previous findings by Kim et al. (2017). It is important to highlight that the peaks in the spectral range of $742-747\mathrm{c}{\mathrm{m}}^{-1}$ observed in MDR and XDR strains may be attributed to thymine (Kim et al., 2017), but were not detected in SENS strains.

In the case of extrapulmonary samples (Fig. 5A), it is likely that a peak corresponding to adenine (Stefancu et al., 2022) in MDR and XDR strains can be identified at a frequency of $733\mathrm{c}{\mathrm{m}}^{-1}$, which slightly shifted relatively to the peak value ( $730\mathrm{c}{\mathrm{m}}^{-1}$) obtained in the work by Stefancu et al. (2022). Concerned SENS strains it could be potentially shifted to $720\mathrm{c}{\mathrm{m}}^{-1}$ or closer to $742-740\mathrm{c}{\mathrm{m}}^{-1}$ (see Fig. 5A).

The second narrow range $980-1\text{,}020\mathrm{c}{\mathrm{m}}^{-1}$, presented in the Figs. 4B and 5B, is also of interest in terms of both DNA methylation markers (Stefancu et al., 2022) in resistance and phenylalanine, which is an amino acid that plays an important role in pathogen metabolism, see the work by Griffin et al. (2012). Additionally, phenylalanine is also involved in the regulation of gene expression and virulence factors in M. tuberculosis (Gouzy, Poquet & Neyrolles, 2014). Since the peaks of phenylalanine located at $1\text{,}002\mathrm{c}{\mathrm{m}}^{-1}$ (Zyubin et al., 2019) and 5-methylcysteine at $1\text{,}005\mathrm{c}{\mathrm{m}}^{-1}$ (Stefancu et al., 2022) are located quite close together, as seen in Fig. 4B, they can overlap for XDR strains. Interestingly, for strains with multiple drug resistance as well as for sensitive strains, these peaks are lower in intensity and significantly shifted to the left. As for extrapulmonary strains, as shown in Fig. 5B, the peak characteristics of phenylalanine are identical for all extrapulmonary strains.

Another demonstrable spectral feature is the Raman spectrum line situated around the frequency $1\text{,}170\mathrm{c}{\mathrm{m}}^{-1}$. It is clearly seen in Fig. 4C, the intensity of this band centred at $1\text{,}169\mathrm{c}{\mathrm{m}}^{-1}$ is close to zero for the sensitive strain, but noticeable for strains with antibiotic resistance.

When considering the SERS for whole cells (Zyubin et al., 2021), this line has been associated with glutathione (GSH) as corresponding to one of the principal indicator lines (Qian & Krimm, 1994) of this biomolecule involved in the host immune response to M. tuberculosis (Allen et al., 2015). Although mycobacteria do not synthesize glutathione themselves, its presence in clinical strains obtained from patients can be of exogenous nature since M. tuberculosis is, in principle, able to import GSH to its cytosol, see the respective discussion in the work by Morris et al. (2013). In addition, it has been mentioned there that GSH is toxic for mycobacteria but one may expect that XDR strains can be safely at its higher concentration. The SERS spectral registration implemented in the present work is able to detect methylation of DNA, and this argues that other intracellular content can affect the Raman spectrum too.

It is also worth noting that M. tuberculosis can synthesize its own functional counterpart of GSH, mycothiol (MSH) protecting the bacterium from a stress action including that of a drug-induced one. Unfortunately, there are no available Raman spectra of this compound at the moment. But the SERS study by Picot et al. (2019) of N-acetyl-l-cysteine, which is one of mycothiol’s components, shows the existence of the band of $\delta {\mathrm{C}}_{\alpha }\mathrm{H}$ mode at the frequency $1\text{,}179\mathrm{c}{\mathrm{m}}^{-1}$ positioned not so far from the discussed one of GSH’s. Thus, one should not discount the possibility of either alternative or combined origin of the drug resistance-indicating behaviour of the considered spectral band. However, any variant leads to the same visible manifestation demanded from the diagnostic point of view.

Ensemble of individual cells

To explore whether the discussed features are general for spectra recorded from bacteria with different drug-resistance statuses, we carried out an additional analysis based on an ensemble of bacteria and TB types. This ensemble includes 77 spectra, among which 28 belong to the SENS class, 31 to the MDR class, and 23 to the XDR class; 36 and 41 spectra were taken from bacteria, which cause pulmonary and extrapulmonary TB types, respectively. For the full set of spectra and their properties, see the section ‘Data Availability’.

As an additional test of the opportunity to use Raman spectra for solving the classification problem, here we use an alternative method of data processing. This approach relies on peak prominence, recently proposed as an effective method for handling spectra with complex backgrounds. It minimizes ambiguity that can arise from varying baseline removal techniques, see e.g., the work by Coca-Lopez (2024). Simultaneously, it is more reliable for searching positions of spectral peaks than their magnitude. Therefore, we applied this method implemented as a standard MATLAB function, determined a set of lines for each spectrum, and excluded as outliers the ones, which have prominence either less than half of the median prominence or a hundred times larger of it.

As the first issue to be investigated, we explored the possibility of distinguishing between different drug-resistance statuses using Raman spectral signatures. As such a signature, we used the number of peaks localized within 10 bins of $100\mathrm{c}{\mathrm{m}}^{-1}$ width centred from $450$ to $1\text{,}450\mathrm{c}{\mathrm{m}}^{-1}$. The agglomerative hierarchical cluster tree analysis of these signatures was carried out using the ward linkage method, i.e., the minimum-variance algorithm operating with Euclidean distances, was applied. The resulting dendrogram giving the plot of the hierarchical binary cluster tree is shown in Fig. 6. One can see that indistinguishable clusters primarily comprise bacteria of a single type. Typically, the genomic variation between bacteria with different drug resistance status is smaller than the variation observed between distinct types of tuberculosis (TB). However, within a given TB type, significant differences may exist between clusters representing distinct drug-resistance statuses. Thus, we show that even this simplified line-counting approach could be used to solve the classification problem. This conclusion also confirms the idea of the significance of the Raman spectrum content preliminary discussed in the work by Lavrova et al. (2023) by the case study of a very limited number of single spectra.

As the second test, we explored the values of three particular Raman spectral lines discussed above. For this goal, we test whether there exists a detected maximum within the range of $10\mathrm{c}{\mathrm{m}}^{-1}$ centred in each of three frequencies considered as markers of drug resistivity. Their mean occurrences averaged over the ensemble are shown in Fig. 7. They indicate that there is a certain clearly expressed correlation between the detection of a maximum in all these ranges and the drug resistivity. Note that the applied method does not measure the peak’s amplitude, i.e., weakly expressed peaks make the same contribution in the resulting number of occurrences. This fact, along with potential shifts of adjacent peaks into the analyzed frequency range, may account for the presence of a maximum observed in approximately half of the drug-sensitive bacteria. However, growing drug resistance leads to more frequent detections. For the line $730\mathrm{c}{\mathrm{m}}^{-1}$, this dependence looks almost linear. The lines $747\mathrm{c}{\mathrm{m}}^{-1}$ and $1\text{,}170\mathrm{c}{\mathrm{m}}^{-1}$ behaves similarly and allow for identifying XDR cases. For the XDR case, the mean probabilities to detect three Raman lines are equal for all of them and higher than for alternatives. Thus, this ensemble-based picture supports the conclusion on their action as markers.

Comparison of data for pathogen population and individual cells

In bacterial populations or conglomerates, characteristic peaks have already been found, such vibrational modes related to the DNA and RNA backbone at $1\text{,}123\mathrm{c}{\mathrm{m}}^{-1}$ were also identified (Dhankhar, 2018) for MDR pulmonary strains at Fig. 4C. The Raman peaks located at shift position $1\text{,}152\mathrm{c}{\mathrm{m}}^{-1}$ represented a cell wall’s N-H, C-H in-plane bendings were also identified (Ullah et al., 2022) for sensitive extrapulmonary strains at Fig. 5C. Separately, it is worth noting that the vibrations corresponding to purines and pyrimidines of DNA at $722\mathrm{c}{\mathrm{m}}^{-1}$ were successfully identified (Zyubin et al., 2019) for XDR pulmonary strains as in Fig. 4A. The cytosine and uracil band modes at $772/776\mathrm{c}{\mathrm{m}}^{-1}$ were discovered in the article (Zyubin et al., 2019), for XDR Pulmonary strains at Fig. 4A.

It should be noted that the peaks of pyrimidine and purine bases, as well as individual bases related to DNA, are manifested mainly in MDR and XDR strains and are completely absent in sensitive strains. This fact requires further and more detailed consideration, such as the identification of metabolic features of such mycobacteria.

However, another marker thymine in DNA was found for individual bacteria at the wavelength $747\mathrm{c}{\mathrm{m}}^{-1}$, which is characteristic for MDR pulmonary strains. Nevertheless, overall, although more detailed information about specific marker peaks can be revealed on individual bacteria, the main information about differences between strains can also be obtained on conglomerates.

In principle, the differences between conglomerates and individual bacteria are based on the heterogeneity of the culture being studied. When studying conglomerates, average values of molecular structure peaks are identified, whereas, in individual bacteria, individual peaks may not be detected or may be greatly shifted depending on the state of the cell itself, in particular on the “age” of the cell. Nevertheless, the overlapping peaks of phenylalanine and 5-methylcytosine found in XDR lung strains that coincide with conglomerates (Zyubin et al., 2019) may serve as markers for strains with extra drug resistance. It should be pointed out that attention to a possible interplay between 5-methylcytosine and drug resistance as involved in the process of methylation of certain cytosine residues, has been drawn a long time ago (see Hattman, Gold & Plotnik, 1972) and considered now as one of the promising methylation markers, which can be detected by Raman spectroscopy (Li et al., 2018; Zhang et al., 2023). Concerning the case of M. tuberculosis, an epigenetic mechanism, which plays a significant role in DNA methylation mediated by 5-methylcytosine has been analysed recently using the methylation motif analysis by Phelan et al. (2018) and with metabolomic tracing by Gong et al. (2021). Thus, our results support this interpretation with results that can be achieved by more simple instrumentation.