Genetic diversity of Babesia bovis studied longitudinally under natural transmission conditions in calves in the state of Rio de Janeiro, Brazil.

Serum and DNA samples from 15 naturally infected calves in Seropédica, Brazil, were obtained quarterly from birth to 12 months of age, in order to longitudinally evaluate their humoral immune response against Babesia bovis and the merozoite surface antigen diversity of B. bovis. Anti-B. bovis IgG antibodies were detected by an indirect fluorescent antibody test (IFAT) and enzyme-linked immunosorbent assay (ELISA). Using DNA amplification, sequencing and phylogenetic analysis, the genetic diversity of B. bovis was assessed based on the genes that encode merozoite surface antigens (MSA-1, MSA-2b and MSA-2c). The serological results demonstrated that up to six months of age, all the calves developed active immunity against B. bovis. Among the 75 DNA samples evaluated, 0, 3 and 5 sequences of the msa-1, msa-2b and msa-2c genes were obtained, respectively. The present study demonstrated that the msa-2b and msa-2c gene sequences amplified from blood DNA of B. bovis-positive calves were genetically diversified. These data emphasize the importance of conducting deeper studies on the genetic diversity of B. bovis in Brazil, in order to design diagnostic antigens and vaccines in the future.


Introduction
Bovine babesiosis is a hemoparasitosis caused by protozoa of the genus Babesia and transmitted by Rhipicephalus microplus ticks. The disease is mainly characterized by hemolytic anemia and fever, with occasional hemoglobinuria and death (Mosqueda et al., 2012). Bovine babesiosis is an important and widespread cattle disease in tropical and subtropical areas worldwide (Bock et al., 2004). The cattle industry is particularly affected by Babesia parasites, since half of the 1.2 billion cattle in the world are at risk of infection and disease (Bock et al., 2004;Grisi et al., 2014).
In Brazil, bovine babesiosis is a serious problem not only in enzootic unstable areas, where most of the calves are not infected with Babesia spp. prior to reduction of their maternal antibodies; but also in enzootic stable regions, due to introduction of cattle from tick-free areas (Trindade et al., 2010).
The well-characterized merozoite surface antigen (MSA) family comprises variable proteins that are expressed on the surface of B. bovis merozoites and sporozoites. MSAs play a crucial role during invasion of the parasite into host erythrocytes (Mosqueda et al., 2002;Yokoyama et al., 2006). This family includes the msa-1 gene and msa-2 loci. While msa-1 is a single genome copy gene, msa-2 comprises four tandem genes, namely: msa-2a1, msa-2a2, msa-2b and msa-2c (Florin-Christensen et al., 2002). Besides contain neutralizing-sensitive epitopes, these antigens are highly immunogenic and therefore have been considered to be candidates for development of vaccines against B. bovis (Mendes et al., 2019). However, the high diversity of these surface antigens is an important obstacle to their use for controlling bovine babesiosis (Genis et al., 2009;Altangerel et al., 2012).
Only live attenuated vaccines are currently available for prevention of clinical babesiosis caused by B. bovis (Shkap et al., 2007). However, clinical babesiosis cases have been increasingly reported among vaccinated cattle due to polymorphism of essential antigenic components of parasites, namely MSAs (Brown et al., 2006).
Few studies regarding the genetic diversity of B. bovis have been conducted in Brazil. For instance, high genetic diversity of B. bovis based on msa-1 sequences has been reported among cattle in the state of Bahia, northeastern Brazil. This finding highlights the importance of conducting extensive studies on this subject before designing immune control strategies in this country (Nagano et al., 2013). Recently, Matos et al. (2017) reported that while the B. bovis msa-1 and msa-2b gene sequences amplified from calves' blood samples from Taiaçu, São Paulo, southeastern Brazil, were genetically distinct, msa-2c sequences were conserved. Moreover, Mendes et al. (2019) targeted the msa-2b and msa-2c genes and reported high genetic diversity of B. bovis in beef cattle sampled in the Pantanal biome, Mato Grosso do Sul, central-western Brazil.
Considering the scarcity of studies investigating the genetic diversity of B. bovis in Brazil, coupled with the fact that msa genes may represent genetic markers for studying the heterogeneity of this protozoa, the aims of the present study were (i) to analyze the genetic diversity of B. bovis based on MSA (MSA-1, MSA-2b and MSA-2c); and (ii) to detect antibodies against B. bovis and Babesia bigemina in cattle sampled in a herd in the state of Rio de Janeiro, southeastern Brazil, through a longitudinal study.

Study design
A longitudinal study was conducted on a dairy herd at the Seropédica Experimental Station (latitude 22° 48' S; longitude 43° 41' W; and altitude 33 m), which is located in the metropolitan microregion of the city of Rio de Janeiro, southeastern Brazil. This experimental study was performed during two periods: the rainy season (between October and March, in which the mean temperature and precipitation are 26 ºC and 1320 mm, respectively) and the dry season (between April and September, in which the mean temperature and precipitation are 23 ºC and 170 mm, respectively). During May 2012 and May 2013, fifteen calves (Bos taurus taurus x Bos taurus indicus) were evaluated every three months from birth until 12 months of age. Calves aged 0 to 2 months were kept in a shed, in individual pens, and received 4 kg of milk per day. They also had access to an area of 0.5 ha of Brachiaria humidicula from the age of 15 days onwards. From 3 to 6 months of age, the calves received 4 kg of milk per day and were kept during the day in an area of 1.5 ha of Brachiaria decumbens and were brought in at night to individual pens. Between the ages of 7 and 12 months, the calves were transferred to an area of 3 ha of B. decumbens and Panicum maximum, where they were kept during the day, and were brought into a collective pen at night. All the animals were dewormed and deticked at monthly intervals using ivermectin (1 mg/50 kg; Merial, Brazil) (Silva et al., 2015).
The first blood sampling was performed after the calves had ingested the colostrum, i.e. not more than one hour after their birth. This was followed by sequential sampling at the ages of 3, 6, 9 and 12 months, thus totaling 75 samples. All the animals were apparently healthy at the time of sample collection. Approximately 10 mL of blood was collected from the coccygeal or jugular vein of each animal, into two types of vacutainer tubes: one containing buffered ethylenediaminetetraacetic acid (EDTA) and the other without EDTA. The samples were kept at 4 °C during transportation to the laboratory. In order to obtain serum samples, the blood samples collected without EDTA were incubated at room temperature for 1 h and then centrifuged at 1000 ×g for 15 minutes. The EDTA blood samples were store at -20 °C until DNA extraction. This project was approved by the university's ethics committee, under the protocol number 017259/14.

Indirect fluorescent antibody test (IFAT)
The serum samples were subjected to a previously described protocol for anti-B. bovis and anti-B. bigemina IgG antibody detection Barci et al., 1994). The antigenic substrates used in the serological assays were prepared in accordance with a previously described protocol . They consisted of blood smears containing erythrocytes parasitized by B. bovis or B. bigemina that had been obtained from splenectomized calves experimentally infected by B. bovis or B. bigemina. The results were analyzed under an epifluorescence microscope (Olympus BX60, Tokyo, Japan) with magnification of 40X. Serum samples that were reactive at dilutions ≥ 1:80 were considered positive.

Enzyme-linked immunosorbent assay (ELISA)
The serum samples were also tested for the presence of IgG antibodies against B. bovis and B. bigemina by an indirect enzyme-linked immunosorbent assay (ELISA), using crude soluble antigens, as previously described by Machado et al. (1997). The cutoff value for the optical density (OD) at 405 nm was determined to be two and a half times (2.5X) the mean value of the negative control serum samples (0.233 OD for B. bovis and 0.275 OD for B. bigemina). Absorbance values at or above this value were considered positive (Machado et al., 1997).
Conventional PCR assays DNA extraction DNA was extracted from 200 µL of each blood sample using the DNeasy® Blood & Tissue kit (Qiagen®, Valencia, California, USA), in accordance with the manufacturer's instructions. The DNA concentration and absorbance ratio (260/280) nm were measured using a spectrophotometer (Nanodrop, Thermo Scientific, USA). The DNA samples were then stored at -20 °C until their use in PCR assays.

PCR amplification and purification of target B. bovis msa gene fragments
Primers previously designed by Tattiyapong et al. (2014) were used in order to amplify B. bovis msa gene fragments (msa-1, msa-2b and msa-2c). For all the msa genes, a common reverse primer (MSAr) was used based on the conserved nature of the GPI-anchor region of these genes.
Two forward primers in different PCR sets were used for msa-1 amplification. One of them (5-TACTTACCTTTTTAATGACAGCCG-3`) targeted some of the Australian msa-1 sequences (DQ028741, DQO28743, DQ028746 and DQ028747), while the other one (5'-ATGGCTACGTTTGCTCTTTTCATTTCAGC-3') targeted the remaining msa-1 gene sequences retrieved from GenBank. PCR amplification of the target genes was performed as previously described by Tattiyapong et al. (2014), with minor modifications. Briefly, 5 µL of target DNA was used as a template in 25 µL reaction mixtures containing 10X PCR buffer, 1.5 mM of MgCl 2 , 0.8 mM of deoxynucleotide triphosphate mixture (Invitrogen, Carlsbad, California, USA), 1.5 U of Taq DNA polymerase (Invitrogen, Carlsbad, California, USA), 0.4 µM of the forward and reverse primers (Integrated Technologies, Coralville, Iowa, USA) and ultrapure sterile water (Promega, Madison, Wisconsin, USA). The PCR cycling conditions were slightly modified, as follows: an initial denaturation step at 95 °C for 5 min was followed by 45 cycles, each consisting of a denaturing step at 95 °C for 30 s, an annealing step at 56 °C for 1 min for msa-1 and msa-2b, or at 58.1 °C for 1 min for msa-2c, and an extension step at 72 °C for 2 min; and then there was a final elongation step at 72 °C for 5 min. A bovine blood sample positive for B. bovis (GenBank accession number KU522551) and ultrapure sterile water (Promega, Madison, Wisconsin, USA) were used as the positive and negative amplification controls, respectively. The PCR products were subjected to electrophoresis on agarose gel (1%) with ethidium bromide staining and were viewed under UV light (Chemic Doc Imaging System, Bio Rad). All the PCR products showing high intensity of the bands of expected sizes were purified using the Silica Bead DNA gel extraction kit (Fermentas, São Paulo, SP, Brazil).

Phylogenetic analysis
Purified amplified DNA fragments from positive samples were subjected to sequence confirmation in an automated sequencer (ABI Prism 310 genetic analyzer; Applied Biosystems/Perkin Elmer) in both directions, using Sanger's method (Sanger et al., 1977). The electropherogram quality was analyzed using the Phred Phrap software (Ewing et al., 1998), in which only nucleotide sequences above 400 bp in size and Phred quality ≥ 20 were used. Additionally, consensus sequences were obtained through analysis on the sense and antisense sequences using the Phred Phrap software. Initially, msa-2b and msa-2c sequences were individually aligned with sequences available in GenBank using Clustal/W (Thompson et al., 1994). Subsequently, sequences were cut down to the same length (the size of the smallest sequence) and finally were manually adjusted in Bioedit v. 7.0.5.3 (Hall, 1999). Phylogenetic analyses based on the maximum likelihood (ML) criterion were inferred by means of RAxML-HPC BlackBox 7.6.3 (Stamatakis et al., 2008) through the CIPRES Science Gateway. The Akaike information criterion, available in Mega 5.05, was applied to identify the most appropriate model for nucleotide substitution.

Sequence analysis
The genetic diversity of the sequences obtained in the present study was analyzed using the DnaSP v5 software (Librado & Rozas, 2009). The similarity values among the nucleotide sequences were assessed by BLASTn, available in the NCBI GenBank database. The nucleotide sequences of the msa-2b and msa-2c gene fragments were converted into amino acid sequences using the Expasy-Translate tool (Artimo et al. 2012). The percentage similarity between the sequences detected in the present study was calculated using the EMBOSS Needle Pairwise Sequences Alignment software (Needleman & Wunsch, 1970).

PCR targeting B. bovis MSAs
Among the 75 DNA blood samples analyzed from the 15 calves, there were six, 24 and 34 PCR-positive results regarding the msa-1, msa-2b and msa-2c genes, respectively. However, because of the low intensity of bands (possibly due to low parasitemia during the blood sample collection), which precluded obtaining sufficient concentration of amplicons for cloning and DNA sequencing, only three msa-2b and five msa-2c sequences were obtained. The sizes of the sequences obtained for the genes msa-2b and msa-2c were 757-762 and 729-744 bp in length, respectively. Genetic diversity and phylogenetic analyses on the MSA-2b and MSA-2c genes The nucleotide diversity per site (Pi), calculated using the DnaSP v5 software, was higher for the msa-2b sequences (0.1439 -SD = 0.06) than for the msa-2c sequences (0.0600 -SD = 0.01). Pairwise comparison of the deduced MSA-2b and MSA-2c amino acid sequences showed that the degree of similarity of the MSA-2b sequences (61.7-94.2%) was lower than that observed for the MSA-2c sequences (84.3-99.1%) ( Table 2). Likewise, the degree of similarity of the MSA-2b sequences (72.6-96.5%) was lower than that of the MSA-2c sequences (89.2-99.5%) ( Table 3). Also, the alignment of the MSA-2c amino acid sequences showed that they were distinct from each other (Figure 1).
The msa-2b sequences analyzed in the present study clustered into six clades (Figure 2), while the Brazilian msa-2b sequences were positioned in three different clades (clades 1, 3 and 6). Two msa-2b sequences amplified in the present study (KX160813, KX160814) were positioned in clade 1, and were grouped with sequences previously reported from Sri Lanka, Thailand, Vietnam, Philippines and Argentina. One msa-2b sequence (KU522556) from Taiaçu, São Paulo, was positioned in clade 3 and clustered with sequences from USA and Thailand. Lastly, and The phylogenetic analysis based on the msa-2c gene showed four clades. Babesia bovis sequences from Brazil were positioned in two clades (clades 3 and 4) (Figure 3). Only one sequence from Seropédica (KX463633) was positioned in clade 3. Clade 4 grouped four sequences from Seropédica, five from Taiaçu (KU522563, KU522565, KU522567, KX420672 and KX420673), and another three sequences previously reported from different regions of Brazil: HM352734-Brazilian, from southeastern Brazil, HM352731-Brazilian, from central-western Brazil and HM352735-Brazilian, from southern Brazil (Ramos et al., 2012), along with several other sequences previously reported from around the world (Figure 3).  The sequences determined in the present study are shown in boldface letters and followed by the letters RJ and the sequences from Taiaçu, SP followed by the letters Tc. The Phylogenetic tree was inferred by using the maximum likelihood method and Model GTR+G+I. The numbers at the nodes correspond to boostrap values higher than 65% accessed with 1.000 replicates. B. bovis MSA-1 sequences were used as an outgroup. The sequences determined in the present study are shown in boldface letters and followed by the letters RJ, and the sequences from Taiaçu, SP followed by the letters Tc. The Phylogenetic tree was inferred by using the maximum likelihood method and Model GTR+G. The numbers at the nodes correspond to boostrap values higher than 70% accessed with 1.000 replicates. B. bovis MSA-2b sequences were used as an outgroup.

Discussion
Using serological and molecular approaches, the current study assessed natural exposure to bovine babesiosis and the genetic diversity of B. bovis in a cattle herd in Seropédica, Rio de Janeiro.
Serological tests (IFAT and ELISA) have been used as indispensable tools for assessing B. bovis epidemiological status in cattle herds in Brazil. Because cross-serological reactivity can occur between B. bovis and B. bigemina, ELISA tests were performed in order to assess exposure to both Babesia species in the present study. The results showed that the calves sampled had been exposed to both Babesia species during their first 12 months of life. In fact, at six months of age, all the calves sampled had already been exposed to B. bovis and had developed active immunity against this parasite. The rearing system within which these calves were maintained allowed early contact with B. bovis-infected ticks, thus stimulating active immunity against the parasite.
Furthermore, it is now accepted that colostral antibodies are not the main source of protection for calves against bovine babesiosis and that innate immune mechanisms are also involved (Bock et al., 2004). A study carried out in the state of Mato Grosso do Sul, central-western Brazil, showed decreased levels of anti-B. bovis colostral antibodies in calves aged 3-4 months. During this period, clinical babesiosis may occur due to low humoral immune response (Madruga et al., 1984).
Given that the rate of exposure to B. bovis in the present study was above 80%, the herd analyzed here could be characterized as presenting enzootic stability for B. bovis. If over 75% of the calves are exposed to B. bovis, a given herd is considered to be endemically stable (Mahoney & Ross, 1972;Trindade et al., 2010;Costa et al., 2015).
In the present study, the diversity of msa gene fragments in the B. bovis population in cattle in a herd in Seropédica, southeastern Brazil, was analyzed. Although a considerable number of DNA samples were positive, only a few fragments were subjected to sequencing, since most of them yielded faint bands. One possible explanation for these findings may be that the parasitemia levels were low during the sampling. Since our main objective was to evaluate the genetic diversity of B. bovis in the same positive animals at different sampling times, the weak intensity of the bands precluded assessment of whether these msa gene fragments had become modified through B. bovis infection.
In the phylogenetic analyses, the Seropédica msa-2b and msa-2c DNA sequences were positioned into multiple clades. However, it needs to be emphasized that only one msa-2c sequence detected was positioned in a separate cluster. Moreover, although the majority of the msa-2c sequences were grouped in the same cluster, the alignment and the similarity analysis showed that these sequences were distinct. On the other hand, while the B. bovis msa-2b sequences recently obtained from calves in Taiaçu (São Paulo, southeastern Brazil) were detected in multiple clades, msa-2c sequences were positioned in a single clade (Matos et al., 2017).
In a previous study targeting the msa-2c gene, it was reported that the Brazilian B. bovis isolates were genetically conserved (Ramos et al., 2012). However, despite the low number of msa-2c sequences amplified, the current study showed that these sequences were distinct from each other. These finding are in agreement with a previously published study that reported high genetic diversity for this genic locus in the B. bovis strains detected in beef cattle in the state of Mato Grosso do Sul, Brazil (Mendes et al., 2019). Additionally, similar results have been described in cattle from Australia (Berens et al., 2005), the Philippines (Tattiyapong et al., 2014) and Israel (Molad et al., 2014). These results suggest that the genetic diversity of B. bovis strains is different at distinct sampling points.
The analyses on msa-2b gene sequences detected in cattle in Seropédica (RJ) and Taiaçu (SP) (Matos et al., 2017) showed that the B. bovis isolates from both localities were positioned into two main genotypic groups. One genotype occurred on both farms. In the same way as in previous studies (Matos et al., 2017;Mendes et al., 2019), the phylogenetic and similarity analyses showed that this fragment was more heterogenic than the msa-2c gene fragments.
One limitation of the present study was the small number of sequences obtained for diversity analysis. This was mainly because the PCR products were not amenable for cloning. The cloning approach would contribute towards providing a more diverse pool of sequences present in each infected animal. In the future, the use of this strategy will enable better clarification of the diversity of B. bovis msa sequences in Brazil.
In summary, the genetic heterogeneity verified herein may result in antigenic variation in B. bovis strains, thus limiting the use of these antigens as recombinant vaccines for use in controlling bovine babesiosis. In addition, the data reported in this study indicate that all the animals turned out to be serologically positive for B. bovis during their first 12 months age, although the genetic diversity of MSAs was found not to be conserved. These data emphasize the importance of conducting deeper studies on the genetic diversity of B. bovis in Brazil, in order to design diagnostic antigens and vaccines in the future.