Comparative genomics analysis and virulence-related factors in novel Aliarcobacter faecis and Aliarcobacter lanthieri species identified as potential opportunistic pathogens

Emerging pathogenic bacteria are an increasing threat to public health. Two recently described species of the genus Aliarcobacter, A. faecis and A. lanthieri, isolated from human or livestock feces, are closely related to Aliarcobacter zoonotic pathogens (A. cryaerophilus, A. skirrowii, and A. butzleri). In this study, comparative genomics analysis was carried out to examine the virulence-related, including virulence, antibiotic, and toxin (VAT) factors in the reference strains of A. faecis and A. lanthieri that may enable them to become potentially opportunistic zoonotic pathogens. Our results showed that the genomes of the reference strains of both species have flagella genes (flaA, flaB, flgG, flhA, flhB, fliI, fliP, motA and cheY1) as motility and export apparatus, as well as genes encoding the Twin-arginine translocation (Tat) (tatA, tatB and tatC), type II (pulE and pulF) and III (fliF, fliN and ylqH) secretory pathways, allowing them to secrete proteins into the periplasm and host cells. Invasion and immune evasion genes (ciaB, iamA, mviN, pldA, irgA and fur2) are found in both species, while adherence genes (cadF and cj1349) are only found in A. lanthieri. Acid (clpB), heat (clpA and clpB), osmotic (mviN), and low-iron (irgA and fur2) stress resistance genes were observed in both species, although urease genes were not found in them. In addition, arcB, gyrA and gyrB were found in both species, mutations of which may mediate the resistance to quaternary ammonium compounds (QACs). Furthermore, 11 VAT genes including six virulence (cadF, ciaB, irgA, mviN, pldA, and tlyA), two antibiotic resistance [tet(O) and tet(W)] and three cytolethal distending toxin (cdtA, cdtB, and cdtC) genes were validated with the PCR assays. A. lanthieri tested positive for all 11 VAT genes. By contrast, A. faecis showed positive for ten genes except for cdtB because no PCR assay for this gene was available for this species. The identification of the virulence, antibiotic-resistance, and toxin genes in the genomes of A. faecis and A. lanthieri reference strains through comparative genomics analysis and PCR assays highlighted the potential zoonotic pathogenicity of these two species. However, it is necessary to extend this study to include more clinical and environmental strains to explore inter-species and strain-level genetic variations in virulence-related genes and assess their potential to be opportunistic pathogens for animals and humans.

To date, only the genomes of A. butzleri, A. cibarius, A. cryaerophilus, and A. thereius have been characterized in detail [11,[15][16][17]. These genomes are featured as low GC content (ca. 27%), with the genome sizes ranging from 1.8 to 2.3 Mb [11,15,17]. Comparative genomics further identified several sets of genes or proteins that may be associated with the pathophysiology of pathogenic Aliarcobacter species. Strains of A. butzleri often carry a full or partial set of the nine virulence determinants that are homologous to genes with known pathogenic mechanisms, including the putative virulent factor mviN [18] or genes associated with adherence (cadF, cj1349, hecA and irgA), invasion (ciaB) or destruction (hecB, tlyA, pldA) of host cell walls [5]. Genes or gene clusters involved in the biosynthesis of lipooligosaccharides and flagella, chemotaxis, and antimicrobial resistance have also been identified in A. butzleri and A. thereius [11,15]. Genome analysis combined with laboratory experiments suggested that A. butzleri, A. cryaerophilus, and A. skirrowii may survive in cold and oligotrophic environments, disinfection regimes, food process procedures, and storage conditions [5]. It was reported that the antimicrobial resistance of pathogenic Aliarcobacter strains might be chromosomally determined and associated with the activity of efflux pumps or the presence of degrading enzymes encoded by genes such as cat (chloramphenicol resistance) [5,11,12,15]. It has also shown that mutations in the quinolone-resistance-determining region of gyrA mediate bacterial susceptibility to fluoroquinolones [19,20]. In addition, exotoxins and endotoxins and toxinantitoxin (TA) systems are critical self-defense mechanisms for bacteria that determine a pathogen's capacity and persistency of pathogenicity [21].
Aliarcobacter lanthieri strain AF1440 T , AF1430, and AF1581 were isolated from pig and dairy cattle manure [22], and A. faecis strain AF1078 T was isolated from a human septic tank [23]. A. lanthieri and A. faecis are phylogenetically closely related and clustered with A. cryaerophilus, A. skirrowii, and A. butzleri, based on the phylogenetic analysis using 16S rRNA and housekeeping (gyrB, rpoB, cpn60, gyrA, and atpA) genes; and equipped with short flagellum for mobility [22,23]. Besides, a recent study isolated and identified A. lanthieri strain R-75363 from the stool culture of an immunocompetent patient who developed persistent abdominal bloating and cramps without fever or diarrhea [24]. Therefore, the focus of this study was to assess the virulence-related factors of these two species through comparative genomics analysis. The objectives of this study were to i) perform whole-genome assembly of the reference strains of A. lanthieri and A. faecis; ii) assess the taxonomic position of A. lanthieri and A. faecis based on genome homology; and iii) identify virulence-, antimicrobial resistance-and toxin-related genes in A. lanthieri and A. faecis. This study provided information on the antibiotic resistance, virulence potential, and general fitness of these two new Aliarcobacter species in natural environments. livestock and human fecal sources, were cultured on modified Agarose Medium (m-AAM) (Oxoid) containing selective antibiotic (cefoperazone, amphotericin-B and teicoplanin) supplements. The plates were incubated at 30 °C under microaerophilic conditions (85% N 2 , 10% CO 2, and 5% O 2 ) for 3 to 6 days as described previously [22,23]. Total genomic DNA was extracted and purified using the Wizard Genomic DNA purification kit (Promega, Madison, WI, USA). The concentration of DNA was determined using the Qubit ™ 2.0 Fluorometer (Life Technologies, Burlington, ON, Canada). Purified DNA was stored at − 20 °C for further use.

Genome sequencing
Library preparation and paired-end whole-genome sequencing of A. faecis AF1078 T and A. lanthieri AF1440 T reference strains were performed at the National Research Council Canada (Saskatoon, Saskatchewan, Canada). In brief, high-molecular-weight genomic DNA was used as input for library preparation using the Illumina TruSeq DNA library preparation kit (Illumina Inc.) to obtain a library with a median insert size of 300 bp. After PCR enrichment, the resultant library was checked on a Bioanalyzer (Agilent Technologies Inc., Mississauga, ON, Canada) and quantified. The libraries were equimolarly pooled and sequenced on an Illumina HiSeq 2500 platform (Illumina, San Diego, CA, USA), generating 2 × 101 bp paired-end reads for each sequenced fragment. Base-calling and primary quality assessments were performed using the Illumina CASAVA pipeline (v1.8.2).
The mate-pair sequencing was performed at the Ottawa Research and Development Centre, Agriculture and Agri-Food Canada (Ottawa, Ontario, Canada). The mate-paired DNA library was prepared using the Nextera ® Mate Pair kit (Illumina, San Diego, CA, USA). DNA fragments with three size ranges (1.8-3.5 Kb, 4.0-7.0 Kb, and 8.0-12.0 Kb) were selected using Pippen SageELF (Sage Scientific, Beverly, MA, USA) and pooled with a mean insert size of 6.1 Kb. The libraries were examined on a Bioanalyzer and then quantified using the KAPA qPCR assay (Wilmington, MA, USA). The sequencing libraries were normalized to 2 nM and then diluted to 6 pM prior to loading on a MiSeq Illumina sequencing platform (Illumina, San Diego, CA, USA), which generated 500 bp mate-paired reads for each sequenced fragment.

Genome assembly and annotation
Automatic trimming (based on a threshold of Q = 25) using Trimmomatic scanning and de novo assembly using SPAdes genome assembler version 3.11.1 [25] with combined Illumina NextSeq500 data set of paired-end and mate-pair reads for each species was performed. The contigs were assembled into scaffolds using Medusa [26], where A. nitrofigilis DSM 7299, A. butzleri RM 4018, and A. cryaerophilus L406 genomes were used as scaffolding references. GapFiller [27] closed scaffold gaps using raw paired-end sequencing data. Finally, both de novo assembled genomes were annotated with Prokka version 1.12 [28]. The genome sequences and annotations were deposited to the JGI IMG/MER under analysis IDs Ga0136198 (A. faecis) and Ga0136182 (A. lanthieri) [29].
An Unweighted Pair Group Method with Arithmetic mean (UPGMA) phylogenetic tree was built using the MASH tool version 2.3 with a sketch size of s = 1000, a k-mer size of k = 21, and 100 random seeds for bootstrap [30]. MASH uses the MinHash technique to assess the pairwise distance between sequences [30]. The bootstrapped phylogenetic tree was visualized using the R package ggtree version 3.2.1 [31]. In this analysis, we included the genomes of strains of nine Aliarcobacter species fetched from the National Center for Biotechnol- ogy  Gene synteny and homology of A. faecis strain AF1078 T and A. lanthieri strain AF1440 T were computed using BLASTp and MCScanX with default parameters (match score ≥ 50, E-value ≤ 10 − 5 , max gaps ≤ 25) [32,33].
TA system is a set of genes encoding a pair of stable toxin and unstable anti-toxin. TAfinder was used to predict type II TA loci in A. faecis strain AF1078 T and A. lanthieri strain AF1440 T with default parameters [36].
The VF and TA genes in the complete genome of A. faecis AF1078 T and A. lanthieri AF1440 T strains were visualized using the circlize package in R [37].
The amplicon sizes of each mono-and multiplex PCR reaction were confirmed by 2.5% agarose gel electrophoresis (Fisher Scientific) using a 100 bp DNA size marker (Life Technologies, Grand Island, NY). The agarose gels were stained in ethidium bromide (0.5 μg/mL), and Alpha Imager (Fisher Scientific) was used for scanning and documentation.

General features of a. faecis and A. lanthieri genomes
The genome of A. faecis AF1078 T (= LMG 28519 T ) reference strain contained 2,327,155 bp in one scaffold, and the genome of A. lanthieri AF1440 T (= LMG 28516 T ) reference strain contained 2,234,737 bp in one scaffold. The overall GC contents of the two genomes were 27.0 and 26.4%, respectively, which were consistent with other Arcobacter sensu lato species (Table 1). The numbers of protein-coding genes predicted in A. faecis and A. lanthieri genomes were 2319 and 2230, respectively ( Table 1). A. lanthieri and A. faecis had a similar percentage (73%) of functionally annotated protein-coding genes with Clusters of Orthologous Groups (COGs) ( Table 1). The phylogenetic tree of whole genome comparison shows that the strains of the two new species clustered with A. vitoriensis and A. cibarius, respectively ( Fig. 1). This confirms previous maximum-likelihood phylogenetic analysis based on 16S rRNA and housekeeping genes [22,23]. Table 2 shows the number of predicted genes associated with the COG functional categories (E-value < 10 − 10 ) in the complete reference genomes of Aliarcobacter species. To avoid errors in gene copy number estimation, draft genomes of A. cibarius and A. thereius were not included. Overall, the genes were assigned to 23 out of 26 COG categories in the genomes of all Aliarcobacter species, as shown in Table 2. Genes encoding bacterial Type II (pulE, pulF) and III (fliF, fliN, ylqH) secretory pathways were identified from genomes of all strains (Table S1). Genes encoding bacterial Type IV (virB4), VI (dotU, vasA, vasK, tssA, virG) and VII (hcp) protein secretion systems were identified in the reference strains of A. faecis but not in those of A. lanthieri or A. butzleri included in this study (Table S1). However, a study found a full Type IV pathway for a particular A. butzleri strain D4963 [39], suggesting significant genetic variation between strains isolated from different geographical locations. Besides, genes (tatA, tatB, and tatC) involved in the Twin-arginine translocation (Tat) secretion pathway were found in A. faecis, A. lanthieri and A. butzleri (Table  S1). Furthermore, the mobilome COG category (code X) was underrepresented in the reference genome of A. lanthieri, carrying only one gene associated with COG2932 from that group ( Table 2). A. faecis and A. butzleri had 11 and seven genes, respectively, related to COG code X ( Table 2). The presence of the mobile genetic elements, such as prophages and transposon, may suggest horizontal gene transfer of potentially antimicrobial resistance and/or adaptation genes.

Detection and comparison of virulence-associated genes of Aliarcobacter species
Twenty-six virulence-associated genes were previously reported in A. butzleri and other Aliarcobacter species [15,35]. This study compared the 26 genes and identified 15 in A. faecis and 20 in A. lanthieri (Table 3, Fig. 2). Besides, additional putative VF and TA genes of A. faecis and A. lanthieri were identified using VFanalyzer and TAfinder (Table 4; Fig. 2) against the VFDB and TADB databases, respectively [34,36]. Other known and putative zoonotic pathogens in the genus Aliarcobacter were also annotated using VFanalyzer (Table 4). The E-values of putative VFs were < 10 − 10 , and in general, the coverages were > 90% (Table 4). Here, we present these genes into functional categories, including motility and export apparatus, invasion and stress resistance, adherence, antimicrobial resistance, TA systems, and general resistance.

Motility and export apparatus
Bacterial flagellum can affect its virulence by determining the physical motility and act as a secretion system for other virulence factors [40]. Flagella genes flaA, flaB, flgG, flhA, flhB, fliI, fliP, cheY1, and motA were found in both A. faecis and A. lanthieri, reference strain genomes (Tables 3 and 4). The flagellum apparatus of pathogenic bacteria is considered a secretion system composed of flagellar proteins, which forms a needle to inject bacterial toxins into the host cell. For example, it was reported that the flagellum of H. pylori is required to colonize the mucosal membrane of the stomach as opposed to penetrating the gastric mucosa [41]. Comparative genomics analysis also claimed that some non-flagellum type III secretion systems were evolved from flagellar secretion systems through a series of genetic deletions, innovations, and recruitments of components from other cellular structures [42].

Invasion and stress resistance
Orthologs of virulence factors ciaB, iamA, and mviN were detected in the genomes of both A. faecis and A. lanthieri reference strains. These genes provide pathogens a competitive advantage to survive in the bacterial community (Table 3). Gene ciaB encodes one of the invasion antigens (Cia proteins), deletion of which resulted in significantly attenuated virulence in C. jejuni [43]. In addition, it has been suggested that flagellum serves as an export apparatus or secretion channel for Cia proteins [43]. Studies showed that mutants of Yersinia and C. jejuni without functional flagellar apparatus lack the ability to secrete Cia proteins in comparison to wild type [43]. Another secretion-associated gene, mviN, encodes peptidoglycan (a.k.a murein) flippase. Murein protects the gram-negative bacterial cell membrane from osmotic stress and serves as an anchor for virulence factors [44,45]. The murein layer is vitally important for bacterial cells' survival and is shown in Table 3. It shows that A. faecis and A. lanthieri carry a single copy of mviN ortholog. Similarly, orthologs of iamA and pldA were found in both A. faecis and A. lanthieri genomes (Table 3). Of these, iamA, an invasion-associated marker gene, was also found in C. jejuni and reported to be associated with diarrhea [46,47]. Previous studies showed that the PLA activity in Legionella spp., E. coli and Mycoplasma hyorhinis was associated with the impairment of host intestine cell membranes through hydrolyzation [48]. The lysis property of PLA also helps bacteria to acquire iron from erythrocytes by penetrating the host cell membranes [49]. Orthologs of waaC and waaF were discovered in A. lanthieri but not in A. faecis (Table 3). These two genes were also virulence determinants involved in the biosynthesis of liposaccharide in A. thereius and other species of the family Campylobacteraceae [50]. The orthologs of irgA and fur2 were identified in both A. faecis and A. lanthieri genomes (Table 3). It was previously suggested that irgA, the enterobactin receptor gene, is induced by low iron, and the regulation depends on the iron-responsive master regulator Fur [51]. In addition, irgA ortholog was described for A. butzleri [15] and to a lesser degree to some Campylobacter species [52].
Furthermore, the urease enzyme secreted by bacteria promotes its own persistence in the stomach, allowing them to quickly migrate into the gastric mucosal epithelial line by chemotaxis, where pH is comparatively higher [53]. Although the urease enzyme gene cluster was found along with some accessory genes (ureB, ureD, ureE, and ureG) in A. butzleri (Table 3) [39], it is not identified in A. faecis and A. lanthieri genomes (Table 3).

Adherence
Adherence mechanisms of bacterial pathovars play a major role in invading the hosts and competing with intestinal commensals [54]. Cell surface adhesion encoding genes represented by orthologs of cadF, hecA, and cj1349 have been considered crucial VFs for pathogenic bacteria [55]. Table 3 indicates that only cadF and cj1349 were are detected, whereas hecA was not found in the genome of either species.
Adhesin encoded by cadF was found to mediate binding to Fibronectin, a protein present on the surface of epithelial cells [56]. Studies have also shown that Campylobacter cells lacking cadF exhibited a 50-90% reduction in adherence to epithelial cells [56]. Both proteins, cadF and cj1349c, were important for C. jejuni to adhere to the outer membrane of chicken cells and increase their virulence [56]. Orthologs of hecA/hecB, previously detected exclusively in a few strains [56], were not found in A. faecis and A. lanthieri genomes (Table 3). Although not all known adherence genes were found in the reference genomes of A. faecis and A. lanthieri (Table 3 and  Table 4), these two species may still exhibit adherence ability linking to pathogenicity.   Table 4). The identified pumps belong to the Resistance Nodulation cell Division (RND) protein superfamily, one of the most studied antiporters found in bacteria [57]. The identified genes were highly similar to acrB, encoding a multidrug efflux pump [58]. AcrB is a well-described antiporter involved in resistance to lipophilic β-lactam antibiotics, such as carbapenems and cephalosporins, fluoroquinolones, tetracyclines (including tigecycline), chloramphenicol, macrolides, trimethoprim, ethidium, rifampicin, and novobiocin [58]. It is of particular interest as previous studies showed that A. butzleri strains exhibited resistance to a variety of antibiotics, where the majority of them belong to β-lactams and some to quinolones and coumarins [15]. As indicated above, RND transporters like AcrB may determine resistance to quinolones and coumarins [58,59]. This class of antibiotics targets bacterial DNA gyrase, type II topoisomerase, which plays an essential role in DNA replication [60]. However, significant data accumulated suggests that the resistance to such antibiotics may be acquired through specific mutations in the DNA gyrase gene [19,20]. According to Vickers [20], resistance to novobiocin (coumarin antibiotic) is acquired through two amino acid residue mutations G(80) K and L(140) R in the B-subunit of DNA gyrase (gyrB) gene in Staphylococcus saprophyticus. Alignment of gyrB gene from novobiocin susceptible strain of Staphylococcus saprophyticus with homologs from A. faecis, A. lanthieri, and A. butzleri showed that these Aliarcobacter species carry A(80) and         R(140) residues in gyrB gene (Fig. 3). This may indicate partial resistance to novobiocin due to gyrB mutations. Subunit-A of DNA gyrase (gyrA) may also define resistance to quinolones. According to a previous study, the mutations of two amino acid residues in gyrA, T(83) and D(87) are enough to gain resistance to a variety of quinolones [61]. The alignment of gyrA genes of five Arcobacter strains with its orthologs in Pseudomonas aeruginosa strain ATCC 27853, susceptible to quinolones, showed that residues T(83) and D(87) marked on P. aeruginosa sequence remain intact for most strains including A. butzleri (Fig. 3). On the other hand, A. lanthieri showed Serine at position 83 instead of Threonine, which still indicates susceptibility to quinolones [61]. As shown in previous studies, A. butzleri is susceptible to a high concentration of quinolones, much higher than   [15,62]. Thus, it is suggested that RND transporters are the main contributors to quinolones resistance in Aliarcobacter species, making Gyrase mutations less significant.

Toxin-antitoxin (TA) systems
TA system is a pair of genes encoding a toxin and its cognate anti-toxin, and it helps bacteria withstand lethal antibiotic exposure or environmental stresses [63]. We identified seven TA systems in A. faecis and three in A. lanthieri ( Table 5). The TetR-type transcriptional regulator is located near a gene encoding a major facilitator superfamily (MFS) efflux transporter (Table 5), showing the resistance to disinfectants of quaternary ammonium compounds (QACs), including benzalkonium chloride (BAC) [64]. Also, hipBA TA systems are present in both species (Table 5). The hipB anti-toxin neutralizes the HipA toxin, a serine/threonine kinase inhibiting cell growth where hipBA modules are found in divergent bacterial genomes, and many are related to the persistence of antibiotic resistance [63].
In A. lanthieri, the AraC-type DNA-binding protein, which regulates the expression of the proteins requiring the sugar L-arabinose, is adjacent to a putative acetyltransferase (Table 5) conserved in most environmental mycobacterial species, such as Mycobacterium smegmatis [65].
On the other hand, in A. faecis, the HigB/HigA TA system was found (Table 5), which regulates VFs pyochelin, pyocyanin, swarming, and biofilm formation in Pseudomonas aeruginosa [66]. Besides, the ParDE TA system was also identified. This TA system helps bacteria resist heat and antibiotics [67]. We also found a TA system in A. faecis related to the OmpR family DNA binding response regulator and a putative gene of acyltransferase ( Table 5). The OmpR protein was found to regulate the expression of a type III secretion system at the transcriptional level in Enterohemorrhagic E. coli [68].

General resistance
Conservatively, general resistance factors are not VFs. They determine overall cell stability as part of the housekeeping processes. We identified five chaperone genes, clpA, clpB, groEL, dnaK and EF-Tu, as general resistance factors in A. faecis, A. lanthieri, and A. butzleri, which previously showed a connection to bacterial virulence.
The main function of chaperones is protein folding, and it might determine cell resistance against abiotic stress [69]. In particular, genes clpA and clpB, encoding members of the Hsp100/Clp ATPases family in chaperones, were found necessarily required for intracellular multiplication and heat tolerance [70]. These chaperones, identified as a part of the Clp proteolytic complex, were first reported in E. coli and later identified in other bacteria, such as Staphylococcus aureus [70]. A study of C. jejuni confirmed that clpB acts in acid resistance and stomach transit [71].
In addition, dnaK and groEL are the significant heat shock genes, helping bacteria to overcome stressful environmental conditions, such as heat and acid environments [72]. Of these genes, dnaK, encoding hsp70, assists in the protein folding process through their substrate binding and ATPase domains [72,73], while groEL, encoding hsp60, provides a protected cavity in a double heptameric ring structure for the folding of newly synthesized proteins [72].
Another general resistance VF is the elongation factor TU (EF-tu), the most abundant protein in bacterial cells [74]. EF-tu is a GTP-transferase that catalyzes the binding of aminoacyl-tRNA to the ribosome during the elongation stage of cell growth [75]. Current data shows that EF-tu can be inhibited by aminoglycoside antibiotics, which induces mistranslation and bacterial death [76].

Validation of in silico identified virulence-related genes using PCR assays
Furthermore, we validated the existence of 11 VAT genes, including six virulence (cadF, ciaB, irgA, mviN, pldA, and tlyA), two antibiotic resistance [tet(O) and tet(W)], and three cytolethal distending toxin (cdtA, cdtB, and cdtC of the cdt operon) [77] genes, in A. faecis and A. lanthieri using species-specific PCR-based assays [38]. A. lanthieri tested positive for all 11 VAT genes. By contrast, A. faecis showed positive for ten genes except for cdtB because no PCR assay for this gene was available for this species [38]. However, our comparative genomics analysis identified all three cdt (cdtA, cdtB, and cdtC) genes in the reference genomes of A. faecis and A. lanthieri strains. To validate our detection of the cdtB in A. faecis AF1078 T genome, we aligned the cdtB gene of A. lanthieri (UnitProt ID: A0A2K9Y5C5) against the protein sequences of A. faecis strain AF1078 T using BLASTp and identified gene 2,690,353,140 as the cdtB gene of A. faecis (identity 78%; E-value = 0) ( Supplementary Fig. S1). Similarly, Campylobacter spp. also showed variable frequency of the cdt genes [78].
Of the Cytolethal Distending Toxin encoded by the cdt operon [77], cdtB is the active subunit, while cdtA and cdtC work as two regulatory subunits that bind to cdtB [77]. The presence of all three genes of the cdt operon may indicate that the A. faecis strain AF1078 T and A. lanthieri strain AF1440 T could potentially be pathogenic; therefore, further in vitro research is warranted to investigate risk assessment analysis associated with human and animal health. In contrary to these results, studies have shown the absence of the cdt genes in A. butzleri [15,79].
In summary, the results of our PCR assays are in congruence with previous studies where a high frequency of cadF, ciaB, mviN, pldA, and tlyA virulence genes was reported in A. butzleri and A. skirrowii strains [80,81]. Similarly, tet(O) and tet(W) antibiotic resistance genes were also detected in both species, which has also been reported in A. cryaerophilus [15]. Our findings indicate that tetracycline resistance is prevalent in the genus Aliarcobacter.

Conclusion
This study provided insights into the virulence-related factors identified in the reference genomes of two new Aliarcobacter species, A. faecis and A. lanthieri, using whole genome sequencing, comparative genomics analysis, and qPCR validation. Our results generally showed genes encoding motility and export apparatus, secretory pathways, abiotic stress resistance, and antimicrobial resistance were found in both A. faecis and A. lanthieri. However, unique genes were also identified for individual species. We acknowledge that further in vitro and in vivo assays are required to evaluate the roles of virulencerelated factors in the pathogenicity of A. faecis and A. lanthieri in human and animal infections.