Bacterial Diversity and Community in Regional Water Microbiota between Different Towns in World’s Longevity Township Jiaoling, China

: Healthy longevity is associated with many factors, however, the potential correlation between longevity and microbiota remains elusive. To address this, we explored environmental microbiota from one of the world’s longevity townships in China. We used 16S rRNA gene high-throughput sequencing to analyze the composition and function of water microbiota. The composition and diversity of water microbiota signiﬁcantly differed between the towns. Lactobacillus , Streptococcus , Bacteroides , Faecalibacterium , and Stenotrophomonas were only dominant in Xinpu, a town with an exceptionally high centenarian population. Several biomarkers were identiﬁed, including Flavobacterium , Acinetobacter , Paracoccus , Lactobacillales , Psychrobacter , Bacteroides , Ruminococcaceae , and Faecalibacterium , and these shown to be responsible for the signiﬁcant differences between towns. The main species contributing to the differences between towns were Cyanobacteria , Cupriavidus and Ralstonia . Based on KEGG pathways showed that the predicted metabolic characteristics of the water microbiota in Xinpu towns were signiﬁcantly different to those of the other towns. The results revealed signiﬁcant differences in the composition and diversity of water microbiota in the longevity township. These ﬁndings provide a foundation for further research on the role of water microbiota in healthy longevity.


Introduction
Increasing human health and longevity is of global interest. The continuous increase in lifespan and the consequent growth of the elderly population have led to an increase in studies on aging. Centenarians are considered a symbol of health and longevity [1,2]. Therefore, many studies have focused on centenarians to identify the factors influencing health and longevity [3][4][5]. Although the mechanism of longevity is complex and there are many controversial debates, researchers have basically reached the consensus that longevity is comprehensively affected by genetics [6][7][8], lifestyle [9], environment [10,11], health care level [12], psychological factors [13,14], and other factors [15][16][17].
In the biomedical field, the research has focused on searching for molecular and biological factors that promote healthy aging and longevity. Considerable attention has been paid to analyzing the role of genetic factors in determining healthy aging and longevity. Siblings of centenarians have a lower risk of suffering from major age-related diseases (diabetes, cancer, and cardiovascular diseases) than controls from the same population [18,19]. Furthermore, the offspring of centenarians are healthier than age-matched controls belonging to the same demographic cohorts [20], suggesting that the biological components of

Sample Collection
In August 2018, 58 water samples were collected from eight towns (TW1, Guangfu, 6 samples; TW2, Nanzai, 4 samples; TW3, Wenfu, 9 samples; TW4, Changtan, 10 samples; TW5, Jiaocheng, 12 samples; TW6, Lanfang, 6 samples; TW7, Sanzhen, 6 samples; and TW8, Xinpu, 5 samples) belonging to the "world's longevity township" of Jiaoling County (Meizhou City, Guangdong Province, China). Water samples were collected from the lake at the front and back of a long-lived family's house or drinking water. In order to avoid potential contamination, we have done the following: (1) When sampling, we first put on clean gloves to remove floating objects on the water surface, and collected water samples with a sterile ribozyme-free sampling tube (Axygen). We unscrewed the lid at 10-15 cm Diversity 2021, 13, 361 3 of 18 underwater and tightened the lid when it was full. (2) Water samples were filtered, DNA extraction and PCR were all performed in a clean workbench after UV sterilization, so as to avoid potential pollution effects caused by environmental factors. After collection, the samples were kept at 4 • C during transportation to the laboratory and then stored at −80 • C until DNA extraction was performed. The distribution of population over 80 years old is provided by the local Jiaoling County government.

DNA Extraction and Polymerase Chain Reaction (PCR) Amplification
The bacterial DNA is extracted by using the Bacteria DNA Kit according to the manufacturer's protocols. First, the 200 mL water sample was filtered with a 0.22 µm sterile filter membrane, and then the side of the filter membrane that had adsorbed microorganisms was placed in a centrifuge tube with grinding beads so that the grinding beads were closely attached to the side with microorganisms. Then, the sample was then vortexed and shaken for 15 min, centrifuged at 13,000× g for 10 min, and 600 µL of the supernatant was transferred to a new 2.0 mL tube. Then, 150 µL of PS buffer and absorption solution were added, respectively, and centrifuged at 13,000× g for 5 min, put the supernatant in a new tube, and added 700 µL of GDP buffer. The HiPure DNA mini-column absorbed the product and eluted it with sterile water. We used Qubit method and 1% agarose gel electrophoresis to detect the concentration and purity of genomic DNA to quantify and report the biomass, and only reached level B (Qubit: 5 µg ≤ m < 10 ug or m ≥ 10 ug; C ≥ 20 ng/µL; AGE: The sample with no degradation, no stickiness, and no impurities) before proceeding to the next experiment. Sample DNA was diluted to 1 ng/L in sterile water. PCR primers were designed to target the V4 and V5 hypervariable regions of the bacterial 16S rRNA gene, which was amplified by PCR (95 • C for 2 min; 25 cycles at 95 • C for 30 s, 55 • C for 30 s, and 72 • C for 30 s; and extension at 72 • C for 5 min) using primers 338F 5 -ACTCCTACGGGAGGCAGCA-3 and 806R 5 -GGACTACHVGGGTWTCTAAT-3 . Phusion high-fidelity PCR Master Mix with GC Buffer (New England Biolabs) was used to ensure the amplification efficiency and accuracy, with the diluted genomic DNA as the template.

16S rRNA Gene Sequencing
The amplicon library was built using the Ion Plus Fragment Library Kit 48 RXNS (Thermofisher, Shanghai, China). The library was sequenced using Ion S5TMXL (Thermofisher, Shanghai, China) at Novogene Bioinformatics Technology Co., Ltd. A single-end sequencing method was used to construct a small-fragment library.

Statistical Analysis
Uparse software (Uparse v7.0.1001, http://www.drive5.com/uparse/ (accessed on 10 September 2018)) [40] was used to cluster all clean reads from the samples. Mothur and the SSUrRNA database [41] from SILVA132 (http://www.arb-silva.de/ (accessed on 10 September 2018)) [42] were used for species annotation. MUSCLE [43] (v3.8.31, http://www.drive5.com/muscle/ (accessed on 10 September 2018)) was used for rapid multiple sequence alignment to determine the phylogenetic relationship of all operational taxonomic units (OTUs). QIIme (v1.9.1) was used to calculate the Shannon index and the Unifrac distance as well as to construct the sample-clustering tree using the unweighted pair group method with arithmetic mean (UPGMA). R (v2.15.3) was used for dilution curve, rank abundance curve, and species accumulation curve calculations, and it was also used to analyze the differences in alpha and beta diversity indices between the groups using the Wilcoxon test in the agricolae package for parametric test. Linear discriminant analysis effect size (LEfSe) software was used to find biomarkers using a default score setting of linear discriminant analysis (LDA) of 4. Based on species abundance, the correlation coefficient (Spearman correlation coefficient or Pearson correlation coefficient) was calculated between the genera, to obtain a correlation coefficient matrix. For edges, graphviz (v2.38.0) was used to draw the network graph. The Tax4Fun function prediction was executed using the nearest-neighbor method based on the minimum 16S rRNA sequence similarity. The specific method involved extracting the 16S rRNA gene sequences from the whole prokaryotic genome of the Kyoto Encyclopedia of Genes and Genomes (KEGG) database and comparing them to those available in the SILVA SSU Ref NR database (BLAST bitscore > 1500) using the BLASTN algorithm to establish a correlation matrix.

Characteristics of the Longevous Populations and Sample Sequencing
Jiaoling County is affiliated to Meizhou City, which is located in the northeast of Guangdong Province. In the eight towns of the county, the proportion of the population over the age of 80 (in a total population over the age of 80 approximately 8365 people) was 6.91%, 10.57%, 9.07%, 9.12%, 22.30%, 11.73%, 8.67%, and 21.64%, respectively. The largest and smallest elderly populations were found in Jiaocheng (TW5 samples) and Guangfu (TW1 samples), respectively ( Figure 1). According to the proportion of the population aged 80-89, 90-99, and >100 to the total population over 80, Xinpu (TW8 samples) had the highest proportion of centenarians, followed by Changtan (TW4 samples), while Nanzhai (TW2 samples) had the lowest (Figure 1). The proportion of centenarians tended to increase from north to south along the Shiku River.
Bacterial genomic DNA was isolated from water samples of the eight towns and the 16S rRNA gene was sequenced using the Ion S5TMXL sequencing platform (by single-end sequencing), generating small-fragment libraries. By shearing and filtering the reads, an average of 85,015 reads per sample was obtained, with an average of 80,136 valid reads after quality control. The efficiency of quality control was 94.31%. The sequences were then clustered into OTUs at 97% identity, yielding 4601 OTUs that were finally annotated using the Silva132 database. Specifically, 100% OTUs were annotated to the bacterial domain level; 95.11% were annotated to the phylum level; 88.46% were annotated to the class level; 78.55% to the order level; 68.14% to the family level; and 46.84% to the genus level.

Composition and Structure of the Water Microbiota
Based on the OTU analysis and clustering, shared and unique OTUs were identified in different groups. Overall, 137 OTUs were shared in the TW1-TW8 samples, and the number of unique OTUs in each group was 399, 98, 246, 139, 900, 138, 217, and 143, respectively (Supplementary Materials Figure S1A). Based on the species annotation, the top ten most abundant bacteria at each classification level (phylum, class, order, family, and genus) were identified in each sample group, and a cumulative bar map of species relative abundance was generated to allow visual assessment of the relative bacterial abundance and their proportions at the different classification levels in each sample. The relative abundance at the phylum level is presented in a histogram ( Figure 2A). The major abundant phyla were Proteobacteria, Cyanobacteria, Bacteroidetes, Firmicutes, Actinobacteria, Acidobacteria, Chloroflexi, Verrucomicrobia, Fusobacteria, and Planctomycetes. Among these, bacteria from Proteobacteria, Bacteroidetes, and Firmicutes were dominant at all taxonomic levels.
Based on the species annotation and abundance at the family level, the 35 most abundant families were identified and clustered at the species and sample levels to generate a heat map (Supplementary Materials Figure S1B). The top ten families in the samples were Burkholderiaceae, unidentified Cyanobacteria, Rhizobiaceae, Chromobacteriaceae, Chitinophagaceae, Pseudomonadaceae, Moraxellaceae, Caulobacteraceae, Xanthomonadaceae, and Pseudonocardiaceae, with different families being dominant in each group. Lactobacillaceae, Xanthomonadaceae, Streptococcaceae, Ruminococcaceae, Bacteroidaceae, Lachnospiraceae, Prevotellaceae, Alteromonadaceae, Aeromonadaceae, and Veillonellaceae were dominant in TW8, corresponding to Xinpu, which is the town with the highest proportion of centenarians. Fusobacteriaceae, Pseudonocardiaceae, Rhodobacteraceae, Micrococcaceae, Moraxellaceae, and Pseudomonadaceae were dominant in TW2, corresponding to Nanzhai, which is the town with the lowest proportion of centenarians.
Based on the species annotation and abundance at the genus level, the 35 most abundant genera were identified and clustered within the sample groups to generate a heat map ( Figure 2C). The top ten genera in the samples were: unidentified Cyanobacteria, Cupriavidus, Ralstonia, Phyllobacterium, Aquitalea, Sediminibacterium, Pseudomonas, Acinetobacter, Caulobacter, and Stenotrophomonas. The dominant genera in TW8 were Lactobacillus, Streptococcus, Bacteroides, Faecalibacterium, Stenotrophomonas, Megamonas, Rheinheimera, Aeromonas, unidentified Lachnospiraceae, and unidentified Prevotellaceae. The dominant genera in TW2 were Fusobacterium, Pseudonocardia, Paracoccus, and Acinetobacter.
To further analyze the phylogenetic relationship of the species at the genus level, representative sequences of the top 100 genera were obtained using multiple sequence alignment. A phylogenetic tree is shown in Figure 2B. The analysis revealed that Ralstonia was the most abundant genus in different groups. For each classification result, the species representing the top ten genera with the highest relative abundance were selected for species-specific classification tree analysis (Supplementary Materials Figure S2). The analysis revealed that Stenotrophomonas sp., a dichloro-diphenyl-trichloroethane (DDT)degrading bacterium, accounted for 1.518% of all species and 2.69% of the selected species. The percentages of S. acidaminiphila, S. koreensis, S. rhizophila, and S. sp PL35a2 S2 among the selected species were 0.01%, 0.29%, 0.67%, and 0.05%, respectively.

Alpha Diversity and Beta Diversity Analyses
Species diversity curves (dilution curves and hierarchical clustering curves) and species accumulation box plots were used to evaluate the differences in species richness and microbial community diversity between the sample groups. The rarefaction curve (Supplementary Materials Figure S3A) tended to be flat, indicating that the amount of sequencing data was sufficient for the analysis. Additional data would generate only a small number of new species (OTUs), and indirectly reflect the abundance of species in the sample. For the rank abundance curve (Supplementary Materials Figure S3B) in the horizontal dimension, the span of the curve gradually increased, indicating a higher richness of species; in the vertical dimension, the curve gradually flattened, indicating a more uniform species distribution. In the box plot of species accumulation (Supplementary Materials Figure S3C), the sample size increased and the box plot position tended to be flat. This meant that the species richness in the environment would not significantly increase with an increasing sample size, indicating that the sample size was sufficient for data analysis. The alpha diversity, assessed by the Shannon index, was significantly different between the TW6 group and the TW1, TW2, and TW3 groups (p < 0.05) ( Figure 3A). The Wilcoxon rank-sum analysis revealed that the alpha diversity in the TW6 group was significantly lower than that in the other groups.
The beta diversity is a comparative analysis of the microbial community composition in different samples. The analysis based on the weighted Unifrac beta-Wilcoxon rank-sum test index ( Figure 3B) revealed statistically significant differences between the TW6 group and the TW1, TW3, TW4, TW5, and TW8 groups (p < 0.05). Significant differences between the TW7 group and the TW1, TW4, TW5, and TW8 groups were also apparent (p < 0.05). To analyze the similarity between samples, the weighted Unifrac distance matrix was used for UPGMA cluster analysis, and the clustering results were integrated and displayed with the relative abundance of each sample at the phylum level ( Figure 3C). The relatively small distance between the TW3 and TW5 groups, and the TW1 and TW7 groups, indicated a similar community structure. Therefore, samples with similar community structures tended to cluster together, while samples with very different communities were further apart. The TW6, TW8, and TW2 groups had the farthest distance and the largest community difference. The communities in each group exhibited certain inter-group differences, indicating the accuracy of grouping. The distance matrix heat map constructed based on the weighted Unifrac distance ( Figure 3D) revealed a big difference between the TW2 and TW8 groups.

Significance Analysis of Species Differences between Groups
The alpha diversity index analysis of differences between groups is used to determine whether the overall community structure in different groups is significantly different. Then, the species responsible for this difference are identified using LEfSe and ternary plot analysis. The identified species constitutes a group biomarker. A histogram of the different species LDA value distribution by LEfSe analysis is shown in Figure 4A, and an evolutionary branch diagram of the different species is shown in Figure 4B. The former revealed biomarker genera with LDA scores greater than four, namely, Phyllobacterium, Flavobacterium, and Cutibacterium in TW1; Acinetobacter, Pseudonocardia, and Paracoccus in TW2; Hydrocarboniphaga in TW3; Sediminibacterium and Aquitalea in TW4; unidentified Cyanobacteria in TW5; Cupriavidus and Caulobacter in TW6; Pseudomonas, Lactobacillales ( Figure 4C) and Psychrobacter in TW7; Bacteroides, Ruminococcaceae ( Figure  4D), Aeromonas, Rheinheimera, and Faecalibacterium in TW8 ( Figure 4).
Next, to identify differences in the dominant species in three groups of samples at each classification level (phylum, class, order, family, genus), the top ten species with the highest average abundance in the three groups of samples at the genus level were used to generate a ternary plot (ternary phase diagram) (Supplementary Materials Figure S4). The dominant genera were: Phyllobacterium in TW1; Acinetobacter and Pseudomonas in TW2; unidentified Cyanobacteria, Ralstonia, and Sediminibacterium in TW3; Aquitalea, Phyllobacterium, Ralstonia, and Sediminibacterium in TW4; unidentified Cyanobacteria and

Significance Analysis of Species Differences between Groups
The alpha diversity index analysis of differences between groups is used to determine whether the overall community structure in different groups is significantly different. Then, the species responsible for this difference are identified using LEfSe and ternary plot analysis. The identified species constitutes a group biomarker. A histogram of the different species LDA value distribution by LEfSe analysis is shown in Figure 4A, and an evolutionary branch diagram of the different species is shown in Figure 4B. The former revealed biomarker genera with LDA scores greater than four, namely, Phyllobacterium, Flavobacterium, and Cutibacterium in TW1; Acinetobacter, Pseudonocardia, and Paracoccus in TW2; Hydrocarboniphaga in TW3; Sediminibacterium and Aquitalea in TW4; unidentified Cyanobacteria in TW5; Cupriavidus and Caulobacter in TW6; Pseudomonas, Lactobacillales ( Figure 4C) and Psychrobacter in TW7; Bacteroides, Ruminococcaceae ( Figure 4D), Aeromonas, Rheinheimera, and Faecalibacterium in TW8 (Figure 4).
To determine individual contributions to the above differences, a similarity percentage (Simper) analysis was performed. The analysis revealed the top ten genera and their abundances, showing their contribution to the difference between two groups (Supplementary Materials Figure S5). The top ten species contributing to the differences between two groups, with their abundance rankings from high to low, were: unidentified Cyanobacteria, Cupriavidus, Ralstonia, Phyllobacterium, Aquitalea, Sediminibacterium, Pseudomonas, Acinetobacter, Caulobacter, and Stenotrophomonas. Next, to identify differences in the dominant species in three groups of samples at each classification level (phylum, class, order, family, genus), the top ten species with the highest average abundance in the three groups of samples at the genus level were used to generate a ternary plot (ternary phase diagram) (Supplementary Materials Figure S4). The dominant genera were: Phyllobacterium in TW1; Acinetobacter and Pseudomonas in TW2; unidentified Cyanobacteria, Ralstonia, and Sediminibacterium in TW3; Aquitalea, Phyllobacterium, Ralstonia, and Sediminibacterium in TW4; unidentified Cyanobacteria and Ralstonia in TW5; Cupriavidus and Caulobacter in TW6; Pseudomonas in TW7; unidentified Cyanobacteria, Stenotrophomonas, and Pseudomonas in TW8. Overall, each group had different dominant genera.
To determine individual contributions to the above differences, a similarity percentage (Simper) analysis was performed. The analysis revealed the top ten genera and their abundances, showing their contribution to the difference between two groups (Supplementary Materials Figure S5). The top ten species contributing to the differences between two groups, with their abundance rankings from high to low, were: unidentified Cyanobacteria, Cupriavidus, Ralstonia, Phyllobacterium, Aquitalea, Sediminibacterium, Pseudomonas, Acinetobacter, Caulobacter, and Stenotrophomonas.

Co-Occurrence Network of the Water Microbiota
The species co-occurrence network diagram enables intuitive visualization of the impact of different environmental factors on microbial adaptability, as well as the dominant species and closely interacting populations of species that occupy a dominant position in a specific environment. These dominant species and populations often play a unique and important role in maintaining a stable structure and function of the environmental microbial community. Lactobacillus and Bifidobacterium are well-known probiotics. In the water microbiota co-occurrence network ( Figure 5), the following bacteria were positively correlated with Lactobacillus: Subdoligranulum, Alloprevotella, Bacteroides, unidentified Prevotellaceae, Faecalibacterium, Qipengyuania, Brachybacterium, unidentified Erysipelotrichaceae, unidentified Acidimicrobiia, Novosphingobium, Dorea, Fusicatenibacter, Empedobacter, Mycobacterium, and Planococcus. In contrast, the following bacteria were negatively correlated with Lactobacillus: Tannerella, unidentified Gammaproteobacteria, Armatimonas, Nannocystis, and Blastopirellula. The following were positively correlated with Bifidobacterium: unidentified Ruminococcaceae, Lacibacter, Pseudomonas, Veillonella, and Prevotella. In contrast, the following were negatively correlated with Bifidobacterium: Corynebacterium, unidentified Rickettsiales, Serratia, unidentified Corynebacteriaceae, Bosea, Bdellovibrio, Brevibacterium, Alcanivorax, Ochrobactrum, and Methyloparacoccus. Hence, in addition to exerting a probiotic effect, Lactobacillus and Bifidobacterium could also exert a synergistic effect via positive and negative regulation of the strains mentioned above.

Predicted Functions of the Water Microbiota
Using the annotation data (Supplementary Materials Figure S6B), the top ten functional annotations for the most abundant bacteria at each annotation level were retrieved, and a histogram of relative functional abundance was generated to visualize the relatively most abundant features at different annotation levels (Supplementary Materials Figure S6A). The identified KEGG pathway annotations mainly focused on cellular processes, environmental information processing, genetic information processing, human diseases, metabolism, and organic systems. Among them, metabolism had the highest proportion, and organic systems had the lowest proportion. Based on the functional annotation and sample abundance information, the top 35 abundance features and their abundances in each sample were used to generate a heat map, and clustered at the level of functional differences. On the horizontal level 1 of the clustering heat map (Figure 6A), the TW8 group had the highest proportion of organic systems and genetic information processing. On the horizontal level 2 of the clustering heat map ( Figure 6B), biosynthesis of other secondary metabolites, metabolism of cofactors and vitamins, immune system, and glycan biosynthesis and metabolism had a higher level of function in the TW8 group than in the others. On the horizontal level 3 of the clustering heat map (Figure 6C), the abundances of exosome, purine metabolism, glycolysis/gluconeogenesis, alanine, aspartate, and glutamate metabolism, amino acid-related enzymes, peptidases, amino sugar and nucleotide sugar metabolism, etc., were only relatively high in the TW8 group. The TW8 group also had unique features in the horizontal level k of the clustering heat map ( Figure 6D). The abundances of K03798, K03723, K02469, K03046, K03657, K05349, K06147, K01955, K01952, K03043, and K02337 were higher in the TW8 group than in the other groups. Overall, the TW8 group showed unique functional characteristics at levels 1, 2, 3, and k.

Predicted Functions of the Water Microbiota
Using the annotation data (Supplementary Materials Figure S6B), the top ten functional annotations for the most abundant bacteria at each annotation level were retrieved, and a histogram of relative functional abundance was generated to visualize the relatively most abundant features at different annotation levels (Supplementary Materials Figure S6A). The identified KEGG pathway annotations mainly focused on cellular processes, environmental information processing, genetic information processing, human diseases, metabolism, and organic systems. Among them, metabolism had the highest proportion, and organic systems had the lowest proportion. Based on the functional annotation and sample abundance information, the top 35 abundance features and their abundances in each sample were used to generate a heat map, and clustered at the level of functional differences. On the horizontal level 1 of the clustering heat map ( Figure 6A), the TW8 group had the highest proportion of organic systems and genetic information processing. On the horizontal level 2 of the clustering heat map (Figure 6B), biosynthesis of other secondary metabolites, metabolism of cofactors and vitamins, immune system, and glycan biosynthesis and metabolism had a higher level of function alanine, aspartate, and glutamate metabolism, amino acid-related enzymes, peptidases, amino sugar and nucleotide sugar metabolism, etc., were only relatively high in the TW8 group. The TW8 group also had unique features in the horizontal level k of the clustering heat map ( Figure 6D). The abundances of K03798, K03723, K02469, K03046, K03657, K05349, K06147, K01955, K01952, K03043, and K02337 were higher in the TW8 group than in the other groups. Overall, the TW8 group showed unique functional characteristics at levels 1, 2, 3, and k.

Discussion
The prediction of local community compositions and functions [44] has improved the understanding of the soil microbial community under current and future changing climatic conditions [45], promoting research on water microbiota. In this study, we aimed to reveal the composition, diversity, and function of the water microbial community in the world's longevity township of Jiaoling, China, and compare the microbiota composition between different towns. We have observed that Xinpu Town and Jiaocheng Town have the largest population over 80 years old. One of the reasons is that Jiaocheng Town is the seat of the county government of Jiaoling County. The urbanization rate is higher than that of other towns, so the population is relatively concentrated. Among them, Xinpu Town has the largest population over a hundred years old. According to our on-site investigation, we found that the environment in this area is closer to the natural environment, with high forest coverage, elegant landscapes and suitable climate. At the same time, as our research results show, we have also observed a high abundance of probiotic Lactobacillus in the water environment of Xinpu Town. We speculate that the reason for the longevity of centenarians may be related to their natural environment. Especially, there is a positive correlation with the abundance of probiotics in the water environment. Of course, these are only speculations based on our existing experimental results, and the probiotics need to be isolated for relevant verification tests in specific circumstances.
Based on the species annotation analysis, the dominant phyla in the samples were Proteobacteria, Cyanobacteria, Bacteroidetes, Firmicutes, Actinabacteria, and Acidobacteria. In terms of Proteobacteria, further analysis showed that it has different categories. Proteobacteria was divided into β-Proteobacteria and γ-Proteobacteria. β-Proteobacteria occupied a dominant position in freshwater and was considered a typical freshwater representative [46,47], while γ-Proteobacteria dominated in sediments. The significant increase in Proteobacteria is directly proportional to the total phosphorus and inorganic nitrogen content, that is, high abundance of Proteobacteria was associated with high nutrient availability [48,49]. Bacteroidetes was another major phylum. Many species in this phylum are related to the gut microbiota of humans [50] and animals [51]. Some species in this phylum have been proposed as indicators for detection of fecal contamination [52]. The relatively high abundance of Bacteroides may be due to untreated wastewater from residential areas. We have also observed that Flavobacteria belonging to Bacteroidetes are heterotrophic bacteria in freshwater [53,54], which play a key role in the carbon cycle (such as biodegradation) [55,56]. Another indicator affected by humans is Faecalibacterium belonging to Firmicutes [57]. However, the relatively high abundance of potential pathogenic bacteria Bacillales was not observed, further indicating that the water resources environment of longevity Township in Jiaoling is good. The dominant genera were Cyanobacteria. Cyanobacteria was a major planktonic bacterial species that can cause FHAB (freshwater harmful algal blooms) [58]. It is an important part of the freshwater ecosystem and plays an important role in nutrient metabolism and the energy cycle [59]. Interestingly, certain species were only dominant in specific areas, which may be closely related to the geographical environment [60][61][62]. This finding might emphasize the influence of the geographical conditions on the composition of the local microbial community.
The environmental microbiota is complex and diverse. This prompted us to question how they affect and are affected by other microorganisms in the local microbiota. Dominant species and populations often play a unique and important role in maintaining the stable structure and function of the microbial community in the environment [63]. Lactobacillus and Bifidobacterium are general probiotics that are widely distributed. For example, large numbers of Lactobacillus have been isolated from the environment (soil and water), food (fermented dairy products, fermented soy products, etc.) [64], and various organisms (human, animals, etc.) [65]. Their presence greatly affects human life and health. Probiotics are currently used to solve problems in various fields, including agriculture, aquaculture, the food processing industry, medical treatment, etc., and, recently, in cosmetics. Molecular methods, such as metagenomics (estimating microbial composition and genome capacity), metatranscriptomics (estimating gene expression), and metaproteomics (estimating protein synthesis), provide insights into the functional profile of the entire soil and water communities. The expression of most functional genes differs among different soil communities, and is driven by the interaction of the climate and soil conditions [66].
In terms of prediction of microbial community function, the high relative abundance of genes related to amino acid metabolism, carbohydrate metabolism, energy metabolism, and xenobiotics biodegradation and metabolism. These functions play an important role in nutrient metabolism and energy circulation. There are some correlations between the abundance of xenobiotic degradation and metabolism genes and the xenobiotic biodegradation and metabolism rate [67,68]. Biodegradation and metabolism genes could be used as indicators of the existence of xenobiotics and their metabolites. It has been reported that Pseudomonas can metabolize chemical pollutants in the environment, such as polycyclic aromatic hydrocarbons [69]. We observed a certain number of Pseudomonas in the water environment, which may indicate the accumulation of polycyclic aromatic hydrocarbons. Polycyclic aromatic hydrocarbons and their biodegradability have been found in many urban rivers [70,71]. Therefore, we must avoid the direct discharge of domestic wastewater into rivers, otherwise it will have a serious impact on human health and microbial ecology.

Conclusions
Our research showed that the spatial changes of the structure and composition of the water microbial community in the world's longevity township of Jiaoling, China, were affected by regional variables. The composition of the bacterial community was mainly Proteobacteria, Cyanobacteria, Bacteroidetes, Firmicutes, Actinabacteria, and Acidobacteria. At the genus level, Cyanobacteria contributed the most to regional differences. We speculate that the reason for the longevity of centenarians may be related to their natural environment. Especially, there is a positive correlation with the abundance of probiotics in the water environment. Of course, these are only speculations based on our existing experimental results, and the probiotics need to be isolated for relevant verification tests in specific circumstances. These data provide a basis for understanding the significant differences in the composition and diversity of water microbial communities in Jiaoling, China. However, the role of these bacteria in health and longevity is unclear, and requires further investigation.  TW1-TW2-TW3, TW2-TW3-TW4,  TW3-TW4-TW5, TW4-TW5-TW6, TW5-TW6-TW7, TW6-TW7-TW8 samples. (Selected top ten species with average abundance at the genus level among the three groups of samples). Figure S5 Simper analysis of contribution of species to differences between groups. (Simper is a decomposition of the Bray-Curtis difference index, used to quantify each species contribution to the difference between two groups. The top ten species and their abundances that contribute to the difference between the two groups). Figure S6

Data Availability Statement:
The authors declare that all the data and the material used in this study are available within this article. All data generated or analysed during this study are available from the corresponding authors upon reasonable request.

Conflicts of Interest:
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.