Molecular characterization of B. anthracis isolates from the anthrax outbreak among cattle in Karnataka, India
BMC Microbiology volume 20, Article number: 232 (2020)
Anthrax, a zoonotic disease is caused by the Gram positive bacterium Bacillus anthracis. During January 2013, an anthrax outbreak among cattle was reported in Gundlupet Taluk, neighboring Bandipur National Park and tiger reserve, India. The present study aims at the molecular identification and characterization of 12 B. anthracis isolates from this outbreak by 16S rRNA gene sequencing, screening B. anthracis specific prophages and chromosomal markers, protective antigen (pag) gene and canonical single nucleotide polymorphism (canSNP) analysis to subtype the isolates into one of the twelve globally identified clonal sub-lineages of B. anthracis.
These isolates had identical 16S rDNA nucleotide sequences with B. anthracis specific dual peaks showing mixed base pair R (G/A) at position 1139 with visual inspection while the automated basecaller software indicated a G. Alternatively the nucleotide A at 1146 position was indicative of the 16S rDNA type 7. Multiple sequence alignment with additional 170 (16S rDNA) sequences of B. cereus sensu lato group from GenBank database revealed 28 new 16S types in addition to eleven 16S types reported earlier. The twelve B. anthracis isolates were found to harbor the four B. anthracis specific prophages (lambdaBa01, lambdaBa02, lambdaBa03, and lambdaBa04) along with its four specific loci markers (dhp 61.183, dhp 77.002, dhp 73.019, and dhp 73.017). The pag gene sequencing identified the isolates as protective antigen (PA) genotype I with phenylalanine-proline-alanine phenotype (FPA phenotype). However, sequence clustering with additional 34 pag sequences from GenBank revealed two additional missense mutations at nucleotide positions 196 bp and 869 bp of the 2294 bp pag sequence among the 5 B. cereus strains with pXO1 like plasmids. The canSNP analysis showed that the isolates belong to A.Br.Aust94 sub-lineage that is distributed geographically in countries of Asia, Africa, Europe and Australia.
The analysis of 16S rDNA sequences reiterated the earlier findings that visual inspection of electropherogram for position 1139 having nucleotide R could be used for B. anthracis identification and not the consensus sequence from base caller. The canSNP results indicated that the anthrax outbreak among cattle was caused by B. anthracis of A.Br.Aust94 sub-lineage.
Bacillus anthracis, the Gram positive, spore forming bacteria is an etiological agent of anthrax and belongs to a monophyletic lineage of Bacillus cereus sensu lato group that currently includes 21 published species . B. anthracis, being genetically closely related to B. cereus and B. thuringiensis is believed to have diverged from its B. cereus ancestor due to the evolutionary acquisition of two virulence plasmids, pXO1 and pXO2  and often considered as a single species . The protective antigen (pag) and capsule (cap) genes, located on these plasmids are used as molecular markers for the routine identification of B. anthracis. However, occurrence of B. cereus and B. thuringiensis strains with either or both of these plasmids pose challenges in unambiguous identification of the species . The 16S rDNA gene sequences of B. cereus group exhibit a relatively high level of sequence similarity (> 99%) but differ by 13 nucleotide positions that were utilized to identify 16S rDNA types associated with B. anthracis . The 16S rDNA gene sequence classified most B. anthracis strains in 16S type 6; however, some were reported in type 7 with B. cereus strains. Later it was reported that the mixed base pair R (G/A) at position 1139 in B. anthracis 16S rRNA gene owing to the presence of multiple rRNA operons in the genome could be used for its differentiation from other closely related B. cereus group species . Further, B. anthracis specific DNA signature sequences (loci A, C, D, and E) and four prophages (lambdaBa01, lambdaBa02, lambdaBa03, lambdaBa04) located in B. anthracis chromosome are also used for definitive discrimination of B. anthracis from other B. cereus group strains [7, 8]. A defined set of 13 canSNPs from different loci in the genome has been used for molecular typing to determine the phylogenetic position and distinguish geographically separated B. anthracis isolates owing to their high level of genetic resolution and discriminatory power . The distribution of these canSNPs identified 12 sub-lineages from the B. anthracis clonal lineages A, B and C representing different geographical locations .
India hosts the world’s largest livestock population where anthrax is endemic with the highest incidence from the southern part of the country [10, 11]. An epidemiological study by National Animal Disease Referral Expert System (NADRES), India describes anthrax as the third important zoonotic disease in the subcontinent . In these endemic regions, anthrax has also been reported in the livestock handlers . To address the same concern, the National Centre for Disease Control (NCDC), India have been taking efforts for strengthening surveillance, early diagnosis and, effective timely containment of this zoonotic disease . In addition, NCDC along with World Health Organisation (WHO) has also developed training programs for veterinary and medical professionals to control the spread of anthrax . The anthrax vaccine (anthrax spore vaccine) produced in India has been subsidized by the government and is distributed at zero cost to the farmers in the endemic regions through veterinary officials .
Recent countrywide surveillance jointly conducted by National Centre for Disease Control (NCDC), Indian Council of Agricultural Research - National Institute of Veterinary Epidemiology and Disease Informatics (ICAR-NIVEDI) and Centers for Disease Control and Prevention (CDC) confirmed anthrax outbreaks in the regions of Odisha and Karnataka . Anthrax associated with infected livestock has later progressed into human cutaneous anthrax due to close proximity with animals. Earlier, human anthrax has been reported due to vector-borne transmission from infected cattle, goat, sheep and buffalo due to handling of infected livestock, slaughtering and even after consumption in few reported cases [18,19,20,21]. These data confirm the endemicity of anthrax in the southern states of the subcontinent, that is, Andhra Pradesh, Karnataka and Tamil Nadu . Previously, twelve B. anthracis isolates harboring both pag and cap genes were isolated from an anthrax outbreak among cattle in Gundlupet Taluk neighboring Bandipur National Park, and tiger reserve, India during January 2013 . The present study aimed to identify and characterize these 12 isolates using molecular tools that included PCR based screening of B. anthracis specific prophages and chromosomal markers, sequencing of 16S rDNA and pag genes. We also employed canSNP analysis technique to delineate the clonal lineages for the set of B. anthracis isolates.
Results and discussion
16S rDNA sequence analysis
B. anthracis genome is characterized by the presence of multiple rRNA operons . Two specific rDNA positions, viz., 1139 and 1148 with nucleotide mismatches among the eleven 16S rRNA operons have been reported earlier for the species identification [5, 6]. Base identification by visual inspection of electropherograms revealed dual peaks representing R (A or G) for the former and W (T or A) for the latter at these nucleotide mismatches (Table 1) . These dual peaks specific to B. anthracis has been used for their identification when the consensus 16S rDNA gene sequence could not be used for differentiation of the closely related B. cereus sensu lato strains i.e., B. cereus, B. anthracis and B. thuringiensis with high level (> 99%) of sequence similarity . The W at position 1148 (1146 in Sacchi et al. ) initially found to be specific for B. anthracis and grouping them in 16S rDNA type 6 . The same position was found among B. cereus and B. thuringiensis strains, whereas the dual peak R at 1139 was found only in B. anthracis strains making this position more reliable for B. anthracis identification . All the 12 B. anthracis strains characterized in this study had identical 16S rDNA nucleotide sequences with A at 1146 position (Table 1) indicative of 16S type 7 shared by B. anthracis, B. cereus and B. thuringiensis strains [5, 6]. Alternatively, the isolates had the B. anthracis specific dual peaks showing mixed base pair R (G/A) at position 1139 with visual inspection (Additional file 1) while the automated basecaller software indicated a G at this position (Table 2) as reported earlier . Conversely, analysis of consensus 16S rDNA sequences ranging from nucleotide positions 352 to 1007 bp is short of the significant number of nucleotides at both 5′ and 3′ regions that otherwise held multiple polymorphic sites, could not discriminate B. anthracis from B. cereus and B. thuringiensis . We further extended our analysis and performed multiple sequence alignment (MSA) with an additional 167 consensus 16S rDNA sequences and 3 sequences with base identification by visual inspection of B. cereus sensu lato to identify the microheterogeneity in the 16S rDNA gene. This included 47 strains including B. anthracis, B. cereus B. thuringiensis and B. mycoides which were representative of eleven 16S types reported by Sacchi et al.  scheme. Also, 33 Bacillus strains (15 inclusivity and 18 exclusivity species) mentioned in the Association of Analytical Communities (AOAC) International panel were included in the analysis . Additionally, 90 strains including B. anthracis, B. cereus, B. thuringiensis and, B. mycoides for which 16S rDNA sequences were available from the GenBank database were included in the MSA. A total of 104 variant sites (V) which included 13 SNP positions reported earlier , 27 new parsimony-informative sites (Pi) and 64 new singletons (S) were observed in the 1402 bp 16S rDNA sequence examined. According to the 13 SNP positions , 50 of the additional 170 strains included in the study exhibited 19 new 16S types not described earlier . Also, it was observed that 28 B. anthracis strains and 13 B. thuringiensis strains from the GenBank database were of 16S types 13 and 7 respectively, otherwise seen exclusively among B. cereus strains reported earlier  (Additional file 2). The position 1139 (1137 in Sacchi et al. ) was uniformly identified with G in most of the B. cereus sensu lato strains, with few exceptions. The exceptions included the twelve isolates from the present study and three B. anthracis strains, Pasteur, A0248 and Sterne that had R (G/A) derived from dual peaks for the mixed base pair observed on visual inspection of the electropherograms  (Table 2). We speculate that the use of consensus 16S rDNA sequence from automated base caller software could have possibly limited the detection of mixed base pair R (G/A) at position 1139 .
Alternatively, when we considered the 27 new Pi sites (positions 14 to 40) along with the 13 previously identified SNP positions for characterizing the population, 28 new 16S types could be identified among the 182 strains (Additional file 2). However, even with the inclusion of 27 additional Pi sites/ SNPs, all the 47 B. cereus sensu lato strains from Sacchi et al.  scheme, 51 strains from GenBank population and the twelve isolates in the present study held on to eleven 16S types reported earlier . About, 19 B. anthracis strains and 10 B. cereus strains were distributed in eight 16S types each viz., 16S types A, C, G, H, I, M, R, T and 16S types J, L, N, U, V, W, X and Y respectively (Additional file 2). Five B. thuringiensis strains were distributed in four 16S types B, D, Z and AA while 5 B. mycoides strains in two 16S types Q and BB (Additional file 2). Interestingly, 8 B. cereus strains and 19 B. thuringiensis strains shared four 16S types E, F, O and P, whereas, 4 B. cereus strains and 2 B. anthracis strains shared two common 16S types, K and S respectively (Additional file 2). Inclusion of additional B. cereus sensu lato strains revealed newer polymorphic nucleotide positions resulting in new 16S rDNA types. Even though 16S rDNA types of B. anthracis strains could solely be identified (16S types A, C, G, H, I, M, R and T), few of them clustered along with B. cereus (16S types E, F, O and P) also. This could result in ambiguous identification of any B. cereus sensu lato strains clustering in 16S rDNA types with multiple species.
Also, the Neighbor joining tree constructed with bootstrap values calculated from 1000 replications grouped the entire Bacillus population in the present study into four major clusters with an average pairwise distance of 0.054–0.182 (Fig. 1). The polymorphisms at nucleotide positions 19 to 35 resulted in Cluster I, an out-group with B. anthracis strains PAK.1, Canadian Bison, K3 and B. cereus strain D17 with 16S types A, G, M, and N respectively. This cluster also included B. thurnigiensis strain Al-Hakam and B. cereus strain NC7401 for which position 13  could not be elucidated due to their shorter 16S rDNA sequences available in the GenBank. The Cluster II mainly consisted of B. mycoides strains and comprised of two new 16S types Q and BB along with 16S type 1 (B. mycoides ATCC 10206) identified earlier . Cluster III comprised majorly B. cereus and B. thuringiensis strains representing the 8 new 16S types of B. cereus, 4 new 16S types of B. thuringiensis and two new 16S types shared among the two species. Cluster IV consisted of B. anthracis strains along with B. cereus and B. thuringiensis strains that shared B. anthracis plasmids and few other strains of these species. This cluster represented 8 new 16S types of B. anthracis and 2 new 16S types shared among B. cereus and B. anthracis. The twelve isolates in the present study also grouped in this cluster along with other B. anthracis strains with 16S type 7 and exhibited 100% similarity in the MSA with no gaps or mismatches in the entire 16S rDNA gene analyzed.
Screening of B. anthracis specific loci and prophages
B. anthracis strains are routinely identified based on the presence of plasmid-borne virulence marker genes (pag and cap) . However, the existence of either or both of these plasmids or their homologs in certain B. cereus sensu lato strains poses challenges in the confirmatory identification of B. anthracis. Assays targeting polymorphism present in chromosome borne housekeeping genes such as rpoB [26, 27], gyrA , gyrB , sspE  and sap  have been used to distinguish B. anthracis from other species of B. cereus senu stricto group. Four putative lambdoid prophages (designated lambdaBa01, lambdaBa02, lambdaBa03 and lambdaBa04) have been reported to be integrated specifically in the 3% of 5.2 Mb B. anthracis chromosome . They are conserved across more than 300 geographically and temporally divergent B. anthracis strains .
The presence of all the four prophages in the chromosome is a unique feature used for definitive discrimination of B. anthracis from certain B. cereus sensu stricto strains that contain some of these lambdoid prophages but with few sequence homologies with B. anthracis prophage genes . B. cereus sensu lato strains such as B. cereus E33L and B. weihenstephanensis KBAB4 have been reported to contain B. anthracis lambda01 prophage in their chromosome  and hence can be determined as non-anthrax causing Bacillus species. Also the BLAST results of these four lambdoid prophages indicated their presence in all the B. anthracis isolates present in the NCBI database and no other Bacillus species were found to contain all of the prophages in their genome (data not shown). However, the B. anthracis isolates in the present study were confirmed for the presence of all the four prophage regions yielding expected amplicon sizes (221, 233, 189 and 157 bp) for lambdaBa01, lambdaBa02, lambdaBa03 and lambdaBa04 markers respectively (Additional file 3A).
Similarly, the presence of four unique “DNA signatures” A, C, D, and E (dhp 61.183, dhp 77.002, dhp 73.019, and dhp 73.017 respectively) within the B. anthracis genome is considered as a powerful tool to differentiate B. anthracis from non-anthrax causing Bacillus species . In the present study, we found that the DNA of all the twelve B. anthracis isolates consistently amplified all the four B. anthracis specific gene markers, dhp 61.183, dhp 77.002, dhp 73.019 and dhp 73.017 with the predicted amplicon size (163 bp, 133 bp, 196 bp, and 241 bp respectively) (Additional file 3B) with no variations within the isolates. These findings corroborated the identification of the B. anthracis isolates in the present study.
Protective antigen gene (pag) sequence analysis and genotyping
Protective antigen (PA), the binding component of anthrax toxin is encoded by pag gene, located on pXO1 plasmid of B. anthracis . B. anthracis strains have been classified into six PA genotypes based on the polymorphism in 7 nucleotide positions observed in the 2294 bp pag gene (GI: NC_001496.1) from 26 diverse B. anthracis strains . These seven nucleotide positions included four synonymous and three missense mutations resulting in four different phenotypes. Therefore, to characterize the PA genotype and phenotype, the pag amplicons from all the 12 isolates were sequenced and analyzed. The PA genotyping helped in elucidating any variation in the gene resulting from gene alteration or gene engineering and hence, aided in understanding the evolution of our B. anthracis isolates . Also, the pag sequences from recently evolved B. anthracis strains and B. anthracis like strains, which were not included in the Price et al.  scheme, were analyzed to evaluate the phylogenetic relationship of our twelve isolates with the global B. anthracis strains in the present scenario. The pag sequences of the twelve isolates, 26 B. anthracis strains from the earlier scheme , 15 B. anthracis strains from the inclusivity panel by AOAC , 5 B. cereus strains with pXO1 like plasmids  and 14 B. anthracis strains from the GenBank were analyzed together for their phylogenetic comparisons. B. anthracis isolates of the present study were found to be of pag genotype I and FPA phenotype (Table 3) representing the Sterne-Ames diversity group as described previously . This pag genotype I have been reported previously in 135 B. anthracis isolates from bioterrorism associated anthrax outbreak in the United States during 2001 , human anthrax cases from Hong Kong  and Lianyungang regions , and cattle related anthrax outbreak from Central Java and Yogyakarta regions of Indonesia . No new pag genotypes other than that reported by Price et al.  scheme were identified among B. anthracis even on including additional strains from recent outbreaks.
However, the sequence clustering on the inclusion of 34 additional pag sequences along with the pag sequence data of 12 isolates revealed two additional polymorphic nucleotide positions among the 5 B. cereus strains with pXO1 like plasmids causing missense mutations at nucleotide positions 196 bp and 869 bp of the 2294-bp pag sequence from Sterne strain (GI: NC_001496.1) (Table 4). No variations among B. anthracis strains were found at these two additional positions of differences. The pag is divided into 4 domains where each domain plays a critical role in toxin action . The transition mutation (T to C) at the nucleotide position 196 resulted in a serine to proline change near the amino-terminal fragment (PA20) of PA domain 1 (amino acids (aa) 1 to 258) that is cleaved due to furin protease before the binding of the PA63 monomer to the host cell receptor . Alternatively, the isoleucine to serine change in aa 289, a part of PA domain 2 (aa 259 to 487) brought about by the transversion mutation (T to G) at 869 bp polymorphic site, could possibly play a role in altering the binding of PA63 monomer and needs further elucidation .
The Neighbour joining tree was constructed based on nucleotide variations at 7 polymorphic sites studied with bootstrap values calculated from 1000 replications grouped the entire B. anthracis group of strains population in the present study into six major clusters with an average pairwise distance of 0.168–0.973 (Table 5). Clusters I, II, III, IV, V and VI were represented by pag types I, II, III, IV, V and VI respectively (Fig. 2). The twelve B. anthracis isolates from the present study grouped with pag type I along with B. anthracis strain Ames Ancestor and Sterne. The five B. cereus strains with pag gene analyzed in the present study were grouped in Cluster V of pag type V excluding the two additional positions of mutations in consideration.
SNPs have played a major role as a genetic marker for detecting and subtyping bacterial pathogens . A small number of SNPs have been used to define genetic groups efficiently even among largely monophyletic clade like B. anthracis. The phylogenetic positions of B. anthracis isolates were evaluated using the 13 canSNPs scheme as proposed earlier (Table 6) . All the B. anthracis isolates from the study were assigned to sub-lineage A.Br.Aust 94, a part of major clonal lineage A of B. anthracis (Fig. 3). Members of this sub-lineage have been identified throughout the five continents viz., Australia, Africa (South Africa, Namibia, Mozambique), America (USA), Europe (Germany, Great Britain, The Netherlands) and Asia (China, Turkey, Georgia, Thailand, India) . The origin of this lineage has been reported in Asia and/or Middle East such as India, China and Turkey  and had possibly spread through the ancient trade route known as the Silk Road that extends from Europe, the Middle East and portions of Asia . The A.Br.Asut94 sub-lineage is also associated with anthrax from humans and non-bovines such as zebra, swine , sheep and dog .
In contrast to other regions of the world that are extensively sampled [47,48,49,50], the Indian sub-continent is still poorly studied . The southern states of the sub-continent are endemic to the disease due to unprotected livestock population and also attributable to warm humid climate, alkaline calcareous soil favoring survival and germination of anthrax spores . The high density of livestock population favors the transmission of disease between more susceptible animals and less susceptible human beings. Multiple animal or human cutaneous anthrax outbreaks from India at Vellore , Chittoor , Ramanathapuram , and Visakhapatnam  have been reported earlier, but without sub-lineage analysis of the B. anthracis isolates. However, the 10 Indian B. anthracis isolates that have been delineated as A.Br.Aust94 sub-lineage by canSNP analysis could not be attributed to any of these outbreaks . Consequently, the present study attributes A.Br.Aust94 lineage of B. anthracis to the anthrax outbreak among cattle in Karnataka state and also reiterates the prevalence of this lineage in India as reported earlier .
In summary, we can conclude that all the twelve B. anthracis isolates shared similar 16S type 7 and pag type I. Our findings also revealed that all the isolates examined belonged to clonal lineage A (worldwide lineage), specifically the sub-lineage A.Br.Aust94 similar to previously characterized strains from the Indian sub-continent. Collectively, the genotyping findings evidenced the fact that the anthrax outbreak cases in Gundlupet region might have originated from a common source of infection.
Bacterial strains, growth conditions and DNA extraction
The bacterial cultures used in the present study isolated and identified as B. anthracis earlier  are listed in Additional file 4. The Bacillus strains were maintained in Brain Heart Infusion (BHI) broth and incubated at 37 °C with shaking. All the media were procured from Hi-Media, India. Antibiotics were procured from Sigma-Aldrich (India). The B. cereus ATCC 14579 strain and B. anthracis BA10 strain from Defence Food Research Laboratory (DFRL), Mysore, India repository were used as the negative and positive controls respectively, for all the PCR experiments. Total genomic DNA (Additional file 4) was extracted from every bacterial culture using Gen Elute Bacterial Genomic DNA Kit (Sigma, India) following the manufacturer’s instructions. The DNA samples thus obtained were diluted to a final concentration of 100 ng/μl for PCR experiments.
16S rDNA sequencing
The 16S rDNA gene sequencing was performed for the 12 B. anthracis isolates using universal bacterial primers, 16S-27F and 16S-1488R (Additional file 5) . All PCR reactions were performed in Bio-Rad S1000 thermal cycler (Bio-Rad, USA) wherein the samples were maintained at 94 °C for 4 min for initial denaturation followed by 30 cycles of denaturation at 94 °C for 30 s, annealing at 56 °C for 45 s, extension at 72 °C for 1 min and a final extension at 72 °C for 8 min. PCR products were then gel purified using Gen Elute gel extraction Kit (Sigma, India) and sequenced in both forward and reverse directions with 16S-27F and 16S-1488R primers using BigDye Terminator v3.1 Cycle sequencing kit (Applied Biosystems, USA). Nucleotide sequences were obtained using Sequencing Analysis Software v5.4 in ABI 3500 Genetic Analyzer (Applied Biosystems, USA). These sequence data have been submitted to the GenBank databases under accession numbers MN421996-MN422007. The 16S rDNA nucleotide sequences from additional 170 strains of B. cereus sensu lato (Additional file 6) were included in the present study for phylogenetic comparisons of the 12 B. anthracis isolates of the present study. This collection included 90 B. cereus sensu lato strains from GenBank (28 B. cereus strains, 28 B. anthracis strains, 30 B. thuringiensis strains, and 4 B. mycoides strains); 33 Bacillus strains (15 inclusivity and 18 exclusivity strains) as per Association of Analytical Communities (AOAC) International panel  and 47 Bacillus strains studied from Sacchi et al.  scheme. A total of 182 (170 from the GenBank and 12 of the isolates under study) 16S rDNA gene sequences were analyzed in this study. Phylogenetic analysis was performed using MEGA software, version X  to form a Neighbor-Joining tree of all the 182 B. cereus sensu lato 16S rDNA gene sequences thus compiled. When applied, bootstrap analysis was computed with 1000 replicates.
Detection for the presence of B. anthracis prophages and B. anthracis specific loci
The presence of B. anthracis specific prophages (lambdaBa01, lambdaBa02, lambdaBa03, lambdaBa04) were examined as described earlier with slight modifications for B. anthracis strains from India (DFR.BHE 1–12) . Briefly, all the PCR reactions were carried out in Bio-Rad S1000 thermal cycler (Bio-Rad, USA) under the following conditions: 94 °C initial denaturation for 4 min, and 72 °C final elongation for 8 min, with 30 cycles of 94 °C for 1 min, 52 °C for 45 s, and 72 °C for 1 min. Additionally, the four B. anthracis specific loci (A, C, D and E) were simultaneously detected in the genomes of all B. anthracis isolates as described earlier . The PCR amplicons were electrophoresed in 2% agarose gel, stained with ethidium bromide and visualized under UV transillumination (G-box, Syngene, India). Primer sequences along with their respective amplicon sizes are listed in Additional file 5.
Protective Antigen (PA) sequencing
The 2158 bp region of the protective antigen gene (pag) (GenBank accession No. M22589) was amplified by PCR using PA-F and PA-R primers (this study). The primers were designed in the present study using Gene runner ver. 5.0 (Additional file 5) to amplify pag from all the twelve B. anthracis isolates along with the DNA of control B. anthracis BA10. The PCR reactions were maintained at 94 °C initial denaturation for 4 min and 72 °C final elongation for 8 min, followed by 30 cycles of 94 °C for 1 min, 56 °C for 1 min and 72 °C for 2 min 30 s. PCR products were then gel purified using Gen Elute gel extraction Kit (Sigma, India). The purified amplicons (PA-F and PA-R) were then cloned and sequenced using the custom synthesized primers (Sigma, India) from the pag sequence (GenBank accession No. M22589) following the strategy as depicted in Fig. 4. Briefly, the purified products were cloned in pTZ57R/T vector using an InsTA clone PCR cloning kit (Thermo Scientific, India) according to manufacturer’s instructions and the recombinant plasmids were extracted using Gen Elute Plasmid Miniprep Kit (Sigma, India). Bidirectional sequencing was then performed on the cloned pag amplicons using M13 F and M13 R primers (Additional file 5). Additionally, the primer walking strategy was followed using the primers PA-F1 and PA-F2 to sequence the intermediate regions of this cloned fragment starting from 2289th and 2826th bp positions respectively (based on GenBank accession No. M22589). The sequencing reads thus obtained using the four sequencing primers were then spliced to omit the overlapping nucleotide regions and form a single consensus sequence. In silico analysis was then carried out to screen the seven significant nucleotide positions according to the Price et al.  scheme for pag genotyping. These sequence data have been submitted to the GenBank databases under accession numbers MN447720-MN447731.
A collection of additional pag sequences from 29 B. anthracis strains and 5 B. cereus strains for which pag gene sequence data was available in NCBI database and were taken for the phylogenetic comparison along with the 12 B. anthracis DFR.BHE strains 1–12 in the present study. A total of these 46 pag gene sequences were determined and in silico screening for seven significant positions of mutations  was carried out. The seven SNPs thus obtained were concatenated into a single consensus nucleotide sequence for each isolate and further used for the maximum likelihood tree construction from an alignment with 1000 replicates using MEGA software, version X .
All the twelve B. anthracis DFR.BHE strains 1–12 were genotyped as described  using the 13 canonical (can) SNP loci (Additional file 5). PCR amplifications were performed under following conditions: 4 min at 94 °C initial denaturation and 72 °C final elongation for 8 min with 30 cycles of 94 °C for 30 s, 56 °C for 30 s and 72 °C for 30 s. The amplicons thus obtained were cloned in pTZ57R/T vector as described previously in section 2.4 and sequenced subsequently using M13F and M13R primers (Additional file 5). The nucleotide sequences thus obtained were further screened in silico for the 13 significant nucleotide positions. The 13 identified canSNPs were concatenated into a single consensus nucleotide sequence for each isolate and were used for the maximum likelihood tree construction from an alignment (gap open cost, 15; gap extension cost, 6.66; end gap cost, free) with 1000 replicates using MEGA software, version X .
Availability of data and materials
All data generated or analyzed during this study are included in this published article and its supplementary information.
Canonical Single Nucleotide Polymorphism
Brain Heart Infusion
Polymersae Chain Reaction
Stakeholder Panel on Agent Detection Assays
Association of Analytical Communities
- pag :
Protective Antigen Gene
Parsimonious informative site
Basic Local Alignment Search Tool
National Center for Biotechnology Information
Next Generation sequencing
National Animal Disease Referral Expert System
- NCDC :
National Centre for Disease Control
Indian Council of Agricultural Research - National Institute of Veterinary Epidemiology and Disease Informatics
Centers for Disease Control and Prevention
Liu Y, Du J, Lai Q, Zeng R, Ye D, Xu J, Shao Z. Proposal of nine novel species of the Bacillus cereus group. Int J Syst Evol Microbiol. 2017;67(8):2499–508.
Keim P, Pearson T, Okinaka RT. Evolution of Bacillus anthracis, causative agent of anthrax. In: Baquero F, Nombela C, Cassell GH, Gutiérrez‐Fuentes JA, editors. Evolutionary Biology of Bacterial and Fungal Pathogens. 2008. p. 523–33.
Helgason E, Økstad OA, Caugant DA, Johansen HA, Fouet A, Mock M, Hegna I, Kolstø AB. Bacillus anthracis, Bacillus cereus, and Bacillus thuringiensis—one species on the basis of genetic evidence. Appl Environ Microbiol. 2000;66(6):2627–30.
Maughan H, Van der Auwera G. Bacillus taxonomy in the genomic era finds phenotypes to be essential though often misleading. Infect Genet Evol. 2011;11(5):789–97.
Sacchi CT, Whitney AM, Mayer LW, Morey R, Steigerwalt A, Boras A, Weyant RS, Popovic T. Sequencing of 16S rRNA gene: a rapid tool for identification of Bacillus anthracis. Emerg Infect Dis. 2002;8(10):1117.
Hakovirta JR, Prezioso S, Hodge D, Pillai SP, Weigel LM. Identification and analysis of informative single nucleotide polymorphisms in 16S rRNA gene sequences of the Bacillus cereus group. J Clin Microbiol. 2016;54(11):2749–56.
Radnedge L, Agron PG, Hill KK, Jackson PJ, Ticknor LO, Keim P, Andersen GL. Genome differences that distinguish Bacillus anthracis from Bacillus cereus and Bacillus thuringiensis. Appl Environ Microbiol. 2003;69(5):2755–64.
Sozhamannan S, Chute MD, McAfee FD, Fouts DE, Akmal A, Galloway DR, Mateczun A, Baillie LW, Read TD. The Bacillus anthracis chromosome contains four conserved, excision-proficient, putative prophages. BMC Microbiol. 2006;6(1):34.
Van Ert MN, Easterday WR, Huynh LY, Okinaka RT, Hugh-Jones ME, Ravel J, Zanecki SR, Pearson T, Simonson TS, U'Ren JM, Kachur SM. Global genetic population structure of Bacillus anthracis. PLoS One. 2007;2(5):e461.
Sonavale KP, Shaikh MR, Kadam MM, Pokharkar VG. Livestock sector in India: a critical analysis. Asian J Agri Ext Eco Socio. 2020;38(1):51–62.
Thappa DM. Cutaneous Anthrax—still a reality in India. Ann Natl Acad Med Sci (India). 2019;55(03):119–23.
Goel AK. Anthrax: a disease of biowarfare and public health importance. World J Clin Cases. 2015;3(1):20.
Asokan GV, Asokan V, Tharyan P. One health national programme across species on zoonoses: a call to the developing world. Infect Ecol Epidemiol. 2011;1(1):8293.
Kumar S, Swain S, Preetha GS, Singh BS, Aggarwal D. Zoonotic diseases in India. Indian J Community Med. 2020;45(Suppl 1):S1.
McKenzie JS, Dahal R, Kakkar M, Debnath N, Rahman M, Dorjee S, Naeem K, Wijayathilaka T, Sharma BK, Maidanwal N, Halimi A. One health research and training and government support for one health in South Asia. Infect Ecol Epidemiol. 2016;6(1):33842.
Sahoo KC, Negi S, Barla D, Badaik G, Sahoo S, Bal M, Padhi AK, Pati S, Bhattacharya D. The landscape of Anthrax prevention and control: stakeholders’ perceptive in Odisha, India. Int J Environ Res Public Health. 2020;17(9):3094.
Govindaraj G, Hiremath J, Reddy GB, Siju SJ, Yogishardhaya R, Prajapati A. ICAR NIVEDI Annual Report 2018–19. 2019. p. 34–5.
Vijaikumar M, Thappa DM, Karthikeyan K. Cutaneous anthrax: an endemic outbreak in South India. J Trop Pediatr. 2002;48(4):225–6.
Rao GRR, Padmaj J, Lalitha MK, Rao PVK, Gopal KV, Kumar HKY, Mohanraj P. An outbreak of cutaneous anthrax in a non-endemic district-Visakhapatnam in Andhra Pradesh. Indian J Dermatol Venereol Leprol. 2005;71(2):102.
Reddy R, Parasadini G, Rao P, Uthappa CK, Murhekar MV. Outbreak of cutaneous anthrax in Musalimadugu village, Chittoor district, Andhra Pradesh, India, July-august 2011. J Infect Dev Ctries. 2012;6(10):695–9.
Balachandrudu B, Bindu SA, Kuma CN, Malakondaiah P. An outbreak of cutaneous anthrax in a tribal area of Visakhapatnam district, Andhra Pradesh. J NTR Univ Health. 2018;7(1):49.
Kingston JJ, Majumder S, Uppalapati SR, Makam SS, Urs RM, Murali HS, Batra HV. Anthrax outbreak among cattle and its detection by extractable antigen 1 (EA1) based sandwich ELISA and immuno PCR. Indian J Microbiol. 2015;55(1):29–34.
Liu Y, Lai Q, Göker M, Meier-Kolthoff JP, Wang M, Sun Y, Wang L, Shao Z. Genomic insights into the taxonomic status of the Bacillus cereus group. Sci Rep. 2015;5(1):1–1.
Ozanich RM, Colburn HA, Victry KD, Bartholomew RA, Arce JS, Heredia-Langner A, Jarman K, Kreuzer HW, Bruckner-Lea CJ. Evaluation of PCR systems for field screening of Bacillus anthracis. Health Secur. 2017;15(1):70–80.
Daffonchio D, Raddadi N, Merabishvili M, Cherif A, Carmagnola L, Brusetti L, Rizzi A, Chanishvili N, Visca P, Sharp R, Borin S. Strategy for identification of Bacillus cereus and Bacillus thuringiensis strains closely related to Bacillus anthracis. Appl Environ Microbiol. 2006;72(2):1295–301.
Ko KS, Kim JM, Kim JW, Jung BY, Kim W, Kim IJ, Kook YH. Identification of Bacillus anthracis by rpoB sequence analysis and multiplex PCR. J Clin Microbiol. 2003;41(7):2908–14.
Qi Y, Patra G, Liang X, Williams LE, Rose S, Redkar RJ, DelVecchio VG. Utilization of the rpoB gene as a specific chromosomal marker for real-time PCR detection of Bacillus anthracis. Appl Environ Microbiol. 2001;67(8):3720–7.
Hurtle W, Bode E, Kulesh DA, Kaplan RS, Garrison J, Bridge D, House M, Frye MS, Loveless B, Norwood D. Detection of the Bacillus anthracis gyrA gene by using a minor groove binder probe. J Clin Microbiol. 2004;42(1):179–85.
Yamada S, Ohashi E, Agata N, Venkateswaran K. Cloning and nucleotide sequence analysis of gyrB of Bacillus cereus, B. thuringiensis, B. mycoides, and B. anthracis and their application to the detection of B. cereus in rice. Appl. Environ. Microbiol. 1999;65(4):1483–90.
Kim K, Seo J, Wheeler K, Park C, Kim D, Park S, Kim W, Chung SI, Leighton T. Rapid genotypic detection of Bacillus anthracis and the Bacillus cereus group by multiplex real-time PCR melting curve analysis. FEMS Immunol Med Microbiol. 2005;43(2):301–10.
Ryu C, Lee K, Yoo C, Seong WK, Oh HB. Sensitive and rapid quantitative detection of anthrax spores isolated from soil samples by real-time PCR. Microbiol Immunol. 2003;47(10):693–9.
Koehler TM. Bacillus anthracis physiology and genetics. Mol Asp Med. 2009;30(6):386–96.
Kolstø AB, Tourasse NJ, Økstad OA. What sets Bacillus anthracis apart from other Bacillus species. Annu Rev Microbiol. 2009;63:451–76.
Vodkin MH, Leppla SH. Cloning of the protective antigen gene of Bacillus anthracis. Cell. 1983;34(2):693–7.
Price LB, Hugh-Jones M, Jackson PJ, Keim P. Genetic diversity in the protective antigen gene of Bacillus anthracis. J Bacteriol. 1999;181(8):2358–62.
Hoffmaster AR, Fitzgerald CC, Ribot E, Mayer LW, Popovic T. Molecular subtyping of Bacillus anthracis and the 2001 bioterrorism-associated anthrax outbreak, United States. Emerg Infect Dis. 2002;8(10):1111.
Shahcheraghi SH, Ayatollahi J. pXO1-and pXO2-like plasmids in Bacillus cereus and B. thuringiensis. Jundishapur J Microbiol. 2013;6(10):e8482.
Cheung DT, Kam KM, Hau KL, Au TK, Marston CK, Gee JE, Popovic T, Van Ert MN, Kenefic L, Keim P, Hoffmaster AR. Characterization of a Bacillus anthracis isolate causing a rare case of fatal anthrax in a 2-year-old boy from Hong Kong. J Clin Microbiol. 2005;43(4):1992–4.
Tan Z, Qi X, Gu L, Bao C, Tang F, Zhu Y. Molecular characterization of Bacillus anthracis directly from patients' eschar and beef in an anthrax outbreak in Jiangsu Province, China, 2012. Am J Trop Med Hyg. 2014;91(3):574–6.
Sanam MU, Asmara W, Wahyuni AE, Wibowo MH. Genetic analysis on protective antigenic gene of Bacillus anthracis isolates of Central Java and Yogyakarta. Jurnal Veteriner. 2015;16(1):15–24.
Petosa C, Collier RJ, Klimpel KR, Leppla SH, Liddington RC. Crystal structure of the anthrax toxin protective antigen. Nature. 1997;385(6619):833.
Abboud N, Casadevall A. Immunogenicity of Bacillus anthracis protective antigen domains and efficacy of elicited antibody responses depend on host genetic background. Clin Vaccine Immunol. 2008;15(7):1115–23.
Kumar S, Stecher G, Li M, Knyaz C, Tamura K. MEGA X: molecular evolutionary genetics analysis across computing platforms. Mol Biol Evol. 2018;35(6):1547–9.
Antwerpen M, Ilin D, Georgieva E, Meyer H, Savov E, Frangoulidis D. MLVA and SNP analysis identified a unique genetic cluster in Bulgarian Bacillus anthracis strains. Eur J Clin Microbiol. 2011;30(7):923–30.
Okutani A, Inoue S, Morikawa S. Comparative genomics and phylogenetic analysis of Bacillus anthracis strains isolated from domestic animals in Japan. Infect Genet Evol. 2019;71:128–39.
Sahin M, Buyuk F, Baillie L, Wölfel R, Kotorashvili A, Rehn A, Antwerpen M, Grass G. The identification of novel single nucleotide polymorphisms to assist in mapping the spread of Bacillus anthracis across the southern Caucasus. Sci Rep. 2018;8(1):1–7.
Simonson TS, Okinaka RT, Wang B, Easterday WR, Huynh L, U'Ren JM, Dukerich M, Zanecki SR, Kenefic LJ, Beaudry J, Schupp JM. Bacillus anthracis in China and its relationship to worldwide lineages. BMC Microbiol. 2009;9(1):71.
Girault G, Blouin Y, Vergnaud G, Derzelle S. High-throughput sequencing of Bacillus anthracis in France: investigating genome diversity and population structure using whole-genome SNP discovery. BMC Genomics. 2014;15(1):288.
Kenefic LJ, Pearson T, Okinaka RT, Schupp JM, Wagner DM, Ravel J, Hoffmaster AR, Trim CP, Chung WK, Beaudry JA, Foster JT. Pre-columbian origins for north American anthrax. PLoS One. 2009;4(3):e4813.
Khmaladze E, Birdsell DN, Naumann AA, Hochhalter CB, Seymour ML, Nottingham R, Beckstrom-Sternberg SM, Beckstrom-Sternberg J, Nikolich MP, Chanturia G, Zhgenti E. Phylogeography of Bacillus anthracis in the country of Georgia shows evidence of population structuring and is dissimilar to other regional genotypes. PLoS One. 2014;9(7):e102651.
Derzelle S, Aguilar-Bultet L, Frey J. Comparative genomics of Bacillus anthracis from the wool industry highlights polymorphisms of lineage a.Br.Vollum. Infect Genet Evol. 2016;46:50–8.
Suma AP, Suresh KP, Gajendragad MR, Kavya BA. Outbreak prediction of anthrax in Karnataka using poisson, negative-binomial and zero truncated models; 2017.
David S, Valentina GO, Lalitha MK. Cutaneous anthrax involving the eyelids. Indian J Med Microbiol. 1999;17(2):92.
Rajasokkappan S, Selvaraju G, Ramkumar P, Rajesh A. Retrospective study of Anthrax outbreaks among animals in Ramanathapuram District. Int J Sci Environ Technol. 2016;5(5):3130–3.
Brosius J, Palmer ML, Kennedy PJ, Noller HF. Complete nucleotide sequence of a 16S ribosomal RNA gene from Escherichia coli. PNAS. 1978;75(10):4801–5.
Authors are grateful to Director, DFRL, Mysore, for providing the necessary facilities to carry out this work. AR is thankful to DST for providing INSPIRE fellowship. AR also acknowledges her colleagues Shreya, Kunal and Monojit for their assistance in completion of the work.
Ethics approval and consent to participate
Consent for publication
The authors declare that they have no competing interest.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Representative electropherograms of the Bacillus anthracis 16S rDNA gene sequences. The electropherogram of 16S rDNA reverse complement read of B. anthracis DFR.BHE 3 strain. Arrows indicate (a) single peak corresponding to position 1148 indicating A (T in reverse read here) (b) dual peak at position 1139 indicating mixed base pair G/A [C/T)in reverse read here] described by Hakovirta et al. 2016 .
Additional nucleotide positions of mutations in the Bacillus cereus sensu lato 16S rDNA gene found in the present study. Numbers in paranthesis refers to polymorphic positions in the 16S rDNA gene and are numbered from 1 to 40. The 28 new 16S types are highlighted with colored boxes. The new genotypes of B. anthracis, B. cereus, B. thuringiensis and B. mycoides are boxed in green, orange, blue and yellow color respectively. The 16S types shared by both B. cereus and B. thuringiensis are boxed in purple color. The genotypes shared by both B. cereus and B. anthracis are boxed in red color.
PCR for the screening of Bacillus anthracis specific prophages and loci in the 12 B. anthracis DFR.BHE strains 1–12. A. PCR with B. anthracis specific prophage genes, (a).lambda01, (b). lambda02, (c). lambda03 and (d). lambda04 and B. PCR with B. anthracis specific loci (a) dhp 61.183 (loci A), (b). dhp 77.002 (loci C), (c). dhp 73.019 (loci D), and (d). dhp 73.017 (loci E). Lane 1–12, B. anthracis DFRL.BHE strains 1–12 along with positive (Lane P: B. anthracis BA10) and negative (Lane N: B. cereus ATCC 14579) controls. Lane M: 100 bp molecular marker.
List of bacterial cultures used in the present study.
List of primers and their amplicons size used in the present study.
List of strains and their accession numbers used for multiple sequence alignment of 16S rDNA gene and pag gene in the present study.
About this article
Cite this article
Roonie, A., Majumder, S., Kingston, J.J. et al. Molecular characterization of B. anthracis isolates from the anthrax outbreak among cattle in Karnataka, India. BMC Microbiol 20, 232 (2020). https://doi.org/10.1186/s12866-020-01917-1
- B. anthracis
- 16S rDNA
- Anthrax detection