Tsc22d3 promotes morphine tolerance in mice through the GPX4 ferroptosis pathway

Differential gene analysis

In this study, using the predefined cutoff values, we identified differentially expressed genes (DEGs) from the gene expression matrix of the GSE7762 dataset, resulting in a total of 500 DEGs (Figure 1A).

Figure 1. (A) Differential gene analysis. A total of 500 DEGs. (B–E) Results of GOKEGG enrichment analysis of DEGs. (B) Biological process analysis. (C) Cellular component analysis. (D) Molecular function analysis. (E) Results of KEGG enrichment analysis.

Functional enrichment analysis

GO and KEGG analysis of DEGs

Functional enrichment analysis, including GO and KEGG analysis, was performed on these differentially expressed genes. According to the GO analysis results, in Biological Process (BP) analysis, they were mainly enriched in negative regulation of brain development, neuronal apoptosis processes, and regulation of nervous system development (Figure 1B). In Cellular Component (CC) analysis, they were mainly enriched in neuronal projection and neuronal cell body (Figure 1C). In Molecular Function (MF) analysis, they were concentrated in the regulation of postsynaptic cytoplasmic calcium ion concentration involving neurotransmitter receptor activity and regulation of postsynaptic membrane potential by G-protein coupled receptor activity (Figure 1D). In KEGG analysis, they were mainly enriched in opioid addiction, dopaminergic synapse, neuroactive ligand-receptor interaction, and MAPK signaling pathway (Figure 1E).

Metascape enrichment analysis

In the Metascape enrichment projects, regulatory aspects of brain development, nervous system development, and the MAPK signaling pathway in mice were visible in the GO enrichment projects (Figure 2A). Additionally, we generated enrichment networks colored by enrichment terms and p-values (Figure 2B–2D), providing a visual representation of the associations and confidence levels of various enrichment projects.

Figure 2. Metascape enrichment analysis. (A) Bar graph of enriched terms across input gene lists, colored by p-values. (B) Network of enriched terms: colored by cluster ID, where nodes that share the same cluster ID are typically close to each other. (C) colored by p-value, where terms containing more genes tend to have a more significant p-value. (D) Protein-protein interaction network. MCODE components identified in the gene lists.

WGCNA

The selection of the soft-thresholding power is a crucial step in WGCNA analysis. The soft-thresholding power for network topology analysis was set to 7 (Figure 3A). A hierarchical clustering tree was constructed for all genes, resulting in 25 modules (Figure 3B). Interaction between important modules was analyzed (Figure 3C), and a heatmap of module-phenotype correlation (Figure 3D) and a scatter plot of the correlation between GS and MM of relevant hub genes (Figure 3E) were generated.

Figure 3. WGCNA. (A) β = 6,0.86. β = 6, 22.89. (B, C) The hierarchical clustering tree of all genes was constructed, and 25 important modules were generated. (D) The heat map of correlation between modules and phenotypes. (E) The scatter map of correlation between GS and MM of related hub genes. (F) The DEGs screened by WGCNA and DEGs were used to obtain Venn map. 32 intersection genes were obtained.

We also created a Venn diagram by intersecting genes selected by WGCNA with DEGs for further analysis (Figure 3F).

Protein-protein interaction (PPI) network construction and analysis

The PPI network of DEGs was constructed using the STRING online database and analyzed using Cytoscape software (Figure 4A). Central genes were identified using five different algorithms (MCC, MNC, DMNC, EPC, EcCentricity) (Figure 4B–4F). The Venn diagram was used to obtain the intersection (Figure 4G), resulting in four core genes (Tsc22d3, Itga6, Ermn, Sult1a1).

Figure 4. Construction and analysis of protein-protein interaction (PPI) networks. (A) Construct the PPI network of DEGs using STRING online database and utilize Cytoscape software for analysis. (B) MCC was used to identify the central gene. (C) MNC was used to identify the central gene. (D) DMNC was used to identify the central gene. (E) EPC was used to identify the central gene. (F) EcCentricity was used to identify the central gene. (G) Core genes (Tsc22d3, Itga6, Ermn, Sult1a1) were obtained by merging using Venn diagrams.

Gene expression heatmap

We visualized the expression levels of core genes in the gene expression matrix of the GSE7762 dataset and created a heatmap (Figure 5A). We observed that core genes (Tsc22d3, Sult1a1) were highly expressed in brain tissue samples injected with morphine, while they were lowly expressed in normal samples, indicating significant differences. Based on these results, we speculate that these core genes may play a regulatory role in morphine tolerance after injection.

Figure 5. Heat map of gene expression. (A) Heat map of the core gene in the gene expression matrix of the data sets GSE7762. (B) CTD analysis. Core genes (Tsc22d3, Itga6, Ermn, Sult1a1) are associated with immune system diseases, drug-related side effects and adverse reactions, chemical and drug-induced liver injury, neurological diseases, heart failure, mental disorders, and neural tube defects.

CTD analysis

In this study, we input the hub gene list into the CTD website to identify diseases related to the core genes, enhancing our understanding of gene-disease associations. We found that the four core genes (Tsc22d3, Itga6, Ermn, Sult1a1) were associated with immune system diseases, drug-related side effects and adverse reactions, chemical and drug-induced liver injury, neurological disorders, heart failure, mental disorders, and neural tube defects (Figure 5B).

Prediction and functional annotation of miRNAs associated with hub genes

In this study, we input the hub gene list into TargetScan to identify relevant miRNAs, enhancing our understanding of gene expression regulation (Table 1). We found that the relevant miRNAs for the Tsc22d3 gene were mmu-miR-196b-5p and mmu-miR-196a-5p; for the Itga6 gene, it was mmu-miR-126a-3p.1; and for the Ermn gene, it was mmu-miR-582-5p, mmu-miR-219a-2-3p, and mmu-miR-496a-3p.2.

Impact of Tsc22d3 on iron death-related pathway protein expression in morphine tolerance

The expression of Tsc22d3 was not significantly different between the CON group and the Non-morphine tolerant group. Compared to the Non-morphine tolerant group, the expression of Tsc22d3 was significantly upregulated in the Morphine tolerant group. Compared to the Morphine tolerant group, the expression of Tsc22d3 was further upregulated in the Morphine tolerant+Tsc22d3_OE group. Compared to the Morphine tolerant group, the expression of Tsc22d3 was significantly downregulated in the Morphine tolerant+Tsc22d3_KO group. The key node proteins (HIF-1alpha, GSH, CCL2, GPX4) in the GPX4 iron death-related pathway did not show significant differences in expression between the CON group and the Non-morphine tolerant group. Compared to the Non-morphine tolerant group, the expression of key node proteins (HIF-1alpha, GSH, GPX4) in the GPX4 iron death-related pathway was significantly downregulated in the Morphine tolerant group. Compared to the Non-morphine tolerant group, the expression of the key node protein (CCL2) in the GPX4 iron death-related pathway was significantly upregulated in the Morphine tolerant group. With the upregulation of Tsc22d3 expression in Morphine tolerant+Tsc22d3_OE, the expression of key node proteins (HIF-1alpha, GSH, GPX4) in the GPX4 iron death-related pathway further decreased significantly. With the upregulation of Tsc22d3 expression in Morphine tolerant+Tsc22d3_OE, the expression of the key node protein (CCL2) in the GPX4 iron death-related pathway further increased significantly. With the downregulation of Tsc22d3 expression in Morphine tolerant+Tsc22d3_KO, the expression of key node proteins (HIF-1alpha, GSH, GPX4) in the GPX4 iron death-related pathway further increased significantly. With the downregulation of Tsc22d3 expression in Morphine tolerant+Tsc22d3_KO, the expression of the key node protein (CCL2) in the GPX4 iron death-related pathway further decreased significantly (Figure 6).

Figure 6. Effect of Tsc22d3 on iron death-related pathway protein expression in morphine tolerance.

Tsc22d3 promotes apoptosis of brain cells in morphine tolerance

Compared to the CON group, the expression of apoptosis proteins (P53, Caspase-3, BCL-2, Bax, SMAC, FAS) in the Non-morphine tolerant group was not significantly different. Compared to the Non-morphine tolerant group, the expression of apoptosis proteins (P53, Caspase-3, Bax, SMAC, FAS) was significantly upregulated in the Morphine tolerant group. Compared to the Non-morphine tolerant group, the expression of the anti-apoptotic protein (BCL-2) was significantly downregulated in the Morphine tolerant group. With the upregulation of Tsc22d3 expression in Morphine tolerant+Tsc22d3_OE, the expression of apoptosis proteins (P53, Caspase-3, Bax, SMAC, FAS) further increased significantly. With the upregulation of Tsc22d3 expression in Morphine tolerant+Tsc22d3_OE, the expression of the anti-apoptotic protein (BCL-2) further decreased in the Morphine tolerant group. With the downregulation of Tsc22d3 expression in Morphine tolerant+Tsc22d3_KO, the expression of apoptosis proteins (P53, Caspase-3, Bax, SMAC, FAS) further decreased significantly. With the downregulation of Tsc22d3 expression in Morphine tolerant+Tsc22d3_KO, the expression of the anti-apoptotic protein (BCL-2) further increased in the Morphine tolerant group (Figure 7).

Figure 7. Tsc22d3 promotes apoptosis of brain cells in morphine tolerance.

Tsc22d3 promotes inflammatory response of brain cells in morphine tolerance

Compared to the CON group, the expression of inflammatory factors (IL6, TNF-alpha, CXCL1, CXCL2) in the Non-morphine tolerant group did not show significant differences. However, when compared to the Non-morphine tolerant group, the expression of inflammatory factors (IL6, TNF-alpha, CXCL1, CXCL2) was significantly upregulated in the Morphine tolerant group. With the upregulation of Tsc22d3 expression in Morphine tolerant+Tsc22d3_OE, the expression of inflammatory factors (IL6, TNF-alpha, CXCL1, CXCL2) further increased significantly. Conversely, with the downregulation of Tsc22d3 expression in Morphine tolerant+Tsc22d3_KO, the expression of inflammatory factors (IL6, TNF-alpha, CXCL1, CXCL2) further decreased significantly (Figure 8).

Figure 8. Tsc22d3 promotes inflammation in brain cells during morphine tolerance.

Research outlook

This study holds promising potential for clinical applications. Further research confirming Tsc22d3 as a key regulatory factor in opioid tolerance could lead to the development of related drugs or treatment strategies to help patients cope with opioid tolerance more effectively.

Morphine tolerance dataset

The morphine tolerance dataset GSE7762 profile was downloaded from the gene expression omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo/), generated from GPL1261 (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE7762), to identify differentially expressed genes (DEGs) in the brain tissue associated with morphine tolerance.

DEGs selection

The R package “limma” was used for probe summarization and background correction of the gene expression matrix from GSE7762. The Benjamini-Hochberg method was employed to adjust the raw p-values. The fold change (FC) was calculated using the false discovery rate (FDR). The cutoff criterion for DEGs was set at P < 0.05. A volcano plot was generated to visualize the DEGs.

Weighted gene co-expression network analysis (WGCNA)

Initially, the Median Absolute Deviation (MAD) for each gene was calculated using the gene expression matrix of GSE7762. The bottom 50% of genes with the smallest MAD values were excluded. The “goodSamplesGenes” method from the R package WGCNA was utilized to remove outlier genes and samples. To categorize genes with similar expression profiles into gene modules, hierarchical average linkage clustering based on TOM-based dissimilarity was performed, with a minimum module size (number of genes) set at 30. The sensitivity was set at 3. For further module analysis, dissimilarity of module characteristic genes was calculated, and a cutting line was chosen for module dendrogram, merging some modules. Additionally, modules with a distance less than 0.25 were merged. It is noteworthy that the grey module is considered a collection of genes that cannot be assigned to any module.

Construction and analysis of protein-protein interaction (PPI) network

The STRING database (http://string-db.org/) aims to collect, score, and integrate all publicly available protein-protein interaction information sources and supplement them through computational predictions. In this study, the list of differentially expressed genes (DEGs) was input into the STRING database to construct a predicted core gene PPI network (confidence > 0.4). Cytoscape software provides biologists with tools for biological network analysis and two-dimensional (2D) visualization. In this study, the PPI network formed by the STRING database was visualized and the core genes were predicted using Cytoscape software. The PPI network was imported into Cytoscape, and the most relevant modules were identified using MCODE. Additionally, five algorithms (MCC, MNC, DMNC, EPC, EcCentricity) were employed to calculate the top ten genes with the best relevance, and the intersection was taken. The resulting core gene list was visualized and exported after visualization.

Functional enrichment analysis

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses are computational methods for assessing gene function and biological pathways. In this study, the list of differentially expressed genes selected by the Venn diagram was input into the KEGG REST API (https://www.kegg.jp/kegg/rest/keggapi.html) to obtain the latest gene annotations for KEGG pathways. This was used as a background for mapping genes to the background set. Enrichment analysis was performed using the R package clusterProfiler (version 3.14.3) to obtain results of gene set enrichment. Additionally, GO annotations for genes from the org.Hs.eg.db package (version 3.1.0) were used as a background, with a minimum gene set size of 5, a maximum gene set size of 5000, and criteria of P-value < 0.05 and FDR < 0.25 considered statistically significant.

Furthermore, the Metascape database provides comprehensive gene list annotation and analysis resources with visualization capabilities. Metascape (http://metascape.org/gp/index.html) was used for functional enrichment analysis of the aforementioned list of differentially expressed genes, and the results were exported.

Gene expression heatmap

The R package heatmap was utilized to create a heatmap of the expression levels of core genes identified in the PPI network in the gene expression matrix of GSE7762. This visualization depicts the expression differences between brain tissues injected with morphine and normal brain tissues.

CTD Analysis

The Comparative Toxicogenomics Database (CTD) integrates a vast amount of data on interactions between chemicals, genes, phenotypes, and diseases, providing significant convenience for studying disease-related environmental exposure factors and potential drug mechanisms. The core genes were input into the CTD website, and the most relevant diseases associated with the core genes were identified. Radar plots depicting the expression differences for each gene were created using Excel.

miRNA

TargetScan (www.targetscan.org) is an online database used for predicting and analyzing miRNA and their target genes. In our study, TargetScan was employed to screen miRNAs regulating central differentially expressed genes (DEGs).

Animal experiments

Fifty C57BL/6J mice (6-8 weeks old) were purchased from Sibeifu (Beijing) Biotechnology Co., Ltd., acclimated for one week, and randomly divided into five groups, each consisting of 10 mice.

Control group (CON): Injected with physiological saline;

Non-morphine tolerant group (Non-morphine tolerant group): Single-dose administration to mice (20 mg/kg, s.c.);

Morphine tolerant group (Morphine tolerant group): Multiple-dose administration to mice (10 ~ 40 mg/kg, three times a day, for 5 days);

Morphine tolerant + Tsc22d3_KO group (Morphine tolerant + Tsc22d3_KO group): Multiple-dose administration to mice (10 ~ 40 mg/kg, three times a day, for 5 days) with Tsc22d3 knocked down using siRNA technology. Plasmids were constructed and injected into mice via tail vein injection.

Morphine tolerant + Tsc22d3_OE group (Morphine tolerant + Tsc22d3_OE group): Multiple-dose administration to mice (10 ~ 40 mg/kg, three times a day, for 5 days) with Tsc22d3 overexpressed. Plasmids were constructed and injected into mice via tail vein injection.

Detection of iron death-related pathway proteins

Brain tissue was isolated, and protein extraction was performed. The tissue chunks were washed 2-3 times with pre-chilled PBS to remove blood stains. They were then cut into small pieces and placed in a homogenizer tube. A 4mm homogenization bead was added, along with lysis buffer (containing various protease inhibitors) at 10 times the tissue volume. Homogenization was performed, and the homogenization tube was taken out and placed on ice for 30 minutes. The supernatant collected after centrifugation at 12000rpm, 4° C, for 10 minutes was considered the total protein solution. The BCA protein assay kit (Invitrogen) was used to determine the protein concentration. Protein denaturation was performed by adding 4:1 protein solution to 5× reducing protein loading buffer, boiling for 15 minutes in a water bath, and storing at -20° C. Different concentration separation gels were prepared based on experimental requirements. The gel was poured after adding TEMED, and a comb was inserted into the gel. Gel electrophoresis was conducted at a constant voltage of 200V for approximately 30 minutes. Once complete, bromophenol blue reached about 1cm from the bottom, and the electrophoresis was terminated for membrane transfer.

For membrane transfer, activated PVDF (0.45um) membrane was used after activation with methanol for 2 minutes. The gel was carefully peeled off and placed on filter paper, and the PVDF membrane was applied to the gel. Three pieces of filter paper were placed on top of the membrane, and any remaining water was removed. The transfer conditions were set as a constant current of 300mA for 30 minutes with the transfer equipment placed in ice water to cool during the process. The transferred membrane was blocked with 5% skimmed milk for 1 hour and then incubated with primary antibodies. Visualization of immune reaction bands was achieved through chemiluminescence using horseradish peroxidase-conjugated IgG secondary antibodies. Band density was measured using a densitometry method.

Data availability

The datasets generated during and/or analyzed during the current study are available from the corresponding author on reasonable request.