Kentrophoros flavus YE was collected on different dates from May 2020 to March 2022 from the intertidal zone near the eastern gate of Yantai University, China (37°28′50″N, 121°27′41″E). Samples were scraped at the chemocline in fine, well-sorted sands with a reduced, sulfide layer at an approximate 15 cm beneath the sediment surface. Wet sands were transported to the laboratory and examined under a stereo microscope (Guiguang GL-6345BI, China). Living cells were observed using a bright field and differential interference contrast microscopy at 100-1000× magnifications (Zeiss Axio Scope A1, Germany). Counts and measurements of the living cells were conducted at a magnification of 1000×.

Approximately 10 cells were isolated from the raw cultures using a micropipette and washed three times with in situ water filtered by a millipore filter of 0.2 µm in pore diameter (Pall Laboratory, US) to remove the irrelevant bacteria in the environmental water. The ciliate cells were then starved for approximately 5 h in the sterilized water, allowing the hosts to digest the potential bacterial food in the cytoplasm. Each individual was transferred into 1.5 ml Eppendorf tubes with a small drop of water for DNA extraction. DNA was extracted using the REDExtract-N-Amp Tissue PCR Kit (Sigma, St. Louis, MO) according to Gong et al. (2016).

The host’s 18S rRNA gene was amplified via polymerase chain reaction (PCR) with the primers 82F (5′-GAAACTGCGAATGGCTC-3′) and 1492R (5′-GGTTACCTTGTTACGACTT-3′) according to Wang et al. (2019). The PCR reaction solution (25 µl) contained 2 µl of 10 µM primers, 1 µl of template, 12.5 µl of MegaFi™ Fidelity 2X PCR MasterMix with MegaFi™ Fidelity DNA Polymerase (ABM, Shanghai, China), and 9.5 µl of nuclease-free H₂O. The cycling parameters were as follows: pre-run of 3 min at 95°C, followed by 30 cycles of 30 s at 95°C, 30 s at 56°C, and 90 s at 72°C, with a final extension of 10 min at 72°C. The associating bacteria’s 16S rRNA gene was amplified via PCR with the primers 27F (5′-AGAGTTTGATCMTGGCTCAG-3′) and 1525R (5′-AAGGAGGTGATCCAGCC-3′). The PCR reaction solution (25 µl) contained 2 µl of 10 µM primers, 1 µl of template, 0.5 µl of dNTP Mix (10 mM of each), 1.25 units of Vazyme Lamp^® DNA Polymerase (Vazyme, Nanjing, China), 2.5 µl 10X of Vazyme Lamp^® buffer with MgCl₂, and 18.75 µl of nuclease-free H₂O. The cycling parameters were as follows: pre-run of 5 min at 94°C, followed by 30 cycles of 90 s at 94°C, 90 s at 53°C, and 3 min at 72°C, with a final extension of 10 min at 72°C. The PCR products were purified using a TIANgel Midi Purification Kit (Tiangen Bio. Co., Shanghai, China) and then inserted into pCE2 TA/Blunt-Zero vectors (5 min TA/Blunt-Zero Cloning Kit, Vazyme, Nanjing, China). The recombinant vectors were then transformed into competent cells Trans1-T1 (TransGen Biotech Co., Ltd., Beijing, China). The recombinants were screened using the blue/white selection. The white colonies were randomly selected and screened using PCR with the primers M13-20 and M13-26. Positive clones with the target inserts were sent for Sanger sequencing (Sangon Biotech Co., Ltd., Shanghai, China). A total of 132 clones were sent, and 46 sequences were successfully obtained, with an average length of 1500 bp.

For the other two ciliate individuals, the PacBio Sequel II platform was selected to complete the collection of the 16S rRNA full-length amplicons. The 16S rRNA gene was amplified using the same primers and PCR cycling parameters as those applied in the clone library construction. The PCR amplicons derived from two individuals were purified separately using the DNA Clean & Concentrator™-5 Kit (Zymo Research, Beijing, China), and each library was constructed according to the size of the amplified fragments (approximately 1500 bp). Sequencing was performed on a PacBio sequel II platform (Frasergen Co., Ltd., Wuhan, China), using a Hifi mode for small-insert (1.6kb) libraries, which allowed the insert being corrected by circular consensus sequences (CCS) for 30-40 times (with CCS≥3) and ensured an accuracy higher than 99.9%. More than 5000 raw CCS were generated from each sample.

Data derived from the clone libraries were examined manually to remove vector sections and low-qualified sequences using Bioedit (Hall et al., 2011). Data generated from PacBio sequencing were pretreated and qualified using the Pacific Biosciences’ SMRT Link v8.0 workflow (https://www.pacb.com/support/software-downloads). Subreads were generated using the lima program integrated in SMRT Link inferring the barcode information of each sample. The subreads were converted to raw CCS using the CCS program with a strict threshold (minPasses≥3, minPredictedAccuracy≥0.99). Cutadapt V3.2 (Martin, 2011) was used to identify and trim the adapters from the raw CCS with a maximum error rate of 0.1. Chimeras and raw CCS shorter than 1300 bp were also filtered by Cutadapt V3.2.

The qualified clone sequences and clean CCS were pretreated following the methods of Fu et al. (2020). Briefly, chimeras were further removed using the vsearch program integrated in QIIME v.1.8.0 (Caporaso et al., 2010), and the singletons were discarded using Mothur v.1.34.4 (Schloss et al., 2009). The resulting data were clustered according to operational taxonomic units (OTUs) at a sequence similarity of 97% using QIIME v.1.8.0 (Caporaso et al., 2010). The OTUs were blasted, aligned, and classified against the SILVA ribosomal RNA gene database (SILVA 128 database, accessed October 2021) (Quast et al., 2012). Analysis of the species composition and abundance was conducted to examine the major taxa of the associated bacteria using the Pacman package (Peng et al., 2021) appended in R (v3.6.2).

For the phylogenetic analyses of the 18S rRNA gene (ciliate host), 34 sequences from the ciliate class Karyorelictea and 11 sequences of the class Heterotrichea were downloaded from the GenBank database and were appended to the new 18S rRNA sequences. Members of class Heterotrichea were selected as outgroups. Two separate phylogenetic analyses of 16S rRNA were performed for the bacterial OTUs “Ca. Kentron” YE and MLAKF. For “Ca. Kentron” YE, the reference sequences were revised from the dataset of Seah et al. (2017). A total of 42 16S rRNA sequences from the class Gammaproteobacteria were appended to the new sequences, and the classes Coxiellaceae and Ectothiorhodospiraceae were used as outgroups (Seah et al., 2017). For the OTU MLAKF, 44 sequences from the classes Muribaculaceae, Bacteroidaceae, Prevotellaceae, and Porphyromonadaceae of the phylum Bacteroidetes were appended. The classes Bacteroidaceae, Prevotellaceae, and Porphyromonadaceae were selected as outgroups (Lagkouvardos et al., 2019).

The assembled 16S and 18S rRNA datasets were first aligned using SINA v.1.2.11 integrated in the SILVA Web server (https://www.arb-silva.de/aligner/) with default parameters (Pruesse et al., 2012). The resulting alignments were refined manually to excise poorly aligned sites in SeaView v4.0 (Galtier et al., 1996). For the dataset of K. flavus YE, “Ca. Kentron” YE, and MLAKF, the final matrix used for the tree analysis comprised 1524, 1580, and 1505 nucleotide positions, respectively. The maximum likelihood (ML) and MrBayes v.3.2.2 (BI) algorithms were used to perform the phylogenetic tree analyses. In all three datasets, the GTR + I + Г model was selected as the best model for both algorithms by jModelTest 2.1.7 (Posada, 2008) using the Akaike Information Criterion (AIC). The ML trees of the rRNA genes were constructed with RAxML V. 8 (Stamatakis, 2014) with 100 random taxon additions to generate the starting trees. Bootstrap support was inferred from 1000 pseudo replicate datasets. The Bayesian inference (BI) approach was implemented by MrBayes v.3.2.2 (Ronquist and Huelsenbeck, 2003). Markov chain Monte Carlo (MCMC) simulations were run for 5,000,000 generations (default temperature parameter and sampling frequency 0.01) with two parallel runs and each run with four simultaneous chains. The mean standard deviation of the split frequencies based on the last 75% of generations was 0.004 (< 0.01) at the end of the analysis. The first 12,500 trees (corresponding to 25% of generations) were removed as burn-in. MEGA v.6 (Tamura et al., 2013) was used to visualize tree topologies and adjust the branches.

Three 16S rRNA-targeted oligonucleotide probes targeting MLAKF were newly designed according to the method of Gong et al. (2014), and they were named C01-612, C01-824, and C01-1125. The specificities of the probes were confirmed using the Test Probe 3.0 tool against the SILVA NR 99 database (https://www.arb-silva.de/search/testprobe/). The non-sense probe NON338, complementary to C01-612/824/1125, was used as a negative control for the hybridization protocol. The specific probe chr4Ca for “Ca. Kentron” (Seah et al., 2017) was also employed in our FISH analysis to virtualize “Ca. Kentron” YE. The primary sequence and specificity of the four probes are shown in Table 1 .

The whole-cell hybridization procedure followed the protocol of Omar et al. (2017). Starved ciliates were transferred on microscopic slides (Frontier FRC-05) and dried at room temperature after fixation using Bouin’s solution (50%, final concentration). Prior to the FISH assay, the slides were washed three times in distilled water and then progressively dehydrated via ethanol gradients (30%, 50%, 80%, and 100%). The concentration of formamide for each probe was estimated in the mathFISH webserver (http://mathfish.cee.wisc.edu/formamide.html) (Yılmaz and Noguera, 2007). The highest concentration (30%) among the three was applied to maximally remove unspecific hybridization (Gong et al., 2014). There was 36 μl of the hybridization buffer (20 mM Tris–HCl at pH 8.0, 0.9 M NaCl, 0.01% SDS = sodium dodecyl sulfate, and 30% formamide) added to the slides together with 4 μl of the targeting fluorescent probes (5 ng/μl final concentration). “Cy3” fluorescence was applied to all MLAKF probes (C01-612/824/1125). For the double label experiment, the MLAKF and ”Ca. Kentron” YE probes were labeled with “Cy3” and “AF488”, respectively. Considering the disparity in the bacterial abundance, the two probes were mixed in various ratios (3.5 μl Cy3-C01-612/824:0.5 μl AF488-chrCa and 2 μl Cy3-C01-824: 2 μl AF488-chrCa) when applied in the experiment. The slides were incubated at 46°C for 3 h. After hybridization, the slides were dipped in a washing buffer (20 mM of Tris-HCl at pH 8.0, 450 mM of NaCl, 0.01% SDS) for 15 min at 48°C and then washed in chilled distilled water three times for 30 s. The specimens on the slides were covered with anti-fade mounting medium (Beyotime, China) and DAPI (50 ng/ml) for 5 min in darkness prior to observations under an epifluorescence microscope (Olympus BX61, Japan) using a SPOT RT3 digital camera (SPOT Imaging Solutions, Sterling Heights, US) and mercury lamp accessories (Olympus U-RFL-T, Japan). Photomicrographs were colored using Image-pro plus 5.1.

The size of the host ciliate was 200-460 × 16-36 μm in vivo, with a length:width ratio of 8-12:1 when fully extended (n=20; SD=71.5 and 7.8 for the length and width, respectively). The cells were highly flattened and had long ribbon shapes ( Figure 1A ). The body was curved ( Figure 1A ). The cells crawled along the bottom of the petri dish and attached to the substrate when stimulated. The ciliate had no cytostome. Nearly the entire dorsal surface of the cell was covered by rod-shaped bacteria that attached to the host erectly with one end ( Figures 1B–E ). Refractile globules (approximately 1 μm in diameter) were found within nearly all the rod-shaped bacteria ( Figures 1C, D ). The ventral side of the ciliates showed no symbiotic bacteria ( Figure 1E ).

The 18S rRNA gene sequence of Kentrophoros flavus YE was deposited in the GenBank database with the accession number OM421748. The acquired 18S rRNA gene fragment was 1566 bp long and had a GC content of 48.98%. The BLAST search on NCBI showed that this 18S rRNA gene fragment had the highest identity (99.68%) with that of the K. flavus isolate QD (the accession number FJ467505) (Gao et al., 2010). The ML and BI trees based on the 18S rRNA gene showed identical topologies ( Figure 2 ) in which the new sequence of K. flavus YE clustered with two species, the K. flavus QD and one isolate of K. gracilis, with high support from the ML (97%) and BI (1.00).

A total of 46 clean 16S rRNA gene sequences were obtained from the clone library technique and were assigned to 17 OTUs ( Figure 3 ). Proteobacteria accounted for a sizable proportion of the assemblages (80.43%, 37 sequences). The most abundant OTU was affiliated with the genus “Candidatus Kentron” (41.30%, 19 sequences). The other OTUs comprised only a few sequences that were closely related to the genera Thalassolituus (10.87%, 5 sequences), Neisseria (8.7%, 4 sequences), Escherichia (4.35%, 2 sequences), Sphingobium (4.35%, 2 sequences), Streptococcus (4.35%, 2 sequences), and “Candidatus Nitrososphaera” (4.35%, 2 sequences).

A total of 6154 and 4729 clean CCS reads were generated from the PacBio sequencing data for the two ciliate individuals, the total of which was clustered to 121 OTUs ( Figure 4 ). Based on the total reads, the phylum Proteobacteria accounted for a considerable proportion of the assemblages (72.12%), followed by Bacteroidetes (10.63%) and Firmicutes (9.78%). The top three abundant bacterial orders were “Ca. Kentron” (53.89%), Bacteroidales (10.63%), and Pseudomonas (6.63%). The measurable OTUs (>2% relative abundance) were affiliated with 7 genera: “Ca. Kentron” (53.89%), Muribaculaceae-like bacteria (the formerly S24–7 group of Bacteroidales) (10.34%), Psychrobacter (5.53%), Cutibacterium (4.73%), Lactobacillus (3.74%), Staphylococcus (3.07%), and Thalassolituus (2.01%). Of the different ciliate individuals, the diversity of the most least abundant OTUs (<10%) varied obviously, for example, Psychrobacter (5.53% vs. 0.01%), Cutibacterium (4.32% vs. 0.41%), Lactobacillus (2.26% vs. 1.48%), Staphylococcus (2.69% vs. 0.38%), and Thalassolituus (1.86% vs. 0.15%). Nevertheless, the OTUs affiliated with “Ca. Kentron” and Muribaculaceae-like bacteria consistently occurred with the highest proportions in both individuals ( Figure 4 ). Sequence of the former was nearly identical to the OTU affiliated with “Ca. Kentron” revealed in the Clone Library (differed by two bases). Accordingly, we named the two consistently occurring OTUs “Ca. Kentron” YE and “Muribaculaceae-like bacteria associated with K. flavus” (MLAKF). For the measurable OTUs in the PacBio data, the sequence length, GC content of the representative 16S rRNA, and the closely related sequences in the NCBI database gene fragments are provided in Table 2 .

The 16S rRNA gene sequences of “Ca. Kentron” YE (derived from the Clone Library and the PacBio platform) and MLAKF were deposited in the GenBank database with accession numbers OM501125, OM501126, and OM468468. The read sequences derived from the PacBio sequencing were also deposited in GenBank, with the bioproject number PRJNA801354.

Phylogenetic analyses based on the 16S rRNA genes were performed to assess the accurate affiliations of the OTU “Ca. Kentron” YE ( Figure 5 ) and the OTU MLAKF ( Figure 6 ). In the phylogenetic trees, “Ca. Kentron” YE fell into the clade of the genus “Ca. Kentron”, belonging to the class Gammaproteobacteria ( Figure 5 ). The BLAST search on the NCBI showed that this 16S rRNA gene fragment had the highest identity (95.43%) with that of “Ca. Kentron” H (accession number LT622015) (Seah et al., 2017). The OTU MLAKF fell into the family Muribaculaceae in phylum Bacteroidetes (100% ML/1.00 BI) ( Figure 6 ). Within Muribaculaceae, MLAKF appeared to be a distinct taxon, as a sister branch to a clade comprising Duncaniella species and several uncultured bacteria (98% ML/1.00 BI) ( Figure 6 ).

The newly designed probes C01-824 and C01-612 for OTU MLAKF have been shown in silico to be highly Muribaculaceae specific in the SILVA NR 99 database (accessed March 2022). C01-824 completely matched the family Muribaculaceae (126/126). C01-612 primarily matched the family Muribaculaceae (151/154), except for one sequence of the family Rikenellaceae (Bacteroidetes) and two sequences of Firmicutes (Bacteria). The third probe, C01-1125, was more MLAKF-specific but matched a smaller range of Muribaculaceae (62/72). It also matched 10 sequences of the family Prolixibacteraceae (Bacteroidetes) and one sequence of Thiomicrospiraceae (Gammaproteobacteria). All matched non-Muribaculaceae families above were not detected in our clone library ( Figure 3 ) and PacBio data ( Figure 4 ), implying their absence in the samples. Compared to the 16S rDNA of “Ca. Kentron” YE, the three probes C01-612, C01-824, and C01-1125 were mismatched by 6, 7, and 7 sites, respectively, excluding the possibility of nonspecific binding to “Ca. Kentron” YE.

The three probes were applied separately to a total of 10 cells of K. flavus (4 cells with C01-612, 3 with C01-824, and 3 with C01-1125). Strong positive signals were observed on all specimens, with rod-shaped structures distinguishable on several slides ( Figures 7C, G, K ; Figures 8C, H, M ). The overlays of the DAPI and probe’s signals showed that the MLAKF bacteria were primarily located on the surface of the ciliate cells ( Figures 7D, H, L ). To distinguish MLAKF and “Ca. Kentron” YE, double-label experiments were conducted on the same ciliate with specific probes labeled with different fluorophores ( Figure 8 ). For two sets of probe combinations, namely C01-612 + chr4Ca and C01-824 + chr4Ca, the red fluorescence (Cy3) of MLAKF ( Figures 8C, H, M ), and the green ones (AF488) of “Ca. Kentron” YE ( Figures 8D, I, N ) were simultaneously observed, both appearing on the host’s surface. Negative control was conducted by using the non-sense probe NON338 ( Supplementary Figure S1 ) in the FISH protocol, in which neither Cy3/AF488 fluorescence nor autofluorescence could be observed. Overlays of the signals on the same cell showed that the MLAKF and “Ca. Kentron” YE were mixed but were distinguishable from each other, with the former appearing more peripheral on the ciliate’s cell surface ( Figures 8E, J, O ).

The newly obtained 18S rRNA gene sequence for Kentrophoros flavus YE was closely clustered with Kentrophoros gracilis and Kentrophoros flavus (97% ML and 1.00 BI in both, Figure 2 ). In a comparison of the primary sequences, K. flavus YE differed from K. flavus and K. gracilis by 4 and 5 sites, respectively. In vivo, the Yantai populations matched the K. flavus population isolate in Qingdao (Xu et al., 2011) in terms of body shape and general appearance. In addition, the living cell characteristics of K. flavus YE were different from the original description of K. gracilis (Xu et al., 2011) in cell size (200–460 × 16–36 μm vs. 250–600 × 30–60 μm) and shape (K. gracilis was more elongated). The cells of K. flavus YE we observed were ciliated on the ventral side and contained ectosymbionts on the dorsal side ( Figure 1E ) corresponding to the “open” type among the four different body involution types recognized by Seah et al. (2017) for the genus Kentrophoros.

The OTU “Ca. Kentron” YE always displayed the highest proportion (41.3% in the clone library data and 53.89% in the PacBio data) among the K. flavus-associating bacteria ( Figures 3 , 4 ). This aligns with the previous study that used the meta-genomic technique to detect the genus “Ca. Kentron” in 13 Kentrophoros morphospecies (Seah et al., 2017). Compared with the most similar relative, “Ca. Kentron” H ( Table 2 ), “Ca. Kentron” YE was different by 70 bases in the primary 16S rRNA sequence, indicating that “Ca. Kentron” YE may represent a new ectosymbiotic species of the genus “Ca. Kentron”.

As seen from the PacBio results ( Figure 4 ), the OTU MLAKF accounted for the second-largest number (10.34%, approximately 1125 reads) among the associated bacteria in two K. flavus individuals. On 10 host cells that were sampled on different dates, the FISH results with three specific probes all produced positive signals of MLAKF on the surface of K. flavus ( Figures 7C, G, K ; Figures 8C, H, M ), suggesting that their presence on the ciliate was relatively stable.

The 16S rRNA gene sequence of MLAKF was 21 bases different from an uncultured bacterium with 98.58% identity ( Table 2 ), and 47 bases were different from the most similar species Duncaniella sp. C9 with 96.83% identity (accession number: MK521457). Combined with its distinct position in the phylogenetic trees ( Figure 6 ), this suggests that MLAKF should represent a new bacterial species in the family Muribaculaceae. The family Muribaculaceae comprises mainly molecular taxa, which are known as prevalent and abundant bacterial components in the mammal gut microbiome (Lagkouvardos et al., 2019). Members of Muribaculaceae have been widely detected in the intestines of mice, humans, ruminants, and homeothermic animals (Van Valkenburgh, 1991; Salzman et al., 2002; Lagkouvardos et al., 2016; Lagkouvardos et al., 2019; Miyake et al., 2020), with only a few species abundant in the natural environment, including manure, wastewater treatment plants, and coral reefs (Lagkouvardos et al., 2019). Our study adds to this list but with another scenario. Muribaculacid bacteria are abundant on the surface of the ciliated protist K. flavus and merge with thiotrophic ectosymbionts. This is the first time a muribaculacid bacterium has been found to be consistently associated with unicellular eukaryotic organisms.

It is worth noting that MLAKF was consistently present in large numbers on the cell surface of the host ( Figure 7 ), emerging with or on an outlayer of the ectosymbiotic “Ca. Kentron” YE ( Figure 8 ). The two symbionts and the ciliate host may have formed a three-partner consortium. First, MLAKF may offer extra organic carbon for the host-symbiont system. Genomic studies of “Ca. Kentron” have revealed that “Ca. Kentron” lacks genes encoding canonical enzymes for autotrophic CO₂ fixation, despite being a food source for Kentrophoros (Seah et al., 2019). It has been suggested that “Ca. Kentron” relies more on the uptake of exogenous organic substrates, heterotrophic carboxylation, and the recycling of host waste for growth (Seah et al., 2019). Extra energy and organic carbon sources should be required for the “Kentron-Kentrophoros” symbiote, and these should come from the environment or microbial activities (Seah et al., 2019). Coincidently, bacteria of Muribaculaceae have been predicted to play important roles in degrading dietary carbohydrates in the animal gut (Lagkouvardos et al., 2019). Genomic studies on a wide range of muribaculacid members have highlighted their prevalent versatility in degrading complex carbohydrates, benzoate resistance, and nitrogen utilization (>90% prevalence) (Ormerod et al., 2016; Lagkouvardos et al., 2019). Therefore, it is possible that MLAKF on the ciliate cells may play similar roles, degrading a wide range of substrate in the sediment (this might also include host waste) for growth, whereby the bio-product may provide accessible organic carbon for “Ca. Kentron”. Simultaneously, MLAKF themselves may offer an extra food source for the ciliate to support the overall growth of the host-symbiont system.

Second, it is known that Kentrophoros likely shuttles between oxic and anoxic zones in marine sediment with the ectosymbionts (Fenchel and Finlay, 1989; Seah et al., 2019). Under anoxic conditions, thiotrophic bacteria would not completely oxidize sulfide to sulfate, leading to less efficient energy generation (Seah et al., 2019). MLAKF may contribute as a compensated food source under anoxic conditions. Members of Muribaculaceae have been revealed to be well equipped with anti-oxygen strategies to support potential microaerobic growth (Ormerod et al., 2016). Nevertheless, the nature of the relationship between the MLAKF bacterium, “Ca. Kentron” YE, and K. flavus remains to be resolved, and single-cell genome sequencing can be a powerful tool for revealing the interactions between ciliated protists and their bacterial symbionts (Martinez-Garcia et al., 2012).

Aside from the thiotrophic bacteria “Ca. Kentron” YE, we detected a high diversity of associated bacteria using 16S rDNA-based approaches ( Figures 3 and 4 ; Table 2 ). This result accords with a previous morphological study of the species K. latum that recorded ectosymbiotic spirochaetes and intracellular prokaryotic symbionts (Raikov, 1971). The existence of particles adhering to thiotrophic bacteria has also been observed (Fenchel, 1969; Fenchel and Finlay, 1989) (e.g., diatoms and other ingested particles within Kentrophoros), which may be another explanation for the high bacterial diversity being detected. Except for the overwhelming OTU “Ca. Kentron” YE, the revealed diversity of the associating bacteria in K. flavus was variable between the clone library and PacBio sequencing ( Figures 3 , 4 ). This might be attributable to the significant difference in sequencing depths between the two techniques (46 clones vs. 10,883 CCS).

In addition to the consistently associated MLAKF bacterium, the comprehensive sequencing by the PacBio platform also detected OTUs related to members in the genera Psychrobacter, Cutibacterium, Lactobacillus, Staphylococcus, and Thalassolituus, all with measurable abundance (2–5.53%) ( Figure 4 ; Table 2 ). Most of these 16S rDNA fragments are most similar to the known species in GenBank ( Table 2 ), indicating that these might represent common taxa. Among these, members of Cutibacterium and Psychrobacter are obligate symbionts in human skin and the anogenital gland of captive pandas (Scholz and Kilian, 2016; Zhou et al., 2021), suggesting their symbiotic potential. However, the sequences might also derive from particles adhering to “Ca. Kentron”, as assumed by previous morphological studies (Fenchel, 1969). Because of the variable result among the different individuals and a lack of FISH analysis, we must refrain from predicting the stable existence or roles of these OTUs.