Abstract

Thermal springs are natural environments present all over the world and their use represents a social-economical resource with an impact on sanus per aquam (SPA) medical and wellness applications. Physical-chemical and microbiological balances characterize these ecological niches and their knowledge is essential to define water properties and support appropriate management. This study is a pilot application of a larger research project, involving metagenomics and aimed to fingerprint springs and map SPA biodiversity. Waters and their deposits were collected in six thermal springs from the Lazio region in Italy. The phylogenic microbial profiles performed by Next Generation Sequencing (NGS) analysis showed a clear separation between different springs. Statistical analyses revealed correlations between the abundance of specific bacteria and environmental variables. Temperature, Sodium and H2S levels appear to play a key role in influencing the microbiota. The extension of this model to other springs will contribute to characterize and map the microbial community in thermal springs, allowing associations with chemical-physical factors. Biodiversity is a still underestimated property of thermal springs and a key element in several SPA applications. The Atlas progress is shedding light on biotic and abiotic components in these ecological niches, opening further perspectives for supporting appropriate use and management of thermal waters.

INTRODUCTION

Since ancient times, thermal springs were known and used for health and social purposes, in natural or artificial pools (Van Tubergen & Van Der Linden 2002). Nowadays, sanus per aquams (SPAs) represent a worldwide resource with a valuable impact for medical and wellness applications, including hydro-balneotherapies, rehabilitation, physical activity and recreational uses (Frosch 2007). However, a more updated and rigorous approach is required to distinguish the different springs, provide evidences for possible therapeutic properties or for the development of bioactive products. Thermal springs are natural water environments present all over the world (Mirete et al. 2016). These habitats comprise ecological niches, characterized by a unique microbial flora active within specific chemical and physical balances (Marsh & Larsen 1953; Yazdi et al. 2015). These unusual or even extreme environments can represent a major source of microbial diversity and a potential for research and development in biotechnology (Trabelsi et al. 2016; Valeriani et al. 2016a; Urbieta et al. 2016). Moreover, knowledge on extremophiles inhabiting hot springs can provide insights into the origin and early evolution of life, as they are considered the direct descendants of the earliest living forms on Earth (Weiss et al. 2016).

Thermal water springs harbor a precious biodiversity that is still not fully classified and exploited, mainly due to technical difficulties in satisfying the specific culture requirements (López-López et al. 2013). In these environments, the biotic components are in a delicate and complex equilibrium with the physical-chemical proprieties (Fierer et al. 2012). The possibility to study the microflora and its modifications can provide a candidate novel indicator to monitor perturbations in the ecological niche. Microbial communities, indeed, represent a constitutive property of thermal waters and play a role in different processes, such as biofilm formation and carbonate precipitation (Starke et al. 2013). This natural constituent is also exploited for several SPA applications such as sludge ripening processes, development of cosmetics or wellness products (Centini et al. 2015).

Analysis of microflora DNA (mfDNA) can provide an overview on the microbial community structure, starting from the detection of each species genome, independently, by knowing culture requirements. Currently, molecular methodologies used for the characterization of extremophiles niches are based on the amplification and sequencing of 16S rRNA genes, that are highly conserved during evolution. Metagenomics and 16S amplicon sequencing represent preferred and consolidated approaches for the characterization of complex matrices, including waters, sewage sludge or sediments (Fierer et al. 2012; Giampaoli et al. 2014; Delforno et al. 2017). By Next Generation Sequencing (NGS), the knowledge on microbial biodiversity in thermal waters can be improved, shedding light on these ecological niches (Inskeep et al. 2013).

In order to evaluate the possibility to characterize the thermal springs by their biodiversity, we applied a massive sequencing approach based 16S amplicon sequencing and set up a map of their distribution. This knowledge can provide additional information on the biological component of the waters, suggesting new hints to understand these ecological niches and their potentialities. Based on these considerations, a national working group of the Italian Society of Hygiene (GSMS-SItI) has proposed the research project Microflora Thermarum Atlas (MTA), to map and characterize biodiversity in Italian thermal springs, considering their natural proprieties and supporting adequate choices for managing surveillance and environmental health in these unique contests. Here, we summarize for the first time the MTA approach and present a pilot application to central Italy, comprising an area around Rome within an historical district where SPA played a role since ancient times.

MATERIALS AND METHODS

Study area and sampling

The thermal spring waters and deposits were collected in the Lazio region (17,242 km2) located in central Italy. In Figure 1, sampling points are indicated with red circles (GPS coordinates are summarized in the supplementary material, S1, available with the online version of this paper). Sampling was done during April–June 2015: water samples were collected directly in the spring, in 1 L borosilicate glass sterile bottles. The thermal spring deposits (3–7 g) covered by microbial mat and in contact with water, were sampled by scraping, with a sterile scalpel. To enhance consistency and reduce local heterogeneity, each sampling was collected in triplicate from three independent points within a 20 cm × 20 cm area and then pooled together. All samples were immediately stored on-site at 4–6 °C and transferred within 3–5 h to the laboratory for processing.

Figure 1

Microflora Thermarum Atlas project: the Italy map shows different thermal spring areas in different regions, while the enlarged detail of the Lazio district indicates the sampling points collected in this study. In the regional map, the thermal spring sites are present in an area characterized by travertine deposits, volcanic units and carbonate rich structures, corresponding to the vestiges of extreme ecological niches where ancient microflora interacted in carbonate and sulfur cycles. Please refer to the online version of this paper to see this figure in color: http://dx.doi.org/10.2166/ws.2017.215.

Figure 1

Microflora Thermarum Atlas project: the Italy map shows different thermal spring areas in different regions, while the enlarged detail of the Lazio district indicates the sampling points collected in this study. In the regional map, the thermal spring sites are present in an area characterized by travertine deposits, volcanic units and carbonate rich structures, corresponding to the vestiges of extreme ecological niches where ancient microflora interacted in carbonate and sulfur cycles. Please refer to the online version of this paper to see this figure in color: http://dx.doi.org/10.2166/ws.2017.215.

Physicochemical water characterization

Water quality parameters, namely temperature, pH, electrical conductivity, free CO2, and H2S were measured in situ using the relevant field meters (Mettler Toledo meters, UK and QRAE II, RAE Systems, USA) while other physical-chemical parameters were determined in the laboratory, applying the APAT 2003 Italian official methods.

DNA extraction

Water samples (1–2 L) were filtered with a 0.47 μm polycarbonate membrane; the filter membrane was transferred upside down into an incubation dish with 2 mL of Lysis buffer (AquaScreen® Fast Extract Kit, Minerva BioLabs, Germany) and incubated at 37 °C for 30 min. Subsequently, the lysis solution was transferred to 2.0 mL tubes with 0.1 mg of glass beads (Sigma Aldrich, USA), vortexed for 1 min and incubated at 56 °C for 15 min. Following purification steps were performed according to AquaScreen® Fast Extract Kit manufacturer's protocol (Minerva BioLabs, Germany). For thermal spring deposits, an aliquot (0.5 g) was extracted following the FastDNA® SPIN Kit protocols (MP Biomedical, USA).

Bacterial community 16S profiling

16S rRNA paired-end sequencing was performed according with ‘16S Metagenomic Sequencing Library Preparation’ protocol (Part# 15044223 rev. A; Illumina, USA). Briefly, the primers containing Illumina adapter and linker sequence and targeting the V1-V2 regions of bacterial 16S rRNA genes were used (Kittelmann et al. 2013). Three libraries with unique tags were generated for each sample as technical replicates. Each amplification reaction had a total volume of 25 μL containing 12.5 μL of KAPA HiFi HotStart ReadyMix (Roche, Pleasanton, CA), 5 μL of each primer (1 μM), and 2 μL template DNA. Reactions were carried out on a Techne®TC-PLUS thermalcycler (VWR International, LLC, Radnor, USA). Thermal cycling conditions were as follows: an initial denaturation at 94 °C for 1 min, and 30 cycles at 94 °C for 20 s, 53 °C for 25 s, and 68 °C for 45 s, with a final extension at 68 °C for 10 min. Following amplification, 5 μL of polymerase chain reaction (PCR) product from each reaction was used for agarose gel (1%) electrophoresis to confirm amplification. The final concentration of cleaned DNA amplicon was determined using the Qubit PicoGreen dsDNA BR assay kit (Invitrogen, Grand Island, NY, USA) and validated on Bioanalyzer DNA 1000 chip (Agilent, USA). Libraries were prepared using the MiSeq Reagent Kit Preparation Guide (Illumina, San Diego, CA, USA). Briefly, the combined sample library was diluted to 4 nM, denatured with 0.2 N fresh NaOH, diluted to 4 pM by addition of Illumina HT1 buffer, and then mixed with an equal volume of 4 pM PhiX (Illumina, San Diego, CA, USA). The library (600 μL) was loaded with read 1, read 2 and index sequencing primers on a 500-cycle (2 × 250 paired ends) reagent cartridge (Illumina), and run on a MiSeq sequencer (Illumina).

Raw sequence data was processed using an in-house pipeline which was built on the Galaxy platform and incorporated various software tools to evaluate the quality of the raw sequence data (e.g. FastQC, http://www.bioinformatics.babraham.ac.uk/projects/fastqc/). All data sets were rigorously screened to remove low quality reads (short reads >200 nt, zero-ambiguous sequences). Demultiplexing was performed to remove PhiX sequences and sort sequences; moreover, to minimize sequencing errors and ensure sequence quality, the reads were trimmed based on the sequence quality score using Btrim (Kong 2011). Operational taxonomic units (OTUs) were clustered at a 97% similarity level and final OTUs were generated based on the clustering results and taxonomic annotation of individual OTUs was based on representative sequences using RDP's 16S Classifier 2.5. Relative abundances of community members were determined with rarefied data and summarized at each taxonomic level (cut off 0.2%). The sequence reads were analyzed, also, in the cloud environment BaseSpace through the 16S Metagenomics app (version 1.0.1; Illumina®): the taxonomic database used was the Illumina-curated version (May 2013 release of the Greengenes Consortium Database; Wang et al. 2013).

Statistical analysis

Microbiota phylogenic distribution of thermal spring waters was analyzed using Metagenassist (Arndt et al. 2012): input files were produced in the format of CSV files, including taxonomic profile file with taxonomic abundance for each sample and further file, containing metadata information (e.g. sodium ions as indicators of salinity and H2S content). Then, water samples and sediments were clustered for salinity and H2S content. Three categories for each cluster (salinity −< 50 mg/L, 51 < T ≤ 200 mg/L and >200 mg/L- and H2S content – sulphurous, slightly sulphurous and not sulphurous) were generated in order to evaluate sample distribution. Multivariate analysis, the PCoA and partial least square-discriminant analysis (PLS-DA) were performed in order to investigate the dissimilarity between groups, using METAGENassist platform (Arndt et al. 2012). We performed feature selection using PLS-DA and 10-fold cross validation to tune algorithm parameters and to check model validity. The putative functional profiles based on the 16S community composition were investigated by automated taxonomic-to-phenotypic mapping using a METAGENassist platform and NCBI microbial taxonomy (Arndt et al. 2012). Agglomerative hierarchical clustering were performed by treating each individual as a separate cluster and then proceeds to combine them until all samples belong to a single cluster. For similarity measure was used Pearson's correlation and the clustering algorithms were set as Ward's linkage (clustering to minimize the sum of squares of any two clusters). The result is presented as a dendrogram in combination with a heatmap.

RESULTS AND DISCUSSION

The Microflora Thermarum Atlas is based on massive sequencing of microflora DNA (mfDNA) and bioinformatic analysis of outputs (Figure 2). The project is aimed to realize a biodiversity map of thermal spring microflora. The final aim is to localize and characterize specific natural properties of the waters, provide information on their geographical distribution and the different SPA applications, including medical therapies or recreational uses in SPA pools. In this study, for the first time, we summarize the Atlas approach and methods, reporting a comprehensive, culture-independent survey of microbial communities in six geothermal springs collected in the Italian region of Lazio (Figure 1). Interestingly, thermal springs location overlapped with several travertine complexes that are predominantly located on the margins of a hydrogeological area belonging to volcanic and carbonate domains. This observation agrees with the role of thermal spring microflora in carbonate precipitation and sulfur oxidizing processes in ancient ecological niches, when chemolithoautotrophic processes occurred and primordial life was supposed to originate, as shown by studies on LUCA, the Last Universal Common Ancestor (Dupraz & Visscher 2005; Weiss et al. 2016).

Figure 2

Flow chart of the project. Key steps from sampling to filling the database and drawing the map. Metagenomics methods are central in this process allowing a comprehensive analysis of the SPA microbiota.

Figure 2

Flow chart of the project. Key steps from sampling to filling the database and drawing the map. Metagenomics methods are central in this process allowing a comprehensive analysis of the SPA microbiota.

The physicochemical characteristics at the water sampling points, are given in Table 1. The temperature values varied, ranging from 20.5 to 54 °C at the sampling point. All waters showed a neutral or slightly acid pH that varied from 6.18 to 6.80, with no main variations for conductivity (2,380–3,640 mg/L) and fixed residue (2,418–3,000 mg/L). Hydrogen sulfide concentration ranged from <0.2 mg/L to 24 mg/L allowing a classification in three categories: sulphurous (TR, 24 mg/L), slightly sulphurous (BU and SU, 2 and 6 mg/L, respectively) and not sulphurous (BA, PA and CA, all with values <0.2 mg/L). Also for sodium ion values, the sampled thermal springs could be aggregated in three groups: <50 (SU, BU), 51–200 (TR, VA and PA) and >200 (CA). Bicarbonates ranged from 763 to 1,455 mg/L.

Table 1

Physical and chemical characteristics of the thermal spring waters (NA = data not available)

Physical/chemical parameter Thermal water spring
 
TR BA PA CA BU SU 
Temperature (°C) 20.5 40.9 39.6 21 54 25 
Conducibility (μS/cm) 2,710 2,380 2,730 NA 2,590 3,640 
pH 6.18 6.3 6.3 6.3 6.8 6.7 
Fixed residue at 180 °C (mg/L) 2,794 2,418 2,768 NA 3,000 2,730 
CO2 free (mg/L) 1,079 380 120 123 871 740 
H2S (mg/L) 24 <0.2 <0.2 <0.2 
Potassium ions (mg/L) 28 28 70 292 42 NA 
Magnesium ions (mg/L) 107 112 89 118 125 187 
Strontium ions (mg/L) 8.4 10 11 8.4 13 NA 
Lithium ions (mg/L) 0.14 0.11 0.48 6.6 0.19 NA 
Bromine ions (mg/L) <0.1 <0.1 <0.1 <0.1 NA 
Sodium ions (mg/L) 126 38 133 1,954 35 40 
Iron ions (mg/L) 0.1 0.6 3.4 <0.1 0.5 0.6 
Bicarbonates (mg/L) 1,455 1,068 763 1,302 840 1,002 
Sulphates (mg/L) 705 1,100 1,325 1,133 1,181 180 
Arsenic (mg/L) 0.22 0.03 0.204 74 0.234 NA 
Silica (mg/L) 22 37 109 23 62 NA 
Barium (mg/L) 0.03 0.06 0.03 <0.01 0.05 NA 
Boron (mg/L) 3.6 0.8 14 50 2.2 NA 
Manganese (mg/L) <0.01 0.19 0.99 <0.01 0.03 NA 
Physical/chemical parameter Thermal water spring
 
TR BA PA CA BU SU 
Temperature (°C) 20.5 40.9 39.6 21 54 25 
Conducibility (μS/cm) 2,710 2,380 2,730 NA 2,590 3,640 
pH 6.18 6.3 6.3 6.3 6.8 6.7 
Fixed residue at 180 °C (mg/L) 2,794 2,418 2,768 NA 3,000 2,730 
CO2 free (mg/L) 1,079 380 120 123 871 740 
H2S (mg/L) 24 <0.2 <0.2 <0.2 
Potassium ions (mg/L) 28 28 70 292 42 NA 
Magnesium ions (mg/L) 107 112 89 118 125 187 
Strontium ions (mg/L) 8.4 10 11 8.4 13 NA 
Lithium ions (mg/L) 0.14 0.11 0.48 6.6 0.19 NA 
Bromine ions (mg/L) <0.1 <0.1 <0.1 <0.1 NA 
Sodium ions (mg/L) 126 38 133 1,954 35 40 
Iron ions (mg/L) 0.1 0.6 3.4 <0.1 0.5 0.6 
Bicarbonates (mg/L) 1,455 1,068 763 1,302 840 1,002 
Sulphates (mg/L) 705 1,100 1,325 1,133 1,181 180 
Arsenic (mg/L) 0.22 0.03 0.204 74 0.234 NA 
Silica (mg/L) 22 37 109 23 62 NA 
Barium (mg/L) 0.03 0.06 0.03 <0.01 0.05 NA 
Boron (mg/L) 3.6 0.8 14 50 2.2 NA 
Manganese (mg/L) <0.01 0.19 0.99 <0.01 0.03 NA 

Sampled waters and sediments harbored an interesting biodiversity as reported in Figure 3 and Table 2. A total of 2,407,415 sequence reads were generated after NGS analysis. The number of sequences for each sample ranged from 30,431 to 352,625, leading to the identification of 970 OTUs defined at 97% identity. Rarefaction curves were calculated for each sample (Figure S2, available with the online version of this paper), showing an adequate and reliable sampling and sequencing effort for describing the bacterial community (Wen et al. 2017). As determined by the Shannon index (biodiversity) and the number of identified species (richness), the complexity of bacterial communities varied strongly between waters (0.884–1.644; 89–529) and sediments (1.607–2.259; 339–1,121). Sediments showed a more heterogeneous structure with an increase in those heterotrophic genera that may probably act in capturing and metabolizing the organic matter, also favoring the precipitation of carbonate and travertine concretions (Pentecost 2005; Schubotz et al. 2015).

Table 2

Summary of NGS analysis after quality assessment of sequences

Sample ID Sample description Number of reads PF Reads PF classified to genus (%) Shannon α-diversity Number of species 
BU Thermal spring water 30,431 73.38 0.884 89 
sed_BU Thermal spring deposit 110,742 96.95 2.259 339 
BA Thermal spring water 145,055 88.77 1.644 305 
sed_BA Thermal spring deposit 352,625 85.07 2.120 1.121 
CA Thermal spring water 167,322 85.05 1.583 436 
sed_CA Thermal spring deposit 240,333 74.35 1.607 715 
PA Thermal spring water 151,974 79.75 1.237 423 
sed_PA Thermal spring deposit 231,81 71.16 1.724 694 
SU Thermal spring water 251,624 83.31 1.032 529 
sed_SU Thermal spring deposit 294,872 88.08 2.171 724 
TR Thermal spring water 105,03 89.48 1.037 277 
sed_TR Thermal spring deposit 325,597 71.66 2.304 757 
Sample ID Sample description Number of reads PF Reads PF classified to genus (%) Shannon α-diversity Number of species 
BU Thermal spring water 30,431 73.38 0.884 89 
sed_BU Thermal spring deposit 110,742 96.95 2.259 339 
BA Thermal spring water 145,055 88.77 1.644 305 
sed_BA Thermal spring deposit 352,625 85.07 2.120 1.121 
CA Thermal spring water 167,322 85.05 1.583 436 
sed_CA Thermal spring deposit 240,333 74.35 1.607 715 
PA Thermal spring water 151,974 79.75 1.237 423 
sed_PA Thermal spring deposit 231,81 71.16 1.724 694 
SU Thermal spring water 251,624 83.31 1.032 529 
sed_SU Thermal spring deposit 294,872 88.08 2.171 724 
TR Thermal spring water 105,03 89.48 1.037 277 
sed_TR Thermal spring deposit 325,597 71.66 2.304 757 

For each sampled point, the number of obtained reads after passing filter (PF), the percentage of genus classification, the Shannon α-Diversity and number of species are indicated for each sample.

Figure 3

Relative genus abundance and thermal spring waters' diversity.

Figure 3

Relative genus abundance and thermal spring waters' diversity.

Sequence analysis provided a total of 12 main bacterial phyla (cut off >0.2%). At this level, waters did not show significant variations with most of the sequences assigned to Proteobacteria (82.3–98.6%) and Firmicutes (0.6–5.6%) phyla. At genus level, waters displayed a different profile with a significant fluctuation between identified bacteria (Figure 3). Thiofaba is the main genus in BU (55.8%), TR waters (45.8%), where H2S values are the highest observed. This genus includes obligate chemolithoautotrophic sulfur-oxidizing bacteria utilizing H2S as electron donor for CO2 reduction (Van Gemerden 1993). The role of sulfur-oxidizing bacteria is crucial as biotic factors can influence carbonate equilibrium, promoting carbonate dissolution (Dupraz & Visscher 2005). In the carbonate equilibrium are involved also sulfate-reducing bacteria, such as Sulfurospirillum (BU = 11%; TR = 13.8% and BA = 25.7%) and Giesbergeria (BA = 11.1%). In CA Thiofaba and also other sulfur-oxidizing bacteria were present: Thiovirga (23.2%) and Sulfurimonas (11.9%). Interestingly, the presence of the moderately halophilic Thiovirga was observed in CA and SU water samples, both rich in sodium ion (CA: >200 mg/L; SU: 40 mg/L).

In some springs, probably, an environmental contamination occurred by sludge and/or fresh waters, as suggested by the identification of bacteria such as Arcobacter (TR: 25.3%; BA: 46.9%) and Hydrogenophaga (SU: 69.5%). This observation is in accordance with the promiscuous structure of these spring environments and the presence of a nearby creek. Conversely, PA microflora differs from other springs, probably due to the specific chemical characteristics, such as the presence of iron (3.4 mg/L) and silica (109 mg/L in Table 1). In particular, Gallionella (37.2%), iron-oxidizing bacteria, and Azoarcus (21.9%), sulfide-oxidizing bacteria, are the dominating genera.

Biodiversity was also related to temperature values, showing a Shannon index decrease with the increase of temperature, in agreement with other observations and with data from the statistical analysis as shown in the supplementary material, Figure S3 (available with the online version of this paper) (Miller et al. 2009). For example, in BU, where the temperature was the highest (54 °C) the calculated Shannon index was 0.888 whether the decrease of temperature as in CA (21 °C) the Shannon index was the highest observed (1.583).

The distribution of the genera seems to be well modeled by PLS discriminant analysis (Figure 4), supporting the hypothesis that some chemical characteristics are closely related to the prevalent microflora and biodiversity composition. The phylogenic grouping showed the clear separation between sulphurous and non sulphurous waters. Surprisingly, samples with different sodium ions content (CA and SU) disclosed a different community structure, supporting the leading role of salinity on microflora composition, as also underlined by metabolic data (Figure 5).

Figure 4

Partial least square-discriminant analysis (PLS-DA) depicts Pearson distance between water samples using phylogeny distribution based on 16S rRNA genes. Samples are, respectively, colored according to environmental categories of salinity (sodium values mg/L) in (a) and H2S content (sulphurous, slightly sulphurous and not sulphurous) in (b).

Figure 4

Partial least square-discriminant analysis (PLS-DA) depicts Pearson distance between water samples using phylogeny distribution based on 16S rRNA genes. Samples are, respectively, colored according to environmental categories of salinity (sodium values mg/L) in (a) and H2S content (sulphurous, slightly sulphurous and not sulphurous) in (b).

Figure 5

The heatmap with summarized putative functional profiles based on the 16S community composition were investigated by automated taxonomic-to-phenotypic mapping using a METAGENassist platform and NCBI microbial taxonomy (Arndt et al. 2012); data: each column represents a pathway and each row represents a water sample. The row Z-score or scaled expression value of each feature is plotted in red–blue color scale. The red color of the tile indicates high abundance and blue indicates low abundance. Please refer to the online version of this paper to see this figure in color: http://dx.doi.org/10.2166/ws.2017.215.

Figure 5

The heatmap with summarized putative functional profiles based on the 16S community composition were investigated by automated taxonomic-to-phenotypic mapping using a METAGENassist platform and NCBI microbial taxonomy (Arndt et al. 2012); data: each column represents a pathway and each row represents a water sample. The row Z-score or scaled expression value of each feature is plotted in red–blue color scale. The red color of the tile indicates high abundance and blue indicates low abundance. Please refer to the online version of this paper to see this figure in color: http://dx.doi.org/10.2166/ws.2017.215.

When analyzing microbial relative abundance at species level, by the principal coordinate analysis (PCoA) method, the overlap of all waters corresponded to the respective sediments (Figure 6). This finding can allow us to suppose a possible link between microbial mat and water, benthonic and planktonic components, as previous studies already underlined (Ward et al. 1998). Moreover, the distribution of sampling points and their hydrogeological setting seemed to play a possible role in defining microflora biodiversity. Indeed, the phylogenic grouping showed a clear separation in respect to the characteristics of the water flow. Specifically, thermal waters segregated following the hydrogeology on which they flow, e.g. carbonate unit as in BA, TR and CA or volcanic unit as in BU, PA or SU (Figure 6). The described microbial communities could be linked to different parameters including fluid hydrodynamics, chemical-physical properties, hydrogeological contexts, as suggested by data reported in previous studies (Di Salvo et al. 2013).

Figure 6

Principal coordinate analysis (PCoA) scatterplot of the normalized relative abundance of all samples. Data are plotted at the species-level classification.

Figure 6

Principal coordinate analysis (PCoA) scatterplot of the normalized relative abundance of all samples. Data are plotted at the species-level classification.

The recurrent presence of specific microbial groups enlightens a kind of biological signature for thermal spring waters, showing an individual bio-fingerprint as well as a common ground with respect to other aquatic environments (Kemp & Aller 2004; Miller et al. 2009). Previous studies have already supposed a strong correlation between the composition of the microbial community and the geochemical properties of hot springs as a bidirectional phenomenon (Inskeep et al. 2013). Therefore, knowledge on microflora composition can support both the acquisition of information on the spring geochemical background as well as on the metabolic properties and biotechnological potentials of that water. Indeed, 16S rRNA gene amplification and sequencing emerged as powerful tools, for elucidating the biodiversity of complex samples and to study of putative metabolic pathways (Inskeep et al. 2013; López-López et al. 2013). The Atlas approach has the advantage of collecting and mapping information from different springs using a common procedure and allowing comparison of data obtained through identical protocols (Figure 2). The availability of a common strategy can support reproducibility and comparison of data, avoiding bias. Indeed, amplification as well as sampling or extraction steps may be critical in selecting microflora species, providing artifacts. All Atlas data are accessible in a single map, because both raw data and final outputs are joined within a single shared database. Presently, the database is analyzed in the Cloud of the University of Rome ‘Foro Italico’ and accessible online to the research network at www.mfATLAS.it, but it is designed to be extended and harbor information also from specific SPA plants and artificial pools, including water management, environmental and epidemiological data. Focus of this study is the acquisition of knowledge on the natural microflora inhabiting thermal waters and its geographical distribution. The characterization of metabolic properties and the identification of new species has interest for classification purposes, but also implications for biotechnology potentials or local SPA applications in the field of wellness. Specific issues involve the possibility to characterize individual properties of the water and its consistency during the time, by joining mfDNA and physical-chemical data into a kind of SPA fingerprint. This can reveal useful for tracing purposes, providing additional information to monitor spring stability even after extraordinary conditions, such as pollutions or earthquakes. Moreover, several SPA-derived products can also be fingerprinted (data not shown), such as extracts, lavenders, cosmetics, sludges and traced through the specific microbial component.

The main limits of this Atlas approach are related to the culture-independent method that would not allow viable isolates. This can be overcome by flanking classical microbiology techniques, even if metagenomics tends to provide a high percentage of unknown species in environmental matrices and only a few will be cultivable on the available media. We observed about 0.3–1% of unclassifiable bacteria kingdom per water sample and 0.03–1% in sediments, but at species level the unclassifiable percentage could exceed 70%. Then, the metagenomics approach draws a general view of the biodiversity framework, without providing an analytical detail of specific microbial contaminations that could be present at very low number, such as pathogens. Indeed, the Atlas metagenomics strategy is not a detection or surveillance strategy, but is aimed to draw the structure of the natural microflora. For identification and quantification of a single microorganism, other molecular methods should be applied, such as real time PCR amplification or the classical culture based techniques. Finally, the NGS equipment is still expensive, as well as reagents, and requires trained personnel for both laboratory and bioinformatic analysis. However, this approach can be easily centralized reducing the final costs for a single analysis. Our laboratory is receiving from different Italian regions several waters, sediments or even mfDNA samples processed following same protocols. The present report on the Lazio region has itself some limitations that may be reduced in the long period, such as the lack of sampling at different seasonalities, a further simplification of protocols, the sequencing and bioinformatics approach (Huson et al. 2011). Nevertheless, observed results are already reliable, showing a consistent stability in the main microbiota components, as already reported for other ecological niches, including soils or biological fluids (Giampaoli et al. 2014; Valeriani et al. 2016b). This is in accordance with the stationary hydrological and physical-chemical conditions of the springs, obviously when excluding major disturbances, such as contaminations, earthquakes or other hydrogeological events. Then, we limited this report to the general approach considering natural springs puddles, without including data from SPA plants by sampling water from therapeutic devices or pools, or simply along the pipeline, where temperature, biofilms, anthropic contamination, materials or treatments can interfere with the natural microflora background, inducing differences that need to be detailed and discussed case by case. Finally, we focused only on analysis of 16S rDNA, considering bacteria taxa and excluding other eukaryotic species belonging to the microflora such as protozoa or metazoans, that can represent important components in biotic balances (Giampaoli et al. 2014; Sanli et al. 2015).

In conclusion, we present an Atlas for mapping and characterizing the microflora of thermal springs and SPA. The project is described in its basic structure and methods, reporting an application to the Lazio region. Thermal waters showed to contain metabolically versatile and specialized bacteria, belonging to an ecological niche with peculiar chemical and physical properties. If salinity appears to be the main factor shaping the structure of the microbial community, it does not drive their diversity. Hydrogen sulfide seems to influence the structure of SPA microbiota, probably acting on sulfur metabolic pathways, as well as on microorganism resistance to redox conditions (Giampaoli et al. 2013). This, as other chemical-physical factors, play a role in establishing biodiversity within a complex and multifactorial network of interactions. Statistical analyses revealed significant correlations between the relative abundance of thermal water bacteria and environmental variables, showing the feasibility of a promising approach open to further investigations. The extension of this model to other Italian and European thermal springs is open and welcome.

CONCLUSIONS

Data from the MTA in the Lazio region are reported, mapping spring microflora biodiversity and presenting the general approach, methods and potentials. Thermal SPA water mfDNA was analyzed by a 16S amplicon NGS approach, showing the complexity of a bacterial community that is metabolically active in a peculiar ecological niche. SPA waters appear as an alive and vital biological fluid, opening to further perspectives for an appropriate management and surveillance of these environments and their exploitation for therapeutic, recreational or biotechnological uses. The proposed study model is aimed to extend the map and characterize microflora biodiversity and geographical distribution through different ongoing projects and collaborations.

ACKNOWLEDGEMENTS

The MTA project is promoted by the GSMS-SItI: National Working Group of the Italian Society of Hygiene, Preventive Medicine and Public Health. Authors thank GeckoBiotech research start up for realizing the online access to the database at www.mfATLAS.it and Dr. Elena Scaramucci for updating references and editing the manuscript.

REFERENCES

REFERENCES
Arndt
,
D.
,
Xia
,
J.
&
Liu
,
Y.
2012
METAGENassist: a comprehensive web server for comparative metagenomics
.
Nucleic Acids Research
40
,
88
95
.
Centini
,
M.
,
Tredici
,
M. R.
,
Biondi
,
N.
,
Buonocore
,
A.
,
Maffei Facino
,
R.
&
Anselmi
,
C.
2015
Thermal mud maturation: organic matter and biological activity
.
International Journal Cosmetic Science
37
(
3
),
339
347
.
Delforno
,
T. P.
,
Lacerda Júnior
,
G. V.
,
Noronha
,
M. F.
,
Sakamoto
,
I. K.
,
Varesche
,
M. B.
&
Oliveira
,
V. M.
2017
Microbial diversity of a full-scale UASB reactor applied to poultry slaughterhouse wastewater treatment: integration of 16S rRNA gene amplicon and shotgun metagenomic sequencing
.
Microbiology open
6
(
3
).
Di Salvo
,
C.
,
Mazza
,
R.
&
Capelli
,
G.
2013
Gli acquiferi in travertino del Lazio: schemi idrogeologici e caratteristiche chimico-fisiche
.
Rendiconti Online – Società Geologica Italiana
27
,
54
76
.
Dupraz
,
C.
&
Visscher
,
P. T.
2005
Microbial lithification in marine stromatolites and hypersaline mats
.
Trends in Microbioliology
13
,
429
438
.
Fierer
,
N.
,
Leff
,
J. W.
,
Adams
,
B. J.
,
Nielsen
,
U. N.
,
Bates
,
S. T.
,
Lauber
,
C. L.
,
Owens
,
S.
,
Gilbert
,
J. A.
,
Wall
,
D. H.
&
Caporaso
,
J. G.
2012
Cross-biome metagenomic analyses of soil microbial communities and their functional attributes
.
Proceedings of the National Academy of Sciences of United States of America
109
(
52
),
21390
21395
.
Giampaoli
,
S.
,
Valeriani
,
F.
,
Gianfranceschi
,
G.
,
Vitali
,
M.
,
Delfini
,
M.
,
Festa
,
M. R.
,
Bottari
,
E.
&
Romano Spica
,
V.
2013
Hydrogen sulfide in thermal spring waters and its action on bacteria of human origin
.
Microchemical Journal
108
,
210
214
.
Giampaoli
,
S.
,
Berti
,
A.
,
Di Maggio
,
R. M.
,
Pilli
,
E.
,
Valentini
,
A.
,
Valeriani
,
F.
,
Gianfranceschi
,
G.
,
Barni
,
F.
,
Ripani
,
L.
&
Romano Spica
,
V.
2014
The environmental biological signature: NGS profiling for forensic comparison of soils
.
Forensic Sciences International
240
,
41
47
.
Huson
,
D. H.
,
Mitra
,
S.
,
Ruscheweyh
,
H. J.
,
Weber
,
N.
&
Schuster
,
S. C.
2011
Integrative analysis of environmental sequences using MEGAN4
.
Genome Research
21
,
1552
1560
.
Inskeep
,
W. P.
,
Jay
,
Z. J.
,
Tringe
,
S. G.
,
Herrgard
,
M.
&
Rusch
,
D. B.
2013
The YNP Metagenome project: environmental parameters responsible for microbial distribution in the yellowstone geothermal ecosystem
.
Frontiers in Microbiology
4
,
67
.
Kittelmann
,
S.
,
Seedorf
,
H.
,
Walters
,
W. A.
,
Clemente
,
J. C.
,
Knight
,
R.
,
Gordon
,
J. I.
&
Janssen
,
P. H.
2013
Simultaneous amplicon sequencing to explore co-occurrence patterns of bacterial, archaeal and eukaryotic microorganisms in rumen microbial communities
.
PLoS ONE
8
,
47879
.
López-López
,
O.
,
Cerdán
,
M. E.
&
González-Siso
,
M. I.
2013
Hot spring metagenomics
.
Life
2
,
308
320
.
Marsh
,
C. L.
&
Larsen
,
D. H.
1953
Characterization of some thermophilic bacteria from the hot springs of Yellowstone National Park
.
Journal Bacteriology
65
,
193
197
.
Mirete
,
S.
,
Morgante
,
V.
&
Gonzalez-Pastor
,
J. E.
2016
Functional metagenomics of extreme environments
.
Current Opinion Biotechnology
38
,
143
149
.
Pentecost
,
A.
2005
Travertine
.
Springer
,
Heidelberg, Germany
.
Sanli
,
K.
,
Bengtsson-Palme
,
J.
,
Nilsson
,
R. H.
,
Kristiansson
,
E.
,
Alm Rosenblad
,
M.
,
Blanck
,
H.
&
Eriksson
,
K. M.
2015
Metagenomic sequencing of marine periphyton: taxonomic and functional insights into biofilm communities
.
Frontiers in Microbiology
6
,
1192
.
Schubotz
,
F.
,
Hays
,
L. E.
,
Meyer-Dombard
,
D. R.
,
Gillespie
,
A.
,
Shock
,
E. L.
&
Summons
,
R. E.
2015
Stable isotope labeling confirms mixotrophic nature of streamer biofilm communities at alkaline hot springs
.
Frontiers in Microbiology
6
,
42
.
Starke
,
V.
,
Kirshtein
,
J.
,
Fogel
,
M. L.
&
Steele
,
A.
2013
Microbial community composition and endolith colonization at an Arctic thermal spring are driven by calcite precipitation
.
Environmental Microbiological Report
5
(
5
),
648
659
.
Urbieta
,
M. S.
,
Donati
,
E. R.
,
Chan
,
K. G.
,
Shahar
,
S.
,
Sin
,
L. L.
&
Goh
,
K. M.
2016
Thermophiles in the genomic era: biodiversity, science, and applications
.
Biotechnology Advances
33
(
6
),
633
647
.
Valeriani
,
F.
,
Biagini
,
T.
,
Giampaoli
,
S.
,
Crognale
,
S.
,
Santoni
,
D.
&
Romano Spica
,
V.
2016a
Draft genome sequence of Tepidimonas taiwanensis strain VT154–175
.
Genome Announcements
4
(
5
),
e00942-16
.
Valeriani
,
F.
,
Protano
,
C.
,
Gianfranceschi
,
G.
,
Cozza
,
P.
,
Campanella
,
V.
,
Liguori
,
G.
,
Vitali
,
M.
,
Divizia
,
M.
&
Romano Spica
,
V.
2016b
Infection control in healthcare settings: perspectives for mfDNA analysis in monitoring sanitation procedures
.
BMC Infectious Diseases
16
,
394
.
Van Gemerden
,
H.
1993
Microbial mats: a joint venture
.
Marine Geology
113
,
3
25
.
Van Tubergen
,
A.
&
van Der Linden
,
S.
2002
A brief history of spa therapy
.
Annals of the Rheumatic Diseases
61
(
3
),
273
275
.
Wang
,
S.
,
Hou
,
W.
,
Dong
,
H.
,
Jiang
,
H.
,
Huang
,
L.
,
Wu
,
G.
,
Zhang
,
C.
,
Song
,
Z.
,
Zhang
,
Y.
,
Ren
,
H.
,
Zhang
,
J.
&
Zhang
,
L.
2013
Control of temperature on microbial community structure in hot springs of the Tibetan plateau
.
PLoS ONE
8
,
e62901
.
Ward
,
D. M.
,
Ferris
,
M. J.
,
Nold
,
S. C.
&
Bateson
,
M. M.
1998
A natural view of microbial biodiversity within hot spring cyanobacterial mat communities
.
Microbiology and Molecular Biology Review
62
(
4
),
1353
1370
.
Weiss
,
M. C.
,
Sousa
,
F. L.
,
Mrnjavac
,
N.
,
Neukirchen
,
S.
,
Roettger
,
M.
,
Nelson-Sathi
,
S.
&
Martin
,
W. F.
2016
The physiology and habitat of the last universal common ancestor
.
Nature Microbiology
1
(
9
),
16116
.
Wen
,
C.
,
Wu
,
L.
,
Qin
,
Y.
,
Van Nostrand
,
J. D.
,
Ning
,
D.
,
Sun
,
B.
&
Zhou
,
J.
2017
Evaluation of the reproducibility of amplicon sequencing with Illumina MiSeq platform
.
PLoS ONE
12
(
4
),
e0176716
.
Yazdi
,
M.
,
Taheri
,
M.
&
Navi
,
P.
2015
Environmental geochemistry and sources of natural arsenic in the Kharaqan hot springs, Qazvin, Iran
.
Environmental Earth Sciences
73
,
5395
5404
.

Supplementary data