A “seed”-based molecular networking strategy for the screening and identification of unknown glucocorticoids in cosmetics

2.1. Chemicals and materials

A total of 36 glucocorticoid reference standards used to establish the standard method were obtained from Shanghai Aladdin Biochemical Technology Co., Ltd. (Shanghai, China) and National Institutes for Food and Drug Control (Beijing, China), and more information of these reference standards has been provided in Table 1. In addition, the reference standard of clobetasol propionate, halobetasol propionate, betamethasone dipropionate, and prednisone acetate were acquired from Shanghai Maclin Biochemical Technology Co., Ltd. (Shanghai, China). Dexamethasone phosphate (Dex-P) were bought from CATO Research Chemicals Inc. (Guangzhou, China). MS grade acetonitrile, methanol, ammonium acetate and formic acid were bought from Shanghai Macklin Biochemical Technology Co., Ltd. (Shanghai, China). Ultrapure water was obtained from a Milli-Q advantage ultrapure water system (Merck Chemical Co., Ltd.). Seven cosmetic samples (labeled C1-C7) were collected from the online and offline markets, other samples came from supervised sampling tasks. Additionally, these cosmetic categories include repair creams, light printing creams, whitening creams, pore unblocking liquids, and skin brightening essences.

2.2. Stock and working solutions

Individual stock solutions of 36 reference standards were prepared with concentration of 1.0 mg/mL using methanol. A mixture working solution with a concentration of 1.0 μg/mL was prepared by diluting the individual standard solutions with acetonitrile. Working solutions of Dex-P were prepared in methanol at a concentration of 1.0 μg/mL. All solutions were stored in the dark at 4°C.

2.3. Sample preparation

Approximately 500 mg of each cosmetic samples were measured and placed in a 20 mL stopper colorimetric tube. Next, 3 mL of a saturated sodium chloride solution was added to the tube containing the cosmetics sample and swirled for 1 min to form a uniform solution. Following this, 10 mL of a 0.1% formic acid solution in methanol was meticulously added to the tube, and the mixture was sonicated for 30 mins. After ultrasonication, the samples were centrifuged at a speed of 4000 revolutions per minute (r/min) for a period of 10 min. Finally, 1 mL of the supernatant liquid was carefully pipetted and filtered through a 0.22-micrometer (μm) filter membrane and stored at a temperature of 4°C until further analysis.

2.4. Chromatographic and mass spectrometric analysis conditions

Ultra-high performance liquid chromatography - quadrupole time-of-flight mass spectrometry/mass spectrometry (UHPLC-Q-TOF-MS/MS) was performed using an Agilent UPLC system equipped with Agilent G6545 Q-TOF-MS/MS system (Agilent Corp, Santa Clara, USA). Chromatographic separation was performed on a Proshell 120 EC-C18 column (2.1 mm × 100 mm, 1.8 μm) with a column temperature of 25°C. The analytical flow rate was 0.35 mL/min. The mobile phase consisted of acetonitrile (A), a solution of 0.1% formic acid, and 2 mmol/L ammonium acetate (B). The percentage of mobile phase A was linearly varied as follows: 0-1 min, 2%, 1-3 mins, 2%-20%, 3-5 mins, 20%-28%, 5-22 mins, 28%-75%, 22-25 mins, 75%-95%, 25-27 mins, 95%-2%, 27-27.5 mins, 95%-2%, 27.5-30 mins, 2%. The injection volume was precisely set at 5 μL for all assays.

For MS analysis, the analyte was ionized via an electrospray ionization (ESI) source operating in positive mode, and data were acquired in a data-dependent acquisition (DDA) mode. The active exclusion mode was configured as follows: excluded after two spectra and released after 9 s. The temperature of the sheath gas was maintained at 250°C, and that of the drying gas was held at 350°C. The flow rate of the sheath gas was set at 10 L/min, and that of the drying gas was 8 L/min. The voltage of spray and capillary were 500 V and 4.0 kV, respectively. The collision voltage was set within the range of 10 to 40 V, and the non-targeted MS/MS mode was adjusted to analyze the top 5 ions.

2.5. Molecular networking analysis

Non-targeted screening of glucocorticoids was performed using FBMN workflow and data mining techniques [19, 20]. The FBMN analysis was conducted based on the cloud-based global natural products social molecular networking platform (GNPS) (https://gnps.ucsd.edu/) [15]. Briefly, the raw data from UPLC-Q-TOF-MS/MS analysis were converted into the mzML format using ProteoWizard MSConvert software, while peak picking was performed, followed by data pre-processing, which included data collection, peak detection, peak deconvolution, and peak alignment; completed by MS-DIAL 4.9.221218 [25]. The detailed parameters of MS-DIAL were as follows: minimum peak height set as 5000, mass slice width set as 0.01, MS1 tolerance set as 0.01 Da, MS2 tolerance set as 0.025 Da, sigma window value set as 0.5, tolerance of retention time (RT) set as 0.05 min. Subsequently, the pre-processed data files were subsequently uploaded to the GNPS webserver for molecular networking construction. It is important to emphasize that the MS data of the 36 standards were integrated with the MS data of seven sets of actual samples when obtaining the data files uploaded to GNPS so as to enable efficient clustering based on 36 known compounds when building molecular networking. For molecular networking parameters, the maximum mass tolerances of the parent ions and fragment ions were set to 0.02 Da and 0.02 Da, respectively, and the minimum cosine score and minimum matched peaks to 0.7 and 5, respectively. Other parameters were set as default. Finally, the molecular networking construction result was visualized using Cytoscape V3.7.1.

2.6. Data mining strategy

The structural identification of glucocorticoids in cosmetics mainly integrated MS-DIAL, MS-FINDER 3.52, and molecular networking analysis. The process began by quickly localizing the seeds of the 36 known reference standards to molecular clusters and individual nodes associated with glucocorticoids. Then, after entering the precursor ion and RT of the unknown node via MS-DIAL, the information of the unknown node was associated with MS-FINDER for structural annotation. A summary of the analytical workflow, including data acquisition, annotation, and interpretation, has been presented in Figure 1.

Close

Figure 1.

Workflow demonstrating the profiling and annotation of glucocorticoids.

Export to PPT

During the structure annotation process, MS-DIAL efficiently identifies precursor peaks with enhanced clarity by performing sophisticated deconvolution of the MS/MS spectra, effectively reducing background noise, and precisely extracting spectra from complex mixtures that contain unknown compounds [26]. The corresponding parameters, such as compound composition (e.g., C, H, O, N, F, and Cl atoms), number of candidate substances, and the search scope of database (e.g., GNPS, PubChem, MassBank, DrugBank, and MINEs), were set in MS-FINDER [26,27]. Through the internal algorithm in MS-FINDER, scores reflecting the similarity between potentially matched compounds and candidate compounds in the database were calculated. The identification of candidate compounds was then determined based on these scores, alongside the fit of the MS/MS spectra, m/z tolerance, elemental composition, chemical valence rules, and elemental ratios [27]. Moreover, the outcomes of candidate compound identification were classified into four levels of confidence, as detailed in Table S1 [28]. To address the potential sources of error caused by using it in data mining and analysis, we use three factors to analyze the results: fragment similarity, mass error, and identification confidence level.

2.7. Molecular docking and molecular dynamics simulations

The ability of glucocorticoids to cause skin atrophy and endocrine disruption was evaluated by molecular docking and dynamics on the active site of hyaluronan synthase 2 (HAS-2) and glucocorticoid receptor (GR) identified using Biovia Discovery studio 4.5 (BIOVIA, 2024), respectively. Structure of HAS-2 of Homo sapiens was obtained from AlphaFold Protein Structure Database (https://www.alphafold.ebi.ac.uk/) [29,30]. The structure of human GR (PDB id, 4UDC) was selected from the Protein Data Bank [23]. The 3D structure of glucocorticoids and energy minimization were constructed using ChemDraw 19.0. Molecular docking of the ligand and protein was applied through using CDOCKER program of Discovery studio 4.5. Hydrogen bonding and other interactions between the ligand and protein were visualized using Discovery Studio 4.5 and PyMoL.

The Dex-P-HAS-2 complex and Dex-P-GR complex with best binding conformations were selected for further MD simulation. The two Dex-P-protein complexes were solvated in orthorhombic cells, where 0.145 mol/L sodium chloride was placed randomly to achieve electrostatic neutrality. The simulation system carried out the energy minimization and equilibrium processes at constant temperature and constant pressure, respectively. After 50 ns dynamics simulation, the Root Mean Squared Deviation (RMSD) values, Root Mean Square Fluctuation (RMSF) values, Radius of gyration (Rg), key residues, and structural changes for Dex-P-protein were calculated, respectively.

3.1. Acquisition of “seed” data by UPLC-Q-TOF-MS/MS

For structure elucidation, known compounds were initially identified as “seed” nodes, which aid in inferring the structures of adjacent nodes [19,20]. A mixture of 36 reference standards, each at a concentration of 1.0 μg/L, was analyzed using UPLC-Q-TOF-MS/MS to obtain the “seed” data. Details of the 36 reference standards have been summarized in Table 1. Total ion chromatograms (TICs) of 36 reference standards have been displayed in Figure S1, and the MS/MS spectra of every standard have been shown in Figures S2-S4.

As indicated in Table 1, several fragment ions appeared frequently among the reference standards, notably at m/z 355 ([C₂₂H₂₆O₄]⁺), 337 ([C₂₂H₂₅O₃]⁺), 319 ([C₂₂H₂₃O₂]⁺), 309 ([C₂₁H₂₅O₂]⁺), and 237 ([C₈H₉N₂O₂]⁺). Chen et al. utilized these characteristic ions, particularly those at m/z 355, 337, and 319, to annotate the structures of betamethasone dibutyrate and betamethasone tributyrate [31]. Approximately half of the glucocorticoids showed a common fragment ion at m/z 279 ([C₂₀H₂₃O]⁺), which Karakka Kal et al. noted as a key characteristic ion of fluorometholone [32]. Similarly, Kim et al. observed this fragment ion in the structural fragments of beclomethasone, clobetasone butyrate, and mometasone furoate [11]. Additionally, fragment ions at m/z 355, 337, 319, 309, 279, and 237 were present in the fragmentation pathways of betamethasone, and dexamethasone [33]. Moreover, most glucocorticoids contained fragment ions at m/z 147, 135, and 121, which corresponded to the identity of [C₁₁H₁₅]⁺, [C₁₀H₁₅]⁺, and [C₈H₉O]⁺, respectively. This observation is consistent with the MS analysis of steroid structures conducted by Jian et al. [10]. Furthermore, Yang et al. found that clobetasol propionate contained the characteristic ions at m/z 147 and 121 [34]. In all, the MS/MS spectra of known glucocorticoids provide valuable structural information about compounds, and the presence of shared fragment ions facilitates the grouping and classification of similar compounds.

3.2. Identification of glucocorticoids using the “seed” molecular networking strategy

The integration of FBMN with UPLC-Q-TOF-MS/MS provides a rapid and efficient approach for the non-targeted screening of unknown compounds [19,20]. To efficiently identify unknown glucocorticoids, data from 36 reference standards were utilized as “seed” to construct molecular networking specific to glucocorticoids with seven cosmetic products (Figure S5). Each node in the molecular network corresponds to a unique feature molecule. As illustrated in Figure 2(a-d), the color of each node indicates the relative abundance of the compound in different groups: orange represents known reference standards used as seeds, grey indicates unidentified compounds, and other colors denote compounds identified using MS-DIAL/MS-FINDER. Consequently, a total of 37,724 positive precursor ions from actual samples, along with 36 reference standards, were organized into 129 clusters (nodes ≥ 2) and 36, 939 single nodes (Figure S6).

Close

Figure 2.

Molecular networking of all components in seven cosmetics and 36 reference standards. (a) A glucocorticoid-related molecular cluster of 27 nodes. (b) Results of 129 molecular networking clusters and partial single node. (c) A glucocorticoid-related molecular cluster of four nodes. (d) Four single node of glucocorticoid. Label of node refers to m/z, and the width of wire refers to the cosine score of fragment ions between neighboring nodes. The red circles represent molecular clusters and single nodes associated with glucocorticoids.

Export to PPT

3.2.1. Known glucocorticoids identified via reference methods

Compounds with similar MS/MS fragmentation patterns tend to form molecular families [20]. Utilizing “seed” compounds as indicators, we rapidly identified a glucocorticoid-related molecular cluster comprising 27 nodes, alongside a single glucocorticoid node among 36,939 individual nodes (Figure 2b-d).

As illustrated in Figure 2(b), betamethasone (m/z 393.2048, RT 9.368 min), paramethasone (m/z 393.2073, RT 9.456 min) and dexamethasone (m/z 393.2073, RT 9.516 min) were included in the “seed” data. From the visualization spectra, it was evident that all three compounds were present in samples C1 and C3. In addition, the “seed” compound dexamethasone acetate (m/z 435.2172, RT 12.766 min) was found in the same molecular clusters, primarily detected in sample C1. Moreover, the single node compound, prednisone (m/z 359.1830, RT 7.896 min), was identified in sample C5 (Figure 2d). The reliability of these identifications was further validated by matching experimental and theoretical fragment ions using MS-FINDER (Figure 3 and S7). According to the confidence levels of identification, these three compounds were considered as level 1 (Table S1) [28]

Close

Figure 3.

Experimental and virtual MS/Mass spectrometry spectra of glucocorticoids. (a) Dex-P, (b) prednisone, (c) betamethasone, (d) paramethasone, (e) dexamethasone, (f) prednylidene.

Export to PPT

3.2.2. Unknown glucocorticoids identified outside the reference methods

The structural analysis of adjacent unknown nodes was conducted based on the known structures of glucocorticoids. By integrating visualized molecular clusters with the application of MS-DIAL/MS-FINDER, similar glucocorticoid structures were identified through the automatic matching of MS/MS spectra of unknown compounds to databases [26,27].

As illustrated in Figure 2(b) and 3(a), an unknown node in sample C3 was adjacent to betamethasone acetate (m/z 435.2152). Both compounds shared fragment ions at m/z 355, 337, 319, 279, and 237, and it could be identified as dexamethasone phosphate (Dex-P, m/z 473.1735, RT 7.014 min). Salem et al. reported fragment ions for Dex-P at m/z 435 ([C₂₂H₂₈O₇P]⁺), 355, 337, and 319, which coincided with our findings, and the identification confidence level of Dex-P was classified as level 2a (Table S1) [35]. Furthermore, the unknown node in sample C2 was adjacent to the “seed” node of betamethasone-17-butyrate-21-propionate (m/z 519.2731). They shared the same fragment ions at m/z 355, 337, 319, and 279, leading to the inference that this unknown compound was betamethasone dipropionate (m/z 505.2592, RT 17.529 min) (Figure S8). Specifically, the identification of betamethasone dipropionate was confirmed by a standard substance, classifying it as level 1 (refer Figure S9, Subpart Figure S9a). Additionally, the node representing 7 alpha-thiospironolactone (m/z 375.1963, RT 13.188 min) in sample C3 shared the same fragment ions at m/z 319, 309, and 147 with the adjacent “seed” node of betamethasone. Lee et al. quantified spironolactone derivatives using the fragment ion at m/z 107.2 ([C₈H₁₁]⁺), which was also observed in 7 alpha-thiospironolactone in our study [36]. This observation further supported the reliability of the 7 alpha-thiospironolactone identification (Level 2a).

Taking the same workflows applied to Dex-P identification, this strategy was also employed for the identification of other compounds. As shown in Figure 2(c), clobetasol propionate (m/z 467.1998, RT 16.839 min), estrone glucuronide (m/z 447.1932, RT 16.839 min), and arisugacin B (m/z 467.1997, RT 16.786 min) were identified within the same molecular cluster and were detected in samples C1 and C4. Sparidans et al. reported the fragment ions for clobetasol propionate as 447 (-HF), 373 (-HF, -C₂H₅COOH), and 355 (-HF, -C₂H₅COOH, -H₂O), aligning perfectly with our findings (Figure S8) [37]. The identification of clobetasol propionate was further supported by studies conducted by Yang et al. and Guo et al. [34,38]. Additionally, validation using the reference standard confirmed the presence of clobetasol propionate (Level 1, Figure S10). Adjacent to the clobetasol propionate node, the node at m/z 447.1932 was identified as estrone glucuronide with a mass error of -19.45 Da, and the node at m/z 467.1997 was identified as arisugacin B with a mass error of -15.635 Da (Figure 3a). Given that the mass errors exceeded the requirement of 10 ppm [39], the confidence level of identification for estrone glucuronide and arisugacin B was considered as level 3.

Furthermore, a single node with an m/z of 373.2009 was proposed as prednylidene (m/z 373.2009, RT 9.516 min), showing characteristic fragments at m/z 355, 337, 319, and 279 (Figure 2d), and this identification was classified as level 2b (Table S1). Another single node with an m/z of 401.1925 was identified as prednisone acetate (RT 11.356 min), with main fragments corresponding to 341.1751 [C₂₁H₂₅O₄]⁺, 323.1652 [C₂₁H₂₃O₃]⁺, 313.1797 [C₂₀H₂₅O₃]⁺, 305.1541 [C₁₈H₂₅O₄]⁺, and 295.1693 [C₂₀H₂₃O₂]⁺. These characteristic ions were consistent with the findings of Guo et al. [38], as well as the reference standards (Level 1, Figure S9b). A single node with an m/z of 485.1902 was identified as halobetasol propionate (RT 16.798 min), and the characteristic ions at m/z 429.1601 and 389.1513 matched those reported by Karakka Kal et al. [32]. Moreover, the fragmentation pattern of halobetasol propionate identified in this study was consistent with the reference standard (Level 1, Figure S11).

The non-targeted screening method based on the “seed” molecular networking technology annotated a total of 14 glucocorticoids in seven samples. Notably, this included Dex-P, prednylidene, and 7 alpha-thiospironolactone, which were not identified by the reference method (Table 2). Previous studies have developed various qualitative and quantitative approaches for determining glucocorticoid concentrations in cosmetics [10,40,41]. For example, Giaccone et al. reported an analytical method capable of detecting 10 glucocorticoids [40], while Kim et al. established a method to identify 43 glucocorticoids [11]. Meng et al. introduced a comprehensive screening method targeting 100 illegal substances, including 40 glucocorticoids [41]. While these methodologies have significantly advanced cosmetic safety monitoring, their primary focus has been on well-known glucocorticoids that are already subject to stringent regulatory oversight [42]. Moreover, conventional methods typically require the use of reference standards for each compound during analysis, which is both resource-intensive and inefficient. Given the structural diversity and complexity of glucocorticoids, the efficient detection of unknown glucocorticoids in cosmetics remains a formidable challenge [12-14]. In contrast to the conventional approaches, the adoption of a “seed” node-centric molecular networking approach can meet the needs of current regulations. This innovative technique eliminates the need for an extensive pre-built standard library, enabling a more simplified and cost-effective workflow for the identification and characterization of unknown glucocorticoids in complex samples.

3.3. Structure deduction of unknown glucocorticoids

Based on the three-point rule of compound-induced cleavage mentioned in the literature, including the preferential location of cleavage, the sequence of transition state carbon ion stability and the source of product ion charge stability [33, 43-45], the proposed fragmentation pathway of Dex-P, 7 alpha-thiospironolactone, and prednylidene were conducted in this work to further verify the accuracy of the glucocorticoid results identified by MS-DIAL/MS-FINDER. Detailed information on the proposed fragmentation pathways have been presented in Figure 4, and Figure S12.

Close

Figure 4.

Mass spectrometry information and structural annotation results for Dex-P and 7 alpha-thiospironolactone. (a) MS/MS spectrometry information of Dex-P, (b) Fragmentation pathway of Dex-P, (c) MS/MS spectrometry information of 7 alpha-thiospironolactone, (d) Fragmentation pathway of 7 alpha-thiospironolactone.

Export to PPT

As seen from the fragmentation pathway of Dex-P (m/z 473.1735) Figure 4(a-b), the first fragment with m/z 453.1659 ([M+H-20]⁺) was associated with the elimination of HF. With subsequent losses of 2H₂O, fragment ions at m/z 435.1566 ([M+H-18]⁺) and then m/z 417.1511 ([M+H-18]⁺) were formed. As reported by Noh et al., loss of HF and H₂O was a common skeletal fragmentation pathway for betamethasone and dexamethasone and their derivatives [44]. Besides, the fragment ion at m/z 337.1815 ([M+H-80]⁺) was connected with the loss of HPO₃ from ion at m/z 417.151. The P-O bond (about 410 kJ/mol) has a higher bond energy than the C-O bond (about 326 kJ/mol), so the loss of HPO₃ was slightly later [46]. Interestingly, the fragment ion at m/z 337.1815 ([M+H-18]⁺) could be also generated by the ion at m/z 355.1894 from another fragmentation pathway. Moreover, the fragment ion at m/z 319.1686 corresponded with the elimination of 2H₂O ([M+H-36]⁺) from m/z 355.1894 and H₂O ([M+H-18]⁺) from ion at m/z 337.1815, respectively. Furthermore, the formation of a fragment ion at m/z 279.1775 ([M+H-40]⁺) was due to the loss of the methyl group and carbonyl group, and the fragment ion at m/z 237.1247 corresponded to the loss of 42 Da attributed to the elimination of C₃H₆.

As presented in Figure 4(c-d), the fragmentation pathway of 7 alpha-thiospironolactone (m/z 375.1963) included the product ions m/z 357.1802, 319.1682, 309.1867, 301.155, 291.1744, 281.1881, 253.123, 191.1067, 189.0881, 147.0801, 135.0809, and 107.0851, respectively. The product resulting from the elimination of H₂O was a signal with m/z 357.1802 ([M+H-18]⁺). According to the results of fragment ion at m/z 357.1802 possessed three different fragmentation pathways. First, the formation of a fragment ion at m/z 319.1682 ([M+H-38]⁺) corresponded with the loss of C₂O and followed by rearrangement. Second, the CH₄S was lost, resulting in a fragment at m/z 309.1867 ([M+H-48]⁺). Third, the fragment ion at m/z 301.155 corresponded to the loss of 56 Da attributed to the elimination of C₃H₄O. The fragment ion at m/z 291.1744 ([M+H-28]⁺) produced from the ion at m/z 319.1682 corresponds to carbanyl group loss. Concurrent with this fragmentation channel, the fragment ions at m/z 191.1067 (-100 Da), 135.0809 (-56 Da), and 107.0851 (-28 Da) resulted from the loss of C₅H₈S, C₆H₁₂, and C₂H₄, respectively. In addition, the fragment ion at m/z 281.1881 had relevance to the elimination of carbanyl group from ion at m/z 319.1682, and the fragment ions at m/z 189.0881 (-92 Da), 147.0801 (-42 Da), and 107.0851 (-40 Da) resulted from the loss of C₆H₄O, C₃H₆, and C₃H₆, respectively. Moreover, the fragment ion at m/z 253.123 ([M+H-48]⁺) was associated with the elimination of OS from the ion at m/z 301.155.

As shown in Figure S12, the prednylidene (m/z 373.2009) fragmentation pathway included the product ions m/z 355.1892, 337.18, 319.1687, 279.1769, 277.1572, 265.1585, 239.1410, 225.1410, 187.1071, and 161.0971, respectively. The fragment ions at m/z 355.1892 and 337.18 were both associated with the loss of H₂O. The m/z 319.1687 [M+H-18]⁺ and 279.1769 [M+H-58]⁺ were connected with the loss of hydroxyl and C₂H₂O₂, respectively. Subsequently, the fragment at m/z 319.1687 formed the fragment ion at m/z 277.1572 and 239.1410, respectively. And then, the fragment ions at m/z 225.1410 and 161.0971 resulted in the loss of CH₂ and C₅H₄, respectively. In addition, the fragment ions at m/z 265.1585, 225.1410, 187.1071, and 161.0971 corresponded to the loss of CH₂, C₃H₄, C₃H₂, and C₂H₂, respectively. These results not only validate the fragmentation pathways proposed for Dex-P, 7 alpha-thiospironolactone, and prednylidene but also corroborate the accuracy of the glucocorticoid identification obtained through MS-DIAL/MS-FINDER.

3.4. Toxicity prediction of identified glucocorticoids based on molecular docking and molecular dynamics simulations

It is worth emphasizing that unknown glucocorticoids found in the non-targeted screening have not been reported on their potential toxicity or adverse effects, probably because they also have few data on their occurrence. To assess their toxicity, we studied the interaction of glucocorticoids with toxic target proteins.

3.4.1. Ability of glucocorticoids to induce skin atrophy and endocrine disorders

Skin atrophy is a natural aspect of the aging process, but the use of glucocorticoids in skin disease treatment can expedite this progression [47,48]. Studies have shown that dexamethasone reduced hyaluronic acid levels by inhibiting HAS-2 activity in fibroblasts and keratinocytes, thereby contributing to skin atrophy [49]. Herein, the study of the interaction between various glucocorticoids and HAS-2 is helpful in understanding the derm-causing capacity of each glucocorticoid. The detailed analysis results are shown in Supplementary Information S1. In Figure 5(a-c) and Table S2, Dex-P showed a strong binding affinity of -399.929 kcal/mol to HAS-2, which is 5.14 to 8.56 times greater than that of dexamethasone and clobetasol propionate. Besides, prednylidene displayed a binding affinity of -62.9058 kcal/mol to HAS-2, surpassing that of dexamethasone (-46.7236 kcal/mol). The binding affinity of 7 alpha-thiospironolactone to HAS-2 was -50.789 kcal/mol, comparable to dexamethasone.

Supplementary Table 2

Download

Close

Figure 5.

Binding energies and interactions of glucocorticoids with HAS-2 and GR. (a) Binding energies (kcal/mol) of all compounds, (b) Two-dimensional interaction diagrams between Dex-P and HAS-2, (c) Three-dimensional interaction diagrams between Dex-P and HAS-2, (d) Two-dimensional interaction diagrams between Dex-P and GR, (e) Three-dimensional interaction diagrams between Dex-P and GR. HAS-2: Hyaluronan synthase 2, GR: Glucocorticoid receptor.

Export to PPT

Typically, endocrine disruptors primarily interfere with the endocrine system of organisms by interacting with 12 classical nuclear receptors, such as androgen receptor, estrogen receptor GR, etc., thereby causing adverse developmental, reproductive, neurological, and immune effects [50]. Additionally, exogenous glucocorticoids can mimic endogenous ones by binding to GR, disrupting the endocrine system [6,51]. Therefore, by studying the interaction between different glucocorticoids and GR, we can elucidate their endocrine-disrupting capabilities. As depicted in Figure 5(d-e), 14 glucocorticoids effectively bound to the ligand-binding domain (LBD) of the GR, like glucocorticoid antagonists such as mifepristone and 3,6-dibromocarbazole [52]. Singh et al. showed that dexamethasone could form hydrogen bonds with the same amino acid residues in the LBD, consistent with our findings [23]. Among these compounds, Dex-P had the strongest binding affinity (-207.562 kcal/mol) (Figure 5a and Table S3). Overall, unknown compounds identified through molecular networking also showed potential for inducing skin atrophy and endocrine disorders, with Dex-P having the greatest potential.

Supplementary Table 3

Download

3.4.2. Molecular mechanism of skin atrophy and endocrine disorders induced by dexamethasone phosphate

To explore the binding mechanism between Dex-P and the target protein, we performed 50 ns molecular dynamics (MD) simulations. In Figure 6(a-b), RMSD values for the Dex-P-HAS-2 complex increased continuously from 0 ns to 30 ns, while for the Dex-P-GR complex, they increased from 0 ns to 25 ns. After this period, RMSD values stabilized, suggesting that the systems had equilibrated. In both protein-ligand complexes, most RMSF were below 3 Å (Figure 6c-d), with residues within 4 Å of the ligand showing values below 2.5 Å (Figure S13). These findings suggested that the Dex-P could stably bind with HAS-2 and GR, respectively. In addition, the Rg value of Dex-P-HAS-2 complex was 26.16 ± 0.40 Å, while that of HAS-2 alone was 27.23 ± 0.71 Å (Figure 6e). The change in Rg value suggested HAS-2 became more compact upon ligand binding. For Dex-P-GR complex, the Rg value was 18.61 ± 0.18 Å compared to 18.47 ± 0.08 Å for GR alone, indicating minimal change in protein compactness and suggesting stable ligand binding (Figure 6f).

Close

Figure 6.

Structural dynamics evolution in molecular systems during 50 ns MD simulations. (a) RMSD plot of HAS-2 and Dex-P-HAS-2 complex, (b) RMSD plot of GR and Dex-P-GR complex, (c) RMSF plot of HAS-2 and Dex-P-HAS-2 complex, (d) RMSF plot of GR and Dex-P-GR complex, (e) Rg plot of HAS-2 and Dex-P-HAS-2 complex, (f) Rg plot of GR and Dex-P-GR complex. MD: Molecular dynamics, RMSD: Root mean square deviation, RMSF: Root mean square fluctuation, Rg: Radius of gyration, Dex-P: Dexamethasone-21-phosphate, HAS-2: Hyaluronan synthase 2, GR: Glucocorticoid receptor.

Export to PPT

Furthermore, we analyzed the frequency of interactions between the ligand and residues in MD simulations and investigated the conformational snapshots of MD trajectory at 0, 30 (25) ns, and 50 ns. In Figure S14, subparts Figures S14a and S14(c-h), ARG315, ARG268, TYR358, ARG259, LYS451, GLU267, and CYS279 occupied prominent positions in Dex-P binding to HAS-2 throughout MD simulation. For the Dex-P-GR complex, high interaction frequency was observed between Dex-P and residues ARG611, LEU608, LEU563, MET604, ARG569, MET560, and GLN570 during MD simulations (Figure S14b). As shown in Figure S15, the number of hydrogen bonds formed between Dex-P and residues increased from 3 to 4, and then to 6 at 0, 25, and 50 ns. Meanwhile, the interactions between Dex-P and ARG611, LEU563, and GLN570 remained. In short, Dex-P can bind to HAS-2 and GR, thereby inducing changes in protein conformation and altering protein activity. Consequently, Dex-P may affect the skin barrier and endocrine systems of organisms. However, the computer prediction is mainly to predict the possible interactions between proteins and ligands, and some experimental studies are needed to verify the theoretical results in the next step of research.

Supplementary Figure 14

Download

Supplementary Figure 15

Download

3.5. Qualitative analysis of dexamethasone phosphate in real samples

Based on the toxicity prediction results, the surreptitious addition of Dex-P, a newly identified risk substance, warrants increased scrutiny during the monitoring of actual samples. For this purpose, after acquiring the corresponding reference standard, we conducted tests on both the Dex-P standard and actual samples under identical chromatographic-mass spectrometric conditions to further confirm the presence of Dex-P. As shown in Figure 7(a-b) and S16, the retention time (RT) of the Dex-P standard was consistent with the detection time in sample C3. Additionally, the Dex-P standard exhibited identical fragmentation patterns as observed in the identification results, including m/z 435.1567, 355.1904, 337.1774, 319.1704, and 237.1250 (Figure S17). These findings indicate that the findings from the results obtained from “seed”-based molecular networking were substantiated by the qualitative analysis based on the Dex-P standard.

Close

Figure 7.

EIC plot of Dex-P in standards and real sample. (a) Dex-P standard; (b) Sample C3; (c) Sample C8; (d) Sample C9; (e) Sample C10; (f) Sample C11. EIC: Extracted ion chromatogram, Dex-P: Dexamethasone-21-phosphate.

Export to PPT

To broaden our screening efforts, we tested 50 additional batches of actual samples to detect the illegal addition of Dex-P, using the validated results as a reference. In the target screening, Dex-P was detected in four samples, with concentrations varying from 203.35 to 27,810.58 ng/g (Figure 7c-f). The National Medical Products Administration (NMPA) clearly stipulated that hormones, including but not limited to anti-androgen substances of steroid structure, glucocorticoids (corticosteroids), estrogens, palamisone, progesterone, substances with androgen effect, etc. were prohibited in cosmetics. In addition, the cosmetic matrices involved in the actual samples of this study include repair creams, light printing creams, pore-unblocking liquids, etc. Woo et al. applied a similar approach to eyebrow tattoos, lipstick tattoos, and hair tint samples [21]. Furthermore, our previous study utilized the FBMN method to screen various actual samples of creams and waters, and screened for antibiotics, antihistamines, phthalates, etc., and these results demonstrated the feasibility of our method for a wide range of cosmetic matrices and repeated regulatory use [53]. Therefore, the use of molecular networking technology for suspicious screening can significantly enhance our efforts to monitor and identify more risk sources in cosmetics in the future.