{'fraction_non_alphanumeric': 0.037993008721443455, 'fraction_numerical': 0.021256311570919107, 'mean_word_length': 4.801310938140107, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 1, 'https://': 2, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Comparative analysis of Lactobacillus gasseri from Chinese subjects reveals a new species-level taxa

Background Lactobacillus gasseri as a probiotic has history of safe consumption is prevalent in infants and adults gut microbiota to maintain gut homeostasis. Results In this study, to explore the genomic diversity and mine potential probiotic characteristics of L. gasseri, 92 strains of L. gasseri were isolated from Chinese human feces and identified based on 16 s rDNA sequencing, after draft genomes sequencing, further average nucleotide identity (ANI) value and phylogenetic analysis reclassified them as L. paragasseri (n = 79) and L. gasseri (n = 13), respectively. Their pan/core-genomes were determined, revealing that L. paragasseri had an open pan-genome. Comparative analysis was carried out to identify genetic features, and the results indicated that 39 strains of L. paragasseri harboured Type II-A CRISPR-Cas system while 12 strains of L. gasseri contained Type I-E and II-A CRISPR-Cas systems. Bacteriocin operons and the number of carbohydrate-active enzymes were significantly different between the two species. Conclusions This is the first time to study pan/core-genome of L. gasseri and L. paragasseri, and compare their genetic diversity, and all the results provided better understating on genetics of the two species.

multiple genomes, including core genome and variable genome. The core genome consists of genes presented in all strains and is generally associated with biological functions and major phenotypic characteristics, reflecting the stability of the species. And variable genome consists of genes that exist only in a single strain or a portion of strains, and is generally related to adaptation to particular environments or to unique biological characteristics, reflecting the characteristics of the species [17]. Pan-genomes of other Lactobacillus species [18], such as Lactobacillus reuteri [19], Lactobacillus paracasei [20], Lactobacillus casei [21] and Lactobacillus salivarius [22] have previously been characterized. The genetic knowledge and diversity of L. gasseri and L. paragasseri is still in its infancy. In addition, previous in silico surveys have reported that Lactobacilli harbour diverse and active CRISPR-Cas systems, which has 6-foldrate occurrence of CRISPR-Cas systems compared with other bacteria [23]. It is necessary to study CRISPR-Cas system to understand the adaptive immune system that protect Lactobacillus from phages and other invasive mobile genetic elements in engineering food microbes, and explore powerful genome engineering tool. Moreover, numerous bacteriocins were isolated from Lactobacillus genus, and these antimicrobials received increased attention as potential alternatives to inhibit spoilage and pathogenic bacteria [24]. A variety of strategies identify bacteriocin culture-based and in silico-based approaches, and to date, bacteriocin screening by in silico-based approaches have been reported in many research investigations [25].
In the current work, strains were isolated from fecal samples collected from different regions in China, and initially identified as L. gasseri by 16S rDNA sequencing. For further investigation, draft genomes of all the strains were sequenced by next generation sequencing (NGS) platform and analysed by bioinformatics to explore the genetic diversity, including subspecies/adjunct species, pan-genome, CRISPR-Cas systems, bacteriocin and carbohydrate utilization enzymes.

Strains and sequencing
Based on 16S rDNA sequencing, 92 L. gasseri strains were isolated from fecal samples obtained from adults and children from different regions in China, with 66 strains being obtained from adults and 26 from children (47 strains were isolated from females, 45 were isolated from males) ( Table 1). The draft genomes of all strains were sequenced using Next Generation Sequencing (NGS) technology and strains were sequenced to a coverage depth no less than the genome 100 ×, and using the genome of L. gasseri ATCC33323 and L. paragasseri K7 as reference sequences.

ANI values
ANI values calculation of Z92 draft genomes was carried out through pairwise comparison at the 95% threshold to further identify their species (Fig. 1). All of the 94 strains were classified into two groups, with 80 strains including L. paragasseri K7 (as type L. paragasseri strain) showing an ANI value range 97-99%, and the other group consisted of 14 strains including the type strain L. gasseri ATCC 33323 (as type L. gasseri strain) with an ANI range 93-94% compared with L. paragasseri. According to a previous report, L. gasseri K7 was reclassified as L. paragasseri based on whole-genome analyses [16], therefore, other 79 strains on the same group with L. paragasseri K7 were preliminarily identified as L. paragasseri, while the remained 13 strains on the other branch with L. gasseri ATCC33323 were identified as L. gasseri.

Phylogenetic analysis
To further verify the results from ANI and evaluate the genetic distance among strains, the phylogenetic relationships between L. paragasseri and L. gasseri were investigated. OrthoMCL was used to cluster orthologous genes and 1282 orthologues proteins were shared by all the 94 genomes. A robust phylogenetic tree based on 1282 orthologues proteins was constructed (Fig. 2). The results indicated that all the 94 strains could be positioned on two branches, in which 80 strains were on the same cluster with L. paragasseri K7 and the other 14 strains were on the cluster with L. gasseri ATCC33323. Surprisingly, all the strains on the cluster with L. gasseri or L. paragasseri were completely consistent with the results from ANI analysis. Therefore, it was confirmed that division of 92 strains isolated from Chinese subjects into two subgroups; 79 strains belong to L. paragasseri, and 13 strains to L. gasseri, is correct. The strains were randomly selected from the fecal samples, suggesting that L. gasseri and L. paragasseri had no preference to either male or female subjects nor region and age. Moreover, the house-keeping genes pheS and groEL were extracted from the genomes and neighbor-joining trees were built. The tree showed that 13 strains of L. gasseri were clustered in a single clade (Fig. 3), which was consistent with phylogenetic data based on orthologous genes. However, there were many branches in the L. paragasseri groups, which indicated a high intraspecies diversity among L. paragasseri and needs further investigation (Fig. 2, Fig. 3).

General genome features and annotation
The general information of the 80 genomes of L. paragasseri strains and 14 genomes of L .gasseri strains are summarized in Table 1. The sequence length of L. paragasseri ranged from 1.87 to 2.14 Mb, with a mean size of 1.97 Mb, and all 14 L. gasseri genomes had an average  Similarly, approximately 85% of L. gasseri ORFs were identified, while 15% were predicted to encode hypothetical proteins. The preference of the two species codons for the start codon were predicted, and the results showed that ATG, TTG and CTG in L. paragasseri with a calculated frequency percentage of 82.6, 10.3 and 7.1%, respectively, and 81.0, 11.7 and 7.4% in L. gasseri, respectively, suggesting that L. paragasseri and L. gasseri had a preference of using ATG as start codon [16].  To further analyse the genome-encoded functional proteins, the COG classification was performed for each draft genome. According to the results of the COG annotation, the genes were divided into 20 groups, and the details are shown in (Additional file 1: Table S1) and (Additional file 2: Table S2). The results revealed that carbohydrate transport and metabolism, defense mechanisms differed in different genomes of L. paragasseri, while L. gasseri showed only difference in defense mechanisms. Notably, due to draft genomes, the possibility of error from missing genes or incorrect copy number is significantly higher [28].

Pan/core-genome analysis
To analyze the overall approximation of the gene repertoire for L. paragasseri and L. gasseri in the human intestine, the pan-genomes of L. paragasseri and L. gasseri were investigated, respectively. The results showed that the pangenome size of all 80 strains of L. paragasseri amounted to 6535 genes while the pan-genome asymptotic curve had not reached a plateau (Fig. 4), suggesting that when more L. paragasseri genomes were considered for the number of novel genes, the pan-genome would continuously increase. Meanwhile, the exponential value of deduced mathematical function is > 0.5 (Fig. 4), these findings indicated an open pan-genome occurrence within the L. paragasseri species. L. paragasseri had a supragenome about 3.3 times larger than the average genome of each strain, indicating L. paragasseri constantly acquired new genes to adapt to the environment during evolution. The pan-genome size of the 14 strains of L. gasseri was 2834 genes, and the exponential value of deduced mathematical function is < 0.5, thus it could not be concluded whether its pan-genome was open or not.
The number of conserved gene families constituting the core genome decreased slightly, and the extrapolation of the curve indicated that the core genome reached a minimum of 1256 genes in L. paragasseri and 1375 genes in L. gasseri, and the curve of L. paragasseri remained relatively constant, even as more genomes Fig. 2 The phylogenetic tree based on orthologous genes. The red area was the L. gasseri cluster and the blue area was the L. paragasseri cluster. The purple circle indicated the strains isolated from infant feces and the gray indicated strains isolated from adults. The pink indicated strains from female and the green represent strains from male subjects were added. The Venn diagram represented the unique and orthologues genes among the 80 L. paragasseri strains. The unique orthologous clusters ranged from 3 to 95 genes for L. paragasseri and ranged from 8 to 125 genes for L. gasseri (Fig. 5). As expected, the core genome included a large number of genes for translation, ribosomal structure, biogenesis and carbohydrate transport and metabolism, in addition to a large number of genes with unknown function (Additional file 5: Figure S1). Identification and characterization of CRISPR in L. paragasseri and L. gasseri The CRISPR-Cas adaptive immunity system provided resistance against invasive bacteriophage or plasmid DNA such as some lytic bacteriophages in engineering food microbes, which consists of CRISPR adjacent to Cas genes. The presence of Cas1 proteins was used to determine the presence or absence of CRISPR-Cas systems, and Cas1 was found among the 39 strains of L. paragasseri, and 13 strains of L. gasseri. The occurrence of Cas1 genes in L. paragasseri and L. gasseri showed differences, in that 12 strains of L. gasseri consisted of two Cas1 genes, and the second Cas1 gene was located in a different region constituting a second putative CRISPR locus. Meanwhile, Cas2 and Cas9 were widespread across the two species, while Cas3, Cas5, Cas6 and Cas7 only occurred in L. gasseri. According to previous method of classification of the CRISPR subtypes, 52 Type II-A systems were detected in all the strains of L. gasseri and 39 strains of L. paragasseri, whereas the Type I-E system only occurred in 12 strains of L. gasseri except FHNFQ57-L4, indicating that subtype II-A was the most prevalence both in L. paragasseri and L. gasseri. The phylogenetic analyses performed with Cas1, Cas2 and Cas9 from the two species showed L. paragasseri was clearly distinct from L. gasseri (Fig. 6). Strikingly, phylogenetic tree based on Cas1 and Cas2 proteins revealed that clusters consisted of only the second Cas1 and Cas2 proteins in Type I-E systems in L. gasseri, and the Cas1 and Cas2 proteins in subtype II-A systems in both L. paragasseri and L. gasseri were clustered in two groups. From this perspective, CRISPR-Cas could be used as an indicator to distinguish L. paragasseri and L. gasseri. Moreover, phylogenetic analysis of Cas9 indicated that the cluster was consistent with Cas1 and Cas2, indicating that co-evolutionary trends happened in CRISPR systems.
The features of all 60 CRISPR loci identified in L. paragasseri and L. gaseseri genomes are summarized in Table S3. The length of DRs were 36 nucleotides (nt) in 36 strains of L. paragasseri except FJSCZD2-L1, FHNFQ53-L2 and FHNXY18-L3, which had DR sequences with 26 nt. The 5′-terminal portion of DRs in L. paragasseri were composed of G (T/C) TTT and the DRs were weakly palindromic. The putative RNA secondary structure of the DRs in L. paragasseri contained two small loops (Fig. 7). The DRs of L. paragasseri shared two variable nucleotides at the 2nd and 29th site (C/T), and the difference affected the RNA secondary structures (Fig. 7). While two CRISPR loci in L. gasseri had different DR sequences and varied in length and content, in which most of them were 28 nt whereas L. gasseri FHNFQ56-L1 and FHNFQ57-L4 had a same DR as L. pargasseri (Additional file 3: Table S3). Further, spacer contents were uncovered for L. paragasseri and L. gaseseri, ranging from 3 to 22 CRISPR spacers (Additional file 3: Table S3). The number of spacers in L. paragasseri and L. gasseri were variable and it provided information about the immunity record.

Distribution of Bacteriocin operons
Identifying bacteriocins in vitro can be a challenging task, however, in silico analysis of genomes for presence of bacteriocin operons could make screening bacteriocin efficient. BAGEL was used to identify the potential bacteriocin operons in the current study. Three hundred twenty-three putative class II bacteriocin and 91 putative class Bacteriolysin (formerly Class III Bacteriocins) operons were identified in all 92 genomes (Additional file 4: Table S4). Class II bacteriocins are small heat stable peptides further subdivided IIa, IIb, IIc and IId based on the structure and activity of the peptides [25]. L. paragassseri genomes contained various bacteriocins including Class IIa (pediocin), Class IIb (gassericin K7B and gassericin T), Class IIc (acidocin B and gassericin A), Class IId (bacteriocin-LS2chaina and bacteriocin-LS2chainb), and Bacteriolysin, whereas all the strains of L. gasseri only encoded bacteriocin-helveticin-J (Bacteriolysin) except L. gasseri FHNFQ57-L4, which contained both bacteriocin-helveticin-J and pediocin operons. Interestingly, gassericin K7B and gassericin T operons co-occured in 43 strains of L. paragasseri, and bacteriocin-LS2chaina and bacteriocin-LS2chainb cooccured in 67 strains of L. paragasseri. Sixteen gassericin A, 31 acidocin B, 69 pediocin and 78 bacteriocinhelveticin-J operons were also predicted in L. paragasseri indicating that helveticin homolog operons were more frequent than other operons. In addition, only one enterolysin A operon was found in L. paragasseri FHNFQ29-L2, FGSYC41-L1 and L. paragasseri FJSWX6-L7 contained a helveticin J operon.
The glycobiome of L. paragasseri and L. gasseri The earliest classifications of lactobacilli were based on their carbohydrate utilization patterns. In the current study, carbohydrate-active enzymes were analyzed by HMMER-3.1 and identified through carbohydrate-active enzyme (Cazy) database. Nineteen glycosyl hydrolase (GH) families, 7 glycosyl transferase (GT) families and 5 carbohydrate esterase (CE) families were predicted for each genome, and the distribution and abundance of GH, GT, CE family genes across the L. paragasseri and L. gasseri were showed by heatmap (Fig. 8).
Twelve strains of L. gasseri formed a single clade using hierarchical clustering method (Fig. 8). Both species of L. gasseri and L. paragasseri appeared to contain a consistent GH65, GH73 and GT51 (murein polymerase) families, while GH42 family (β-galactosidase and α-L- The sequence of consensus DR sequences in L. gasseri strains. The height of the letters indicates the frequency of the corresponding base at that position. c-e Predicted RNA secondary structures of CRISPR DR in L. paragasseri. f-g Predicted RNA secondary structures of the CRISPR DR in L. gasseri arabinopyranosidase) was only found in L. paragasseri. Additionally, the gene number of GT8 (α-transferase) family in L. gasseri was less than that in L. paragasseri. The results revealed that carbohydrate utilization patterns of L. gasseri differed from L. paragasseri. Carbohydrate-active enzymes abundance in L. paragasseri showed high diversity, but the difference was not a result of gender and age difference, and may be associated with diet habits of the host individual. Diversity do not correlate with gender and age and could be result of sugar diet habits of the host individual.

Discussion
NGS technologies have made sequencing easier to get high-quality bacterial genomes, and provides the possibility to better understand the genomic diversity within some genus [29]. In this study, genome sequences for 92 strains from human feces, which were preliminary identified as L. gasseri by 16S rDNA sequencing, combined with two publicly available genomes L. gasseri ATCC33323 and L. paragasseri K7, were further analysed. ANI values of 94 draft genomes were calculated Fig. 8 The distribution and number of GH, CE and GT family genes. Gene copy number was indicated by color ranging from green (absent) to red. The strain number in red and in black indicated L. gasseri and L. paragasseri, respectively through pairwise comparison at the 95% threshold, together with phylogenetic analysis based on orthologous genes and house-keeping genes (pheS and groEL) were performed to ensure the species affiliations and eliminated the mislabeled genomes only using ANI [30]. Seventy-nine strains were determined as L. paragasseri, and the remained 13 (14%) strains were L. gasseri, revealing that the most (86%) of isolates initially identified as L. gasseri by 16S rDNA sequencing were L. paragasseri. The current results were highly in line with previous publication by Tanizawa and colleagues [16], in which they reported that a large portion of genomes currently labelled as L. gasseri in the public database should be re-classified as L. paragasseri based on whole-genome sequence analyses as well. All those results indicated that L. gasseri and L. paragesseri are sister taxon with high similarity but not the same species, and the cultivable "L. gasseri" isolated from environment actually contained both L. gasseri and L. paragasseri species, which might be the reason for the high intraspecies diversity among "L. gasseri" exhibited. Meanwhile, groEL, a robust singlegene phylogenetic marker for Lactobacillus species identification [31], could serve as a marker to distinguish L. paragasseri and L. gasseri. Our current results provide a basis for distinguishing the two species by genotype. L. gasseri and L. paragasseri had no preference to colonize the female or male subjects, and the strains distribution had no trend on age neither infants nor adults. Nevertheless, a high intraspecies diversity in L. paragasseri may be caused by diet habits, health condition and others, which needs further research.
In general, the genome size of L. paragasseri and L. gasseri were smaller than other Lactobacillus species, which had an average size 1.96 Mb, while other Lactobacillus had a genome approximately 3.0 Mb, such as L. paracasei [20], L. casei [21], Lactobacillus rhamnosus [32]. Additionally, the G + C contents in L. paragasseri (34.9%) and L. gasseri (34.82%) were lower than that in other Lactobacillus species. For instance, average G + C contents were 38.96% in L. reuteri [19], 46.1-46.6% in L. casei, 46.5% in L. paracasei [20], and 46.5-46.8% in L. rhamnosus [33], and the average G + C content among lactobacilli genera is estimated at 42.4%. As previously found in bifidobacterial genomes, that the preferred start codon was ATG, also analysis of start codons in L. paragasseri and L. gasseri showed that they preferably used ATG as start codon [34].
Pan-genomes of L. paragasseri and L. gasseri were analysed, and the pan-genome size of the 80 strains among L. paragasseri and 14 strains of L. gasseri plus currently genome public strains of L. gasseri ATCC33323 and L. paragasseri K7 were 6535 and 2834 genes, respectively, and the core genomes were 1256 and 1375 genes, respectively, suggesting that open pan-genome within the L. paragasseri species and its pan-genome will increase if more L. paragasseri genomes were considered for the number of novel gene families and an open pan-genome implies that gene-exchanging within a species is higher [28]. But it could not conclude whether the pan-genome of L. gasseri was open or not because of limited number of sequenced genomes.
It has been reported that lactic acid bacteria are enriched resource for Type II CRISPR systems [35] and some previous studies on L. gasseri CRISPR-Cas reported that L. gassseri harboured type II-A CRISPR-Cas system with diversity in spacer content, and confirmed functionality [36]. However, the former results on "L. gasseri" might not be the real L. gasseri, as L. paragasseri was distinguished from L. gasseri recently, which might be mixed in the previous research. In the current result, L. gasseri and L. paragasseri were distinguished and separately, then were loaded for CRISAP-Cas analysis, respectively. The results showed that 39 of 79 L. paragasseri strains carried Type II systems and all the strains of L. gasseri harboured Type II and Type I CRISPR-Cas system (except FHNFQ57-L4), implying that both L. paragasseri and L. gasseri are main candidates for gene editing and cleavage of lytic bacteriophages in food industry. In the current study we found that Cas1, Cas2 and Cas9 were widespread across both L. paragasseri and L. gasseri species, and the L. gasseri species had a second Cas1 and Cas2, while the second Cas1 and Cas2 were clustered in a single clade through phylogenetic analyses. Similarity, the Cas9 gene was different between the two species, suggesting that CRISPR-Cas could provide a unique basis for resolution at the species-level [37], and the CRISPR-Cas systems may contribute to the evolutionary segregation [33].
It has been reported that L. gasseri produces variety of bacteriocin to inhibit some pathogens. Screening bacteriocin in vitro was complex and difficult while in silico analysis could make it rapid, generally using BAGEL to identify the potential bacteriocin operons. In the current study, most of the L. gasseri strains only had a single bacteriocin operon (Bacteriocin_helveticin_J), while L. paragasseri showed a variety of bacteriocin operons belonging to class II such as gassericin K7B, gassericin T and gassericin A. With the current results, although bacteriocin was not separated and verified in vitro, we presume that the strains with the high-yielding of bacteriocin, which was commonly known as L. gasseri, should actually be L. paragaseri rather than L. gasseri. For instance, previously L. gasseri LA39 was reported to produce gassericin A [38] and L. gasseri SBT2055 [39] could produce gassericin T, according to our results, they might belong to L. paragasseri species instead of L. gasseri. To confirm our hypothesis, more L. gasseri strains should be isolated and screened for bacteriocin to verify.
In order to investigate L. paragasseri and L. gasseri carbohydrate-utilization capabilities, carbohydrate-active enzymes were predicted for all the strains and these families have predicted substrates and functional properties for each strain. Analyzing the abundance of Cazy revealed that carbohydrate utilization patterns of L. gasseri significantly distinguished with L. paragasseri in genotype, which provided foundation for fermentation experiment with unique carbon sources. Moreover, 10.83% core genes had predicted function of carbohydrate transport and metabolism, which is the reason for strains diversity and separation.

Conclusion
Ninety two strains isolated from Chinese subjects were initially identified as L. gasseri by 16S rDNA sequencing, while based on whole-genome analyses they were reclassified. According to ANI values and phylogenetic analysis based on both orthologous genes and housekeeping genes, 13 strains and 79 strains were reclassified as L. gasseri and L. paragasseri, respectively, which revealed a new species-level taxa from Chinese subjects. Pan-genome structure for L. paragasseri was open, meanwhile, L. paragasseri had a supragenome about 3.3 times larger than the average genome size of individual strains. After species reclassification, genetic features CRISPR-Cas systems, bacteriocin, and carbohydrateactive enzymes were analysed, revealing differences in the genomic characteristics of L. paragasseri and L. gasseri strains isolated from human feces and mine potential probiotic characteristics in the two species. To our knowledge, this is the first study to investigate pan/coregenome of L. gasseri and L. paragasseri, compared the genetic features between the two species.

Methods
Isolation of strains, genome sequencing and data assembly Ninety-two strains isolated from adult and infant feces from different regions in China were listed in Table 1. Strains were selected in Lactobacillus selective medium (LBS) [4] and incubated at 37°C in an anaerobic atmosphere (10% H 2 , 10% CO 2 , and 80% N 2 ) in a anaerobic workstation (AW400TG, Electrotek Scientific Ltd., West Yorkshire, UK) for 18-24 h and 16S rRNA genes were sequenced for species identification. All the identified L. gasseri strains were stocked at -80°C in 25% glycerol [40]. Draft genomes of all the 92 L. gasseri strains were sequenced via Illumina Hiseq× 10 platform (Majorbio BioTech Co, Shanghai, China), which generated 2 × 150 bp paired-end libraries and construct a paired-end library with an average read length of about 400 bp. It used double-end sequencing, which single-ended sequencing reads were 150 bp. The reads were assembled by SOAPde-novo and local inner gaps were filled by using the software GapCloser [41]. Two publicly available genomes (L. gasseri ATCC33323 [26] and L. gasseri K7 [27]) from National Centre for Biotechnology Information (https://www.ncbi.nlm.nih.gov/) were used for comparison and the latter one has recently been reclassified as L. paragasseri [16].

Phylogenetic analyses
All the genomic DNA were translated to protein sequences by EMBOSS-6.6.0 [44]. OrthoMCL1.4 was used to cluster orthologous genes and extracted all orthologous proteins sequences of 94 strains. All orthologous proteins were aligned using MAFFT-7.313 software [45] and phylogenetic trees were constructed using the python script (https://github.com/jvollme/fasta2phylip) and the supertree was modified using Evolgenius (http:// www.evolgenius.info/evolview/). The house-keeping genes, pheS [46] and groEL [47], were extracted from the genomes using BLAST (Version 2.2.31+) [48], and the multiple alignments were carried out through Cluster-W (default parameters), and the single gene neighborjoining trees were built by MEGA 6.0 [49], with bootstrap by a self-test of 1000 resampling.

Pan/core-genome analysis
Pan-genome computation for L. paragasseri and L. gasseri genomes was performed using the PGAP-1.2.1, which analysed multiple genomes based on protein sequences, nucleotide sequences and annotation information, and performed the analysis according to the Heap's law pan-genome model [17,52]. The ORF content of each genome was organised in functional gene clusters via the Gene Family method and a pan-genome profile was then built.

Bacteriocin identification
The bacteriocin mining tool BAGEL3 was used to mine genomes for putative bacteriocin operons [55]. To determine the bacteriocins pre-identified by BAGEL3, BLASTP was secondly used to search each putative bacteriocin peptide against those pre-identified bacteriocins from BAGEL screening, and only the consistent results from both analysis were recognized as truly identified bacteriocin.
The L. gasseri glycobiome Analysis of the families of carbohydrate-active enzymes was carried out by using HMMER-3.1 (http://hmmer. org/) and with below a threshold cutoff of 1e-05. The copy number of the verified enzymes were summarized in a heatmap with hierarchical clustering method and Pearson distance [35].