Salivary metabolomics profile of patients with recurrent aphthous ulcer as revealed by liquid chromatography-tandem mass spectrometry.

Objective We compared the salivary nontargeted metabolite profiles between patients with recurrent aphthous ulcer (RAU) and healthy individuals to investigate the metabolic alterations associated with RAU. Methods Saliva samples were collected from 45 patients with RAU and 49 healthy individuals, and the salivary metabolites were quantified using liquid chromatography-tandem mass spectrometry. The metabolomic profiles were then analyzed using multivariate and univariate statistical methods, and enrichment of the metabolites in various biological pathways was assessed. Results In total, 206 significant differentiating metabolites (Wilcoxon test, false discovery rate [FDR] of <0.05) were identified between patients with RAU and healthy individuals. These metabolites were implicated in tryptophan metabolism, steroid hormone biosynthesis, and other metabolic pathways. Two commonly circulating steroids, estrone sulfate and dehydroepiandrosterone sulfate, were significantly lower in the saliva of patients with RAU (Wilcoxon test, FDR < 0.05, power > 0.9). Principal component analysis and partial least-squares discriminant analysis revealed metabolic perturbations involving RAU, and receiver operating characteristic curve analysis with several metabolites showed good diagnostic ability for RAU. Conclusions The results of this study indicate that patients with RAU are characterized by metabolic imbalances. Psychogenic factors, endocrinopathies, and immunosuppression may contribute to the onset of RAU.

Introduction

Recurrent aphthous ulcer (RAU) is one of the most common diseases affecting the oral mucous membrane. The incidence of RAU varies but is generally considered to be close to 20% in any given population.[1-3] RAU always occurs in the mucous membrane, including that of the lips, tongue, and soft palate; it occurs less often in the isometric mucosa of the hard palate. The oral ulcers can be round or oval, exhibit circumscribed margins and yellow or gray floors, and are surrounded by a 1-mm-wide inflamed region. Pressing the center of these ulcers leads to a burning sensation.[4,5]

The exact pathogenesis of RAU remains unknown. Various causes of RAU have been reported in several studies and include the presence of certain oral microbial communities,[6,7] immunological factors,[4] endocrinopathies, and psychological and hereditary factors.[1,8] The general consensus in the medical community is that low immunity and imbalanced T-cell subgroups are involved in the development of RAU,[4,9,10] but the exact pathogenesis is unclear. Oral dysbacteriosis or adventitious infections may facilitate the onset of RAU; in particular, the microflora can produce virulence to and inhibit the proliferation of oral epithelial cells.[6] Babaee et al.[11] found that the oxidant/antioxidant status of blood and saliva was markedly different between patients with and without RAU. Psychogenic factors have also been proven to be associated with the onset of RAU. Mental stress may change the levels of certain hormones, such as cortisone, which affects the immune status and increases the risk of RAU.[8,9,12] However, no single potential pathogenic factor can independently explain the mechanism of RAU, and the various conclusions from different studies remain under debate.[13]

No metabolomics research about RAU has been reported to date. Therefore, to investigate the metabolic alterations in the saliva of patients with RAU, we employed liquid chromatography–tandem mass spectrometry (LC-MS/MS) to detect salivary metabolites in both patients with RAU and healthy individuals. The relative content of salivary metabolites was analyzed using univariate and multivariate statistical methods to examine the different expression of salivary metabolites between the two groups.

Materials and methods

Ethics statements

Written informed consent was obtained from all patients in accordance with the Declaration of Helsinki. This study was approved by the Institutional Review Board on Bioethics and Biosafety of BGI (BGI-IRB, Shenzhen 518083, China) with approval number BGI-IRB 15079.

Sample collection

All participants in this study were volunteers either with or without RAU. Each volunteer filled out a questionnaire about his or her lifestyle habits and other background information. Saliva samples (3 mL) were collected between 7:00 and 10:00 AM. After collection, the samples were stored in a freezer at −80°C.

Metabolite extraction

All saliva samples were stored at −80°C until sample preparation. The samples were first placed in a −20°C environment for 30 minutes and then thawed at 4°C in a refrigerator until the ice disappeared. A 100-µL sample was mixed with 300 µL of methanol and vortexed for 1 minute, and centrifugation was then performed at 14,000 × g and 4°C for 20 minutes. Finally, 250 µL of supernatant of each sample was preserved for LC-MS/MS detection. Quality control (QC) samples, which were a mixture of equal volume taken from each real sample, also underwent LC-MS/MS analysis for quality assurance of the experiment.

Ultra-performance LC-MS

We employed an ultra-performance LC (UPLC) system (2777C; Waters, Milford, MA, USA) with an ACQUITY UPLC BEH C18 column (100 × 2.1 mm, 1.7 µm) to separate the salivary metabolites. Mobile phase A was water and 0.1% formic acid. Mobile phase B was methanol and 0.1% formic acid. The gradient elution was as follows, with a flow rate of 0.4 mL/minute: 0% B for the first 2 minutes, 100% B from 2 to 12 minutes, maintenance at 100% B for the next 2 minutes, and a return to 0% B in the final 1 minute. The SYNAPT G2XS QTOF (Waters) was coupled with the UPLC system. Electrospray ionization was used as the ion source, and data-independent acquisition (MSE) data were acquired in positive mode with a scan range of 50 to 1200 Da. The scan time was 0.2 s. The sampling cone was set at 40 V, and the capillary voltage was 1.3 kV. The source temperature and desolation temperature were 120°C and 500°C, respectively. The desolation and cone gas flow rates were maintained at 800 and 50 L/h, respectively. The collision energy was set at 20 to 40 eV.

Data processing

LC-MS/MS raw data were imported into the commercial software Progenesis QI (version 2.0; Nonlinear Dynamics, Newcastle upon Tyne, UK) for peak alignment and peak-picking. The exported result was imported to metaX,[14] an R package software developed in-house for subsequent data processing and statistical and biomarker analysis.

Features with >80% missing values in the real samples or >50% missing values in the QC samples were deleted first. The k-nearest neighbors method was then performed for missing value imputation.[15,16] Next, QC-based robust local regression signal correction, which is an accepted correction method in the metabolomics field, was employed for QC based on the QC samples to correct baseline drift and eliminate the batch effect.[15,17] Data were normalized by the projected quasi-Newton method.[15] After data normalization, we filtered out the low-quality ions with a relative standard deviation (RSD) of >30% in the QC samples. Student’s t test and the Wilcoxon rank sum test were both used to screen significant features between the RAU and control groups.[18] The p-values were adjusted using the Benjamini–Hochberg method for multiple-hypothesis testing, and the false discovery rate (FDR) was limited to 0.05. The metabolomics data were also analyzed by multivariate statistical methods, including principal component analysis (PCA) and partial least-square discrimination analysis (PLS-DA).[15,19] Before implementation of PLS-DA, the Pareto scaling method was used to scale the data.[15] A receiver operation characteristic (ROC) curve was calculated for each feature to select prominent discriminatory variables between the two groups.[20,21] Hierarchical cluster analysis was performed with the selected features to observe relative intensity changes.[22] Permutational multivariate analysis of variance (PERMANOVA) analysis was performed to evaluate the correlation between the phenotypes and metabolomic profiles.[23]

Identification was conducted using Progenesis QI (version 2.0). Exact mass and isotope similarities were matched with the HMDB database (http://www.hmdb.ca/) to obtain putative metabolites. For more convincing results, the deconvolution experimental MS/MS spectra of each feature were compared to the theoretical fragment of the spectral library by MetFrag (http://msbi.ipb-halle.de/MetFrag).

Results

Participant characteristics

The volunteers comprised 47 men and 47 women ranging in age from 21 to 64 years (mean age, 27.8 years). The experimental group comprised 45 volunteer patients with RAU, and the control group comprised 49 healthy volunteers without RAU. In total, 94 saliva samples were used in this study. PERMANOVA showed that the onset time (p = 0.0004) and immediate family members with RAU (IFMRAU) (p = 0.0172) were significant factors explaining the variation in the examined salivary metabolic samples (Table 1). The onset of RAU was correlated with the presence of an IFMRAU onset history (χ2 test, p < 0.05, power = 0.824). This indicates that the occurrence of RAU may be associated with genetic factors.

Table 1.

PERMANOVA analysis of phenotype and metabolomics profile

Phenotype	R2	p-value
Age	0.01303	0.2791
Sex	0.013426	0.2631
OT	0.076922	0.0004
CDD	0.020249	0.1145
LMS	0.02043	0.1073
PSF	0.021954	0.0939
LC	0.004766	0.8079
IFMRAU	0.035742	0.0172
PAF	0.010792	0.3737
LI	0.010829	0.3864

PERMANOVA: permutational multivariate analysis of variance; OT: onset time (days); CDD: chronic digestive diseases, yes or no; LMS: long-term mental stress, yes or no; PSF: preference for spicy food, yes or no; LC: long-term constipation, yes or no; IFMRAU: immediate family members with recurrent aphthous ulcer, yes or no; PAF: preference for acidic foods, yes or no; LI: long-term insomnia, yes or no.

QC analysis

In total, 7859 features were extracted after data preprocessing. After normalization, 7218 features with an RSD of <30% constituted 91.84% of all QC samples, and the correlation coefficient among these QC samples was >0.95. These results indicate that there was no obvious batch influence on the acquired data after preprocessing. We also analyzed the data by PCA to observe the consistency of the QC samples. The PCA score plots representing the QC samples could be gathered together relative to the other samples (Figure 1(a)). These results indicate that our experiment was stable and that the data we obtained were eligible for subsequent statistical analysis.

Figure 1.

PCA and PLS-DA score plots. QC of PCA and multivariable statistical analysis. (a) PCA score plot includes QC samples. The QC samples are clustered tightly, indicating high experimental quality. (b) PCA score plot between RAU and control groups. (c) PLS-DA score plot. (d) R2 and Q2 of PLS-DA with 200 permutation tests. The left points represent Q2 and R2 of 200 permutation tests, and the right two points represent Q2 and R2 of the real PLS-DA model. Both values are higher than those in the permutations tests, respectively, indicating the model’s robustness. PCA, principal component analysis; PLS-DA, partial least-square discrimination analysis; QC, quality control; RAU, recurrent aphthous ulcer

Univariate analysis

We performed a univariate analysis using Student’s t test and the Wilcoxon rank sum test on the metabolomic profile to screen out metabolites showing statistically significant associations with RAU. Features with an FDR of <0.05 were considered to be significantly expressed between the RAU and control groups. Next, we screened the features with an intensity fold change of >1.20 or <0.83. A large number of variances (total of 1063 features) met the above criteria and were selected as those that differentiated patients with RAU from healthy controls.

Multivariate analysis

To identify potential biomarkers from the massive metabolite data, we performed a multivariate statistical analysis after the univariate analysis. First, we used the unsupervised method (PCA). The PCA scoring plot showed that the salivary metabolomic profiles between the RAU and control groups were separated into two clusters (Figure 1(b)), allowing for a totally visual judgment regarding the difference in the composition of small-molecule metabolites in the saliva of patients with RAU and control participants. We further built a predictive model by PLS-DA, a supervised discriminant analysis method with which to distinguish the two groups. Figure 1(c) shows the plot of the PLS-DA scores, demonstrating obvious separation between the RAU and control groups. R2 and Q2, parameters of the PLS-DA model, reflect the goodness and predictive abilities of the PLS-DA model. In this model, R2 and Q2 are 0.91 and 0.76, respectively. These values are high, basically explaining the huge model variance for the two groups. To further validate the reliability of the model, a permutation test with 200 repeats was performed (Figure 1(d)). The p-value of R2 was 0.01 and that of Q2 was 0.00, proving that this model was not overfitting. The variable importance in projection (VIP) calculated from the PLS-DA model reflects the influence of features in constructing the PLS-DA model.

Screening and identification of differentiating metabolites

To screen out significant differentiating metabolites from the massive numbers of small-molecule metabolites, we used the VIP of the first two principal components of the PLS-DA model, then combined the fold change and q-value as filtering factors. The screening conditions were as follows: VIP, ≥1; fold change, >1.20 or <0.83; and q-value, <0.05. We identified 939 differentiating features. Among the identified features, 638 ions were identified. With only the mass accuracy search, it appeared that one feature corresponded to many candidates. Therefore, we further confirmed the identification by combining isotope similarity and matched MS/MS fragments to select candidates with the highest score. We set the threshold of isotope similarity at >60, combined with a total score (considering overall mass error, isotope similarity, and fragment matching) of >30. Using this method, we finally identified 206 differentiating metabolites.

Metabolites involved in metabolic pathways

Differentiating metabolites were mapped to the metabolic pathways in the Kyoto Encyclopedia of Genes and Genomes (KEGG, http://www.genome.jp/kegg/) to help understand the metabolites involved in the main biochemical metabolic pathways. We performed a pathway impact analysis by MetaboAnalysis (http://www.metaboanalyst.ca/). Tryptophan metabolism and steroid hormone biosynthesis were enriched. The intensity of tryptamine, formyl-5-hydroxykynurenamine, 5-methoxytryptamine, and indoleacetaldehyde enriched in tryptophan metabolism changed prominently. We also found that 5-methoxytryptophan, the downstream metabolite of tryptamine, was significantly increased in patients with RUA. Estrone sulfate (E1S), dehydroepiandrosterone sulfate (DHEAS), and 17β-estradiol 3-sulfate were enriched in steroid hormone biosynthesis. Specific information, including statistical comparisons of the p-values and power values of these metabolites, are shown in Table 2, and their relative intensities in the two groups are shown in the boxplot in Figure 2. We performed an ROC analysis with the combination of tryptamine, formyl-5-hydroxykynurenamine, 5-methoxytryptamine, indoleacetaldehyde, and 5-methoxytryptophan. The area under the ROC curve was 0.989 (Figure 3(a)), indicating that the diagnostic capability of these different metabolites was very good. The ROC analysis of E1S, DHEAS, and 17β-estradiol 3-sulfate showed an area under the curve of 0.813, which was also good (Figure 3(b)).

Table 2.

Eight differentiating metabolites enriched in tryptophan metabolism and steroid hormone biosynthesis

Metabolites	RT (minutes)	m/z	Ratio	p-value	Power
Tryptamine	3.94	183.0886	1.379766	5.84 × 10−16	1
Formyl-5-hydroxykynurenamine	4.51	439.157	0.562064	0.000897	1
5-Methoxytryptamine	4.60	229.0746	1.964598	5.43 × 10−12	1
Indoleacetaldehyde	4.68	319.1463	0.45484	1.23 × 10−11	0.999694
5-Methoxytryptamine	3.58	217.095	1.483523	2.59 × 10−17	1
Estrone sulfate	4.39	368.1492	0.527014	4.25 × 10−7	1
17β-Estradiol 3-sulfate	4.63	353.1402	0.637748	1.89 × 10−6	1
Dehydroepiandrosterone sulfate	6.24	351.163	0.653309	0.003535	0.999789

RT: retention time; m/z: mass-to-charge ratio; Ratio: fold change in patients relative to controls.

p-value: Wilcoxon test after Benjamini–Hochberg correction.

Power: Statistical power at α = 0.05, two-tailed.

Figure 2.

Box plots of (a) the five metabolites enriched in tryptophan metabolism and (b) the three metabolites enriched in steroid hormone biosynthesis. ***q-value < 0.001, **q-value < 0.01.

Figure 3.

ROC curves of metabolites enriched in tryptophan metabolism and steroid hormone biosynthesis. (a) ROC curve of the combination of tryptamine, formyl-5-hydroxykynurenamine, 5-methoxytryptamine, indoleacetaldehyde, and 5-methoxytryptophan. The area under the curve is 0.989. (b) ROC curve of the combination of estrone sulfate, 17β-estradiol 3-sulfate, and dehydroepiandrosterone sulfate. ROC, receiver operating characteristic.

Discussion

We compared the salivary metabolomic profiles between patients with RAU and healthy individuals. Hundreds of differentiating metabolites were discovered by our salivary metabolomic method, exhibiting the metabolic alterations association with RAU. These differentiating metabolites were found to be enriched in various important biological function pathways, including tryptophan metabolism and steroid hormone metabolism.

The PERMANOVA results showed that the presence of IFMRAU has a significant impact on the salivary metabolite composition. Alterations in the salivary metabolite composition and the onset of RAU were correlated with the presence of an IFMRAU onset history, indicating that genetic factors play an important role in the onset and development of RAU. Many gene polymorphisms are reportedly correlated with RAU, including those of the interleukin-6, -1, and -10 genes.[1,24,25]

Tryptamine, an endogenous metabolite that can induce the release of serotonin,[26] was significantly increased in the saliva of patients with RAU, suggesting that imbalanced tryptophan metabolism may be associated with the incidence of oral ulcers. Recent studies have shown that elevated salivary serotonin is positively correlated with detrimental psychological factors including depression and stress.[27] Additionally, depression is reportedly an important psychological factor in the pathogenesis of RAU,[8] and an increased tryptamine level may play a role in the occurrence of RAU caused by negative emotions. However, the exact molecular mechanism of how the salivary tryptamine and tryptophan metabolism pathways are involved in the onset of RAU remains unclear and requires further study.

The salivary levels of DHEAS and E1S, which are two of the most abundant steroids in the human circulation,[28] showed significant alterations between patients with RAU and healthy individuals. This supports existing evidence that hormone imbalances increase the risk of RAU onset.[1,4] The salivary level of DHEAS was significantly lower in patients with RAU than controls. A low salivary DHEAS level is reportedly more common in individuals with depression,[29,30] as mentioned above, and depression is an important risk factor for the onset of RAU.[8] Additionally, DHEAS has been shown to be a negatively predictor of salivary immunity, in that a decreased DHEAS level may be associated with a decline in salivary immunity; thus, DHEAS may be correlated with an increased risk of RAU onset.[4,9,10]

In summary, we compared the metabolomic profiles between patients with RAU and healthy individuals and observed salivary metabolic alterations in patients with RAU. Tryptophan and hormone metabolism imbalances have been shown to be correlated with detrimental psychological factors including depression and stress, and a decline in salivary immunity may also play a role in the pathogenesis of RAU. Further studies involving a larger population are recommended.