Brucellosis is a bacterial zoonotic disease caused by species belonging to the genus Brucella and results in high economic impact (McDermott et al., 2013). Brucella spp. may be transmitted to humans resulting in a severe disease requiring a specific and long-term antibiotic treatment with significant burden to public health systems (Franc et al., 2018). The genus Brucella currently contains 12 validly published species (Olsen and Palmer, 2014; Whatmore and Foster, 2021; Occhialini et al., 2022). Brucella melitensis, Brucella abortus, B. suis, B. ovis, B. neotomae, and B. canis are often referred to as “classical” Brucella species in the literature. These have been identified between 1887 and 1968 and were differentiated on the basis of phenotypic traits and host preference (Buddle, 1956; Stoenner and Lackman, 1957; Carmichael and Bruner, 1968; Moreno et al., 2002). Since 2007, the wider use of genetic methods of identification and differentiation has led to the identification of B. ceti, B. pinnipedialis, B. microti, B. inopinata, B. papionis, and B. vulpis, which also preserve the tradition of naming the species in accordance with their original host (with the exception of B. inopinata isolated from a breast implant; Foster et al., 2007; Scholz et al., 2008, 2010; Whatmore et al., 2014; Scholz et al., 2016). The use of molecular methods made it possible to identify other potential new species recovered not only from mammals but also from amphibians, reptiles, and fish (Soler-Llorens et al., 2016; Al Dahouk et al., 2017; Eisenberg et al., 2017; Mühldorfer et al., 2017; Eisenberg et al., 2020), leading to further expansion of our knowledge of the genus Brucella. Based on whole-genome comparisons, a merge of the Ochrobactrum and Brucella genus was recently proposed (Hordt et al., 2020; Leclercq et al., 2020). Nonetheless, the greatest impact on public health and livestock infections around the globe have to date only been caused by B. melitensis, B. abortus and B. suis (Godfroid et al., 2011; Franc et al., 2018). An increase of human cases due to B. canis is currently suspected (Hensel et al., 2018; Zhou et al., 2020).

The majority of brucellosis cases are reported in the Mediterranean countries, South and Central America, Africa, Asia, Arabian Peninsula, Indian subcontinent, Eastern Europe, and the Middle East (Pappas et al., 2006; Nicoletti, 2010). Over the past decade, there has been a decline in incidence in many regions previously considered to be highly endemic, but also new reservoirs have been identified, such as in Africa and the Middle East, possibly resulting from a better implementation of diagnostic methods (Wang and Jiang, 2020). The real incidence of the disease in the mentioned regions is most probably largely underestimated because registration is based on passively collected data (Dean et al., 2012).

In Kazakhstan, brucellosis remains a major livestock and public health problem. More than 1,300 cases of human brucellosis are registered annually corresponding to an incidence of 7.6 per 100,000 inhabitants. Seropositivity to Brucella antigens in cattle and small cattle is 0.6% and 0.4%, respectively (Charypkhan and Ruegg, 2022).

The high zoonotic potential, re-emergence of the infection in previously disease-free regions, and identification of new reservoirs underscore the need for modern molecular epidemiology approaches such as genotyping to trace source reservoirs and paths of introduction. For Brucella genotyping at subspecies level, two methods are most widely used, Multilocus Sequence Typing (MLST) and Multiple Loci VNTR polymorphisms (Variable Number of Tandem Repeats, MLVA; Whatmore et al., 2016; Vergnaud et al., 2018). Robust phylogenetic relationships can be obtained from nucleotide sequencing data owing to the strictly clonal evolution of classical Brucella spp. (Whatmore and Foster, 2021). The first MLST scheme included nine loci, seven housekeeping genes, the outer membrane protein gene omp25 and int-hyp. Six B. abortus MLST9 sequence types (STs) were initially described (Whatmore et al., 2007). Twenty-six B. abortus MLST9 STs are recorded in the current version of the Brucella spp. MLST database.¹ Including 12 housekeeping genes in the genotyping resulted in the more discriminatory MLST21. Thirty MLST21 STs were initially described in 172 B. abortus strains defining three clades (A, B, C including C1 and C2; Whatmore et al., 2016) and 43 STs are recorded in the current version of the Brucella spp. MLST database (see text footnote 1). The rare clades A and B include strains originating almost exclusively from Africa and allow defining the most ancestral nodes within the B. abortus phylogeny. Clades C1 and C2 are found on five continents and their presence in Africa seem to result from livestock importation (Whatmore et al., 2016; Vergnaud et al., 2018). Several MLVA genotyping schemes have been proposed for Brucella, one most commonly used scheme is MLVA16 (Scholz and Vergnaud, 2013). MLVA16 combines two panels of markers: one VNTR panel with a low discriminatory ability allows determining the species and main genetic lineages, and the other VNTR panel with a high discriminatory ability but low phylogenetic value allows differentiating strains in local outbreaks (Scholz and Vergnaud, 2013; Whatmore and Foster, 2021). The current version (Brucella v4_6_5) of the Brucella MLVA data hosted by MLVAbank² contains in vitro data derived from more than 1,400 B. abortus strains. The eight VNTR loci with the low discriminatory power (MLVA8) cluster B. abortus into 28 genotypes represented by at least two entries. MLVA alone is not suitable for phylogenetic reconstructions because of the high homoplasy at VNTR loci (Keim et al., 2004), but interestingly, MLVA using the low discriminatory markers (MLVA8, MLVA10 or MLVA11) and MLST clustering are highly congruent allowing to indirectly deduce phylogenetic signal from MLVA data (Vergnaud et al., 2018). Consequently, MLVA may constitute a first-line assay with low phylogenetic resolution until whole-genome sequencing (WGS) becomes sufficiently low-cost to be directly used as first line assay. MLVA genotypes can be readily deduced from WGS data with sequencing reads longer than 200 bp so that the highly discriminatory MLVA loci can also constitute a strain identity control tool (Vergnaud et al., 2018; Holzer et al., 2021; Pelerito et al., 2021).

The availability of whole-genome sequence (WGS) data opened the way to whole-genome and core genome MLST assays (wgMLST and cgMLST, respectively) with much higher resolution and phylogenetic value than these classical genotyping tools (Janowicz et al., 2018; Abdel-Glil et al., 2022). Genome-scale MLST assays as well as genome-wide SNP-genotyping confirmed the presence of the four major clusters A, B, C1, C2 in a collection of B. abortus with strains collected at a global scale (Whatmore and Foster, 2021; Abdel-Glil et al., 2022).

The current knowledge on the genetic diversity of B. abortus circulating in Kazakhstan is limited to MLVA-typing data. Interestingly, the MLVA genotyping investigations demonstrated low genetic diversity among strains circulating in Kazakhstan (Shevtsov et al., 2015; Daugaliyeva et al., 2018). Inclusion in the analysis of strains deposited since 1935 made it possible to describe the predominance of B. abortus clade C1, with genetic proximity of the majority of strains to Russian, Chinese and European strains. A few strains were assigned to B. abortus clade C2 (Shevtsova et al., 2016; Daugaliyeva et al., 2018). In order to better understand the population structure and origins of B. abortus in Kazakhstan, we selected 49 representative strains from Kazakhstan for whole-genome sequencing. We also selected 36 B. abortus strains representing the genetic diversity previously uncovered by MLVA within the historical collection maintained by Inrae, Nouzilly, France (Vergnaud et al., 2018) to determine phylogenic relations among B. abortus strains and to evaluate the usefulness of genomic data in the epidemiological control of the infection in a highly endemic region, such as Kazakhstan.

Brucellosis in cattle remains a major problem for cattle breeders in Kazakhstan. Genetic monitoring of genotypic dynamics is not included in the infection control strategy at the state level, but is implemented in the country through scientific grants and is mainly based on MLVA technology. This is the first study to describe the genetic diversity of B. abortus isolated between 1947 and 2015 in the territory of Kazakhstan using genome-wide sequence data that confirmed the presence of the two major subclades C1 and C2. Clades A and B were absent, in agreement with the previously published MLVA clustering analysis (Shevtsova et al., 2016; Daugaliyeva et al., 2018) and well established very strong geographic association of clades A and B with East and West Africa, respectively, (Whatmore et al., 2016; Vergnaud et al., 2018).

In this study, in five among the 10 MLVA16 genotypes represented by at least two strains, unrelated wgSNP genotypes were observed. The high percentage of homoplasia in the VNTR analysis is primarily due to the targeted selection of strains with identical MLVA16 from outbreaks which are unrelated by epidemiological data, and to the high level of homoplasia associated with the most variable VNTR loci. Four out of five cases of MLVA16 homoplasia involved combinations of the “historical” and currently circulating strains. Therefore, MLVA16 homoplasia seems higher in strains collected over an extended period of time. Eight wgSNP clusters were identical in MLVA16 or differed only in one or two hypervariable loci (Bruce07 or Bruce09), generally supporting the MLVA clustering. Previous studies of B. abortus strains using whole-genome SNP analysis and MLVA demonstrated similar results in outbreak differentiation and detection of imported strains, with a conclusion of a need to use genome-wide SNPs for reliable phylogenetic inference (Garofolo et al., 2017; Allen et al., 2020; Suárez-Esquivel et al., 2020; Holzer et al., 2021). Comparable strain clustering for MLVA11 and wgSNP allows MLVA in combination with epidemiological data to be considered as a first choice of methods to select strains for WGS in highly endemic regions, but only WGS analyzes will provide sufficiently precise phylogenetic information and might progressively become the first choice methods if global sequencing costs keep decreasing.

The identification of very closely related wgSNP genotypes in different regions with a difference of isolation by several years indicates a long-term circulation of B. abortus genotypes in Kazakhstan. Also, identification of the same wgSNP genotype cluster in unrelated outbreaks in the same region, such as Kostanay, indicates circulation of infected animals between farms. Thus, uncontrolled cattle trade and movement, as well as keeping animals from different farms on the same pastures, is postulated to be among main factors in the spread of the cattle brucellosis infection in Kazakhstan (Syzdykov et al., 2014; Charypkhan et al., 2019).

The Bayesian phylodynamic approach suggested that B. abortus lineages currently circulating in Kazakhstan were introduced in the 19th-20th centuries from Europe, mainly from Russia (North Caucasia). It may be interesting to note that evolution of subclades C1 and C2 showed similarities in this respect. The import of clade C2 into the United States was dated to the end of the 17^th century (Kamath et al., 2016), with closest neighbor lineages corresponding to strains isolated in Western Europe.

The introduction of B. abortus clade C1 in Kazakhstan might have happened during periods of human migration and because of the importation of numerous livestock for breeding with native breeds of cattle. Migration processes to Kazakhstan from the territory of the Russian Empire, and later the USSR, began with the accession of the northern territories of Kazakhstan in the Russian Empire in the 1730s and continued until the 1970s (Qazaqstan_tarihy, 2018; Presidential_Library_RF, 2022). The imperial period was characterized by mass migration of peasants, whose migration had an impact on the livestock breeding system of the nomadic people. Migrations had a wave-like character and were associated with the abolition of serfdom and the resettlement of free peasants from the European part of Russia since 1861. In 1889 the law on resettlement provided land plots and loans to peasant settlers. Crop failure and famine in European Russia and the Stolypin reforms of 1906–1911 further stimulated these migrations (Rather and Abdullah, 2018). In 1897 Russians made up 12 percent of the total population of Kazakhstan, i.e., 600,000 inhabitants. From the end of the 19th century to 1916, about 1,400,000 European Russians arrived in the Kazakh steppes, making up 40% of the population of the steppe regions of Kazakhstan (Bell-Fialloff, 2016). There are no exact data on imported livestock for that period, however, there are some data indicating that people moved with their livestock, and original Kalmyk-cattle-breeds appeared in Kazakhstan along with immigrants from Voronezh, Stavropol, Astrakhan, and other provinces of southern Russia (Belmont, 1960). During this period, many Kazakhs were forced to rebuild their traditional way of life with the transition to agriculture and a semi-nomadic lifestyle due to the seizure of land in favor of settlers for farming and taking rangelands to graze Russian riches’ cattle (Sailaukyzy et al., 2018). In the structure of Kazakhs herds, the number of cattle had increased because of greater demand for cattle meat. At the end of the 19th and the beginning of the 20th century, selective transformation of aboriginal cattle began, for which a massive import of cattle from various regions of Russia was carried out. For example, to create the Alatau breed from 1904 to 1940, Swiss and Kostroma cattle were brought to Kazakhstan and Kyrgyzstan from the Smolensk, Sumy and Kostroma regions (Soldatov, 2001). The Aulieatinskaya cattle breed has been developed since 1885 by crossing the Dutch Black Pied breed with aboriginal cattle and subsequent improvement in the 1930s of the 20th century by East Friesian bulls (Dmitriev and Ernst, 1989). All these processes might have constituted opportunities for a wide distribution of European strains of B. abortus in Kazakhstan, but it is important to note that essentially subclade C1 strains contaminate Kazakhstan whereas Western Europe is associated with subclade C2. Subclade C1 is common in North Caucasus (Kovalev et al., 2021). At the same time, the disastrous large-scale collectivization carried out in 1929–1933 resulted in a four-fold reduction in the local cattle population (Zhumasultanova, 2021), which could lead to a reduction in the historically circulating strains of B. abortus. The observed reduction in the genetic diversity of B. abortus strains in Kazakhstan in the 21st century is possibly associated with the successful implementation of epidemic-control activities similar to methods implemented in the 1970–1980s. During 4 years, from 1981 to 1985, the incidence of brucellosis in cattle decreased from 3.5 to 2.2% and a reduction in the most affected areas was as high as 30% (Shablov et al., 1988). Another factor for consideration is a two-fold reduction in the cattle population from 1991 to 2000.

The isolation of two B. abortus subclade C2 strains Kaz025 and Kaz027 genetically closest to the B. abortus 104 M vaccine strain is intriguing. We could not find records of the use of this vaccine in Kazakhstan. The B. abortus 104 M strain was first isolated from an aborted fetus of cattle in the central European Russia in 1929, selected and proposed as a vaccine strain for human vaccination by Kh. S. Kotlyarova in 1950 (Kotlyarova, 1950; Shumilov et al., 1983). Despite proven immunogenicity, strain 104 M did not find wide application in the USSR and was used in experimental vaccination of cattle and in limited production trials on small cattle. In this connection, the probability of importation of a vaccinated animal is low. We are more inclined toward the introduction of genetically close pathogenic strains. Expanding the results of genome-wide data on B. abortus strains isolated in the central European Russia will improve understanding of the origin of strains Kaz025 and Kaz027.

The new data of various origins contributed in this report strengthen the strong association of clades A and B with East and West Africa, respectively. The topology of the observed phylogeny combined with human history is pointing to East Africa as current most parsimonious scenario for the origin of B. abortus. The WGS data analysis of B. abortus strains from Kazakhstan shows that currently circulating lineages were introduced only recently in Kazakhstan, most of them during the 19^th or 20^th century. The closest currently known lineages are present in Caucasia, in agreement with the history of recent population migrations. This recent introduction is reminiscent of the situation described in Costa Rica, in which WGS data analysis allowed to identify five independent introductions responsible for the current population structure of B. abortus in Costa Rica (Suárez-Esquivel et al., 2020). Costa Rica was contaminated by B. abortus subclade C2 strains, imported from neighboring countries, whereas Kazakhstan was contaminated essentially by subclade C1 strains. This reflects the progressive spread of Brucella worldwide.

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: Bioproject PRJNA892249 and PRJNA901374.

GV editing of the original manuscript. GV, AC, and AShe conceptualization, designed the experiments, and analyzed the data. AShe and AC wrote the first draft of the manuscript. KB, MZ, ES, AShu, AA, DK, and YR genotyping and wrote sections of the manuscript. TK bacteriological researches and wrote sections of the manuscript. All authors contributed to the article and approved the submitted version.

This study was supported by the Ministry of Education and Science of the Republic of Kazakhstan (grant no. AP08052352) and by the french «Agence Nationale de la Recherche» grant ASTRID430 Maturation ANR-14-ASMA-0002-02.

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.

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.