Phylogenetic relationships of Mycobacterium tuberculosis isolates in Poland: The emergence of Beijing genotype among multidrug-resistant cases

Introduction The epidemiological situation of tuberculosis (TB) in Poland urges for its continuous and scrupulous monitoring. The objective of this study was to explore the genetic diversity of multidrug-resistant (MDR) and drug-susceptible (DS) Mycobacterium tuberculosis isolates from Poland with a combination of spoligotyping and high-resolution mycobacterial interspersed repetitive unit-variable number tandem repeat (MIRU-VNTR) analysis. The results were placed in the Northern and Eastern Europe context. Methods The study included 89 (39 MDR and 50 DS) M. tuberculosis isolates collected from as many patients between 2018 and 2021 in Poland. The analysis was done using spoligotyping, and MIRU-VNTR typing at 24 standard loci. The data were compared to those available on Poland and neighbors and global M. tuberculosis datasets. Results The main identified families were Beijing (28.1%) and Haarlem (16.8%) while 34.8% of isolates were in the heterogeneous L4-unclassified group. Although the Beijing family was the most prevalent (61.5%) among MDR-TB cases, it accounted for only 2% of DS isolates. Among foreign-born patients, a higher ratio of MDR isolates were observed when compared with those who Poland-born (64.3% vs. 40%). Furthermore, all patients from the Former Soviet Union (FSU) countries were infected with MDR-TB. Discussion Whereas DS M. tuberculosis population in Poland is dominated by L4 isolates, MDR isolates are mostly of the Beijing genotype. The rise in the prevalence of the Beijing isolates in Poland, coupled with high proportion of the Beijing genotype among foreign-born TB patients may reflect an ongoing transmission of this family, imported to Poland mainly from FSU countries.

1 Introduction 2 Materials and methods

Bacterial isolates
The study included a convenience sample of an 89 (39 MDR and 50 DS) M. tuberculosis isolates, deposited in the culture collection of the (i) Department of Internal Medicine, Pulmonary Diseases and Allergy of the Medical University of Warsaw (n=31), (ii) Mazovian Centre for the Treatment of Lung Diseases and Tuberculosis in Otwock (n=23) and (iii) Mycobacterium tuberculosis Laboratory, Krakow (n=35) ( Table S2). The sample covered all MDR isolates retrieved during the study period. Drugsusceptible isolates were selected based on the availability of medical records. The isolates were recovered from different pulmonary TB patients (70 males, 19 females; age range, 19 to 90 years; mean age, 51 ± 16.3 years), diagnosed with TB between 2018 and 2021. The study sample represented 20.7% and 0.3% of all bacteriologicallyconfirmed MDR-and DS-TB cases respectively, reported in Poland during the survey period. However, a moderate sample size is a notable limitation of this study.
Primary isolation, culturing, and species identification were performed with standard mycobacteriological methods (Clinical & Laboratory Standards Institute, 2018).
All personal data were anonymized, therefore informed consent of patients was not needed (Medical University of Warsaw Bioethics Committee decision no. AKBE/22/2019). All experimental protocols and methods were carried out in accordance with the guidelines and recommendations of the Medical University of Warsaw.
Apart from country of birth and place of isolation no other epidemiological data were available for this study.

Drug susceptibility testing
Drug susceptibility testing (DST) was performed using either the standard 1% proportion method on the Löwenstein-Jensen medium or BACTEC MGIT system (Becton Dickinson, USA), following the WHO protocols (World Health Organization, 2018). The M. tuberculosis H37Rv reference strain was used as a quality control. Drug resistance profiles were categorized in accordance with WHO updated criteria (World Health Organization, 2020).

DNA isolation
Extraction of genomic M. tuberculosis DNA was done using the cetyl-trimethyl ammonium bromide (CTAB) method, as described elsewhere (de Almeida et al., 2013). The purified DNA was dissolved in TE buffer and quantified with the NanoDropTM 2000 Spectrophotometer (ThermoFisher Scientific, USA). The DNA samples were diluted to the required concentration (ca. 10 ng/µL) and stored at -20°C until used.

Molecular typing
Spoligotyping was performed using commercial kits (Ocimum Biosolutions, India) and following the published protocol (Kamerbeek et al., 1997). All profiles were assessed by two independent researchers. Spoligotype International types (SITs) of M. tuberculosis were assigned according to SITVIT2 (http:// www.pasteur-guadeloupe.fr:8081/SITVIT2/). MIRU-VNTR analysis was done at 24 standard loci, essentially as described previously (Supply et al., 2006), except that the reaction components were adjusted as listed in Table S1. DNA fragments were visualized and analyzed using capillary electrophoresis system (Qiaxcel, Qiagen, USA) (Qiagen, 2016). VNTR alleles were considered as discrete variables.
For both, spoligotyping and MIRU-VNTR analysis, M. tuberculosis H37Rv and M. bovis BCG reference strains were used as quality controls.
The adopted hierarchy of clades, from the highest to lowest rank, states as follows: lineage, family, SIT (or MLVA type), spoligotype (or MIRU-VNTR type). A spoligotyping cluster was defined as two or more isolates sharing identical spoligotypes. The same criteria, i.e. exact match at 24 MIRU-VNTR loci was applied for a MIRU-VNTR cluster.
The discriminatory power of spoligotyping and MIRU-VNTR at all loci was calculated with the Hunter and Gaston discriminatory index (HGDI), using the following formula: where N is the total number of isolates, and nj is the number of isolates representing each type (Hunter and Gaston, 1988).

Dendrogram and minimum spanning tree construction
Dendrogram was constructed based on MIRU-VNTR typing data, using MIRU-VNTRplus software and UPGMA algorithm (Allix-Beǵuec et al., 2008) and a dataset of 186 MIRU-VNTR profiles, deposited in the database as a reference.
Minimum spanning trees (MST) were drawn based on MIRU-VNTRplus software (Allix-Beǵuec et al., 2008). Single-locus variants (SLVs) within MSTs were defined as isolates which differed from the ancestral isolate at one of the MIRU-VNTR locus. Those SLVs formed MST clusters.

Family assignment
Based on the position of the isolates on MIRU-VNTR phylogenetic tree ( Figure S1) and single VNTR-locus signatures of LAM sublineages, essentially described elsewhere (Mokrousov et al., 2014) (Figure S1), the isolates were assigned to the final genetic families (Table S2).
Since T family is known as ill-defined, polyphyletic, and heterogeneous group (Coll et al., 2014), its designation was not used. Isolates grouped within L4 lineage, with undefined families were designated as L4-unclassified.

Genetic diversity of Mycobacterium tuberculosis isolates
In total, 40 spoligotypes were identified, split into 10 clusters (n=59, 66.3%, 2-19 isolates per cluster) and 30 (33.7%) unique patterns (Tables 1, 2). Fourteen (15.7%) isolates had patterns not deposited previously in the SITVIT2 database, and thus were called orphan types. Six of them formed two new STs, A and B, with 4 and 2 isolates, respectively. Other 8 orphan types were labeled as Orphan C-H (Table 1).
Four MIRU-VNTR clusters were described, shared by 2 to 9 isolates, totaling 15 (16.8%) isolates ( Table 2). The largest cluster (n=9) consisted of 6 MDR and 3 pre-extensively drug resistant (pre-XDR) isolates, all but one isolated from Polish patients in Masovian voivodeship (one was from Ukraine, isolated in Lesser Poland). Second cluster included two pre-XDR isolates from polish patients, yet obtained in different voivodeships (Masovian and Lesser Poland). In the third cluster there was one MDR and one pre-XDR isolate, from patients born in Poland or Georgia, retrieved in the same voivodeship (Lesser Poland). Fourth cluster comprised two susceptible isolates from polish patients, isolated in Masovian voivodeship. The remaining 74 (82.2%) isolates had unique patterns. MIRU-VNTR-based phylogenetic tree is depicted in Figure S1.
Among patients who were foreign-born, a higher ratio of MDR isolates were observed when compared with those who were born in Poland, i.e. 64.3%; (9/14) vs. 40% (30/75). Furthermore, all patients from the FSU had MDR-TB.
At the family level, L4-unslassified (26/45; 57.8%) and Beijing (9/30; 30%) were the most common families among Polish DS-and MDR-TB patients, respectively. Approximately half (8/14; 57.1%) of the foreign-born patients were infected with Beijing strains, compared to one-fourth (22.7%; 17/75) of Polish natives. Single isolates of CAS-Delhi (n=1) and CAS (n=1) families were detected in an India-born and Nepal-born patient, respectively. Importantly, the isolates from the two patients were the only representatives of L1 and L3 in the study sample.

Discussion
The rise and expansion of DR-TB, galvanized by the increased migration flows, represents one of the main challenges of TB control in Europe. However, molecular epidemiological studies of TB in eastern European countries, including Poland, are seriously lacking. The present study provides an important insight into the genetic diversity of DS and MDR M. tuberculosis strains circulating in Poland.

Phylogeny of TB in Poland
Among the whole population studied, L4 predominated (69.7%), with a prevalence of 94% and 38.5% for DS-and MDR-TB isolates, respectively. This is consistent with previous studies from Poland, in which the L4 lineage accounted for a similar proportion of isolates (69.5% and 71%) (Krawczyk et al., 2011;Kruczak et al., 2019). However, the previously reported L4 prevalence among MDR-TB isolates was 2-fold higher, and reached 65.2% (Jagielski et al., 2010). L4 has been repeatedly identified as a major lineage in countries neighboring Poland, such as Latvia (67.6%) (Pole et al., 2020) and Belarus (51.5%) (Zalutskaya et al., 2013). At the family level, MDR was found among all families of L4, with somewhat higher prevalence among LAM and L4-unclassified SIT53 compared to other molecular types. This finding is unsurprising, as previous studies have reported the associations of SIT53 and different LAM genotypes with increased drug resistance and MDR-TB (Sheen et al., 2013;Ioannidis et al., 2017).
The second largest phylogenetic group included L2 isolates, all of the Beijing family, and comprised 28.1% of the analyzed isolates. This lineage was found with a much higher prevalence (61.5%) among MDR-TB cases compared to DS-TB (2%) (P<0.05). According to previous reports, the frequency of the Beijing family isolates was strikingly lower in Poland 10-15 years ago, and accounted for 3.9% and 8% among DR population (Sajduda et al., 2006;Jagielski et al., 2010). On the contrary, high proportion of this family among TB is characteristic for FSU countries, such as Ukraine (33%) (Dymova et al., 2011) and Russia, Kaliningrad (45.6%) (Mokrousov et al., 2009).
In this study, almost all (96%) Beijing isolates were MDR. The Beijing strains are epidemiologically important since they have been associated with an increased acquisition of drug resistance, enhanced virulence, and high transmissibility (Mokrousov et al., 2006;Lasunskaia et al., 2010;Barbier and Wirth, 2016). Furthermore, as the predominant genotype in East Asia, the Beijing family has been emerging in various areas of the world and is often associated with disease outbreaks (Luo et al., 2014;Erie -20. et al., 2017;Perdigão et al., 2020). The increasing prevalence rate of the Beijing family in Poland might be attributed to its importation from FSU and its local active circulation and/or withincountry transmission.
Single M. tuberculosis isolates of L1 and L3 were DS and were recovered from patients born in India and Nepal. This is consistent with studies from Czech Republic, Germany, and Slovakia, where the two lineages (L1 and L3) accounted for less than 2.5% of the study sample (Roetzer et al., 2011;Kozińska et al., 2019;Pinkovaé t al., 2019). Strains of the L1 and L3 are characteristic mainly for the Asian continent and usually are imported to European countries with migration flows. Phylogenetic analysis using MST revealed the presence of four MST clusters. The first was formed by nearly 75% of the Beijing isolates, mostly of 94-32 or 94-15 MLVA type ( Figure 1A). Importantly, the 94-32 MLVA type, also known as Central-Asian Russian cluster, is indistinguishable by MIRU-VNTR from the Central Asia Outbreak (CAO) cluster (Merker et al., 2020). The isolates within this cluster were obtained from patients from different cities, had two different SITs (1 and 265) and different resistance profiles (MDR or pre-XDR). All these findings may indicate a local evolution of Beijing isolates or transmission in a distant past, rather than a recent transmission event. Type 94-15 differs at one MIRU-VNTR locus when compared to type 94-32, one of the major Beijing types, circulating in Russia . Interestingly, type 94-15 was not found in the largest global Beijing dataset  nor in the recent large study in six provinces of Northwestern Russia (Vyazovaya et al., 2023). Thus, the 94-15 type might be a newly formed local clone. The second MST cluster was formed by two MDR-TB Beijing isolates of different MLVA types, with either Georgian or Ukrainian origin. The third MST cluster comprised two Haarlem isolates, of different SITs (50 and 207). The fourth MST cluster was formed within the L4-unslassified group, with three isolates of "?-15" MLVA type ( Figure 1D). The isolates within this cluster were all DS, harbored identical SIT, and were recovered from patients living in the same city. Apart from the domicile, no other epidemiological links were revealed to support a direct transmission between the patients within MST clusters.
Concerning the origins of the isolates, the L4 strains circulating in Poland were only isolated from Polish patients. On the other hand, the L2 Beijing genotype was more frequently found among foreign-born patients (P<0.05). Overall, the high prevalence of the Beijing isolates in Poland evidenced in this study, coupled with a high proportion of the Beijing genotype among foreign-born TB patients may reflect an ongoing circulation of this genotype imported to Poland mainly from FSU countries (Niemann et al., 2010;Konstantynovska et al., 2019). Of note, one patient born in Poland was infected with type MLVA 1071-32, which has been recently described as an emerging resistant Beijing variant that is circulating in Siberia. This variant has been rarely found in the European part of Russia and only sporadically in Balkan countries (Mokrousov et al., 2021).
Interestingly, a patient from Moldova was infected with MDR L4 Ural strain of SIT262. This spoligotype was recently shown to be associated with pre-XDR-TB and belongs to the Ural family circulating in Moldova since the 1990s (Sinkov et al., 2018).

Mycobacterium tuberculosis phylogeography in Poland and Western/Central Europe
To place our results in a wider phylogeographic context, the main genotype families from this study were compared with those observed previously in Poland and neighboring countries (Sajduda et al., 2006;Konstantynovska et al., 2019;Pinkováet al., 2019;Merker et al., 2020;Yang et al., 2022) or deposited in the SITVIT2 database ( Figure 3). In total, there are 434 isolates from Poland in the SITVIT2 database. All were deposited between 2000 and 2004 and were mostly resistant to isoniazid and streptomycin, or had no detailed drug resistance profile.
In previous studies from Poland, the Beijing prevalence among DR-TB population was remarkably lower than in this study (61.5%), and varied from 3.9% in 2006 up to 20.6% in 2015 (Sajduda et al., 2006;Jagielski et al., 2010;Kozinśka and Augustynowicz-Kopec, 2015). As depicted in Figure 3, Beijing is one of the most common genotypes of Distribution of genotype families in multidrug-resistant and susceptible isolates.
M. tuberculosis in the FSU countries. Data of patient origin might suggest that migrants from FSU are significant source of Beijing isolates and thus MDR-TB in Poland.
The prevalence of the mainly DS-TB Haarlem strains in Poland has slightly decreased between 2010 and the present, i.e. from 36.4% % to 16.8% (Krawczyk et al., 2011). A similar trend was observed in Estonia, between 1999 and 2015 (Pole et al., 2020). It can be only speculated that infections with Haarlem isolates are associated with better treatment outcomes, and thus are easier to eradicate, if lowtransmissible (Reiling et al., 2013). Haarlem is a "European" genotype with a clear geographically increasing East-to-West gradient (Smit et al., 2014;Pichat et al., 2016;Vluggen et al., 2017;Pole et al., 2020), and low prevalence in FSU countries.
During the past 20 years, the proportion of LAM isolates in Poland remained at a similar level. Most of LAM isolates in Poland belong to the LAM-RUS branch that was described in Russia and Kazakhstan and in Latvia and Estonia, too. Interestingly, the XDR isolates of SIT266 (LAM family), highly-prevalent in Belarus, were not found in this study (Zalutskaya et al., 2013). This SIT remains also at a low prevalence in Russia (Mokrousov et al., 2021). This might be largely explained by the lack of mass migration from Belarus to Russia and to Poland.

Conclusions
This work describes the recent genetic diversity of MDR and DS M. tuberculosis strains circulating in Poland, assessed with a combination of spoligotyping and 24-loci MIRU-VNTR-typing. There are two major findings from the study. First, the populations of DS-and MDR-TB isolates differed significantly. The DS-TB was dominated by L4-unclassified isolates, whereas the MDR-TB isolates were mostly of the Beijing genotype. Second, the rise in the prevalence of Beijing isolates in Poland, along with high proportion of Beijing genotype among foreign-born TB patients may reflect an ongoing and successful transmission of this family, imported to Poland mainly from FSU countries. An ongoing local circulation of the Beijing family is supported by its genetic structure, i.e. most of the Beijing isolates were SLVs. On the contrary, significant distance of LAM, Haarlem, and L4-unclassified strains may be due to divergent sources of some recently imported isolates but can also be explained by their long-term local evolution in Poland.

Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding authors.

Author contributions
ZB: performed the analysis, wrote the paper, MM: performed the analysis; AB: performed the analysis; MP: collected the isolates; JK: collected the patients clinical data; RK: collected the patients clinical data; PS: conceived and designed the analysis; IM: performed the analysis; wrote the paper; TJ: conceived and designed the analysis; wrote the paper. All authors contributed to the article and approved the submitted version.

Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher's note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fcimb.2023.1161905/ full#supplementary-material SUPPLEMENTARY FIGURE 1 Dendrogram based on the 24 VNTR loci for 89 M. tuberculosis isolates from Poland (yellow boxes, with either black or red borders, representing drug susceptible and multidrug-resistant isolates, respectively) and 186 reference isolates from SITVIT2 database (grey boxes). Only families identified in polish isolates are labeled.

SUPPLEMENTARY TABLE 2
Data on 89 isolates included in this study.