MLVA16 Typing of Portuguese Human and Animal Brucella melitensis and Brucella abortus Isolates

To investigate the epidemiological relationship of isolates from different Portuguese geographical regions and to assess the diversity among isolates, the MLVA16Orsay assay (panels 1, 2A and 2B) was performed with a collection of 126 Brucella melitensis (46 human and 80 animal isolates) and 157 B. abortus field isolates, seven vaccine strains and the representative reference strains of each species. The MLVA16Orsay showed a similar high discriminatory power (HGDI 0.972 and 0.902) for both species but panel 1 and 2A markers displayed higher diversity (HGDI 0.693) in B. abortus compared to B. melitensis isolates (HGDI 0.342). The B. melitensis population belong to the “Americas” (17%) and “East Mediterranean” (83%) groups. No isolate belonged to the “West Mediterranean” group. Eighty-five percent of the human isolates (39 in 46) fit in the “East-Mediterranean” group where a single lineage known as MLVA11 genotype 116 is responsible for the vast majority of Brucella infections in humans. B. abortus isolates formed a consistent group with bv1 and bv3 isolates in different clusters. Four MLVA11 genotypes were observed for the first time in isolates from S. Jorge and Terceira islands from Azores. From the collection of isolates analysed in this study we conclude that MLVA16Orsay provided a clear view of Brucella spp. population, confirming epidemiological linkage in outbreak investigations. In particular, it suggests recent and ongoing colonisation of Portugal with one B. melitensis lineage usually associated with East Mediterranean countries.

Brucellosis in humans is largely dependent on the animal reservoir and can occur by direct contact with infected animals or through consumption of contaminated dairy products. B. melitensis and B. abortus are highly pathogenic and a frequent cause of human infection. B. abortus is the most disseminated species worldwide and B. melitensis produces the most severe infection in humans and predominates in the Mediterranean Basin [7].
Although several countries have implemented eradication programs which resulted in the reduction or elimination of the disease, the infection remains enzootic in many regions of the world [8]. In the scope of eradication programmes Brucella species and biovars are currently typed by microbiological, serological and chemical tests. For many years, the phenotypic characterization was one of the few tools available to provide epidemiological data but does not allow tracing back the source of origin, hindering the application of measures to prevent additional spread and disease elimination from the primary source [9,10]. In the countries where brucellosis has been eradicated or strictly controlled, continued surveillance is essential to prevent the re-emergence of the disease [8].
Genetic typing of the bacterium is often used for molecular epidemiological purposes and can be useful to observe population structures. In the last decades, attempts made to develop epidemiological markers for Brucella have shown limited success.
Most of the tests lacked reproducibility and showed a moderate capability to differentiate strains [9,11]. The genetic conservation within Brucella has resulted in past difficulties in establishing the true relationships between some classical Brucella species and biovars and in identifying molecular markers for certain groups. In recent years, the sequencing of multiple genetic loci (Multilocus Sequence Typing; MLST) has rapidly gained acceptance as a tool for the characterization of microbial populations. MLST has provided the first genus-wide and robust view of the population structure of this highly clonal genus [10]. Over the past few years, polymorphic tandem repeat loci have been identified by the analysis of published genome sequences of the strains B. melitensis biovar 1 16 M, B. suis biovar 1 1330 and B. abortus biovar 1 9-941 [12][13][14][15][16][17]. One of the applications based on these markers is the multiple-locus variable number tandem repeats analysis (MLVA) with a set of 16 repeat loci (MLVA16 Orsay ) [9,13]. MLVA has a much higher discriminatory power than MLST and produces highly relevant clustering. A number of reference strains have been typed by both MLVA and MLST [6]. These linking strains make it possible to predict with high confidence a MLST type from the MLVA type for most new isolates. Only those strains showing weak clustering to such a linking strain will need to be analysed by the more demanding MLST or draft whole genome sequencing analysis to become new linking strains.
In Portugal, brucellosis is an endemic disease where B. melitensis biovars 1 and 3 and B. abortus biovars 1 and 3 are the prevailing species (LNIV, National Laboratory for Veterinary Research, NRL for Brucellosis, unpublished data). Human infections due to B. melitensis pose a significant public health problem with 80 confirmed cases in 2009 (Portuguese Health Authority, unpublished data). The archipelago of Azores is officially free of B. melitensis but B. abortus is present in three of the nine islands, S. Miguel, S. Jorge and Terceira. Animals move between the archipelago islands and from the islands to mainland Portugal [18]. In this report we present the results of a study performed by Brucella MLVA16 Orsay on 283 human and animal field isolates of B. melitensis and B. abortus. The aims of this study were to evaluate MLVA16 Orsay and its MLVA11 subset in order to assess the diversity among strains for epidemiological purposes in human brucellosis and to estimate the epidemiological relationship of isolates from different geographical origins in Portugal.

MLVA 16 typing and clustering of B. melitensis and B. abortus populations
The MLVA16 Orsay assay (panels 1, 2A and 2B) was performed with a collection of 126 B. melitensis (46 human and 80 animals) and 157 B. abortus isolates, 7 vaccine strains and the representative reference strains of each species (Table 1, Figure 1). All loci could be amplified in all samples, except for three B. abortus isolates, one negative for bruce09, another for bruce16 and a third for bruce30. The combination of the three panels (1, 2A and 2B) provided a similar high discriminatory power in the B. melitensis and B. abortus populations (HGDI of 0.972 and 0.902, respectively) (  Figure 2. One hundred fifty seven bovine B. abortus field strains were studied. MLVA11 discriminated seven genotypes (not including the reference and vaccine strains), four of which were not previously observed (181 to 184). MLVA16 Orsay distinguished a total of 40 genotypes with 20 singleton genotypes ( Figure 3). The isolates fell into two clusters, A1 and A2 respectively. Cluster A1 corresponded to 29 MLVA16 Orsay genotypes (1 to 29 containing the bv1 (544, 2308 and 99 W), bv2 (86/8/59) and bv4 (292) strains, the vaccine strains S19 and RB51 and all B. abortus bv1 field isolates. The vast majority of A1 isolates belonged to MLVA11 genotype 82, 181, 182 or 183. MLVA11 genotype 82 was observed only in mainland Portugal and S. Miguel island. MLVA11 genotype 181 seems to be specific to Terceira. MLVA11 genotype 182 was observed in mainland Portugal, S. Miguel, Terceira and S. Jorge islands. Cluster A2 grouped all isolates from bv3 (genotypes 30 to 40). These isolates were from mainland and shared MLVA11 genotype 72. B. abortus bv1 and bv3 isolates were clearly differentiated by bruce11 (4 and 3 repeats respectively) and bruce55 (3 and 1 repeats respectively).

Epidemiological relationship between human and animal B. melitensis isolates
Forty six B. melitensis strains isolated from Portuguese patients and 80 B. melitensis strains from cattle, sheep and goats from the different Portuguese regions ( Figure 1) were studied. The human isolates represented 21 MLVA16 Orsay genotypes, thirteen of which included only human isolates and eight were common to human and animals ( Figure 2). Seven human isolates belonged to cluster M1 and 39 to M2. The M1 isolates were associated with three MLVA11 genotypes, 130, 147, 186. All the other 39 human isolates showed MLVA11 genotype 116. Cluster analysis grouped three isolates, with MLVA16 Orsay genotype 10, obtained from three different cultures from one patient (strains LNIV-860Bm3-09, LNIV-861Bm3-09 and LNIV-862Bm3-09). In some cases an association with isolation year and geographic origins was observed. Sixteen human isolates in cluster M2 shared the genotype 53 (10 isolates were from an outbreak in Beira Litoral (BL) in 2005 and the other six were isolated in France from Portuguese emigrants, in the same year). Additionally, three human strains (one from BL and two isolated in France) clustered as genotype 52, which differs from genotype 53 by the loss of two repeat units at locus bruce16. Six strains (one human, four goat and one bovine isolates) from a recent outbreak (2011) in BL, as well as a sheep strain isolated in 2005 in Algarve (ALG) displayed the genotype 33. The genotype 60 clustered a human strain isolated in France in 2006 and the human strains LNIV-917Bm3-11 and LNIV-923Bm3-11, isolated from the same patient before and three months after treatment with rifampicin and doxycycline.
A minimum spanning tree (MST) analysis was performed using all isolates from this study and typing data from more than one thousand isolates previously published in Brucella2010 database hosted at http://mlva.u-psud.fr and by other authors [19][20][21][22].
Epidemiological relationship between B. abortus isolates from different geographical origins A total of 157 bovine B. abortus isolates from the seven Portuguese mainland regions ( Figure 1) and S. Miguel, S. Jorge and Terceira islands from Azores were evaluated based on MLVA16 Orsay data ( Figure 4). Isolates from Ribatejo e Oeste (RO), as well as most of the isolates from Alentejo (A) clustered in genotype 12, irrespective of the isolation year. Further, in Beira Interior (BI) two strains (LNIV-243Ba1-05 and LNIV-246Ba1-05) from an outbreak in 2005 and four strains (LNIV-DNA181Ba1-08, LNIV-DNA182Ba1-08, LNIV-DNA183Ba1-08 and LNIV-DNA184Ba1-08), isolated from a neighbouring herd in 2008 displayed the genotype 10. Strains from Azores clustered in 18 genotypes. S. Miguel isolates were grouped in 11 genotypes which were also found in isolates from mainland ( Figure 3) and S. Jorge isolates clustered apart from most of the other Azores and mainland strains (genotypes 4, 5 and 6). One strain from Terceira (LNIV-399Ba1-06, genotype 7) showed 13 repeats in bruce12, one more than the other isolates. Eight milk isolates, also from Terceira, presented 20 repeats in bruce19 (genotypes 18, 19 and 20), one less than the other strains isolated from tissues. Additionally, two isolates of RB51 (1 from S. Miguel and 1 from Trás-os-Montes) displayed the genotype 22, like the RB51 vaccine strain and its parental strain 2308. These strains had three unit repeats in bruce43, one more compared to all the other bv1 isolates.
A MST analysis using all B. abortus isolates was performed with previously published typing data recovered from the compilation in Brucella2010 database ( Figure 5).

Discussion
Control of brucellosis, particularly in the final stages of an eradication program, requires an accurate surveillance and highly discriminatory methods to characterize an outbreak strain and determine the infection source and transmission routes. Correct application of bacterial typing will increase the efficacy of control measures aimed to contain or interrupt the outbreak. For epidemiological trace-back purposes, the array of phenotypic tests used in the classification of Brucella species and biovars is not accurate enough particularly in geographical regions with a predominance of a small number of biovars, like Portugal, where most isolates are B. abortus bv1 and bv3 and B. melitensis bv1 and bv3.
MLVA16 Orsay and its MLVA11 subset were used to assess the diversity among strains for epidemiological purposes in human brucellosis and to estimate the epidemiological relationship between B. melitensis (n = 126) and B. abortus (n = 157) isolates collected during the last 11 years from different geographical areas in Portugal. The combination of the three panels provided a similar high discriminatory power for both species ( MLVA11, B. melitensis bv3 isolates are separated from bv3 reference strain (Ether), differing from bv2 reference strain (63/ 9) in one locus (bruce12). Similar results were obtained by other authors and methods [1,23], questioning the phylogenetic robustness of some biovars established by classical microbiological typing. MLVA16 Orsay was able to identify the epidemiological relationship between strains. Genotype 53, assigned to 16 human isolates (10 from an outbreak occurred in 2005 in Beira Litoral and six isolated from Portuguese emigrants, in the same year, in France) probably represents an outbreak. Additionally, three human isolates (one from Beira Litoral and two from France, both isolated in 2005) clustered together as genotype 52, which differs from genotype 53 by the loss of two repeat units at the highly variable locus bruce16, suggesting a close relatedness [24,25]. The six strains from the recent outbreak (2011) in Beira Litoral (one human, four goats and one bovine isolate) share genotype 33, clearly showing the epidemiological relation between isolates from human and animal origin. A sheep strain isolated in 2005 in Algarve displayed the same genotype 33. These strains had the same MLVA16 Orsay profile as a strain isolated in France in 1996 (presumably from a Portuguese patient) and described by Al Dahouk et al. [9], suggesting a common source of infection and confirming the persistence of genotypes in certain areas for long periods. A previous isolate (LNIV-917Bm3-11) and the strain LNIV-923Bm3-11, recovered from the same patient after antibiotic treatment presented genotype 60. These results confirm the relapse even after three months of treatment with rifampicin and doxycycline. The antibiotic minimal inhibitory concentration (MIC) was determined by agar dilution method, as recommended by the Clinical Laboratory Standards Institute [26] and both strains exhibited the same susceptibility to rifampicin and doxycycline (MIC,1 mg/mL).
B. abortus strains formed a consistent group with bv1 and bv3 isolates completely separated in clusters A1 and A2. As noticed by other authors [27], B. abortus bv3 isolates are found in a separate group from the biovar reference strain (Tulya). When comparing MLVA11 results with typing data published in Brucella2010 database, seven genotypes were defined. Genotypes 181 to 184 were observed for the first time in isolates from S. Jorge and Terceira islands. All B. abortus bv3 isolates presented the same MLVA11 genotype 72 previously described by Le Flèche et al. [13].
The epidemiological relationship between isolates from different Portuguese geographical origins was evaluated by MLVA16 Orsay . S. Jorge isolates clustered in three genotypes (4, 5 and 6), separately from most of the other strains from Azores and mainland, which differ from each other only in one repeat in bruce16 or in bruce30. The one repeat difference in the highly variable microsatellites from panel 2B reveals them as the same strain, presumably with a recent common ancestor. These genotypes were also found in two isolates from S. Miguel, one from Terceira and two from Alentejo showing the importance to monitor national movement of strains by MLVA16 Orsay . Further, Bruce42 distinguished these isolates from others from Azores and mainland. Seventeen out of 34 B. abortus isolates from S. Miguel  clustered in two genotypes (12 and 14) along with strains from Alentejo, Ribatejo e Oeste and Beira Interior. These strains showing identical genotypic traits are clones with a common origin, revealing linkage of genotypes with the geographic origin [25]. These results also confirm the report by Batista and Nunes [18] indicating that there is a high cattle movement out of S. Miguel to mainland. In Beira Interior, two strains from an outbreak in 2005 and four strains isolated from a neighbouring herd in 2008 displayed the same genotype 10, again confirming the epidemiological linkage identified in outbreak investigations.
Considering the vaccine strains, three Rev1 field isolates from two sheep and one goat (genotypes 19 and 20) shared panel 1 genotype 47 described by García-Yoldi et al. [28] and differed from strain 16 M in bruce07, bruce09 and bruce18. Also two isolates from Azores without island identification, isolated in 2001, presented the same genotype as the vaccine strain S19 (genotype 17), which is not unexpected, since vaccination with S19 was used in Azores (S. Miguel and Terceira) until 1984 [29]. Two RB51 field isolates (one from S. Miguel and one from Trás-os-Montes) displayed the same genotype as the RB51 vaccine strain and its parental strain 2308 (genotype 22), which is in accordance with previous observations [13,27]. These results confirm the genetic stability of Rev1, S19 and RB51 vaccine strains [13,28]. Also, the RB51 strains differed from all other B. abortus isolates in bruce43.
The set of 16 loci used in MLVA16 Orsay showed to be extremely discriminant and highly efficient in distinguishing strains within a local outbreak. Though able to correctly differentiate regional geographically related B. abortus isolates from Azores, it could not distinguish B. melitensis strains from closely related geographical regions from mainland. In fact, data on animal movement in Portugal show that all mainland regions are connected.
MLVA16 Orsay provides a clear view of the Brucella population, confirming epidemiological linkage in outbreak investigations. In areas with low prevalence of brucellosis, the epidemiological answers given by MLVA16 Orsay are very important to prevent the expansion and achieve the eradication of the disease. In particular, our results suggest a recent and ongoing colonisation of Portugal with one lineage (MLVA11 genotype 116) usually associated with East Mediterranean countries. Draft whole genome sequencing of selected isolates may help decipher sub lineages within Portugal lineage 116 and provide estimates to date this contamination event.
The data produced in this work can be queried in the Brucella MLVA database release (starting from the Brucella2012 release at http://mlva.u-psud.fr).

Brucella isolates
Two hundred and eighty-three new isolates were investigated: 126 B. melitensis and 157 B. abortus isolates. Fifteen reference strains were included as controls: the type strains of the three biovars of B. melitensis (biovars 1, 2 and 3) and biovars 1-6 and 9 of B. abortus reference strains; B. abortus biovar 1 strains 2308 and 99 W; S19 and RB51 vaccine strains as well as three Rev1 and two RB51 vaccine field isolates. All isolates and strains used in this study are summarized in Table 1 and the geographical origin of isolates is shown in Figure 1. All isolates were biotyped according to standard bacteriological procedures [30] [6]. Amplifications were performed in a MyCycler thermal Cycler (Biorad, France). An initial denaturation step at 96uC for 5 min was followed by 30 cycles of denaturation at 96uC for 30 s, primer annealing at 60uC for 30 s and elongation at 70uC for 1 min. The final extension step was performed at 70uC for 5 min. Five microliter of amplification products were loaded on a 3% standard agarose gel to analyse panel 2A and 2B loci (tandem repeats with a unit length shorter than 8 bp) and on a 2% standard agarose gel for panel 1 loci (tandem repeats with a unit length larger than 10 bp) and run under a voltage of 8 V/cm for 120 min. Depending on the tandem repeat unit length, a 20 bp (20 bp PCR Molecular Ruler, Biorad, France) or 100 bp ladder was used as molecular size markers (100 bp DNA ladder, Invitrogen).

Data analysis
The 16 loci have been classified in three panels, named panel 1, composed of 8 minisatellite loci (bruce06, bruce08, bruce11, bruce12, bruce42, bruce43, bruce45 and bruce55), panel 2A (bruce18, bruce19 and bruce21) and panel 2B (bruce04, bruce07, bruce09, bruce16 and bruce30) composed of three and five microsatellite markers, respectively. The total number of repeat units at each locus was determined by the correlation with the amplicon size according to previously published tables (Le Flèche et al. 2006 version 3.3: http://mlva.u-psud.fr Brucella support web site for MLVA typing). The copy number for each locus and each panel was managed as a character dataset using BioNumerics version 6.6 (Applied Maths, Belgium). Cluster analysis was based melitensis isolates (excluding vaccine strains) using the MLVA11 data. Colour codes are associated with the main B. melitensis MLVA groups. The analysis was performed using all isolates from this study and typing data from more than one thousand isolates previously published [19][20][21][22], Brucella2010 database hosted at http://mlva.u-psud. fr]. Portuguese isolates are shown in yellow. The B. melitensis population is included in ''American'' (17%) and ''East-Mediterranean'' (83%) groups. Thirty nine human isolates (85%; 39 out of 46) are MLVA11 genotype 116 usually associated with the ''East-Mediterranean'' B. melitensis strains and seven isolates (15%) are included in the ''American'' group. The numbers represent the MLVA11 genotypes found in this study. The size of the shapes indicates the number of strains described in each genotype. Each of the circles showing yellow and dark blue or red colors included the Portuguese genotype (yellow) and the genotypes found in ''East-Mediterranean'' (dark blue) or ''American'' (red) B. melitensis isolates. doi:10.1371/journal.pone.0042514.g003  on the categorical coefficient and unweighted pair group method using arithmetic averages (UPGMA). Three different character data sets were defined (panel 1, 2A and 2B markers) as combined sets using the composite data set tool provided by BioNumerics. Different weightings were given to each panel according to Al Dahouk et al. [9]: panel 1 markers got an individual weight of 2 (total weight 16), panel 2A markers got an individual weight of 1 (total weight 3) and markers of panel 2B got an individual weight of 0.2 (total weight 1). The genetic diversity of each locus and respective panels was determined using the Hunter & Gaston diversity index (HGDI) [31] via the online tool V-DICE available at the HPA website (http://www.hpa-bioinformatics.org.uk/cgibin/DICI/DICI.pl).