Identification of Metabolism-Related Proteins as Biomarkers of Insulin Resistance and Potential Mechanisms of m6A Modification

Background: Insulin resistance (IR) is a major contributing factor to the pathogenesis of metabolic syndrome and type 2 diabetes mellitus (T2D). Adipocyte metabolism is known to play a crucial role in IR. Therefore, the aims of this study were to identify metabolism-related proteins that could be used as potential biomarkers of IR and to investigate the role of N6-methyladenosine (m6A) modification in the pathogenesis of this condition. Methods: RNA-seq data on human adipose tissue were retrieved from the Gene Expression Omnibus database. The differentially expressed genes of metabolism-related proteins (MP-DEGs) were screened using protein annotation databases. Biological function and pathway annotations of the MP-DEGs were performed through Gene Ontology and Kyoto Encyclopedia of Genes and Genomes pathway analyses. Key MP-DEGs were screened, and a protein–protein interaction (PPI) network was constructed using STRING, Cytoscape, MCODE, and CytoHubba. LASSO regression analysis was used to select primary hub genes, and their clinical performance was assessed using receiver operating characteristic (ROC) curves. The expression of key MP-DEGs and their relationship with m6A modification were further verified in adipose tissue samples collected from healthy individuals and patients with IR. Results: In total, 69 MP-DEGs were screened and annotated to be enriched in pathways related to hormone metabolism, low-density lipoprotein particle and carboxylic acid transmembrane transporter activity, insulin signaling, and AMPK signaling. The MP-DEG PPI network comprised 69 nodes and 72 edges, from which 10 hub genes (FASN, GCK, FGR, FBP1, GYS2, PNPLA3, MOGAT1, SLC27A2, PNPLA3, and ELOVL6) were identified. FASN was chosen as the key gene because it had the highest maximal clique centrality (MCC) score. GCK, FBP1, and FGR were selected as primary genes by LASSO analysis. According to the ROC curves, GCK, FBP1, FGR, and FASN could be used as potential biomarkers to detect IR with good sensitivity and accuracy (AUC = 0.80, 95% CI: 0.67–0.94; AUC = 0.86, 95% CI: 0.74–0.94; AUC = 0.83, 95% CI: 0.64–0.92; AUC = 0.78, 95% CI: 0.64–0.92). The expression of FASN, GCK, FBP1, and FGR was significantly correlated with that of IGF2BP3, FTO, EIF3A, WTAP, METTL16, and LRPPRC (p < 0.05). In validation clinical samples, the FASN was moderately effective for detecting IR (AUC = 0.78, 95% CI: 0.69–0.80), and its expression was positively correlated with the methylation levels of FASN (r = 0.359, p = 0.001). Conclusion: Metabolism-related proteins play critical roles in IR. Moreover, FASN and GCK are potential biomarkers of IR and may be involved in the development of T2D via their m6A modification. These findings offer reliable biomarkers for the early detection of T2D and promising therapeutic targets.


Introduction
Insulin resistance (IR) is a significant pathophysiological basis and necessary stage of type 2 diabetes (T2D) as well as a critical pathogenic component of the metabolic syndrome [1]. The condition may develop as early as 13 years prior to the clinical diagnosis of T2D [2]. Therefore, the precise identification of IR and the provision of appropriate interventions for high-risk groups could effectively delay or prevent its development.
The pathogenesis of IR involves the attenuation or blocking of insulin signal transition as well as impaired glycolipid, lipid, and energy metabolism [2,3]. Under normal body conditions, elevated glucose levels in the body stimulate insulin secretion by islet β-cells. The insulin molecules then circulate to the target tissue via the blood system and bind to the insulin receptors on the tissue cell membrane. Next, a series of specific signals trigger the cells to produce biological effects that promote glucose intake and use, thereby reducing the blood glucose level. Importantly, the phosphoinositide 3-kinase/protein kinase B (PI3K/Akt) signaling pathway is a key pathway for insulin signal transduction. T2D occurs when there is insufficient insulin secretion to compensate for IR, thereby disrupting signal transduction in pathways involved in the uptake and usage of insulin by target tissues [2,4]. The effector organs of IR are the liver, skeletal muscles, and adipose tissues. In the liver, IR is characterized by increased glucose production and impaired glucose uptake [5]. By contrast, reduced glucose uptake is the defining feature of IR in adipose tissue and skeletal muscles [4,6], the latter of which are the primary sites of glucose disposal. Additionally, IR occurs in adipocytes, whereas muscle IR-KO mice are typically healthy.
Adipocytes play a crucial role in the regulation of mammalian metabolic homeostasis. These cells have a significant impact on energy expenditure, glucose and lipid balance, and IR development [5,6]. Adipocyte dysfunction causes the dysregulation of adipokines, which may influence inflammatory responses, glucose and lipid homeostasis, and other metabolic illnesses, including diabetes [6][7][8]. Adiponectin and leptin, which are mainly expressed in adipose tissue, can promote energy consumption and metabolism, reduce fat accumulation, and improve insulin sensitivity [9,10]. Although fibroblast growth factor 21 (FGF21) is known to regulate glucose metabolism [11], the molecular mechanisms underlying IR are not fully understood. Identifying biomarkers for the early and accurate diagnosis and prevention of T2D and conducting further research focusing on the molecular mechanism of IR are necessary.
Bioinformatics analysis of transcriptional profiles using RNA sequencing (RNA-seq) methods is an innovative technique for investigating disease pathophysiology, identifying disease biomarkers, and proposing therapeutic targets [12,13]. According to several transcriptomic profiling studies, the gene expression profiles in peripheral blood cells of patients with T2D vary considerably, indicating that the proteins encoded by differentially expressed mRNAs might play a role in the etiology of the disease [14,15]. The integration of bioinformatics methodologies and expression profiling techniques may help in determining novel biomarkers and discovering the molecular pathogenesis of T2D.
Recently, the N 6 -methyladenosine (m 6 A) modification of mRNA has been proposed to play a vital role in the progression of T2D [16][17][18]. The aims of this study were to identify metabolism-related proteins that could function as potential biomarkers of IR and to elucidate the role of m 6 A modification in the pathogenesis of this condition. The results of this study may contribute reliable biomarkers for the early detection and prevention of T2D.

Design of the Study
First, the RNA-seq dataset GSE174475 from the Gene Expression Omnibus (GEO) database was used to identify differentially expressed genes (DEGs) between patients with IR and matched healthy individuals (controls). Second, DEGs coding for metabolismrelated proteins (MP-DEGs) were selected. Functional and pathway enrichment analyses of the MP-DEGs were then carried out using the Gene Ontology (GO) and Kyoto Encyclo- pedia of Genes and Genomes (KEGG) databases. Third, a protein-protein interaction (PPI) network of the MP-DEGs was constructed and used to select functional modules and hub genes. Subsequently, 10 hub genes were selected for further analysis. Key hub genes were further identified using the least absolute shrinkage and selection operator (LASSO) regression analysis method, and their clinical value for the diagnosis of IR was evaluated using receiver operating characteristic (ROC) curves. Furthermore, the relationships between candidate gene expression and m 6 A modification were explored. Finally, adipose tissues from patients with IR and healthy individuals were collected to validate their expression of the key MP-DEGs and their correlation with m 6 A modification. A flowchart of this study is shown in Figure 1.

Data Acquisition and Processing
The RNA-seq dataset GSE174475, which contains data on pre-processed adipose tissue samples from 27 patients with IR and 16 healthy individuals (controls), was downloaded from the GEO database [19]. The BioBase package was used to normalize the data.

Selection of MP-DEGs
The limma package in R is a powerful tool for assessing DEGs [20]. We used this tool to examine the DEGs between the patient and control tissues in the GSE174475 dataset (detection criteria: p < 0.05 and |log2 fold change| ≥ 1).
The list of genes coding for metabolism-related proteins was downloaded from the Human Protein Atlas (HPA) protein annotation database [21]. The DEGs screened using limma were then intersected with the genes on the HPA list to screen out the MP-DEGs, and these were then examined for differences in expression between the two tissue groups.

Functional and Pathway Enrichment Analyses of the MP-DEGs
The GO and KEGG enrichment analyses of the MP-DEGs were conducted using the clusterProfiler package in Bioconductor. We first identified the GO biological processes (BPs), cellular components (CCs), and molecular functions (MFs) of the MP-DEGs using the human genome as a background reference (cutoff values: p < 0.05 and count = 2). After classifying the MP-DEGs into upregulated and downregulated groups, KEGG pathway enrichment analysis was conducted for each group, using a p-value of less than 0.05 and a count of 2 or higher as cutoff values. Finally, the findings of the enrichment study were displayed as dot plots and bar plots using the clusterProfiler package.

Construction of the Protein-Protein Interaction Network of the MP-DEGs and Gene Expression Analysis
PPIs in functional protein association networks were evaluated with the STRING database, using core factors as question proteins. A PPI network of the MP-DEGs was built on the basis of STRING11.5 [22]. Cytoscape was used to visualize the network, with a confidence score of greater than 0.4, and the disconnected genes were concealed. Text mining data were used to filter out interactions that have been confirmed in published literature. To identify functional clusters of genes in the PPI network, Molecular Complex Detection (MCODE) in Cytoscape was used, with a degree cutoff of 2, node score cutoff of 0.2, k-core of 2, and maximum depth of 100 [23]. Modules with proven scores greater than 5 were eliminated. The top 10 node genes were identified using 10 different analytical methods in the CytoHubba plug-in in Cytoscape, and the hub genes were then screened using the intersection of the results. In CytoHubba, the maximal clique centrality (MCC) technique is the most accurate one [24].

LASSO Regression Analysis
LASSO analysis was performed using the glmnet package in R, and the most advantageous cost of the penalty parameter was determined via 10-fold cross-validation [25]. First, the glmnet package was used to create a model with "binomial" as the family setting. Next, the important variables were selected according to the coefficients in the abovementioned model. Then, "deviance" was used as a minimum target parameter in the cross-verification fitting model. Finally, the target variables were selected on the basis of the variable coefficient and the optimal variables of the optimal lambda.

Clinical Sample Collection and Processing
Human cervical adipose tissues from 42 patients with IR and 52 controls were obtained from Xuanwu Hospital, Capital Medical University (Beijing, China) from December 2021 to December 2022. Peripheral venous blood samples of 10 mL were collected. Peripheral venous blood samples of 5 mL were used to measure fasting plasma glucose (FPG), triglycerides (TG), total cholesterol (TC), high-density lipoprotein cholesterol (HDL-C), and low-density lipoprotein cholesterol (LDL-C) by an automated chemistry analyzer (7600, Tokyo, Hitachi, Japan), and plasma FIns levels by aγ-automatic radioimmunoassay counter (XH-6020, Hitachi, Japan). All tissue samples and clinical data were collected with each patient's consent. The study protocol was approved by the Ethics Committee of Capital Medical University (Approval No. 2017SY24).
All adipose tissue samples were stored in a refrigerator at −80 • C. Total RNA was extracted from each sample using the SV Total RNA Isolation System (Z3100, Promega, Madison, WI, USA) according to the manufacturer's instructions. The reverse transcriptionquantitative polymerase chain reaction (RT-qPCR) was performed on a QuantStudio5 system (Applied Biosystems, Thermo Fisher, Waltham, MA, USA) using the RT2 SYBR Green ROX FAST Mastermix (Qiagen, Germantown, MD, USA). The gene-specific primers used are shown in Supplementary Table S1. GAPDH was used as an endogenous reference to normalize the relative amount of each target gene. The expression level was calculated using the 2 −∆∆CT method.

Grouping of Study Subjects
The study subjects of this experiment were divided into an insulin-resistant-positive group and a normal control group. The insulin resistance index was evaluated by the homeostasis model assessment (HOMA) and the HOMA-IR (homeostasis model assessment insulin resistance) index. HOMA-IR value >2.5 is IR positive [26]. HOMA-IR = FIns×PG 22.5 . The IR group was IR positive, and individuals were excluded from the study based on exclusion criteria. The exclusion criteria were as follows: patients with T2D; patients with coronary heart disease and myocardial disease; patients with malignant tumors; patients with acute and chronic infectious diseases; patients with severe liver and kidney diseases; those who use drugs to lower glucose; and those on hormones, psychotropics, and other prescription medications. The normal control group was IR negative, and individuals were excluded from the study based on exclusion criteria. HOMA-IR value <2.5 is IR negative. The exclusion criteria were as follows: patients with T2D; patients with coronary heart disease and myocardial disease; patients with malignant tumors; patients with acute and chronic infectious diseases; patients with severe liver and kidney diseases; and those who use drugs to lower glucose; those on hormones, psychotropics, and other prescription medications.

Gene-Specific MeRIP-qPCR
The EpiQuik TM CUT&RUN m 6 A RNA Enrichment Kit (A-P-9018, EpiGentek, Farmingdale, NY, USA) was used to perform the m 6 A-based methylated RNA immunoprecipitation (MeRIP) assay according to the manufacturer's instructions. Next, using an mRNA purification kit (Z5300; Promega, USA), a portion of the total RNA was converted into mRNA. The enriched RNA was immediately quantified using qPCR after enrichment of the methylated RNA with an m 6 A antibody. In brief, A/G immunomagnetic beads were first mixed with m 6 A antibodies, following which the mixture was added to the mRNAs. Next, the m 6 A immunomagnetic beads were enriched using a magnetic frame, and the enriched RNA-antibody complexes were digested with protease. Finally, routine qPCR was performed for quantification of the mRNAs. The primer sequences used for the qPCR are listed in Supplementary Table S1.

Statistical Analysis
DEGs were examined using the limma package in R (version 4.0). The PPI network and hub genes were analyzed using Cytoscape 8.0. GO and KEGG enrichment analyses of the MP-DEGs were conducted using the clusterProfiler package in Bioconductor. LASSO regression analysis was performed using the glmnet package in R. Correlations between FASN, GCK, FBP1, and FGR and m 6 A-related genes were analyzed using Pearson's correlation coefficient for normal data and Spearman's test for non-normal data. The average mRNA and methylation levels are expressed as the mean ± standard deviation. Differences between groups were compared using the t-test. Statistical significance was set at a p-value of less than 0.05.

Identification of Differentially Expressed Genes
DEG analysis was performed to assess differences between the IR and control groups of the GSE174475 dataset. Basic information of samples for GSE174475 used in the study is presented in Supplementary Table S2. The genes in this dataset had the same quartile, Max, and Min values (Figure 2A). In total, 344 genes were identified as DEGs (Supplementary File S1). CDC20, ABCC3, MATK, SPOCD1, and MFSD12 were the top five upregulated genes, whereas NAALAD2, PRTG, Clorf185, NCAM2, and COL6A6 were the top five downregulated ones, with the lowest p-values ( Figure 2B). The heatmap showed that the DEGs were consistently significantly upregulated or downregulated in the IR group relative to their expression levels in the control group ( Figure 2C).

Selection of the MP-DEGs
Annotated metabolism-related genes from the HPA database were used to identify the genes encoding metabolism-related proteins among the DEGs between the IR and control groups. After intersecting the HPA-sourced list with the DEGs identified from the GSE174475 dataset, 69 MP-DEGs were filtered out, of which 15 were upregulated and 54 were downregulated ( Figure 3A and Supplementary File S2). ABCC3, MATK, CHI3L1, SDS, ACP5, IL4I1, PLA2G7, SLC29A3, ATP6V0D2, and SLC9A7 were the top 10 upregulated genes, whereas SLC6A2, CA3, AQP6, APOB, SLC2A4, GCK, SLC27A2, MOGAT1, SLC16A9, and FASN were the top 10 downregulated ones ( Figure 3B and Supplementary File S2). Among the top 15 significantly downregulated genes in IR, SLC6A2, CA3, AQP6, APOB, and SLC2A4 had the smallest p-values according to the heatmap of the MP-DEGs ( Figure 3C).

Functional and Pathway Enrichment of the MP-DEGs
The screened MP-DEGs were mainly enriched in the GO BP categories regulation of hormone levels and hormone metabolic process, the CC category low-density lipoprotein particles, and the MF categories carboxylic acid transmembrane transporter activity and organic acid transmembrane transporter activity (Figure 4). The MP-DEGs were shown using the cnetplot function in the clusterProfiler package. SLC27A2, FASN, ELOVL6, SLC16A9, SLC6A12, and SDS were enriched in multiple terms of the BP category ( Figure 5 and Supplementary File S3). KEGG pathway analysis was used to determine the pathways enriched in the upregulated and downregulated genes. The upregulated genes were not enriched in any pathway. By contrast, the downregulated genes were enriched in the insulin signaling, IR, and AMPK signaling pathways ( Figure 6).

Protein-Protein Interaction Network Construction and Module Analysis
The STRING database was used to create a PPI network from the 69 MP−DEGs, following which the network was used to explore the interactions between the metabolismrelated proteins. Cytoscape was used to visualize the entire PPI network ( Figure 7A). The PPI network comprised 72 edges and 69 nodes. The MCODE plug-in in Cytoscape was used to create the functional modules. The findings revealed that only one module, comprising 9 genes and 34 edges, had an established score of greater than 5 ( Figure 7B). The top 10 hub genes were screened using the 10 topological approaches of the CytoHubba plug-in in Cytoscape. Ten genes were identified using all ten methods, namely, FASN, GCK, FGR, FBP1, GYS2, PNPLA3, MOGAT1, SLC27A2, PNPLA3, and ELOVL6 ( Figure 7C and Supplementary File S4). FASN was selected as the key gene because it had the highest MCC score. The first node gene that interacted with FASN was determined using CytoHubba. Of the 10 genes screened, 1 was upregulated and 9 were downregulated ( Figure 7D).

Selection and Validation of the Primary Hub Genes
The 10 hub genes identified from the PPI network were submitted to LASSO regression analysis to identify the primary hub genes, whereupon GCK, FBP1, and FGR were selected in the model (Figure 8). These three genes, as well as FASN (another key gene in the PPI network), were studied further using ROC curve analysis. GCK, FBP1, and FGR were highly accurate in diagnosing IR (AUC = 0.80, 95% CI = 0.67-0.94; AUC = 0.86, 95% CI = 0.74-0.94; AUC = 0.83, 95% CI = 0.64-0.92), whereas FASN showed moderate accuracy in this regard (AUC = 0.78, 95% CI = 0.64-0.92; Figure 8). These results indicate that FASN, GCK, FBP1, and FGR can be applied as reliable biomarkers of IR, with the expression of FASN, GCK, and FBP1 being downregulated and that of FGR being upregulated in patients with this condition.

FASN Is Lower in IR and Associated with FPG and HbA1c
To verify the FASN and GCK diagnostic values for IR, the mRNA levels of FASN and GCK for adipose tissue were detected in clinical samples from 42 patients with IR and 52 controls. Clinical information of the samples is presented in Supplementary Table S3. The FASN and GCK mRNA expression was lower in patients with IR than in the controls ( Figure 9A, p < 0.05). Moreover, FASN mRNA expression was negatively correlated with FPG, HbA1c, Homa-IR, and FIns (r = −0.309, p = 0.002; r = −0.296, p = 0.004; r = −0.279, p = 0.006, r = −0.325, p = 0.001) in all participants (Table 1). GCK mRNA expression was negatively correlated with FPG and Homa-IR (r = −0.220, p = 0.017; r = −0.230, p = 0.022) in all participants ( Table 1). The strength of the association between the mRNA levels of FASN, GCK, and IR was evaluated by univariate and multivariate regression models. The univariate logistic model showed that FASN and GCK were associated with IR. Subsequently, several multivariate logistic regression models were established by adjusting for covariates. The multivariate logistic models showed that FASN and GCK were also associated with IR, indicating that the association of FASN and GCK with IR was independent (Supplementary Table S4). We further evaluated the ability of the FASN and GCK to predict IR as a biomarker by ROC analysis. The results showed that the FASN was moderately effective for the detection of IR, with the area under the curve being 0.78 and the cutoff value being 0.33% (Table 2 and Supplementary Figure S1). These results indicated that FASN could be used as a reliable biomarker of IR.

Correlations of Candidate Gene Expression with m 6 A Modification in Insulin Resistance
m 6 A modification plays a vital role in the progression of T2D. To explore its involvement in IR, the GSE174475 dataset was used to investigate the relationship between the expression of FASN, GCK, FBP1, and FGR and that of 15 m 6 A-related genes in IR [27]. The findings revealed that the expression of these four key genes was significantly correlated with that of IGF2BP3, FTO, EIF3A, WTAP, METTL16, and LRPPRC (p < 0.05; Figure 10). Specifically, FASN expression was significantly positively correlated with FTO, EIF3A, and LRPPRC expression (p < 0.05) and negatively correlated with IGF2BP3 and WTAP expression (p < 0.05). The expression of GCK was significantly positively correlated with that of EIF3A and LRPPRC (p < 0.05) and negatively correlated with that of WTAP (p < 0.05). FGR expression was significantly positively correlated with IGF2BP3 and WTAP expression (p < 0.05) and negatively correlated with FTO, EIF3A, METTL16, and LRPPRC expression (p < 0.05). The expression of FBP1 was significantly positively correlated with that of IGF2BP3 and WTAP (p < 0.05) and negatively correlated with that of FTO, EIF3A, METTL16, and LRPPRC (p < 0.05). These results suggest that the changes in expression of the candidate genes may be associated with m 6 A modification. Based on the results on the ability of the FASN to predict IR as a biomarker and the correlation with m 6 A modification, meRIP-qPCR was conducted to assess FASN methylation levels in clinical samples. The results showed that the methylation levels of FASN in adipose tissue of patients with IR were decreased ( Figure 9B, p < 0.05). Additionally, we found FASN mRNA expression was positively correlated with the methylation levels of FASN (r = 0.359, p = 0.001) in all participants ( Figure 9C). These data suggested that the downregulation of FASN may be caused by the decrease in the methylation of FASN in IR.

Discussion
To the best of our knowledge, this is the first study to identify metabolism-related hub genes associated with IR and the role of m 6 A modification in the pathogenesis of this condition. Additionally, we attempted to validate these essential genes in clinical adipose tissue samples. The findings offer a novel perspective for understanding the molecular basis of IR in adipocytes.
The adipose tissue transcriptome was selected because adipocytes are essential regulators of metabolic homeostasis and play an indispensable regulatory role in IR. Our study had a relatively large transcriptome sample size, being based on the GSE174475 RNA-seq dataset, which includes tissue samples from 27 patients with IR and 16 healthy individuals. Compared with other datasets, GSE174475 comprises transcriptomes in adipocytes rather than in whole adipose tissue, thus eliminating the effect of many non-adipocyte cells in adipose tissue.
Of the 344 DEGs identified in the GSE174475 dataset, 69 were identified in the GeneCards database to encode metabolism-related proteins. GO functional enrichment analysis revealed that these MP-DEGs were enriched in the BP categories regulation of hormone levels and hormone metabolic process, the CC category low-density lipoprotein particle, and the MF categories carboxylic acid transmembrane transporter activity and organic acid transmembrane transporter activity. SLC27A2, FASN, ELOVL6, SLC16A9, SLC6A12, and SDS were all enriched in various biological processes according to the GO enrichment circle map. These genes might be more significant than others in IR. KEGG enrichment analysis showed that the downregulated MP-DEGs were primarily enriched in the insulin signaling, IR, and AMPK signaling pathways. Studies have shown that AMPK plays a major role in regulating cellular energy balance [28]. Therefore, the results indicate that the development of IR and T2D may be regulated by these metabolism-related proteins through their interference with the insulin signaling pathway or energy balance.
Using the PPI network of MP-DEGs and the top 10 hub genes assessed by 10 topological approaches in CytoHubba, we discovered that FASN, GCK, FGR, FBP1, GYS2, PNPLA3, MOGAT1, SLC27A2, PNPLA3, and ELOVL6 were all simultaneously included in the functional gene modules created by MCODE. These 10 hub genes were highly correlated with IR. The logistic LASSO model can be used to select a greater and more accountable set of predictors from the regression's massive and underlying multicollinearity set of variables [25]. Through LASSO regression analysis, the 10 hub genes were reduced to three key genes, namely, GCK, FBP1, and FGR. FASN was found in the functional module, being a downregulated gene with the highest MCC score (the most accurate method). These four genes were included in the ROC curve analysis, whereupon all four of them showed reasonable diagnostic value for IR, indicating that they are potential biomarkers for the early detection and prevention of T2D. Consistent with our study results, FBP1 has been previously suggested as a biomarker for the early detection of T2D [29]. GCK encodes glucokinase, which catalyzes the phosphorylation of glucose to glucose-6-phosphate. This enzyme plays an indispensable role in the regulation of glucose metabolism and insulin secretion in pancreatic β-cells [30]. Additionally, it was noted in a review that polymorphisms in GCK are risk predictors of T2D [31]. A study on building specific biosignatures via machine learning showed that GCK methylation was the major feature in two study models and could serve as a biomarker of T2D [32]. Our present study showed that low GCK expression in IR was related to early diagnosis.
In our validation of the clinical samples, we found the lower expression of FASN has a certain diagnostic value for IR. Many studies have indicated that FASN can be used as a biomarker for the diagnosis and prognosis of various cancers, such as triplenegative breast cancer, gastric cancer, and glandular cancer [33][34][35][36][37][38]. Decreased FASN expression in the subcutaneous adipose tissue of obese individuals and mice with IR or T2D has been reported [39][40][41]. Similarly, a study suggests that FASN can be used as a biomarker for IR, and extracellular FASN has the potential to become an alternative biomarker for the condition [42]. Moreover, studies have shown a correlation of FASN mRNA expression with its protein levels in visceral and subcutaneous adipose tissues as well as with impaired insulin sensitivity [40]. The expression and activity of FASN are influenced by insulin [42][43][44], and our results also showed a negative correlation between FASN and insulin. Therefore, we speculate that insulin decreases FASN mRNA expression by regulating transcription. FASN encodes fatty acid synthase, which plays an important role in lipid metabolism as the primary enzyme involved in the de novo synthesis of fatty acids [43,45]. A study by Sun et al. revealed that FTO regulates hepatic lipogenesis via FTO-dependent m 6 A demethylation in FASN mRNA [46]. Additionally, it was recently demonstrated that METTL3 inhibits hepatic insulin sensitivity via the N6-methylation of FASN mRNA and promotes fatty acid metabolism [18]. N6-methyladenosine (m 6 A) is a methylation modification formed by the catalysis of a methyltransferase at the sixth nitrogen atom, which regulates gene expression by modulating RNA degradation, splicing, export, and translation [47]. The m 6 A modification is dynamically regulated by an enzyme system, which mainly includes three types of proteins: writers that catalyze m 6 A methylation, erasers that catalyze demethylation, and readers that recognize methylation [48]. Changes in the m 6 A level on a specific gene do not always correlate with changes in the total m 6 A level. This is because the m 6 A modification on a specific gene is dynamically regulated by different methyltransferases ("writers"), demethylases ("erasers"), and binding proteins ("readers"), and each gene has specificity [49]. Our study revealed that FASN expression was downregulated in IR and that FASN mRNA expression was positively correlated with the methylation levels of FASN. The downregulation of FASN may be caused by the decrease in the methylation of FASN in IR. These results suggest that the m 6 A modification of FASN may contribute to the development of IR and T2D. The genes of FASN that serve as potential biomarkers of IR should be studied in a larger population, and the role of m 6 A modification in the pathogenesis should be validated in vivo and in vitro.
The current study had several limitations. First, the clinical samples were designed for a case-control study. Second, the validation sample was from a single-center case-control study. Cohort and nested case-control studies with more samples are required to verify the predictive value of these biomarkers. And body composition data such as body fat and muscle mass should also be analyzed in future studies. Moreover, in vivo and in vitro molecular mechanism studies should be carried out to validate how FASN is regulated by m 6 A modification in IR.

Conclusions
In this study, 69 MP-DEGs were screened from the gene expression profiles of human adipose tissues with IR. Four key MP-DEGs were found to distinguish the IR condition. Moreover, the expression of FASN and GCK in clinical samples was validated. Additionally, we found the downregulation of FASN may be caused by the decrease in the methylation of FASN in IR. These results will help toward the identification of reliable biomarkers for the early detection and prevention of T2D.
Supplementary Materials: The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/nu15081839/s1, File S1: Design of the study of 343 differentially expressed genes (DEGs) between the insulin resistance (IR) and control groups in the GSE174475 dataset; File S2: 69 differentially expressed genes of metabolism-related proteins (MP-DEGs); File S3: Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses of the 69 differentially expressed genes of metabolism-related proteins (MP-DEGs); File S4: Top 10 differentially expressed genes of metabolism-related proteins (MP-DEGs) identified using 10 topological analysis methods of CytoHubba; Table S1: The Primer sequences used in the study; Table S2: Basic information of samples for GSE174475 used in the study; Table S3: Clinical information of samples used in the study; Table S4: Logistic regression models analysis of insulin resistance; Figure S1:. ROC curves for FASN and GCK.  Informed Consent Statement: Informed consent was obtained from all subjects involved in the study.

Data Availability Statement:
The data presented in this study are not publicly available due to confidentiality reasons. The datasets are available from the corresponding author upon reasonable request.