*2.5. Statistics and Clustering Analysis of (GTG)5 and ERIC Rep-PCR Patterns*

Patterns obtained from (GTG)5 and ERIC rep-PCR runs were imported into BioNumerics tool, package version 8. The cluster analysis was performed calculating Dice coefficients using a tolerance of 1% and optimization of 1%. Clustering based on Dice coefficients was performed using the unpaired group method with arithmetic average (UPGMA) method. For each fingerprint assay, the cluster cut-off analysis was applied to define the most reliable clusters. The cophenetic correlation index was used as a statistical tool to evaluate the quality of branches. For each of the genotyping experiments, the discriminatory index (D), expressed as the probability that two strains consecutively taken from a sample would be placed into different clades, was calculated as described by Hunter and Gaston [27].

Correlation analysis was performed considering the most abundant clades (Clade 5, 6 and 7), in order to obtain statistically significant results. The Spearman correlation index was calculated using the Hmisc v4.4–2 R package [28] implemented in R v 4.0.3 [29]. The correlation plot was obtained using the Corrplot v 0.84 package [30].

#### **3. Results**

#### *3.1. Kombucha Characteristics and Isolated Strains*

Within 3 days of cultivation in 600 mL beakers, the two kombucha tea samples, GTK and BTK, demonstrated a thin exopolysaccharide layer that became thicker with time (from 2–3 mm to 10 mm at the end of the cultivation period). Optical microscopy observations showed a high number of free bacterial and yeast cells, as well as aggregates of cells within the matrix, making cell counting uncertain and not enlightening even from the first day of fermentation.

Titratable acidity reached a maximum of 12 g/L at the end of fermentation in the black kombucha trial (GTK) and 6 g/L in the green one (BTK). pH dropped from approximately 3.7 to 2.75 for both samples, as a result of acid formation. Although no inhibition compounds were determined in this study, the difference in the final amount of acetic acid of the GTK and BTK samples could be due to the occurrence of more antibacterial compounds in GTK [31]. Ethanol was nearly 0 at the beginning of the cultivation time and reached 0.28% (*v*/*v*) and 0.30% (*v*/*v*) at the final time, for samples GTK and BTK, respectively. The low value of observed ethanol is in agreement with the gradual increase in acetic acid during cultivation time.

The isolation was performed from preculture (P1 and P2) and samples at 0, 6, and 12 days of cultivation. A total of 122 isolates were selected based on the colony and/or cellular morphologies. Microscopically, presumptive acetic acid bacteria cells appeared in single or grouped pairs or short chains, Gram-negative, catalase-positive, and KOH positive. A total of 46 isolates, which oxidized both ethanol and acetic acid, as confirmed by the ethanol chalk–ethanol test on GYC and the colorimetric assay on ACB, and produced a well-defined cellulose layer and/or formed both a cellulose layer and cellulose particles dispersed in the surrounding liquid, were considered in this study (Table 1).


In bold are reported the strains deposited in UNIMORE Microbial

> **Table 1.** Features of 46 isolates from Kombucha

 tea used in this study. Site of isolation and RFLP analysis.

#### *Appl. Sci.* **2021**, *11*, 1595




kombucha; GTK: green tea kombucha.

Qualitative cellulose tests confirmed their ability to produce BC. Significant macroscopic diversity in the cellulosic material was observed, coupled with high variability in terms of amount of BC, as reported in Figure 1. The strains with a BC production level of at least higher than 4 g/L were deposited in the UMCC culture collection with the code reported in Table 1.

**Figure 1.** Bacterial cellulose yields of selected strains isolated from Kombucha tea. Each value is the mean of three parallel replicates ± standard deviation.

#### *3.2. Diversity of Acetic Acid Bacteria from Kombucha Tea*

The RFLP analysis on 16S rRNA, using *Rsa*I and *Alu*I restriction endonucleases, was conducted as a preliminary test to investigate the diversity of isolates recovered from kombucha tea samples. The mapping analysis based on the size of the fragments obtained by *Rsa*I and *Alu*I showed high homogeneous results, indicating the isolates as members of the same species. This evidence was in agreement with a previous study [16] which grouped *Komagataeibacter* strains according to *Rsa*I and *Alu*I restriction enzymes and observed similarities between the pattern of *K. xylinus* DSMZ 2004 with the *Komagataeibacter* strains investigated in this study. The low variability of the fingerprinting pattern is also consistent with previous studies on the high similarity of the 16S rRNA gene among the *Komagataeibacter* genus and especially within the *K. xylinus* species [32,33].

From the genotyping of strains performed by using (GTG)5 and ERIC fingerprinting techniques, two dendrograms were generated (Figure 2). Based on the cluster cut-off analysis, profiles obtained using the amplification of interspersed tandem repeats GTG (Figure 2) were capable of discriminating between the strains considering a minimum percentage of similarity of 27.9%. A total of five major clades were created with a discrimination power of 0.42, calculated using the D index (Simpson index) [27]. In order to

improve the discrimination power of the analysis, we considered a minimum percentage of similarity of 70%. At this percentage of similarity, the resulted discrimination power was 0.94, defining 20 different biotypes. Among the detected clades, eight of them included just one isolate while seven were formed by two isolates. It was possible to distinguish two major clades, of which, one included nine isolates (clade 18) and another included five isolates (clade 19). The two reference strains, represented by K1G4 (=UMCC 2947) and K2G30 (=UMCC 2756) clustered into two different clades, represented by clades 19 and 7, respectively. The similarity percentage among the two references was 40.2%.

**Figure 2.** Dendrograms obtained from UPGMA analysis, using Dice's coefficient, of the digitized patterns obtained from (GTG)5 rep-PCR. The similarity threshold for biotypes discrimination was 70%. The cophenetic coefficient is represented by numbers and dots colored as red–orange–yellow–green, based on the branch quality.

The phylogenetic tree obtained from the digitized pattern profile using ERIC rep-PCR is represented in Figure 3. Based on the clustering cut-off analysis, ERIC was able

to discriminate using a minimum similarity percentage of 40%, identifying five major clades, as in the case of (GTG)5. The clustering cut-off analysis discriminated a total of 19 biotypes, with a discrimination power of 0.95, slightly higher compared to (GTG)5. A total of six clades were represented by just one isolate, while the remaining isolates were well distributed among the detected clades. Three clades were represented by five isolates (clade 7, 16, and 17) and two clusters by four isolates (clade 8 and 9). The remaining clades were represented by three and two isolates. The number of detected clades was not so different compared to (GTG)5 and, also in this case, the reference strains (K1G4 and K2G30) were found to be separated into two different groups (clade 15 for K1G4 and clade 12 for K2G30). The similarity of Dice indexes based on the patterns obtained from ERIC rep-PCR was 56.7%; quite similar to that obtained using (GTG)5 rep-PCR (40.2%).

**Figure 3.** Dendrograms obtained from UPGMA analysis, using Dice's coefficient, of the digitized patterns obtained from ERIC rep-PCR. The similarity threshold for biotypes discrimination was 70%. The cophenetic coefficient is represented by numbers and dots colored as red–orange–yellow–green, based on the branch quality.

#### *3.3. Improving the Fingerprinting of Komagataeibacter Strains*

In order to improve the sensitivity of the intraspecific discrimination, the data obtained from both rep-PCR assays were combined. The clustering analysis is reported in Figure 4. Based on the cluster cut-off analysis, the combination of (GTG)5 and ERIC rep-PCR allowed for discriminating the isolates, grouping them into 22 clades, considering a similarity percentage threshold of 94%. The discrimination power, calculated as the Simpson index, was 0.94, mining a high intraspecies diversity. Most of the isolates were grouped into three clusters (clade 5, 6, and 7) representing 43.5% of the isolates (*n* = 20). The remaining isolates (*n* = 26) were dispersed in 19 clades, of which, most of them consisted of just one isolate (*n* = 12). The reference strains were clustered into two different clades, as observed in both the case of (GTG)5 and ERIC rep-PCR. Clade 10, in which K1G4 was clustered with K1G6, and clade 19, where K2G30 clustered alone, had a similarity percentage of Dice index of 84.9%, higher considering the separated results from (GTG)5 and ERIC rep-PCR. This result is in agreement with the meaning of the combined analysis, since the two reference strains belong to the species of *K. xylinus*, as previously stated [14,34].

**Figure 4.** Dendrogram obtained combining the digitalized patterns from both (GTG)5 rep-PCR and ERIC rep-PCR. The dendrogram was drawn from UPGMA analysis using Dice's coefficient. The discrimination of biotypes was performed considering a similarity threshold of 94%. The cophenetic coefficient is represented by numbers and dots colored as red–orange–yellow–green, based on the branch quality. Cluster colors were defined in agreement with the Figure 1.

### *3.4. Correlation between Strains and Cellulose Yield*

In order to understand putative correlation of the BC yield obtained from isolates with time and environmental factors, such as sampling fraction (pellicle or liquid) and clustering methods applied in this study, a Spearman correlation index was calculated. At the beginning all clades that resulted from the combined fingerprinting approach (Figure 4) were compared with all variables (Figure 5a). All factors seemed to strongly influence the BC yields. Only the time of sampling seemed to slightly influence BC production, but with a low value (R2 = −0.38, *p* = 0.009). In order to understand how the sampling time could affect BC yield, the three most abundant clades, represented by clade 5, 6 and 7 were evaluated by calculating the Spearman coefficient, considering both time and sampling fraction. Regarding clade 5 (*n* = 6) (Figure 5b), the time seems to correlate positively (R2 = 0.84, *p* = 0.03). Indeed, the amount of BC (Figure 1) was higher for isolates sampled at day 12 (K2A46, K2A47, K2G46 and K1A34) and lower for isolates sampled at day 0 and from the preculture (K1G5 and K1G3 respectively). In clade 6 (*n* = 8) (Figure 5c), BC yield seems to be not statistically correlated by the sampling time (R2 = −0.70, *p* > 0.05). Regarding clade 7 (*n* = 6) (Figure 5d), the sampling time negatively affects the BC amount (R2 = −0.92, *p* = 0.009).

**Figure 5.** Correlation plot of the most abundant clades. Spearman correlation among time, environmental variables and obtained yields. (**a**) All retrieved clades were compared with time and environmental variables. Yield obtained within (**b**) clade 5, (**c**) clade 6 were compared with time and fraction variables. For clade 7 (**d**), includes only isolates sampled from liquid fraction.

#### **4. Discussion**

Bacterial cellulose is receiving a great deal of attention due to several features that make it suitable for many biotechnological applications, such as biomedical devices, food packaging components and engineering materials [1,2,35].

It is well known from the literature that kombucha tea represents one of the most abundant reservoirs of BC producer bacteria, in particular, strains belonging to the *Komagataeibacter* genus [10,12].

In this study, kombucha tea was chosen as a selective source for recovering cellulose producing bacteria. A total of 46 isolates were selected as the main components of the liquid and pellicle fractions of GTK and BTK kombucha tea samples and their phenotypic and genotypic variability was evaluated. RFLP of 16S rRNA gene and an analysis of the interspecific variability of isolates was conducted, and the results were discussed, considering the available genome sequences of strains K1G4 and K2G30.

The data obtained by RFLP analysis showed a very low variability in terms of the number of different species isolated from kombucha tea. This result is in agreement with previous studies reflecting both the low species diversity in kombucha tea [10], as well as the difficulty of taxonomic assignment of *Acetobacteraceae* by using 16S rRNA gene, due to its similarity in high phylogenetic related organisms, as extensively reported in the literature [36].

Considering the analysis of the intraspecific variability within bacteria, usually, genotyping is performed using the amplification of highly conserved tandem repeat regions that are widely dispersed in the genome. The amplification of repeated elements of genomic DNA, such as (GTG)5 and ERIC sequences, offers some advantages in terms of costs and rapidity, and is well developed for several Gram-negative and Gram-positive bacterial species, such as *Lactobacilli* [37] and *Enterococci* [38] and acetic acid bacteria [39].

In this study, the clustering analysis using the digitized patterns of (GTG)5 elements revealed the clustering of the 46 isolates into five large clusters considering a cluster cut-off value of 50% in similarity. The discrimination power of this analysis was 0.94 considering a similarity threshold of 70%. At this similarity value, a total of 20 different biotypes were detected. Similar results were obtained using the amplification of ERIC repeated elements. Considering the cluster cut-off analysis, the comparison of the digitized patterns of ERIC rep-PCR showed a clustering of the 46 isolates with a similarity percentage of 50%, the same result obtained using (GTG)5 rep-PCR. The discrimination power of ERIC fingerprint was 0.95, very similar to the (GTG)5 rep-PCR, considering a similarity percentage of 70. Based on these results, a total of 19 different clades were identified, compared to the 20 obtained using (GTG)5 rep-PCR. In both the fingerprint approaches, K1G4 and K2G30, previously identified as *K. xylinus* [14,34], were considered as reference strains and they were clustered into two different clades. Among the two reference strains, the percentage of similarity was 44.2% in the case of (GTG)5 rep-PCR and 56.7% for ERIC rep-PCR. These results support the hypothesis that most of the strains analyzed belong to the *K. xylinus* species and that the discrimination was performed at the strain level.

In order to improve the genotyping analysis, a combined model of the two fingerprint assays was performed (Figure 3). The resulting analysis showed a clusterization of the isolates into seven clades based on the cluster cut-off value of 88% of pattern similarity. The discrimination power of this analysis was 0.94 at 94% of similarity. The number of biotypes defined by the combined model was 22; meaning that an improvement of the discrimination of the isolates was achieved compared to the single results obtained by the two approaches. Additionally, concerning the reference strains (K1G4 and K2G30) in this case, we observed a grouping into two different clades, with a similarity percentage of 84.9%, meaning that all the isolates with equal or major than that value possibly belong to the species of *K. xylinus*. Based on this consideration, all of the isolates could be assigned to the *K. xylinus* except one isolate, the K2G44, clustering out of the group, with a similarity percentage of 78.8%. The data were in agreement with previous studies of the fingerprinting analysis of acetic acid bacteria [39,40], in which the differences of strains among the same species (as *Komagataeibacter xylinus*) were about 80%. In this case, combining the two methods, the detection of different strains was improved.

The correlation analysis deriving from the combined strain clustering allowed us to link the BC yield with isolation time and environmental variables. A possible explanation

of the obtained results may lie in the genetic plasticity features of acetic acid bacteria. It is well known in the literature that acetic acid bacteria are able to adapt themselves to environmental stressors, such as temperature and high concentration of organic acids, by acquiring genetic material from the environment or losing it. Indeed, acetic acid bacteria possess a large number of transposons genes and plasmids [41]. One example is provided by *Acetobacter* species which strains adapt themselves to environmental conditions along cultivation and fermentation cycles [42,43]. In *Komagataeibacter* species the genetic plasticity has not yet been evaluated. However, we believe this aspect should be elucidated in the light to better understand the variability of BC yield within strains of *Komagataeibacter* genus.

The kombucha tea samples analyzed in this study allowed us to recover a pool of acetic acid bacteria strains, which have been characterized and deposited into the UMCC culture collection as promising candidates for BC production.

**Author Contributions:** M.G. conceptualized the research, provided resources and supervised the work; D.M. performed experiments; S.L.C. performed genotyping analysis and wrote part of the original draft; L.D.V. reviewed and edited the manuscript; K.A. wrote part of the original draft; M.B. participated in the investigation, reviewed and edited the manuscript. All authors have read and agreed to the published version of the manuscript.

**Funding:** Part of this research was granted by The Research Doctorate School in Food and Agricultural Science, Technology and Biotechnology (STEBA).

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Acknowledgments:** The National Institute of Applied Sciences and Technology (INSAT, Tunisia), is acknowledged for providing Kombucha samples. A special thanks to The Research Doctorate School in Food and Agricultural Science, Technology and Biotechnology (STEBA) is acknowledged for assigning the "*Michela Stanca*" award. The JRU MIRRI-IT is greatly acknowledged for technical and scientific support.

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

#### **References**

