Comparative Pathogenomics of Aeromonas veronii from Pigs in South Africa: Dominance of the Novel ST657 Clone

The pathogenomics of carbapenem-resistant Aeromonas veronii (A. veronii) isolates recovered from pigs in KwaZulu-Natal, South Africa, was explored by whole genome sequencing on the Illumina MiSeq platform. Genomic functional annotation revealed a vast array of similar central networks (metabolic, cellular, and biochemical). The pan-genome analysis showed that the isolates formed a total of 4349 orthologous gene clusters, 4296 of which were shared; no unique clusters were observed. All the isolates had similar resistance phenotypes, which corroborated their chromosomally mediated resistome (blaCPHA3 and blaOXA-12) and belonged to a novel sequence type, ST657 (a satellite clone). Isolates in the same sub-clades clustered according to their clonal lineages and host. Mobilome analysis revealed the presence of chromosome-borne insertion sequence families. The estimated pathogenicity score (Pscore ≈ 0.60) indicated their potential pathogenicity in humans. Furthermore, these isolates carried several virulence factors (adherence factors, toxins, and immune evasion), in different permutations and combinations, indicating a differential ability to establish infection. Phylogenomic and metadata analyses revealed a predilection for water environments and aquatic animals, with more recent reports in humans and food animals across geographies, making A. veronii a potential One Health indicator bacterium.


Introduction
Aeromonas spp. are Gram-negative, rod-shaped, non-sporulating, non-acid-fast, and facultative anaerobic bacteria that have been recognized as emerging nosocomial pathogens [1]. They belong to the family Aeromonadaceae, class Gammaproteobacteria that encompasses three genera, viz., Tolumonas, Oceanimonas and Aeromonas [2,3]. Aeromonas share many similar biochemical characteristics with Enterobacterales but are easily differentiated, with Aeromonas being oxidase-positive [4]. The Aeromonas genus has a complex, dynamic taxonomy due to the expanding number of new species and its high intraand interspecies genetic variability [5,6]. This genus currently comprises 36 species, with A. dhakensis,

Assessment of Pathogenic Potential and Prediction of Putative Virulence Factors
The pathogenic potential of the A. veronii isolates was assessed using the PathogenFinder web service under the automated mode [39]. All isolates were subjected to the pathogenicity prediction using Fasta formatted genome data. Virulence determinants (sequences and functions), corresponding to four major bacterial virulence factors (adherence, motility, secretion system and toxin) associated with A. veronii, were searched for in the pathogenic bacteria database, VFanalyzer [40] and validated by blasting assembled genomes to a pseudomolecule generated by concatenating a set of target genes using the NCBI in-house BLASTN tool. The known A. veronii B565 (4,551,783 bp, Accession number: NC_015424) was used as the reference genome.

Global Phylogenomic Relationship and Metadata Analysis
A phylogenomic analysis was undertaken to compare the genomes of the study isolates with all available A. veronii genomes downloaded from GenBank and PATRIC via CSI Phylogeny-1.4 [41]. The genome of Tolumonas auensis DSM 9187 (GenBank accession number: CP001616) of the Aeromonadaceae family was used to root the tree, facilitating the configuration of the phylogenetic distance between the isolates on the branches. A bootstrapped indicator with 100 replicates was applied to identify recombined regions and provide the phylogenetic accuracy in groups with little homoplasy. Figtree (http://tree.bio.ed.ac.uk/software/figtree/) was used to visualize and edit the phylogenetic tree. The phylogeny was visualized alongside annotations for isolate metadata (WGS in silico molecular typing, source, and country) using Phandango [42], to provide a comprehensive analysis of the generated phylogenomic tree.

Data Availability
The raw read sequences and assembled whole-genome contigs have been deposited in GenBank. The data is available under project number PRJNA564235.

Prevalence and Phenotypic Resistance Profiles of A. Veronii Isolates
Five (5) carbapenem-resistant (CR) A. veronii isolates were obtained from 5 different pig samples, giving a carriage rate of 1.5% (5/345) among the sampled pig population. The isolates exhibited similar resistance phenotypes, all of them displaying 100% resistance to ampicillin, amoxicillin and imipenem, except for piperacillin-tazobactam to which only one isolate showed resistance (Table 1).

Genome Confirmation, Statistics, Annotation, and Visualization
The global Pathogenwatch platform confirmed the phenotypic identity of the A. veronii genomes. The genomic features of the sequences, in terms of size, GC content, number of contigs, N50, L50, number of RNAs and number of protein-coding sequences are listed in Supplementary Materials Table S1. The genome size of the A. veronii isolates ranged from 4.64 Mb to 4.77 Mb, with GC content of 58.2-58.5 and a coverage range between 99% to 102%. Comparative visualization analysis, via the GView server (Figure 1), showed a similarity of DNA synteny with >98% coverage and identity from the AVNIH2 reference in all the A. veronii genomes.   (Figure 1), showed a similarity of DNA synteny with >98% coverage and identity from the AVNIH2 reference in all the A. veronii genomes.

Figure 1.
Comparative visualization of the isolates (n = 5) with the A. veronii reference strain (AVNIH2, Accession number: LRBO00000000). The map was constructed using the GView online server (https://server.gview.ca/). The concentric circles represent comparisons between AVNIH2 and, starting with the inner circle, genome assemblies from A. veronii genomes (isolate ID: A5, A31, A34, A86 and A136). Colour codes are given for each isolate with a synteny identity, ranging from >98%.
Pan-genome and ortholog analysis revealed that all the isolates formed a total of 4349 clusters and shared 4296 orthologous clusters (98.8%) ( Figure 2). The singletons ranged from 10 to 42 gene clusters. The isolates shared a similar range of orthologous clusters (n range = 26-30) with each other; however, no unique orthologue gene cluster (n = 0) was found ( Figure 3).

Defence Systems (CRISPRCas Cluster and Restriction-Modification (R-M) System), Antibiotic Resistome, Mobilome and Genetic Environment Analysis
The isolates shared a similar resistome (same resistance and efflux genes), with the mobilome comprised of variable insertion sequence families (Table 2 and Figure 4), while plasmids, integrons and intact prophages were absent. However, isolates contained the same incomplete prophage (PHAGE_Escher_PA28 [Accession number: NC_041935]). The bla CPHA3 and bla OXA-12 genes were located on the chromosome, with >97% homology to Aeromonas veronii strain AVNIH1 (Accession no.  (Table S3). The isolates contained 4-5 clustered, regularly interspaced, short palindromic repeat (CRISPR) arrays with no Cas element ( Table 2). None of the strains harboured the restriction-modification (R-M) system.

Molecular Typing, Global Epidemiological and Phylogenomic Analysis
In silico determination of the sequence types (STs) of the isolates, using the Aeromonas MLST scheme resulted in an undescribed ST comprising new alleles for gltA_473, groL_445, gyrB_465, metG_470, ppsA_512 and recA_519. Allele sequences were submitted for curation and assigned to the  (Table 2). BURST (Based Upon Related Sequence Types) analysis identified the novel ST657 as a satellite clone (more distantly) with no single-, double nor triple-locus variants of global curated A. veronii STs.
The metadata analyses of the five isolates, together with the 49 global strains, showed the diversity of A. veronii sequence types (n = 29), country of origin (n = 13) and sources. The most prevalent sequence types were ST23 (n = 9), ST657 (n = 5, study isolates), ST91 (n = 4), ST166 (n = 3), ST485 (n = 3) and ST50 (n = 2) ( Figure 5 and Table S4). The A. veronii isolates were from different sources, predominantly animal hosts (n = 38; mostly from fish sources [n = 25], pigs [n = 6] and cow [n = 4]), followed by humans (n = 10) and the environment (n = 6). The USA (n = 17), China (n = 10), Greece (n = 9), South Africa (n = 5) and India were the countries with the highest number of isolates deposited on the databases.  The phylogenomic analyses via the WGS SNP tree differentiated the global strains into three clusters (I, II and III) (Figures 5 and 6). The study isolates were found in Cluster I with a 100% branch conservation and close lineage to two international stains, KLG7 (UK) and A8-AHP (India) from the environment and fish source, respectively. Cluster II contained strains (n = 14) which were mostly of animal origin (n = 11), except for three strains which were from humans (AVNIH2, VBF557 and FC951). Cluster III was the largest group and contained strains from diverse sources (humans, animals, and environment) (Figures 5 and 6). Six highly conserved genetic subclades (A-E) were identified, which depicted a clustering of isolates mainly according to their sequence types (clonal lineages) and sources ( Figure 6). Moreover, the results of the global phylogenetic tree demonstrated the complexity and diversity of A. veronii regarding geography, source, and clonality, with many unresolved clusters. Figure 6. A circular cladogram of global A. veronii genomes depicting the association between isolates in three clusters (I, II and III). The KwaZulu-Natal (South Africa) isolates belonged to the smallest cluster I and were mainly related to international strains KLG7 (UK) and A8-AHP (India). Sub-clades (A-E) depicted a clustering of isolates, mainly according to sequence types/sources.

Pathogenic Potential and Putative Virulence Factors
The mean pathogenicity score of 0.60 indicated the potential pathogenicity of A. veronii in humans and was found to match 30 pathogenic families. The whole virulome analysis predicted a total of 200 putative virulence-encoding genes belonging to six major virulence factor classes of Aeromonas, namely, adherence factor (lateral flagella, mannose-sensitive hemagglutinin pilus, polar flagella, tap type IV pili and type I fimbriae), secretion system (T2SS and T3SS), toxins (aerolysin/cytotoxic enterotoxin and hemolysin), anaerobic respiration (nitrate reductase), antiphagoctyosis (capsular polysaccharide) and immune evasion (capsule and LOS) with minor Figure 6. A circular cladogram of global A. veronii genomes depicting the association between isolates in three clusters (I, II and III). The KwaZulu-Natal (South Africa) isolates belonged to the smallest cluster I and were mainly related to international strains KLG7 (UK) and A8-AHP (India). Sub-clades (A-E) depicted a clustering of isolates, mainly according to sequence types/sources.

Pathogenic Potential and Putative Virulence Factors
The mean pathogenicity score of 0.60 indicated the potential pathogenicity of A. veronii in humans and was found to match 30 pathogenic families. The whole virulome analysis predicted a total of 200 putative virulence-encoding genes belonging to six major virulence factor classes of Aeromonas, namely, adherence factor (lateral flagella, mannose-sensitive hemagglutinin pilus, polar flagella, tap type IV pili and type I fimbriae), secretion system (T2SS and T3SS), toxins (aerolysin/cytotoxic enterotoxin and hemolysin), anaerobic respiration (nitrate reductase), antiphagoctyosis (capsular polysaccharide) and immune evasion (capsule and LOS) with minor differences (Figure 7 and Table S5). A total of 195 conserved virulence factors were observed across the isolates. differences (Figure 7 and Tables S5). A total of 195 conserved virulence factors were observed across the isolates.

Discussion
Aeromonas spp. are important human pathogens and colonizers and have also been increasingly isolated in animals (food-animals, wildlife, and companion animals), and the environment (soil, water, and air) globally. They, thus, have potential as One Health indicator bacteria for monitoring the spread of antibiotic-resistant bacteria between humans, animals and the environment [4,[43][44][45][46][47]. However, information on this pathogen in food animals, using high-throughput technologies such as whole genome sequencing in Africa, is lacking. In this study, we describe for the first time in Africa, the comparative genomics of five A. veronii isolates recovered from pigs in South Africa. We also show the phylogenomic relationship between this novel strain and all globally deposited A. veronii genomes with complete metadata (country, sources and sequence types), as their incidence and geographic spread is vital to understand the evolution of this emerging pathogen which is on a global rise.
Analyses of the genomic data revealed a high degree of genomic synteny (>98%), suggesting a close association between the study isolates ( Figure 1). All the isolates shared orthologous clusters (98.8%) with no unique gene cluster (Figure 2), indicating a relatively large set of core functions with low variable sections as well as a vast array of similar central networks (Figure 3) which are crucial for their survival in the microbial community [48]. The high degree of similarity between the isolates predicted by the different analyses also corroborated the novel clonal lineage (ST657), where the isolates possessed the same genetic make-up with low variation. For instance, mobilome analysis of the ST657 highlighted the lack of plasmids, integrons and intact prophages in this lineage. However,

Discussion
Aeromonas spp. are important human pathogens and colonizers and have also been increasingly isolated in animals (food-animals, wildlife, and companion animals), and the environment (soil, water, and air) globally. They, thus, have potential as One Health indicator bacteria for monitoring the spread of antibiotic-resistant bacteria between humans, animals and the environment [4,[43][44][45][46][47]. However, information on this pathogen in food animals, using high-throughput technologies such as whole genome sequencing in Africa, is lacking. In this study, we describe for the first time in Africa, the comparative genomics of five A. veronii isolates recovered from pigs in South Africa. We also show the phylogenomic relationship between this novel strain and all globally deposited A. veronii genomes with complete metadata (country, sources and sequence types), as their incidence and geographic spread is vital to understand the evolution of this emerging pathogen which is on a global rise.
Analyses of the genomic data revealed a high degree of genomic synteny (>98%), suggesting a close association between the study isolates ( Figure 1). All the isolates shared orthologous clusters (98.8%) with no unique gene cluster (Figure 2), indicating a relatively large set of core functions with low variable sections as well as a vast array of similar central networks (Figure 3) which are crucial for their survival in the microbial community [48]. The high degree of similarity between the isolates predicted by the different analyses also corroborated the novel clonal lineage (ST657), where the isolates possessed the same genetic make-up with low variation. For instance, mobilome analysis of the ST657 highlighted the lack of plasmids, integrons and intact prophages in this lineage. However, variability in chromosome-borne insertion sequence (IS) families was observed (Table 2 and Figure 4). More so, the isolates contained the same incomplete prophage (Escher_PA28) which did not harbor any resistance and virulence genes. A similar scattered IS pattern was previous reported in A. veronii by Tekedar et al. [20]. Further analyses of this clonal lineage (ST657) depicted a unique satellite-variant, implying that it was distantly related to global sequence types, hence does not have any ancestral linkage with the STs found in the Aeromonas MLST database.
Bacteria often accrue defence systems which offer protection against foreign DNA invasion and viral predation [49,50]. Regarding the CRISPR-Cas defence systems consisting of two main components, CRISPR array and associated genes (Cas), the A. veronii genomes encoded the CRISPR elements with no Cas system (Table 2). This shows that the genomes of the isolates contained short repeat clusters which have been implicated in bacterial adaptation strategies, ranging from immune evasion and tissue tropism to the modulation of environmental stress tolerance. This finding was similar to a study by Tekedar et al. [20], which confirmed the presence of CRISPR elements in all compared A. veronii isolates (n = 41), with only a few (n = 4) harbouring the Cas elements. This implies that the complete CRISPR-Cas systems are less prevalent in A. veronii, probably because of the lack of nucleotide biosynthesis capacity [51]. Moreover, the Cas systems were predominantly found in human isolates, suggesting the need for further studies to understand the CRISPR-Cas-mediated host interactions.
A. veronii has been reported in different food animals and products including pigs, chicken, cattle, sheep, buffalo and fish [10,46,52]. The consumption of undercooked/raw meat or meat products is an important route for human infection with Aeromonas spp. [4,10]. The carriage rate of carbapenem-resistant A. veronii isolates was 1.5%, which was comparable to the prevalence rate in food-producing animals in Europe (<1%) and in the lower range of the resistance reported in both Africa (2-26%) and Asia (1-15%) [53]. Although the overall prevalence of carbapenem-resistant A. veronii in food animals appears to be low, the transmission of these pathogens from food animals to their derived products could be a threat to consumers, supporting the transmission of resistant bacteria and their determinants between commensal and pathogenic microorganisms with unknown, but potentially severe, consequences for human health [53,54].
The resistance phenotypes corroborated the presence of bla CPHA3 and bla OXA-12 conferring resistance to imipenem and penicillin (ampicillin and amoxicillin) (Tables 1 and 2). The Aeromonads, including veronii and A. hydrophilia, have been reported to harbour conserved resistance genes on their chromosome, conferring intrinsic resistance against these antibiotics [20,55,56]. More so, they often harbour genes that code for the production of β-lactamases such as class B metallo-β-lactamase, class C cephalosporinase, and class D penicillinases [3,5,41,42]. However, Aeromonas spp. are reported to be susceptible to monobactams, third-and fourth-generation cephalosporins, aminoglycosides, and fluoroquinolones, as found in this study [3].
The pathogenic potential (P score ), with the probability ranging from 0 to 1, is used to predict the ability of bacteria to cause infection in humans [39,48]. This theoretical estimation of the pathogenic potential, using trained algorithms to differentiate between pathogenic or commensal strains, predicted a relatively higher average probability (P score ≈ 0.60), suggesting that the clone (ST657) ( Table 2) could be potentially pathogenic to humans. However, it originated from a non-human source, highlighting the role of Aeromonas spp. from food animal sources as potential human pathogens [4]. The global emergence of Aeromonas spp. in all One Health settings (human-animal-environment), makes it a potential One Health indicator organism.
Several insertion sequences were found in our isolates (Table 2 and Figure 4). Insertion elements are significant in the evolution of Aeromonas genome [57], contributing to its resistome and virulome through the incorporation of additional genes, genome reduction and rearrangement, gene decay and inactivation, and expansion of flanking regions [58,59]. The virulome analysis revealed the possession of a battery of determinants which play a significant role in their survival and pathogenesis, comparable to previously reported studies on Aeromonas spp. [5,20,55,60] and supporting the pathogenic potential of this pathogen. The ST657 clone contained an array of six putative virulence factor classes, which were mostly conserved across the isolates, suggesting that A. veronii relies heavily on these factors for host invasion, immune evasion, tissue damage, and competition in diverse ecological niches (Figure 7 and Table S5). Adherence factors were the most prevalent putative virulence determinants followed by secretion systems and toxins (Figure 7 and Table S5) in contrast to virulence factors possessed by A. veronii isolated from fish samples, where the secretion system and its components were found to be the most predominant [20]. This observation could imply that predominating virulence factors may be host-specific. Interestingly, differences in the virulome of the ST657 lineage were evident. Some of the isolates lacked specific genes within the sub-components of adherence factors, secretion systems and immune evasion virulence factors. For example, the T2SS sub-component of the secretion system, which is known for exporting hydrolytic enzymes and aids in the gut colonization [61,62], lacked the exeH gene in isolate A31. The expression of this putative virulome probably confers a competitive advantage, contributing to its success as a pathogen [48]. Moreover, the genomic detection of these virulence genes could aid in identifying targets for the development of novel vaccines for this emerging pathogen [63,64].
Global epidemiological comparison of deposited A. veronii genomes revealed the diverse nature of this pathogen regarding its host, clonal lineages, and geographical distribution ( Figure 5 and Table S4). This diversity implies that the A. veronii can serve as good One Health indicator pathogen to understand and track the geographical spread of antibiotic resistance. Phylogenomic analysis depicted the clustering of group members from disparate geographies. Interestingly, the study isolates were closely related to strains from the UK and India, from different lineages (Figure 7), but not close enough to suggest import into South Africa from other countries. The small number of global deposited strains with insufficient metadata made it challenging to make much inference from the tree analysis on the transmission dynamics of this species, as there were many unresolved clusters. It is thus recommended that more studies be carried out in all the sectors (human, animal and environment) to harness the ability of genomics and bioinformatic analysis in making useful predictions about the dynamics of emerging pathogens in the One Health context.

Conclusions
The comparative genomics of A. veronii revealed the clonal dominance of the novel strain, ST657, isolated from South Africa. The genomic data presented lends useful insights into the pan-genome, resistome, defense system, virulome, pathogenic potential, clonal lineages, global dissemination, and phylogenetic relationship of this pathogen. To the best of our knowledge, this is the first comprehensive genomic analysis of A. veronii isolates in Africa and presents this species as a potential One Health indicator.