Relationship between dental and periodontal health status and the salivary microbiome: bacterial diversity, co-occurrence networks and predictive models

1 Vol.:(0123456789) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports Relationship between dental and periodontal health status and the salivary microbiome: bacterial diversity, co‑occurrence networks and predictive models M. Relvas1, A. Regueira‑Iglesias2, C. Balsa‑Castro2, F. Salazar1, J. J. Pacheco1, C. Cabral1, C. Henriques1 & I. Tomás2* The present study used 16S rRNA gene amplicon sequencing to assess the impact on salivary microbiome of different grades of dental and periodontal disease and the combination of both (hereinafter referred to as oral disease), in terms of bacterial diversity, co‑occurrence network patterns and predictive models. Our scale of overall oral health was used to produce a convenience sample of 81 patients from 270 who were initially recruited. Saliva samples were collected from each participant. Sequencing was performed in Illumina MiSeq with 2 × 300 bp reads, while the raw reads were processed according to the Mothur pipeline. The statistical analysis of the 16S rDNA sequencing data at the species level was conducted using the phyloseq, DESeq2, Microbiome, SpiecEasi, igraph, MixOmics packages. The simultaneous presence of dental and periodontal pathology has a potentiating effect on the richness and diversity of the salivary microbiota. The structure of the bacterial community in oral health differs from that present in dental, periodontal or oral disease, especially in high grades. Supragingival dental parameters influence the microbiota’s abundance more than subgingival periodontal parameters, with the former making a greater contribution to the impact that oral health has on the salivary microbiome. The possible keystone OTUs are different in the oral health and disease, and even these vary between dental and periodontal disease: half of them belongs to the core microbiome and are independent of the abundance parameters. The salivary microbiome, involving a considerable number of OTUs, shows an excellent discriminatory potential for distinguishing different grades of dental, periodontal or oral disease; considering the number of predictive OTUs, the best model is that which predicts the combined dental and periodontal status. Oral diseases are a major public health concern, having a negative effect on individuals, communities and society at large1. Recent data from the Global Burden of Disease Study (GBD) revealed that around 3.5 billion people worldwide have suffered from an untreated dental condition, mainly in the form of caries, periodontitis, tooth loss and edentulism2. It is also well-known that oral diseases such as periodontitis have significant associations with chronic systemic conditions like cardiovascular disease and diabetes, with poor oral health having an adverse impact on the development and control of these types of disorder3–5. Several sets of diagnostic criteria have been developed over the last few decades for use in the classification of dental caries. These include: the Decayed, Missing and Filled Teeth (DMFT) index6; the Dundee Selectable Threshold Method (DSTM)7; the International Caries Detection and Assessment System (ICDAS)8; and the Caries Assessment Spectrum and Treatment (CAST) tool9. The Community Periodontal Index of Treatment Needs (CPITN)10 and the recent International Classification of Periodontal and Peri-implant Diseases and Conditions OPEN 1Institute of Research and Advanced Training in Health Sciences and Tecnologies (IINFACTS), IUCS-Cespu-Instituto Universitário de Ciencias da Saúde, Gandra, Paredes, Portugal. 2Oral Sciences Research Group, Special Needs Unit, Department of Surgery and Medical-Surgical Specialties, School of Medicine and Dentistry, Health Research Institute of Santiago (IDIS), Universidade de Santiago de Compostela, Galicia, 15872 Santiago de Compostela, Spain. *email: [email protected] 2 Vol:.(1234567890) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ are also used to categorise periodontal status, with the latter distinguishing forms of the disease based on a multidimensional staging and grading system11,12. In recent years, scales have been designed based on a combination of different dental and/or periodontal parameters. These include: the Total Dental Index (TDI)13; the Modified Total Dental Index (MTDI)14,15; the Asymptotic Dental Score (ADS)16; and the Brief Oral Health Status Examination (BOHSE)17. In 2013, our group developed a new scale of overall oral health, combining parameters relating to caries, gingivitis and periodontitis. The aim of this clinical scale, which is applicable to the adult population, is to obtain a single numerical value that synthesises the patient’s oral health status, allowing the establishment of an individualised grade of dental or periodontal health18. A symbiotic relationship between resident oral microbiota and a host is essential for the maintenance of oral health19,20. Indeed, several studies using 16S rRNA gene next-generation sequencing (NGS) have shown that local compositional changes in the bacterial communities, which are associated with dental caries and periodontitis, are far more complex than previously believed21,22. Assuming saliva may reflect the presence of oral disease23, investigations using 16S rRNA gene NGS have been conducted to identify which salivary microbiota profiles are related to oral health24, dental caries25 and periodontitis26. However, very few 16S-based studies have analysed the impact of various dental and periodontal parameters simultaneously to demonstrate which of them have the greatest repercussions for the salivary microbiota27–29. As a consequence, our study used 16S rRNA gene NGS to examine the bacterial diversity and the co-occurrence network patterns of the salivary microbiota in patients who have been classified clinically by our self-designed and previously validated scale of overall oral health, including dental and periodontal parameters. Also, another objective was to evaluate the diagnostic potential of salivary microbiome to discriminate between different clinical conditions. Materials and methods Selection of the study group. The participants in our study were recruited between January, 2016 and July, 2017 from a group of patients aged between 25 and 65 who had visited the Dental Clinic of Instituto Superior de Ciências da Saúde Norte, Cooperativa de Ensino Superior, Politécnico e Universitário (CESPU) (Gandra, Paredes, Portugal) for a dental check-up. The following exclusion criteria were applied: (1) fewer than 18 teeth; (2) the presence of removable dental prostheses or orthodontic devices; (3) antibiotic treatment or the routine use of oral antiseptics in the previous 3months; and (4) the incidence of any etiological factors (systemic diseases, medication, radiotherapy etc.) that could trigger changes in the production and/or composition of saliva30. After applying these criteria, 270 patients remained in our sample. The research was conducted based on the principles of the Declaration of Helsinki (revised in 2000) concerning studies involving human experimentation31. The research protocol was approved by the Instituto Superior de Ciências da Saúde–Norte, CESPU (registration number 35/CE-IUCS/2019), and all the participants provided their written informed consent to their involvement in the study. Application of our oral health scale. Examination of the oral cavity (excluding the third molars) was performed by a single trained examiner using a conventional dental mirror, a probe, and a CP11 manual periodontal probe calibrated at 3, 6, 8 and 11mm. Six sites per tooth were examined, taking into account the dental and periodontal variables included in our new scale. The findings are expressed in terms of the number of sites revealing specific oral conditions, and were evaluated according to tested and reproducible indices. The following four clinical parameters were utilised: (1) the number of sites with supragingival plaque, based on the O’Leary Plaque Index32; (2) the number of sites with caries and their severity (1 = affecting the enamel, 2 = affecting the enamel and dentine, and 3 = affecting the enamel, dentine and pulp); (3) the number of sites with gingival inflammation, based on the Ainamo and Bay Gingival Index33; and (4) the number of sites with periodontal pockets ≥ 4mm and their severity (average depth of the pockets). The participants’ dental and periodontal health grades (DG and PG, respectively) corresponded to the grades assigned to at least two of the three variables analysed in each of these categories. If there were differences between the grades allocated to each of the variables in a category, the parameters for "number of caries" and "number of periodontal pockets ≥ 4mm" took precedence. If the same grade was allocated to two variables in a category, but the third variable’s grade was two levels higher, the value assigned to it was one grade higher than that of the matching variables. The oral health grade (oral grade—OG) was determined by the category (dental or periodontal) with the highest ranking, enabling patients to be classified based on the score for their dental and periodontal health and the combination of both conditions (Table1). Our scale of overall oral health was used to select a convenience sample of 81 patients from the 270 initially recruited. These 81 participants were given the following OGs: 0 for 17 of them; 1 for 25; 2 for 28; and 3 for 11. In relation to the subscales, 47 patients had a PG of 0 and different DGs, and 46 had a DG of 0 and different PGs. Collection of saliva samples and the 16S rDNA gene amplicon sequencing. Unstimulated saliva samples (2–3mL) were collected from each participant using the spitting method34. None of the subjects had eaten or brushed their teeth for at least 1 h before they provided their sample. Total DNA was extracted from the saliva using a commercial kit (QIAamp DNA Mini kit; Qiagen, GmbH, Hilden, Germany) and according to the manufacturer’s instructions, albeit with minor modifications, including the addition of a lysozyme treatment (20mg/ml at 37 ºC for 30min). The isolated DNA was eluted in 200µl of distilled and apyrogenic water, and its quality and concentration were assessed using a Nanodrop spectrophotometer (ND-2000 Spectrophotometer, Wilmington, USA). DNA samples with spectrophotometer ratios (Abs 260/280) between 1.5 and 2.0 were considered to be acceptable for inclusion in the study35. We used the 3 Vol.:(0123456789) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ ZymoBIOMICS Microbial Community DNA Standard (Catalog No. D6306, Zymo Research, Irvine, CA, USA), which is a mix of genomic DNA isolated from pure cultures of eight bacterial and two fungal strains, as both the bacterial mock community and the positive control for the downstream procedures. Mock DNAs were amplified and sequenced in the same way as all the other samples used in the experiment. A polymerase chain reaction (PCR) amplification of the 16S rDNA gene was performed with the KAPA HiFi HotStart DNA Polymerase KAPA HiFi HotStart ReadyMixPCR Kit (Cat. No. KK2602, 7958935001; Kapa Biosystems, F. Hoffmann-La Roche Ltd, Basel, Switzerland). The V3–V4 hypervariable region was amplified as previously described36 using the following primers in a limited-cycle PCR: V3-V4-Forward (5′-CCT ACG GGNGGC WGC AG-3). V3-V4-Reverse (5′-GAC TAC HVGGG TAT CTA ATC C-3). A set of modified primers, V3-V4-F and V3-V4-R, were also used. This set contained a 1–3bp "heterogeneity spacer" that we designed to mitigate the issues caused by low sequence diversity amplicons. Each PCR amplification was carried out on a total volume of 10µl, which comprised 4µl of DNA, 0.2-μM from each forward and reverse primer, and a Kapa ready mix (Kapa Biosystems). The PCR conditions were modified by conducting: (1) an initial denaturation at 95°C for 3min; (2) 25 three-step cycles of 95°C for 30s, 55°C for 30s and 72°C for 30s; and (3) a final 5min extension at 72ºC. Water, up to a total volume of 50μl, was added after the first PCR step. The reactions were purified using AMPure XP beads (Beckman Coulter, Brea, CA, USA) with a 0.9X (V3-V4 amplicon) ratio, according to the manufacturer’s instructions. PCR products were eluted from the magnetic beads with 32μl of Buffer EB (Qiagen N.V, Hilden; Germany), with 30μl of the eluate transferred to a fresh 96-well plate. The primers described above contain overhangs that enable the addition of full-length Nextera adapters. Barcodes are available for multiplex sequencing in a second PCR step, which produces sequencing-ready libraries. To this end, 5μl of the first amplification was used as a template for the second PCR, with Nextera XT v2 adaptor primers added up to a final volume of 50μl. The PCR mix and thermal profile employed for the first PCR were also used for the second, but only for eight cycles. After the second PCR, 25μl of the final product was purified and normalised with the SequalPrep normalisation kit (Invitrogen, Carlsbad, CA, USA), according to the manufacturer’s protocol. Libraries were eluted in a 20-μl volume and pooled for sequencing. Final pools were quantified with a qPCR using the Kapa library quantification kit for Illumina Platforms (Kapa Biosystems) on an ABI 7900HT real-time cycler (Applied Biosystems, Foster City, CA, USA). Sequencing using v3 chemistry with a loading concentration of 18pM was performed in Illumina MiSeq (Illumina Inc., San Diego, CA, USA) with 2 × 300bp reads. In all cases, 10% of the PhIX control libraries were spiked to increase the diversity of the sequenced samples. In parallel, negative control tests of the sample-collection buffer, DNA-extraction and PCR-amplification steps were conducted routinely under the same conditions and using reagents. One such non-template control was subjected to the library preparation and then sequenced. As expected, this yielded very few reads (220 per sample), in contrast to an average of 188,713 reads/library in the sample-derived collections. The sequences obtained were deposited in the SRA database under accession number PRJNA623352. Bioinformatic processing. The raw reads were obtained in the fastq format and processed according to the Mothur pipeline proposed by Schloss etal.37, albeit with some modifications. Four samples with a very low number of raw sequences were removed (fewer than 436 per sample). The final raw-read total was 15,285,797 (mean ± standard deviation per sample = 188,713 ± 37,473, median = 189,918, maximum–minimum number of sequences in a sample = 282,936–119,882). The Needleman-Wunsch alignment algorithm was used to obtain the contigs from the raw sequences38. The quality filtering criteria applied were contigs with: (1) errors in the primers; (2) ambiguous bases; (3) more than eight homopolymers; (4) mismatches; and (5) overlap of less than 25bp. Lengths below 400bp or greater than 535bp were discarded. The unique candidate sequences were aligned against the SILVA-based bacterial reference alignment39, again using the Needleman-Wunsch algorithm38. A total of 52,544 potential chimeras were detected with the Table 1. Scale of overall oral health, involving grades of dental and periodontal health. Grade 0 Grade 1 Grade 2 Grade 3 Grades of dental health Supragingival plaque (O’Leary index) 0 1–56 57–112 > 112 Caries 0 1–4 5–8 ≥ 9 Severity of the caries (median) 0 1 2 3 Grades of periodontal health Gingival inflammation (Ainamo and Bay index) 0 1–56 57–112 > 112 Periodontal pockets ≥ 4mm 0 1–56 57–112 > 112 Severity of the pockets (mean) < 4mm 4–4.9mm 5–5.9mm ≥ 6mm 4 Vol:.(1234567890) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ VSEARCH tool40, representing 4.06% of the total unique sequences. The final total of high-quality reads was 4,269,754. The opticlust algorithm and the Matthews Correlation Coefficient Metric (MCC) were used to group the sequences into operational taxonomic units (OTUs) with a similarity of 97%41. A total of 3788 OTUs were obtained from all the samples, and the native Bayesian classifier proposed by Wang etal.42 was employed to categorise them. The Human Oral Microbiome Database, version 15.2, was used for this purpose43 and includes 1015 different consensus strains. A classification threshold of 80% was employed for the taxonomic designations at all levels. The non-classified level was assigned the term "unclassified". The mock D6306 sample was processed using the same protocol to evaluate the quality of the sequencing and the bioinformatic pipeline. First, a reference database was constructed from the 58 different 16S rRNA sequences provided by the manufacturer. These corresponded to all the variants of the 16S rRNA genes of the eight species included in the sample. The mock sample was then extracted from the pipeline and analysed separately against the reference database, producing errors of 0.0000073% in the bp from a total of 57,350 sequences with a length of 400bp. Fourteen OTUs were obtained from the mock sample, with eight of them highly abundant (between 4241 and 14,971 counts) and present in very similar numbers to the theoretical percentages indicated by the manufacturer. Only six OTUs were present in very low abundances (one or two counts), which could be due to sequencing errors. As a consequence, the protocol applied was considered to be of high quality. Statistical analysis. The statistical analysis of the 16S rDNA sequencing data at the species level was performed according to the protocol recently proposed by McMurdie and Holmes44, using implementations in R such as the phyloseq, DESeq2 and microbiome packages (versions 1.34.0, 1.30.0 and 1.12.0, respectively)45–47. An independent filter had previously excluded from the statistical analysis the OTUs with an abundance of ≤ 2 counts and a presence in ≤ 3 samples (3.89%)48, leaving a final total of 403 OTUs(the final number of high-quality reads was 4,235,218). The impact of the different oral-health scale grades (DG, PG and OG) on the salivary microbiome was investigated in relation to: (1) four indicators of alpha diversity and the structure of the bacterial community; and (2) the composition of the core microbiota and the testing of differential abundances. The individualised impact of the dental subscale was analysed in the participants with a PG of 0 while the periodontal subscale was assessed in those with a DG of 0. Four different comparative analyses of the grades for the dental and periodontal subscales and scale of overall oral health were performed: grade 0 vs 123; grade 0 vs 1; grade 0 vs 23; and grade 1 vs 23. Alpha diversity indicators and structure of the salivary bacteria community. The phyloseq and microbiome packages were used in the following ways to obtain data relating to the four alpha diversity indicators45,47: (1) taxa richness: the absolute count data ("observed") and the coverage index; the latter specifies how many of the more abundant OTUs are required to achieve a particular proportion of the occupied ecosystem (95%); and (2) diversity and evenness of the taxa present in the samples: the Shannon and Pielou indices49,50. The Student t-test was used to conduct different comparative analyses. A principal coordinate analysis (PCoA) was employed to visualise the clustering of the salivary samples in relation to their respective dental, periodontal or oral health grades. The "ape" package, version 5.351, was used to rotate the phylogenetic tree, while the weighted unifrac algorithm52 was applied in the phyloseq package to obtain the phylogenetic distance matrix. A non-parametric permutational multivariate analysis of variance (PERMANOVA)53 was used to measure the multivariate community-level differences between the groups. These analyses were performed using the vegan package, version 2.5.654. Composition of the core salivary microbiome and testing the differential abundances. We used the microbiome package47 to identify the central taxa present with a prevalence rate of 100% for each grade, or combination of grades from the dental and periodontal subscales and the oral scale. The DESeq2 package46 was used to identify the bacteria with the most significant changes in differential abundance at the species level for the different conditions. Improvements to the stability and dispersion of the counts (variance) were required before it was possible to calculate the differential abundances for the different species present in the samples being compared. To this end, we used the estimate SizeFactors function in DESeq246 to transform the stabilisation of the variance. The differential abundances were measured with the log2foldchange value, and the different conditions were compared using the Wald test with the Benjamini–Hochberg correction (Q parameter = 0.1, FDR < 10%). The differential-abundance measurements were statistically significant if the adjusted p-value was < 0.05 (−log10 adjusted p-value = 1.3). Co‑occurrence networks in the salivary microbiome. After filtering out the OTUs with > 10 counts in the total samples of each clinical group, we conducted an analysis of the co-occurrence networks for OG0, OG23, DG23 and PG23. The SparCC method was used to generate the networks55, as this has been shown to enable researchers to detect the linear relationships in a compositional dataset to a high degree of precision56. Default parameters and the SpiecEasi package (version 1.1.0) were employed to run SparCC57, and the correlation matrix we obtained was filtered using an absolute correlation score greater than or equal to 0.5. The networks were subsequently visualized using the igraph package (version 1.2.4.1)58, where each node represents one OTU and each edge the correlations between the OTU abundances. A set of measures was calculated to describe the topology of the resulting networks: the number of nodes and edges; the density; the average number of neighbours; the characteristic path length; the clustering coefficient, also known as transitivity; the centralization; the modularity; the number of subnetworks; and the number of modules59. 5 Vol.:(0123456789) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ We calculated the Betweenness Centrality (BC)60 to measure the relative importance of each taxon within the network (how influential a taxon is within a network). This determines the fraction of the shortest paths through one particular bacterial taxon to another. The BC of a taxon in a network reflects the importance of the control the taxon exerts over the interactions of other taxons in the same network60. In line with Banerjee etal.61, a combined score based on a high degree value and a high BC value was used as a threshold for defining the hub or keystone OTUs in the microbial communities. Diagnostic value of salivary microbiome for discriminating the clinical condition. We conducted a supervised classification in the form of a sparse partial least-squares–discriminant analysis (sPLS-DA) to facilitate the categorization of the different dental, periodontal or oral grades and identify the OTUs that best distinguished the three groups within each subscale and scale (grade 0 vs others, grade 1 vs others, grade 2, 3 vs others). The sPLS-DA was performed using the MixOmics package (version 1.6.3), which is dedicated to the integrative examination of “omics” data62. The balanced error rate (BER) was calculated for different types of distance (centroid, maximum and mahalanobis) to enable us to identify the recommended number of components for use in relation to the latent information in the OTU table. Receiver operating characteristic (ROC) curves were constructed with the true positivity rate (sensitivity) as a function of the false positivity rate 1-specificity, and area under the curve (AUC) values were used to distinguish between each clinical grade. The component with the highest discrimination value was chosen in the dental and periodontal subscales and the oral health scale. It should be noted that simulations with an AUC value equal to or higher than 0.70 are generally considered to be acceptable predictive models63. Results The 81 participants in the convenience sample had the following OG scores: 0 for 17 of them; 1 for 25; 2 for 28; and 3 for 11. Focusing on the group with the highest grade of oral pathology (OG23), we can differentiate two main subgroups: patients with DG0_PG23 (n = 14) and patients with PG0_DG23 (n = 20). In relation to the subscales, the results were as follows: (1) 47 patients had a PG of 0 and DGs between 0 and 3 (17 had a DG of 0, nine a DG of 1, 11 a DG of 2, and 10 a DG of 3); and (2) 46 had a DG of 0 and PGs between 0 and 3 (17 had a PG of 0, 14 a PG of 1; 14 a PG of 2; and one a PG of 3). The four patients excluded due to the low number of raw sequences obtained were: two of OG1 and two of OG2 (two of DG0, one of DG1 and one of DG2; two of PG0, one of PG1 and one of PG2). Impact of the dental and periodontal subscales and the scale of overall oral health on the sali‑ vary microbiome: alpha diversity indicators and the structure of the bacterial community. As seen in Table2, worsening dental or periodontal health revealed a trend of increasing alpha diversity in both the dental and periodontal subscales. In our comparisons, however, significant differences in the number of OTUs were only observed in the former: DG0 vs DG123 (209.82 ± 37.01 vs 240.57 ± 35.02; p = 0.009); and DG0 vs DG23 (209.82 ± 37.01 vs 243.55 ± 33.09, p = 0.006). Meanwhile, in our oral health scale comparisons, grade increments were linked to progressive increases in bacterial richness, especially in OG0 vs OG23 (no. observed OTUs = 212.52 ± 38.44 vs 240.24 ± 32.77, p = 0.015; coverage index 95% = 65.17 ± 15.59 vs 78.24 ± 20.42, p = 0.013). This was also the case for the Shannon Index values (OG0 = 3.67 ± 0.23, and OG23 = 3.83 ± 0.26, p = 0.026). The PCoA revealed a grouping of the salivary samples taken from the participants with DGs, PGs and OGs of 0. This was in contrast to the picture for the other grades, whose compositional distributions were more diverse (Fig.1). This visual observation was confirmed by the PERMANOVA test, which produced significant results for the comparison of grades 0 and 123 (dental subscale, p = 0.0009; periodontal subscale, p = 0.0229; oral scale, p = 0.0008). These findings were mainly at the expense of the contrast between grades 0 and 23 (dental subscale, p = 0.0005; periodontal subscale, p = 0.0287; oral scale, p = 0.0008). Focusing on the group with the highest grade of oral pathology (OG23), Permanova’s test revealed that the structure of the salivary microbiome was different depending on the predominance of dental pathology (PG0_DG23) or periodontal pathology (DG0_PG23) (p = 0.027). Impact of the dental and periodontal subscales and the scale of overall oral health on the salivary microbiome: composition of the core microbiome and testing differential abun‑ dance. Composition of the core microbiome. Figure2 uses Venn diagrams to portray the core microbiome’s number of OTUs for the different DGs, PGs and OGs. The core microbiome associated with the participants’ dental and periodontal health contained 57 species, representing 14.14% of the total number of OTUs and 63.06% of the total abundance. There were only nine taxa in DG0 and eight in PG0 (the specific core of grade 0), exemplifying abundances of 7.80% and 1.34%, respectively. Of these specific core species, five were common to both the dental and periodontal health conditions: Neisseria macacae, Butyrivibrio sp. HMT_455, Campylobac‑ ter concisus, Porphyromonas catoniae and Corynebacterium durum. There were 66 species in the core microbiome associated with the most severe dental disease (DG23) and 73 with the most serious periodontal disease (PG23), representing 16.37% and 18.11% of the total salivary microbiota, and 67.14% and 67.98% of the total abundance, respectively. There were only eight and 10 taxa present in DG23 and PG23 (the specific core of grade 23), exemplifying 2.54% and 2.46% of the abundance, respectively. Of these specific core species, only Porphyromonas endodontalis was common to both pathological conditions. There were 35 taxa common to both the dental and periodontal subscales, regardless of the grade (non-specific core), representing abundances of 50.83% and 52.75%, respectively. Of these non-specific core species, 25 were common to both subscales, with the most abundant being (abundance > 1%): Granulicatella adiacens, Haemophi‑ lus parainfluenzae, Leptotrichia sp., Porphyromonas pasteri, Prevotella sp., Prevotella melaninogenica, Prevotella 6 Vol:.(1234567890) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ salivae, Rothia mucilaginosa, Streptococcus sp., Streptococcus. oralis_subsp._dentisani_clade_058, Streptococcus salivarius, Veillonella sp. and Veillonella parvula. Species belonging to the salivary core microbiome are shown in the supplementary data 1. Testing differential abundance. The results for the dental health subscale revealed differential abundances for the different grades in 102 species (25.31% of the total OTUs), 39 of which were core species and 63 non-core (9.67% and 15.63% of the total OTUs, respectively). There were 36 species associated with oral health and 66 with some grade of dental pathology. In the comparison of DG0 vs DG23, there were differential abundances in the main dental-health related OTUs (with > 1% levels in DG0; log2foldchange values ranged from 5.86 to 1.41): Porphyromonas pasteri, Fusobacterium periodonticum, Veillonella parvula, Alloprevotella sp. HMT 473, Alloprevotella tannerae and Neisseria subflava. Those related to high grades of dental pathology (> 0.5% levels in DG23; log2foldchange values ranged from −1.98 to −1.35) were: Alloprevotella sp. HMT 308, Streptococcus parasanguinis clade 411, Atopobium sp., Fusobacterium nucleatum_subsp._vincentii, Megasphaera micronuci‑ formis and Alloprevotella rava (Table3). The results for the periodontal health subscale revealed differential abundances for the different grades in 27 species (6.69% of the total OTUs), eight of which were core species and 19 non-core (1.98% and 4.71% of the total OTUs, respectively). Twelve species were associated with oral health and 15 with some grade of periodontal pathology. In the comparison of PG0 vs PG23, there were differential abundances in the two main periodontalhealth related OTUs (with > 1% levels in PG0): Haemophilus parainfluenzae and Capnocytophaga leadbetteri, with log2foldchange values of 1.16 and 2.31, respectively. Those associated with high grades of periodontal pathology (with > 0.1% levels in PG23; log2foldchange values ranged from −3.07 to −2.13) were: Tannerella forsythia, Mycoplasma faucium, Fretibacterium sp., and Bacteroidetes [G-5] bacterium HMT 511 (Table4). The results for the scale of overall oral health revealed differential abundances for the different grades in 88 species (21.83% of the total OTUs), 22 of which were core species and 66 non-core (5.45% and 16.37% of the total OTUs, respectively). There were 22 species associated with oral health and 66 with some grade of oral pathology, including dental and periodontal conditions. In the comparison of OG0 vs OG23, the differential abundances in the main oral-health related OTUs (with > 1% levels in OG0; log2foldchange values ranged from Table 2. Alpha diversity indicators in the different grades of the dental and periodontal subscales, as well as in the scale of overall oral health. The following four comparisons were performed in both the dental and periodontal subscales as well as the scale of overall oral health: grade 0 vs. grade 123; grade 0 vs. grade 1; grade 0 vs. grade 23; grade 1 vs. grade 23. No. observed OTUs Coverage index (95%) Shannon index Pielou index Dental health grade (DG) DG0 209.82 ± 37.01 65.11 ± 15.64 3.66 ± 0.23 0.68 ± 0.03 DG123 240.57 ± 35.02 74.64 ± 20.54 3.79 ± 0.24 0.69 ± 0.03 DG1 233.12 ± 40.84 71.37 ± 20.97 3.79 ± 0.26 0.69 ± 0.04 DG23 243.55 ± 33.09 75.95 ± 20.76 3.79 ± 0.23 0.69 ± 0.03 Significant grade comparisons (p-value < 0.05) 0 vs. 123 – – – – – – – 0 vs. 23 – – – – – – – Periodontal health grade (PG) PG0 208.52 ± 35.75 64.88 ± 15.36 3.66 ± 0.23 0.68 ± 0.03 PG123 219.48 ± 32.32 71.51 ± 18.80 3.76 ± 0.27 0.70 ± 0.04 PG1 214.61 ± 39.62 66.07 ± 18.30 3.69 ± 0.24 0.68 ± 0.03 PG23 224.00 ± 24.38 76.57 ± 18.45 3.83 ± 0.28 0.70 ± 0.04 Significant grade comparisons (p-value < 0.05) – – – – – – – – – – – – – – – – Oral health grade (OG) OG0 212.52 ± 38.44 65.17 ± 15.59 3.67 ± 0.23 0.68 ± 0.03 OG123 234.85 ± 36.47 75.28 ± 20.39 3. 80 ± 0.26 0.69 ± 0.04 OG1 226.17 ± 41.02 70.52 ± 19.86 3.76 ± 0.26 0.69 ± 0.03 OG23 240.24 ± 32.77 78.24 ± 20.42 3.83 ± 0.26 0.70 ± 0.04 Significant grade comparisons (p-value < 0.05) 0 vs. 123 0 vs. 123 0 vs. 123 – – – – – 0 vs. 23 0 vs. 23 0 vs. 23 – – – – – 7 Vol.:(0123456789) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ 1.71 to 0.88) were: Porphyromonas pasteri, Haemophilus parainfluenzae, Veillonella parvula, and Capnocytophaga leadbetteri. Those associated with high grades of oral pathology (> 1% levels in OG23; log2foldchange values ranged from −2.86 to −1.33) were: Alloprevotella sp. HMT 308, Fusobacterium nucleatum_subsp._vincentii and Porphyromonas gingivalis (Table5). The relative abundances of total OTUs in dental, periodontal, oral health grades, and analyses of differential abundance between grades are detailed in the supplementary datas 2–4. Impact of the dental and periodontal subscales and the scale of overall oral health on the sali‑ vary microbiome: co‑occurrence network patterns and discriminatory potential of the salivary microbiome. Co‑occurrence network patterns. Table6 shows the topological parameters of the co-occur- rence networks in the different clinical groups. In general terms, the most important differences were detected between OG0 and OG23 followed by OG0 and DG23. Figure1. Principal Coordinate Analysis (PCoA), including PERMANOVA test values in the comparison between different grades of dental, periodontal and oral health (graphic made using the phyloseq package, version 1.34.0). 8 Vol:.(1234567890) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ The network coverage in OG0 was higher than in OG23 (65.75% and 43.42%, respectively), as was the number of edges (1867 and 436, respectively). There was also a better balance between the number of positive and negative correlations in OG0 with respect to OG23 (56.5% and 43.5% and 83.1% and 16.9%, respectively). The network in OG0 had a higher density, average number of neighbours, and higher centralization values than the network in OG23 (0.05, 15.09 and 0.16 vs 0.02, 5.98 and 0.08); in contrast, the characteristic path lengths and modularity scores were lower (2.86 and 0.36 in OG0 vs 4.27 and 0.62 in OG23). The OG0 network had fewer subnetworks and a higher number of modules (two and 28 vs five and 19 in OG23). All of these differences (except for the modularity values) were also observed between the OG0 and DG23 networks, albeit to a lesser extent. In contrast, all of the other parameters had similar values in the comparison between OG0 and PG23. In the OG0 network, the three main hubs or keystone OTUs, based on their combined scores in the three main modules, were Porphyromonas pasteri, Porphyromonas endodontalis and Prevotella salivae. Although all three species were part of the core microbiome, only Porphyromonas pasteri was present in copious numbers (relative abundance of 8.07%) and differentially abundant in relation to the oral health status. In the OG23 network, the main hubs or keystone taxa found were Fusobacterium periodonticum, Treponema socranskii and Prevotella sp. HTM 305. Only Fusobacterium periodonticum was abundant (relative abundance of 2.95%) in the core microbiome, but none of the three species were differentially abundant. In the DG23 network, the main hubs or keystone taxa detected were Tannerella forsythia, Fusobacterium nucleatum_subsp. vincentii and Prevotella oris. Only Fusobacterium nucleatum_subsp. vincentii was part of the core microbiome and differentially abundant in relation to the dental disease status. In the PG23 network, the main hubs or keystone taxa detected were Granu‑ licatella adiacens, Porphyromonas endodontalis and Campylobacter gracilis. Of these three species, Granulicatella adiacens and Porphyromonas endodontalis were abundant (relative abundances of 1.92% and 1.45%, respectively) and belonged to the core microbiome, but none of the three taxa were differentially abundant (Figs.3, 4, 5, 6). Although there was a predominance of positive correlations in the keystone OTUs of the OG23 network, there were similar numbers of positive and negative associations in the keystone OTUs of the DG23 network. Discriminatory potential of the salivary microbiome. With regard to the diagnostic value of salivary microbiome, the model obtained to distinguish the grades of the dental subscale was composed of 60 OTUs and the Figure2. Number of OTUs of the salivary core microbiome present in different grades of dental, periodontal and oral health, as well as their respective relative abundance values. 9 Vol.:(0123456789) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ Table 3. Species that presented differential abundance among different grades of dental health, and exhibited an abundance > 0.1% (green square, in favour of health; red square, in favour of pathology). The species that belonged to the core are colored orange. No OF OTUS SEICEPS SUNEG DENTAL HEALTH GRADES ABUNDANCE (%) DIFFERENTIAL ABUNDANCE 0 1231 23 0 vs 123 0 vs 1 0 vs 23 1 vs 23 Otu0002 Porphyromonas pasteri 7.95 4.19 5.47 3.67 Otu0006 Fusobacterium periodoncum 6.10 2.76 4.22 2.17 aluvrap allenollieV 7100utO 2.11 0.80 1.20 0.65 Otu0018 Alloprevotella sp._HMT_473 2.27 0.74 0.56 0.81 Otu0019 Alloprevotella sp._HMT_308 0.54 1.61 1.16 1.79 Otu0021 Alloprevotella tannerae 1.25 0.83 1.41 0.60 Otu0022 Streptococcus parasanguinis_clade_411 0.50 1.67 1.29 1.83 92.0 silarorev alletoverP 8200utO 1.37 1.14 1.46 deifissalcnu muibopotA 0300utO 0.36 0.98 0.55 1.16 Otu0032 Fusobacterium nucleatum_subsp._vinceni 0.35 1.06 1.39 0.92 deifissalcnu muiretcabirO 7300utO 0.65 0.50 0.93 0.33 Otu0045 Campylobacter rectus 0.75 0.33 0.26 0.35 deifissalcnu allemeG 6400utO 0.38 0.77 0.09 1.04 00.0 81.050.039.0 acaitnarua alletoverP 9400utO 90.0 20.070.084.1avalfbus airessieN 1500utO Otu0053 Stomatobaculum sp._HMT_097 0.38 0.22 0.37 0.15 47.0 40.045.061.0 deifissalcnu aihtaenS 4500utO Otu0055 Capnocytophaga spugena 0.55 0.18 0.26 0.14 30.0 44.051.056.0 iihahs alletoverP 7500utO Otu0058 Absconditabacteria_(SR1)_[G-1] bacterium_HMT_875 0.85 0.11 0.12 0.11 223_TMH_.ps alleyegreB 9500utO 0.36 0.21 0.32 0.16 01.0 91.031.043.0 muvrap muiretcabirO 3600utO Otu0065 Megasphaera micronuciformis 0.28 0.96 0.80 1.03 60.0mugnol mulucabotamotS 0700utO 0.34 0.34 0.34 ibrom allenotaC 2700utO 0.37 0.19 0.26 0.16 43.0 10.052.053.0 602_TMH_.ps alleyegreB 6700utO snecsergin alletoverP 0800utO 0.11 0.41 0.44 0.39 Otu0085 Capnocytophaga unclassified 0.60 0.24 0.35 0.19 Otu0088 Alloprevotella rava 0.23 0.56 0.25 0.68 Otu0089 Lachnoanaerobaculum unclassified 0.20 0.10 0.16 0.07 Otu0090 Absconditabacteria_(SR1)_[G-1] bacterium_HMT_345 0.26 0.09 0.13 0.07 Otu0102 Lachnospiraceae_[G-3] [G-3]_bacterium_HMT_100 0.29 0.06 0.16 0.03 Otu0105 Saccharibacteria_(TM7)_[G-1] unclassified 0.09 0.28 0.23 0.30 deifissalcnu sulihpomeaH 6010utO 0.52 0.21 0.06 0.26 Otu0110 Lachnospiraceae_[G-2] bacterium_HMT_096 0.17 0.04 0.03 0.04 Otu0113 Absconditabacteria_(SR1)_[G-1] bacterium_HMT_8740.170.070.06 0.08 41.0 80.021.020.0muicuaf amsalpocyM 4110utO 63.053.070.0 alocitned alletoverP 6110utO 0.34 10.0 iedaw aihcirtotpeL 0210utO 0.18 0.27 0.15 Otu0124 Saccharibacteria_(TM7)_[G-1] bacterium_HMT_346 0.06 0.16 0.17 0.16 50.0 silicarg retcabolypmaC 2310utO 0.18 0.32 0.12 Otu0143 Bacteroidetes_[G-5] bacterium_HMT_511 0.01 0.10 0.23 0.05 50.0 40.050.011.0airea aihtoR 6410utO 11.0 61.031.030.0 siditiruelp alletoverP 8510utO 61.0 40.031.010.0 iihahs aihcirtotpeL 2810utO 30.0 40.030.061.0032_TMH_.ps amenoperT 6810utO Otu0211 Campylobacter concisus 0.11 0.03 0.02 0.03 20.0 41.050.010.0 729_TMH_.ps amenoperT 9320utO Otu0244 Peptostreptococcaceae_[XI][G-7] unclassified0.120.010.02 0.00 41.0 70.021.000.0 719_TMH_.ps allenollieV 6420utO 40.0 02.080.010.0 xarovirahccasitlum alletoverP 7420utO Otu0251 Capnocytophaga sp._HMT_336 0.03 0.09 0.22 0.04 10.0 22.070.000.0 muluvrap muibopotA 6620utO 30.0 71.070.000.0 614_TMH_.ps muibopotA 9620utO 16 Vol:.(1234567890) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ periodontal condition (healthy or pathology). Interestingly, more than half of the non-specific core species in the present series were the same as those previously identified by Takeshita in ≥ 75% of Japanese adults28; some of the most predominant species with an abundance ≥ 1% were: Granulicatella adiacens, Haemophilus parain‑ fluenzae, Prevotella sp. (P. melaninogenica), Rothia mucilaginosa, Streptococcus sp. (S. salivarius), and Veillonella sp. (V. parvula). These results confirm that several bacterial taxa in the salivary microbiome could be present in ethnically diverse populations, regardless of their dental and periodontal health statuses. In terms of our differential abundance findings, the number of taxa present in the salivary microbiota at significantly different levels for the subscale or overall scale grades did not exceed 25%, and mainly included non-core species (5–16%); for the oral health scale, the figure was 22%. This number of taxa was higher for the dental than the periodontal subscale (25% vs 7% of the total OTUs). This is the first time that this finding has been reported in the literature. If it is assumed that the salivary microbiome comprises a mix of bacterial communities originating from various sites in the oral cavity67, our observations provided evidence that the relative abundances of the most predominant bacteria in saliva are not strongly related to the grade of overall oral health. Figure6. Main module of the co-occurrence network associated with grade 0 of the dental health subscale (node = 101; degree = 2663) (graphic made using the igraph package, version 1.2.4.1). 17 Vol.:(0123456789) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ Figure7. Potential of the salivary microbiome to discriminate the different grades of dental and periodontal health and the combination of both: ROC curves and AUC values (graphics made using the MixOmics package, version 1.6.3). 18 Vol:.(1234567890) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ This corroborates the view that the dominant source of the salivary microbiome is probably the bacterial communities on the mucosal surfaces, especially those from the dorsum of the tongue67,84,85. When we compared our two subscales, we detected that parameters like supragingival plaque levels and the number and severity of caries had a greater differential abundance impact on the salivary microbiota than periodontal parameters such as levels of gingival inflammation and the number and severity of the periodontal pockets. Consequently, it is clear that supragingival dental conditions affect the abundance of salivary microbiota more than subgingival periodontal conditions, with the former having a greater influence on oral health. Nevertheless, in the present series, specific bacteria were involved in different oral conditions. In our view, for a species to be evaluated as a possible bacterial biomarker associated with a given oral condition, it must be: predominant in most individuals (core species), especially in those with the oral condition with which there is an intended association; and present in high levels of abundance. Applying these criteria to our 16S-based data, the main bacteria associated with dental health (2–3-fold increment in relative abundance in DG0) were: Porphyromonas pasteri, Fusobacterium periodonticum and Veillonella parvula; for periodontal health (1.69-fold increase in the relative abundance in PG0), the species was Haemophilus parainfluenzae. All of these bacteria continued to be associated with overall oral health (1.5–2-fold increase in relative abundance in OG0). The exception was Fusobacterium periodonticum, probably because it has been found to be more abundant in patients with periodontitis27. These outcomes are in accordance with the findings of most 16S-based studies concerning the identification of the above-mentioned taxa as core species present in > 70% of people28,68,72,79. Other authors have, like us, also detected differential abundances of these species, which provide support for dental or periodontal health27,68–70,76. In our study, the main bacteria associated with high grades of dental pathology (2.5–3.5-fold increment in relative abundance in DG23) were: Streptococcus parasanguinis clade 411, Alloprevotella sp. HMT 308, Atopo‑ bium sp., and Fusobacterium nucleatum_subsp._vincentii. Although most of these bacteria have been described as core species in other studies28, there is no link to dental pathology in relation to differential abundance. In the periodontal subscale, Tannerella forsythia was the only species associated with high periodontal pathology grades (6-fold increment in relative abundance in PG23). This bacterium has also been identified in the core microbiome in our previous studies27,74. However, although it was abundant in the subgingival samples taken from participants with periodontitis75, these levels in the salivary microbiota were not high (< 0.3%) compared to those of the species mentioned previously. With respect to the other two species that compose the red complex, Porphyromonas gingivalis and Treponema denticola, although they did not present differential abundance, these taxa were ten and six times more abundant, respectively, in the highest grades of periodontal pathology than in periodontal health. Only Fusobacterium nucleatum_subsp._vincentii continued to be the taxon associated with the highest levels of oral pathology, with a 3-fold increase in OG23. Interestingly, although not a core species, Porphyromonas gingivalis was nine times more abundant for the highest oral pathology grades, revealing its role not only in periodontal but also in dental pathology, as reported in the literature86,87. Impact of the dental and periodontal subscales and the scale of overall oral health on the salivary microbiome: co‑occurrence network patterns and discriminatory potential of the sali‑ vary microbiome. The analysis of networks representing microbe-microbe interactions has increased our knowledge of how microbes potentially interact in their environment88. Despite the importance of this type of investigation, our review found that co-occurrence results were reported in only a few papers that evaluated the salivary microbiota of patients with and without caries, the periodontally healthy and those with periodontitis, and adults with different oral health conditions25,69. It is important to notethat any comparison of the co-occurrence results in the different studies in the literature should be undertaken with caution: the findings obtained may well be affected by methodological differences such as the different correlation values employed as cut-off points89 or the use of different keystone OTU definitions61. The positive and negative correlations detected in a co-occurrence network can describe the tendency of different species to co-occur in a variety of circumstances. Consequently, two species exhibiting a significantly positive correlation in relation to abundance could be a direct indication of a shared preference for a particular combination of environmental conditions; alternatively, it might be a true ecological interaction where, for example, two species can grow better due to metabolite exchanges. In contrast, a negative correlation between species may indicate a competitive interaction for nutrients or differences in physiological requirements to the extent that they never occupy the same niche90. As co-occurrence networks are undirected and weighted, there is no directionality of the interactions to the degree that the underlying cause of such patterns cannot be explained90. In our study, like the findings of previous authors25,69, the topological characteristics of the salivary microbiome networks differed between the statuses of oral health and the presence of high grades of oral or dental disease. The health-associated network presented as a bacterial community with more interconnections between its members and a greater balance between the coactive and competitive interactions. The disease-associated network was less dense and synergistic exchanges predominated, suggesting that the antagonism between the oral bacteria was not a major driving force in the formation of the disease-associated microbial community26. However, in contrast, high grades of periodontal disease did not affect the characteristics of the salivary microbiome co-occurrence network. One of the most useful features of a co-occurrence network analysis is that hubs or keystone OTUs, which are highly associated taxa in a microbiome, can be identified91. However, in the salivary microbiome studies in the literature that included co-occurrence network analyses, attempts to detect keystone OTUs were uncommon26. In the present study, seven of the 12 keystone OTUs identified (58.3%) were part of the salivary core microbiome. In this sense, it has been suggested that the contribution of keystone taxa will be greater if they are part 19 Vol.:(0123456789) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ of the core microbiome and consistently present, highlighting the importance of such taxa for microbiome functioning92. On the other hand, eight and 10 of the 12 keystone OTUs identified (66.6% and 83.3%, respectively) were low-abundance species (< 1%) with no associated differential abundance results. This confirms that keystone OTUs could have an impact on microbiome functioning, irrespective of the abundance parameters61. Our results revealed that the keystone OTUS identified in the main modules of the networks varied between the clinical groups. This suggests that the possible drivers of the salivary microbiome’s structure and functioning differ in relation to oral health and disease, and even vary in different types of disease. The exception was the species Porphyromonas endodontalis, which was identified as a keystone OTU in the OGO and PG23 networks, making it particularly important. Given our findings, and despite the well-known limitations, we believe it is vital to conduct identifications at the species level as part of a co-occurrence network analysis. This is because we identified species in the present networks that were associated with different clinical conditions but belong to the same genus: Porphyromonas in OG0 and PG23, Prevotella in OG0, OG23 and DG23, and Fusobacterium in OG23 and DG23. However, evidence derived from further experimentation is required before hub taxa in inferred networks can be classified as keystone93,94. Our review of the relevant literature identified only a few studies that attempted to evaluate the diagnostic accuracy of the salivary microbiome for diagnosing oral diseases like periodontitis72,76. Moreover, to the best of our knowledge, our research is the first to assess the potential of the salivary microbiome for distinguishing different grades of dental and periodontal disease, or a combination thereof. The results we obtained revealed the excellent discriminatory power (AUC values > 0.88) of this approach. In line with the differential abundance results, a higher number of predictive variables (140 OTUs) was required to discriminate the periodontal condition than the dental condition (60 OTUs). Interestingly, the best model of the salivary microbiome, i.e., the one associated with the lowest number of predictor variables (30 OTUs), was the model for a combination of dental and periodontal disease, evidencing the impact of both of these conditions on saliva. Focusing on the group with the highest grade of oral pathology (OG23), the heterogeneity of this group (based mainly on the presence of two subgroups, DG0_PG23 and PG0_DG23) did not affect the discriminatory power of the salivary microbiome. The main limitation of this study is clinical, due to the difficulty in recruiting patients with the highest level of periodontal pathology. This meant that grade 2 and grade 3 patients had to be placed in the same group, despite our preference to conduct an analysis for each grade individually. On the other hand, we are aware that relative abundance measurements derived from 16S rRNA gene amplicon sequencing do not accurately reflect absolute taxon concentrations. This would require, at least, one broad-range qPCR determination for the calculation of inferred bacterial concentrations, which are a reasonable proxy of species-specific qPCR values95. Although this was not an objective in the present study, another interesting avenue for future research would be to evaluate the impact that smoking has on the salivary microbiota of patients with different grades of dental, periodontal or oral health. In conclusion, the simultaneous presence of dental and periodontal pathology has a potentiating impact on the richness and diversity of the salivary microbiota. The structure of the bacterial community in oral health differs from that present in dental, periodontal or oral disease, especially in high grades. The non-specific microbiome core comprises a greater number of species present in higher abundance than the specific core of a particular dental or periodontal condition (health or pathology). The number of taxa in the salivary microbiota with differential abundances between the DGs, PGs or OGs represents, at most, a quarter of the bacterial community and are mainly non-core species. Supragingival dental parameters condition the abundance of the microbiota more than subgingival periodontal parameters, with the former contributing more to the impact that oral health has on the salivary microbiome. The oral health-associated network has a bacterial community with more interconnections between its members, and a greater balance between its synergistic and competitive interactions, than the network associated with high grades of oral or dental disease; meanwhile, high grades of periodontal disease do not condition the characteristics of the salivary microbiome co-occurrence network. The possible keystone OTUs are different in relation to oral health and disease, and even these vary in different types of disease: half of them belongs to the core microbiome and are independent of the abundance parameters. The salivary microbiome, involving a considerable number of OTUs, shows an excellent discriminatory capability for distinguishing different grades of dental, periodontal or oral disease; considering the number of predictive OTUs, the worst model is that which predicts the periodontal status and the best model is that which predicts the combined dental and periodontal status. Received: 25 July 2020; Accepted: 8 December 2020 References 1. Peres, M. A. et al. Oral diseases: A global public health challenge. Lancet 394, 249–260. https ://doi.org/10.1016/S0140 -6736(19)31146 -8 (2019). 2. Kassebaum, N. J. et al. Global, regional, and national prevalence, incidence, and disability-adjusted life years for oral conditions for 195 countries, 1990–2015: A systematic analysis for the global burden of diseases, injuries, and risk factors. J. Dent. Res. 96, 380–387. https ://doi.org/10.1177/00220 34517 69356 6 (2017). 3. Carrizales-Sepúlveda, E. F., Ordaz-Farías, A., Vera-Pineda, R. & Flores-Ramírez, R. Periodontal disease, systemic inflammation and the risk of cardiovascular disease. Heart Lung Circ. 27, 1327–1334. https ://doi.org/10.1016/j.hlc.2018.05.102 (2018). 4. Preshaw, P. M. & Bissett, S. M. Periodontitis and diabetes. Br. Dent. J. 227, 577–584. https ://doi.org/10.1038/s4141 5-019-0794-5 (2019). 5. Falcao, A. & Bullón, P. A review of the influence of periodontal treatment in systemic diseases. Periodontol 2000. 79, 117–128. https ://doi.org/10.1111/prd.12249 (2019). 20 Vol:.(1234567890) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ 6. Klein, H., Palmer, C.E. & Knutson, J.W. Studies on dental caries. I. Dental status and dental needs of elementary school children.Pub. Health Rep.53, 751–765(1938). 7. Fyffe, H. E., Deery, C., Nugent, Z. J., Nuttall, N. M. & Pitts, N. B. Effect of diagnostic threshold on the validity and reliability of epidemiological caries diagnosis using the Dundee selectable threshold method for caries diagnosis (DSTM). Community Dent. Oral Epidemiol. 28, 42–51. https ://doi.org/10.1034/j.1600-0528.2000.28010 6.x (2000). 8. Ismail, A. I. et al. The international caries detection and assessment system (ICDAS): An integrated system for measuring dental caries. Community Dent. Oral Epidemiol. 35, 170–178. https ://doi.org/10.1111/j.1600-0528.2007.00347 .x (2007). 9. Leal, S. C., Ribeiro, A. P. D. & Frencken, J. E. Caries assessment spectrum and treatment (CAST): A novel epidemiological instrument. Caries Res. 51, 500–506. https ://doi.org/10.1159/00047 9042 (2017). 10. Ainamo, J. et al. Development of the World Health Organization (WHO) community periodontal index of treatment needs (CPITN). Int. Dent. J. 32, 281–291 (1982). 11. Caton, J. G. et al. A new classification scheme for periodontal and peri-implant diseases and conditions—Introduction and key changes from the 1999 classification. J. Clin. Periodontol. 45, S1–S8. https ://doi.org/10.1111/jcpe.12935 (2018). 12. Tonetti, M. S., Greenwell, H. & Kornman, K. S. Staging and grading of periodontitis: Framework and proposal of a new classification and case definition. J. Periodontol. 89, S159–S172. https ://doi.org/10.1002/JPER.18-0006 (2018). 13. Mattila, K. J. et al. Association between dental health and acute myocardial infarction. BMJ 298, 779–781. htt ps ://doi.org/10.1136/ bmj.298.6676.779 (1989). 14. Beck, J., Garcia, R., Heiss, G., Vokonas, P. S. & Offenbacher, S. Periodontal disease and cardiovascular disease. J. Periodontol. 67(Suppl 10S), 1123–1137. https ://doi.org/10.1902/jop.1996.67.10s.1123 (1996). 15. Mattila, K. J. Dental infections as a risk factor for acute myocardial infarction. Eur. Heart J. 14(Suppl K), 51–53 (1993). 16. Janket, S. J. et al. Asymptotic dental score and prevalent coronary heart disease. Circulation 109, 1095–1100. https ://doi. org/10.1161/01.CIR.00001 18497 .44961 .1E (2004). 17. Kayser-Jones, J., Bird, W. F., Paul, S. M., Long, L. & Schell, E. S. An instrument to assess the oral health status of nursing home residents. Gerontologist 35, 814–824. https ://doi.org/10.1093/geron t/35.6.814 (1995). 18. Relvas, M., Diz, P., Seoane, J. & Tomás, I. Oral health scales: Design of an oral health scale of infectious potential. Med. Oral Patol. Oral Cir. Bucal 18, e664–e670. https ://doi.org/10.4317/medor al.18427 (2013). 19. Lynge Pedersen, A. M. & Belstrøm, D. The role of natural salivary defences in maintaining a healthy oral microbiota. J. Dent. 80, S3–S12. https ://doi.org/10.1016/j.jdent .2018.08.010 (2019). 20. Marsh, P. D., Head, D. A. & Devine, D. A. Ecological approaches to oral biofilms: Control without killing. Caries Res. 49(suppl 1), 46–54. https ://doi.org/10.1159/00037 7732 (2015). 21. Simon-Soro, A., Belda-Ferre, P., Cabrera-Rubio, R., Alcaraz, L. D. & Mira, A. A tissue-dependent hypothesis of dental caries. Caries Res. 47, 591–600. https ://doi.org/10.1159/00035 1663 (2013). 22. Camelo-Castillo, A. J. et al. Subgingival microbiota in health compared to periodontitis and the influence of smoking. Front. Microbiol. 6, 119. https ://doi.org/10.3389/fmicb .2015.00119 (2015). 23. Kaczor-Urbanowicz, K. E. et al. Saliva diagnostics—Current views and directions. Exp. Biol. Med. 242, 459–472. https ://doi. org/10.1177/15353 70216 68155 0 (2017). 24. Zaura, E. et al. On the ecosystemic network of saliva in healthy young adults. ISME J. 11, 1218–1231. https ://doi.org/10.1038/i smej .2016.199 (2017). 25. Zhou, J. et al. Exploration of human salivary microbiomes–insights into the novel characteristics of microbial community structure in caries and caries-free subjects. PLoS ONE 11, e0147039. https ://doi.org/10.1371/journ al.pone.01470 39 (2016). 26. Chen, C. et al. Oral microbiota of periodontal health and disease and their changes after nonsurgical periodontal therapy. ISME J. 12, 1210–1224. https ://doi.org/10.1038/s4139 6-017-0037-1 (2018). 27. Belstrom, D., Paster, B. J., Fiehn, N. E., Bardow, A. & Holmstrup, P. Salivary bacterial fingerprints of established oral disease revealed by the human oral microbe identification using next generation sequencing (HOMINGS) technique. J. Oral Microbiol. 8, 30170. https ://doi.org/10.3402/jom.v8.30170 (2016). 28. Takeshita, T. et al. Bacterial diversity in saliva and oral health-related conditions: The Hisayama study. Sci. Rep. 6, 22164. https :// doi.org/10.1038/srep2 2164 (2016). 29. Chen, T. et al. High-throughput sequencing analyses of oral microbial diversity in healthy people and patients with dental caries and periodontal disease. Mol. Med. Rep. 16, 127–132. https ://doi.org/10.3892/mmr.2017.6593 (2017). 30. Tschoppe, P., Wolgin, M., Pischon, N. & Kielbassa, A. M. Etiologic factors of hyposalivation and consequences for oral health. Quintessence Int. 41, 321–333 (2010). 31. World Medical Association. World medical association declaration of Helsinki: Ethical principles for medical research involving human subjects. JAMA 310, 2191–2194. https ://doi.org/10.1001/jama.2013.28105 3 (2013). 32. O’Leary, T. J., Drake, R. B. & Naylor, J. E. The plaque control record. J. Periodontol. 43, 38. h ttps ://doi.org/10.1902/jop.1972.43.1.38 (1972). 33. Ainamo, J. & Bay, I. Problems and proposals for recording gingivitis and plaque. Int. Dent. J. 25, 229–235 (1975). 34. Navazesh, M. & Christensen, C. M. A comparison of whole mouth resting and stimulated salivary measurement procedures. J. Dent. Res. 61, 1158–1162. https ://doi.o rg/10.1177/00220 34582 06101 00901 (1982). 35. Mesquita, B. et al. Helicobacter pylori identification: A diagnostic/confirmatory method for evaluation. Curr. Microbiol. 69, 245–251. https ://doi.org/10.1007/s0028 4-014-0578-8 (2014). 36. Willis, J. R. et al. Citizen science charts two major “stomatotypes” in the oral microbiome of adolescents and reveals links with habits and drinking water composition. Microbiome 6, 218. https ://doi.org/10.1186/s4016 8-018-0592-3 (2018). 37. Schloss, P. D. et al. Introducing mothur: Open-source, platform-independent, community-supported software for describing and comparing microbial communities. Appl. Environ. Microbiol. 75, 7537–7541. https ://doi.org/10.1128/AEM.01541 -09 (2009). 38. Needleman, S. B. & Wunsch, C. D. A general method applicable to the search for similarities in the amino acid sequence of two proteins. J. Mol. Biol. 48, 443–453. https ://doi.org/10.1016/0022-2836(70)90057 -4 (1970). 39. Pruesse, E. et al. SILVA: A comprehensive online resource for quality checked and aligned ribosomal RNA sequence data compatible with ARB. Nucleic Acids Res. 35, 7188–7196. https ://doi.org/10.1093/nar/gkm86 4 (2007). 40. Rognes, T., Flouri, T., Nichols, B., Quince, C. & Mahe, F. VSEARCH: A versatile open source tool for metagenomics. PeerJ 4, e2584. https ://doi.org/10.7717/peerj .2584 (2016). 41. Westcott, S. L. & Schloss, P. D. OptiClust, an improved method for assigning amplicon-based sequence data to operational taxonomic units. mSphere 2, e00073-17. https ://doi.org/10.1128/mSphe reDir ect.00073 (2017). 42. Wang, Q., Garrity, G. M., Tiedje, J. M. & Cole, J. R. Naive Bayesian classifier for rapid assignment of rRNA sequences into the new bacterial taxonomy. Appl. Environ. Microbiol. 73, 5261–5267. https ://doi.org/10.1128/AEM.00062 -07 (2007). 43. Chen, T. et al. The human oral microbiome database: A web accessible resource for investigating oral microbe taxonomic and genomic information. Database 2010, baq013. https ://doi.org/10.1093/datab ase/baq01 3 (2010). 44. McMurdie, P. J. & Holmes, S. Waste not, want not: Why rarefying microbiome data is inadmissible. PLoS Comput. Biol. 10, e1003531. https ://doi.org/10.1371/journ al.pcbi.10035 31 (2014). 45. McMurdie, P. J. & Holmes, S. phyloseq: An R package for reproducible interactive analysis and graphics of microbiome census data. PLoS ONE 8, e61217. https ://doi.org/10.1371/journ al.pone.00612 17 (2013). 21 Vol.:(0123456789) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ 46. Love, M. I., Huber, W. & Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550. https ://doi.org/10.1186/s1305 9-014-0550-8 (2014). 47. Lahti, L. & Shetty, S. Tools for microbiome analysis in R. Microbiome package version1.9.15.ht t p://micro biom e .gi thu b .co m/micro biome (2017). 48. Bourgon, R., Gentleman, R. & Huber, W. Independent filtering increases detection power for high-throughput experiments. Proc. Natl. Acad. Sci. USA. 107, 9546–9551. https ://doi.org/10.1073/pnas.09140 05107 (2010). 49. Spellerberg, I. F. & Fedor, P. J. A tribute to Claude Shannon (1916–2001) and a plea for more rigorous use of species richness, species diversity and the ‘Shannon–Wiener’ Index. Global Ecol. Biogeogr. 12, 177–179. https ://doi.org/10.1046/j.1466-822X.2003.00015 .x (2003). 50. Pielou, E. C. The measurement of diversity in different types of biological collections. J. Theor. Biol. 13, 131–144. https ://doi. org/10.1016/0022-5193(66)90013 -0 (1966). 51. Paradis, E. & Schliep, K. ape 5.0: An environment for modern phylogenetics and evolutionary analyses in R. Bioinformatics 35, 526–528. https ://doi.org/10.1093/bioin forma tics/bty63 3 (2019). 52. Lozupone, C. A., Hamady, M., Kelley, S. T. & Knight, R. Quantitative and qualitative beta diversity measures lead to different insights into factors that structure microbial communities. Appl. Environ. Microbiol. 73, 1576–1585. https ://doi.org/10.1128/ AEM.01996 -06 (2007). 53. Anderson, M. A new method for non-parametric multivariate analysis of variance. Austral. Ecol. 26, 32–46. https ://doi.org/10.11 11/j.1442-9993.2001.01070 .pp.x (2001). 54. Oksanen, J. et al. Vegan: community ecology package. R package version 2.5-6. https ://cran.r-pr oje ct.org,htt ps ://gi thu b.com/vegan devs/vegan (2019). 55. Friedman, J. & Alm, E. J. Inferring correlation networks from genomic survey data. PLoS Comput. Biol. 8, e1002687. https ://doi. org/10.1371/journ al.pcbi.10026 87 (2012). 56. Weiss, S. et al. Correlation detection strategies in microbial data sets vary widely in sensitivity and precision. ISME J. 10, 1669–1681. https ://doi.org/10.1038/ismej .2015.235 (2016). 57. Kurtz, Z., Mueller, C., Miraldi, E. & Bonneau, R. SpiecEasi: sparse inverse covariance for ecological statistical inference. R package version 1.1.0. https ://githu b.com/zdk12 3/Spiec Easi (2020). 58. Csardi, G. & Nepusz, T. The igraph software package for complex network research. InterJournal Complex Syst. 1695, 1–9 (2005). 59. Golbeck, J. Chapter3 network structure and measures. In Analyzing the social web (ed. Morgan, K.) 25–44 (2013). 60. Martín-González, A. M., Dalsgaard, B. & Olesen, J. M. Centrality measures and the importance of generalist species in pollination networks. Ecol. Complex. 7, 36–43. https ://doi.org/10.1016/j.ecoco m.2009.03.008 (2010). 61. Banerjee, S., Schlaeppi, K. & van der Heijden, M. G. A. Keystone taxa as drivers of microbiome structure and functioning. Nat. Rev. Microbiol. 16, 567–576. https ://doi.org/10.1038/s4157 9-018-0024-1 (2018). 62. Lê Cao, K. A., Boitard, S. & Besse, P. Sparse PLS discriminant analysis: Biologically relevant feature selection and graphical displays for multiclass problems. BMC Bioinformatics. 12, 253–253. https ://doi.org/10.1186/1471-2105-12-253 (2011). 63. Hosmer, D. J., Lemeshow, S. & Sturdivant, R. Applied Logistic Regression (Wiley, New York, 2013). 64. Boers, S. A., Jansen, R. & Hays, J. P. Understanding and overcoming the pitfalls and biases of next-generation sequencing (NGS) methods for use in the routine clinical microbiological diagnostic laboratory. Eur. J. Clin. Microbiol. Infect. Dis. 38, 1059–1070. https ://doi.org/10.1007/s1009 6-019-03520 -3 (2019). 65. Robinson, C. K., Brotman, R. M. & Ravel, J. Intricacies of assessing the human microbiome in epidemiologic studies. Ann. Epide‑ miol. 26, 311–321. https ://doi.org/10.1016/j.annep idem.2016.04.005 (2016). 66. de la Cuesta-Zuluaga, J. & Escobar, J. S. Considerations for optimizing microbiome analysis using a marker gene. Front. Nutr. 3, 26. https ://doi.org/10.3389/fnut.2016.00026 (2016). 67. Gomar-Vercher, S., Simon-Soro, A., Montiel-Company, J. M., Almerich-Silla, J. M. & Mira, A. Stimulated and unstimulated saliva samples have significantly different bacterial profiles. PLoS ONE 13, e0198021. https ://doi.org/10.1371/journ al.pone.01980 21 (2018). 68. Yang, F. et al. Saliva microbiomes distinguish caries-active from healthy human populations. ISME J. 6, 1–10. https ://doi. org/10.1038/ismej .2011.71 (2012). 69. Zhou, J. et al. Influences of pH and iron concentration on the salivary microbiome in individual humans with and without caries. Appl. Environ. Microbiol. 83, e02412-16. https ://doi.org/10.1128/AEM.02412 -16 (2017). 70. Rupf, S. et al. Comparison of initial oral microbiomes of young adults with and without cavitated dentin caries lesions using an insitu biofilm model. Sci. Rep. 8, 14010. https ://doi.org/10.1038/s4159 8-018-32361 -x (2018). 71. Armitage, G. C. Development of a classification system for periodontal diseases and conditions. Ann. Periodontol. 4, 1–6. https :// doi.org/10.1902/annal s.1999.4.1.1 (1999). 72. Acharya, A. et al. Species-level salivary microbial indicators of well-resolved periodontitis: A preliminary investigation. Front. Cell. Infect. Microbiol. 9, 347. https ://doi.org/10.3389/fcimb .2019.00347 (2019). 73. Damgaard, C. et al. Porphyromonas gingivalis in saliva associates with chronic and aggressive periodontitis. J. Oral Microbiol. 11, 1653123. https ://doi.org/10.1080/20002 297.2019.16531 23 (2019). 74. Belstrom, D. et al. Microbial profile comparisons of saliva, pooled and site-specific subgingival samples in periodontitis patients. PLoS ONE 12, e0182992. https ://doi.org/10.1371/journ al.pone.01829 92 (2017). 75. Li, Y. et al. Oral microbiome in Chinese patients with aggressive periodontitis and their family members. J. Clin. Periodontol. 42, 1015–1023. https ://doi.org/10.1111/jcpe.12463 (2015). 76. Chen, H. et al. A Filifactor alocis-centered co-occurrence group associates with periodontitis across different oral habitats. Sci. Rep. 5, 9053. https ://doi.org/10.1038/srep0 9053 (2015). 77. Huse, S. M., Ye, Y., Zhou, Y. & Fodor, A. A. A core human microbiome as viewed through 16S rRNA sequence clusters. PLoS ONE 7, e34242. https ://doi.org/10.1371/journ al.pone.00342 42 (2012). 78. Ling, Z., Liu, X., Wang, Y., Li, L. & Xiang, C. Pyrosequencing analysis of the salivary microbiota of healthy Chinese children and adults. Microb. Ecol. 65, 487–495. https ://doi.org/10.1007/s0024 8-012-0123-x (2013). 79. De Filippis, F. et al. The same microbiota and a potentially discriminant metabolome in the saliva of omnivore, ovo-lacto-vegetarian and Vegan individuals. PLoS ONE 9, e112373. https ://doi.org/10.1371/journ al.pone.01123 73 (2014). 80. Ding, T. & Schloss, P. D. Dynamics and associations of microbial community types across the human body. Nature 509, 357–360. https ://doi.org/10.1038/natur e1317 8 (2014). 81. Noecker, C., McNally, C. P., Eng, A. & Borenstein, E. High-resolution characterization of the human microbiome. Transl. Res. 179, 7–23. https ://doi.org/10.1016/j.trsl.2016.07.012 (2017). 82. Cao, Y., Fanning, S., Proos, S., Jordan, K. & Srikumar, S. A review on the applications of next generation sequencing technologies as applied to food-related microbiome studies. Front. Microbiol. 8, 1829. https ://doi.org/10.3389/fmicb .2017.01829 (2017). 83. Callahan, B. J., McMurdie, P. J. & Holmes, S. P. Exact sequence variants should replace operational taxonomic units in marker-gene data analysis. ISME J. 11, 2639–2643. https ://doi.org/10.1038/ismej .2017.119 (2017). 84. Segata, N. et al. Composition of the adult digestive tract bacterial microbiome based on seven mouth surfaces, tonsils, throat and stool samples. Genome Biol. 13, R42-r42. https ://doi.org/10.1186/gb-2012-13-6-r42 (2012). 85. Hall, M. W. et al. Inter-personal diversity and temporal dynamics of dental, tongue, and salivary microbiota in the healthy oral cavity. NPJ Biofilms Microbiomes 3, 2. https ://doi.org/10.1038/s4152 2-016-0011-0 (2017). 22 Vol:.(1234567890) Scientific Reports | (2021) 11:929 | https://doi.org/10.1038/s41598-020-79875-x www.nature.com/scientificreports/ 86. Tezal, M., Scannapieco, F. A., Wactawski-Wende, J., Grossi, S. G. & Genco, R. J. Supragingival plaque may modify the effects of subgingival bacteria on attachment loss. J. Periodontol. 77, 808–813. https ://doi.org/10.1902/jop.2006.05033 2 (2006). 87. Martin, F. E., Nadkarni, M. A., Jacques, N. A. & Hunter, N. Quantitative microbiological study of human carious dentine by culture and real-time PCR: Association of anaerobes with histopathological changes in chronic pulpitis. J. Clin. Microbiol. 40, 1698–1704. https ://doi.org/10.1128/jcm.40.5.1698-1704.2002 (2002). 88. Jiang, D. et al. Microbiome multi-omics network analysis: Statistical considerations, limitations and opportunities. Front. Genet. 10, 995. https ://doi.org/10.3389/fgene .2019.00995 (2019). 89. Lupatini, M. et al. Network topology reveals high connectance levels and few key microbial genera within soils. Front. Environ. Sci. 2, 10. https ://doi.org/10.3389/fenvs .2014.00010 (2014). 90. Das, P., Ji, B., Kovatcheva-Datchary, P., Bäckhed, F. & Nielsen, J. Invitro co-cultures of human gut bacterial species as predicted from co-occurrence network analysis. PLoS ONE 13, e0195161. https ://doi.org/10.1371/journ al.pone.01951 61 (2018). 91. Manirajan, B. A. et al. Diversity, specificity, co-occurrence and hub taxa of the bacterial-fungal pollen microbiome. FEMS Microbiol. Ecol. 94. https ://doi.org/10.1093/femse c/fiy11 2 (2018). 92. Shetty, S. A., Hugenholtz, F., Lahti, L., Smidt, H. & de Vos, W. M. Intestinal microbiome landscaping: Insight in community assemblage and implications for microbial modulation strategies. FEMS Microbiol. Rev. 41, 182–199. https ://doi.org/10.1093/femsr e/ fuw04 5 (2017). 93. Banerjee, S., Schlaeppi, K. & van der Heijden, M. G. A. Reply to “can we predict microbial keystones?”. Nat. Rev. Microbiol. 17, 194. https ://doi.org/10.1038/s4157 9-018-0133-x (2019). 94. Röttjers, L. & Faust, K. Can we predict keystones?. Nat. Rev. Microbiol. 17, 193. https ://doi.org/10.1038/s4157 9-018-0132-y (2019). 95. Tettamanti Boshier, F. A. et al. Complementing 16S rRNA gene amplicon sequencing with total bacterial load to infer absolute species concentrations in the vaginal microbiome. mSystems 5, e00777-19. https ://doi.org/10.1128/mSyst ems.00777 -19 (2020). Acknowledgements The authors would like to thank the Genomics Unit at the Centre de Regulació Genòmica (CRG, Barcelona, Spain) for assistance with the 16S rDNA gene amplicon sequencing. Author contributions Selection of patients and collection of salivary samples: M.R, F.S., J.J.P., C.C.; Extraction of DNA from salivary samples: C.H.; Bioinformatic and statistical analysis: C.B.C., A.R.I., I.T.; Study of the literature: A.R.I., I.T.; Preparation and review of the manuscript: C.B.C., A.R.I., I.T. Funding This investigation was supported by the Instituto de Salud Carlos III (General Division of Evaluation and Research Promotion, Madrid, Spain) and co-financed by FEDER (“A way of making Europe”) under grant ISCIII/PI17/01722, and the CESPU under grants MVOS2016 and MVOS-PT-IINFACTS-2019. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. Competing interests The authors declare no competing interests. Additional information Supplementary Information The online version contains supplementary material available at https ://doi. org/10.1038/s4159 8-020-79875 -x. Correspondence and requests for materials should be addressed to I.T. Reprints and permissions information is available at www.nature.com/reprints. Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/. © The Author(s) 2021