*Gupta RS, Lo B and Son J (2019) Corrigendum: Phylogenomics and Comparative Genomic Studies Robustly Support Division of the Genus Mycobacterium into an Emended Genus Mycobacterium and Four Novel Genera. Front. Microbiol. 10:714. doi: 10.3389/fmicb.2019.00714* **Phylogenomics and Comparative Genomic Studies Robustly Support Division of the Genus Mycobacterium into an Emended Genus Mycobacterium and Four Novel Genera** by Gupta, R. S., Lo, B., and Son, J. (2018). Front. Microbiol. 9:67. doi: 10.3389/fmicb.2018.00067

In the original article, there was an error. Based on the branching position of Mycobacterium vulneris (van Ingen et al., 2009) in different phylogenomic trees and on multiple identified molecular signatures that this species shared with a clade of rapid growing mycobacteria, we proposed a reclassification of M. vulneris, into a new genus, Mycolicibacterium, corresponding to a clade of rapid-growing mycobacteria. However, it was noted in our article that the branching of M. vulneris, which is a slow-growing species with rapid-growing mycobacteria, was anomalous.

In a Frontiers commentary, Tortoli (2018) indicated that the genome sequence of M. vulneris, originally available in the NCBI genome database (accession CCBG00000000; Croce et al., 2014), was mislabeled and very likely corresponded to Mycobacterium porcinum (a rapid grower). Tortoli (2018) also reported the sequencing of the type strain of M. vulneris, DSM 45247T and this genome sequence (accession NCXM01000000) showed the branching of M. vulneris within the slow-growing group of mycobacteria, belonging to the genus Mycobacterium.

Our own analysis with this new genome sequence also confirms the branching of M. vulneris within the delimited genus Mycobacterium, encompassing different slow-growing mycobacteria. As a result, the transfer of M. vulneris into the genus Mycolicibacterium as proposed in Table 11 of our article was incorrect as a direct result of the mislabeling of the available genome sequence for this species. To correct this error, we propose that the species Mycolicibacterium vulneris (Gupta et al., 2018) should be reinstated to its previous basonym Mycobacterium vulneris (van Ingen et al., 2009) and as part of the genus Mycobacterium (Gupta et al., 2018).

The authors apologize for this error and state that this does not change the scientific conclusions of the article in any way.

# REFERENCES

Copyright © 2019 Gupta, Lo and Son. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# Revisiting the Taxonomy of the Genus Arcobacter: Getting Order From the Chaos

Alba Pérez-Cataluña1 , Nuria Salas-Massó1 , Ana L. Diéguez2 , Sabela Balboa2 , Alberto Lema2 , Jesús L. Romalde2 \* and Maria J. Figueras1 \*

1 Departament de Ciències Mèdiques Bàsiques, Facultat de Medicina, Institut d'Investigació Sanitària Pere Virgili, Universitat Rovira i Virgili, Reus, Spain, 2 Departamento de Microbiología y Parasitología, CIBUS-Facultad de Biología, Universidade de Santiago de Compostela, Santiago de Compostela, Spain

#### Edited by:

Martha E. Trujillo, Universidad de Salamanca, Spain

#### Reviewed by:

John Phillip Bowman, University of Tasmania, Australia Javier Pascual, Deutsche Sammlung von Mikroorganismen und Zellkulturen (DSMZ), Germany

#### \*Correspondence:

Jesus L. Romalde<EMAIL_ADDRESS>María J. Figueras<EMAIL_ADDRESS>
#### Specialty section:

This article was submitted to Evolutionary and Genomic Microbiology, a section of the journal Frontiers in Microbiology

Received: 15 February 2018 Accepted: 14 August 2018 Published: 04 September 2018

#### Citation:

Pérez-Cataluña A, Salas-Massó N, Diéguez AL, Balboa S, Lema A, Romalde JL and Figueras MJ (2018) Revisiting the Taxonomy of the Genus Arcobacter: Getting Order From the Chaos. Front. Microbiol. 9:2077. doi: 10.3389/fmicb.2018.02077 Since the description of the genus Arcobacter in 1991, a total of 27 species have been described, although some species have shown 16S rRNA similarities below 95%, which is the cut-off that usually separates species that belong to different genera. The objective of the present study was to reassess the taxonomy of the genus Arcobacter using information derived from the core genome (286 genes), a Multilocus Sequence Analysis (MLSA) with 13 housekeeping genes, as well as different genomic indexes like Average Nucleotide Identity (ANI), in silico DNA–DNA hybridization (isDDH), Average Amino-acid Identity (AAI), Percentage of Conserved Proteins (POCPs), and Relative Synonymous Codon Usage (RSCU). The study included a total of 39 strains that represent all the 27 species included in the genus Arcobacter together with 13 strains that are potentially new species, and the analysis of 57 genomes. The different phylogenetic analyses showed that the Arcobacter species grouped into four clusters. In addition, A. lekithochrous and the candidatus species 'A. aquaticus' appeared, as did A. nitrofigilis, the type species of the genus, in separate branches. Furthermore, the genomic indices ANI and isDDH not only confirmed that all the species were welldefined, but also the coherence of the clusters. The AAI and POCP values showed intra-cluster ranges above the respective cut-off values of 60% and 50% described for species belonging to the same genus. Phenotypic analysis showed that certain test combinations could allow the differentiation of the four clusters and the three orphan species established by the phylogenetic and genomic analyses. The origin of the strains showed that each of the clusters embraced species recovered from a common or related environment. The results obtained enable the division of the current genus Arcobacter in at least seven different genera, for which the names Arcobacter, Aliiarcobacter gen. nov., Pseudoarcobacter gen. nov., Haloarcobacter gen. nov., Malacobacter gen. nov., Poseidonibacter gen. nov., and Candidate 'Arcomarinus' gen. nov. are proposed.

Keywords: Arcobacter, Aliiarcobacter gen. nov., Pseudoarcobacter gen. nov., Haloarcobacter gen. nov., Malacobacter gen. nov., Poseidonibacter gen. nov., taxonomic criteria

# INTRODUCTION

fmicb-09-02077 September 1, 2018 Time: 10:25 # 2

The genus Arcobacter was created by Vandamme et al. (1991) to accommodate Gram-negative, curved-shaped bacteria belonging to two species Campylobacter cryaerophila (now Arcobacter cryaerophilus) and Campylobacter nitrofigilis (now A. nitrofigilis), considered atypical campylobacters due to their ability to grow at lower temperatures (15◦C–30◦C) and without microaerophilic conditions (Vandamme et al., 1991). The latter species was selected as the type species for the new genus (Vandamme et al., 1991). One year later the genus was enlarged with the addition of two new species, A. skirrowii with an animal origin being isolated from aborted ovine, porcine and bovine fetuses, and from lambs with diarrhea, and A. butzleri, which was recovered from cases of human and animal diarrhea (Vandamme et al., 1992). Another two new species were incorporated into the genus in 2005. A. halophilus was isolated from water from a hypersaline lagoon in Hawaii (Donachie et al., 2005), and A. cibarius was isolated from broiled carcasses in Belgium (Houf et al., 2005). These species were assigned to the genus Arcobacter on the basis of the 16S rRNA gene similarity (94% and 95% for A. nitrofigilis with A. halophilus and A. cibarius, respectively). However, these values are equal, or even below, the cut-off of 95% for genus definition (Rosselló-Mora and Amann, 2001; Yarza et al., 2008, 2014; Tindall et al., 2010).

From 2009 onward, new species were being described yearby-year, reaching a total number of 27 in 2017. In some of these descriptions, the similarity of the 16S rRNA gene was the decisive character for taxonomic assignation at genus level, although phylogeny based on housekeeping genes (rpoB first and then gyrB and hsp60) was also included as additional, more discriminatory tools for the species (Collado et al., 2009a, 2011; De Smet et al., 2011). Using this approach, A. molluscorum, A. ellisii, A. defluvii, or A. bivalviorum were defined, among others (Collado et al., 2009a, 2011; Figueras et al., 2011a,b; Levican et al., 2012), which showed 16S rRNA similarities ranging from 91.1 to 94.7%, not supporting their common affiliation. On the other hand, the most closely related species, which showed a similarity of 99.1% were A. ellisii and A. defluvii (Collado et al., 2011), giving evidence for the first time of the poor resolution of the 16S rRNA gene for separating closely related species in the genus Arcobacter. However, the phylogenetic analysis based on the concatenated sequences of gyrB, rpoB, and cpn60 genes, together with the DNA–DNA hybridization results, clearly supported the existence of these two differentiated taxa (Figueras et al., 2011a). Also in 2011, A. trophiarum was discovered from the intestinal tract of healthy fattening pigs, which interestingly showed the closest similarities (≥97.4%) with the other species also recovered from humans or animals, i.e., A. cryaerophilus, A. thereius, A. cibarius, or A. skirrowii (De Smet et al., 2011; Figueras et al., 2014; Van den Abeele et al., 2014).

In 2013, the species A. cloacae and A. suis were described, using a Multilocus Sequence Analysis (MLSA) approach including five housekeeping genes (Levican et al., 2013) for the first time. Simultaneously, and due to the highest 16S rRNA gene similarity with A. marinus (95.5%), the species A. anaerophilus was incorporated to the genus (Sasi-Jyothsna et al., 2013). However, this species showed atypical characteristics, including lack of motility and obligate anaerobic metabolism, which led to the original description of the genus Arcobacter being emended (Sasi-Jyothsna et al., 2013). The most recently described species from shellfish are A. lekithochrous, A. haliotis, and A. canalis (Diéguez et al., 2017; Tanaka et al., 2017; Pérez-Cataluña et al., 2018a). The first one included several isolates recovered from scallop larvae and from tank seawater of a Norwegian hatchery (Diéguez et al., 2017), the second species came from an abalone of Japan (Tanaka et al., 2017) and the third from oysters submerged in a water channel contaminated with wastewater (Pérez-Cataluña et al., 2018a). However, Diéguez et al. (2018) evidenced that the species A. haliotis is a later heterotypic synonym of A. lekithochrous. Additionally, the low 16S rRNA gene similarity of A. lekithochrous with the known Arcobacter species (91.0–94.8%) found in the A. lekithochrous description made Diéguez et al. (2017) suggest that certain species might belong to other genera and recommend that a profound revision of the genus might clarify the taxonomy.

On the other hand, adding 2.5% NaCl to the enrichment medium and subculturing on marine agar, Salas-Massó et al. (2016) recognized seven potential new species from water and shellfish (mussels and/or oysters), and recovered new isolates of A. halophilus and A. marinus of which only the type strains had been known. In addition, during the characterization of the most recently described species A. canalis (Pérez-Cataluña et al., 2018a) and when trying to define the seven mentioned new species, we observed that the Arcobacter species formed several different clusters distant enough to suspect they might correspond to different genera, in agreement with Diéguez et al. (2017).

There are clear criteria for describing new bacterial species (Tindall et al., 2010; Figueras et al., 2011a,b). However, the description of a genus is usually based on a cut-off of <95% similarity in the 16S rRNA gene sequence, and a G+C (% mol) content differing by more than 10% (Rosselló-Mora and Amann, 2001; Yarza et al., 2008; Tindall et al., 2010; Yarza et al., 2014). Nowadays, genomic data like the Average Nucleotide Identity (ANI) and the in silico DNA–DNA hybridization (isDDH) are used to define bacterial species, although have not yet been fully explored for delineating genera (Konstantinidis and Tiedje, 2005; Goris et al., 2007; Richter and Rosselló-Móra, 2009; Qin et al., 2014; Chun et al., 2018).

A percentage of Average Amino-acid Identity (AAI) ranging from 60 to 80% between the compared genomes of species or strains and a Percentage of Conserved Proteins (POCPs) above 50% has been proposed if they are to belong to the same genus (Konstantinidis and Tiedje, 2005; Qin et al., 2014). Finally, the Relative Synonymous Codon Usage (RSCU) has also been used by some authors to infer evolutionary and ecological links among bacterial species (Ma et al., 2015; Farooqi et al., 2016).

Very recently, Waite et al. (2017) carried out a comparative genomic analysis of the class Epsilonproteobacteria. Using 16S and 23S rRNA, 120 single-copy marker proteins and AAI analysis they proposed its reclassification as the new phylum Epsilonbacteraeota. In that study, Waite et al. (2017) also proposed a reclassification of the genus Arcobacter as a new

Family Arcobacteraceae, within the class Campylobacteria, order Campylobacterales. One weakness of this study, specifically regarding the genus Arcobacter, is that only seven validated species were included in the analysis. The new family therefore comprised only the genus Arcobacter. However, these findings also support the need for a clarification of the taxonomy of the current genus Arcobacter.

The rise of genome sequencing has dramatically changed the landscape of systematics of prokaryotes, improving different aspects such as the identification of species, the functional characterization for resolving taxonomic groups, and the resolution of the phylogeny of higher taxa (Whitman, 2015). It seems clear that the incorporation of genomics into the taxonomy will boost its credibility providing reproducible, reliable, highly informative means to infer phylogenetic relationships among prokaryotes, and avoiding unreliable methods and subjective difficult-to-replicate data (Chun and Rainey, 2014; Chun et al., 2018).

Within this modern taxonomy context, the objective of the present study was to reassess the taxonomy of the known and newly recognized Arcobacter species by using a MLSA of 13 housekeeping genes, the whole genome sequences and the derived genomic analysis. The latter analysis included ANI, isDDH, AAI, POCP, and RSCU of all Arcobacter type strains. In addition, phylogenies based on 16S and 23S rRNA gene sequences were also performed with comparative purposes. The new taxonomic criteria were stable when including whole genome sequences of a second strain of each species or of unassigned sequences obtained from the public databases.

#### MATERIALS AND METHODS

#### Bacterial Strains

All 27 valid species included in the genus Arcobacter have been studied. They are represented by 39 strains, and 13 strains that are potentially new species (**Table 1**). Furthermore, 50 genomes of Arcobacter strains identified at species level were investigated, 39 of which were obtained in our laboratory (27 from known species and 13 from potentially new species) and the others from the public databases1,2 . Five genomes that had been deposited as Arcobacter sp. in the databases were also included in the study. If there was more than one strain of a known Arcobacter species, two representative genomes for each species were included in the analysis. The only exceptions were: A. acticola (Park et al., 2016) and A. pacificus (Zhang et al., 2015), whose taxonomic positions were only inferred by the phylogenetic analysis of the 16S rRNA gene sequences published in their species descriptions, together with a MLSA of three housekeeping genes (atpA, gyrB, and rpoB) for A. pacificus (Zhang et al., 2015; Park et al., 2016). The strains considered potentially new species, and named hereafter as 'candidate species,' had been recognized with an MLSA analysis of five housekeeping genes (atpA, gyrA, gyrB, hsp60, and rpoB) (data not shown).

Culturing for genome sequencing was carried out either on blood agar (DIFCO, Madrid, Spain) or marine agar (Scharlau, Sentmenat, Spain) at 30◦C in aerobiosis for 24–72 h, depending on the requirements. DNA was extracted using Easy-DNATM gDNA Purification kit (Invitrogen, Madrid, Spain) following the manufacturer's instructions. The integrity of the DNA was evaluated by electrophoresis of 10 µl of the sample in a 1.5% agarose gel. The total amount of DNA was quantified using QubitTM with the dsDNA Broad Range Assay kit (Invitrogen). Paired-end libraries were constructed with 50 ng of DNA using Nextera DNA Library Preparation Kit (Illumina, Lisbon, Portugal) and sequenced with MiSeq platform (Illumina). Sequencing generated 2 × 300 bp paired-end reads. Clean reads were assembled with SPAdes (Nurk et al., 2013) and the CGE assembler (Larsen et al., 2012) in order to select the better assembly. Before depositing the genomes in the NCBI database, FASTA files were screened for eukaryotic and prokaryotic sequences using BLASTn, and for adaptors with VecScreen standalone software3 . The five housekeeping genes used in the first MLSA analysis (atpA, gyrA, gyrB, hsp60, and rpoB) were extracted from each genome and compared with the Sanger sequences of these genes obtained originally for the identification of the strain. The existence of a single and identical copy of these genes confirmed that the genomes were not contaminated and belonged to the correct strain. Finally, contigs were deleted if they had less than 200 bp. The genomes were deposited in the GenBank database and **Table 1** lists the accession numbers.

The 55 genomes were annotated with a local installation of Prokka v1.2 (Seemann, 2014) using an e-value of 1e-06. The annotation was performed with Prokka, with the prediction tools Prodigal v2.6 (Hyatt et al., 2010) and ARAGORN v1.2 (Laslett and Canback, 2004). The prediction tool Barrnap v0.64 included in Prokka v1.2 was used for the annotation of rRNA genes. Coding sequences (CDS) were annotated, combining the Rapid Annotation Subsystems Technology (RAST) (Overbeek et al., 2014) using the classic RAST scheme and the Annotation Tools of PATRIC server (Wattam et al., 2017). The characteristics of each genome (i.e., N50, number of contigs, number of CDS, G+C content) were obtained from NCBI annotations.

#### Analysis of Housekeeping Genes, Ribosomal Genes, and Core Genome

Thirteen housekeeping genes (atpA, atpD, dnaA, dnaJ, dnaK, ftsZ, gyrA, hsp60, radA, recA, rpoB, rpoD, and tsf) were obtained from the genomes using BLASTn search. Sequence similarities of housekeeping genes were determined using the MegAlign program (DNASTAR R , Madison, WI, United States). Genes were aligned using ClustalW (Larkin et al., 2007) and phylogenies based on individual genes and on the concatenated sequences was constructed with MEGA version 6.0 (Tamura et al., 2013) using the Neighbor-Joining (NJ) and Maximum-Likelihood (ML) algorithms.

1https://www.ncbi.nlm.nih.gov/genome/

2https://gold.jgi.doe.gov/

3 ftp://ftp.ncbi.nlm.nih.gov/blast/demo/

4http://www.vicbioinformatics.com/software.barrnap.shtml

The phylogenetic analysis of the core genome was assessed with the Roary software (Page et al., 2015) using 80% as cut-off for the BLASTp search. The core genome alignment was extracted with the latter software and the phylogeny was inferred using SplitsTree version 4.14.2 as described in Sawabe et al. (2007) using SplitsTree version 4.14.2, with a neighbor net drawing and Jukes-Cantor correction (Bandelt and Dress, 1992; Huson and Bryant, 2005).

Furthermore, the 16S and 23S rRNA genes of each genome were obtained using RNammer (Lagesen et al., 2007). In some cases, 16S rRNA gene sequences were obtained in our laboratories by Sanger sequencing or from the GenBank. The similarity of the 16S rRNA genes was calculated using MegAlign version 7.0.0 (DNASTAR R , Madison, WI, United States). Phylogenetic trees were reconstructed with MEGA version 6.0 (Tamura et al., 2013) also using the NJ and ML algorithms. Alignments obtained for both genes were visually analyzed in order to localize signature sequences for strains or groups of strains.

## Genomic Indices

In order to ensure the correct assignation at species level of each analyzed genome, the ANI and the isDDH were calculated between all the genomes (Konstantinidis and Tiedje, 2005; Richter and Rosselló-Móra, 2009; Qin et al., 2014). The ANIb was calculated using JSpeciesWS (Richter et al., 2016), the resulting matrix was clustered and visualized using ggplot2 2.2.1 package (Wickham, 2009) and the isDDH was calculated with the GGDC software using results obtained with the formula 2 (Meier-Kolthoff et al., 2013). Two other indices (AAI and POCP) described for genus classification (Konstantinidis and Tiedje, 2005; Luo et al., 2014; Qin et al., 2014) were calculated among the genomes that corresponded to the type strains of the accepted species and the reference strains of the candidate species. The AAI was calculated with the Lycoming College Newman Lab AAIr Calculator5 using the Sequence-Based Comparison Tools output file from RAST (Overbeek et al., 2014). The POCP was determined as described by Qin et al. (2014) using the following parameters to consider a peptide as a conserved protein: an e-value lower than 1e-5 and an identity percentage higher than 40% from an aligned region higher than 50%.

Finally, the RSCU was computed using the Codon Adaptation Index (CAI) developed by Sharp and Li (1987) through the CAIcal web-server (Puigbò et al., 2008). Statistical differences in the RSCU were assessed by a multinomial regression approach using the R software environment (R Core Team, 2015). The principal component analysis (PCA) was performed by the R software environment (R Core Team, 2015, and visualized using ggplot2 2.2.1 and ggfortify 0.4.4 (Wickham, 2009; Horikoshi and Tang, 2015; Tang et al., 2016) or pca3d 0.10 (Weiner, 2017) packages.

#### Phenotypic Analysis and Metabolic Inference

Phenotypic characterization of each described species was obtained from this study, from the original descriptions or from the summary published by On et al. (2017). For the potentially new Arcobacter species, the phenotype was characterized following the recommended minimal standards described for new taxa of the family Campylobacteraceae (Ursing et al., 1994; On et al., 2017) and with complementary tests used in the description of other Arcobacter species (Levican et al., 2013).

Inference of the metabolic routes from the genome sequences was performed with the software package Traitar (Microbial Trait Analyzer) (Weimann et al., 2016), using the protein coding genes files obtained with Prokka v1.2 (Seemann, 2014). Traitar software is based on phenotypic data extracted from the Global Infectious Disease and Epidemiology Online Network (GIDEON) and Bergey's Systematic Bacteriology. The software uses two prediction models: the phypat classifier, which predicts the presence/absence of proteins found in the phenotype of 234 bacterial species; and the phypat+PGL classifier, which uses the same information as the phypat combined with the information of the acquisition and loss of protein families and phenotypes during evolutive events. A total of 67 traits available within the software, related to oxygen requirement, enzymatic activities, proteolysis, antibiotic resistance, morphology and motility and the use of different carbon sources, were tested and the combined results of the two predictors were analyzed using a heat map.

# RESULTS AND DISCUSSION

#### Strains and Genomes

All the 27 species currently included in the genus Arcobacter and 13 candidate species have been investigated in the present study, which has analyzed 55 genomes, 16 of them from the public databases and 39 sequenced in this study (**Tables 1**, **2**). It was not possible to analyze the genomes from A. acticola and A. pacificus because we were unable to get the type strains of the species. The contigs obtained and the N50 values complied with the recently proposed minimal standards for the use of genomes in taxonomic studies (Chun et al., 2018). The genome size ranged from 1.81 Mb for A. skirrowii F28 to 3.60 Mb for A. lekithochrous CECT 8942T (**Table 2**). The G+C content ranged from 26.1% in A. molluscorum CECT 7696T to 34.9% in 'A. aquaticus' W112-28. The G+C values agree with the range from 24.6% (which corresponded to the type strain of A. anaerophilus) to 31% indicated for the genus Arcobacter in the recent emended description by Sasi-Jyothsna et al. (2013). Interestingly, 26 genomes (47.3%) showed the presence of Clustered Regularly Interspaced Short Palindromic Repeats (CRISPRs) and CRISPRassociated genes, related with the immune response of the bacteria.

#### Taxonomic and Phylogenetic Analysis

Similarities in the 16S rRNA gene sequences among type and representative strains of the different Arcobacter species (all the 27 species currently included in the genus and the 13 new candidate species) showed a wide range of values (**Supplementary Tables S1, S4**). They ranged from 90.8% (observed between A. anaerophilus and A. faecis) to 99.9% (between A. butzleri and 'A. lacus'). The lower range of

5http://lycofs01.lycoming.edu/∼newman/AAI/

aGenome not available; bGenome sequenced in this study; cGenome obtained from NCBI database; dGenome obtained from JGI Gold atabase; ePNC means PobleNou Channel, which is a freshwater channel heavily (geometric mean of E. coli counts 4.1 × 104 c.f.u./100ml) contaminated with wastewater where shellfish were exposed for 72h (Salas-Massó et al., 2016, 2018). fThis strain was obtained from F.J. García from the Laboratorio Central de Veterinaria de Algete, MAGRAMA, Madrid, Spain; gThese strains were recovered at the Faculty of Pharmacy, University of the Basque Country (UPV-EHU), Vitoria-Gasteiz, Spain, by R. Alonso, I. Martinez-Malaxetxebarria and A. Fernández-Astorga.

similarity (90.8%) is due to the fact that those species, as occurred with others, were assigned within the genus based on the premise that 16S rRNA gene similarity was higher with any type strain of Arcobacter than with other taxa. However, in some cases being below the 95% cut-off value for genus delimitation (Rosselló-Mora and Amann, 2001; Yarza et al., 2008; Tindall et al., 2010; Figueras et al., 2011a,b). It is interesting to point out that 16S rRNA gene sequence similarities among A. nitrofigilis, the type species of the genus, and the other described species ranged from 93.2% (with A. thereius) to 95.9% (with A. venerupis). Furthermore, A. nitrofigilis showed higher similarities than the threshold value of 95% with only seven species (A. acticola, 'A. caeni,' A. cloacae, A. defluvii, A. ellisii, A. suis, and A. venerupis) out of the 27 accepted species. In any case, from the analysis of the similarities in the 16S rRNA gene sequences among the Arcobacter species it is clear that this gene has limited value and that other approaches available in the genomic era of taxonomy are needed for their study.

Phylogenetic analysis based on the core genome made up of 286 genes (**Figure 1** and **Supplementary Table S5**) and also on the concatenated sequences of 13 housekeeping genes of the representative Arcobacter strains (**Figure 2**) revealed that the Arcobacter species could be grouped into 4 major monophyletic clusters. Cluster 1, comprised seven validated species: A. butzleri, A. cibarius, A. cryaerophilus, A. lanthieri, A. skirrowii, A. thereius, and A. trophiarum, together with A. faecis (species described but not validated yet) and five candidate taxa 'A. hispanicus,' 'A. lacus,' 'A. miroungae,' 'A. porcinus,' and 'A. vitoriensis' (**Figure 1**). Cluster 2 embraced the species A. aquimarinus, A. cloacae, A. defluvii, A. ellisii, A. suis, and A. venerupis, as well as the non-validated A. acticola and the candidatus 'A. caeni.' Cluster 3 included five species, A. canalis, A. halophilus, A. marinus, A. molluscorum,

and A. mytili, together with two candidates, 'A. neptunis' and 'A. viscosus.' Finally, Cluster 4 included the species A. anaerophilus, A. bivalviorum, and A. ebronensis, as well as the candidates 'A. mediterraneus,' 'A. ponticus,' and 'A. salis.' The split decomposition network analysis of the core genome showed that the species A. lekithochrous CECT 8942T and A. nitrofigilis DSM 7299T appeared as orphan species. Furthermore, with this analysis the candidatus 'A. aquaticus' W112-28 also appeared in a separate branch near to A. nitrofigilis DSM 7299T . On the other hand, both analyses, MLSA and core genome, confirmed the existence of two sub-clusters in Cluster 1 (again A. butzleri and 'A. lacus' were located in the most distant branch within the cluster), and also two subgroups could be observed in Cluster 4, one comprising the species A. anaerophilus and A. ebronensis, and the other including the rest of species within this cluster (**Figures 1**, **2**). All the clusters and sub-clusters showed a similarity in the concatenated sequences of the 13 housekeeping genes higher than 85% (**Figure 2**).

Phylogenies based on the 16S and 23S rRNA gene sequences, undertaken with the NJ and ML approacheserealso constructed with comparative purposes. 16S rRNA based tree showed also the four major clusters although less defined (**Supplementary Figure S1A**). Species within Cluster 1, showed 16S rRNA gene sequence similarities ranging from 96.1 to 99.9%. Cluster 2 yielded similarities among species for the 16S rRNA gene between 96.7 and 99.6%, whereas within Cluster 3 ranged between 93.0 and 99.1%. Finally, Cluster 4 included species with a range of 16S rRNA sequence similarity from 94.0 to 99.5%. With the exception of Cluster 3, similarity values within the clusters (>94–95%) were within the classical boundaries for genus assignation in bacterial taxonomy (Rosselló-Mora and Amann, 2001; Yarza et al., 2008, 2014; Tindall et al., 2010; Figueras et al., 2011a,b). Our results agree with those from a recent study by Yarza et al. (2014), who investigated 568 taxa and described a threshold in 16S rRNA sequence identity of 94.5% for genus delineation.

Similar groups and topology, with only minor differences, were obtained when the 23S rRNA gene sequences were used to analyze the phylogeny of the genus (**Supplementary Figure S2**). In this analysis, the recently described species A. acticola, and A. pacificus could not be included because of the unavailability of the type strains and/or whole genome sequences. The same four major clusters formed in the 23S rRNA gene phylogenetic tree, and the species A. lekithochrous and A. nitrofigilis appeared also as orphan species (**Supplementary Figure S2**). Within Cluster 1 two subgroups could also be obtained, differentiating the species A. butzleri and 'A. lacus' from the rest of the species. Similarly, the species A. anaerophilus and A. ebronensis formed a differentiated subgroup in Cluster 4.

The visual analysis of the alignments obtained with the sequences of the 16S and 23S rRNA genes allowed the localization of signature motifs, especially in the 16S rRNA gene, for the different clusters established in the phylogenetic analysis. In these sequences, a total of 16 locations were found, presenting nucleotide combinations characteristic for the clusters (**Supplementary Figure S3**). Some of these motifs were located in helix regions as interactions with proteins of the ribosomal 30S subunit, such as helix 21 (region V4) or helix 28/44 (region V9), and therefore had a considerable level of protection against mutations (Adilakshmi et al., 2008; Kitahara et al., 2012). There are some studies on the presence of signature regions with taxonomic/phylogenetic implications in the ribosomal genes (Martínez-Murcia et al., 1992, 2007; Ue et al., 2011;Reháková ˇ et al., 2014; Martínez-Murcia and Lamy, 2015). Some regions with signature motifs detected in the present study have also shown implications for phylogenetic analysis in cyanobacteria, including regions H15, H17, H21, H22-H23, H41, and H44 (Reháková et al., 2014 ˇ ). A tree was also constructed weighting such positions (**Supplementary Figure S1B**), which allowed a better definition of the main clusters observed with the whole 16S rRNA sequences although, as expected, differentiation among species within each cluster was lower. Two sub-clusters were observed in Cluster 1, where the species A. butzleri and 'A. lacus' grouped into a well-differentiated branch with respect to the other species in the cluster (**Supplementary Figure S1B**). In this analysis, A. pacificus was clearly located in the Cluster 3, whereas in Cluster 4, A anaerophilus was the borderline species, while A. ebronensis and 'A. mediterraneus' were located in an independent branch (**Supplementary Figure S1B**). Therefore, the signature motifs described here might be a new tool for identification of the different clusters and/or genus.

#### Genomic Indices

The results of the calculations of the ANI and the isDDH among the 36 studied genomes are given in the **Supplementary Table S2** and **Supplementary Figure S4**. The results of the ANI and isDDH calculations showed that the genomes grouped into the same clusters observed by the analyses of the MLSA of the 13 housekeeping and core genes (**Figures 1**, **2**). Ranges of ANI within each cluster were from 75.2 to 95.4%, whereas isDDH values were between 19.5 and 65.4% (**Figure 2** and **Table 3**). These results confirm the phylogenetic analysis for the 13 new candidate species because all of them showed ANI and isDDH values of <96% and <70%, respectively, which are the cut-off values proposed for the delineation of new species (Konstantinidis and Tiedje, 2005; Goris et al., 2007; Richter and Rosselló-Móra, 2009; Figueras et al., 2017). As discussed in other studies, the ANI and isDDH indices provided reliable information for the delineation of Arcobacter species and are also included in the minimal guidelines to define species using genomes (Whiteduck-Léveillée et al., 2015, 2016; Figueras et al., 2017; Chun et al., 2018). Although those indices are not considered useful for delimiting genera, each of the four clusters showed values that ranged between 75.2 and 81.8% as their lowest ANI, which might be the suitable range for separating different, closely related genera. These values are relatively similar to those reported by Qin et al. (2014) that found 68–82% interspecies ANI values among the genera that they studied. Values of ANI obtained for the candidate species 'A. aquaticus' were lower than the other results, from 70.0% with A. cryaerophilus LMG 24291T to 71.9% with A. bivalviorum CECT 7835T and more in line with the Qin et al. (2014) results of 68% (**Supplementary Table S2**). In the case of the isDDH the lower values among species in the same cluster ranged between 19.5 and 24.8%, and again these might be the levels associated to different genera.

With the aim of confirming if the clusters observed might represent different genera, as suggested by the phylogenetic analyses, the similarity indices AAI and POCP were also calculated (**Supplementary Table S3**). In agreement with the 60–80% AAI that have been described for species belonging to the same genus (Konstantinidis and Tiedje, 2005) all our clusters showed lower ranges of between 67.6 to 80.3% (**Table 3**). All the clusters also complied with the POCP proposed for genus separation above 50% (Luo et al., 2014; Qin et al., 2014) because as shown in **Table 3** all clusters showed the lowest values from 67.0 to 75.4%.

It is widely known that synonymous codon usage varies among organisms and that it is related to differences in G+C content, replication strand skew, or gene expression (Suzuki et al., 2008; Farooqi et al., 2016). The interaction of these factors may vary among species depending on their evolutionary process (Ma et al., 2015). It has also been suggested that the extent of codon usage bias plays a role in the adaptation of prokaryotic organisms to their environments and lifestyles (Botzman and Margalit, 2011). To analyze the overall codon usage trends of the Arcobacter species, the frequencies of the different codons were obtained from the whole genomes and the RSCU was computed using the CAI, which is a useful tool for estimating codon usage bias (Ma et al., 2015; Farooqi et al., 2016). A first finding was that all the Arcobacter species presented a preferential use of the codons finishing in A or T (**Supplementary Figure S5**), which might be expected due to their low G+C% content. The characteristic pattern showed by A. aquaticus is noteworthy (**Supplementary Figure S5**), which supports its differentiation from the other species in Cluster 3 as well as its unique taxonomy. Such difference was the only statistically significant (p < 0.05) in the multinomial regression analysis carried out.

Next, the codon usage trends were analyzed by PCA to reveal possible evolutionary relationships. Interestingly, different groups of strains could be observed in the threedimensional graphic (**Figure 3**), which correlated with those clusters established in the different phylogenetic analyses, as shown above. As reported previously for different species of Mycoplasma (Marenda et al., 2005; Ma et al., 2015), PCA provides an additional pathway to investigate the evolutionary direction of the Arcobacter species. In addition, similarities in the synonymous codon usage patterns might reflect similar lifestyles (pathogenic vs. non-pathogenic) and adaptation to certain environments (marine water, shellfish, etc.).

# Metabolic Inference and Phenotypic Analysis

Phylogenetic and genomic analysis confirmed the existence of four clusters among the validated and candidate Arcobacter

species, which comply with the cut-off values established for the differentiation of independent genera. A thorough phenotypic analysis was therefore carried out to determine if the description of new taxa at genus level was possible or if such clusters were only clades or genomovars within the genus Arcobacter. In fact, this is what has occurred in a recent polyphasic study of 52 A. cryaerophilus strains (including genome information) in which, despite four different genomospecies

**184**

being recognized, the phenotypic characterization did not allow their differentiation into separate species and were therefore considered genomovars (Pérez-Cataluña et al., 2018a).

Phenotypic inference using Traitar confirmed the lack of reaction of Arcobacter species to most of the tests commonly used for bacterial identification (**Supplementary Figure S6**). Thus, all the type and representative strains rendered negative results, regardless of the predictor employed, for use as the sole carbon source of sugars (D-Mannitol, D-Mannose, Salicin, or Trehalose, among others) and carboxylic acids (Citrate or Malonate). Such results have been previously reported in the original descriptions of the species (see review of On et al., 2017). On the other hand, there was some incongruence between results from Traitar and those obtained by classical characterization for some tests, including growth on MacConkey agar or urea hydrolysis (data not shown). A possible explanation is related with the macroaccuracy of the predictors employed in the Traitar analysis (82.6–85.5%), as reported in the original description of the microbial trait analyzer (Weimann et al., 2016). The fact that

fmicb-09-02077 September 1, 2018 Time: 10:25 # 11

some of the Arcobacter species studied are halophilic cannot be ignored, since some of the media usually employed in the wet-lab characterization are developed for non-halophilic microorganisms.

The heat maps built from the combined results of both predictors in the Traitar analysis revealed the existence of similarity groups regarding the metabolic characteristics of the Arcobacter type strains (**Supplementary Figure S6**). In most case, clustering of strains supported the groups obtained with genomic tools, although some incongruence was also observed, such as for A. butzleri (better related here to A. defluvii, A. ellisii or A. cloacae), A. mytili (closest Traitar species 'A. caeni') or A. venerupis (forming a branch with A. ebronensis and 'A. ponticus'). In any case, Traitar might be helpful as a first-step method for phenotypic inference, although further verification should be made, especially in environmental bacterial species with special growth requirements (i.e., halophilic conditions).

A deep review of the characteristics reported in the original descriptions of the Arcobacter species, together with results obtained in our respective laboratories, allowing phenotypic traits to differentiate the clusters established by the phylogenetic and genomic analyses (**Table 4**). Growth at 37◦C in microaerophilic condition, the halophilic character, the ability to grow in presence of glycine, safranin, oxgall, or triphenyltetrazolium chloride (TTC), the presence of some enzymatic activities, such as catalase, urease or indoxyl acetate hydrolysis, and resistance to cefoperazone among others, were the main differentiating traits. Most of these characters are included in the minimal standards for describing new species in the families Campylobacteraceae and Helicobacteraceae (On et al., 2017), and they should, therefore, also be maintained for the new family Arcobacteraceae proposed by Waite et al. (2017), once this taxonomical change is validated. The phenotypic differentiation proposed in **Table 4** enabled to further describe the new genera that corresponded to the different clusters of Arcobacter species determined in the present study.

#### Stability of the Genomic-Based Clustering

In order to test the stability of the new taxonomical scheme proposed, we analyzed the whole genome sequences using second strains from each species or from unassigned sequences obtained from the public databases. That analysis is shown in **Supplementary Figure S7** and included 55 genomes. These new phylogenetic analyses of the core genome also using a Split network showed that the four clusters were maintained, but the two clusters (Clusters 3 and 4) that include species able to grow in media containing 2.5% NaCl appeared in the right place (**Supplementary Figure S7**). The genome of Arcobacter sp. LPB0137 obtained from the NCBI database grouped with the species A. lekithochrous CECT 8942T , while the genomes Arcobacter sp. LA11 and CAB grouped together in a separate branch near to Cluster 4. Interestingly, the ANI and isDDH values of 91.4% and 45.8% between strain F2176, previously identified as A. nitrofigilis (Figueras et al., 2008), and the type strains of this species along with the phylogenetic position (**Supplementary Figure S7**), revealed that this strain belonged to another potentially new species. Furthermore, strains L and AF1028, deposited at the NCBI database as Arcobacter sp. were identified as A. defluvii and A. faecis, respectively, because they clustered with the type strains of those species (**Supplementary Figure S7**). This was also confirmed by the ANI and isDDH results being above 96% and 70%, respectively.

Collado and Figueras (2011), in their review about the epidemiology and clinical significance of the genus Arcobacter, reported that these bacteria should be considered quite atypical within the class Epsilonproteobacteria because of the great diversity of hosts and habitats from which they have been isolated. In order to show if the clusters obtained have a relationship with their ecological habitat, the origin of each strain is also given in **Supplementary Figure S7**. Despite the fact that only two strains from each species were included in the analysis, each of the clusters embraced species that had been recovered from common or related origins. Cluster 1 included by strains isolated from humans and animals, from wastewater and from broiler skin (A. cibarius CECT 7203T ). The fact that some strains isolated from wastewater that was contaminated by humans or animal excreta, gives evidence of the relationship of these sources. This finding agrees with the high abundance of Arcobacter in wastewater and in water contaminated with fecal pollution (Collado et al., 2008, 2010). Among the species of Cluster 1, both by metagenomics analysis or direct plating without enrichment (Fisher et al., 2014; Levican et al., 2016), the species A. cryaerophilus was the prevalent species in wastewater, while the species A. butzleri is normally predominant in studies that investigate water and food samples of animal origin, such as different types of meats using an enrichment step (Collado et al., 2009b; Collado and Figueras, 2011; Hsu and Lee, 2015; and references therein). So far, only the species A. cryaerophilus, A. thereius, A. trophiarum, A. cibarius or A. skirrowii have been recovered from humans or animals (De Smet et al., 2011; Figueras et al., 2014; Van den Abeele et al., 2014) and all these species are as commented in the same cluster.

Cluster 2 included strains from different origins but was dominated by species that came from wastewater, shellfish or food products. In this sense, A. defluvii CECT 7697T and 'A. caeni' RW17-10 were isolated from wastewater, while the strain A. defluvii L was recovered from a microbial fuel cell. Strains of A. defluvii have also been recovered from shellfish in other studies (Levican et al., 2014; Salas-Massó et al., 2016). The strain A. suis CECT 7833T was isolated from pork meat, but other isolates have also been obtained from buffalo milk in Italy (Levican et al., 2013; Giacometti et al., 2015). The other five strains in the cluster were isolated from shellfish, wastewater and seawater (**Table 1** and **Supplementary Figure S7**). The other two clusters (Clusters 3 and 4) included strains isolated from seawater shellfish giving evidence of the marine origin of these clusters. The orphan species (A. nitrofigilis DSM7299T , A. lekithochrous CECT 8942T , and 'A. aquaticus' W112-28) also corresponded

TABLE 4 | Differential phenotypic traits among the different clusters of Arcobacter species obtained on the basis of the characteristics of the type and representative strains of the species included in each group.

+, positive result; −, negative result; V, variable result in all the species of the cluster; aWith the exception of A. skirrowii; bWith the exception of A. pacificus; cA. lekithochrous needs sea salts to grow; dWith the exception of A. ebronensis; eWith the exception of A. molluscorum; fWith the exception of A. cibarius; gWith the exception of A. anaerophilus. ND, not determined.

to strains isolated from marine environments and their phylogenetic position was close to the two marine clusters (3 and 4).

As indicated in the review by Collado and Figueras (2011), there are many uncultured or not-yet-described species of Arcobacter, which have been recognized on the basis of nearly full-length 16S rRNA gene sequences, and which probably outnumber those species that were already known at that time. Their hosts and/or habitats are very diverse and include cod larvae, cyanobacterial mats, activated sludge, tidal and marine sediments, estuarine and river water, plankton, coral, tubeworms, snails, etc. (Collado et al., 2011; and references therein). In the near future new species can be expected to emerge that will reinforce the value of the different genera proposed in this study.

#### CONCLUSION

Genomic information obtained through next-generation sequencing leads to great advances in the systematics of prokaryotes (Whitman, 2015), not only to the general understanding of prokaryotic biology but also for the resolution of the phylogeny of taxa higher than species. Single gene phylogeny, including 16S rRNA gene, has often limitations that analysis of complete genome sequences can overcome. The study aims to use this modern taxonomy approach to clarify the relationships of the diverse Arcobacter species.

The results obtained in the present study confirmed the opinion of some authors on the need for a clarification of the taxonomy of the genus Arcobacter. The phylogenetic analyses derived from the MLSA of 13 genes and of the core genome as well as the existence of signature regions in the 16S rRNA gene have shown, together with the genomic indexes ANI (75.2–81.8%), isDDH (19.5–24.8%), AAI (67.6–80.3%), and POCP (67.0–75.4%), to be useful tools for delimiting several genomic and phylogenetic groups within this genus. The intragenus ranges and cut-off values established here might also be helpful for future taxonomic studies in other bacterial groups.

Such genomic variability, together with the determination of combinations of differentiating phenotypic traits allowed the division of the current genus Arcobacter in at least six different genera for which the names Aliiarcobacter gen. nov., Pseudoarcobacter gen. nov., Haloarcobacter gen. nov., Malacobacter gen. nov., and Poseidonibacter gen. nov. are proposed. In addition, the candidate species 'A. aquaticus' also constitutes a new genus for which the name Candidate 'Arcomarinus' gen. nov. is proposed, although such proposal should be formulated in parallel to the formal description of the species.

According to Tindall et al. (2010) "the type strain of a genus is the most important reference organism to which a novel species has to be compared." In the case of the genus Arcobacter, the type species has rarely been isolated (Collado et al., 2009b; Toh et al., 2011; Levican et al., 2016; Salas-Massó et al., 2016) and in fact, all the analyses show that A. nitrofigilis is an orphan species and the only representative of the genus Arcobacter, for which an emended description is provided.

The other genera are described here while taking into account the species validated at the time of writing but with the confidence that the formal description of the candidate species would fit in such descriptions. Thus, the genus Aliiarcobacter gen. nov. is described comprising seven species Aliiarcobacter cryaerophilus comb. nov., A. butzleri comb. nov., A. skirrowii comb. nov., A. cibarius comb. nov., A. thereius comb. nov., A. trophiarum comb. nov., A. lanthieri comb. nov., and A. faecis comb. nov. On the other hand, the genus Pseudoarcobacter gen. nov. includes

the species Pseudoarcobacter defluvii comb. nov., P. ellisii comb. nov., P. venerupis comb. nov., P. cloacae comb. nov., P. suis comb. nov., P. aquimarinus comb. nov., and P. acticola comb. nov. Four species, Malacobacter halophilus comb. nov., M. mytili comb. nov., M. marinus comb. nov., M. molluscorum comb. nov., and M. pacificus comb. nov. are compiled in the new genus Malacobacter gen. nov., whereas the genus Haloarcobacter gen. nov. comprises three species Haloarcobacter bivalviorum comb. nov., H. anaerophilus comb. nov., and H. ebronensis comb. nov. Finally, the genus Poseidonibacter gen. nov. has a unique species Poseidonibacter lekithochrous comb. nov.

#### Emended Description of the Genus Arcobacter Vandamme et al., 1991 emend. Vandamme et al., 1992 and Sasi-Jyothsna et al., 2013

Arcobacter (Ar'co.bac.ter. L. n. arcus, bow; Gr. n. bacter, rod; M. L. masc. n. Arcobacter, bow-shaped rod).

Cells are Gram-negative, curved rods 0.2–0.9 µm in diameter and 1–3 µm long. Coccoid bodies are found in old cultures but are not rapidly produced under aerobic conditions. Motile with a rapid corkscrew motion. Each cell possesses a single polar flagellum. Does not swarm. Chemoorganotrophic. Utilizes organic and amino acids as carbon sources, but not carbohydrates. Respiratory metabolism with oxygen as the terminal electron acceptor; anaerobic growth with aspartate and fumarate, but not with nitrate. Nitrate usually reduced to nitrite. Requires NaCl for growth. Grows at temperatures of 10◦C– 35◦C but not at 42◦C. Catalase, oxidase, urease, and nitrogenase positive. Phosphatase, sulfatase and indole negative. Does not hydrolyze esculin, casein, DNA, gelatine, hippurate or starch. Fluorescent pigments are not produced. Unable to grow with glycine (1% wt/vol), safranin (0.05% wt/vol), oxgall (1% wt/vol), or 2,3,5-triphenyltetrazolium chloride (0.04%, wt/vol). Positive for the hydrolysis of indoxyl acetate. Poly-β-hydroxybutyrate not produced.

The base composition of the DNA is 28.1–28.4% G+C as determined from the genomes.

The type species is Arcobacter nitrofigilis.

#### Description of Aliiarcobacter gen. nov.

Aliiarcobacter (A.li.i.ar.co.bac'ter, L. pronoun alius other, another; N.L. masc. n. Arcobacter a bacterial generic name; N.L. masc. n. Aliiarcobacter the other Arcobacter).

Cells are Gram-negative, curved rods 0.2–0.5 µm in diameter and 1–3 µm long. Motile by single polar flagellum. Does not swarm. Chemoorganotrophic. Oxidase and catalase positive. No growth occur at 4% NaCl. Growth occurs at 15◦C–42◦C. Carbohydrates are not fermented. Nitrate usually reduced to nitrite. Positive for the hydrolysis of indoxyl acetate and negative for urease. Growth does not occur in the presence 2,3,5 triphenyltetrazolium chloride (0.04%, wt/vol) or glycine (1% wt/vol). Some species may grow in the presence of safranin (0.05% wt/vol) or oxgall (1% wt/vol). Fluorescent pigments are not produced. Some species are sensitive to cefoperazone (64 mg/l). Range of DNA G+C content is 26.4–29.4 mol%.

The type species is Aliiarcobacter cryaerophilus.

#### Description of Aliiarcobacter cryaerophilus comb. nov.

Basonym: Campylobacter cryaerophila Neill et al., 1985.

Other synonym: Arcobacter cryaerophilus Vandamme et al., 1991.

The description is the same given by Neill et al. (1985). The type strain is A169/BT (= NCTC 1185T = ATCC 43158T ).

#### Description of Aliiarcobacter butzleri comb. nov.

Basonym: Campylobacter butzleri Kiehlbauch et al., 1991.

Other synonym: Arcobacter butzleri Vandamme et al., 1992.

The description is the same given by Vandamme et al. (1992). The type strain is LMG 10828T (= CDC D2686T = ATCC 49616T ).

## Description of Aliiarcobacter skirrowii comb. nov.

Basonym: Arcobacter skirrowii Vandamme et al., 1992.

The description is the same given by Vandamme et al. (1992). The type strain is Skirrow 449/80T (= LMG 6621T = CCUG 10374T ).

#### Description of Aliiarcobacter cibarius comb. nov.

Basonym: Arcobacter cibarius Houf et al., 2005. The description is the same given by Houf et al. (2005). The

type strain is LMG 21996T (= CCUG 48482T ).

#### Description of Aliiarcobacter thereius comb. nov.

Basonym: Arcobacter thereius Houf et al., 2009.

The description is the same given by Houf et al. (2009). The type strain is LMG 24486T (= CCUG 56902T ).

#### Description of Aliiarcobacter trophiarum comb. nov.

Basonym: Arcobacter trophiarum De Smet et al., 2011.

The description is the same given by De Smet et al. (2011). The type strain is 64T (= LMG 25534T = CCUG 59229T ).

### Description of Aliiarcobacter lanthieri comb. nov.

Basonym: Arcobacter lanthieri Whiteduck-Léveillée et al., 2015.

The description is the same given by Whiteduck-Léveillée et al. (2015). The type strain is AF1440T (= LMG 28516T = CCUG 66485T ).

#### Description of Aliiarcobacter faecis comb. nov.

Basonym: Arcobacter faecis Whiteduck-Léveillée et al., 2016.

The description is the same given by Whiteduck-Léveillée et al. (2016). The type strain is AF1078T (= LMG 28519T = CCUG 66484T ).

# Description of Pseudoarcobacter gen.

#### nov.

Pseudoarcobacter (Pseu.do.ar.co.bac'ter, Gr. adj. pseudes, false; N.L. masc. n. Arcobacter a bacterial generic name; N.L. masc. n. Pseudoarcobacter, false Arcobacter).

Gram-negative, cells are rod shaped and motile. Cell size 0.2–0.9 µm in diameter and 0.4–2.2 µm long. Some species may present cells up to 10 µm in length. Oxidase and catalase positive. No growth occurs at 4% NaCl. Growth occurs at 15– 37◦C, but not at 42◦C. Carbohydrates are not fermented. Reduce nitrate to nitrite. Positive for the hydrolysis of indoxyl acetate. Some species may hydrolyze urea. Growth does not occur in the presence 2,3,5-triphenyltetrazolium chloride (0.04%, wt/vol) or glycine (1% wt/vol). Some species may grow in the presence of safranin (0.05% wt/vol) or oxgall (1% wt/vol). Sensitive to cefoperazone (64 mg/l). Range of DNA G+C content is 26.3– 28.0 mol%.

The type species is Pseudoarcobacter defluvii.

### Description of Pseudoarcobacter defluvii comb. nov.

Basonym: Arcobacter defluvii Collado et al., 2011.

The description is the same given by Collado et al. (2011). The type strain is SW28-11T (= CECT 7697T = LMG 25694T ).

#### Description of Pseudoarcobacter ellisii comb. nov.

Basonym: Arcobacter ellisii Figueras et al., 2011b.

The description is the same given by Figueras et al. (2011b). The type strain is F79-6T (= CECT 7837T = LMG 26155T ).

#### Description of Pseudoarcobacter venerupis comb. nov.

Basonym: Arcobacter venerupis Levican et al., 2012.

The description is the same given by Levican et al. (2012). The type strain is F67-11T (= CECT 7836T = LMG 26156T ).

# Description of Pseudoarcobacter cloacae comb. nov.

Basonym: Arcobacter cloacae Levican et al., 2013.

The description is the same given by Levican et al. (2013). The type strain is SW28-13T (= CECT 7834T = LMG 26153T )

# Description of Pseudoarcobacter suis

# comb. nov.

Basonym: Arcobacter suis Levican et al., 2013.

The description is the same given by Levican et al. (2013). The type strain is F41T (= CECT 7833T = LMG 26152T ).

#### Description of Pseudoarcobacter aquimarinus comb. nov.

Basonym: Arcobacter aquimarinus Levican et al., 2015.

The description is the same given by Levican et al. (2015). The type strain is W63T (= CECT 8442T = LMG 27923T ).

#### Description of Pseudoarcobacter acticola comb. nov.

Basonym: Arcobacter acticola Park et al., 2016.

The description is the same given by Park et al. (2016). The type strain is AR-13T (= KCTC 52212T = NBRC 112272T ).

### Description of Malacobacter gen. nov.

Malacobacter (Ma.la.co.bac'ter; Gr. n. malaco, soft, with soft boy, mollusc; Gr. n. bacter, rod; N.L. masc. n. Malacobacter, bacteria isolated from molluscs).

Gram-negative, cells are rod shaped and motile. Cell size 0.1– 0.6 µm wide and 0.5–3.6 µm long. Oxidase positive and catalase variable among species. Halophilic, no growth can be obtained without NaCl and capable to grow up to 4% NaCl. Growth occurs at 15◦C–37◦C. Does not grow at 37◦C in microaerophilic conditions nor at 42◦C in anaerobiosis. Carbohydrates are not fermented. Does not reduce nitrate to nitrite. Negative for the hydrolysis of urea. Some species may hydrolyze indoxyl acetate. Growth does not occur in the presence of oxgall (1% wt/vol) or 2,3,5-triphenyltetrazolium chloride (0.04%, wt/vol). Some species may grow in the presence of glycine (1% wt/vol) or safranin (0.05% wt/vol). Sensitive to cefoperazone (64 mg/l) variable among species. Range of DNA G+C content is 26.1– 27.3 mol%.

The type species is Malacobacter halophilus.

# Description of Malacobacter halophilus comb. nov.

Basonym: Arcobacter halophilus Donachie et al., 2005.

The description is the same given by Donachie et al. (2005). The type strain is LA31BT (= ATCC BAA-1022T = CIP 108450T ).

## Description of Malacobacter mytili comb. nov.

Basonym: Arcobacter mytili Collado et al., 2009a.

The description is the same given by Collado et al. (2009a). The type strain is F2075T (= CECT 7386T = LMG 24559T ).

#### Description of Malacobacter marinus comb. nov.

Basonym: Arcobacter marinus Kim et al., 2010.

The description is the same given by Kim et al. (2010), with the exception of variable result among strains for the hydrolysis of the indoxyl-acetate under microaerobic conditions (Salas-Massó et al., 2016). The type strain is CL-S1T (= KCCM 90072T = JCM 15502T ).

# Description of Malacobacter canalis comb. nov.

Basonym: Arcobacter canalis Pérez-Cataluña et al., 2018b.

The description is the same given by Pérez-Cataluña et al. (2018b). The type strain is F138-33T (= CECT 8984T = LMG 29148T ).

# Description of Malacobacter molluscorum comb. nov.

Basonym: Arcobacter molluscorum Figueras et al., 2011a.

The description is the same given by Figueras et al. (2011a). The type strain is F98-3T (= CECT 7696T = LMG 25693T ).

# Description of Malacobacter pacificus comb. nov.

Basonym: Arcobacter pacificus Zhang et al., 2015.

The description is the same given by Zhang et al. (2015). The type strain is SW028T (= DSM 25018T = JCM 17857T = LMG 26638T ).

#### Description of Haloarcobacter gen. nov.

Haloarcobacter (Ha.lo.ar.co.bac'ter, Gr. n. halo, salt; N.L. masc. n. Arcobacter, a bacterial generic name; N.L. masc. n. Haloarcobacter, Arcobacter salt loving).

Gram-negative, cells are rod shaped and motile. Cell size 0.1– 0.5 µm in diameter and 0.9–2.5 µm in length. Oxidase positive and catalase variable among species. Halophilic, growth can be obtained within the range of 0.5% (variable among species) and up to 4% NaCl. Growth occurs at 15–42◦C. Growth at 37◦C in microaerophilic conditions or at 42◦C in anaerobiosis variable among species. Carbohydrates are not fermented. Some species may reduce nitrate to nitrite. Negative for the hydrolysis of urea (with the exception of H. ebronensis). Some species may hydrolyze indoxyl acetate. Growth does not occur in the presence of oxgall (1% wt/vol) (with the exception of H. molluscorum) or 2,3,5-triphenyltetrazolium chloride (0.04%, wt/vol). No growth on CCDA. Some species may grow in the presence of glycine (1% wt/vol) or safranin (0.05% wt/vol). Sensitive to cefoperazone (64 mg/l). Range of DNA G+C content is 27.3–29.9 mol%.

The type species is Haloarcobacter bivalviorum.

### Description of Haloarcobacter bivalviorum comb. nov.

Basonym: Arcobacter bivalviorum Levican et al., 2012. The description is the same given by Levican et al. (2012). The type strain is F4T (= CECT 7835T = LMG 26154T ).

# Description of Haloarcobacter anaerophilus comb. nov.

Basonym: Arcobacter anaerophilus Sasi-Jyothsna et al., 2013.

The description is the same given by Sasi-Jyothsna et al. (2013). The type strain is JC84T (= KCTC 15071T = MTCC 10956T = DSM 24636T ).

# Description of Haloarcobacter ebronensis comb. nov.

Basonym: Arcobacter ebronensis Levican et al., 2015.

The description is the same given by Levican et al. (2015). The type strain is F128-2T (= CECT 8441T = LMG 27922T ).

#### Description of Poseidonibacter gen. nov.

Poseidonibacter (Po.se.i.do.ni.bac'ter, Gr. n. Poseidon, God of the sea; Gr. n. bacter, rod; N.L. masc. n. Poseidonibacter referring to the marine habitat of this bacteria).

Gram-negative, cells are rod shaped and motile. Oxidase and catalase positive. Halophilic, no growth can be obtained without seawater or the addition of combined marine salts to the medium. Growth occurs at 15◦C–25◦C, but not at 37◦C or 42◦C. Range of pH for growth is 6–8. Carbohydrates are not fermented. Reduce nitrate to nitrite. Negative for the hydrolysis of indoxyl acetate and urea. Growth occurs in the presence of safranin (0.05% wt/vol), and 2,3,5-triphenyltetrazolium chloride (0.04%, wt/vol), but not in the presence of glycine (1% wt/vol) sensitive to cefoperazone (30 µg). Possess ubiquinone MK-6 as a respiratory quinone. DNA G+C content is 28.7 mol%.

The type species is Poseidonibacter lekithochrous.

# Description of Poseidonibacter lekithochrous comb. nov.

Basonym: Arcobacter lekithochrous Diéguez et al., 2017.

The description is the same given by Diéguez et al. (2017). The type strain is LFT1.7T (= CECT 8942T = DSM 100870T ).

# AUTHOR CONTRIBUTIONS

MF and JR designed the work. AP-C, NS-M, and AD performed the phenotypic and phylogenetic experiments. AP-C and SB carried out the genome sequencing and analysis. AP-C, AL, and JR performed the bioinformatic work. JR, MF, AP-C, and AD wrote the paper.

# FUNDING

This work was supported in part by Grants JPIW2013-69095- C03-03 from the Ministerio de Economía y Competitividad (MINECO), AQUAVALENS of the Seventh Framework Program (FP7/2007-2013) grant agreement 311846 from the European Union and AGL2013-42628-R and AGL2016-77539-R (AEI/FEDER UE) from the Agencia Estatal de Investigación (Spain).

#### ACKNOWLEDGMENTS

The authors thank Dr. F. J. García (Laboratorio Central de Veterinaria de Algete, MAGRAMA, Madrid, Spain) and Drs. R. Alonso, I. Martinez-Malaxetxebarria, and A. Fernandez-Astorga [Faculty of Pharmacy, University of the Basque Country (UPV-EHU), Vitoria-Gasteiz, Spain], for kindly providing some of the

Arcobacter strains. AP-C thanks Institut d'Investigació Sanitària Pere Virgili (IISPV) for her Ph.D. fellowship and NS-M thanks the Universitat Rovira i Virgili (URV), the Institut de Recerca i Tecnologia Agroalimentària (IRTA) and the Banco Santander for her Ph.D. fellowship.

#### REFERENCES

#### SUPPLEMENTARY MATERIAL

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb. 2018.02077/full#supplementary-material

cryaerophilus as a species complex that embraces four genomovars. Front. Microbiol 9:805. doi: 10.3389/fmicb.2018.00805

Demequinaceae fam. nov. Int. J. Syst. Evol. Microbiol. 61, 1322–1329. doi: 10. 1099/ijs.0.024299-0

pig and dairy cattle manure. Int. J. Syst. Evol. Microbiol. 65, 2709–2716. doi: 10.1099/ijs.0.000318

**Conflict of Interest Statement:** The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Copyright © 2018 Pérez-Cataluña, Salas-Massó, Diéguez, Balboa, Lema, Romalde and Figueras. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# Corrigendum: Revisiting the Taxonomy of the Genus Arcobacter: Getting Order From the Caos

Alba Pérez-Cataluña1 , Nuria Salas-Massó1 , Ana L. Diéguez 2 , Sabela Balboa2 , Alberto Lema2 , Jesús L. Romalde2 \* and Maria José Figueras 1 \*

1 Departament de Ciències Mèdiques Bàsiques, Facultat de Medicina, Institut d'Investigació Sanitària Pere Virgili, Universitat Rovira i Virgili, Reus, Spain, 2 Departamento de Microbiología y Parasitología, CIBUS-Facultad de Biología, Universidade de Santiago de Compostela, Santiago de Compostela, Spain

Keywords: Arcobacter, Aliarcobacter gen. nov., Pseudoarcobacter gen. nov., Haloarcobacter gen. nov., Malacobacter gen. nov., Poseidonibacter gen. nov., taxonomic criteria

#### **A Corrigendum on**

#### Edited and reviewed by:

Yasir Muhammad, King Abdulaziz University, Saudi Arabia

#### \*Correspondence:

Jesús L. Romalde<EMAIL_ADDRESS>Maria José Figueras<EMAIL_ADDRESS>
#### Specialty section:

This article was submitted to Evolutionary and Genomic Microbiology, a section of the journal Frontiers in Microbiology

Received: 02 November 2018 Accepted: 03 December 2018 Published: 21 December 2018

#### Citation:

Pérez-Cataluña A, Salas-Massó N, Diéguez AL, Balboa S, Lema A, Romalde JL and Figueras MJ (2018) Corrigendum: Revisiting the Taxonomy of the Genus Arcobacter: Getting Order From the Caos. Front. Microbiol. 9:3123. doi: 10.3389/fmicb.2018.03123 **Revisiting the Taxonomy of the Genus Arcobacter: Getting Order From the Chaos.**

by Pérez-Cataluña, A., Salas-Massó, N., Diéguez, A. L., Balboa, S., Lema, A., Romalde, J. L., et al. (2018). Front. Microbiol. 9:2077. doi: 10.3389/fmicb.2018.02077

In the original article, there was an error. The NCTC number for the strain type Arcobacter cryaerophilus, was incorrect. A correction has been made to the **Conclusion, Description of Aliiarcobacter cryaerophilus comb. nov.:**

The description is the same given by Neill et al. (1985). The type strain is A169/BT (= NCTC 11885T = ATCC 43158T ).

There were also typographical errors of some Basonym descriptions.

A correction has been made to the Conclusion, Description of Pseudoarcobacter ellisii comb. nov., Malacobacter mytili comb. nov., Malacobacter canalis comb. nov., Malacobacter molluscorum comb. nov and Malacobacter pacificus comb. nov.:

Basonym: Arcobacter ellisii Figueras et al., 2011b

Basonym: Arcobacter mytili Collado et al., 2009

Basonym: Arcobacter canalis Pérez-Cataluña et al., 2018

Basonym: Arcobacter molluscorum Figueras et al., 2011a

and, Basonym: Arcobacter pacificus Zhang et al., 2016.

After indication of the nomenclature editors for validation of new bacterial names, the term Aliiarcobacter has been corrected to Aliarcobacter throughout the article.

In the original article, the reference for Zhang et al., 2015 was incorrectly written as Zhang, Z., Yu, C., Wang, X., Yu, S., and Zhang, X. H. (2015), It should be Zhang, Z., Yu, C., Wang, X., Yu, S., and Zhang, X. H. (2016).

The authors apologize for these errors and state that they do not change the scientific conclusions of the article in any way. The original article has been updated.

# REFERENCES

sp. nov. Int. J. Syst. Bacteriol. 35, 342–356. doi: 10.1099/00207713-3 5-3-342

Copyright © 2018 Pérez-Cataluña, Salas-Massó, Diéguez, Balboa, Lema, Romalde and Figueras. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# Corrigendum (2): Revisiting the Taxonomy of the Genus Arcobacter: Getting Order From the Chaos

Alba Pérez-Cataluña1 , Nuria Salas-Massó1 , Ana L. Diéguez 2 , Sabela Balboa2 , Alberto Lema2 , Jesús L. Romalde2 \* and María J. Figueras 1 \*

1 Departament de Ciències Mèdiques Bàsiques, Facultat de Medicina, Institut d'Investigació Sanitària Pere Virgili, Universitat Rovira i Virgili, Reus, Spain, 2 Departamento de Microbiología y Parasitología, CIBUS-Facultad de Biología, Universidade de Santiago de Compostela, Santiago de Compostela, Spain

Keywords: Arcobacter, Aliiarcobacter gen. nov., Pseudoarcobacter gen. nov., Haloarcobacter gen. nov., Malacobacter gen. nov., Poseidonibacter gen. nov., taxonomic criteria

#### Edited by:

**A Corrigendum on**

Martin G. Klotz, Washington State University, United States

#### Reviewed by:

David W. Ussery, University of Arkansas for Medical Sciences, United States Iain Sutcliffe, Northumbria University, United Kingdom George Garrity, Michigan State University, United States

#### \*Correspondence:

Jesús L. Romalde<EMAIL_ADDRESS>María J. Figueras<EMAIL_ADDRESS>
#### Specialty section:

This article was submitted to Evolutionary and Genomic Microbiology, a section of the journal Frontiers in Microbiology

Received: 01 August 2019 Accepted: 17 September 2019 Published: 01 October 2019

#### Citation:

Pérez-Cataluña A, Salas-Massó N, Diéguez AL, Balboa S, Lema A, Romalde JL and Figueras MJ (2019) Corrigendum (2): Revisiting the Taxonomy of the Genus Arcobacter: Getting Order From the Chaos. Front. Microbiol. 10:2253. doi: 10.3389/fmicb.2019.02253 **Revisiting the Taxonomy of the Genus Arcobacter: Getting Order From the Chaos**

by Pérez-Cataluña, A., Salas-Massó, N., Diéguez, A. L., Balboa, S., Lema, A., Romalde, J. L., et al. (2018) Front. Microbiol. 9:2077. doi: 10.3389/fmicb.2018.02077

The original paper, which contained a description of five new genera with 25 species (all new combinations) contained many errors that prevented the proposed names from being included in a Validation List in the International Journal of Systematic and Evolutionary Microbiology. The List Editors could correct part of the problems, so that the generic names Pseudarcobacter, Malaciobacter. Halarcobacter, and Poseidonibacter and the species assigned as "comb. nov." to these genera could be validly published [Oren, A. and Garrity, G. M. (2019). List of new names and new combinations previously effectively, but not validly, published. Validation List no. 185. Int. J. Syst. Evol. Microbiol. 69, 5–9]. The corrections made are explained in footnotes to the list. However, because of the nature of some of the changes required, the List Editors could not make the corrections for the proposed genus "Aliiarcobacter" and the eight proposed new combinations in the Validation List. The corrigendum published [**Revisiting the Taxonomy of the Genus Arcobacter: Getting Order From the Caos (sic),** by Pérez-Cataluña, A., Salas-Massó, N., Diéguez, A. L., Balboa, S., Lema, A., Romalde, J. L., et al. (2018). Front. Microbiol. 9:3123. doi: 10.3389/fmicb.2018.03123] failed to correct the remaining errors and introduced new problems. This necessitated a new Corrigendum in order to effectively publish the names Aliarcobacter and eight "comb. nov." species to be submitted subsequently for List validation in the International Journal of Systematic and Evolutionary Microbiology.

#### **Description of Aliarcobacter gen. nov.**

Aliarcobacter (A.li.ar.co.bac'ter. L. pronoun alius other, another; N.L. masc. n. Arcobacter a bacterial generic name; N.L. masc. n. Aliarcobacter the other Arcobacter). Cells are Gram-negative, curved rods 0.2–0.5 mm in diameter and 1–3 mm long. Motile by single polar flagellum. Does not swarm. Chemoorganotrophic. Oxidase and catalase positive. No growth occur at 4% NaCl. Growth occurs at 15 ◦ C−42 ◦ C. Carbohydrates are not fermented. Nitrate usually reduced to nitrite. Positive for the hydrolysis of indoxyl acetate and negative for urease. Growth does not occur in the presence 2,3,5-triphenyltetrazolium chloride (0.04%, wt/vol) or glycine (1% wt/vol). Some species may grow in the presence of safranin (0.05% wt/vol) or oxgall (1% wt/vol). Fluorescent pigments are not produced. Some species are sensitive to cefoperazone (64 mg/l). Range of DNA GC+C content is 26.4–29.4 mol%. The type species is Aliarcobacter cryaerophilus.

#### **Description of Aliarcobacter cryaerophilus comb. nov.**

Basonym: Campylobacter cryaerophila Neill et al. 1985. The description is the same given by Neill et al. (1985). The type strain is A169/BT (= NCTC 11885T = ATCC 43158T ).

#### **Description of Aliarcobacter butzleri comb. nov**.

Basonym: Campylobacter butzleri Kiehlbauch et al. 1991. The description is the same given by Vandamme et al. (1992). The type strain is LMG 10828T (= CDC D2686T = ATCC 49616T ).

#### **Description of Aliarcobacter skirrowii comb. nov.**

Basonym: Arcobacter skirrowii Vandamme et al. 1992. The description is the same given by Vandamme et al. (1992). The type strain is Skirrow 449/80T (= LMG 6621T = CCUG 10374T ).

#### **Description of Aliarcobacter cibarius comb. nov**.

Basonym: Arcobacter cibarius Houf et al. 2005. The description is the same given by Houf et al. (2005). The type strain is LMG 21996T (= CCUG 48482T ).

#### **Description of Aliarcobacter thereius comb. nov.**

Basonym: Arcobacter thereius Houf et al. 2009.

The description is the same given by Houf et al. (2009). The type strain is LMG 24486T (= CCUG 56902T ).

#### **Description of Aliarcobacter trophiarum comb. nov.**

Basonym: Arcobacter trophiarum De Smet et al. 2011.

#### REFERENCES

The description is the same given by De Smet et al. (2011). The type strain is 64T (= LMG 25534T = CCUG 59229T ).

#### **Description of Aliarcobacter lanthieri comb. nov.**

Basonym: Arcobacter lanthieri Whiteduck-Léveillée et al. 2015. The description is the same given by Whiteduck-Léveillée et al. (2015). The type strain is AF1440T (= LMG 28516T = CCUG 66485T ).

#### **Description of Aliarcobacter faecis comb. nov.**

Basonym: Arcobacter faecis Whiteduck-Léveillée et al. 2016. The description is the same given by Whiteduck-Léveillée et al. (2016). The type strain is AF1078T (= LMG 28519T = CCUG 66484T ).

The first Corrigendum stated that the original article had been updated. This action was reversed with the publication of the Erratum. The authors apologize for this error and state that this does not change the scientific conclusions of the article in any way.

#### ACKNOWLEDGMENTS

The authors are indebted to Prof. Aharon Oren (Jerusalem, Israel) for expert guidance on following the rules on bacterial nomenclature as provided by The International Code of Nomenclature of Prokaryotes (ICNP).

an aerotolerant bacterium isolated from veterinary specimens. Int. J. Syst. Bacteriol. 42, 344–356. doi: 10.1099/00207713-42-3-344

Copyright © 2019 Pérez-Cataluña, Salas-Massó, Diéguez, Balboa, Lema, Romalde and Figueras. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# Species Delimitation, Phylogenetic Relationships, and Temporal Divergence Model in the Genus Aeromonas

#### J. G. Lorén1 , Maribel Farfán1,2 \* and M. C. Fusté1,2

1 Departament de Biologia, Sanitat i Medi Ambient, Secció de Microbiologia, Facultat de Farmàcia i Ciències de l'Alimentació, Universitat de Barcelona, Barcelona, Spain, 2 Institut de Recerca de la Biodiversitat, Universitat de Barcelona, Barcelona, Spain

#### Edited by:

Jesus L. Romalde, Universidade de Santiago de Compostela, Spain

#### Reviewed by:

Tomochika Fujisawa, Institut Pasteur, France David R. Arahal, Universitat de València, Spain Antony T. Vincent, Institut de Biologie Intégrative et des Systèmes, Canada

> \*Correspondence: Maribel Farfán<EMAIL_ADDRESS>
#### Specialty section:

This article was submitted to Evolutionary and Genomic Microbiology, a section of the journal Frontiers in Microbiology

Received: 12 January 2018 Accepted: 05 April 2018 Published: 20 April 2018

#### Citation:

Lorén JG, Farfán M and Fusté MC (2018) Species Delimitation, Phylogenetic Relationships, and Temporal Divergence Model in the Genus Aeromonas. Front. Microbiol. 9:770. doi: 10.3389/fmicb.2018.00770 The definition of species boundaries constitutes an important challenge in biodiversity studies. In this work we applied the Generalized Mixed Yule Coalescent (GMYC) method, which determines a divergence threshold to delimit species in a phylogenetic tree. Based on the tree branching pattern, the analysis fixes the transition threshold between speciation and the coalescent process associated with the intra-species diversification. This approach has been widely used to delineate eukaryote species and establish their diversification process from sequence data. Nevertheless, there are few examples in which this analysis has been applied to a bacterial population. Although the GMYC method was originally designed to assume a constant (Yule) model of diversification at between-species level, it was later evaluated simulating other conditions. Our aim was therefore to determine the species delineation in Aeromonas using the GMYC method and asses which model best explains the speciation process in this bacterial genus. The application of the GMYC method allowed us to clearly delineate the Aeromonas species boundaries, even in the controversial groups, such as the A. veronii or A. media species complexes.

Keywords: species delimitation, GMYC, diversification model, Aeromonas, mdh, recA

# INTRODUCTION

As in other organisms, diversification in bacteria leads to entities that we call species, which group together organisms that have evolved separately from others. Species can be differentiated from populations of individual strains by a maximum likelihood method that determines the point of transition of the evolutionary processes from the level of species (speciation and extinction) to the population (coalescence). The entities determined in this way maintain the biological properties and the levels of sequence divergence of the traditionally defined species.

Several methods have been proposed for species delimitation in bacteria. The classical approach uses a sequence divergence of 3 per cent pairwise distance in 16S rRNA (Schloss and Handelsman, 2006) or the 1 percent recommended by Acinas et al. (2004) as a threshold for separating species. However, rates of substitution may vary among lineages, as can the levels of variation within and between species, which makes it difficult to assume a universal threshold. Other widely applied methods for species delimitation are based on multilocus sequences, or more recently whole

genome sequencing of several individuals has been used to determine the average nucleotide identity (ANI) (Richter and Roselló-Mora, 2009) and in silico DNA-DNA hybridization (isDDH) (Meier-Kolthoff et al., 2013) values to separate species. When comparing these methods with traditional phenotypic approaches, discrepancies can sometimes arise.

The study of evolutionary patterns in prokaryotes is hampered by the lack of a reliable fossil record, limited morphological differentiation and frequently complex taxonomic relationships. The few studies in this field suggest a constant rate of diversification in free-living or symbiotic bacteria (Martin et al., 2004; Vinuesa et al., 2005; Sanglas et al., 2017), while in some pathogens, such as Borrelia burgdorferi, it seems to follow a radiation pattern (Morlon et al., 2012). Differences in the cladogenesis of B. burgdorferi could be explained by its association with vertebrate and arthropod hosts, which could restrict the gene flow between populations. Understanding the evolution of prokaryote biological diversity therefore remains a significant challenge for biologists (Barraclough et al., 2009).

The main objective of our research in the last years has been the study of diversification in the genus Aeromonas, a γ-Proteobacteria that comprises a group of Gram-negative, rod-shaped bacteria, which are found in aquatic environments worldwide and are members of the microbiota (as well as primary or secondary pathogens) of fish, amphibians and other animals (Martin-Carnahan and Joseph, 2005; Janda and Abbott, 2010). Aeromonas is an ideal genus to study the diversification processes in bacteria because its species, a combination of free-living bacteria and host-associated strains, can be isolated from a wide variety of habitats. We began by investigating the rate and pattern of cladogenesis in this bacterial genus from a collection of strains including the type strains of all the Aeromonas species, using the sequences of five housekeeping genes (Lorén et al., 2014). However, the frequently high intra-specific diversity in bacteria is not always adequately reflected by the type strain.

We therefore performed a second analysis using molecular data from the sequences of two housekeeping genes (mdh and recA) obtained from 150 strains belonging to 27 species of Aeromonas (Sanglas et al., 2017). Phylogenies that mix variation between and within species often need to be reduced to trees with only one sequence per species to avoid invalidating the results (Fontaneto et al., 2012). To fulfill these conditions we constructed a tree with the consensus sequence for each species and used the BEAST program to obtain the species tree (Sanglas et al., 2017). In both cases the results of the analysis allowed us to determine that the process of speciation in Aeromonas follows a constant model of diversification.

The genus Aeromonas, includes 30 taxonomic species, many of them described recently. However, some descriptions were based on only one strain, or using techniques such as the 16S rRNA gene sequences, which have been demonstrated to have a low discriminatory power for delimiting species in Aeromonas (Figueras et al., 2000). Consequently, the improper characterization of some Aeromonas species has led to several reclassifications, as in the case of A. aquariorum (Beaz-Hidalgo et al., 2013), A. ichthiosmia (Huys et al., 2001), A. culicicola (Huys et al., 2005), among others. In addition, the genus includes several species complexes, such as A. hydrophila (Martin-Carnahan and Joseph, 2005), A. veronii (Silver et al., 2011) or A. media (Talagrand-Reboul et al., 2017), which are constituted by groups of strains with high intra-specific heterogeneity, giving rise to controversies about the species delineation.

In this study, we aimed to determine the species delineation and corroborate the diversification process in Aeromonas by applying a sequence-based method that allows the use of several strains per species, the Generalized Mixed Yule Coalescent (GMYC) method. This approach has been widely used to delineate eukaryote species and establish the diversification process from sequence data (Pons et al., 2006; Fontaneto et al., 2007; Lahaye et al., 2008). Nevertheless, there are only a few previous examples in which this analysis has been applied to a bacterial population (Barraclough et al., 2009; Powell, 2012). The method determines a divergence threshold to delimit species in a phylogenetic tree. Based on the tree branching pattern, the analysis fixes the transition threshold between speciation and the coalescent process associated with the intra-species diversification. The pruned tree provided by the threshold allows the use of classical methods of diversification analysis to determine the speciation process followed by the species in the analyzed population.

It is presupposed that branch lengths between species are derived from speciation and extinction rates, whereas branch lengths within a species reflect a coalescence process at the population level (Pons et al., 2006). The method determines the locations of ancestral nodes that define putative species and applies a likelihood ratio test to assess the fit of the branch lengths to a mixed lineage birth-population coalescent model (Pons et al., 2006). As the analysis to determine the transition threshold is based on the tree branching pattern, the tree has to accomplish certain requirements, for example, be ultrametric, fully resolved (without multifurcations) and reliable (with a high support in the nodes).

Our aim was therefore to determine the species delineation in Aeromonas based on the GMYC method. Although initially this method assumed diversification at between-species level follows a constant (Yule) model (Pons et al., 2006; Fontaneto et al., 2007), it was later evaluated simulating other diversification models (Fujisawa and Barraclough, 2013). Despite previous results indicate that Aeromonas follows a constant (Yule) model of diversification (Lorén et al., 2014; Sanglas et al., 2017), we also determined which model best explains the speciation process in this bacterial genus.

# MATERIALS AND METHODS

#### Gene Sequences

A collection of 147 Aeromonas strains, representative of the 30 species recognized up to August 2017, was selected for the study. Two housekeeping genes (mdh and recA) were chosen for the analysis; for each strain, the full-length sequences for both genes were previously obtained (Sanglas et al., 2017). In the case of the more recently accepted Aeromonas species: A. lacus, A. aquatica and A. finlandiensis, sequences were obtained from genomes

available in the GenBank1 . The strains and sequences are listed in Supplementary Table S1, including the GenBank accession numbers.

#### Data Sets

Phylogenetic reconstruction was carried out from the concatenated sequences of mdh and recA genes. For each gene, the translated sequences were aligned using the ClustalW program implemented in MEGA6 (Tamura et al., 2013) and translated back to obtain the nucleotide alignments. Both alignments were concatenated with the DAMBE program (v5.3.10; Xia, 2013) and later the sequences were checked for the presence of incongruences or gaps. Divergent and ambiguously aligned blocks were also removed using the Gblocks program (Castresana, 2000).

#### Sequence Analysis

To ascertain if the sequences used allowed a good species discrimination we determined the intra- and inter-specific distances among the different Aeromonas species. Data were graphically depicted with box plots obtained using the R package ggplot2 (Wickham, 2016). Sequences were also analyzed with different distance models that consider the influence of saturation, base frequency, transition–transversions or multiple substitutions. All the analyses were conducted with the function dist.dna implemented in the R package ape (Paradis et al., 2004). To determine the polymorphic sites along the sequences we construct a dots plot graph with the R package phyclust (Chen, 2010, 2011). The R package NbClust (Charrad et al., 2014) was used to determine the number of clusters from the p distances of the alienated sequences (Supplementary Table S2).

#### Phylogenetic Reconstruction

Bayesian phylogenetic trees were reconstructed with the BEAST program (v1.8.1; Drummond and Rambaut, 2007; Drummond et al., 2012) from the data sets. The model of evolution for the each gene was determined using the jModelTest 2 program (Darriba et al., 2012). The general time-reversible model with discrete gamma distribution and invariant sites (GTR + G + I) was selected as the best-fit model of nucleotide substitution. The Bayesian analysis was performed using a GTR model with four gamma categories, a Yule process of speciation, and a constant clock model of evolution as the tree priors, as well as other default parameters. We used the divergence time between Escherichia coli and Salmonella enterica estimated by Ochman and Wilson as the calibration point (Ochman and Wilson, 1987a,b). Accordingly, we calibrated the divergence of Aeromonas with a normally distributed prior with a mean of 140 Ma and a standard deviation of 10 Ma. We performed three independent Markov Chain Monte Carlo (MCMC) runs 20 million generations, sampling every 2,000 generations. Posterior distributions for parameter estimates and likelihood scores to approximate convergence were visualized with the Tracer program (v1.6.0; Rambaut et al., 2014). Visual inspection of traces within and across runs, as well as the effective sample sizes (ESS) of each parameter (>200), allowed us to confirm that the analysis was adequately sampled. A maximum clade credibility (MCC) tree was chosen by TreeAnnotator (v1.8.1; Drummond et al., 2012) from the combined output of the three MCMC runs using the LogCombiner program2 after the removal of the initial trees (20–25%) as burn-in. The MCC tree was visualized with the program FigTree (v1.4.2)3 .

#### GMYC Analysis

The GMYC analysis uses the information contained in a tree to delimit species and determine the diversification model. To analyze the influence of the priors chosen, we constructed 3 different trees from the same DNA sequence alignment, changing the model used to express the branching pattern of the tree (the Yule or a coalescent model) and the rate of molecular evolution, considering a constant or a relaxed clock, as recommended by Michonneau (unpublished). The topology of the trees would be almost identical in all cases; however, the branch lengths would vary.

The analysis was conducted using the function gmyc available in the splits package implemented in R (Ezard et al., 2014), with the single threshold option as recommended by Fujisawa and Barraclough (2013). The method compares the null (no threshold) and alternative hypothesis (one threshold) and infers the number of genetic entities. Branching events between species are modeled with a Yule model (Barraclough and Nee, 2001), while branching events within species are adjusted to a neutral coalescent process (Hudson, 1991).

FIGURE 1 | Intra- and inter-specific distances in Aeromonas species. Box plots showing the intra- (turquoise) and inter-specific (salmon) p distances obtained from the data set. The ends of the boxes correspond to the first and third quartiles. The ends of the horizontal lines indicate the highest and lowest values. The black vertical line dividing the boxes represents the median of the data.

1http://www.ncbi.nlm.nih.gov/genbank/

2http://beast.bio.ed.ac.uk/logcombiner 3http://beast.bio.ed.ac.uk/figtree

Results are evaluated by a likelihood ratio test (LRT) between the null hypothesis, which considers that all strains analyzed constitute a single species, and the alternative hypothesis (GMYC model) that assumes the existence of a species-delimiting threshold. The LRT significance was calculated using a chi-square test with 2 degrees of freedom.

The gmyc function gives a likelihood score for the model that considers all sequences belonging to the same species, and a likelihood score considering the sequences split in different species. The output also lists how many clusters and entities are associated with the highest likelihood score with the corresponding confidence intervals (CI) and the estimated threshold time when there was a transition between the speciation- and coalescent-level events. The R package also contains functions that plot (1) the number of lineages-throughtime (LTT), with the inferred position of the threshold (red vertical line); (2) the likelihood profile through time; (3) the tree with the clusters highlighted in red. Additionally, the "support" for the delineated species can be plotted, indicating whether the results are reliable or not.

The use of only one indirect calibration point in the construction of the phylogeny, due to the absence of more reliable calibration data, can be a source of uncertainty. As the GMYC method relies on relative rather than absolute branch lengths, we repeated the analysis with the same alignment but removing the outgroup, obtaining the corresponding tree. Following the suggestion of Fujisawa and Barraclough (2013), we scaled branch lengths in the tree to have a root age of 1.0 before running the analysis.

#### Diversification Analyses

A standard lineage-through-time (LTT) plot was constructed using the R package ape (Paradis et al., 2004) to graphically visualize and evaluate the temporal pattern of lineage diversification in Aeromonas.

We used the birth–death likelihood (BDL) tests implemented in laser (Rabosky, 2009) to detect the diversification model and the speciation and extinction rates (λ and µ) from the pruned tree obtained by cutting the original tree at the threshold given by the GMYC analysis. To test the null hypothesis of no-rate change versus variable-rate change in diversification, we have applied the ML approach of Rabosky, the test 1AICRC (Rabosky, 2009). This statistic is calculated as: 1AICRC = AICRC – AICRV, where AICRC is the Akaike information criterion (AIC) score for the best fitting

rate-constant diversification model, and AICRV is the AIC for the best fitting variable-rate diversification model. Thus, a positive value for AICRV indicates that the data are best approximated with a rate-variable model, while a negative AICRV value suggests a rate-constant model of diversification. We tested eight different models, two of which were rate-constant (pure-birth or Yule and birth-death) and six were rate-variable (DDL, DDX, Yule 2-,3-,4 and 5- rates) (Lorén et al., 2014).

We calculated the gamma (γ) statistic (Pybus and Harvey, 2000; Fordyce, 2010) and its significance by simulating 5,000 phylogenies, as described previously (Lorén et al., 2014). This statistic compares the relative node positions in a phylogeny with those expected under a constant diversification rate model, in which the statistic follows a standard normal distribution. Positive γ values evidence that nodes are closer to the tips than expected under the constant rate model. When γ is negative, the internal nodes are closer to the root than expected under a constant model, indicating a decrease in diversification through time. In addition, we compared the observed empirical gamma value with a distribution of the gamma statistics obtained by simulation.

## RESULTS

#### Sequence Analysis

The analysis involved 147 Aeromonas strains. The number of total positions analyzed was 1,979 bp. All positions containing gaps and missing data were eliminated in the construction of the phylogenetic tree.

The intra- and inter-specific distances determined from the sequences (Supplementary Figure S1 and Supplementary Table S3), with the exception of the conflicting A. veronii group and the A. bestiarum/A. piscicola cluster, allowed a good discrimination of the species with no overlap among the distances (**Figure 1**).

We also determined the sequence evolution of the concatenated sequences of mdh and recA. As can be seen in **Figure 2**, we analyzed our sequences for the influence of saturation (**Figure 2A**), base frequencies (**Figure 2B**), transitions and transversions (**Figure 2C**), and the gamma correction for multiple substitutions (**Figure 2D**). The sequences showed no influence of saturation, base frequency bias or the transitiontransversion ratio; only the heterogeneity in the inter-site substitution rate could have some effect (**Figure 2D**).

The nucleotide substitution analysis (**Figure 3**) showed a clear diversity among the sequences of the strains belonging to the same species. Nevertheless, the different species exhibited segregating site patterns that allowed them to be separated.

#### Phylogenetic Reconstruction

The best-fit models of sequence evolution were implemented according to the Akaike Information Criterion (AIC) scores for the substitution models evaluated, using jModeltest. The general time reversible (GTR) model was selected as the best model of evolution for the concatenated sequences using a discrete gamma distribution and a fraction of invariable sites (GTR + G + I). **Figure 4** shows the Aeromonas Bayesian phylogeny, all the strains belonging to the same species clustering together in the same group. The posterior values obtained for each node were close to 1 for the majority of the main clades. The figure also shows the clusters obtained with the p distance analysis (Supplementary Table S2), which fully matched those established from the phylogeny.

**Figure 5** depicts the phylogeny obtained with the corresponding lineage-through-time plot (LTT plot). The LTT plot (**Figure 5B**) clearly shows two different linear relationships with a sudden change in the slopes at approximately −25 Ma. This breakpoint roughly matches the region of the tree between −25 Ma and the present (**Figure 5A**), which accumulates the majority of the nodes (75%). We applied a linear regression

bottom indicates the divergence time in millions of years (Ma, Mega annum).

model to the LTT plot points between the crown age (−235.6 Ma) and the breakpoint (−25 Ma) (segment I), and those obtained from the breakpoint to the present (segment II). The fit of the two segments to a linear model was very good, with R squared values of 0.9898 and 0.9849 for segments I and II, respectively. The two regression lines intersected at −26.5 Ma, and the slopes for segment I and II were 0.0178 and 0.0628 (**Figure 5B**). Using the Chow test (Chow, 1960), we verified that the two linear regression slopes differed significantly (F test = 1715.8, P value = 1.12e−100).

#### Generalized Mixed Yule Coalescent

The existence of distinct Aeromonas lineages was confirmed by the branch length analysis. The resulting LTT plot showed a steep upturn in branching rates toward the present, marking the transition from the between- to within-species rate of lineage branching (**Figure 6**). The GMYC model fitted a transition in branching rate occurring at −26.6 Ma. The support in the majority of the GMYC clusters was high, confirming the reliability of the results (**Figure 6**). The number of ML entities determined was 31 with a confident interval (CI) of 24–36 (**Table 1**). The λ values for the diversification and coalescent processes were 0.0092 and 0.2779, respectively. The scaling parameter p for the diversification process was close to 1 (1.154), which indicates a constant speciation rate model with no extinction (Yule model), while the value for the coalescent process was clearly below 1 (6.186e−08) indicating a deficit of recent coalescent events (Fujisawa and Barraclough, 2013). No substantial differences were detected between Bayesian models using a coalescent or Yule prior when constructing the tree (**Table 1**).

The calibration point used to date the phylogeny also had no influence. Identical results were obtained when the dated tree was compared with the one obtained by repeating the analysis with the same alignment but removing the outgroup.

#### Diversification Analysis

As suggested by Rabosky (2006), we calculated the significance of 1AICRC for the set of analyzed models by using the Yule model

to simulate 5,000 phylogenies of the same size and diversification rate as those obtained from our data, and determined the P value from the resulting distributions. As can be seen in **Table 2**, we cannot reject the null hypothesis of a Yule model to a level of significance of α = 0.05, which means that the diversification in Aeromonas is constant.

For further corroboration, we determined the gamma statistic of Pybus and Harvey, a powerful tool principally used for comparing models of decreasing speciation rate through time and a constant rate of diversification. The estimated γ value from the chronogram (pruned tree) was 0.798 with a P value of 0.425, indicating that we cannot reject the null hypothesis of constant diversification for our phylogeny. Data obtained from the simulation also corroborate a constant diversification process in Aeromonas (**Figure 7**).

# DISCUSSION

The choice of appropriate genes is crucial for the reconstruction of reliable molecular phylogenies. They should be housekeeping genes, evolving at a constant rate (good molecular clocks), and have unsaturated synonymous and non-synonymous substitutions. Unlike the majority of phylogenetic studies, in our analysis we checked our sequences to ensure that the genes used in this work fulfilled the above conditions and allowed a good species separation.

The phylogeny constructed from the concatenated sequences corroborates the monophyletic origin of this group of bacteria. In the chronogram obtained, the majority of the nodes were strongly supported, with posterior values close to 1. In addition, the main clade distribution was in

#### TABLE 1 | GMYC analysis.

fmicb-09-00770 April 18, 2018 Time: 17:26 # 9

aML, maximum likelihood; CI, confidence interval. bdiv, diversification process. ccoal, coalescent process.

agreement with previously published phylogenies (Roger et al., 2012; Colston et al., 2014; Lorén et al., 2014; Sanglas et al., 2017). We obtained a perfect clustering of the strains belonging to the same species, including those considered synonymous.

The analyses based on the LTT plot detected an increase in the diversification rates from a certain point in time. As **Figure 5** illustrates, a high percentage of nodes accumulate close to the present due to the greater similarity of the sequences belonging to the same species.

The LTT plot constructed from the phylogenetic tree clearly shows two different linear relationships with a sudden change in the slopes at −26.5 Ma. The fit of the two segments to a linear model was very good, with high R squared values. The Chow test verified that the two linear regression slopes differed significantly (P value = 1.12e−100), which indicates a change in the rate of lineage branching.

In this work we carried out a quantitative analysis of sequence data to delimit putative Aeromonas species by detecting shifts in the rate of lineage branching. The branch-length analysis is based on a probabilistic model that distinguishes between species diversification and coalescent processes. This approach compensates for undefined species limits by including confidence intervals when allocating species-defining nodes and does not rely on population limits; thus, species represented by a single strain can be included (Pons et al., 2006).

Some studies recommend the use of BEAST for the tree construction as input for the GMYC analysis (Monaghan et al., 2009; Tang et al., 2014). From a set of given data (sequences)

and a substitution model, the Bayesian Inference performs a probability analysis in search of the best set of trees that maximize the posterior probability. From an initial tree, it uses the Markov Chains Monte Carlo to evaluate the posterior probability of the different states proposed. After generating a large number of trees (frequently around 10,000), it uses the subsequent probability to ascertain how many times a node is repeated in each tree and this is done for each node. Subsequently, the Bayes theorem is used to build a single phylogenetic tree (maximum clade credibility, MCC). The same method provides an estimate of the reliability of the results obtained through the posterior probability values of each node, without the need to evaluate it a posteriori.

Although Monaghan et al. (2009) recommended the use of the coalescent tree prior instead of the Yule prior for constructing the tree for the GMYC analysis, in our work, in agreement with Talavera et al. (2013), we obtained identical results with the different BEAST options.

Barraclough et al. (2009) used the GMYC method to delimit bacterial species based on 16S rRNA sequences, concluding that

a AIC, Akaike Information Criterium. b1AICRC test (Rabosky, 2006). cP value obtained by simulation of 5,000 phylogenies.

the 16S rRNA is a rather conservative molecule for surveying bacterial diversity and cannot be used to distinguish putative species. They suggest that studies using sequences of multiple genes or whole genomes of individuals would be more useful. Accordingly, we applied the GMYC method to delimit species in the genus Aeromonas, and the results corroborate the species previously described with other methods (phenotypic analysis, MLSA, genome sequence analysis) for this bacterial genus. The resulting tree exhibited the branching rate pattern of the speciesto-population transition. Only the controversial A. veronii and A. media species complexes, the A. piscicola/A. bestiarum group, and two of the most recently proposed species were not entirely resolved.

The A. veronii group includes two biovars, A. veronii Veronii and A. veronii Sobria, as well as two other Aeromonas species considered as synonymous, A. culicicola (Huys et al., 2005) and A. ichthiosmia (Huys et al., 2001). When we analyzed the A. veronii species complex, four groups of strains were revealed: a main cluster consisting of 9 strains, including the two biovars, A. veronii Veronii and A. veronii Sobria; a second group with two strains, including A. ichthiosmia; a third comprising three strains, including A. culicicola; and the fourth with two strains (**Figure 8**).

Considering the GMYC results, it would be interesting to review the inclusion of A. culicicola and A. ichthiosmia in the already controversial A. veronii group (Huys et al., 2001; Esteve et al., 2003). In addition, Colston et al. (2014) determined the ANI and isDDH values from the genomes of different Aeromonas species. When the authors compared the type strain of A. culicicola and A. ichthiosmia with all the strains defined as A. veronii, the isDDH values obtained were below the threshold of 70% (69.1–69.6% for A. culicicola and 67.4–68.2% for A. ichthiosmia) and ANI values close to 96%, suggesting these species are different from A. veronii.

Aeromonas media described by Allen et al. (1983) traditionally comprised two hybridization groups: A. media HG5A and HG5B (Hänninen and Siitonen, 1995), with the strains LMG 13459 (CDC 0862-83) and CECT 4232 (ATCC 33907) representing both groups (Popoff et al., 1981). This was later corroborated by Küpfer et al. (2006) and more recently by Talagrand-Reboul et al. (2017). Our results reveal the existence of two groups of strains outside the GMYC threshold that determine the

species in Aeromonas (**Figure 8**), which is in accordance with Talagrand-Reboul et al. (2017), who analyzed the relationships of 40 A. media strains with the recently proposed species of A. rivipollensis. These authors also identified two groups of strains inside the A. media cluster with ANIb (ANI calculator and EzGenome) and isDDH values between groups of ≤94.6% and ≤61.2%, respectively. These values, which are below the considered threshold for isolates of the same species, support the existence of two different species within this group.

Another controversial GMYC clade grouped together A. bestiarum/A. piscicola (**Figure 8**). Colston et al. (2014) determined 0.04 substitutions per site for this pair of strains, which was clearly lower than those determined for other Aeromonas species, such as A. eucrenophila/A. tecta (1.0) or A. schubertii/A. diversa (0.9). Furthermore, the ANI or isDDH values for A. bestiarum/A. piscicola were close to the threshold to be considered as members of the same species.

Recently, three new Aeromonas species have been formally accepted: A. finlandensis, A. lacus and A. aquatica (Beaz-Hidalgo et al., 2015). Nevertheless, our analysis showed that A. lacus and A. aquatica should probably not be considered new species (**Figure 6**). In the case of A. lacus, the results obtained are in agreement with the ANI values calculated from the genomes in the original description of these species. When A. jandaei was compared with A. lacus, the ANI value obtained was 95.38%, similar to the species cutoff of 96%.

The GMYC analysis also reveals that diversification in Aeromonas is a constant process and follows a Yule model. This result is independent of the number of sequences per species used, and agrees with a study of macroevolutionary dynamic models in a huge number of eukaryote and prokaryote taxa (Maruvka et al., 2013). Prokaryote populations are less affected by extinction and founder effects than larger less abundant organisms (Lynch and Conery, 2003; Butterfield, 2011). Bacterial species may be described as metapopulations extending over time and evolving separately from other species (Cohan, 2002a,b; Achtman and Wagner, 2008). Although microbial dispersion can be limited by geographical barriers and the resulting physical isolation may influence microbial evolution, due to their small size bacteria generally have an unrestricted capacity for dispersion (Finlay and Clarke, 1999) via a variety of passive mechanisms and over long distances (Fierer, 2008).

## REFERENCES

# CONCLUSION

The GMYC method clearly delineated the species boundaries in Aeromonas, even in the controversial groups, such as the A. veronii or A. media species complexes. Some of these taxonomic uncertainties are due to the use of inappropriate methods for defining species [biochemical tests with miniaturized methods (API) or 16S rRNA sequences]. It would therefore be interesting to revise bacterial species separation using the emerging new analyses based on genome sequences (ANI, isDDH). The method presented here, widely used to delimit species in eukaryotes, offers an alternative that is easy to use and fully practicable for all laboratories. In addition, the GMYC method has promising potential for phylogenetic community ecology studies.

# AUTHOR CONTRIBUTIONS

JGL and MCF conceived and designed the research. JGL and MF performed the computations. JGL, MCF, and MF analyzed the data and discussed the results. All authors contributed to and revised the final manuscript version. They all approved the final version to be published.

# FUNDING

This work was supported by a grant from the University of Barcelona for Open Access Publishing in scientific journals.

# ACKNOWLEDGMENTS

We thank Dr. Ariadna Sanglas for her invaluable help and contribution of sequences to this study.

# SUPPLEMENTARY MATERIAL

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb. 2018.00770/full#supplementary-material

**Conflict of Interest Statement:** The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Copyright © 2018 Lorén, Farfán and Fusté. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

# Genome Data Provides High Support for Generic Boundaries in Burkholderia Sensu Lato

Chrizelle W. Beukes1 , Marike Palmer1 , Puseletso Manyaka1 , Wai Y. Chan1 , Juanita R. Avontuur1 , Elritha van Zyl1 , Marcel Huntemann2 , Alicia Clum2 , Manoj Pillay2 , Krishnaveni Palaniappan2 , Neha Varghese2 , Natalia Mikhailova2 , Dimitrios Stamatis2 , T. B. K. Reddy2 , Chris Daum2 , Nicole Shapiro2 , Victor Markowitz2 , Natalia Ivanova2 , Nikos Kyrpides2 , Tanja Woyke2 , Jochen Blom3 , William B.Whitman4 , Stephanus N. Venter1 \* and Emma T. Steenkamp1

1 Department of Microbiology and Plant Pathology, Forestry and Agricultural Biotechnology Institute, University of Pretoria, Pretoria, South Africa, 2 DOE Joint Genome Institute, Walnut Creek, CA, United States, 3 Bioinformatics and Systems Biology, Justus-Liebig-University Giessen, Giessen, Germany, 4 Department of Microbiology, University of Georgia, Athens, GA, United States

#### Edited by:

Sabela Balboa Méndez, Universidade de Santiago de Compostela, Spain

#### Reviewed by:

Paulina Estrada De Los Santos, Instituto Politécnico Nacional, Mexico Radhey S. Gupta, McMaster University, Canada

> \*Correspondence: Stephanus N. Venter<EMAIL_ADDRESS>
#### Specialty section:

This article was submitted to Evolutionary and Genomic Microbiology, a section of the journal Frontiers in Microbiology

Received: 20 April 2017 Accepted: 07 June 2017 Published: 26 June 2017

#### Citation:

Beukes CW, Palmer M, Manyaka P, Chan WY, Avontuur JR, van Zyl E, Huntemann M, Clum A, Pillay M, Palaniappan K, Varghese N, Mikhailova N, Stamatis D, Reddy TBK, Daum C, Shapiro N, Markowitz V, Ivanova N, Kyrpides N, Woyke T, Blom J, Whitman WB, Venter SN and Steenkamp ET (2017) Genome Data Provides High Support for Generic Boundaries in Burkholderia Sensu Lato. Front. Microbiol. 8:1154. doi: 10.3389/fmicb.2017.01154 Although the taxonomy of Burkholderia has been extensively scrutinized, significant uncertainty remains regarding the generic boundaries and composition of this large and heterogeneous taxon. Here we used the amino acid and nucleotide sequences of 106 conserved proteins from 92 species to infer robust maximum likelihood phylogenies with which to investigate the generic structure of Burkholderia sensu lato. These data unambiguously supported five distinct lineages, of which four correspond to Burkholderia sensu stricto and the newly introduced genera Paraburkholderia, Caballeronia, and Robbsia. The fifth lineage was represented by P. rhizoxinica. Based on these findings, we propose 13 new combinations for those species previously described as members of Burkholderia but that form part of Caballeronia. These findings also suggest revision of the taxonomic status of P. rhizoxinica as it is does not form part of any of the genera currently recognized in Burkholderia sensu lato. From a phylogenetic point of view, Burkholderia sensu stricto has a sister relationship with the Caballeronia+Paraburkholderia clade. Also, the lineages represented by P. rhizoxinica and R. andropogonis, respectively, emerged prior to the radiation of the Burkholderia sensu stricto+Caballeronia+Paraburkholderia clade. Our findings therefore constitute a solid framework, not only for supporting current and future taxonomic decisions, but also for studying the evolution of this assemblage of medically, industrially and agriculturally important species.

Keywords: Burkholderia, Paraburkholderia, Caballeronia, phylogenomics, Robbsia andropogonis, Burkholderia rhizoxinica

# INTRODUCTION

The genus Burkholderia was originally introduced to accommodate an assemblage of seven Pseudomonas species (Yabuuchi et al., 1992), two of which were later transferred to Ralstonia (Gillis et al., 1995; Yabuuchi et al., 1995). Since then, the number of Burkholderia species has grown substantially, to about 108 in 2015 (Estrada-de los Santos et al., 2016), spanning a range of

human, animal and plant pathogens, as well as numerous strains with significant biotechnological potential (Depoorter et al., 2016; Estrada-de los Santos et al., 2016). The latter includes the so-called plant beneficial and environmental (PBE) species (Suárez-Moreno et al., 2012), many of which are plant-associated (e.g., those with plant growth promoting activities, the symbiotic diazotrophs and free-living species with diazotrophic, bioremedial and antibiotic activities) (Depoorter et al., 2016; Estrada-de los Santos et al., 2016). Because of this heterogeneity, new genera [e.g., 'Caballeronia' (Gyaneshwar et al., 2011) and 'Paraburkholderia' (Sawana et al., 2014)] has been introduced to accommodate most of the PBE species (Oren and Garrity, 2015a,b, 2017), while retaining the pathogens in Burkholderia sensu stricto. Most recently, a third genus, Robbsia was introduced to accommodate the phytopathogen previously referred to as B. andropogonis (Lopes-Santos et al., 2017).

Overall, the taxonomy of Burkholderia sensu lato remains in significant flux (Estrada-de los Santos et al., 2016). With their review of the group, Estrada-de los Santos et al. (2016) recognized two monophyletic groups [Groups A and B; A consists of Caballeronia and Paraburkholderia as circumscribed by Gyaneshwar et al. (2011) and Sawana et al. (2014), respectively, while Group B includes most of the notable human, animal and plant pathogens, as well as the so-called "B. cepacia complex"]. They showed that B. andropogonis (now Robbsia andropogonis) is separated into its own group, and they designated two socalled "Transition Groups" (i.e., 1 and 2; neither were supported as monophyletic and both contained mainly environmental species). Since then, Dobritsa and Samadpour (2016) have proposed the transfer of species in Transition Group 2 to a new genus. However, to complicate the issue, this new genus was named "Caballeronia" although its proposed usage is not synonymous with the one previously proposed by Gyaneshwar et al. (2011) for accommodating all the PBE isolates.

The proposals to split Burkholderia sensu lato were based almost entirely on evidence from 16S ribosomal RNA (rRNA) phylogenetic trees with limited and in some cases no statistical support (Gyaneshwar et al., 2011; Dobritsa and Samadpour, 2016; Eberl and Vandamme, 2016). Even phylogenies based on conventional multilocus sequence analysis (MLSA) using the combined sequence information for 4–7 genes (Gevers et al., 2005) produced phylogenies in which the major groups were not supported as monophyletic (Estrada-de los Santos et al., 2013). Also, the most comprehensive phylogenetic hypothesis to date (based on 21 conserved gene sequences) lacked sufficient representation across this diverse assemblage (Sawana et al., 2014). Thus, uncertainties remain regarding the genomic and evolutionary coherence of Burkholderia sensu lato and its lineages. This, in turn, blurs the boundaries of the Burkholderia sensu lato genera currently recognized and also casts doubt on the appropriateness and legitimacy of their taxonomic circumscriptions.

In this study, we aimed to resolve the relationships within Burkholderia sensu lato, particularly those pertaining to Paraburkholderia and Caballeronia, by making use of whole genome sequence data. For this purpose, we utilized all of the sequences for type strains (or appropriate representatives) available in the public domain. To increase representation of the so-called environmental species, we also determined the sequences for eight additional taxa via Phase III of the GEBA (Genomic Encyclopedia of Bacterial and Archaeal type strains) project (Whitman et al., 2015). These included the rhizobial species P. aspalathi (Mavengere et al., 2014) and P. diazotrophica (Sheu et al., 2013), and the soil bacteria P. hospita (Goris et al., 2002), P. phenazinium (Viallard et al., 1998), P. sartisoli (Vanlaere et al., 2008), P. terricola (Goris et al., 2002), as well as the plantassociated diazotrophic species P. caballeronis (Martínez-Aguilar et al., 2013) and P. tropica (Reis et al., 2004).

### MATERIALS AND METHODS

#### Whole-Genome Sequencing of Eight Paraburkholderia Type Strains

The eight type strains (P. aspalathi LMG 27731T , P. hospita LMG 20598T , P. diazotrophica LMG 206031T , P. phenazinium LMG 2247T , P. sartisoli LMG 24000T , P. terricola LMG 20594T , P. tropica LMG 22274T and P. caballeronis LMG 26416T ) were obtained from the Belgian Coordinated Collections of Microorganisms (University of Gent, Belgium). Routine growth of these bacteria in the laboratory and extraction of high quality genomic DNA were completed as described previously (Steenkamp et al., 2015). Whole genome sequencing was performed by the Joint Genome Institute (JGI) following standard protocols1 and using the Illumina HiSeq-2500 1TB platform with an Illumina 300 base pair (bp) insert standard shotgun library.

All raw sequences were filtered using BBDuk (Bushnell, 2017), which removes known Illumina artifacts, and PhiX. Reads with more than one "N" or with quality scores (before trimming) averaging less than 8 or reads shorter than 51 bp (after trimming) were discarded. The remaining reads were mapped to masked versions of human, cat and dog reference sequences using BBMap (Bushnell, 2017) and discarded if identity values exceeded 93%. The remaining reads were then assembled into contigs using Velvet version 1.2.07 (Zerbino and Birney, 2008) (the settings used were velveth: 63 –shortPaired and velvetg: –very clean yes –exportFiltered yes –min contig lgth 500 –scaffolding no –cov cutoff 10). The Velvet contigs were then used to generate 1–3 kbp simulated paired end reads using wgsim version 0.3.02 (the settings used were –e 0 –1 100 –2 100 –r 0 R 0 –X 0). We then assembled the quality filtered Illumina reads with the simulated read pairs using Allpaths-LG version r46652 (Gnerre and MacCallum, 2011) (the settings used were PrepareAllpathsInputs: PHRED 64 = 0 PLOIDY = 1 FRAG COVERAGE = 125 JUMP COVERAGE = 25 LONG JUMP COV = 50 and RunAllpathsLG: THREADS = 8 RUN = std shredpairs TARGETS = standard VAPI WARN ONLY = True OVERWRITE = True).

The standard JGI microbial genome annotation pipeline (Huntemann et al., 2015) was used to predict and annotate

1http://www.jgi.doe.gov

2https://github.com/lh3/wgsim

genes in each of the eight assembled genomes. For this purpose, we specifically used the Prodigal algorithm to identify proteincoding genes (Hyatt et al., 2010). Additional annotation was performed using JGI's Integrated Microbial Genomes (IMG) system (Markowitz et al., 2014).

### Sequence Datasets and Multiple Alignments

Protein-coding gene datasets were generated for the eight bacteria sequenced here, as well as all the Burkholderia sensu lato type strains (or suitable conspecific strains) for which whole genome sequences were available (Supplementary Table S1). This was achieved by using the EDGAR (Efficient Database framework for comparative Genome Analyses using BLAST score Ratios) server3 (Blom et al., 2016) to identify single-copy orthologous genes shared among all of the genomes examined. The respective amino acid and nucleotide sequences for each gene dataset were then batch-aligned using the Multiple Sequence Comparison by Log-Expectation (MUSCLE) (Edgar, 2004) iteration-based alignment tool implemented in CLC Main Workbench 7.6 (CLC Bio).