A Survey of Antimicrobial Resistance Determinants in Category A Select Agents, Exempt Strains, and Near-Neighbor Species

A dramatic increase in global antimicrobial resistance (AMR) has been well documented. Of particular concern is the dearth of information regarding the spectrum and prevalence of AMR within Category A Select Agents. Here, we performed a survey of horizontally and vertically transferred AMR determinants among Category A agents and their near neighbors. Microarrays provided broad spectrum screening of 127 Francisella spp., Yersinia spp., and Bacillus spp. strains for the presence/absence of 500+ AMR genes (or families of genes). Detecting a broad variety of AMR genes in each genus, microarray analysis also picked up the presence of an engineered plasmid in a Y. pestis strain. High resolution melt analysis (HRMA) was also used to assess the presence of quinolone resistance-associated mutations in 100 of these strains. Though HRMA was able to detect resistance-causing point mutations in B. anthracis strains, it was not capable of discriminating these point mutations from other nucleotide substitutions (e.g., arising from sequence differences in near neighbors). Though these technologies are well-established, to our knowledge, this is the largest survey of Category A agents and their near-neighbor species for genes covering multiple mechanisms of AMR.


Introduction
Ever-increasing levels of antimicrobial resistance (AMR) are considered a top global health concern and pose a "threat to the success and continuation of clinical medicine as we know it." [1,2]. The issue of AMR is also challenging for biosecurity reasons. While morbidity/mortality and ease of Horizontally transferred resistance determinants also pose a significant threat. MDR Yersinia pestis strains isolated in Madagascar in the 1990s harbored plasmids with up to eight AMR determinants, at least four of which confer resistance to recommended therapeutics or prophylactic treatments [13][14][15][16][17]. Fortunately, the strains responsible for more recent outbreaks in Madagascar, the Democratic Republic of Congo, and Uganda do not carry these plasmids and are susceptible to the standard therapeutics [18,19]. Even so, AMR determinants present on mobile elements have the potential to spread widely; a number of reports have described the means of transferring naturally occurring plasmids harboring AMR determinants to and from these strains [14,20,21]. Furthermore, ever increasing numbers of cloning vectors are being described that are capable of being transferred to and from Y. pestis and other Category A agents, such as Bacillus anthracis and Francisella tularensis [22][23][24][25][26]. In most cases, common AMR markers are used for selection and could potentially be suggestive of engineered strains.
In spite of the large number of publications describing the pathophysiology of these Category A agents, there still remains a dearth of information regarding the full spectrum and prevalence of AMR and their molecular mechanisms within these agents and their near neighbors. For this reason, we performed a survey of 62 Bacillus spp., 23 Francisella spp., and 42 Yersinia spp. strains for both horizontally and vertically transferred AMR determinants; included amongst these strains were 46 Category A Select Agents. Microarrays were used for broad spectrum screening of all 127 strains for the presence/absence of 500+ AMR genes and gene families, while high resolution melt analysis (HRMA) was used to assess the presence of quinolone resistance-associated mutations in 100 of these strains. To our knowledge, this is the largest survey of Category A agents and their near-neighbor species for genes covering multiple mechanisms of AMR.

Results and Discussion
Phenotypic antimicrobial susceptibility testing (AST) is still considered the gold standard for determination of resistance profiles. However, ASTs cannot be used to trace the spread of AMR strains nor the mobile genetic elements responsible for transmission of resistance. Phenotypic testing further cannot distinguish between various resistance mechanisms with overlapping specificities and are generally inadequate for epidemiological and forensic tracking. While not always predictive of phenotype, molecular tests can nonetheless determine the genetic basis and overall mechanisms for resistance, as well as the potential for transmission to other strains or species via horizontal gene transfer.
Here, we genetically characterized 127 archived Category A agents and their near neighbors with respect to AMR. To accomplish this, we determined the presence/absence of >500 AMR determinants using the ARDM v.3.1. In parallel, we used a combination of targeted PCR and HRMA to identify mutations in gyrA, parC, and parE that confer resistance to fluoroquinolones used as first-line therapeutics. By combining the broad screening capability of microarrays with the fine discriminatory capabilities of HRMA, we detected markers acquired through both horizontal transfer and mutation. The dataset generated here provides a baseline measurement of the spectrum and prevalence of various resistance mechanisms in both Category A select agents and their near neighbors. These baseline data may be helpful in development of predictive algorithms designed to assess risk associated with outbreaks with naturally occurring or engineered threats, as well as in guiding development of new therapeutics.

Testing for Horizontally-Acquired Genes-Microarray Analysis
Based on the ARDM v.2 [27] and v.3 [28], the ARDM v.3.1 contains probes directed against emerging AMR determinants of global importance (e.g., bla NDM ), as well as those derived from several Select Agent species (including F. tularensis, B. anthracis, Y. pestis), their near neighbors, and more distantly related species in the same genera. In total, the ARDM v.3.1 content (4480 probes per two sub-arrays) covers 512 AMR determinants, conferring resistance to 19 categories of antimicrobials ( Figure 1a). Microarray content was derived from Actinobacteria, Firmicutes, Bacteroidetes, and Proteobacteria phyla (Figure 1b). Pertinent to the current study, the following determinants were derived from Select Agents and their near neighbors: 7 genes from Francisella spp. (blaA, blaB, tet, bcr/cfl1, bcr/cflA2, and bcr/cflA3 from F. tularensis, and bla from F. philomiragia); 6 genes from Yersinia spp. (qacE, yegB, emrA, emrB, and bacA from Y. pestis CO92, and vatF/sat from Y. enterocolitica); and 30 genes from a wide variety of Bacillus species (   Twenty F. tularensis and three F. philomiragia strains were characterized on the ARDM v.3.1. Not surprisingly, all F. tularensis strains were positive for the six F. tularensis-derived genes and the three F. philomiragia strains were positive for the single F. philomiragia-derived gene, bla ( Table 2). One of the F. philomiragia strains, Jensen O#319-029, was also positive for the presence of bcr/cflA1. The reference sequence for bcr/cflA1 was compared with the published sequence for this strain (CP009343.1) and the other two F. philomiragia strains via BLAST. The bcr/cflA1 reference sequence had 91.6% sequence identity with Jensen O#319-029 while sequence identity with Jensen O#319L and Jensen O#319-036 was only~65%. Similar BLAST comparisons showed~80% or lower sequence identities between the three F. philomiragia strains and the other F. tularensis-derived genes on the microarray. Though not Int. J. Mol. Sci. 2020, 21, 1669 5 of 31 PCR-confirmed, these results support previous observations [27,29] that the specificity of the ARDM can enable detection of closely related genes (or families of genes), but that a sequence identity of >90% is required.  [30], seven Y. enterocolitica strains, and nine strains from six species related more distantly to both Y. pestis/Y. pseudotuberculosis and Y. enterocolitica. Not surprisingly, all of the Y. pestis (17 isolates) and Y. pseudotuberculosis (9 isolates) strains were positive for the same five Y. pestis-derived genes: bacA, qacE, emrA, emrB, and yegB. Importantly, Y. pestis KIM(pCD1Ap+) was also positive for the presence of bla TEM ; KIM(pCD1Ap+) is a KIM5 strain harboring plasmid pCD1 into which bla TEM was inserted into the yadA gene (i.e., as a selectable marker) [31]. Therefore, the ARDM was able to accurately detect the presence of a non-native AMR gene engineered into Y. pestis.
All Y. enterocolitica strains except strain 8081 were positive for the Y. enterocolitica-derived gene, vatF/sat (Table 3); strain 8081 tested negative for this gene. A BLAST search of Y. enterocolitica 8081 whole genome sequence for vatF/sat yielded a single sequence with only 79% identity, indicating that the microarray's negative results were indeed accurate. None of the more distantly related Yersinia spp. strains were positive for any of the 512 resistance determinants on the ARDM v.3.1, including those derived from Yersinia spp. (Table 3). DNA preparations from 62 different strains from the B. cereus group were characterized on the microarray, including 41 B. anthracis, 13 B. cereus, and eight B. thuringiensis strains. Not surprisingly, all of the B. anthracis strains were positive for the same sixteen B. anthracisor B. cereus-derived genes (Table 4). A plasmid-borne gene encoding a chloramphenicol acetyltransferase originally isolated from Staphylococcus aureus, cat(pSCS5), was also detected in a number of B. anthracis samples. The cat(pSCS5) gene has no identical allele in B. anthracis. However, six of the ARDM's 10 probes for this gene have a 25-nt sequence with 92% identity to the B. anthracis-derived cat gene represented on the ARDM v.3.1 ( Figure 2); this latter cat gene was observed in all B. anthracis and most B. cereus and B. thuringiensis strains. These cat(pSCS5)-positive calls therefore represent false-positive results. To minimize the potential for future false-positives, redesign of the ARDM-or more specifically, the group of cat(pSCS5) probes-is needed to improve specificity and eliminate redundant sequences.     Four B. cereus and five B. thuringiensis were also positive for the same 16 genes observed in all of the B. anthracis strains. Of the remaining samples that did not harbor all 16 genes, the number of positive genes ranged from two (B. thuringiensis var. israelensis) to 15 genes (B. cereus D17 and 03BB108). Perhaps not surprisingly, the number of genes detected was generally correlated to the phylogenetic similarity between strains ( Figure 3), with fewer genes detected in near-neighbor strains more distantly related to B. anthracis. In general, the microarray-based survey of AMR determinants did not reveal any significant surprises, with Select Agent species harboring AMR genes expected to be present. Similarly, nearneighbor species harbored the expected genes (e.g., F. philomiragia-derived bla in all three F. philomiragia strains, vatF/sat in all but one Y. enterocolitica strains). The microarray was, however, able to detect a non-native AMR determinant (blaTEM) in a single Y. pestis strain (KIM (pCD1Ap)); this gene was engineered into a plasmid in the parent strain as a selectable marker in knock-out experiments. The ability to detect such a non-native gene in a Select Agent species-within <24-h and without a priori knowledge of its potential presence-illustrates the value of a broad-based screening method such as microarray analysis.

HRMA
High-resolution melt analysis (HRMA) combines the finely detailed discriminatory capabilities of sequence analysis with the rapidity, relevance, and user-friendliness of PCR. Originally developed for genotyping studies, HRMA determines the melting profile of double stranded DNA sequences using a DNA intercalating dye. The shape and position of the melting curve yield information about the composition of a given DNA sequence from which mutational data can be derived; substitution of a G/C with A/T will reduce the Tm while an A/TG/C change will result in a similar increase in Tm. In this manner, HRMA can provide the capability to distinguish allelic variants which play an In general, the microarray-based survey of AMR determinants did not reveal any significant surprises, with Select Agent species harboring AMR genes expected to be present. Similarly, near-neighbor species harbored the expected genes (e.g., F. philomiragia-derived bla in all three F. philomiragia strains, vatF/sat in all but one Y. enterocolitica strains). The microarray was, however, able to detect a non-native AMR determinant (bla TEM ) in a single Y. pestis strain (KIM (pCD1Ap)); this gene was engineered into a plasmid in the parent strain as a selectable marker in knock-out experiments. The ability to detect such a non-native gene in a Select Agent species-within <24-h and without a priori knowledge of its potential presence-illustrates the value of a broad-based screening method such as microarray analysis.

HRMA
High-resolution melt analysis (HRMA) combines the finely detailed discriminatory capabilities of sequence analysis with the rapidity, relevance, and user-friendliness of PCR. Originally developed for genotyping studies, HRMA determines the melting profile of double stranded DNA sequences using a DNA intercalating dye. The shape and position of the melting curve yield information about the composition of a given DNA sequence from which mutational data can be derived; substitution of a G/C with A/T will reduce the T m while an A/T→G/C change will result in a similar increase in T m . In this manner, HRMA can provide the capability to distinguish allelic variants which play an important role in resistance arising from mutations. As these mutations are not detectable with the ARDM alone, HRMA is a valuable approach that is highly complementary to the ARDM.
High level resistance to ciprofloxacin-a first-or second-line therapeutic for tularemia, plague, and anthrax-often arises from point mutations in the QRDRs of gyrA, gyrB, parC, parE. Loveless [34] previously demonstrated that HRMA could be used to detect point mutations in QRDRs of F. tularensis, Y. pestis, and B. anthracis. Here, we applied their strategy to assess sequence differences in 10 Francisella spp. strains (gyrA, parE), 34 Yersinia spp. strains (gyrA only), and 60 Bacillus spp. strains (gyrA, parC) to determine the practicality of this method for rapid screening.

Francisella spp.
QRDRs in gyrA and parC were assessed by HRMA in seven F. tularensis and three F. philomiragia strains using PCR primers designed to bracket mutations observed to result in ciprofloxacin resistance [9,34,35]. In the gyrA assay, all F. tularensis strains clustered together (Table 5, Figure 4); these results are consistent with the phenotypic susceptibility of all strains for which ASTs have been reported [36]. Formation of a separate cluster by all three F. philomiragia, with T m s nearly a degree higher, was not surprising. There are three A/T→C/G and one C/G→A/T changes within the gyrA target sequence; these changes would be expected to-and did-cause a significant change in T m . Differences within the primer sequences were also likely responsible for the significantly higher C t s observed with the F. philomiragia samples (average C t of 32.7 versus 13.8 for F. tularensis, p < 0.025).  Comparison of reference gyrA sequence with published sequences from F. philomiragia strains (see Table 8). All F. tularensis strains tested are identical to the reference sequence within the stretch interrogated. Primer sequences are underlined. Arrows indicate nucleotide positions where changes are demonstrated to result in quinolone resistance [34,35]. Sequence differences predicted to result in increased or decreased Tms are highlighted in green and yellow, respectively.   Table 9). All F. tularensis strains tested are identical to the reference sequence within the stretch interrogated. Primer sequences are underlined. Arrows indicate nucleotide positions where changes are demonstrated to result in quinolone resistance [34,35]. Sequence differences predicted to result in increased or decreased T m s are highlighted in green and yellow, respectively. The HRMA primers used here for the parE assay were designed to detect a 5-bp deletion previously noted by Loveless [34]. All of the strains-both F. tularensis and F. philomiragia-clustered together, in spite of second T m peak observed with F. philomiragia O#319-036. Clustering together of the F. tularensis and F. philomiragia strains was not unexpected, though the clustering was not as tight as observed with gyrA. Though the 5-bp region interrogated had a single nucleotide difference (T→A) not predicted to cause a change in T m , there were multiple differences within the primer regions that might affect the T m s of half the amplicons in each PCR reaction (highlighted in Figure 5b). As with the gyrA amplifications, C t s for the three F. philomiragia strains were significantly greater than those of the F. tularensis strains (p < 0.005), likely due to these changes in the primed regions. the F. tularensis and F. philomiragia strains was not unexpected, though the clustering was not as tight as observed with gyrA. Though the 5-bp region interrogated had a single nucleotide difference (TA) not predicted to cause a change in Tm, there were multiple differences within the primer regions that might affect the Tms of half the amplicons in each PCR reaction (highlighted in Figure 5b). As with the gyrA amplifications, Cts for the three F. philomiragia strains were significantly greater than those of the F. tularensis strains (p < 0.005), likely due to these changes in the primed regions. Comparison of reference parE sequence with published sequences for the F. philomiragia strains (see Table 8). All F. tularensis strains tested are 98% identical to the reference sequence within this locus. Primer sequences are underlined. The sequence within the red box indicates the 5-nt deletion of interest [34]; mutations/substitutions predicted to increase or decrease amplicon Tms are highlighted in green or yellow, respectively.

Yersinia spp.
Multiple studies have reported ciprofloxacin resistance arising in Y. pestis from mutations in the gyrA QRDR [34,37-40], but to our knowledge, analogous mutations in parC or parE have not been documented. For this reason, only gyrA-specific HRMA assays were used with Yersinia spp. The site bracketed in the Y. pestis assay encompassed four mutations causing phenotypic ciprofloxacin resistance, with an additional two mutations found within the primer sequences [37,38].  Table 9). All F. tularensis strains tested are 98% identical to the reference sequence within this locus. Primer sequences are underlined. The sequence within the red box indicates the 5-nt deletion of interest [34]; mutations/substitutions predicted to increase or decrease amplicon T m s are highlighted in green or yellow, respectively.

Yersinia spp.
Multiple studies have reported ciprofloxacin resistance arising in Y. pestis from mutations in the gyrA QRDR [34,[37][38][39][40], but to our knowledge, analogous mutations in parC or parE have not been documented. For this reason, only gyrA-specific HRMA assays were used with Yersinia spp. The site bracketed in the Y. pestis assay encompassed four mutations causing phenotypic ciprofloxacin resistance, with an additional two mutations found within the primer sequences [37,38].
HRMA was performed on DNA from thirty-five Yersinia spp. strains, including 10 Y. pestis, eight Y. pseudotuberculosis, and seven Y. enterocolitica isolates. In many cases, two peaks (or a peak and a shoulder) were observed ( Figure 6); T m s of the lower peaks ranged between 79.03 • C and 79.95 • C and of the higher peaks between 78.83 • C and 80.25 • C. All of the Y. pestis strains formed a single tight cluster with average T m of 79.83; five Y. pseudotuberculosis, four Y. enterocolitica, and six more distantly related strains clustered with these Y. pestis strains (  Based on the published sequences available for the Yersinia spp. strains (Figure 7), the basis for principal component 1 used in the clustering analysis is unclear. Though all Y. pestis strains were placed in the same cluster, several Y. pseudotuberculosis strains with sequences identical to the reference sequence were placed in cluster 2. Similarly, amplicons from other strains having predicted C/GA/T changes that would decrease Tm (typically within the primer regions) were split evenly between clusters 1 and 2. The predictability in clustering could potentially be improved through analysis of additional replicates. Furthermore, inclusion of quinolone-resistant strains with documented QRDR mutations (to which we did not have access) would help demonstrate the true efficacy in using HRMA to detect both species-specific and sequence-specific changes predictive of  Based on the published sequences available for the Yersinia spp. strains (Figure 7), the basis for principal component 1 used in the clustering analysis is unclear. Though all Y. pestis strains were placed in the same cluster, several Y. pseudotuberculosis strains with sequences identical to the reference sequence were placed in cluster 2. Similarly, amplicons from other strains having predicted C/G→A/T changes that would decrease T m (typically within the primer regions) were split evenly between clusters 1 and 2. The predictability in clustering could potentially be improved through analysis of additional replicates. Furthermore, inclusion of quinolone-resistant strains with documented QRDR mutations (to which we did not have access) would help demonstrate the true efficacy in using HRMA to detect both species-specific and sequence-specific changes predictive of phenotype. As is, based on the current knowledge of known point mutations causing fluoroquinolone resistance and the location of sequence differences, none of the tested strains would be predicted to be resistant. Indeed, for those strains for which phenotypic susceptibility has been determined (Y. pestis Nairobi, PBM19, Pestoides F, Harbin35, Java9, CO92, Pestoides G, Angola, and Nicholisk 41), none are resistant [41]. At present, however, the utility of HRMA to detect such differences is not clear.  Table 10). The Y. pestis reference sequence is shown at the top with primer sequences underlined. Mutations at the positions shown in red have been demonstrated to cause quinolone resistance [37,38]. Sequence substitutions predicted to cause a decrease in Tm are highlighted in yellow. Y. entero = Y. enterocolitica; Y. krist = Y. kristensenii; Y. pseudo = Y. pseudotuberculosis.

Bacillus spp.
HRMA was performed on DNA preparations from 37 B. anthracis strains. This sample set included DNA preparations from both 18 wild-type strains and 19 plasmid-cured strains selected for development of high-level ciprofloxacin resistance [6] to assess whether documented mutations within gyrA and parC QRDRs could be detected.
In HRMA for gyrA, all B. anthracis wild-type strains produced Tms in the range of 75.9 °C to 76.2 °C, whereas Tms for B. anthracis gyrA mutants had Tms at least 0.5 °C lower ( Figure 8, Table 6). Principal component analysis yielded three clusters: wild-type strains (Cluster 1 = wt), Cluster 2 with most gyrA-mutant strains, and Cluster 3 comprising two samples with significantly lower Tms (mutant strains S3-5 and S3-15). Interestingly, strain S3-16 has two mutations in gyrA (C254→T and A266→C) that should thermodynamically balance each other; Loveless previously observed that this strain clustered with wild-type B. anthracis [34]. However, in our hands, this strain was found in Cluster 2 with other C/G→A/T-mutated B. anthracis strains (Table 6). Unpublished sequences of the targeted loci (S. Sozhamannan, C. Bernhards personal communication, manuscript in preparation) support clustering of the remaining mutants not previously described [6,34] Figure 7. Alignment of interrogated gyrA sequences with published reference sequences from all non-pestis Yersinia spp. strains tested (see Table 11). The Y. pestis reference sequence is shown at the top with primer sequences underlined. Mutations at the positions shown in red have been demonstrated to cause quinolone resistance [37,38]. Sequence substitutions predicted to cause a decrease in T m are highlighted in yellow. Y. entero = Y. enterocolitica; Y. krist = Y. kristensenii; Y. pseudo = Y. pseudotuberculosis.

Bacillus spp.
HRMA was performed on DNA preparations from 37 B. anthracis strains. This sample set included DNA preparations from both 18 wild-type strains and 19 plasmid-cured strains selected for development of high-level ciprofloxacin resistance [6] to assess whether documented mutations within gyrA and parC QRDRs could be detected.
In HRMA for gyrA, all B. anthracis wild-type strains produced T m s in the range of 75.9 • C to 76.2 • C, whereas T m s for B. anthracis gyrA mutants had T m s at least 0.5 • C lower ( Figure 8, Table 7). Principal component analysis yielded three clusters: wild-type strains (Cluster 1 = wt), Cluster 2 with most gyrA-mutant strains, and Cluster 3 comprising two samples with significantly lower T m s (mutant strains S3-5 and S3-15). Interestingly, strain S3-16 has two mutations in gyrA (C254→T and A266→C) that should thermodynamically balance each other; Loveless previously observed that this strain clustered with wild-type B. anthracis [34]. However, in our hands, this strain was found in Cluster 2 with other C/G→A/T-mutated B. anthracis strains (Table 7). Unpublished sequences of the targeted loci (S. Sozhamannan, C. Bernhards personal communication, manuscript in preparation) support clustering of the remaining mutants not previously described [6,34] away from wild-type B. anthracis. Melt curves for wild-type strains, Sterne and ∆ANR, and gyrA-mutated ∆ANR strains, S1-1, S1-2, S2-1, S2-2, S2-3, S3-1, S3-3, and S3-4. (b) Melt curves for the remainder of wild-type and gyrA-mutated ∆ANR strains (see Table 6). curves for the remainder of wild-type and gyrA-mutated ∆ANR strains (see Table 7).
HRMA and clustering of the same set of B. anthracis samples based on parC sequences were complicated by higher PCR failure rates and a narrower temperature range between T m s of wild-type and mutant amplicons. Furthermore, T m s for control strains could shift from day to day by as much as 0.5 • C, making it difficult to compare samples tested on different days, other than to indicate that they clustered with (or separate from) specific control samples (see footnote to Table 7). , whereas only a single non-wild-type cluster was observed on Day 2 (2B). 4 Substitutions shown were described by [6] or [34].
Most of the B. anthracis strains with documented parC mutations clustered separately from wild-type B. anthracis ( Figure 9; Table 7, Clusters 2A, 3A, 3B). Decreased T m s observed with mutant strains were consistent with C242→T changes. However, in contrast to previous observations [34], two strains that should have wild-type parC genes (S1-1 and S1-2 [6]) clustered separately from wild-type strains, with strain S1-2 as the sole member of Cluster 3A. From the melt curve, it is unclear why this sample clusters separately from all other strains tested. However, as with the gyrA results, unpublished sequences of the remaining mutants (S2-1 through S3-18) support their exclusion from the wild-type cluster (S. Sozhamannan, C. Bernhards, personal communication; manuscript in preparation).
strains that should have wild-type parC genes (S1-1 and S1-2 [6]) clustered separately from wild-type strains, with strain S1-2 as the sole member of Cluster 3A. From the melt curve, it is unclear why this sample clusters separately from all other strains tested. However, as with the gyrA results, unpublished sequences of the remaining mutants (S2-1 through S3-18) support their exclusion from the wild-type cluster (S. Sozhamannan, C. Bernhards, personal communication; manuscript in preparation).
Twelve B. cereus and eight B. thuringiensis strains were also analyzed by HRMA to determine whether their gyrA and parC sequences were sufficiently similar to B. anthracis to allow robust amplification, for species discrimination, and for detection of QRDR mutants.
Almost two-thirds (14/19) of the near neighbor strains clustered closely with wild-type B. anthracis (cluster 1) in gyrA assays (Table 8). Two additional clusters were comprised of five strains (cluster 2) and a single strain (cluster 3) with T m s that were 1-1.5 • C higher than the wild-type B. anthracis cluster ( Figure 10). Cluster 2 represented four strains distantly related to B. anthracis (B. thuringiensis kurstaki, morrisoni, israelensis, and B. cereus ATCC 10876) and B. cereus E33L. Although closely related to B. anthracis, [42], E33L's published sequence has, like other strains clustering with it, a T→C substitution within the interrogated sequence ( Figure 10); this substitution-though in a location different from the others-would be expected to cause an increase in T m . The published sequence for B. cereus 4342 shows two T→C substitutions within the interrogated sequence, which may explain its larger increase in T m and formation of its own outlying cluster. it, a TC substitution within the interrogated sequence ( Figure 10); this substitution-though in a location different from the others-would be expected to cause an increase in Tm. The published sequence for B. cereus 4342 shows two TC substitutions within the interrogated sequence, which may explain its larger increase in Tm and formation of its own outlying cluster.    Table 12). Primer sequences are underlined; nucleotide substitutions at the locations shown in red have been documented to cause resistance. Substitutions predicted to increase or decrease amplicon Tms are highlighted in green or yellow, respectively.   Table 13). Primer sequences are underlined; nucleotide substitutions at the locations shown in red have been documented to cause resistance. Substitutions predicted to increase or decrease amplicon T m s are highlighted in green or yellow, respectively.
In parC HRMA assays, roughly half of both B. cereus and B. thuringiensis strains clustered with wild-type B. anthracis controls (cluster 1) and most of the remainders formed a separate cluster (cluster 2). There were two outliers: B. cereus 03BB108 clustered with B. anthracis mutant S3-5 (cluster 3) and B. thuringiensis kurstaki formed its own cluster (cluster 4) ( Table 8). Comparisons of published sequences ( Figure 11) support our hypothesis that strains with sequences identical to the reference should-and do-cluster with B. anthracis wild-type, while those with a T→C change in the interrogated region should have higher T m s and cluster separately (cluster 2; Figure 11). Similarly, B. thuringiensis kurstaki, having a predicted C→T substitution, should produce an amplicon with a lower T m , which we observed here (cluster 4). It is unclear, however, why B. cereus 03BB108 (predicted T→C substitution) produced an amplicon with lower T m and did not fall within cluster 2. For this strain-and all of these near neighbors-it is possible that various substitutions within the primed regions affected the melt curves in ways that could not always be predicted.  Table 12). Primer sequences are underlined; nucleotide substitutions at the locations shown in red have been documented to cause resistance. Substitutions predicted to increase or decrease amplicon Tms are highlighted in green or yellow, respectively.  Table 12 for accession numbers). Primer sequences are underlined; mutations at the locations shown in red have been documented to cause fluoroquinolone resistance in B. anthracis. Substitutions predicted to increase or decrease amplicon Tms are highlighted in yellow or green, respectively.
In parC HRMA assays, roughly half of both B. cereus and B. thuringiensis strains clustered with wild-type B. anthracis controls (cluster 1) and most of the remainders formed a separate cluster (cluster 2). There were two outliers: B. cereus 03BB108 clustered with B. anthracis mutant S3-5 (cluster 3) and B. thuringiensis kurstaki formed its own cluster (cluster 4) ( Table 7). Comparisons of published sequences ( Figure 11) support our hypothesis that strains with sequences identical to the reference should-and do-cluster with B. anthracis wild-type, while those with a TC change in the interrogated region should have higher Tms and cluster separately (cluster 2; Figure 11). Similarly, B. thuringiensis kurstaki, having a predicted CT substitution, should produce an amplicon with a lower Tm, which we observed here (cluster 4). It is unclear, however, why B. cereus 03BB108 (predicted TC substitution) produced an amplicon with lower Tm and did not fall within cluster 2. For this strain-and all of these near neighbors-it is possible that various substitutions within the primed regions affected the melt curves in ways that could not always be predicted.
We   Table 13 for accession numbers). Primer sequences are underlined; mutations at the locations shown in red have been documented to cause fluoroquinolone resistance in B. anthracis. Substitutions predicted to increase or decrease amplicon T m s are highlighted in yellow or green, respectively.
We did not have access to Y. pestis or F. tularensis strains with known mutations causing ciprofloxacin resistance in this study. Therefore, we were unable to assess whether HRMA could be used to discriminate between wild-type and resistant strains within these two genera. Fortunately, we were clearly able to distinguish between wild-type and mutant strains of B. anthracis in a readily observed and relatively predictable fashion.
Furthermore, in spite of numerous failures in parC amplifications, clustering of most near neighbor species for gyrA, parC, and parE in all three genera could be predicted through comparisons of the sequences targeted. F. tularensis could be easily discriminated from F. philomiragia based on gyrA, but not parE melt curves, potentially due to the short (5-nt) parE sequence targeted. Similarly, as expected, all Y. pestis and wild-type B. anthracis clustered with most, but not all, near neighbor strains having identical target sequences, though some anomalous clustering was observed (e.g., Y. pseudotuberculosis strains AMC TB4 and IP32983, B. cereus 03BB108). It is certainly possible that the unpredictable clustering may simply have been due to noisy data [43]. More likely, these anomalous results likely represent contaminants within the samples or errors introduced through use of whole genome amplicons in these HRMA studies. While ϕ29 DNA polymerase is much less prone to errors than other DNA polymerases, any errors introduced in the targeted QRDRs (or even potentially within the primer regions) could result in anomalous melt curves and subsequent clustering. A final possibility is that these anomalies may be due to spontaneous mutations or genetic drift. Assessment of this possibility would require additional amplicon sequencing, which was not performed here.
The ARDM v.3.1 has a total of four subarrays (two pairs of identical subarrays, 2240 probes/ subarray). The subarrays are based on two previous versions-ARDM v.2 (238 AMR determinants) [27,53,54], and ARDM v.3 (542 AMR determinants) [28]-and was re-designed to remove non-relevant or redundant content and emphasize AMR-related content from Select Agents. Each microarray slide could therefore analyze two samples (each sample = two subarrays). Each determinant on the microarray was represented by 8-10 probes, typically corresponding to 4-5 duplicate pairs of probes. Samples were called positive for a gene if signals from half of the gene's probes exceeded the 95% threshold (mean signal from lowest 2128 probes + 3 standard deviations, SD) or if ≥70% of the probes had signals exceeding two lower stringency thresholds [29]. Using a combination of these thresholds has allowed us to detect families of closely related genes using probes with minor mismatches.

High Resolution Melt Analysis (HRMA)
The HRMA protocol, primers, and conditions used here were based on those described by Loveless [34] and are shown in Table 14. The analyses were performed using Type-it HRMA PCR kit (cat#206542, Qiagen, Germantown, MD, USA). The total volume of HRMA reactions was 25 µL and the reaction mixes were prepared according to the manufacturer instructions. Final concentrations of the primers (Eurofins/Operon, Louisville, KY, USA) in the reaction mix were 0.7 µM; template DNA was added to achieve the final concentration of approximately 100 pg/µL. The HRMA PCRs were performed on Rotor-Gene Q MDx (Qiagen, Hilden, Germany) using the following cycling conditions: 5 min at 95 • C, followed by 45 cycles of (10 s at 95 • C, 30 s at 55 • C), ending with melt curve generation (65 • C to 95 • C at 0.1 • C increments, 2 s each). Samples with C t > 35 cycles or without sigmoidal amplification curves were called negative or invalid, respectively. Melt curves of valid amplifications were analyzed and clustering performed using ScreenClust software with unsupervised cluster analysis [55]. Table 14. Primers, probes for HRMA, derived from Loveless [34].

Conclusions
Both PCR and microarrays are valuable tools for the tracking the genetic underpinnings of AMR resistance. Here, we used two complementary technologies-microarray analysis and HRMA-to survey 127 Select Agents, exempt strains, and near-neighbor species for a broad variety of resistance mechanisms acquired through both horizontal transfer and gene mutations. To our knowledge, this is the largest survey of Category A agents, exempt strains, and near-neighbor species for genes covering multiple mechanisms of AMR.
Both of the methods used here had benefits and limitations. The broad detection capabilities of the microarray analysis enabled detection of the genes that we expected to find within the Select Agent species for all three genera. We further detected a non-native gene within a Select Agent strain (KIM (pCD1Ap)+) and detected families of similar genes in near-neighbor species; detection of these related genes was highly dependent on the phylogenetic distance (and hence, sequence conservation) between the near neighbor species and the species from which the reference sequence was derived.
However, each microarray analysis is only as good as its content, and while its content is relatively comprehensive (500+ AMR genes), the ARDM v.3.1 is certainly not all-encompassing and is unable to detect small mutations that can affect specificity. HRMA, on the other hand, was able to rapidly discriminate between ciprofloxacin-resistant B. anthracis strains with point mutations within QRDRs of gyrA and parC, though several strains clustered anomalously (albeit separately from wild-type strains). One limitation of HRMA is that each reaction must be designed for each target sequence and reactions are not easily multiplexable. Furthermore, changes in DNA sequence detectable by HRMA do not always represent a change in amino acid sequence, and near-neighbor species may be erroneously identified as resistant mutants.
Provided suitable sensitivity can be achieved and data processing simplified, this is where rapidly evolving technologies such as next generation sequencing (NGS) may find tremendous application [56]. Several excellent studies have recently demonstrated use of targeted pyrosequencing and next generation sequencing approaches to rapid detect and identify both horizontally and vertically transmitted AMR determinants within Select Agent strains [57][58][59]. While these recent studies emphasized a limited number of AMR genes, we hope that these technologies will prove useful in performing similarly broad-based surveys of Select Agents and other pathogenic species.  Acknowledgments: The authors would like to express their gratitude to Eric I. Thompson, DBPAO, for the many hours spent enabling us to receive and work with these DNA preparations.

Conflicts of Interest:
The authors declare no conflict of interest. DTRA was active in both the research design and project management, and in preparation of the manuscript. This manuscript is approved for public release: distribution is unlimited.