Genomic Characterization of Lactiplantibacillus plantarum Strains Possessing Differential Antiviral Immunomodulatory Activities

: Lactiplantibacillus plantarum strains are used in the food industry for their probiotic prop ‐ erties. Some of these bacteria have immunomodulatory effects on the host and are able to improve resistance against different pathogens, including viruses. However, to date, the bacterial genes in ‐ volved in the immunomodulatory effect are not known. In this work, the complete genomes of L. plantarum MPL16, CRL1506, CRL681 and TL2766 were used to perform comparative genomics with the aim of identifying the genes involved in their differential immunomodulatory effects. L. planta ‐ rum WCFS1, a strain with proven probiotic activity, was also used for comparisons. The analysis of the genes involved in the metabolic pathways of the five strains did not reveal differences in the metabolism of amino acids, lipids, nucleotides, cofactors and vitamins, nor in the genes associated with energy metabolism or the biosynthesis of lipoproteins and teichoic acids. However, differences were found between the five strains when considering carbohydrate metabolism pathways, partic ‐ ularly in the presence/absence of glycosylhydrolases and glycosyltransferases. In addition, a great variability was detected in the predicted surface proteins of each L. plantarum strain. These results suggest that the surface molecules expressed in the different strains of L. plantarum could be in ‐ volved in their differential ability to modulate the innate antiviral immune response.


Introduction
Lactiplantibacillus plantarum has one of the largest genomes known among lactic acid bacteria [1] and is therefore able to survive in a wide range of environmental niches including the gastrointestinal tract, easily colonizing the intestine of humans and other mammals [2,3]. In addition, many strains of this bacterial species have shown to possess beneficial properties for the host, including their ability to beneficially modulate the immune system [4]. Due to these properties, L. plantarum is considered one of the most widely used bacterial strains in the food industry both as a starter culture and as a probiotic [5].
Previously, we performed in vitro studies in porcine intestinal epithelial (PIE) cells and in vivo studies in mice with different L. plantarum strains (MPL16, CRL1506, CRL681 and TL2766) and we demonstrated that those strains possess a differential ability to modulate the respiratory and intestinal innate antiviral immune responses [6]. Quantitative and qualitative differences were found in the expression of genes with immunological functions in PIE cells when they were treated with the distinct L. plantarum strains before the challenge with the Toll-like receptor 3 (TLR3) agonist polyinosinic:polycytidylic acid [poly(I:C)]. The transcriptional profile performed in vitro on PIE cells challenged with poly(I:C) allowed the selection of a new immunobiotic strain, L. plantarum MPL16, with the ability to stimulate in vivo antiviral immunity not only in the intestinal mucosa, but also in the respiratory tract, when administered orally [6]. On the other hand, our studies showed that the immunomodulatory strain L. plantarum CRL1506 only stimulates antiviral immunity at the intestinal level, whereas the CRL681 and TL2766 strains did not present immunomodulatory effects in the context of the response mediated by the activation of TLR3 [6]. Therefore, having strains of L. plantarum with different abilities to modulate intestinal and respiratory antiviral immunity provides an interesting opportunity to investigate the bacterial molecule(s) that could be involved in their differential immunomodulatory activities.
With the rapid development of next-generation sequencing techniques, the complete genomes of multiple strains of L. plantarum have been sequenced. This scientific-technological advance has provided a better understanding of the relationships between genotypes and functions [7,8]. Comparative genomics allows efficient analysis of intra-strain differences in the genomes of multiple L. plantarum strains, leading to better functional analysis at the molecular level [7][8][9]. These studies revealed that L. plantarum exhibits a high conservation of genes involved in the synthesis or degradation of proteins and lipids, whereas it presents a high diversity of genes involved in the transport and catabolism of sugars [10,11] as well as in the bacterial surface-expressed molecules [1,3,12], which determines the ability of each strain to interact with the host.
Considering that the complete genomes of L. plantarum MPL16, CRL1506, CRL681 and TL2766 has been sequenced, the aim of this work was to carry out comparative and functional genomics studies in order to identify the bacterial gene(s) responsible for the distinct immunomodulatory properties of each strain. For this purpose, special emphasis was given to the study of genes with potential immunomodulatory function, including surface and secreted proteins, adhesion factors, and genes involved in the biosynthesis of exopolysaccharides (EPS), lipoproteins, and teichoic acids.

General Genomic Characteristics of L. plantarum
The complete genomes of L. plantarum MPL16, CRL1506 and CRL681 were sequenced with the Illumina MiSeq platform and uploaded to the NCBI repository. Genomes are available under the accession numbers indicated in Table 1. For the comparative genomic analysis, L. plantarum WCFS1, with recognized probiotic activity [13], and L. plantarum TL2966, which has not possess immunomodulatory properties [14,15] were also incorporated. The general genomic characteristics of the L. plantarum evaluated in this work are shown in Table 1. The mean whole genome size and GC content of the strains were 3.26 Mb and 44.2%, respectively. This is in agreement with recent results obtained from the genomic analysis of 127 strains of L. plantarum which revealed a size and GC content of 3.32 Mb and 44.5%, respectively [1]. These data are also in line with the study reported by Mao et al. [11], describing that within the complete genomes of 114 strains of L. plantarum, genomes have a size between 2.94 and 3.90 Mb and a GC content of 44.1 to 46.5%. Among the five strains selected for this work, L. plantarum CRL681 presented the highest number of protein-coding genes (3081) whereas the lowest number (2918) was found for the CRL1506 strain.
Using the BlastKOALA tool, the functional characterization of L. plantarum genomes was carried out according to the KEGG database. As observed in Figure 1, no remarkable differences were observed in the number of genes associated with the metabolism pathways of carbohydrates (216 ± 3), amino acids (136 ± 2), lipids (40 ± 3), nucleotides (66 ± 2), glycans (136 ± 2), cofactors and vitamins (70 ± 3) or energy metabolism (68 ± 2) among the MPL16, CRL1506, CRL681, WCFS1 and TL2766 strains. L. plantarum CRL1506 presented a lower number of genes associated with the metabolism of terpenoids and polyketides as well as for the biosynthesis of other secondary metabolites compared with the other strains ( Figure 1). There were also no significant differences in the numbers of genes associated with "genetic information processing", "environmental information processing" or "cellular processes" when the five strains of L. plantarum were compared ( Figure 2). The only exceptions were the numbers of genes related to "folding, sorting and degradation" and "transduction signal" for the CRL1506 strain which were lower than those found in the other lactobacilli. The MPL16, CRL1506, CRL681, WCFS1 and TL2766 strains were then compared by a Venn diagram using the orthologous genes ( Figure 3). The analysis revealed that the five strains of L. plantarum share 2282 genes, whereas 194, 53, 106, 69 and 62 unique genes were detected for MPL16, CRL1506, CRL681, WCFS1 and TL2766, respectively. Among the unique genes for L. plantarum MPL16, two bacterial proteins with immunoglobulin-like (Ig-like) domains (OG0003093 and OG0002455), a PTS sugar transporter (OG0003210), an N-acetyl transferase of the GNAT family (OG0003311), two glycosyltransferases (GT) (OG0003314, OG0003357), two glycosylhydrolases (GH) (OG0003274, OG0003356), one lysophospholipase (OG0003265), the stress response membrane protein GlsB/YeaQ/YmgE (OG0003255), and several proteins associated with the type I-E CRISPR system (OG0003362, OG0003363, OG0003364, OG0003365, OG0003366, OG0003367, OG0003368, OG0003369) were detected (Table S1).
Among the seven genes shared by the MPL16 and CRL1506 strains, a polysaccharide pyruvyl transferase (OG0002887) was found, whereas among the 15 genes shared by the CRL1506 and CRL681 strains, a family 2 glycosyltransferase (OG0002860) was noticed ( Figure 3, Table S1).
Considering that the Venn diagram analysis revealed differences in the content of GH and GT, the difference between the strains with respect to the number of genes for GH and GT was further investigated in detail. Previous genomic analyzes found that L. plantarum has six main families of enzymes involved in carbohydrate metabolism: GH, GT, carbohydrate esterases, carbohydrate-binding enzymes, auxiliary active enzymes, and polysaccharide lyases [11]. The study also found that the most abundant sugar metabolism genes in the L. plantarum strains belonged to the GH family, followed by the GT family. Bearing in mind that these enzymes affect the ability of the strains to adapt to their environments, we then proceeded to study the number of genes that code for these enzymes. As observed in Figure 4, genes for several GT and GH families were detected in the genomes of L. plantarum MPL16, CRL1506, CRL681, WCFS1 and TL2766. The most abundant GT families in L. plantarum were GT2 and GT4. The highest number of genes for GT2 was found in the CRL1506 and CRL681 strains, whereas L. plantarum WCSF1 had the highest number of genes for GT4. The WCSF1 strain was also the one with the highest number of genes for GT111 and the only one with GT32 ( Figure 4). On the other hand, the analysis of GH gene numbers showed that GH1 and GH13 were the most abundant in the L. plantarum strains studied. The highest number for GH1 was found in the WCSF1 strain, whereas L. plantarum MPL16 had the lowest number of genes for GH13 ( Figure 4). The MPL16 strain also presented a lower number of genes for the GH1, GH25, and GH65 families compared with the other strains evaluated and it did not contain genes for GH31 or GH38 families. These results show a variability between the strains when considering genes associated with carbohydrate metabolism, in line with previous genomic studies [10,11].
Unique genes or genes shared between specific strains detected in our comparative genomic analysis are involved in bacterial cell division, membrane biosynthesis, cell wall peptidoglycan biosynthesis and catabolism, EPS biosynthesis, or they encode for proteins that are secreted or expressed on the cell surface. These results suggest that the cell wall and surface molecules expressed in the different strains of L. plantarum could be involved in their differential ability to modulate the innate antiviral immune response.

Study of "Probiotic Markers" Genes
Some studies have proposed a set of genes involved in stress resistance, active metabolism in the host, adhesion, and immunomodulation as genes with probiotic functions [13,16,17]. Based on such studies, Carpi et al. [1] recently created an updated list of "probiotic marker" genes including genes responsible for resistance to stress (acid, osmotic, oxidative, temperature), bile salt hydrolase activity, adhesion capacity, and intestinal persistence. A total of 75 probiotic marker genes have been reported, of which approximately 70% correspond to genes located in the core/soft core genome (Table S2), whereas 12 genes are located in the shell/cloud genome ( Table 2). We next analyzed the presence/absence of these probiotic marker genes in the genomes of the five strains evaluated in this work. As expected, most of the core genome probiotic marker genes were found in L. plantarum MPL16, CRL1506, CRL681, TL2766 and WCFS1 ( Figure 5). The celB (cellobiose-specific PTS system component) and gpmB (histidine phosphatase family protein) genes were not detected in the WCFS1 and CRL681 genomes, respectively. The analysis also revealed that L. plantarum MPL16 does not possess the celB, eno2 (enolase 2), oppA1 (OppA1 oligopeptide-binding protein), oppA2 (OppA2 oligopeptidebinding protein), nor treA (trehalose-6-phosphate hydrolase) genes belonging to the core genome of probiotic genes proposed for this bacterial species ( Figure 5). It should be noted that the treA and celB genes have been associated with persistence in the gastrointestinal tract, whereas eno2, oppA1 and oppA2 are involved in resistance to bile salts [1].
In addition, when the analysis was carried out considering the probiotic marker genes of the cloud genome, an absence of several of them was found in the L. plantarum CRL1506, CRL681, TL2766 and WCFS1 genomes ( Figure 6). The genes bshA (bile salt hydrolase), clpP (Clp-ATP-dependent protease), gbuB (glycine/carnitine transport protein GbuB), gla2 (aquaporin), oppA4 (oligopeptide-binding protein OppA4), and xylA (xylose isomerase) were not detected in any of the four genomes. The study also showed that L. plantarum MPL16 does not have the bshA, clpP, gla2, oppA4 and xylA genes. The absence of cbh (coliglycin hydrolase) and oppA3 (OppA3 oligopeptide binding protein) was also observed in the MPL16 strain, whereas it presented the gbuB gene, which is absent in L. plantarum CRL1506, CRL681, TL2766 and WCFS1 ( Figure 6). The clpP gene was associated with resistance to acid stress, whereas gla2 was associated with resistance to osmotic stress [1]. On the other hand, bshA, cbh, oppA3 and oppA4 are involved in resistance to bile salts [1]. The study of the presence/absence of the probiotic marker genes proposed for L. plantarum [1] revealed a notable difference between the MPL16 strain and the other strains evaluated. However, the analysis was not able to clearly discriminate the other immunomodulatory strains from the non-immunomodulatory strains evaluated in this work. Considering these results, a comparative study focused on the molecules of the bacterial surface that could be responsible for the interaction of L. plantarum strains with the host's immunological receptors was carried out.

Study of Genes Associated with the Synthesis of Exopolysaccharides
It was reported that L. plantarum might have four clusters for EPS biosynthesis [12]. To analyze whether the MPL16, CRL1506, CRL681 and TL2766 strains had these EPS biosynthesis-related genes in their genomes, a sequence comparison was made using the genes described in L. plantarum WCFS1 as a reference ( Figure 7). The genes corresponding to the eps1 cluster (or cps1) were not found in the genomes of any of the L. plantarum studied, except for the eps1C and eps1D genes that were found in the CRL1506 strain and the eps1H gene that was observed in L. plantarum CRL681. Figure 7. Analysis of the presence/absence of genes involved in the biosynthesis of exopolysaccharides (EPS) in Lactiplantibacillus plantarum. The genomes of L. plantarum strains TL2766, MPL16, CRL681 and CRL1506 were studied and compared with the probiotic strain WCFS1 as reference. Four clusters for EPS biosynthesis (eps1, eps2, eps3 and eps4) and the rfb cluster involved in the incorporation of rhamnose to EPS are shown.
The first five genes of the eps2 (or cps2) cluster have been reported to be highly conserved among L. plantarum strains, including the eps2A, esp2B, and eps2C genes (also called wzd, wze, and wzh, respectively) that are involved in regulation of EPS biosynthesis dependent on tyrosine phosphoregulation [12,20]. In line with these previous studies, it was observed that L. plantarum MPL16, CRL1506, CRL681 and TL2766 possess the three genes eps2A, esp2B and eps2C (Figure 7). The other two conserved genes are eps2D, which encodes an epimerase enzyme that catalyzes the epimerization of UDP-acetylglucosamine to UDP-N-acetylgalactosamine, common precursors for EPS biosynthesis, and eps2E, which encodes a glycosyltransferase that mediates the initiation of synthesis of the repeated unit of the EPS [21]. Interestingly, these two genes were detected in L. plantarum CRL1506 and TL2766 but not in strains MPL16 or CRL681 (Figure 7). The other components of the eps2 cluster comprising eps2F-J have been shown to possess very limited sequence homology between the corresponding regions from different strains [12]. The eps2G, eps2H and eps2J genes of L. plantarum WCFS1 were not detected in any of the strains studied whereas eps2F was only observed in CRL1506 and eps2I in CRL681 and TL2766 (Figure 7).
Sequence and homology analysis revealed that eps3 is relatively conserved among different strains of L. plantarum, but that its genetic arrangement was different from that of eps2 [12]. It was observed that L. plantarum CRL1506, CRL681 and TL2766 possess the complete eps3 cluster, whereas in strain MPL16 the EPS biosynthesis membrane proteins cps3F (lp_1222) and cps3G (lp_1224) were not detected, nor was the O-acetyltransferase cps3I (lp_1226) (Figure 7). The eps3A gene (lp_2115) was also not observed in the MPL16 genome. On the other hand, it has been reported that the organization of the eps4 (or cps4) cluster in L. plantarum WCFS1 and NUC116 is similar to that of eps2 with esp4A-J genes [12]. The bioinformatic analysis carried out in this work showed that L. plantarum CRL1506, CRL681 and TL2766 have the complete cluster, except for the cps4A gene that codes for the EPS chain length regulatory protein (lp_2108), which is absent in the TL2766 strain. None of the eps4 cluster genes were detected in MPL16 ( Figure 7).
It has been described that some strains of L. plantarum such as ST-III possess a group of conserved rmlACBD genes also called rfb ( Figure 7). These genes produce a rhamnose precursor that is essential for the biosynthesis of some EPS [22]. Consistently, rhamnose has been detected as a component of surface EPS in L. plantarum WCFS1 [20]. Five of the four genes of this cluster were not detected in the genome of L. plantarum CRL1506, whereas the lp_1187 gene was not detected in the MPL16 strain. Both L. plantarum CRL681 and TL2766 displayed the complete rfb cluster in their genomes ( Figure 7).

Study of Extracellular Proteins
The extracellular proteins, which together constitute the secretome, are involved in variable processes such as adherence to the host, recognition, degradation and absorption of nutrients, as well as signal transduction. It was suggested that extracellular and surfaceexposed proteins play important roles in the various interactions that lactobacilli establish with their environments [23,24]. Therefore, we studied the genes that code for extracellular proteins as well as for their transport systems in L. plantarum MPL16, CRL1506, CRL681 and TL2766.
The Sec system is the main modulator of protein transport across the cell membrane in Gram-positive bacteria [27]. Most of the genes of the Sec system were identified in the genomes of L. plantarum MPL16, CRL1506, CRL681 and TL2766 when the comparative analysis was performed with the WCFS1 strain ( Figure 8). All the genomes analyzed presented the genes for secA (ATP-dependent motor protein), secYEG (membrane channelforming complex) and two yidC genes (membrane protein insertases). The Sec system also has two proteins encoded by the ftsY and ffh genes, which bind to the cell membrane and interact directly with the secYEG channel, forming a binding pocket for the signal sequence that allows the preprotein to be recognized and presented to the secYEG channel [27]. Both genes, ftsY and ffh, were detected in the genomes of MPL16, CRL1506, CRL681 and TL2766 (Figure 8). It was established that this system also has three genes (secD, secF and YajC) that make up the trimeric complex secDF-YajC, which directs the motor protein SecA [27]. The secD or secF genes were not detected in the genomes of the L. plantarum in this study, which is in line with previous reports for the WCFS1 strain [26] as well as for other probiotic strains such as L. plantarum NCU116 [12]. On the other hand, it has been shown that YajC can bind to YidC and thus control the movement of the polypeptide chain through the secYEG channel, suggesting that YajC fulfills the function of preprotein translocation in L. plantarum NCU116 and WCFS1 [12] as well as in the strains of this work. The FPE pathway allows uptake of exogenous DNA across the plasma membrane [28] whereas holins are small integral membrane proteins associated with permeabilizing lytic enzymes found in most lactobacilli [29]. Recently it was suggested that these two protein translocation systems could complement the Sec system in protein secretion function in L. plantarum [12]. Studies described that L. plantarum WCFS1 [26] and NCU116 [12] possess genes for the main components of the FPE system including the comGA-GC operon and the comC homologue, but to lack the comGD and comGF genes. Similarly, the L. plantarum strains analyzed in this study had the comGA-GC and comC genes ( Figure 8) but not the comGD and comGF genes. Single copies of the comE and comF genes were identified in other species of lactobacilli, which are also involved in the DNA uptake process [30]. Similar to the WCFS1 strain, L. plantarum CRL1506, CRL681 and TL2766 presented the comEA, comEB, comEC, comFC and comFA genes ( Figure 8). However, comEC and comFA were not detected in the genome of L. plantarum MPL16. The comEC gene codes for a channel protein involved in DNA binding and uptake whereas the comFA gene codes for a membrane-associated DNA-dependent ATPase. The genes encoding holins were found in the genomes of MPL16, CRL1506, CRL681 and TL2766 (Figure 8). One holin gene copy was found in each strain (Table S3) as described for L. plantarum WCFS1 [26]. It should be noted that this bacterial species may have more than one copy of this gene, as has been reported for L. plantarum NCU116, which has three copies [12].
The extracellular proteins in the L. plantarum strains evaluated in this work were identified by bioinformatic analysis taking into account their subcellular location (secreted or expressed on the surface), their type of anchorage and the secretory mechanism (proteins with signal peptides), following the guidelines previously described for L. plantarum WCFS1 [26] and NCU116 [12]. The WCFS1 genome was used as reference ( Figure 9). We studied 304 proteins which were classified into eight groups, including: C-terminal anchored proteins (with and without cleavage sites), N-terminal anchored proteins (with and without cleavage sites), lipid-anchored proteins, proteins anchored to the wall by the LPxTG domain, secreted proteins (with cleavage sites) and proteins secreted by minor pathways (bacteriocins). Contrary to what was observed for the genes of protein secretion systems (Figure 8), notable differences were detected between the strains studied when analyzing the genes for extracellular proteins (Figure 9).
These results reveal that the surface of each of the strains of L. plantarum studied in this work owns a great variation in the surface/secreted proteins.

Study of Adhesion Factors
Adhesion factors are believed to help probiotic bacteria persist in the host's gastrointestinal tract and compete against pathogens as well as stimulate the immune system. Therefore, the property of adhesion to cells and molecules of the intestinal mucosa has been associated with health-promoting effects. We next investigated the presence of genes involved in the adhesion to the gastrointestinal tract in the genomes of L. plantarum MPL16, CRL1506, CRL681, TL2766 and WCFS1. Proteins with multiple copies of the MucBP or Mub_B2 mucin-binding domains have been identified in L. plantarum strains [12]. Genomic analysis revealed that the five strains evaluated in this work presented several proteins with multiple MucBP domains (PF06458) (Figure 10). All the strains presented the protein MucBP1 protein. L. plantarum CRL1506, CRL681, TL2766 and WCFS1 but not MPL16 had the genes for MucBP2 protein, whereas L. plantarum CRL681, TL2766, MPL16 and WCFS1 but not CRL1506 contained the MucBP3 protein. MucBP4 protein was detected in all strains except L. plantarum CRL681 (Figure 10). The MucBP5 protein, which in addition to contain MucBP domains has LPxTG cell wall anchor domains (IPR019931), was found in WCFS1, MPL16, CRL1506 and TL2766 strains. Proteins simultaneously containing MucBP and DUF285 domains (PF03382) were also found. One copy of the gene encoding the MucBP-DUF1 protein was detected in the genomes of TL2766 and WCFS1, whereas L. plantarum MPL16 had two genes for this putative adhesion molecule ( Figure 10). Analysis of their sequences showed that they are not copies of each other (data not shown). On the other hand, the five strains presented the MucBP-DUF2 protein. Two proteins with alternating mucin-binding domains (PF17965) and MubB2-like domains (PF17966) were also observed. One of these proteins also had MGB domains (PF17883) and was designated here as Muc-MubB2-MBG protein. This protein was found in the genomes of L. plantarum CRL1506, TL2766, MPL16 and WCFS1 but not in the CRL681 strain ( Figure 10). The other protein also possessed YGX-like MBG domains and was designated as Muc-MubB2-YGX protein. It was found in the genomes of L. plantarum CRL681, TL2766, MPL16 and WCFS1 and partially in the genome of the CRL1506 strain ( Figure 10). A protein with a MucBP domain and several MubB2 domains was also detected in all strains except in L. plantarum CRL681.
Interestingly, some of the proteins with known adhesion-mediating domains were found exclusively in one of the L. plantarum strains studied. A protein with alternating MubB2 and mucin-binding domains but without MBG or YGX domains (Muc-MubB2 protein) was detected only in L. plantarum CRL1506 (Figure 10). It was also observed that only the CRL1506 strain had a protein with a YceG-type domain (PF02618), which is an endolytic murein transglycosylase (also known as MltG, YrrL or YceG) and functions as a peptidoglycan terminase by endolytically cleaving strands of nascent peptidoglycan to finish their elongation. It has been proposed that this protein is involved in the adhesion of Lacticaseibacillus paracasei under conditions of heat stress [31]. Only L. plantarum CRL681 had a protein with two MubB2 domains and a protein with a SEC10/PgrA surface exclusion domain ( Figure 10). The SEC10/PgrA domain has been described in adhesins from Gram-positive bacteria such as SpyAD from group A streptococci and has been implicated in adhesion to the cell surface [32]. Proteins containing this domain have also been described in strains of Limosilactobacillus reuteri and it has been postulated that they mediate adhesion to the intestinal mucosa [33].
On the other hand, it was observed that only L. plantarum MPL16 has a protein with a SplA domain (PF03217) that also had a domain with mannosyl-endo-beta-N-acetylglucosaminidase activity (PF01832) (Figure 10). SplA domains have been associated with the ability of Lactobacillus acidophilus and Levilactobacillus brevis to adhere to intestinal epithelial cells (IECs) [34]. Additionally, only the MPL16 strain had a protein with a bacterial domain similar to immunoglobulin (Ig-like Big_3, PF07523).
In addition, the analysis of the secretome in L. plantarum WCFS1 identified other adhesion factors and three of them contained domains to adhere to collagen, chitin, and fibronectin [24]. These include lp_1229, which encodes the Msa protein. The deletion of the msa gene results in the loss of the ability of L. plantarum WCFS1 to agglutinate with yeast [35]. Recent studies comparing the mannose-specific binding ability of L. plantarum WCFS1 and 299v showed that L. plantarum 299v is superior to WCFS1 [36]. Variations in the promoter region of the msa gene produced higher levels of expression of the Msa adhesin in L. plantarum 299v than in WCFS1. On the other hand, it was reported that the adhesion capacity of the highly adhesive strains L. plantarum AR326 and AR269 depends on a protein that is 100% homologous with the glyceraldehyde-3-phosphate dehydrogenase (GAPDH) of L. plantarum WCFS1 [37]. The results indicated that GAPDH mediates the adhesion of these lactobacilli to IECs. Studies carried out with L. plantarum KLDS1.0391 showed that the luxS gene plays an important role in the bacteria's resistance to the gastrointestinal environment as well as in its adhesion capacity [38]. Therefore, the presence of the genes msa, gapdh and luxS as well as collagen binding proteins (cbp) and fibronectin (fbp) in the genomes of L. plantarum MPL16, CRL1506, CRL681, TL2766 and WCFS1 were also investigated ( Figure 11). Figure 11. Analysis of the presence/absence of adhesion factor genes in Lactiplantibacillus plantarum. The genomes of L. plantarum strains TL2766, MPL16, CRL681 and CRL1506 were studied and compared with the probiotic strain WCFS1 as reference.
Of the two fibronectin-binding proteins detected in the L. plantarum WCSF1 genome, only fbp2 was present in all the strains, whereas the gene for fbp1 was absent in the genome of the MPL16 strain. Two collagen-binding proteins were also observed. Of these, cbp1 was present in all the studied strains, whereas the cbp2 gene was not detected in L. plantarum MPL16 (Figure 11). In addition, the luxS, gapdh, aad (alpha-acetolactate decarboxylase) genes and a putative cell surface adhesin were detected in all L. plantarum strains (Figure 11). The gene for Cna protein type B (cnaB) was detected in L. plantarum CRL1506, CRL681 and TL2766 but not in MPL16 whereas msa was absent in the genome of L. plantarum CRL681.

Study of Lipoproteins and Teichoic Acids
Some studies have established that the composition of the lipoproteins expressed on bacterial surfaces conditions their interaction with the host's immune system, since they are considered important microbial-associated molecular patterns (MAMPs) capable of being recognized by receptors such as TLR2. In this sense, it has been shown that the lgt gene (lp_0755), which encodes prolipoprotein diacylglyceryl transferase, an enzyme responsible for the lipidation of lipoprotein precursors, is important for the immunomodulatory properties of the probiotic strain L. plantarum WCFS1 [39]. As can be seen in Figure  12, the lgt gene was detected in all the L. plantarum strains evaluated in this work. Species belonging to the Lactobacillus group have been reported to be producers of diacylated lipoproteins [40]. However, recent studies with L. plantarum WCFS1 have suggested that this species could also contain triacylated lipoproteins [39]. It has been reported that the genome of L. plantarum WCFS1 has three acyltransferases (lp_0856, lp_0925 and lp_1181) that could be involved in lipoprotein triacylation. Bioinformatic analysis revealed that the acyltransferases lp_0856 and lp_0925 are present in the genomes of MPL16, CRL1506, CRL681, TL2766 and WCFS1 strains ( Figure 12). The acyltransferase lp_1181 is also present in all strains except in L. plantarum MPL16. Teichoic acids (TA), and in particular, lipoteichoic acids (LTA), have also been shown to play an important role in the immunomodulatory activity of L. plantarum WCFS1 both in vitro [41,42] and in vivo [43]. These studies demonstrated that the deletion of the dltX-D (or dltD) gene, involved in the incorporation of D-alanine into LTA, alters the immunomodulatory properties of L. plantarum WCFS1. Our own studies using a dltD mutant of L. plantarum CRL1506 revealed the importance of LTA in the immunomodulatory activity of the strain in the context of intestinal inflammation mediated by TLR3 activation [44]. As expected, the dltD gene was found in all five L. plantarum genomes studied (Figure 12).

Discussion
In this study, the complete genomes of L. plantarum MPL16, CRL1506, CRL681, TL2766 and WCFS1 were compared with the aim of identifying the genes involved in their differential immunomodulatory properties. The comparative analysis of the genes involved in the metabolic pathways of the five strains under study did not reveal differences in the metabolism of amino acids, lipids, nucleotides, cofactors and vitamins or in the genes associated with energy metabolism. However, differences were found between MPL16, CRL1506, CRL681, WCFS1 and TL2766 strains when considering carbohydrate metabolism pathways, particularly in the presence/absence of GH and GT. These findings are in line with previous studies carried out for other L. plantarum strains. Molenaar et al. [45] investigated 20 L. plantarum strains and reported that genes involved in protein and lipid synthesis or degradation were largely conserved, whereas genes involved in sugar transport and catabolism were highly variable among the strains. In another study, 24 different phenotypes were found when 185 strains of L. plantarum were studied to determine their fermentation and growth characteristics [3], indicating differences among the strains in their ability to metabolize sugars. These studies were complemented and extended by subsequent genomic studies in which 108 [10], 127 [1] and 133 [11] different strains of L. plantarum were used. They confirmed the wide variety of phenotypes that can be found in this species when carbohydrate metabolism is considered. Genomic studies focused on the carbohydrate-active enzymes of L. plantarum found no significant differences in the number of genes in the families of carbohydrate esterases, carbohydrate-binding enzymes, auxiliary active enzymes, and polysaccharide lyases for 114 L. plantarum isolated from different niches [11]. The same study found variations in the numbers of GH and GT genes in L. plantarum strains. The GH family includes the most abundant enzymes, which are necessary for the degradation of glycosidic bonds between carbohydrates. The GH encoded in all L. plantarum strains were members of the GH1, GH109, GH13, GH31 and GH25 families. On the other hand, within the GT family, GT2 (necessary for cellulose synthesis) and GT4 (necessary for sucrose synthesis) were identified in all L. plantarum strains.
One study showed a lack of association between gene content and habitat when various strains of L. plantarum were studied [46]. However, more recent works identified that certain origins were weakly associated with two phylogenetic groups based on SNP analysis: groups that were designated as G2 and G3 [10]. G2 contained L. plantarum strains isolated mainly from meat samples whereas G3 contained strains from plants. In groups G1, G4 and G5, strains from plants, milk and dairy products as well as intestine isolates were located. For the five defined phylogenetic groups, no critical differences in genome size, number of coding sequences, or G + C content was found. However, the study found variations in the presence of some genes in the genomes of the five groups [10]. Strains in G1 and G2 had a higher capacity to metabolize various carbohydrates such as glucose, fructose, galactose, and lactose, but genes related to these metabolic processes were rare in groups G4 and G5. Thus, the differences between L. plantarum MPL16, CRL1506, CRL681, TL2766 and WCFS1 in relation to carbohydrate metabolism could be associated with the different origin of each of the strains (Table 1). It should be noted that the differences in the metabolic capacities of the MPL16, CRL1506, CRL681 and TL2766 strains, and therefore their ability to persist more efficiently in the gastrointestinal tract, could not be associated with their different immunomodulatory capacities.
In an attempt to explain the differences observed between the L. plantarum strains in this study, a group of genes proposed as "probiotic markers" in a recent study published by Carpi et al. [1] was first analyzed. The study performed a comprehensive pan-genomic analysis of L. plantarum considering only closed complete genomes available in open access repositories. The ultimate goal was to generate information that would serve as a reference to help in the characterization and identification of genes involved in the beneficial properties of potential probiotic strains. The comparative genomics study using the complete genomes of 127 L. plantarum strains identified 1,436 genes in the core genome, 414 in the soft core, 1858 genes in the shell, and 13,203 genes in the cloud genome, respectively, out of a total of 16,911 genes. Interestingly, this work detected 75 "probiotic marker" genes in the genomes of L. plantarum [1] based on previous works that proposed that genes associated with stress resistance, active metabolism in the host, adhesion and intestinal persistence are involved in the beneficial properties of lactobacilli [13,16,17]. Of the 75 probiotic marker genes, 70% corresponded to genes located in the core and softcore genome. Most of these genes were present in all the L. plantarum genomes used for the study [1]. On the other hand, 12 genes (clpP1, bshA, oppA3, oppA4, srtA, xylA, gbuB, gla2, glf, dps, glpF1 and cbh/bsh) were located in the shell and cloud genome, indicating that they were only present in some strains. In fact, genes such as bshA and oppA4 were found in only one strain of L. plantarum each. Interestingly, the study did not correlate a gene or group of genes with documented probiotic properties for the L. plantarum strains. For example, the list of strains that possessed any of these 12 genes did not include L. plantarum WCFS1 and strains with documented probiotic activity such as L. plantarum Zhang-LL [47] or NCU116 [48][49][50][51] possessed 2 and 3 of those 12 genes, respectively. As a first attempt to find genes that could explain the differential immunomodulatory activity of the MPL16, CRL1506, CRL681, TL2766 and WCFS1 strains, the study of the presence/absence of the probiotic marker genes proposed for L. plantarum was carried out [1]. Although the results revealed a marked difference between L. plantarum MPL16 and the other strains, it was not possible to clearly discriminate the immunomodulatory strains from the non-immunomodulatory ones. Considering these results, a comparative study focused on the molecules of the bacterial surface that could be responsible for interacting with the host's immunological receptors: EPS, adhesion factors, and LTA was carried out.
Studies performed by our research group demonstrated that the EPS produced by some immunobiotic bacteria have important roles in the immunomodulatory effect in the context of TLR3 activation in PIE cells. We showed that EPS produced by Lactobacillus delbrueckii TUA4408L [52], L. delbrueckii OLL1073R-1 [53], and Streptococcus thermophilus ST538 [44] participate in the ability of immunobiotic strains to increase the production of type I interferons (IFNs) and antiviral factors in PIE cells after poly(I:C) challenge. It was reported that L. delbrueckii TUA4408L and its EPS are capable of increasing the activation of the transcription factors IRF-3 and NF-κB, leading to an increase in the expression of IFN-β, MxA and RNAseL and a decrease in the replication of the rotavirus in PIE cells [52]. It was also described that the EPS of the strain OLL1073R-1 increases the expression of IFN-α, IFN-β, MxA and RNAseL in PIE cells challenged with poly(I:C) and that this effect depends on its ability to regulate the expression of the A20 protein [53]. The capacity of EPS-producing S. thermophilus ST538 to modulate the TLR3-mediated innate antiviral response in PIE cells was also demonstrated [44], by a successful development of two mutant strains lacking the epsB or epsC genes, which were unable to produce the macromolecule. Studies in PIE cells demonstrated that EPS from S. thermophilus ST538 was able to significantly enhance IFN-β, IL-6 and CXCL10 expression in response to TLR3 stimulation, whereas in S. thermophilus ΔepsB and ΔepsC mutant strains this effect was significantly decreased [44]. On the other hand, some studies have shown that EPS produced by L. plantarum strains are involved in their immunomodulatory properties. L. plantarum NCU116, isolated from a traditional Chinese sauerkraut called "paocai", has been shown to possess probiotic functions including immunomodulation ability [54,55]. This strain produces a heteropolysaccharide called EPS116, which has been shown to attenuate intestinal inflammation through the STAT3 signaling pathway and inhibit cancer cell proliferation through the TLR2-dependent mechanism [56]. EPS isolated from L. plantarum JLK0142 and designated EPS0142 has also been described as possessing immunomodulatory properties [57]. In vivo studies in mice immunocompromised by chemotherapy showed that oral administration of EPS0142 increases intestinal IgA levels, serum TNF-α and IL-2 levels, and the proliferative capacity of spleen lymphocytes [57].
Considering this background, studies were carried out to compare the EPS clusters in the genomes of the MPL16, CRL1506, CRL681, TL2766 and WCFS1 strains. Interestingly, L. plantarum MPL16, which is the strain with the most remarkable immunomodulatory activity among the lactobacilli evaluated here, was the strain that did not present any complete EPS cluster of the 4 described for the species. In line with these findings, it has been shown that the MPL16 strain does not have the ability to produce EPS [58]. No differences could be recognized between the EPS clusters that could explain the different immunomodulatory properties of the CRL1506, CRL681 and TL2766 strains.
Adhesion factors are considered key genes in the probiotic activity of many lactobacilli [12,16]. Particularly, in L. plantarum WCFS1, various proteins with adhesion function have been described, including mucin binding proteins, collagen binding proteins, chitin binding proteins, fibronectin binding proteins [24], agglutination and adhesion protein Msa [35,36] as well as GAPDH [57], which are essential for the strain to optimally fulfill its probiotic functions. L. plantarum NCU116 is another strain with probiotic properties for which exhaustive studies of its adhesion factors have been carried out. In its genome, 17 proteins that can be directly associated with adherence were identified through analysis of the composition of their domains. Five of these proteins possess a bacterial immunoglobulin-like domain (Ig-like) including PF17961, PF02368 and PF07523 according to the Pfam database [12]. These proteins are involved in extracellular binding and adhesion processes [59]. Two collagen-binding domain proteins were also found in strain NCU116, which can mediate adherence to host IECs. In addition, 8 proteins with multiple copies of the MucBP (PF06458) or Mub_B2 (PF17966) mucin-binding domains were identified. They play an important role in the adhesion of lactobacilli to the gastrointestinal tract [12]. The presence of a 1343 amino acid protein in L. plantarum NCU116 with MucBP and Mub B2 domains (PF17966) was also reported, which was identical (98.91% homology) to the adhesion protein Msa of L. plantarum WCFS1, which has been shown to participate in bacterial adherence [13,35].
When proteins with adhesion functions in L. plantarum MPL16, CRL1506, CRL681, TL2766 and WCFS1 were evaluated, a great variability between strains was observed (Figure 13). L. plantarum MPL16 differed from the other strains under study due to the absence of the fbp1, cbp2 and cnaB genes; and for being the only one to present genes for a protein with a SplA domain (PF03217) and for a protein with a bacterial domain similar to immunoglobulin (Ig-like Big_3, PF07523). Figure 13. Schematic representation of glycosyltransferases, glycosylhydrolases, proteins involved in resistance to passage through the gastrointestinal tract and bacterial surface molecules detected in the Lactiplantibacillus plantarum genomes. The genomes of L. plantarum strains TL2766, MPL16, CRL681 and CRL1506 were studied and compared with the probiotic strain WCFS1 as reference.
As mentioned above, the SplA domains have been associated with the ability of L. acidophilus and L. brevis to adhere to IECs [34]. On the other hand, bacterial Ig-like domains, which are frequently observed in cell surface proteins, have been shown to mediate cell-to-cell recognition via surface receptors [60,61]. L. acidophilus NCFM possesses a bacterial immunoglobulin-like domain-containing protein SlaP, which has been designated IgdA [60]. This protein is considered to be conserved among lactobacilli that form the S layer [62]. In silico analyzes revealed that IgdA and the corresponding orthologs were unique to host-adapted lactobacilli, whereas the immunoglobulin domain was specific to vertebrate-adapted species [60]. It was further observed that IgdA mutants of L. acidophilus possessed a visibly altered cell surface, which contributed to their increased salt sensitivity, altered immunogenicity profile, and severely reduced adhesion to intestinal Caco-2 cells, extracellular matrices, and mucin in vitro [60]. L. plantarum MPL16 was the only strain in this study in which the SlaP protein and a protein with bacterial Ig-like domains were detected; however, the data obtained in this work do not allow us to predict whether these proteins would be involved in its differential immunomodulatory activity. The development of mutants of the MPL16 strain that do not express these genes could contribute to revealing their specific role in the capacity of the strain to modulate the intestinal and respiratory antiviral immune response.
Bacterial lipoproteins are widely recognized MAMPs, which interact with TLR2. Lipoproteins are conjugated with two or three acyl (di-or tri-acyl) chains, which are essential for their proper anchoring in the cell membrane, as well as for interaction with TLR2. After export across the cell membrane, lipoprotein precursors (pre-prolipoproteins) undergo lipid modification, which is catalyzed by three conserved enzymes. Prolipoprotein diacylglyceryl transferase (Lgt) transfers a diacylglyceryl moiety to the cysteine residue, which is then directly cleaved by lipoprotein signal peptidase (Lsp) from the N-terminal signal sequence of the lipid-modified cysteine residue. Finally, lipoprotein N-acyl transferase (Lnt) adds a third acyl chain to the free amino group of the lipidated cysteine. Di-and triacyl lipoproteins produced by bacteria are differentially recognized by different TLR2 heterodimers: TLR2/6 and TLR1/2, respectively [39].
Information on the role of lipoproteins in the interaction of probiotic bacteria with the host is scarce. By deleting the lgt gene (lp_0755), which encodes prolipoprotein diacylglyceryl transferase, an enzyme responsible for the lipidation of lipoprotein precursors, the functions of the set of lipoproteins in the physiology of the probiotic strain L. plantarum WCFS1 was investigated [39]. These authors compared the ability of Δlgt mutant and wild type strain to stimulate TLR2 signaling and regulate the inflammatory response. The deletion of the lgt gene in the WCFS1 strain significantly reduced its ability to stimulate TLR1/2, whereas its ability to stimulate TLR2/6 signaling was not affected. Using human blood mononuclear cells, it was also shown that the Δlgt mutant induced a proinflammatory profile compared with the wild-type strain, characterized by increased production of IL-12, TNF-α, IL-1β and IL-8 and lower levels of the anti-inflammatory cytokine IL-10 [39]. It was then suggested that the loss of lipoproteins on the surface of the probiotic bacteria WCFS1, which were found in higher concentration in the culture medium, reduced its anti-inflammatory capacity. These results are in contrast to what occurs in pathogens such as Streptococcus agalactiae [63] and Staphylococcus aureus [64], in which the deletion of the lgt gene significantly increased its proinflammatory capacity.
Results of TLR1/2 and TLR2/6 activation showed that triacylated lipoproteins are dominant in L. plantarum WCFS1 [39]. This contrasts with previous reports describing the absence of orthologous genes for lnt, which is responsible for transferring a third acyl chain to the N-terminal cysteine of lipoproteins. Therefore, Gram-positive bacteria with low GC content such as species belonging to the Lactobacillus group have been considered as producers of diacylated lipoproteins [40]. The L. plantarum WCFS1 genome has been found to possess three acyltransferases (lp_0856, lp_0925, and lp_1181) that could be involved in a similar role to lnt in lipoprotein triacylation [39]. In this work, we revealed that the three acyltransferases are present in the genomes of CRL1506, CRL681 and TL2766 strains but that L. plantarum MPL16 does not have the acyltransferase lp_1181. Whether this difference impacts the composition of the surface-expressed lipoproteins of the MPL16 strain and consequently on its interaction with the host's immune system is an interesting point for future research.
LTA are the most proinflammatory components of the Gram-positive bacterial envelope and this effect has been shown to be highly dependent on D-alanine substitution of the TA backbone [65]. Although the potency of immunomodulatory capacity differs among bacterial strains [66], it has been shown that LTA purified from L. plantarum WCFS1 can induce a potent pro-inflammatory response in immune cells in vitro [41,42], and that this response was dependent on D-alanine substitution in the LTA backbone. Indeed, the absence of D-alanine in LTA changed the ability of L. plantarum WCFS1 and its purified LTA to modulate immune responses in vitro towards a more anti-inflammatory cytokine profile [42]. In vivo studies compared the effects of L. plantarum WCFS1 and its negative mutant for D-Ala in LTA (ΔdltX-D) [43] on the distribution of populations of intestinal and systemic T cells and dendritic cells (DCs) in mice. It was shown that the distribution of not only pro-, but also anti-inflammatory T cell and DCs populations depends on the functionality of the dltX-D gene and LTA D-alanine in the cell envelope of L. plantarum WCFS1 [65]. Not only the proinflammatory immune responses were reduced in the absence of D-Ala substitution, but also anti-inflammatory responses, such as the generation of regulatory T cells. In line with these studies, our group demonstrated that a mutant strain of L. plantarum CRL1506 (Δdlt) has altered immunomodulatory capacity in the context of intestinal inflammation mediated by TLR3 activation compared with the wild-type strain [44]. L. plantarum CRL1506 (Δdlt), similar to the wild strain, was able to increase IFN-β expression levels in PIE cells stimulated with poly(I:C), but did not induce changes in IL-6 levels or MCP-1. The CRL1506 (Δdlt) strain was also unable to increase IL-10 levels or reduce CD3 + NK1.1 + CD8αα + intraepithelial lymphocytes or TNF-α, IL-6, or IL-15 levels in the intestine of mice challenged with the TLR3 agonist. As expected, the dltD gene was found in the five genomes of L. plantarum studied (MPL16, CRL1506, CRL681 and TL2766), and therefore could not be associated with their different immunomodulatory capacities.
Interestingly, L. plantarum MPL16 was the strain that presented the least number of genes associated with resistance to passage through the gastrointestinal tract. Compared with CRL1506, CRL681, and TL2766 strains, the genome of L. plantarum MPL16 did not have the cbh, treA, celB, eno2, oppA1, oppA2, and oppA3 genes and was the only one to have the gbuB gene ( Figure 13). These results allow us to speculate that the MPL16 strain would be less resistant to passage through the gastrointestinal tract. These findings also allow us to speculate that the action of gastric juice and digestive enzymes could more easily expose surface molecules of L. plantarum MPL16 that could interact with the immunological receptors of the intestinal mucosa.
In summary, the comparison of genes potentially involved in the in vitro and in vivo differential effects of the L. plantarum strains in this work revealed no obvious differences between the bacteria evaluated. The comparative and functional genomic studies carried out do not allow proposing bacterial genes involved in the immunomodulatory effect, since it was not possible to correlate a group of specific genes with the differential immunomodulatory capacity of the L. plantarum strains. Proteomic studies are ongoing to evaluate the possible role of protein expression patterns with the immunomodulatory differences among strains.

Microorganisms
L. plantarum CRL1506, L. plantarum CRL681 and L. plantarum MPL16 belong to CER-ELA Culture Collection. L. plantarum TL2766 was kindly provided by the Culture Collection of the Nutritional Immunology Group of the International Center for Agricultural and Nutritional Immunology of Tohoku University (Sendai, Japan). The published genomic sequence of L. plantarum WCFS1 was obtained from the GenBank database (ID number AL935263.2) [67].

Assembly and Annotation of Nucleotide Sequences
The quality of each sequence was controlled with FASTQC [68]. Sequences with low quality or with contaminations were eliminated with PrinSeq lite v0.20.4 [69]. The resulting readings were assembled using Bruijn graphs through the SPAdes v3.11.1 tool [70]. The obtained assemblies were uploaded to public access repositories (DDBJ/ENA/GenBank).
Genes were predicted with the Prokaryotic Genome Annotation Pipeline (PGAP) v4.8 using the stand-alone configuration. A virtual machine service in the cloud, provided by Microsoft Azure ® , was used because the scorer requires a high computing resource. The configuration used was D8v3 consisting of 32 GB of RAM, 8 Intel Xeon ® E5-2673 v3 (Haswell) 2.4 GHz processors and 200 GB of storage.

Bioinformatic Tools for Comparative Genomics Studies
The functional characterization of individual genes from the genomes of L. plantarum CRL1506, CRL681, MPL16, TL2766 and WCFS1 was performed with the BlastKOALA tool [71].
Amino acid sequences encoded in L. plantarum genomes were compared for ortholog group inference with OrthoFinder v2.5.2. The Venn Diagram was made with the Venn Diagram tool (https://bioinformatics.psb.ugent.be/webtools/Venn/; accessed on 15 th November 2021). The dbCAN2 server, which allows automated notation of active carbohydrate enzymes, was used for the analysis of glycosylhydrolases and glycosyltransferases [72].
The GenBank database was used to obtain the amino acid sequences of the genes belonging to the potential probiotic markers, including: proteins involved in the biosynthesis of EPS, extracellular proteins, adhesion factors, lipoproteins and teichoic acids. These genes were searched in the genomes of L. plantarum with BLAST [73] in the standalone mode.
The phylogenetic tree was constructed from protein sequences possessing a MucBP domain. Using a multiple alignment program (MUSCLE) [74], the sequences were aligned and the tree was built from the maximum likelihood estimation statistical test [75], both tools are available in MEGAX [76].
Gene presence/absence graphs were performed with Python 3.7 using the Pandas and Seaborn libraries.

Conflicts of Interest:
The authors declare no conflict of interest.