Sensitive Metabolites and Single Nuclear Polymorphisms of OCT1 on Imatinib Meyslate in Patients with Gastrointestinal Stromal Tumors (GISTs)

Imatinib mesylate (IM) is highly ecacious in the treatment of gastrointestinal stromal tumors (GISTs). Therapeutic drug monitoring (TDM) and pharmacogenetic screening are armed for better management of IM therapy. The goal of this study was to explore the most sensitive metabolites and the role of organic cation transporter 1 (OCT1) on single nuclear polymorphisms (SNP) of IM treatment, which could supply a greater mechanistic understanding of therapeutic effect or resistance of IM. A total of 40 human serum samples from patients with GISTs were collected for TDM. Basing on the results of TDM, the untargeted metabolomic analysis was to determine characteristics of the serum. Principal component analysis (PCA), orthogonal partial least-squares discrimination analysis (OPLS-DA), and heat map were used for multivariate analysis. In addition, KEGG database were used to identify the pathways of metabolites. Differential metabolites were identied based on a heat map with a t-test threshold (p < 0.05), fold-change (FC) threshold (FC > 1.5 or FC < 2/3) and variable importance (VIP) in the projection (VIP>1). The potential differential metabolites were included D-sphingosine, 1-sphingosine phosphate, phytosphingosine and phosphoethanolamine, which were belonged to the metabolism pathway of sphingolipid. Because of the function of OCT1 transporters was related with drug in the liver, the OCT1 1386C>A (rs 622342) was statistically signicant (P<0.05) with IM plasma concentration. Thus, this study demonstrated that sphingolipid metabolism should be considered as the potential pathway of IM treated GISTs, which could bring us a clue to new mechanism for IM treatment of patients with GISTs. GISTs:Gastrointestinal tumors; TDM:Therapeutic drug monitoring; OCT1:The role of organic cation 1; SNP:Single nuclear polymorphisms; PCA:Principal component analysis; OPLS-DA:Orthogonal partial least-squares discrimination analysis; KEGG:Kyoto Encyclopedia of Genes and VIP:Variable Importance for the TKIs:Tyrosinase inhibitors; KIT:Tyrosine kinase; PDGFRα:Platelet-derived growth factor receptor α; CYP450:Cytochrome OATP1:Organic transporter protein 1; NCCN:The National Comprehensive Cancer Network; cluster analysis; ionization; range; ; SPK1:Sphingosine kinase-1;


Background
Gastrointestinal stromal tumors (GISTs) is the most common mesenchymal tumor of Gastrointestinal tract, and has the potential of malignant differentiation of recurrence and metastasis [1]. Surgical resection is the main treating means of GISTs [2], while oral targeted small molecule drugs are mainly used for adjuvant therapy of GISTs with condition of high recurrence and metastasis rate [3]. Imatinib mesylate (IM) is the rst generation of tyrosinase inhibitors (TKIs) that have been applied clinically, and it has made remarkable achievements as an oral small-molecule targeted drug for the rst-line treatment of completely resected and high-risk GIST patients as well as unresectable, relapsed, and metastatic or advanced GIST patients [4]. The problem of IM in the treatment of GIST is that nearly 85% of patients develop secondary drug resistance after 1-3 years of IM treatment [5]. At present, the mechanism of IM resistance mainly focused the following three aspects: i) acquired drug resistance caused by secondary mutations of receptor tyrosine kinase (KIT) or platelet-derived growth factor receptor (PDGFRα). Since IM is an inhibitor for targeting KIT and PDGFRα receptors, mutations in these two sites would be directly determined the sensitivity of IM therapy [6,7]. ii) IM resistance caused by the in uence of enzymes related to drug metabolism. Because IM is mainly metabolized by liver with the oral bioavailability of 98%, effects of liver metabolism enzymes, such as cytochrome P450 (CYP450) [8], and drug transport enzymes in the liver, such as coding of P-glycoprotein ABCB1 [9], coding of breast cancer resistance protein ABCG2 [10], coding organic cation transporter 1 (OCT1) SLC22A1 [4] and coding organic anion transporter protein 1 (OATP1) SLCO1A2 [11], and other sites gene polymorphism would affect IM concentration in serum actually. The secondary resistance to IM was related with insu cient serum concentration. iii) The last mechanism of IM resistance to GISTs patients was veri ed due to activation of other alternative survival pathways [12].
As an oral small-molecule targeted drug, IM has huge intra-individual differences (about 75%) and inter-individual differences (60%) [13], and it was pointed out that the determination of plasma drug concentration of IM is of great clinical signi cance for patients with metastatic, relapsed and unresectable GIST in the National Comprehensive Cancer Network on the Diagnosis and Treatment of GISTs (NCCN, 2020 edition) [14]. Plasma drug concentration refers to the drug concentration in plasma, including free drug concentration and drug concentration bound to plasma protein [15]. At present, the grain concentration C min of 1100µg/L was considered as an effective therapeutic concentration in the world consensus [16]. The metabolism correlated gene polymorphism has a particularly signi cant in uence on IM serum concentration reported by our team [4].
Combined with the research results on the mechanism of secondary drug resistance of IM, it was speculated that it might be related to the metabolic signaling pathway of IM. Therefore, we grouped the GISTs patients according to the serum concentration of IM and analysis the metabolites of the different groups, which was aimed to study the potential IM resistance mechanism with liver metabolism. In this study, we conducted the nontargeted metabolomics assay, which is one of the main means of metabolomic research and could be comprehensively and systematically analysis the metabolites of all small molecular weight ( 1500 Da) [17]. The classically discriminative pattern analysis, were principal component analysis (PCA), hierarchical cluster analysis (HCA) and orthogonal partial least squares discriminant analysis (OPLS-DA) [18]. After a series of screening of differential compounds, the metabolic pathway analysis with the help of bioinformatic tools was consequently carried out [19]. The utility of bioinformatic tools could automate the highlighting of distinct metabolite features from different groups of samples, which contributed great efforts to identify the early subclinical markers [20].
Then the gene polymorphisms of OCT1, which is correlated with the function of drug transportation in the liver [21], were carried out in this study. As a member of the solute transporter superfamily [22], OCT1 transporter is mainly expressed in the basolateral membrane of hepatocytes and plays an important role in drug uptake into the liver. Its main function is to transport drugs into hepatocytes for metabolism or to pump drugs into intestinal epithelial cells [23]. The OCT1 transporter encoding gene (SLC22A1) showed signi cant polymorphism [24,25]. A sequencing analysis of 2171 unrelated individuals from 67 world populations reported 29 variants leading to amino acid substitution, which could be divided into 30 haplotype variants consisting of 16 alleles, most of which had signi cant effects on OCT1 transporter activity [26]. Polymorphic expression of the encoding OCT1 gene can cause partial loss of the function of transporters, resulting in reduced or loss of transport e ciency, thus affecting the blood concentration of IM and further affecting the therapeutic e cacy, and more likely to lead to drug resistance of IM caused by insu cient effective dose [27,28]. In view of this, combined with previous studies and polymorphism functional variability locus retrieval results, we selected OCT1 1759 G > A (rs 6935207), OCT1 201 C > G (rs5812592), OCT1 1386 C > A(rs 622342), OCT1 1022 C > T (rs2282143) and OCT11222 A > G (rs628031) to study their correlation with the plasma concentration of IM in GIST patients, so as to provide potential evidence for subclinical markers of GIST patients to take IM with OCT1 SNP. Figure 1 showed the therapeutic drug monitoring (TDM) and pharmacogenetic screening, especially untargeted metabolomic analysis after IM treatment.

Ethical Considerations
This study was conducted in accordance with the principles of the Declaration of Helsinki, and the study protocol was approved by the Research Ethics Committee of the Sichuan Cancer Hospital, (Approval No. SCCHE-02-2020-044). The written informed consent was obtained to all individuals before participating in this study.

Study Participants
All plasma samples were collected from 29 GIST patients who took IM regularly and 11 GIST patients who had never taken IM collected from 2019 to 2020 in Sichuan Cancer Hospital. The inclusion criteria for the GIST group were depicted as: (i) male or female, aged ≥ 18; (ii) histopathological diagnosis of GIST; (iii) Continuous regular oral administration of IM 400 mg/d for more than 1-month, good compliance, no missed or less drug phenomenon; (iv) the patient's tumor growth site was unlimited, regardless of whether they underwent surgical resection. The exclusion criteria were commanded as: (i) pregnant women, very poor mental and physical conditions; (ii) taking IM and other targeted drugs at the same time; (iii) combined with other organ failure and infection; (iv) poor compliance, expected survival < 3 months; (v) communication and cognitive impairment.

Plasma Sample Collection
Patients taking IM need to take the medicine regularly for more than 1 month to maintain steady state of plasma concentration. Before collection of blood, patients would be instructed to take IM about 12:00. Then blood was collected about 10:00 the next morning. The collected peripheral venous blood was placed in an anticoagulation tube. At the same time, the basic information of the GIST patients was clearly speci ed, including name, gender, age, whether to take IM, etc. All samples were divided into GIST group (D), who had never been taken IM; and GIST group (N) who took IM regularly. According to the results of IM plasma concentration, the GIST group N was classi ed into 3 groups, that were, IM high concentration group (A), IM medium concentration group (B) and IM low concentration group (C). Each sample was taken from 1 mL of whole blood and stored at -80°C for non-targeted metabolomics analysis, and multiple freeze-thaw cycles should be avoided. The remaining whole blood was centrifuged at 5000 rpm for 5 min to obtain plasma, and then stored at -20°C until IM plasma concentration analysis.

Sample Preparation for Plasma drug Concentration Testing
For plasma IM concentration analysis, we added 1000 µL of organic treating agent (ORG-1) to 400 µL of each plasma sample for protein precipitation, then vortexed for 1 min and centrifuged at 145000 rpm for 8 min (TG-16, Sichuan Shuke Instrument Co., Ltd., China). The supernatants were transferred into sampling bottles, and the plasma concentration of IM was detected by Two-dimensional liquid chromatograph(2D-LC-UV) FLC-2701, MLC2420, Hunan Demeter Instrument Co., Ltd., China .
The operating column temperature was 40°C with detection wavelength at 264nm.

Sample Preparation for Untargeted Metabolic Pro ling
Plasma samples for the untargeted metabolic pro ling were prepared with 300 µL of methanol, which was contained 1 ppm 2chlorophenylalanine. Then it was added for protein precipitation, vortexed for 2 min and incubated in a − 20°C freezer for 0.5 h.
Consequently, samples were vortexed for 2 min and centrifuged at 12000 r/min for 10 min at 4°C. The supernatant was transferred to a centrifuge tube and then centrifuged at 12000 r/min for 15 min at 4°C. The supernatant was transferred into sampling bottles for QTOF/MS-6545 mass spectrometer and 1290 In nity LC ultra-high performance liquid chromatograph were supplied by Aglient)analysis.

HUPLC-MS Instrumentation and Operating Condition
The HUPLC-MS system consisted of a ABSCIEX Q TRAP 3200 mass spectrometer (ABSCIEX, Germany) equipped with an electrospray ionization (ESI) interface coupled to 1290 In nity LC ultra-high performance liquid chromatograph were supplied by Aglient with a quaternary gradient pump (Agilent 1,260 in nity). The gradient elution system consists of solvent A (water containing 0.04% acetic acid) and solvent B (acetonitrile containing 0.04% acetic acid), and the control. The column temperature was maintained at 40°C with the ow rate at 0.35 mL/min. Subsequently, the mass spectrometer was in positive and negative ion modes for accurate mass measurement, respectively. ProteoWizard software was conducted to analyze the whole data.

Pharmacogenetic Analysis
The peripheral venous blood of the patients is collected as the similar assay as depicted in the Sect. 5.3. The ve candidate SNPs (OCT1 1795 G > A(rs 6935207), OCT1 201 C > G(rs58812592), OCT1 1386 C > A(rs 622342), OCT1 1022 C > T(rs2282143), OCT11222 A > G(rs628031)) were analyzed via rst generation gene sequencing. The experimental instruments, reagent consumables and the primers used for each SNP were supplied in Supplementary data S4 and S5. The speci c process of PCR was as usual report [29,30].

Data Analysis
The metabolite's data were converted to mzML format using ProteoWizard software. Peak areas were corrected using the "SVR" method and peaks with > 50% missingness in each set of samples were ltered. The screened peaks were corrected, by searching laboratory self-built databases, integrating public libraries, which contained with the Human Metabolome Database (www.hmdb.ca), KEGG network (www.genome.jp), PubChem (www.pubchem.ncbi.nlm.nih.gov), Small Molecule Pathway Database (www.smpdb.ca), and LIPID MAPS Lipidomics Gateway (www.lipidmaps.org). Finally, statistical analysis was performed by the R program. Statistical analysis was divided into univariate statistical analysis and multivariate statistical analysis, univariate statistical analysis included Student'-test and fold difference analysis, and multivariate statistical analysis included principal component (PCA) analysis, partial least squares discrimination (PLS-DA) analysis and orthogonal partial least squares discrimination (OPLS-DA) analysis.
The pharmacogenetic data were performed using SPSS 25.0 software. The Shapiro-Wilk method was used for normality test. If it conformed to the normal distribution (P > 0.05), the data were represented by mean (M) ± standard deviation (SD), and used independent sample t-test for comparison between the two groups. If not (P < 0.05), the data were represented by median and interquartile range (IQR), and used the Mann-Whitney U test for comparison between the two groups. The cutoff for statistical signi cance was set at P < 0.05 (two-sided). Fisher's exact probability method was used to test whether the frequency distribution of genotypes at each point complies with the Hardy-Weinberg equilibrium.

Patient Characteristics and Groups
This study included 29 GISTs patients who regularly took IM orally and 11 GISTs patients who did not take IM. The included sample situation was displayed in Table 1. The average age of GISTs patients with IM or without IM were respectively 57.30 ± 8.42 and 54.12 ± 6.78, which demonstrated no signi cant difference. Similarly, the height, weight and sex were also recorded. The GISTs patients were divided into 4 groups, according to with or without IM treatment. Speci cally, the rst three groups (A, B and C) were received the IM 400mg/qd treatment. In detail, group A represented the IM plasma concentration 2000 ng/mL (n = 9); group B represented the IM plasma concentration ranging from 2000 ng/mL to 1100 ng/mL (n = 10) and group C represented the IM plasma concentration 1100ng/mL (n = 10). The last group D (n = 11) was a representative of IM missed patients. As shown in Supplementary S1, the difference between the highest IM concentration to the lowest IM concentration was more than 12 folds. This consequence re ected the huge individual differences when GISTs patients receiving the treatment of IM.
3.2. In uencing factors of gender, age and surgical operation with mean plasma concentration of IM According to our previous literature analysis [4], the gender, age and surgical operation could be some controversial factors to in uence the mean plasma concentration of IM. Herein, in this study, we conducted those factors to validate again.
Firstly, the plasma concentration of IM was detected. Consequently, the Shapiro-Wilk test for normality [31], Mann-Whitney U test for evaluating treatment effects in randomized elements [32], and independent t test for verifying statistically signi cant difference [33] were carried out to analyze the factors, which might be related with the mean plasma concentration of IM. The speci c statistical data analysis of gender, age and surgical operation with C min of IM were supplied in the Supplementary data S1 to S3. As displayed in the Fig. 2A to C, the C min of IM in the female patients was higher than that of male patients, which was accordance with the report from a team of United Kingdom [34] (n = 93) and Switzerland [35] (n = 2478). The potential elements should be drawing to the mean body weight and clearance rate. As long-term results from the report of United Kingdom [34], when plasma concentration of IM was normalized for body weight, the differences in mean concentrations were no longer apparent. Herein, this consequence suggested that the higher plasma concentrations of IM in women could be partially explained by the lower in body weight, when compared with men. Furthermore, the clearance rate of women was higher than man from the report of Switzerland [35]. Due to these two reasons, the plasma concentration of IM was higher in women group. The age was demonstrated no signi cantly difference in our study, which was different with the Gotta et al. [35].
The Gotta et al. de ned "young" as < 30 years of age and "elders" as people up to 70 years of age. The young group exhibited low IM concentrations compared with those of the elders (p < 0.05). However, rstly, the median age in our study was 57.30 ± 8.42. Speci cally, the age of included GIST patients was from 40 to 78 years old. There were none of samples below 30 years old. Then, the sample capacity of our study was insu cient, which was limit by the COVID-19 [36]. As listed in the Fig. 2B, the surgical operation played important role in the plasma concentration of IM. The difference between gastric surgery and nongastric surgery group was statistically signi cant. Moreover, partial gastric excision and total gastric removal were statistically signi cant in the gastric surgery group. Speci cally, the concentration of IM in the total gastric surgery was lower than that of in partial gastric resection patients. Similarly, Yoo et al. and Hompland et al. [37,38] reported that C min was signi cantly lower in patients who had previously undergone major gastrectomy than in those with previous wedge gastric resection or without gastric surgery. Obviously, decreased absorption of IM might be caused in part by the lack of gastric acid secretion in GISTs patients who had undergone major gastric resection. The gastric acid secretion was extremely important for IM absorption, because of IM tablets dissolved rapidly at pH 5.5 or less [39].

Correlation among mean plasma concentration of IM min with body surface area and some essential serum index
Linear mixed model analyses could display similar trends, although sometimes without reaching statistical signi cance. The body surface area was demonstrated by several reports [37,[40][41][42], which should be related with plasma concentration of IM.
As displayed in the Fig. 2D to I, the body surface area was exhibited negative correlation with plasma concentration of IM (P = 0.043, r 2 = 0.129). In yet other words, the higher dosage of IM should be given to the patients with smaller body surface area. Combined with ours' and others' consequence about the in uence of body surface area, the xed dosage of IM treated on the GISTs patients might not be suitable.
The laboratory data of ALT, AST, creatinine (CRE), total bilirubin (TBil), and albumin was also taken into consideration in this study. As demonstrated in the Fig. 2E to I, there was no signi cant correlation between ALT (P = 0.747, r 2 = 0.003), AST (P = 0.667, r 2 = 0.006), CRE (P = 0.905, r 2 = 0.019), Tbil (P = 0.326, r 2 = 8.468 e − 4 ) and albumin (P = 0.583, r 2 = 0.017) with plasma concentration of IM. The consequence of ALT and AST was consistent with the report of Yoo et al [37], while the result of CRE and Tbil were inconsistent with [40][41][42]. The reason should be associated with at least two hands. On the one hand, the number of samples among these above-mentioned studies was ranged from 25 to 89 GISTs patients. As we all known, the more included eligible samples, the statistics of data was more reliable. However, among these studies, the included samples were not very satis ed. On the other hand, the different ethnic population should not be ignored. The association of albumin with IM min was controversial. It could be suggested that more IM would be bound to albumin in patients with higher albumin, resulting in higher total IM min from the demonstration of Yoo [37], when the albumin exhibited active correlation with IM min .
Sometimes, the albumin was without reaching statistical signi cance with IM min [38,43]. The actual relationship of albumin with plasma concentration of IM min was still unclear. Herein, not only large scale and multi-center trial of plasma concentration of IM related with in uencing factors needed to be applicable, but also the mechanism behind the association should be explored, which would be promoted the rational usage of IM treating on patients with GISTs.

Principal Component Analysis (PCA)
PCA was belonged to unsupervised analysis assay. The advantage of PCA was that it could reduce the number of highly correlated metabolic features to a smaller set of principal components. This superiority made the PCA scores plots provide a visual description of the pattern described by the model that can be used for the identi cation of batch effects [44]. The results of PCA were illustrated in the Fig. 3.
The IM treated groups (A, B and C) could be distinguished from none-treated group D, especially in the ESI + mode, no matter in 2D or 3D plot. In the mode of ESI-, IM treated groups (A and B) were also tightly centralization when compared with group D.
These consequences encouraged us for further statistical analysis.

Orthogonal partial least squares discriminant analysis (OPLS-DA)
The problem of insensitivity of variables with small correlation could be overcame by OPLS-DA model, which was the de ciency of PCA. Herein, we constructed the OPLS-DA model for the next data management. Firstly, we evaluate reliable and predictive ability of OPLS-DA model, which was supplied in the Supplementary S2. The parameters, R 2 X, R 2 Y and Q 2 were all meet the standard. When the R 2 X, R 2 Y and Q 2 were closer to 1, the more suitable and reliable for this model [45]. Meanwhile, the Q 2 > 0.5 is considered as an effective model, Q 2 > 0.9 belongs to an excellent model [46].
As exhibited in the Fig. 4, no matter in the mode of ESI-or ESI+, the Q 2 was always more than 0.8, which indicated the models were all successful and biochemical changes between groups were clear. The abscissa of OPLS-DA plot represents the scores of major components. Therefore, the difference between groups could be judged from the direction of the abscissa. There was a clear separation between the IM treated groups (A, B and C) and group D. Concurrently, the samples from both groups tended to cluster in a concentrated manner, with a high degree of aggregation.

Differential metabolites screening
To distinguish the most important metabolites between the groups, we used the above analysis methods to digit the most potential differential metabolites.

Volcano analysis
After the OPLS-DA analysis, the speci c numbers of statistical metabolites were conducted via volcano analysis. The variable importance in projection (VIP) of the OPLS-DA model can be used to preliminarily screen out the metabolites that differ among different groups. As exhibited in Fig. 5, the bigger the plot is, the value of VIP is higher. The value of VIP represents the in uence intensity of corresponding metabolite differences in the classi cation discrimination of each group of samples in the model. Generally, the value of VIP ≥ 1 metabolite is considered to be signi cant metabolite. The red plots represent the signi cant upregulation metabolites, while the green plots suggest the signi cant down-regulation metabolites. The rest, grey plots, recommend the insigni cant metabolites.
Also, according to the Fig. 5, the total numbers of metabolites were 1020 and 2238 in the mode of ESI-and ESI+, respectively. This result suggested that the ESI-was less sensitive than ESI + when collecting IM metabolites. Herein, the mode of ESI + was carried out as the next analytical method.
It can be seen from Fig. 5B, D and F that there were totally 168 and 195 down-regulation metabolites, when compared with IM none-treatment group, respectively. Meanwhile, three were 125, 113 and 84 signi cant up-regulation metabolites, when compared with IM none-treatment group.

Cluster heat map analysis
Heat maps are one of the most widely used biological informatic graphic expression. Heat map is the representation of the metabolomics dataset with hierarchical clustering analysis (HCA) [47]. The graphical depictions of the metabolites are altered signi cantly across the different groups in the heat map in different colors. The heat map takes differential metabolites as the vertical axis, while the samples of model group or control group as the horizontal axis. The content of each differential metabolite in each sample could be directly seen according to the color depth. The darker the color is, the higher the content is [48]. As illustrated in Fig. 6A to C, the metabolites in the model group (A, B and C) and control group exhibited different behaviors (colors), which suggested these differential metabolites with good discriminating ability.

Pearson correlation coe cient
Pearson correlation coe cient was used to evaluate whether the relationship between the two groups was linear or nonlinear, or in other words, the strength of correlations [49]. Notably, 0 < Pearson correlation coe cient ≤ + 1 represents a positive correlation, while − 1 ≤ Pearson correlation coe cient < 0 represents a negative correlation. The greater the absolute value, the stronger the correlation [50]. We screened the top 50 differential compounds with the highest VIP value for Pearson correlation analysis. The each of top 50 compounds was calculated via Pearson correlation coe cient. Then a matrix graph was obtained, as displayed in the Fig. 6D-F, which displayed that red color represents a positive correlation between the two groups, while green represents a negative correlation between the two groups. The consequence of correlation would supply the information for narrowing the huge number of metabolites in order to dig the most potential differential metabolites.

Fold change
After qualitative analysis of differential metabolites, combined with the grouping situation of speci c samples, the difference fold changes (FC) in the quantitative information of metabolites in each group were compared. The FC value represented the difference between the two groups and was calculated as the average of each individual peak area: (Mean value of peak area obtained from IM patients)/(mean value of peak area obtained from none IM treated patients). If the FC value was greater than As summarized in Table 2, there were glucose (P0145, P0918, P0925), securinine (P0635), carnitine (Carnitine C18:1(P0077), Carnitine C18:2 (P0079), Carnitine C18:3 (P0082) and mixture of those three (P0851, P0889), D-sphingosine (P0881), 1sphingosine phosphate (P0927), phytosphingosine (P0883) and phosphoethanolamine (P0448) could be speculated into corresponding metabolism pathway, which we would like analyze with the help of KEGG pathway enrichment method as follow.
In this section, the distribution and probability density of top 20 differential samples in the Table 2 were inspected by violin box [53]. The thin black line extending from violin box represents the 95% con dence interval, the black bar in the middle of violin box represents the median value, and the outer shape represents the distribution density of the samples. As depicted in the Fig. 7A to C, the differential metabolites in the group D to A, B and C, exhibited different distribution patterns, which indicated the capacity of distinguish of these metabolites in the IM none treating or treating groups.

KEGG pathway enrichment analysis
KEGG (http://www.kegg.jp) is an encyclopedia of not only genes and genomes but also a professional tool for metabolites and nonmetabolic pathways [54,55]. Pathways were considered signi cantly enriched, if p < 0.05, impact number of metabolite hits in the pathway > 1 [56]. The greater number of metabolite hits in the pathway with the lower p value, the more matching the signi cant pathway. The most in uenced metabolic pathway was set as a pathway in uence cut off value > 0.1 to lter for less important pathways. As displayed in Fig. 7D-F, the signal pathway as described included sphingolipid metabolism, drug metabolism (CYP 450), butanoate, caffeine metabolism and so on. The sphingolipid metabolism was with greatest commonness among those three groups, which impact numbers were above 1.75 with p value below 0.05. What follows is the metabolism map of sphingolipid metabolism, as described in Fig. 8G. The key metabolites such as D-sphingosine (P0881), 1sphingosine phosphate (P0927), phytosphingosine (P0883) and phosphoethanolamine (P0448) were belonged to the metabolism pathway of sphingolipid, which was also existed in the Table 2. 3.8 Genotype of SNP Five SNP (OCT1 1795 G > A(rs 6935207), OCT1 201 C > G(rs58812592), OCT1 1386 C > A(rs 622342), OCT1 1022 C > T(rs2282143), OCT11222 A > G(rs628031)) of 29 patients were performed by rst-generation gene sequencing, and the typing results were shown in Fig. 8. As displayed in the Fig. 8, the results of genetic test demonstrated that the genotypes at OCT1 201C > G (rs58812592) were all CC (wild-type). The other four candidate SNPs of the OCT1 were with mutations. he genotypic distribution of the four SNPs were summarized in the Table 3. Herein, only the remaining four SNPs of genotypes were analyzed by Hardy-Weinberg equilibrium. Consequently, none of the genotypes showed deviation from the Hardy-Weinberg equilibrium [57] (P > 0.05). Then the relationship between genotypic SNP and IM plasma concentration were conducted through Mann-Whitney U test, which exhibited in the Table 4. Two genotypes at OCT1 1386C > A (rs 622342) was statistically signi cant (P < 0.05) with IM plasma concentration. The consequence of 1022 C > T(rs2282143)and 1222 A > G rs628031 were similar with the report from Francis [58], which included the 73 patients of IM treated chronic myeloid leukemia patients.

Discussion
There is a long journey for IM treating on the patients with GISTs. However, the IM-resistance is the obstacle for trapping the further development of IM. Speci cally, more than 50% of patients with advanced and metastatic forms of the disease develop secondary resistance to IM-based therapy within 2 years after initiation of treatment [59]. It is widely known that IM is a successfully orally targeted drug, which metabolized by CYP 450 [60]. The interpatient differences of pharmacokinetic (PK) variability have been estimated to be about 60% in IM steady-state trough concentrations in patients with GISTs [61].
Basing on our previous literature report [4], surgery operation, some serum index and SNP distribution of OCT1 could play some potential role in disturbing the steady-state plasma concentration of IM. What's more, with the aim of digging the most potential metabolite of IM treating on GISTs, we divided the blood samples into three parts, those were, detecting for IM plasma concentration, IM plasma metabolites and SNP distribution of OCT1, respectively. The ultimate purpose of this study is to nd out the key factors affecting the clinical e cacy of IM and dig out the sensitive metabolites of IM in the treatment of GIST, so as to provide evidence for the later study of the mechanism of IM resistance.
Firstly, the gastric surgery operation and body surface area were associated with IM steady-state concentration. The consequence was in accordance with some reports [37,38]. With the help of untargeted metabolomics assay, we found nearly 2000 metabolites in the mode of ESI + from the HUPLC-MS detector. Through a series of metabolomics analysis, the sphingolipids metabolism was selected as the most potential pathway, which could be probably used to discriminate IM from treating or not. Beyond that, the sphingolipids metabolism pathway could bring us a clue to new mechanism for IM-resistance of patients with GISTs. .Therefore, the dynamic balance of CER/S1P determines whether tumor cells go to apoptosis or proliferation. Sphingosine kinase-1 (SPK1) generates S1P by phosphorylating SP, a product of CER, and is a key enzyme for S1P synthesis, as well as regulating CER /S1P homeostasis [67].Studies have shown that the overexpression of SPK1 can promote the expression of PI3K/Akt signaling pathway in triple-negative breast cancer, thus promoting the proliferation and metastasis of triple-negative breast cancer cells, which is the potential mechanism of tumor drug resistance [68].When activated by phosphorylation, SPK1 accelerates the synthesis of S1P on the cell membrane. As the second messenger, the synthetic S1P participates in the S1P receptor response in an autocrine or paracrine manner, including stimulating a series of tumor biological functions, such as inhibition of apoptosis, promotion of autophagy, promotion of proliferation, promotion of migration and invasion, and promotion of angiogenesis and in ammation, thus promoting the further occurrence of tumor drug resistance [69].The combination of SPK1 and bortezomib can also reverse the sensitivity of IM-resistant chronic myelogenous leukemia (K562) cells by down-regulating Mcl-1 expression [70].Although the role of sphingolipid metabolism pathway on the GISTs patients still remained to be unclear, we could speculate the potential association of sphingolipid metabolism pathway with IM on the treatment of patients with GISTs according to results in this study.
After the nding of sphingolipid metabolism pathway, the SNP distribution of OCT1 was veri ed, because of the function of OCT1 was considered to be the main in ux transporter involved in IM uptake into chronic myelogenous leukemia cells [71]. Lower OCT1 activity in diagnostic chronic myeloid leukemia blood mononuclear cells, the poorer molecular response was occured to IM [28]. The molecular response of IM was usual correlated with plasma concentration of IM [28]. Studies have also shown that drug resistance is easily caused when the effective blood drug concentration is lower for a long time [72]. Due to the actual role of OCT1 in GIST cells was still unclear, we chose 5 SNP of OCT1 in this study according to some reports. The consequence was demonstrated that only 1386C > A(rs 622342)was signi cantly difference in IM plasma concentration. The weakness of this study lies in insu cient samples. The more samples, the stronger the statistical reliability was performed in the study of SNP. Herein, this consequence should be recommended as a preliminary experiment for later research.

Conclusion
In summary, we performed TDM of IM with an untargeted metabolomics investigation to improve our understanding of IM metabolism on GISTs patients, which found D-sphingosine, 1-sphingosine phosphate, phytosphingosine and phosphoethanolamine were selected as candidate IM treating biomarkers for GISTs, and sphingolipid pathway should be a potential pathway to in uence the effect of IM in the treatment of GISTs. The 1386C > A(rs 622342)of OCT1 was signi cantly related with IM plasma concentration, while it also needed to be more samples to con rm.

Availability of data and materials
All data generated or analyzed during this study are included in this published article.

Declarations
Ethics approval and consent to participate All participants provided informed consent to the protocol (Approval No. SCCHE-02-2020-044) was approved by the Review Ethic Board at the Sichuan Cancer Hospital.

Consent for publication
Not applicable.

Con icts of interest
The authors have no con icts of interest to declare.   OPLS-DA score plots of each group in both negative (ESI-) and positive (ESI+) modes. A and B were represented that OPLS-DA score plot of group D compared with group A under the ESI-or ESI+ mode. C and D were suggested that group D compared with group B under the ESI-or ESI+ mode. E and F were indicated that Group D compared with group C under the ESI-or ESI+ mode.

Figure 5
Volcano plots of each group in both negative (ESI-) and positive (ESI+) modes. A and B were represented that volcano plot of group D compared with group A under the ESI-or ESI+ mode. C and D were suggested that group D compared with group B under the ESI-or ESI+ mode. E and F were indicated that Group D compared with group C under the ESI-or ESI+ mode. The abscissa represents the multiplier value of the metabolite difference, which chose Log base 2. The ordinate represents the t-test P value of univariate analysis, which chose Log base 10.
Page 24/26  The violin box, KEGG enrichment and sphingolipid metabolism map analysis. A, B and C were belonged to the violin box analysis of group D compared to group A, B and C. The vertical position of each histogram represents the raw intensity. The thin black line extending from violin box represents the 95% con dence interval, the black bar in the middle of violin box represents the median value, and the outer shape represents the distribution density of the samples. D, E and F exhibited the metabolism view of pathway impact analysis obtained from differential metabolites in group D comparison with group A to C.
The color and size of each circle is based on p-values (green: higher p-values and red: lower p-values) and pathway impact values (the larger the circle the higher the impact score) calculated from the topological analysis, respectively. Pathways were considered signi cantly enriched, if P < 0.05, impact 0.1 and number of metabolite hits in the pathway >1 [54]. G represented the sphingolipid pathway map, which should be one of most potential metabolic pathways to distinguish IM treatment or not.