*Article* **The New Human** *Babesia* **sp. FR1 Is a European Member of the** *Babesia* **sp. MO1 Clade**

**Claire Bonsergent 1,\*, Marie-Charlotte de Carné 2 , Nathalie de la Cotte <sup>1</sup> , François Moussel <sup>3</sup> , Véronique Perronne <sup>2</sup> and Laurence Malandrin 1,\***


**Abstract:** In Europe, *Babesia divergens* is responsible for most of the severe cases of human babesiosis. In the present study, we describe a case of babesiosis in a splenectomized patient in France and report a detailed molecular characterization of the etiological agent, named *Babesia* sp. FR1, as well as of closely related *Babesia divergens*, *Babesia capreoli* and *Babesia* sp. MO1-like parasites. The analysis of the conserved 18S rRNA gene was supplemented with the analysis of more discriminant markers involved in the red blood cell invasion process: *rap-1a* (rhoptry-associated-protein 1) and *ama-1* (apical-membrane-antigen 1). The *rap-1a* and *ama-1* phylogenetic analyses were congruent, placing *Babesia* sp. FR1, the new European etiological agent, in the American cluster of *Babesia* sp. MO1-like parasites. Based on two additional markers, our analysis confirms the clear separation of *B. divergens* and *B. capreoli*. *Babesia* sp. MO1-like parasites should also be considered as a separate species, with the rabbit as its natural host, differing from those of *B. divergens* (cattle) and *B. capreoli* (roe deer). The natural host of *Babesia* sp. FR1 remains to be discovered.

**Keywords:** *Babesia divergens*; *Babesia* sp. MO1; *Babesia capreoli*; *rap-1a*; *ama-1*; phylogeny

#### **1. Introduction**

Babesiosis is a tick-borne disease affecting a wide range of vertebrates worldwide. Symptoms of this disease are caused by the intraerythrocytic development of Protozoa of the genus *Babesia*, causing fever, jaundice, hemoglobinuria and anemia, possibly leading to death, depending on the *Babesia* species and the host. About one hundred species of *Babesia* have been described and transmission of the parasite between hosts occurs almost exclusively through Ixodid tick bites [1].

Even though humans are not natural hosts of *Babesia*, human infections caused by several different species of *Babesia* have been reported worldwide. *Babesia microti, B. duncani* (WA1) [2] and to a lesser extent *B. divergens*-like (*Babesia* sp. MO1 clade) [3] have been reported to cause disease in humans in the USA. The most prevalent species is *B. microti* responsible of infections that follow a relatively benign course [4]. In Asia, a few cases have recently been reported, caused by *B. divergens*-, *B. venatorum*- or *B. crassa*-like strains [5–7].

In Europe, the first case of human babesiosis was described in 1957 in Croatia [8,9]. In 1997, a review on human babesiosis in Europe reported 24 cases in splenectomized (20/24) and non-splenectomized (4/24) patients, 46% of which were fatal even in nonsplenectomized patients (2/4) [10]. At that time, the molecular diagnosis of the parasite species was lacking and cases were attributed to *B. divergens* based on morphological and/or serological grounds. A few years later, molecular analysis revealed a new etiological babesiosis agent, *Babesia* sp. EU1, which was found to be responsible for human cases in Austria, Italy [11], Germany [12] and Sweden [13]. Human babesiosis cases due to *B.*

**Citation:** Bonsergent, C.;

de Carné, M.-C.; de la Cotte, N.; Moussel, F.; Perronne, V.; Malandrin, L. The New Human *Babesia* sp. FR1 Is a European Member of the *Babesia* sp. MO1 Clade. *Pathogens* **2021**, *10*, 1433. https://doi.org/10.3390/pathogens 10111433

Academic Editors: Estrella Montero, Jeremy Gray, Cheryl Ann Lobo and Luis Miguel González

Received: 1 October 2021 Accepted: 30 October 2021 Published: 4 November 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

*microti* have been reported in Europe but usually they are imported cases from the USA [14], with only one autochtonous case reported to date in Germany [15]. Severe sporadic cases are usually attributed to *B. divergens* [16–24]. However, molecular confirmation of the species is not always undertaken [25–27]. Serological analysis and morphology on smears are not sufficient to ascertain *B. divergens* as the etiological agent. Even for specialists, the morphological distinction of *B. divergens* from *Babesia* sp. EU1 on smears is impossible [11,12]. Confirmed cases of babesiosis due to *B. divergens* can remain serologically negative [28–30], and serology can be confusing due to dot-like reactivity patterns of most human positive sera, concentrated at the apical pole of the parasite [31]. This reactivity pattern was confirmed with a serum from a clinically and molecularly confirmed human *B. divergens* case in Finland [18,31].

The phylogenetic group including *B. divergens* gathers different named or as yet unnamed species that are very closely related, and we will refer to this group as *B. divergens*like. *B. divergens* is indeed closely related and can be confused with *B. capreoli*, a parasite frequently found in roe deer in Europe, due to their high 18S rRNA sequence relatedness [32,33]. However, the conservation of three base differences in this gene between isolates of *B. divergens* (pathogen of cattle/humans) and isolates of *B. capreoli* (pathogen of roe deer), linked to different in vitro host ranges, allowed the delineation of these two species [33]. *Babesia* sp. MO1 also belongs to this phylogenetic group, and is responsible for a small number of severe or fatal human babesiosis in splenectomized patients in the USA [3,34–36]. Cottontail rabbits are the natural hosts of *Babesia* sp. MO1 [37,38]. In vitro cultivation features as well as in vivo experimental infections demonstrated the incapacity of this genetic variant to infect cattle, and, combined with 18S rRNA sequence differences, led to its species differentiation from *B. divergens* and the provisional name *Babesia* sp. MO1 [39–41].

In splenectomized patients, babesiosis due to *B. divergens* is fulminant with symptoms that appear within 1–3 weeks post infection, with persistent high fevers and headaches, followed by severe intravascular hemolysis, hemoglobinuria, and jaundice. Babesiosis in splenectomized patients is often fatal in Europe, as diagnosis and therefore adequate treatment are often delayed due to uncharacteristic flu-like symptoms and the infrequency of cases [42,43]. Severe symptoms and fatal cases also occur in non-splenectomized patients with known or unknown predispositions such as splenic dysfunction or a rudimentary spleen [18,19]. In immunocompetent patients, *B. divergens* infection is associated with flu-like symptoms shortly after a tick bite [29] or may remain asymptomatic [44].

In the present study, we describe an unusually mild babesiosis in an asplenic patient in France, originally suspected to be caused by *B. divergens*. Intrigued by the unusual course of infection, we carried out the molecular characterization of the responsible agent. As the 18S rRNA gene is rather conserved within the *B. divergens* taxonomic group, and therefore not sufficiently informative, we supplemented the molecular description with two additional and more variable markers: the apical membrane antigen 1 (*ama-1*) and the rhoptryassociated-protein-1a (*rap-1a*) genes. Molecular characterization and polymorphism of these two genes were also analyzed for different members of the *Babesia divergens*-like phylogenetic group, including the phylogenetically closely related *B. capreoli* and *Babesia* sp. AR1 identical to *Babesia* sp. MO1 but from a patient in Arkansas [35], and compared to available sequences of these genes for *B. divergens*.

#### **2. Results**

#### *2.1. Babesia* sp*. FR1: Report of the Clinical Case*

A 56-year-old man came into the emergency room with a suspected meningitidis syndrome. He was Caucasian and his only notable antecedent was a splenectomy in 2001 following a skiing accident (pneumococcal vaccine administered in November 2016, no *Haemophilus* nor meningococcal vaccines).

The patient stayed on the Île de Ré from 4 August 2017 to 24 August, then from 29 August to 3 September, at a house located at the edge of a forest. He also stayed in Béthune from 25 August to 28 August. He had a fever for 2 weeks associated with headaches. He then developed severe asthenia, sweats, tachypnea, myalgia, and elbow, shoulder, and knee arthralgia.

The first blood test (5 September) revealed thrombopenia: 99 giga/L, CRP 65.6 mg/L, ASAT 82U/L, ALAT 73U/L. On 9 September 2017, he developed vomiting, photophobia and a stiff neck, which led to the patient being transferred to hospital (11 September). Nothing specific was revealed by non-injected brain CT. Lumbar puncture was normal, and culture was sterile. Blood tests revealed the following: platelets 71 giga/L, leukocytes 8.40 giga/L (PNN 7.056 giga/L, lymphocytes 0.670 giga/L), hemoglobin 14.4g/dL, ASAT 57 U/L, ALAT 61 U/L, GGT 128 U/L, PAL 188 U/L, normal kidney function.

On 12 September, when admitted to the infectious disease unit, clinical examination showed fever, asthenia, and non-significant axillary lymph nodes. The same day, a blood smear showed red blood cells with *Babesia* corpuscles inside, reaching a parasitemia of 3.7% (Figure 1). Blood analysis revealed thrombopenia (60 giga/L) and hemolysis signs without anemia (Hb 14.4 g/dL, LDH 890 U/L). Lyme, HCV and HBV serologies were all negative, and protein electrophoresis was normal.

**Figure 1.** Blood smears of *Babesia* sp. FR1 used to diagnose the *Babesia divergens*-like infection of the patient. Human red blood cells infected with dividing pear shaped merozoites are visible, as well as rounded trophozoites. Bar = 5 µm.

Antimicrobial therapy was undertaken on the same day with Atovaquone (750 mg/12h) and Azithromycine (500 mg on day 1 then 250 mg per day). The patient rapidly felt better with apyrexia and disappearance of all symptoms. On 14 September, parasites were still detected on the blood smear and cytolysis was persistent.

Diagnosis of a *B. divergens*-like infection was confirmed by serology (IFAT with *B. divergens* antigen) with a titer of 1:1024 and by PCR on the 18S rRNA gene as described in materials and methods. Sequencing of the amplified 18S rRNA gene portion confirmed that the responsible agent was closely related to *B. divergens*, the most commonly responsible agent of human babesiosis in France, but different.

Control of the patient's infectious status was performed 16 months later. Serology using the same antigen remained positive with a titer reduced to 1:128. PCR was negative.

#### *2.2. Analysis of 18S rRNA Sequences and Position of Babesia* sp*. FR1 within the B. Divergens-Like Phylogenetic Group*

A 1641 bp sequence was obtained for *Babesia* sp. FR1, covering the positions that are discriminant among members of the *B. divergens*-like phylogenetic group: nucleotide positions 631, 663, 819, and 1637. The sequence is highly similar (99.95%) to published *Babesia* sp. MO1-like and *B. capreoli* sequences, with only one nucleotide modification at position 819 and 663 respectively. It is also related to *B. divergens* (99.9%) with two nucleotide substitutions at positions 631 and 1637 (Table 1).



likely vector according to [36]. <sup>b</sup> corresponding to position described in [33]. <sup>c</sup> not determined.

> The sequence of 1643 bp from the Arkansas case [35] obtained in this study (named *Babesia* sp. AR1) was 100% identical to the first *Babesia* sp. MO1 case from Missouri (Gen-Bank AY048113) [3], to the Kentucky case (GenBank AY887131) [34] and to the cottontail rabbit isolates [38]. They differ from *B. divergens* by three mutations at positions 631, 819 and 1637 (99.8% identity), and from *B. capreoli* by two mutations at positions 663 and 819 (99.9% identity) (Table 1).

#### *2.3. Major Differences in Ama-1 and Rap-1a Genes within the B. divergens-Like Phylogenetic Group*

Before analyzing the detailed sequence polymorphism of *ama-1* and *rap-1a* genes between the *B. divergens*-like phylogenetic group members, some major differences in gene sequences appeared on the alignment (Figure 2).

Complete *ama-1* sequences (sizes between 1803 and 1857 bp) were obtained in this study for four clonal lines of *B. capreoli*, for *Babesia* sp. FR1 and for *Babesia* sp. AR1, and were compared to *B. divergens ama-1* sequences [45]. An 18 bp sequence located between bases 553 and 570 was absent only in *Babesia* sp. AR1, corresponding probably to a deletion, which did not modify the translation frame (Supplementary Figure S1). *Babesia* sp. FR1 *ama-1* sequence differs from all the others by an insertion of a 36 bp sequence located between bases 1496 and 1531. The inserted sequence is highly similar to an upstream 36 bp sequence (differing by two nucleotides) and seems therefore to correspond to a gene conversion event of this small gene portion.


**Figure 2.** Major differences in *ama-1* and *rap-1a* genes between members of the *B. divergens*-like phylogenetic group. (**a**) Copy number, and presence (Y) or absence (N) of insertion/deletion in *ama-1* and *rap-1a* genes; (**b**) Schematic representation of the major differences and their positions in gene sequence. The colors used in the part (**a**) correspond to the colors used in the graphical representation of the corresponding deletions/insertions in the part (**b**). The deleted or absent regions are dashed.

> Regarding the *rap-1a* gene, partial sequences (sizes between 1203 and 1236 bp) were obtained for four clonal lines of *B. capreoli*, for *Babesia* sp. FR1 and for *Babesia* sp. AR1, and were compared to *B. divergens rap-1a* sequences [46]. Comparison of *rap-1a* sequences from *B. divergens*, *Babesia* sp. AR1, and *Babesia* sp. FR1 highlighted a 3 bp deletion (nucleotides 1012 to 1014) in *Babesia* sp. AR1 only and a 33 bp deletion located between bases 1034 and 1066 in the 12 *rap-1a B. divergens* sequences performed, which did not modify the translation frame. As superposed chromatograms were observed for *B. capreoli* at the 30 end of the *rap-1a* gene, the presence of multiple copies of this gene was suspected and confirmed by cloning/sequencing and subsequent specific amplifications of each gene copy from each of the four isolates. Two *rap-1a* gene copies were observed, which we named *rap-1a1* and *rap-1a2*. The *rap-1a1* copy was characterized by the absence of the two deletions, as found for *Babesia* sp. FR1. The *rap-1a2* copy contained the 33 bp deletion only and therefore resembled *B. divergens rap-1a*.

> These major deletions/insertions were not included in the sequence identity calculations, nor in the phylogenetic analyses, as they represent one-time, usually non-reversible events, with different evolutionary tempo and mode compared to substitutions.

#### *2.4. Intraspecific Sequence Diversity of Rap-1a and Ama-1 within B. divergens and B. capreoli*

The genetic variability within the *B. divergens* and *B. capreoli* strains have been analyzed previously for the 18S rRNA gene and no variations were found [33].

Intraspecific genetic variability of *B. divergens ama-1* and *rap-1a* genes was previously analyzed in studies performed at our lab and was found to be very low [45,46]. Sequence identities higher than 99.5% were highlighted for both genes, when comparing sequences of the same set of nine and twelve French isolates, for *ama-1* and *rap-1a* respectively (Table 2 and Table S1). The *ama-1* sequences showed between 99.9 to 100% conserved sites; two similar nucleotide substitutions were noted in *ama-1* sequences of 1505B F14, 3601B E2 and Rouen87 F5 isolates, compared to *ama-1* sequences of the other six clonal lines. The *rap-1a* sequences showed sequence identities between 99.6 and 100%, corresponding to a pairwise maximum of six nucleotide substitutions.


**Table 2.** Genetic variability of 18S rRNA, *ama-1* and *rap-1a* genes within *B. divergens* and *B. capreoli*.

In the case of *B. capreoli*, we analyzed the genetic polymorphism of *ama-1* and *rap-1a* for four isolates collected and cultivated at our laboratory from previous studies [33,47]. The *ama-1* sequences showed sequence identities between 99.7 to 100%. The 2770 F6 and CVD08 005 *ama-1* sequences were identical and up to nine polymorphic sites were identified resulting in six different nucleotide substitutions between *ama-1* sequences of 2704C and 2801 F10 isolates, and were compared to the other two identical *ama-1* sequences.

The two copies of the *rap-1a* gene (*rap-1a1* and *rap-1a2*) were identified in all four *B. capreoli* isolates. Sequence variability of each gene copy was low (less than seven nucleotide substitutions), and identities ranged between 99.5–100% and between 99.6–100% among the *rap-1a1* and the *rap-1a2* sequences, respectively. Sequence identities between *rap-1a1* and *rap-1a2* copies ranged between 96.1 and 96.5%. Most substitutions specific to each gene copy (39 positions) were non silent (39 substitutions resulting in 28 amino acid modifications), with a majority of substitutions on the first (nine substitutions) and/or second codon position (16 substitutions).

#### *2.5. Genetic Variability within the B. divergens-Like Phylogenetic Group*

As explained above, sequence identities were calculated without the regions corresponding to deletions/insertions and are presented as a contingency table including all three analyzed genes (Table 3). In general, the *ama-1* gene seemed to be more conserved than the *rap-1a* gene as the percentage of sequence identities ranged between 94.3 to 98.7% for *ama-1* and between 86.6 to 98.7% for *rap-1a*. For both genes, the lowest sequence identities were evidenced between *B. divergens* and all other analyzed *Babesia* within the group. The highest sequence identities for *ama-1* and *rap-1a* were obtained between *Babesia* sp. FR1 and *Babesia* sp. AR1 sequences (98.7% identities for both genes). *B. capreoli* was found to be more closely related to *Babesia* sp. AR1 and *Babesia* sp. FR1 than to *B. divergens*.

It was not possible to determine if one of the two copies of *B. capreoli rap-1a* was more related to the unique *rap-1a* gene sequence of other members of the phylogenetic group, as sequence identity values were highly similar.

**Table 3.** Contingency table for 18S rRNA, *rap-1a* (partial cds) and *ama-1* genes. For the 18S rRNA sequences, the number of nucleotide differences between 18S rRNA gene sequences are indicated in red, instead of the identity percentages. For *ama-1* and *rap-1a* genes, percentage of identities are indicated in green and blue, respectively. Sequence identities with each *rap-1a* copy (*rap-1a1* and *rap-1a2*) are indicated. The identities are calculated excluding the deletions and duplication indicated in Figure 2. Accession numbers of sequences used to perform the analysis are indicated in supplementary Table S1.


## *2.6. Phylogenetic Analysis*

The phylogenetic analyses based on 18S rRNA, *ama-1* and *rap-1a* genes were concordant and confirmed the placement of the *Babesia* sp. FR1 into the *B. divergens*-like phylogenetic group, with strong bootstrap values of 100 (Figures 3–5). According to the 18S rRNA phylogenetic analysis, and despite the high level of conservation of this marker, two sister groups were supported by good bootstrap values, and *Babesia* sp. FR1 clustered with *B. capreoli*, *Babesia* sp. AR1 and *Babesia* sp. MO1 (bootstrap of 73), and not with *B. divergens* (forming the second cluster supported by a bootstrap value of 89) (Figure 3). The separation of these two clusters was also well supported in the phylogenetic analysis with *ama-1* and *rap-1a* as markers (Figures 4 and 5). The *B. divergens* clade was supported by bootstrap values of 99 and 100 (*ama-1* and *rap-1a* respectively). The *B. capreoli*/*Babesia* sp. AR1/*Babesia* sp. FR1 clade was also well supported by bootstrap values of 100 (*ama-1* and *rap-1a*), but splits on the one hand into a subclade with *B. capreoli* (bootstraps of 100 and 83) and on the other hand into a second subclade with *Babesia* sp. FR1 and *Babesia* sp. AR1 (bootstraps of 100 and 99). The two *Babesia capreoli rap-1a* copies clustered into two sister groups with strong support (100).

**Figure 3.** Maximum likelihood unrooted phylogenetic tree of *Babesia* from the *Babesia divergens*-like phylogenetic group based on partial 18S rRNA sequences (1189 bp in the final data set). Branch support/bootstrap values are indicated at each node. *Babesia* sp. FR1 sequence obtained in this study is emphasized in red. (**a**) Scale bar indicates nucleotide substitution rate per site. (**b**) Topology of the tree allowing a better visualization of the bootstrap values; hosts of *Babesia* isolates are indicated.

**Figure 4.** Maximum likelihood unrooted phylogenetic tree of *Babesia* from the *Babesia divergens*-like phylogenetic group based on partial *ama-1* gene sequences (1728 bp in the final data set). Branch support/bootstrap values are indicated at each node. *Babesia* sp. FR1 sequence obtained in this study is emphasized in red. (**a**) Scale bar indicates nucleotide substitution rate per site. (**b**) Topology of the tree allowing a better visualization of the bootstrap values.

**Figure 5.** Maximum likelihood unrooted phylogenetic tree of *Babesia* from the *Babesia divergens*-like phylogenetic group based on partial *rap-1a* gene sequences (1138 bp in the final data set). Branch support/bootstrap values are indicated at each node. *Babesia* sp. FR1 sequence obtained in this study is emphasized in red. Scale bar indicates nucleotide substitution rate per site.

#### **3. Discussion**

Most human babesiosis cases are recorded in North America and are mainly due to *Babesia microti*, sporadically to *B. duncani* (*Babesia* sp. WA1) and to *Babesia* sp. MO1-like parasites. Sporadic cases were reported in Asia, Africa, and South America, with diverse and often partially characterized etiological agents [48]. In Europe, human babesiosis is rare and *B. divergens* is the main causal agent [43,49]. The most impacted countries are France, Ireland and Great Britain, and in France, Western regions and Normandy are most affected [10], due to substantial farming of bovines, the natural host of *B. divergens* [50].

The patient was most probably bitten by a tick on the Île de Ré, even if he had no recollection of a tick bite. This is the most probable place of tick acquisition by the patient, as it is close to a forest, where the abundance of the potential vector *I. ricinus* is high, increasing the risk of contracting tick-borne pathogens [51]. The patient was asplenic, which is also a major risk factor for severe or fatal babesiosis [10,43,49]. However, the symptoms in this patient developed slowly (two to three weeks between the onset of symptoms and admission to hospital) despite aspleny, while *B. divergens'* course of infection in such cases is usually fulminant [42]. Biological diagnosis of the provisionally named *Babesia* sp. FR1 was based on a blood smear, which led to the administration of antibabesial therapy (Atovaquone and Azithromycine) as soon as practicable. This treatment was effective, as symptoms rapidly disappeared, and parasite clearance was attested 16 months later by a negative PCR, correlated with a reduction of the serology titer.

Molecular characterization of *Babesia* sp. FR1 required a deeper analysis. Sequence and phylogenetic analysis of 18S rRNA revealed that it was genetically close but different from typical *B. divergens* isolates infecting cattle or humans (two polymorphic sites at positions 831 and 1637) [33], but that it closely resembled the American *B. divergens*-like parasite *Babesia* sp. MO1. Despite the genetic difference with *B. divergens*, the infection could be diagnosed using *B. divergens*-specific serological tools (IFAT), confirming anyway a close relationship with *B. divergens.*

We decided to explore new markers to improve knowledge on the *B. divergens*-like phylogenetic clade and to correctly position this new isolate within this species complex. We chose *rap-1a* and *ama-1* for two reasons. First, both genes code for proteins involved in the process of red blood cell invasion by the parasite [52], and as host range/specificity is an important biological feature in the description of this intra-erythrocytic obligatory parasite, they represent markers of interest. Second, we know from previous studies that both genes were well conserved among *B. divergens* isolates from cattle or humans [45,46]. Their interspecies divergence remained to be determined.

Regardless of the marker used, the sequences of *B. divergens* (cattle as a natural host) are grouped in a cluster well-separated from the other two clusters corresponding to *B. capreoli* and *Babesia* sp. MO1/AR1/FR1. Phylogenies based on more discriminant markers (*rap-1a* or *ama-1*) placed *Babesia* sp. FR1 in the cluster formed by isolates responsible for cases of human babesiosis in the USA represented by *Babesia* sp. AR1. This cluster is separated from the cluster of *B. capreoli* sequences, and from the cluster of *B. divergens* sequences. We can therefore conclude that *Babesia* sp. FR1 is not a *B. divergens*.

The phylogenetic group containing *Babesia* sp. MO1 is sometimes referred to as the *B. divergens* US lineage [53,54]. However, *Babesia* sp. FR1, which clusters with *Babesia* sp. MO1 and *Babesia* sp. AR1, was clearly acquired locally and is an autochthonous case as the patient did not travel to the USA in the months before the onset of the symptoms. Therefore, we not only confirm in the present study that *Babesia* sp. MO1-like sequences form a well-supported taxon, but we also highlight that the geographical distribution of this group is not restricted to the USA, but extends to Europe as it includes *Babesia* sp. FR1. The three clusters within *Babesia divergens*-like, i.e., *B. divergens*, *B. capreoli,* and *Babesia* sp. MO1 like might be associated with their natural host rather than with geographic distribution. Humans are only incidental hosts for parasites belonging to the *B. divergens*-like group, the natural hosts being cattle for *B. divergens*, roe deer for *B. capreoli* and rabbits for *Babesia* sp. MO1. The natural host for *Babesia* sp. FR1 has not been characterized but could also well be a Laporidae, especially in the Île de Ré context, a highly touristic and populated island where cattle and cervids are rare or absent, due to limited forested areas dominated by resinous trees (mainly maritime pine trees) and typical local productions (vineyards and salt marshes). The European rabbit (*Oryctolagus cuniculus*) is highly abundant on this island where it has been pullulating since the 2000s, and could therefore be the potential natural host for *Babesia* sp. FR1.

In this study, we did not include isolates described as *B. divergens* in sika deer described in Japan or in humans in China [6,53,54]. Our goal in this study was to characterize and correctly place the new *Babesia* sp. FR1 isolate in the phylogenetic group of *Babesia divergens*like, among biologically well-characterized isolates, i.e., whose host range has been studied and whose parasites have been cultured [33,38–41,55]. The isolates described in Japan from sika deer and named *B. divergens* are not included in the *B. divergens*-like phylogenetic group because they differ at the 18S rRNA sequence from all other members of this group by at least six conserved substitutions all of which are different from those described within this group. These isolates form a sister group to *B. divergens*-like. The name *B. divergens* should be reserved for isolates from cattle or humans whose 18S rRNA sequences match the many descriptions already published [18–23,29,30,33,49]. There is no evidence that isolates

from sika deer are capable of infecting either cattle, gerbils (*B. divergens* experimental host), or humans.

The sequences named as "*B. divergens*" and described in humans in China [6] actually match those described for *B. capreoli*, with the characteristic differences at positions 631 and 663 [33]. This information raises the possibility of human infections by *B. capreoli*, a species that has never been molecularly characterized as responsible for symptomatic cases of human babesiosis in Europe and described as not growing in vitro in human red blood cells [33]. However, short-term asymptomatic carriage of parasites could occur in humans in geographic areas with high parasite and vector prevalences. We unsuccessfully attempted to obtain DNA from these parasites to include them in our study.

In the present study, we characterized *rap-1a* genes in the *B. divergens*-like phylogenetic group. The *rap-1a* genes belong to a multigene family, and multiple copies have already been demonstrated in a few *Babesia* species: two copies in *B. bovis* and *B. canis*, four to five in *B. ovis*, at least seven copies in *Babesia* sp. Xinjiang, eleven in *B. bigemina* and twelve in the *B. motasi*-like group members [52,56–61]. The multiple copies of *rap-1a* are usually different, allowing their differentiation, except in the case of *B. divergens* where the presence of two identical copies was highlighted when its genome was sequenced [52]. It is highly probable that two identical and therefore undistinguishable copies of *rap-1a* exist in all the *B. divergens* isolates characterized. We cannot exclude the presence of two identical copies also for members of the *Babesia* sp. MO1-like clade. In all *B. capreoli* isolates analyzed, two different but closely related *rap-1a* copies, named *rap-1a1* and *rap-1a2* (sequence identities of about 96%) were identified. Each copy is equally different from either *B. divergens*, *Babesia* sp. FR1 or *Babesia* sp. AR1 *rap-1a* genes, and the two copies place as sister groups to each other, indicating a gene duplication that occurred after *B. capreoli* speciation. While genetic divergence occurred between the two *rap-1a* copies of *B. capreoli*, it was not the case between *B. divergens rap-1a* copies. *B. divergens rap-1a* gene has been previously characterized and its genetic variability among cattle and human isolates was found to be limited [51,62,63]. A greater sequence diversity among *B. capreoli* isolates compared to *B. divergens* was also highlighted in the case of the merozoite surface antigen *Bc37/41* compared to *Bd37*, and a greater selection pressure was hypothesized [64]. This could also explain the sequence divergence between the two copies of *rap-1a* in *B. capreoli* and not between the two *B. divergens* copies of *rap-1*. But a more recent event of the *rap-1a* gene duplication in *B. divergens* could also explain the difference in sequence divergence between the two copies. Whether the last common ancestor of members of the *B. divergens*-like group possesses one or two copies of *rap-1a*, or when and how many times *rap-1a* gene duplication occurred in the speciation process is difficult to evaluate. Despite the presence of all motives that characterized RAP-1 family members, *rap-1b* sequence identity with the other two copies of *rap-1a* in the *B. divergens* genome is extremely low (45%) and was not amplified with the primers used.

In conclusion, we describe here a case of human babesiosis in Europe (France) due to a *Babesia* isolate more closely related to the American *Babesia* sp. MO1 and AR1 than to *B. divergens*. Using two discriminant molecular markers, our study confirms the existence of three phylogenetic clades within the *B. divergens*-like group that would deserve the rank of species as their phylogenetic classification corresponds with their natural hosts; *B. divergens* natural host is cattle, *B. capreoli* infects mainly roe deer, and *Babesia* sp. MO1-like parasites probably infect Laporidae.

In Europe, diagnosis of human babesiosis is complicated due to its infrequency which often leads to a delayed detection and treatment. This delayed treatment promotes the development of a fulminant manifestation of this parasitic disease, in particular for asplenic or immunocompromised patients. Thus, the infection results often in the death of the patient, which could probably have been avoided by a more precocious diagnosis and treatment [23,26]. In Europe, the serological and molecular tools developed to diagnose *B. divergens* infections should principally be adequate to detect *Babesia* sp. FR1 infections. However, it needs to be taken into account that atypical immunofluorescence patterns (dots

and weak fluorescence of the parasite surface) may lead to a negative conclusion when carrying out an immunofluorescence test. Prophylactic treatments are advised, such as wearing long clothes and performing skin examination for tick detection after exposure to high-risk environments.

#### **4. Materials and Methods**

#### *4.1. Babesia Isolates and DNA Origins*

A preliminary identification of *Babesia* sp. FR1 responsible for the mild form of babesiosis was performed on blood smears stained with May-Grünwald Giemsa. Further diagnosis was carried out by serology (IFAT with *B. divergens* antigen) [31], as well as 18S rRNA gene amplification [33] and sequencing from blood DNA extracted using the Nucleospin Blood kit according to the manufacturer's instructions (Macherey-Nagel, Düren, Germany).

*B. capreoli* isolates were collected and characterized in previous studies performed at our lab [33,47]. We included in our analysis four in vitro cultivated isolates from roe deer blood samples or spleen, from three different regions of France (Supplementary Table S1). *B. divergens* isolates were also cultivated in vitro from acute piroplasmosis cases in cows (11 isolates) or in humans (one isolate) [55].

DNA from cultured *B. capreoli* clonal lines (2704C, 2770 F6, 2801 F10 and CVD08 005) was extracted as previously mentioned.

DNA from one American case of human babesiosis was kindly provided by Mayo Medical Laboratory, Rochester, USA. This fatal case occurred in an asplenic patient, in Arkansas in 2015, with a possible acquisition through transfusion [35]. Due to its geographical origin, we named this isolate *Babesia* sp. AR1. It was characterized as a *Babesia* sp. MO1-like babesiosis etiological agent.

As all these isolates are very closely related to *B. divergens*, we have used throughout the manuscript the terminology *B. divergens*-like phylogenetic group to qualify the following isolates, species or clades: *B. divergens*, *B. capreoli*, *Babesia* sp. MO1, *Babesia* sp. AR1, and the new *Babesia* isolate, named *Babesia* sp. FR1.

#### *4.2. Comparison of 18S rDNA Sequences within the B. divergens Taxonomic Group*

*Babesia* sp. AR1 18S rRNA sequences (560 bp) from the previously mentioned isolates were kindly provided by Mayo Medical Laboratory, Rochester, USA. Published *B. capreoli* as well as *B. divergens* 18S rRNA sequences were used as a comparison [33]. Their origin and accession numbers are described in Supplementary Table S1. The partial 18S rRNA sequence of *Babesia* sp. FR1 was obtained using the same primers [33] (Table 4). With the aim of obtaining sequences of comparable sizes, the 18S rRNA partial sequence of *Babesia* sp. AR1 was also amplified with these primers and sequenced. The alignment was done using the ClustalW program as implemented in the Geneious R6 software (https://www.geneious.com accessed on 1 October 2021).

#### *4.3. Amplification of Ama-1 (Apical Membrane Antigen-1) and Rap-1a (Rhoptry Associated protein-1) Genes for B. capreoli, Babesia* sp*. AR1, and Babesia* sp*. FR1*

PCR was performed to amplify the *ama-1* and *rap-1a* genes of *B. capreoli*, *Babesia sp.* AR1, and *Babesia* sp. FR1 isolates using ama1-S1/ama1-R3 and rap1-fw/rap1-rev primers respectively (Table 4). Reactions were carried out in 30 µL reaction mixtures containing 1 X GoTaq buffer, 4 mM MgCl2, 0.2 mM of each dNTP (Eurobio Scientific, Les Ulis, France), 1 unit GoTaq G2 Flexi DNA Polymerase (Promega, Madison, WI, USA), 0.5 µM of each primer and 1 µL of DNA template. The amplification conditions comprised 5 min at 95 ◦C followed by 40 cycles of 30 s at 95 ◦C, 30 s at the temperatures indicated in Table 4, 1 min 30 s at 72 ◦C, and a final extension at 72 ◦C for 5 min. The amplified fragments were purified with the ExoSAP-IT reagent according to the manufacturer's instructions (Affymetrix, Santa Clara, CA, USA) and sequencing was performed on both strands (Eurofins Genomics, Ebersberg, Germany) using the same primers for *rap-1a* gene, or using primers distributed

## along the sequence for *ama-1* gene (Table 4). Sequences were then assembled using the Geneious R6 software.


**Table 4.** Description of the primers used in this study for gene amplification as well as sequencing.

For the four isolates of *B. capreoli* (2704C, 2770 F6, 2801 F10 and CVD08 005), preliminary sequencing results of *rap-1a* gene highlighted superposed chromatograms at the gene 3 0 end, suggesting the presence of multiple copies of this gene, a frequent feature for this gene. The PCR products were therefore cloned in the pGEM-T easy vector according to the manufacturer's instructions (Promega, Madison, WI, USA), to determine the number and sequences of the putative different copies of the *rap-1a* gene. *Escherichia coli* strain BL21 was transformed with the plasmid constructions and colonies with the expected inserts were selected by direct colony PCR using vector primers: T7 and SP6. Recombinant plasmids were then isolated using the Nucleospin Plasmid kit (Macherey-Nagel, Düren, Germany) and both strands of the inserts were sequenced using vector primers (Eurofins, Genomics, Ebersberg, Germany). Then primers were designed to selectively amplify the different *rap-1a* gene copies of *B. capreoli* isolates (Table 4). Primers rap1-a1-fw or rap1-a2-fw were associated with primer rap1-rev to amplify the 30 part of respective *rap-1a* copies. Primers rap1-a1-rev or rap1-a2-rev were associated with primer rap1-fw to amplify the 50 part of *rap-1a* copies (same reaction and cycling conditions as above). PCR products were purified and sequenced as already mentioned.

#### *4.4. Comparison of Ama-1 and Rap-1a DNA Sequences for B. divergens-Like Phylogenetic Group Members*

The resulting *ama-1* and *rap-1a* DNA sequences of *B. capreoli*, *Babesia* sp. AR1 and *Babesia* sp. FR1 were aligned along with the published *ama-1* and *rap-1a* DNA sequences of *B. divergens* isolates using the ClustalW program as implemented in the Geneious R6 software.

#### *4.5. Phylogenetic Analysis*

Phylogenetic relationships within the *Babesia divergens*-like group were inferred using published sequences available in GenBank (Supplementary Table S1) and sequence data produced in the present study (18S rRNA, *ama-1* and *rap-1a* sequences). Sequences were aligned using Muscle as implemented in MEGA version X [65]. Phylogenetic analyses used a trimmed alignment of 1189 bp with complete deletion option for the 18S rRNA gene, 1728 bp and 1138 bp with complete deletion option for the *ama-1* and *rap-1a* coding sequences respectively. Maximum likelihood phylogenetic trees were produced using MEGA-X, with 1000 bootstrap replications based on the Tamura 3-parameter model [66] for the 18S rRNA gene and the *ama-1* gene, and based on the Kimura 2-parameter model [67] model for *rap-1a* gene. For the *ama-1* and *rap-1a* sequences, the three codon positions were included. The appropriate model of nucleotide substitution for ML analysis was selected based on the Bayesian Information Criterion (BIC) computed by MEGA-X.

#### *4.6. Genbank Deposition*

Nucleotidic sequences obtained in this study were submitted to GenBank with accession numbers MZ825347 and OK086051 for the 18S rRNA sequence of *Babesia* sp. FR1 and *Babesia* sp. AR1 respectively, MZ836259 and MZ836261 for the *ama-1* sequences of *Babesia* sp. AR1 and *Babesia* sp. FR1 respectively, MZ836258 and MZ836260 for the *rap-1a* sequences of *Babesia* sp. AR1 and *Babesia* sp. FR1 respectively. Accession numbers for *rap-1a* and *ama-1* sequences of *B. capreoli* are presented in Supplementary Table S1.

**Supplementary Materials:** The following are available online at https://www.mdpi.com/article/ 10.3390/pathogens10111433/s1, Table S1: Description of *B. divergens* (A) and *B. capreoli* (B) clonal lines used in the study and Genbank accession numbers. Figure S1: (A) Major differences in *ama-1* gene and AMA-1 protein sequences. (B) Major differences in *rap-1a* gene and RAP-1A protein sequences.

**Author Contributions:** Conceptualization, L.M. and M.-C.d.C.; methodology, C.B. and F.M.; validation, C.B., N.d.l.C. and L.M.; formal analysis, C.B.; investigation, C.B., V.P. and M.-C.d.C.; resources, V.P., F.M., M.-C.d.C. and L.M.; data curation, C.B.; writing—original draft preparation, C.B., M.-C.d.C. and L.M.; writing—review and editing, C.B., M.-C.d.C., N.d.l.C., F.M. and V.P..; visualization, F.M., C.B. and L.M.; supervision, L.M.; project administration, L.M.; funding acquisition, C.B., N.d.l.C. and L.M. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Institutional Review Board Statement:** Ethical review and approval were waived for this study, because this case report did not involve human subject research.

**Informed Consent Statement:** Written informed consent has been obtained from the patient to publish this paper.

**Acknowledgments:** The authors wish to thank M. Jouglin for technical support.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

