Comparative genomics of Aeromonas veronii: Identification of a pathotype impacting aquaculture globally

Aeromonas veronii is a gram-negative species abundant in aquatic environments that causes disease in humans as well as terrestrial and aquatic animals. In the current study, 41 publicly available A. veronii genomes were compared to investigate distribution of putative virulence genes, global dissemination of pathotypes, and potential mechanisms of virulence. The complete genome of A. veronii strain ML09-123 from an outbreak of motile aeromonas septicemia in farm-raised catfish in the southeastern United States was included. Dissemination of A. veronii strain types was discovered in dispersed geographical locations. Isolate ML09-123 is highly similar to Chinese isolate TH0426, suggesting the two strains have a common origin and may represent a pathotype impacting aquaculture in both countries. Virulence of strain ML09-123 in catfish in a dose-dependent manner was confirmed experimentally. Subsystem category disposition showed the majority of genomes exhibit similar distribution of genomic elements. The type I secretion system (T1SS), type II secretion system (T2SS), type 4 pilus (T4P), and flagellum core elements are conserved in all A. veronii genomes, whereas the type III secretion system (T3SS), type V secretion system (T5SS), type VI secretion system (T6SS), and tight adherence (TAD) system demonstrate variable dispersal. Distribution of mobile elements is dependent on host and geographic origin, suggesting this species has undergone considerable genetic exchange. The data presented here lends insight into the genomic variation of A. veronii and identifies a pathotype impacting aquaculture globally.


Introduction
Aeromonas veronii is a Gram-negative, rod-shaped, mesophilic, motile bacterium in the Aeromonadaceae family. It is widespread in aquatic environments across the globe [1]. Mesophilic creation and de novo assembly were conducted using CLC Workbench and Sequencher version 5.4 (Gene Codes Corporation). For annotation, the draft genome was submitted to RAST [22] and NCBI's Prokaryotic Genome Automatic Annotation Pipeline (PGAAP).

Pan-Core genome
Core-and pan-genome analyses were performed using EDGAR 2.0 (Blom et al. 2016) employing a generic orthology criterion based on BLAST score ratio values (SRV). BLAST scores were normalized in relation to the best hit possible. Based on the distribution of SRVs in the dataset, cutoffs were estimated as described [23]. The number of genes in the core and in the pan-genome were calculated for up to 500 random unique combinations of every number of genomes from 2 to 41. Then the number of core genes or pan-genes was plotted as a function of the number of compared genomes. An exponential decay function (core genome) or a Heaps' power law function (pan-genome) was used to extrapolate the development of the size of the core or pan-genome.

Phylogenetic tree and ANI calculation
To establish taxonomic positions of evaluated A. veronii genomes, a phylogenetic tree was constructed based on the complete core genome of 41 A. veronii genomes (Table 1). Gene sets of the core genome were aligned using MUSCLE [24], resulting in an alignment of 117,137 genes in total. The alignments were concatenated, yielding a multiple alignment of 2,898,457 bp per genome, 118,836,737 bp in total. This concatenated alignment was then used to compute a Kimura distance matrix, and it was used as input for the Neighbor-Joining algorithm as implemented in PHYLP [25]. A second phylogeny was constructed based on this concatenated alignment using the approximate maximum likelihood method of FastTree [26], which automatically generates local support values calculated using the Shimodaira-Hasegawa test [27]. Finally, RaxML [28] was used to construct a maximum likelihood phylogeny. To confirm this result, the optimal model for the dataset was estimated using the ModelFinder of IQ-TREE [29] using the extended model selection option. Based on the results, the GTR model with 10 rate categories for the model of rate heterogeneity was selected. Phylogenetic trees were individually created with RaxML as described above for all 2857 gene sets of the core genome. The resulting 2857 phylogenetic trees were checked for taxa showing outlier phylogenies in a significant number of individual gene trees using the rogue taxa check of RaxML [30]. No such taxa were found. The pylogeny based on individual substitution final tree was computed based on the concatenated alignment of all genes using the described models with partitioning for every codon position. The resulting tree was bootstrapped using the rapid bootstrapping method of RAxML with 100 iterations.

Subsystems coverage
The genomes of ML09-123 and 40 other A. veronii genomes downloaded from NCBI were submitted for annotation to the Rapid Annotations using Subsystems Technology (RAST) for subsystem categorization and comparison [22]. For annotation purposes, the following criteria were chosen for the annotation pipeline: classic RAST for annotation, RAST gene caller for ORF identification, and Figfam (version release70 with automatic fix errors and fix frameshifts options).

Secretion systems
RAST annotated protein files were submitted to MacSyFinder [11,31], which employs systems modeling and similarity searches to identify protein secretion systems and related appendages. Analysis was conducted using the default settings using an ordered/unordered replicon dataset with linear/circular topology. All available protein secretions systems were searched with a maximum e-value of 1.0, a maximum independent e-value of 0.001, and minimal profile coverage of 0.5.

Prophages
The presence of prophages in A. veronii genomes was determined using PHASTER (PHAge Search Tool Enhanced Release) [36]. Nucleotide sequences from all genomes were concatenated to serve as an input file prior to submission to the PHASTER server. Computed results were arranged into three categories: score > 90 was considered intact phage element; a score between 70-90 was deemed questionable; and score < 70 is considered incomplete phage region. Additionally, nucleotide sequences of all the identified phage regions were concatenated and aligned in MAUVE [37] to identify conserved phage regions.

Integron identification
Identification of integrons and their components in the genomes was performed by Integron-Finder [38].

CRISPR (Clustered regularly interspaced short palindromic repeats) and Cas (CRISPR associated genes) elements analysis
CRISPR-Cas systems and their distribution in A. veronii genomes were determined by CRISPRfinder [39] and MacSyFinder [31]. The following criteria were used for the identification of Cas elements: maximal e-value, 1.0; independent e-value, 0.001; minimal profile coverage, 0.5. Three different categories were obtained (mandatory, accessory, forbidden) defined as follows: mandatory elements are identifiable and ubiquitous; accessory elements could be essential but not identifiable. If the identified element is partly homologous, it is considered a forbidden element.

Putative antibiotic resistance-related genes identification
Because the majority of the genomes are not closed, contig files were concatenated, and concatenated nucleotide files were submitted to CARD (Comprehensive Antibiotic Resistance Database) [40] to perform BLASTN search. Results were organized based on an E-value <1 � 10 −10 , and duplicates were removed.

Virulence in catfish
All fish disease challenges were conducted in compliance with protocol #17-288 approved by the Mississippi State University Institutional Animal Care and Use Committee (IACUC). The approved protocol included humane endpoints, and when morbid fish met established criteria, they were immediately euthanized by immersion in tricaine methane sulfonate (MS-222). Criteria for euthanasia were loss of balance, hanging at the water surface, or non-responsiveness to external stimuli. Some of the fish died during the study as a result of the experimentally induced systemic infection due to its rapid progression. All personnel on this experiment received IACUCapproved trained in animal care and welfare by the University Laboratory Animal Veterinarian. The relative virulence of A. hydrophila ML09-119 [21] and A. veronii ML09-123 were tested using a channel catfish challenge model. Briefly, 12 month-old channel catfish (18.2 ± 0.53 cm, 82.3 ± 6.39 g) reared indoors for disease research at Mississippi State University were stocked into twenty-seven 40-liter flow-through tanks (10 fish/tank) and acclimated for four days. Tanks were assigned randomly to two treatment groups (A. veronii and A. hydrophila) and four doses (1 x 10 4 , 1 x 10 5 , 1 x 10 6 , 1 x 10 7 ) with three replicate tanks for each dose. A negative control group was also included. Water temperature was maintained at 32˚C (±2) during the experiments. Fish were fed to satiation twice daily and monitored three times daily for morbidity and mortality. On the day of challenge, fish were anesthetized by immersion in tricaine methane sulfonate (MS-222), and each pre-determined dose was administered intraperitoneally (IP) in 0.1 ml volume. Fish were monitored twice daily, and mortalities were recorded for a total of seven days. The experiment was terminated when no fish mortality was observed for three consecutive days. Mean percent mortalities for each dose and treatment group were calculated and arcsine transformed. The mean transformed percent mortality for fish challenged with four doses for A. hydrophila and A. veronii was compared against sham challenged fish using the Student's t-test (p < 0.05).

Genome features
A. veronii strains used in this study and their genome features are listed in Table 1. G+C ratios range from 58.1-58.9%. Only one strain (AVNIH1) carries a plasmid.

Pan/Core genome analyses
Pan/core genome analysis revealed 8710 total genes in the pan-genome for 41 A. veronii genomes and 2855 genes in the core genome. The extrapolation of the core genome size ( Fig  1) predicted a core genome size of 2791, which is close to the actual core genome of 2857. Extrapolation of the pan-genome showed that it is open with a growth factor γ of 0.240, which is similar to the previously reported 0.260 for A. veronii [41].

Phylogenetic tree and ANI calculation
The phylogenetic relationships of 41 A. veronii genomes were built from the complete core genome of the genomes. The phylogeny was analyzed with the distance matrix-based neighbor joining approach as well as with the maximum likelihood methods FastTree and RaxML (S1 Fig and Fig 2). Results demonstrated A. veronii to be a diverse species with many unresolved clades. Four highly conserved genetic subgroups were identified: 1) U.S. dairy cattle isolates and Greece surface sediment isolates (strain pamvotica), 2) U.S. strain ML09-123 and China strain TH0426 (catfish isolates), 3) U.S. human isolates (strains CECT 4257, CCM 4359, AER 397) and China pond sediment isolate B565, and 4) U.S. surface water and Germany surface water isolates. The remaining A. veronii isolates did not fit into these or other subgroups. These findings were confirmed and supported by Average Nucleotide Identity (ANI); ANI values in conserved branches (genetic subgroups) were above 99.91% (Fig 3).

Subsystems coverage
The subsystems predicted by RAST annotation are listed in Fig 4. SEED subsystem categorization analysis predicted the Chinese isolate strain TH0426 carries the most elements (3559). The most abundant system is "amino acid and derivatives" biosynthesis and utilization. The second most abundant system is "carbohydrates", followed by "cofactors, vitamins, prosthetic groups, and pigments". Interestingly, U.S. catfish isolate strain ML09-123 and Chinese catfish isolate strain TH0426 carry the most abundant elements in "phages, prophages, transposable elements, and plasmids".

Secretion systems, flagellum, TAD, and T4P
All the evaluated genomes encode T1SS, T2SS, T4P, and flagellum core components. Some genes from these systems were absent in some strains, but this could be attributed to errors in annotation or assembly. Strains ML09-123, TH0426, RU31B, X12, and Hm21 do not encode TAD and T5SS components (Fig 5). The majority of A. veronii genomes encode T3SS. Exceptions are Sri Lanka fish isolate Ae52, pond sediment isolate B565 from China, six human isolates (, AVNIH2, AVNIH1, AMC35, CECT4257, CCM 4359, and AER397), and unknown source isolate LMG 13067 from the U.S. (Fig 5). For the T5aSS and T5bSS systems, only 17 dairy cattle isolates from the U.S. and strain pamvotica from Greece encode both systems, whereas the remaining genomes either carry one gene from each system or none (Fig 5). For the T4SS, 13 dairy cattle isolates encode eight of the T4SS-Type T accessory genes, whereas Comparative genomics of Aeromonas veronii four strains (TTU2014-108AME, TTU2014-108ASC, TTU2014-115AME, and TTU2014-115ASC) do not encode this system. Strains CIP107763, AVNIH1, and pamvotica encode some of the T4SS-type T components (S1 Table). Strain AVNIH1 possesses a large plasmid that encodes several T4SS-type F elements. The genomes from strains VBF557 and pamvotica encode some of the T4SS-type F elements (S1 Table).
A. veronii genomes either carry all the mandatory genes of T6SSi or one gene, tssH, from this system. The presence of tssH in the genome of strains that do not have a T6SSi operon is not surprising because tssH is the most frequent T6SSi gene located outside of the T6SSi locus [11], so TssH may contribute to other functions besides T6SSi. T6SSi shows similar distribution as the T3SS with the following exceptions: strains A29, AER39, CB51 and CIP 107763 have T3SS but do not encode T6SSi system. Only 12 out of 41 genomes encode the TAD system (Fig 5).

Insertion elements
All evaluated A. veronii genomes carry insertion elements. The number of insertion elements per genome varies from 16 to 163. Strain AVNIH1 has the highest number of insertion elements, whereas strain A29 carried the least number of insertions. Insertion families and subgroups IS1595-subtype-ISPna2, IS51, IS3, IS2, IS481, IS4, IS903, and ISL3 are encoded by a majority of the evaluated A. veronii genomes (S2 Table).   Table).

CRISPR-Cas elements
The genomes of strains A29, LMG 13067, AVNIH2, AVNIH1, and AMC35 have CRISPR regions. All of the other A. veronii genomes carry questionable CRISPR regions. Cas systems are divided into three categories (Type-I, Type-II, and Type-III), with these three categories further subdivided into 10 different subcategories (Type I-A to F, Type II-A and B, Type III-A and B). Only a handful of A. veronii genomes encode these elements (TH0426, LMG 13067, AVNIH1, and AVNIH2), while the remaining genomes do not carry Cas elements (S4 Table).

Integrons
Only two strains (AVNIH1 and Ae52) encode complete integron regions. The other A. veronii genomes either carry a cluster of attC sites lacking integron-integrases (CALIN) or none of the integron regions (S5 Table).

Antimicrobial resistance genes
Resistome genes and components were categorized into four categories: antibiotic efflux, antibiotic inactivation, antibiotic target alteration, and antibiotic target replacement. All of the Comparative genomics of Aeromonas veronii evaluated A. veronii genomes had relatively similar resistome profiles. The AVNIH1 genome possessed the largest number of antimicrobial resistance elements. Antimicrobial resistance elements encoded by channel catfish isolate ML09-123 and Chinese channel catfish isolate TH0426 were consistent with the profiles from most of the evaluated A. veronii genomes. Antimicrobial inactivation genes imiH, imiS, and cphA derivates were conserved across the 41 A. veronii genomes (Fig 7 and S6 Table).

Discussion
In the current study, we sequenced the genome of A. veronii strain ML09-123 from diseased catfish in the U.S. aquaculture industry. Comparative genomics analysis suggests that ML09-123 and Chinese isolate TH0426 have a common origin. By the injection route of exposure, A. veronii strain ML09-123 is virulent in catfish and is a potential pathogen that could impact U. S. aquaculture. However, experimental infection by injection bypasses the normal mucosal barriers, so this trial does not address the potential necessity for A. veronii to require predisposing conditions to cause infection in catfish under standard aquaculture conditions. ANI calculations can be used for an accurate digital investigation of bacterial systematics [42]. As such, ANI was used to evaluate where the A. veronii ML09-123 genome falls with respect to the publicly available A. veronii genomes. Four discrete groups were identified with a high degree of genomic homology (>99%). Interestingly, members within the groups have disparate geographic origins. Core-genome based phylogenetic analysis was consistent with and supported genetic groupings based on ANI calculations.
Pan/core genome analysis revealed a relatively low proportion of core-genome genes relative to pan-genome genes (30.9%), indicating a substantial amount of gene acquisition in A. veronii clades and strains. However, the core genome size (2272 genes) relative to the average number of protein coding genes for A. veronii (3927 genes) is 70.6%. Thus, the species has a relatively large set of core functions with highly variable sections of genome that potentially enable different environmental adaptations.
Based on SEED subsystem analysis, the Chinese catfish isolate TH0426 and U.S. catfish isolate ML09-123 carry the most abundant elements in the category "phages, prophages, transposable elements and plasmids" (Fig 4). Thus, bacteriophages appear to have contributed significantly to gene acquisition in these two strains.
Secretion systems play an important role in bacterial metabolism and pathogenesis, including host invasion, immune evasion, tissue damage, and bacterial competition [43]. Comparative secretion systems analysis revealed secretion system element distribution correlates well with the four A. veronii genetic subgroups based on ANI and core genome comparison (Fig 5). Previously, epidemic Aeromonas hydrophila genomes'secretion systems were evaluated [44,45], and the presence of each secretion system varied based on the geographic location.

Comparative genomics of Aeromonas veronii
The T1SS is widespread in gram-negative bacteria and is involved in one-step transportation of unfolded substrates directly into extracellular space by bypassing the periplasm [46].  Table. The T1SS is comprised of three proteins that span the cell envelope (outer membrane protein [omp], ATP-binding cassette [ABC], and the membrane fusion protein [mfp]). All evaluated A. veronii genomes possess the T1SS, but the number of mandatory elements varies depending on isolate. Because this system is often responsible for secreting tissue destructive enzymes and toxins, it may contribute to host invasion and competition. Presence of T1SS increases virulence in Serrattia marcescens [47] and Vibrio cholera [48] (Fig 5).
Similarly, all 41 A. veronii genomes encode the T2SS system, which is known for exporting hydrolytic enzymes that are important virulence factors of A. hdyrophila [49,50]. T2SS helps A. veronii colonize the gut of leeches; in particular, export of hemolysin by T2SS is an important step for initial colonization of the leech gut [51] (Fig 5). T2SS is one of the most well-conserved secretion systems within A. veronii, suggesting A. veronii strains rely heavily on the T2SS. T2SS components are similar to Tad pilus and T4P, especially in their structural components [52,53]. T4P can contribute to motility, adhesion, and signaling [54]. All the mandatory genes of T4P are carried by the evaluated A. veronii genomes.
Flagellum is evolutionarily similar to T3SS. We identified that all the A. veronii genomes encode the mandatory genes for flagellum. Interestingly, eleven of the A. veronii genomes do not encode T3SS (Fig 5). Seven of the nine human isolates lack T3SS, suggesting this system may not be critical for mammalian virulence, whereas Silver et al. reported the importance of T3SS in A. veronii strain HM21R virulence in mice [55]. T3SS regulatory network [56] and effector proteins [57] have been investigated in A. hydrophila. In Aeromonas salmonicida, T3SS is linked to immunosuppression in fish [58]. A. veronii uses T3SS to protect itself from leech immune cells [55]. SctN is a component of the T3SS; interestingly, the gene encoding this protein is present in one to four copies in all the evaluated A. veronii genomes. SctN contributes to preparing substrates for export [59] and regulates transfer of effector proteins across host cell membranes [60].
T4SS facilitates exchange of mobile genetic elements between bacteria. In Aeromonas culicicola, T4SS enables conjugative genetic material transfer between bacteria [61]. In our analyses, some A. veronii strains encode components of two different types of T4SS, type T and F (S1 Table). Presence of T4SS in some of the evaluated A. veronii genomes correlates with the number of insertion elements present in corresponding genomes, which suggests that T4SS might be affecting the genome structure of some A. veronii strains by promoting gene uptake. Comparative genomics of Aeromonas veronii T5SS is one of the most widespread, simplest, and abundant systems in bacterial genomes. It is divided into five subtypes (T5aSS-T5eSS) [11], and it utilizes a unique mechanism for secreting substrates, including toxins and receptor proteins [43]. It is known for secreting a large number of proteins, even more than T2SS, and some of its substrates are well-known to contribute to virulence [62,63]. In our analysis, all the cattle isolates and strain pamvotica encode mandatory genes of T5aSS and T5bSS, but the other evaluated genomes carry only one or none of these mandatory genes (Fig 5). Therefore, even though this system is known for contributing to adherence, biofilm formation, virulence, toxin export, and mediating cell-to cell interactions in other pathogens, T5SS may not be important in A. veronii virulence.
There are three type of T6SS: T6SSi, T6SSii, and T6SSiii. The majority of the A. veronii genomes encode the T6SSi system, which is well-known for contributing to virulence and the ability to interact with environment, host, and competitor bacteria [64]. In A. hydrophila, T6SS contributes to virulence and translocates the effector protein Hcp into eukaryotic cells, which can modulate activation of host immune cells [65,66]. Interestingly, only one T6SSi gene (tssH, also known as clpV) is conserved across all the evaluated A. veronii genomes (except strain X11). TssH is an ATPase implicated in the recycling of other T6SS components, TssB/TssC, which are part of the tubular sheath [67,68]. It was proposed that TssH is used to provide energy to contract the tubular sheath [69]. However, because tssH is present in A. veronii genomes that lack genes encoding core T6SS components, TssH may have another important role. In our analysis, we could not identify a pattern for presence of this system based on genetic subgroups, geographic location, or host. However, presence of T6SS in A. hydrophila genomes varies based on the geographic location [70].
The Tad system is essential for colonization, biofilm formation, and virulence of some species [71]. For example, tight-adherence genes are required for virulence of Actinobacillus actinomycetemcomitans [72]. The majority of A. veronii genomes (30 out of 41) do not encode this system. Interestingly, all the genomes carry at least two, some as many as four, copies of tadZ, which encodes one of the core components of TAD system. While functionality of this component remains unclear, it is thought to be peripherally linked with the inner membrane [71].
Aeromonas spp. are important pathogens in humans and other organisms. Surface polysaccharides, iron binding systems, exotoxins and other extracellular enzymes, secretion systems, adhesins, and flagella are some of the main factors contributing to Aeromonas virulence [73]. Interestingly, A. hydrophila and A. veronii-Aeromonas sobria groups contain more virulence genes than the Aeromonas caviae-Aeromonas media group, suggesting that A. hydrophila and A. veronii-A. sobria may be more adapted to a pathogenic lifestyle [74]. A. veronii isolates from diseased Gibel carp have some variation in the virulence genes they carry, but both encode and express several secreted enzymes [75]. Our virulence database search showed that secretion systems and their components are the most prevalent putative virulence elements in the A. veronii genomes, followed by adherence and immune evasion (Fig 6).
Insertion elements are important in bacterial genome evolution, contributing by gene inactivation and structural changes to the genome [76]. Insertion elements can affect virulence by causing expansion of flanking regions, gene inactivation and decay, genome rearrangement and reduction, and incorporation of additional genes [77]. In our analysis, two genetic subgroups (U.S. and China catfish isolates and U.S. cattle and Greece surface sediment isolates) encode similar insertion families, whereas the remaining A. veronii isolates demonstrate a more scattered IS distribution pattern. T4SS contributes to DNA uptake or exchange from other bacteria. In our analysis, strains pamvotica and AVNIH1 encode either all the accessory T4SS genes or some of them, and these two strains encode the largest number of insertion elements (S2 Table). Thus, in A. veronii there appears to be a correlation between number of insertion elements and presence of T4SS.
Bacteriophages can mediate horizontal gene transfer, including virulence genes, from one bacterial strain to another. Specifically, prophage regions can provide specific mechanisms for bacterial attachment, invasion, and survival capabilities [78]. In our panel of isolates, all A. veronii genomes have phage elements except two (strains CECT4486 and CCM7244). U.S. catfish isolate ML09-123 and China catfish isolate TH0426 have different numbers of bacteriophage elements, but they share one intact phage region. We did not identify a conserved phage element present in all 41 A. veronii genomes (S3 Table).
CRISPR-Cas systems consist of two main components, CRISPR array and CRISPR associated genes (Cas), which are separated from each other by spacers [79]. They provide protection against viral predation and foreign DNA invasion (lysogenic bacteriophages, plasmids, or transposons) [80]. Our results showed that all the evaluated A. veronii genomes have either confirmed or questionable CRISPR elements, but only strains TH0426, LMG 13067, AVNIH1, X11 and AVNIH2 encode Cas elements along with CRISPR regions ( Table 2).
Integrons are a major type of genetic element resposible for spread of antibiotic resistance genes. They consist of two main components: 1) integron-integrase (intI) and its promoter region (P intI ) and 2) attachment site of the integron (attI) and constitutive promoter (Pc) that integrates gene cassettes at attI. For a functional integron, only these core elements (intI gene and attI integration site) are required [81]. Gene cassettes typically have an open reading frame (ORF) with flanking attC genes, which mediate integration into integrons. Complete integron elements include integrase and at least one attC region. Incomplete integrons include ln0 elements that consist of integrase without attC sites and clusters of attC lacking integronintegrase (CALIN) [38, [81][82][83]. Strain AVNIH1 has one large plasmid that carries three integrons. On the other hand, strain Ae52 encodes a complete integron in its chromosome. Many of the A. veronii genomes either carry only one, two, or none of the CALIN elements in their genome.
Strain AVNIH1 has more putative antimicrobial resistance genes than other A. veronii genomes we analyzed, which could be the result of encoding three integrons. Our antimicrobial resistance analysis showed a majority of the putative antimicrobial resistance elements are shared by all the A. veronii genomes with some exceptions (Fig 7) (S6 Table).
Virulent A. hydrophila (vAh) is an important cause of MAS outbreaks in catfish aquaculture in the southeast U.S. By the experimental injection route of exposure, A. veronii strain ML09-123 had similar virulence as vAh strain ML09-119. One important difference is that ML09-119 did not cause any mortality at doses of 10 4 and 10 5 CFU/fish, while ML09-123 caused a dose response that correlated with mortalities at the intermediate doses (Fig 8). This suggests the two pathogens use different mechanisms to regulate virulence gene expression. A quorum sensing regulatory mechanism for virulence genes could explain the mortality pattern for vAh. N-acylhomoserine lactone (AHL)-mediated quorum sensing mechanism is important for regulating virulence of A. hydrophila SSU [84]. Virulence of A. hydrophila and A. salmonicida towards burbot (Lota lota L.) larvae is also controlled by quorum sensing [85].
In conclusion, our comparative A. veronii genomes analyses shows that A. veronii strain ML09-123 is highly similar to China strain TH0426. This observation was confirmed by ANI, core genome comparison, and subsystem category distribution (particularly secretion systems, insertion elements, and putative antibiotic resistance genes). Our results confirm that A. veronii is a potential pathogen in Ictalurus catfish and that its virulence by the experimental injection route is similar to vAh. Thus, strains ML09-123 and TH0426 appear to be inherited from a common ancestor strain and affect aquaculture in both China and the U.S. We also found that pan genome is still open and the core genome may reduce as further strains are added to the repertoire. However, approximately 30% of A. veronii genomes show considerable variation, particularly in putative virulence genes, which contributes to a relatively large pangenome for the species.
Supporting information S1