Skip to main content

Targeting the inward rectifier potassium channel 5.1 in thyroid cancer: artificial intelligence-facilitated molecular docking for drug discovery

Abstract

Background

Recurrent and metastatic thyroid cancer is more invasive and can transform to dedifferentiated thyroid cancer, thus leading to a severe decline in the 10-year survival. The thyroid-stimulating hormone receptor (TSHR) plays an important role in differentiation process. We aim to find a therapeutic target in redifferentiation strategies for thyroid cancer.

Methods

Our study integrated the differentially expressed genes acquired from the Gene Expression Omnibus database by comparing TSHR expression levels in the Cancer Genome Atlas database. We conducted functional enrichment analysis and verified the expression of these genes by RT-PCR in 68 pairs of thyroid tumor and paratumor tissues. Artificial intelligence-enabled virtual screening was combined with the VirtualFlow platform for deep docking.

Results

We identified five genes (KCNJ16, SLC26A4, TG, TPO, and SYT1) as potential cancer treatment targets. TSHR and KCNJ16 were downregulated in the thyroid tumor tissues, compared with paired normal tissues. In addition, KCNJ16 was lower in the vascular/capsular invasion group. Enrichment analyses revealed that KCNJ16 may play a significant role in cell growth and differentiation. The inward rectifier potassium channel 5.1 (Kir5.1, encoded by KCNJ16) emerged as an interesting target in thyroid cancer. Artificial intelligence-facilitated molecular docking identified Z2087256678_2, Z2211139111_1, Z2211139111_2, and PV-000592319198_1 (-7.3 kcal/mol) as the most potent commercially available molecular targeting Kir5.1.

Conclusion

This study may provide greater insights into the differentiation features associated with TSHR expression in thyroid cancer, and Kir5.1 may be a potential therapeutic target in the redifferentiation strategies for recurrent and metastatic thyroid cancer.

Peer Review reports

Introduction

Thyroid cancer is projected to become the fourth leading type of cancer with a substantial rise in its global incidence [1]. According to the histological classification, follicular thyroid cell-derived tumors can be classified as papillary thyroid carcinoma (PTC), follicular thyroid carcinoma (FTC), poorly differentiated thyroid carcinoma, and anaplastic thyroid carcinoma (ATC) [2]. Approximately 30% of the well-differentiated thyroid cancers (WDTC) progress to dedifferentiated thyroid cancers (DeTCs) and develop recurrence and metastasis [3]. Hence, the dedifferentiation status of DeTC and ATC indicates a more invasive approach and worse prognosis. Radioactive iodine (RAI) is the dominant therapy for WDTC patients with recurrent and metastatic diseases. The loss of differentiation features represents the reduced expression of thyroid-specific genes, such as the thyroid-stimulating hormone (TSH) receptor (TSHR), sodium/iodide symporter (NIS), thyroglobulin, thyroperoxidase (TPO), and paired box transcription factor 8 [4]. Such patients display resistance to RAI therapy (RAIR) and cannot benefit from RAI treatment. Thus, restoring differentiation may be a promising strategy for overcoming cancer invasion and RAIR.

TSHR is a representative thyroid-specific gene that reveals the differentiation status of the thyroid cells. The expression of TSHR consistently induces the formation of colonies with epithelial morphology and thyroglobulin expression [5]. TSH-TSHR signaling induces thyrocyte proliferation and the expression of other thyroid-specific genes. TSHR is required for optimal NIS expression and its localization to the plasma membrane, thus ensuring a response to RAI therapy [6, 7]. In the absence of functional TSHR, the thyroid gland develops to a normal size, whereas the expression of thyroperoxidase and sodium/iodide symporter gets significantly reduced [8]. Evidence from numerous studies on the resected thyroid tissue suggest that TSHR is more persistently expressed than other differentiation markers (including NIS and thyroglobulin ), thereby indicating it is a prominent differentiation marker [9].

As a marker of thyroid differentiation, TSHR expression may affect the function of cancer cells that have undergone malignant transformation, such as epithelial-mesenchymal transition (EMT) and dedifferentiation. Promoting TSHR expression may facilitate the transformation of DeTC to WDTC. Taken together, we analyzed the Cancer Genome Atlas (TCGA) database and Gene Expression Omnibus (GEO) database, grouped by TSHR expression level, screened WDTC and ATC, finally identified two differentially expressed genes (DEGs), namely KCNJ16 and SYT1. Subsequently, we verified the expression change of TSHR, KCNJ16 and SYT1 in 68 pairs of normal and tumor tissues. During exploration, we found the significant role of KCNJ16 in thyroid differentiation and cancer associated signaling pathways. Therefore, we targeted the protein encoded by KCNJ16 and used artificial intelligence facilitated molecular docking to search potential molecular for drug discovery.

Materials and methods

Patients sample collection

Seventy pairs of primary tumors and adjacent normal tissues were collected from patients diagnosed with thyroid tumor between January 2021 and September 2021 at the affiliated Tongji Hospital of Tongji Medical College. We recorded their clinical and pathological information, including demographics, historical types, metastasis and invasion status, BRAF V600E status, and tumor stage. Of the 70 patients, six, 36, and 26 patients were diagnosed with benign thyroid tumors, PTC, and papillary thyroid microcarcinoma (PTMC). Moreover, two patients had medullary thyroid carcinoma (MTC). We excluded these patients from the subsequent analysis because the sample size was considerably small and MTC lacked TSHR expression. All patients signed an informed consent form. All methods were approved by the Research Ethics Committee of the Tongji Medical College, and this study was conducted in accordance with declaration of Helsinki.

Gene expression profiles and data preprocessing

Transcriptomic data of WDTC were collected from TCGA database. We downloaded HTseq-counts of 510 patients with differentiated thyroid carcinomas (DTC) from TCGA website [10], which were ordered by the TSHR mRNA expression levels. We classified 100 samples of the lowest TSHR expression into the “low-TSHR” group; contrarily, 100 samples of the highest TSHR expression were classified into the “high-TSHR” group. Eventually, we identified significant DEGs using the EdgeR package in R (The R Project for Statistical Computing). DEGs related to KCNJ16 were available by comparing the “high-KCNJ16” and “low-KCNJ16” group. The dividing line represented the median expression level of KCNJ16. The criteria were set as P = 0.05 and |log2 fold change (FC)| = 2.

We obtained RNA sequencing data of ATC from the National Center for Biotechnology Information GEO. mRNA expression data were extracted and merged from GSE29265, GSE33630, GSE65144, and GSE76039 (Additional Table 1). We performed normalization and batch-normalization using an Sva package of the R Project. Considering the data scale, we selected 10 samples of the lowest and highest expression of TSHR mRNA as the “low-TSHR” and “high-TSHR” groups, respectively. DEGs were identified using the Limma package, which processed the data of both groups. We used criteria similar to those described above. Gene Ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG) and Gene Set Enrichment Analyses (GSEA) of DEGs were performed using the clusterProfiler package of the R Project. The threshold for statistical significance was set at a P-value of 0.05.

In addition, we performed a gene expression profiling interactive analysis (GEPIA) [11], and used cBioPortal for Cancer Genomics [12] to determine the correlation between the two targeted genes. Gene expression in different phenotypes of thyroid cancer was from UALCAN. TIMER and GSCA database were used to analyze immune infiltration, Survival curve of target genes in thyroid cancer was from Kaplan-Meier Plotter. MEXPRESS supplied analysis of DNA methylation status.

Total RNA extraction and real-time polymerase chain reaction

Fresh surgical specimens were immersed in RANwait (Biosharp, Anhui, China), a storage reagent to protect RNA from degradation. Total RNA from thyroid tissues was extracted using the TRIzol reagent (CWBIO, Jiangsu, China) according to the manufacturer’s instructions. We converted RNA (2ug) to cDNA using the Hifair III 1st Strand cDNA Synthesis SuperMix for quantitative PCR (qPCR) (YEASEN, Shanghai, China). Real-time PCR was performed using the ChamQ Universal SYBR qPCR Master Mix (Vazyme, Nanjing, China). The primers were synthesized by Tsingke (Beijing, China; Additional Table 1). Upon pairing the tumor tissues and their adjacent paratumor tissues from different patients, we calculated the 2−ΔCt (ΔCT = CT(target gene) - CT(GAPDH)) of the target genes to compare the expression levels of the tumor tissues and paired normal tissues. Quantitative data were obtained from three replicate experiments.

Construction of the plasmid and stable transfected cell lines

The nucleotide sequence of the KCNJ16-CDS was extracted from the CCDS database [13]. The CDS of KCNJ16 was provided in Additional File 1. KCNJ16-CDS was inserted into the over-expression vector PLVX-EF1α-IRES. The lentiviral vector was packaged using the VSVG and Gag plasmids. PTC cell line K1 and FTC cell line FTC133 were transfected with the lentivirus, following which the stably transfected cell lines were screened using puromycin.

Artificial intelligence-enabled deep docking

The REAL library of enamine can be accessed through a graphical interface [14]. We used the DD docking protocol as described in part A [15] and narrowed down the library running VirtualFlow for the final docking. The PDB file of Kir5.1, from the A.I. system AlphaFold, was used to predict its 3-D structure [16, 17]. We performed receptor preparation using AutoDock Vina v.1.2.0, and the grid box parameters were set as the x center = 25.073, y center = 15.936, and z center = -41.275. The results were visualized using PyMOL v2.5.3. The programs were used with the majority of the default docking parameters. Runtime environment: ubuntu 20.04; slurm; open shserver.

Statistical analysis

Statistical calculations and the visualization of gene expression were performed using GraphPad Prism (Version8.0.2). Initially, we performed normality and lognormality tests. We performed paired t-tests to compare targeted gene expression in the tumor tissues and corresponding normal tissues. Spearman’s correlation and linear regression analyses depicted the associations between the two genes. To describe the distribution of TSHR and KCNJ16 expression in different clinicopathological features, we performed the Student’s t-test and Mann-Whitney test for two groups of normally distributed variables and non-normally distributed variables, respectively. The Kruskal-Wallis test was performed to compare the multigroup data. Statistical significance was set at P < 0.05.

Results

DEGs screened by the TSHR expression level

We analyzed TSHR expression in the clinicopathological features of thyroid cancer (Table 1). The age and invasion status, including the capsular, vascular, and extrathyroidal extensions, were closely related to TSHR expression in tumor tissues. Patients aged ≥ 55 years and those with cancer extension revealed decreased TSHR mRNA expression. TSHR, as a differentiation indicator, had varying expression in different phenotypes of thyroid cancer. The correlation between TSHR expression and the two prognostic indicators indicated the role of TSHR in thyroid cancer. Therefore, the differentially expressed genes were analyzed according to the expression level of TSHR. DEGs were compared between the “high-TSHR” and “low-TSHR” groups.

Figure 1 A depicted the flow chart of the DEG screening. Data were collected from 510 patients with DTC and 52 patients with ATC from TCGA and GEO database, respectively. Finally, five genes were differentially expressed simultaneously by integrating 50 DEGs from TCGA with 27 DEGs from GEO, namely KCNJ16, SLC26A4, TG, TPO, and SYT1 (Fig. 1B and D). KCNJ16 and SYT1 were our target genes because SLC26A4, TG, and TPO were the classic genes in thyroid and the downstream of the TSH-TSHR signaling pathway [2]. As in Fig. 2, we analyzed the association between KCNJ16, SYT1, and TSHR using GEPIA and cBioPortal. The mRNA expression of KCNJ16 revealed a positive correlation with TSHR (R = 0.51, P < 0.05), which was consistent with a previous result that KCNJ16 was upregulated in the “high-TSHR” group. SYT1 expression still displayed a negative correlation with TSHR (R = -0.43, P < 0.05). Binary interactions between TSHR and SYT1 was observed in the UniProt database.

Fig. 1
figure 1

DEGs screened by the TSHR expression level

(A) The flow chart of DEG screening process from TCGA database and GEO database. (B) The expression regulation of 50 DEGs between “high-TSHR” group and “low-TSHR” group in thyroid cancer patients from TCGA database. (C) The expression regulation of 27 DEGs between “high-TSHR” group and “low-TSHR” group in ATC patients From GEO database. (D) Venn diagram of the two DEGs datasets. DEG: differentially expressed gene, TSHR: thyroid stimulating hormone receptor, ATC: anaplastic thyroid cancer

Fig. 2
figure 2

The primary analysis of correlation between TSHR and DEGs.

(A) The correlation analysis between KCNJ16-TSHR mRNA expression and SYT1-TSHR mRNA expression on the cBioPortal website (http://www.cbioportal.org/). (B) The correlation analysis between KCNJ16-TSHR mRNA expression and SYT1-TSHR mRNA expression on the GEPIA website (http://gepia.cancer-pku.cn/). (C) The binary interaction between TSHR and SYT1 on the Uniport protein database (https://www.uniprot.org/)

Table 1 Association of TSHR mRNA expression with different variables of thyroid tumor patients

Expression analysis and prognostic value of DEGs in thyroid cancer

We used TCGA database in UALCAN to investigate DEGs expression in thyroid cancer. TSHR was decreased in primary tumor and in any stages or histological subtypes when compared to normal group, indicating it was a strong biomarker of thyroid cancer (Fig. 3A). The down-regulation of KCNJ16 was significant in stage 4 and tall thyroid papillary carcinoma (Fig. 3B). KCNJ16 expression change might come from DNA methylation, cause its expression in thyroid cancer was negatively correlated with DNA promoter methylation and CpG hypermethylation appeared in metastatic group and thyroid papillary carcinoma - tall cell ( > = 50% tall cell features) (Additional Fig. 1A and B). Tall cell morphology was an aggressive variant of papillary thyroid carcinoma that has been associated with poor outcomes [18]. SYT1 expression increased in primary tumor, almost in all stages when compared to normal group. And SYT1 was also observed in tall thyroid papillary carcinoma (Fig. 3C).

Fig. 3
figure 3

Expression of TSHR, KCNJ16 and SYT1 in thyroid cancer phenotypes from the UALCAN.

(A) Box plots of TSHR transcripts in sample types, stages nodal metastasis status and histological subtypes of thyroid cancer from TCGA database. (B) Box plots of KCNJ16 transcripts in thyroid cancer. (C) Box plots of SYT1 transcripts in thyroid cancer. N0: No regional lymph node metastasis, N1: Metastases in 1 to 3 lymph nodes, PTC: papillary thyroid carcinoma, FTC: follicular thyroid carcinoma, *: P < 0.05

Then we applied KM plotter to investigate their prognostic value in thyroid cancer, TIMER and GSCA to analyze immune cell infiltration and survival. From Fig. 4, high expression of TSHR (HR = 0.21, logrank P = 0.017) and low expression of SYT1 (HR = 3.55, logrank P = 0.00058) were significantly correlated with favorable RFS. KCNJ16 supported good prognosis of OS (HR = 0.21, logrank P = 0.00075) and RFS (HR = 0.41, logrank P = 0.018). The tumor immune microenvironment of thyroid cancer was composed of cancer cells, tumor-associated fibroblasts and immune cells. Immune cells have complicated roles in tumor growth, differentiation, angiogenesis and migration [19]. Changes in TSHR, KCNJ16 and SYT1 expression level were accompanied by the alteration of immune cells infiltration. But the cumulative survival curves showed that the infiltration of B cell, CD8 + T cell, CD4 + T cell, macrophage, neutrophil and dendritic cell didn’t cause a significant change in survival time of thyroid cancer patients (Additional Fig. 2A). But from the GSCA database, TSHR and KCNJ16 shared a similar correlation with immune infiltrates. And SYT1 expression was positively correlated with the infiltration of iTreg, nTreg and macrophage (Additional Fig. 2B).

Fig. 4
figure 4

Survival curve of TSHR, KCNJ16 and SYT1 in thyroid cancer patients from Kaplan-Meier plotter

(A) The OS probability of high expression and low expression of TSHR, KCNJ16 and SYT1 in thyroid cancer patients. (B) The RFS probability of high expression and low expression of TSHR, KCNJ16 and SYT1 in thyroid cancer patients. Logrank P indicates statistical significance. OS: overall survival, RFS: recurrence-free survival, HR: hazard ratio

DEGs expression in clinical samples and their role in thyroid cancer

RT-qPCR experiments were performed to verify the mRNA expression levels of the identified DEGs. We analyzed mRNA expression levels in 68 thyroid tumor samples and adjacent normal tissues. The expression of TSHR and KCNJ16 was lower in the tumor tissues than that in the paired normal tissues, whereas the expression of SYT1 was upregulated in the tumor tissues (Fig. 5A C). The correlations were nearly consistent with the prior analysis. KCNJ16 mRNA expression was positively correlated with TSHR mRNA expression (Pearson r = 0.7734, P < 0.05). However, our results did not demonstrate an association between SYT1 and TSHR expression (Pearson r = 0.045, P = 0.7283) (Fig. 5D).

Fig. 5
figure 5

Target genes expression in clinical samples and DEGs enriched pathways

(A) Individual expression and median expression of TSHR in paratumor and tumor tissue. (B) Individual expression and median expression of KCNJ16 in paratumor and tumor tissue. (C) Individual expression and median expression of SYT1 in paratumor and tumor tissue. (D) The correlation analysis of KCNJ16-TSHR expression and SYT1-TSHR expression of clinical samples. (E) The enriched pathways of DEGs in “high-KCNJ16” group of thyroid cancer patients from TCGA database. (F) The GSEA analysis of listed pathways in ATC patients from GEO database. The running enrichment score from − 0.25 to 0.75. DEG: differentially expressed gene, GSEA: gene set enrichment analysis, ATC: anaplastic thyroid cancer. *: P < 0.05

Subsequently, we explored KCNJ16 expression in different clinical and pathological features of thyroid cancer (Table 2). Conspicuously, KCNJ16 expression was decreased in patients with capsular and vascular invasion or an extrathyroidal extension invading the larynx, trachea, or esophagus, in the absence of significant changes among other clinicopathologic features. Simultaneously, we observed the downregulation of TSHR expression in patients with capsular/vascular invasion (Table 1). Both differences were statistically significant.

Table 2 Association of KCNJ16 mRNA expression with different variables of thyroid tumor patients

Considering the coherent correlation between KCNJ16 and TSHR expression and their tendency to decline in tumor tissues, we focused on the potential role of KCNJ16 in thyroid cancer. Gene expression quantification of thyroid cancer in TCGA database was extracted and divided into the “high-KCNJ16” and “low-KCNJ16” groups by the median expression. We screened 614 DEGs, with 52 upregulated genes and 562 downregulated genes. Functional and pathway enrichments were implemented (Fig. 5E). The results depicted the enriched biological processes related to cell growth and differentiation, including cell growth (GO:0016049), epithelial cell proliferation (GO:0050673), columnar/cuboidal epithelial cell differentiation (GO:0002065). Moreover, DEGs were enriched in EMT (GO:0001837) and cell-cell adhesion (GO:0022408). KEGG pathways included cell adhesion molecules and ECM-receptor interactions (GO and KEGG pathways were listed in the Additional Tables 2 and 3, respectively)[20,21,22].

To further explore KCNJ16 and these enriched pathway in ATC, the sample data of ATC were downloaded from the GEO database. DEGs were received from “high-KCNJ16” group when compared to “low-KCNJ16” group. KCNJ16 associated DEGs in ATC enriched in some similar pathways. And GSEA analysis was performed to observe up-regulation or down-regulation of those interested pathways (Fig. 5F). Pathways of “thyroid gland development”, “positive regulation of epithelial cell differentiation” and “negative regulation of epithelial to mesenchymal transition” were upregulated in “high-KCNJ16” group, indicating that KCNJ16 might promote differentiation and inhibit migration. Although “regulation of MAP kinase activity” pathway was up-regulated, “ERK1 and ERK2 cascade” and NF-κB signaling pathway were down-regulated. Normalized enrichment score and P-value of the interested pathways were listed in additional Table 4. These signaling pathways were involved not only in the normal physiological function and cell differentiation of the thyroid, but also in the invasive changes of thyroid cancer [23,24,25]. Therefore, the role of KCNJ16 and enriched pathways in thyroid cancer need to be further investigated.

Kir5.1 (encoded by KCNJ16) reportedly forms a heterotetramer with Kir4.1 (encoded by KCNJ10); however, the expression of KCNJ10 in thyroid tissues was previously lower than that of KCNJ16. KCNJ15 and KCNJ16 were highly expressed in thyroid tissues (Additional Fig. 3). We hypothesized that, Kir4.2/Kir5.1 was the functional heterotetramer in thyroid tissues, instead of the classical Kir4.1/Kir5.1. This hypothesis was supported by the observation that KCNJ10 was seldom expressed in thyroid samples, whereas KCNJ15 and KCNJ16 were synchronously expressed in every sample (Additional Fig. 4A and B). Thus, we concluded that KCNJ16 might affect the differentiation and development of thyroid cancer alone or by the Kir4.2/Kir5.1 heterotetramer.

To determine the interaction between KCNJ16 and TSHR, over-expression vectors of KCNJ16 were constructed and transfected into the PTC cell line K1, FTC cell line FTC133 and ATC cell line CAL62. In K1, FTC133 and CAL62, KCNJ16 was significantly over-expressed following transfection (Additional Fig. 4C). But no significant increase in TSHR expression was observed in K1 and FTC133 when KCNJ16 was over-expressed. Only in CAL62, we observed a nearly eight-fold increase in TSHR expression (Additional Fig. 4D). Then we used flow cytometry to detect TSHR expression on membrane. In three cell lines, membrane TSHR expression was not detected, even in CAL62 with KCNJ16 over-expression (Additional Fig. 4E).

Using the deep docking platform and VFVS to screen ligands against Kir5.1

First, we narrowed down a subset library of Real to approximately 1% using the deep neural network (Fig. 6A). Eventually we considered 2 072 compounds for the final docking. Every ligand tried nine poses using two scoring methods (qvina02 and smina), and we ranked the maximum affinity value. The A.I. system predicted the 3-D structure of Kir5.1 (Fig. 6B). The PDB file of Kir5.1 was processed with Autodock vina to add polar hydrogen and set the grid box for the interaction pocket (Fig. 6C). Figure 6D depicted an example of the superimposition of the docking compounds-Kir5.1 from the crystal structure. The corresponding scores of the ligands were listed in Table 3 (leading 10) and Additional Table 5 (leading 100). Consequently, we selected the leading 10 candidates to analyze the docking conformations in the best pose out of nine (Fig. 6E). Maximum candidates formed hydrogen bonds at ASP101, ASP113, and ASP105, which belong to the extracellular domain near TM1. Contrarily, hydrogen bonds usually interacted with GLU141, GLU142, and VAL145, which belong to the extracellular domain near TM2, thereby supporting the hypothesis that a series of three amino acids in this region were important for channel gating. Overall, 11 residues, including ASP113, ASP101, ASP105, GLU141, GLU142, HIS116, THR109, VAL145, ILE108, ASN114, and PRO106, were identified important for Kir5.1 binding.

Fig. 6
figure 6

Virtual screening for ligands against Kir5.1

(A) A schematic of the virtual screening workflow. (B) Predicted crystal structure of Kir5.1. Model confidence indicates the degree of similarity with the true structure. (C) Visualization of the target docking position and size of Kir5.1 (inside the box), processed by Autodock Vina. (D) Diagram of interaction position between docking compound and Kir5.1. Structure domain in purple is the box section of Kir5.1 in (C). Indigo structure is the screened compound. (E) Docking poses of the leading 10 candidate compounds. Yellow dotted lines represent hydrogen bonds. Indigo structures are the candidate compounds. Green structures are hydrogen-bonded amino acids of Kir5.1. Kir5.1: protein encoded by KCNJ16.

Table 3 The top 10 ligands against Kir5.1

Discussion

In this study, we focused on TSHR and TSHR related DEGs, KCNJ16 and SYT1. In thyroid cancer, TSHR expression was reduced in primary tumor and in any stages or histological subtypes, indicating its role as a biomarker. The two DEGs showed opposite tendency in the regulation of gene expression in thyroid cancer. The mRNA expression level of SYT1 was higher in tumor tissues, while the expression of KCNJ16 was positively correlated with that of TSHR and was decreased in tumor tissues. TSHR and KCNJ16 demonstrated a relationship with capsular/vascular invasion, and they showed similar function in immune infiltration. Furthermore, GO and KEGG enrichment analyses demonstrated that KCNJ16 may play a significant role in cell growth and differentiation, EMT, MAPK signaling pathway and NF-κB signaling pathway. Targeting Kir5.1 (encoded by KCNJ16), Z2087256678_2, Z2211139111_1, Z2211139111_2, and PV-000592319198_1 (-7.3 kcal/mol) were identified as the most potent commercially available molecular. Other molecular namely, Z2211137992_1, Z2211137992_2 (-7.2 kcal/mol), Z2717222271_1 (-7.1 kcal/mol), Z2211137992_4 (-7.0 kcal/mol), Z1411107749_2 (-7.0 kcal/mol), and PV-000377841410_1 (-6.9 kcal/mol) also revealed potent interactions with Kir5.1. Herein, we reported on Kir5.1 as a potential biomarker for more invasive thyroid cancer, and redifferentiation therapeutic targets for recurrence/metastases or PDTC and ATC patients. And we successfully developed and characterized 10 potent Kir5.1 interaction compounds through AI-assisted virtual screening.

TSHR is differentially expressed in the thyroid cancer cells from normal tissues or benign diseases [26,27,28]. Generally, our results were in accordance with other reports stating TSHR was expressed at lower levels in cancer tissues. A prior study reported that approximately half of the samples (44.2%) presented a lower staining intensity in lymph node metastases than that in their corresponding primary tumors [29]. However, our results revealed that lymph node metastases did not affect TSHR expression in the primary sites, thus indicating its expression was heterogeneously distributed. Recently, thyrocytes from subcutaneous loci belonging to distant metastases in PTC displayed lower degrees of differentiation (assessed by the TDS score, which considers TSHR expression as one of the indicators) and preferentially upregulated dedifferentiation-related transcription factors and pathways. Unsurprisingly, these samples were obtained from patients who had progressed to the RAI-refractory (RAIR) status [30]. Furthermore, the normal expression and function of TSHR was essential for the character of thyroid cancer cells and the prognosis of patients with an RAIR status.

Potassium channel regulation may influence tumor progression via multiple pathways, such as cell adhesion or migration, angiogenesis and apoptosis. Most Kir subunits are overexpressed in multiple types of cancers. Kir2.1 (KCNJ2) was overexpressed in 44.23% of small-cell lung cancer (SCLC) tissues and is correlated with the clinical stage and chemotherapy response in patients with SCLC [31]. Kir3.1 was another protein associated with cancer, overexpressed in non-small cell lung cancers (NSCLCs), 40% of primary human breast cancers, and breast cancer cell lines [32, 33]. A functional and pathway enrichment analysis revealed that KCNJ16 was associated with epithelial cell proliferation and differentiation, thyroid hormone generation and metabolism, epithelial to mesenchymal transition, and cell-cell adhesion. In other words, KCNJ16 may be a regulator of thyroid differentiation and cancer migration. Consistent with our findings, KCNJ16 has been reported in other cancers. Its downregulation was observed in a hepatocellular carcinoma group and PDAC [34]. And from the Human Protein Atlas, the protein expression score of KCNJ16 was similar in thyroid gland and kidney. Interesting thing was that KCNJ16 was a favorable prognostic marker in renal cancer and thyroid cancer (data not shown). KCNJ16 had the diagnostic, prognostic, and therapeutic value in clear cell renal cell carcinoma [35]. Thus, KCNJ16 may play a protective role in cell differentiation and cancer progression.

KCNJ16 expression was positively correlated with that of TSHR in our study. Particularly, the expression significantly declined in patients with vascular/capsular invasion. Moreover, the over-expression of KCNJ16 upregulated the expression of TSHR in ATC cell line CAL62. But the CT value of TSHR mRNA expression was over 30, indicating its low expression. In sharp contrast to the ubiquitous expression of TSHR in DTC clinical samples, TSHR expression appears to be absent in most DTC-derived cell lines [36]. Many literatures have reported that thyroid cancer cell lines lost TSHR expression. The absence of this protein in commonly used DTC cell lines could be attributed to gene mutation and culture-specific selective pressures in vitro and has been well documented in the literature [37,38,39]. And in our study, we did not detect TSHR membrane expression in cell lines, even in CAL62 with escalated TSHR mRNA expression. Both low mRNA expression and complex translation process from RNA to protein could result in no-detection. The regulation role of KCNJ16 will need more investigation in primary thyroid cell, cause cell lines have lost most differentiation features.

The complex regulatory mechanisms of Kir5.1 in cancer progression remain largely unknown. In a previous study, the potassium channel subunits, namely KCNQ1 and KCNE2, formed a thyroid-stimulating hormone-stimulated, constitutively active, thyrocyte K + channel. KCNQ1 and KCNE2 are required for normal thyroid hormone biosynthesis, and this K + voltage affects the NIS [40, 41]. So according to KCNQ1-KCNE2 K + channel was necessary for iodine uptake, the Kir may also play the similarly role in thyroid tissues. However, this hypothesis requires further investigation. The Kir4.1, Kir4.2, Kir5.1 mRNA levels had significantly increased during mouse thyroid gland development [42]. Targeting Kir5.1 may be a promising strategy for patients with an RAIR status. Samantha et al. performed a novel study comprising a wet-bench high-throughput screening of 80 475 compounds as the potent inhibitors of Kir4.1/kir5.1 heterotetramer [43]. On an average, a novel potent drug costs billions of dollars and takes more than decades to develop, particularly comprising time-consuming wet-bench experiments. Therefore, computer-assisted drug discovery can reduce the cost of ultra-large screening.

Our study had some limitations. First, our sample capacity was insufficient, specifically owing to the lack of DeTC and ATC specimens. We could not confirm DEGs expression in DeTC and ATC, despite considering ATC while selecting GEO datasets. Moreover, we verified the expression of DEGs by RT-PCR, which could only reflect the relative quantification of DEGs, but not their location or distribution in cells, thus necessitating further studies. In addition, we explored the possible mechanisms between the targeted genes using the GO and KEGG enrichment analysis. The role of these genes in thyroid differentiation and cancer progression requires further verification. Even with the accelerated A.I., our initial library contained only approximately 200 000 compounds. Our library was relatively small than the largest commercial library, i.e., Real comprising 1.3 billion compounds. While A.I. is a promising tool for identifying potential derivatives of a target and improving drug development efficiency, the actual effects of drugs and their applicability to various diseases must be evaluated through basic and clinical trials. Therefore, the conclusions of this study can only be considered suggestive, and generalization significance is lacking. Validation of the candidates from our screening through wet-bench experiments is still necessary.

In summary, we identified five common DEGs related to TSHR mRNA expression in PTC and ATC using an integrated microarray analysis. TSHR expression was closely associated with the clinicopathological characteristics, particularly the age, historical type, stage, and invasion status. KCNJ16 expression was consistent with TSHR expression, and KCNJ16 may affect thyroid differentiation and malignant transformation by regulating TSHR expression. A.I. assistant virtual screening demonstrated that Z2087256678_2, Z2211139111_1, Z2211139111_2, PV-000592319198_1, Z2211137992_1, Z2211137992_2, Z2717222271_1, Z2211137992_4, Z1411107749_2, and PV-000377841410_1 had potent interactions with Kir5.1. This study may provide greater insights into the dedifferentiation features associated with TSHR expression in thyroid cancer and Kir5.1, which may be a potential therapeutic target in redifferentiation strategies for recurrent and metastatic thyroid cancer, thus enabling the adoption of RAI therapy and other treatments.

Data Availability

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

Abbreviations

TSHR:

Thyroid stimulating hormone receptor

TC:

Thyroid cancer

PTC:

Papillary thyroid carcinoma

FTC:

Follicular thyroid carcinoma

ATC:

Anaplastic thyroid carcinoma

WDTC:

Well-differentiated thyroid cancers

DeTCs:

Dedifferentiated thyroid cancers

RAI:

Radioactive iodine

NIS:

Sodium/iodide symporter

TG:

Thyroglobulin

TPO:

Thyroperoxidase

RAIR:

Resistance to RAI therapy

EMT:

Epithelial-mesenchymal transition

TCGA:

The Cancer Genome Atlas

GEO:

Gene Expression Omnibus

DEGs:

Differentially expressed genes

MTC:

Medullary thyroid carcinoma

GO:

Gene Ontology

KEGG:

Kyoto Encyclopedia of Genes and Genomes

References

  1. Kim J, Gosnell JE, Roman SA. Geographic influences in the global rise of thyroid cancer. Nat Rev Endocrinol. 2020;16(1):17–29.

    Article  PubMed  Google Scholar 

  2. Xing M. Molecular pathogenesis and mechanisms of thyroid cancer. Nat Rev Cancer. 2013;13(3):184–99.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  3. Barreno LRQ, Mello JBH, Barros-Filho MC, Francisco AL, Chulam TC, Pinto CAL et al. Characterization of BRAF mutation in patients older than 45 years with well-differentiated thyroid carcinoma. Braz J Otorhinolaryngol. 2020.

  4. Bruno R, Ferretti E, Tosi E, Arturi F, Giannasio P, Mattei T, et al. Modulation of thyroid-specific gene expression in normal and nodular human thyroid tissues from adults: an in vivo effect of thyrotropin. J Clin Endocrinol Metab. 2005;90(10):5692–7.

    Article  CAS  PubMed  Google Scholar 

  5. Ludgate M, Gire V, Crisp M, Ajjan R, Weetman A, Ivan M, et al. Contrasting effects of activating mutations of GalphaS and the thyrotropin receptor on proliferation and differentiation of thyroid follicular cells. Oncogene. 1999;18(34):4798–807.

    Article  CAS  PubMed  Google Scholar 

  6. Xu S, Chen G, Peng W, Renko K, Derwahl M. Oestrogen action on thyroid progenitor cells: relevant for the pathogenesis of thyroid nodules? J Endocrinol. 2013;218(1):125–33.

    Article  CAS  PubMed  Google Scholar 

  7. Morgan SJ, Neumann S, Marcus-Samuels B, Gershengorn MC. Thyrotropin and insulin-like Growth factor 1 receptor crosstalk Upregulates Sodium-Iodide Symporter expression in primary cultures of human thyrocytes. Thyroid. 2016;26(12):1794–803.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  8. Postiglione MP, Parlato R, Rodriguez-Mallon A, Rosica A, Mithbaokar P, Maresca M, et al. Role of the thyroid-stimulating hormone receptor signaling in development and differentiation of the thyroid gland. Proc Natl Acad Sci U S A. 2002;99(24):15462–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. Rowe CW, Paul JW, Gedye C, Tolosa JM, Bendinelli C, McGrath S, et al. Targeting the TSH receptor in thyroid cancer. Endocr Relat Cancer. 2017;24(6):R191–r202.

    Article  CAS  PubMed  Google Scholar 

  10. Atlas TCG. https://portal.gdc.cancer.gov/. Accessed June 2021.

  11. Analysis GEPI. http://gepia.cancer-pku.cn/. Accessed October 2021.

  12. Genomics cfC. http://www.cbioportal.org/. Accessed October 2021.

  13. database C. https://www.ncbi.nlm.nih.gov/projects/CCDS/CcdsBrowse.cgi. Accessed September 2022.

  14. VirtualFlow. https://virtual-flow.org/. Accessed September 2022.

  15. GitHub. https://github.com/jamesgleave/DD_protocol. Accessed June 2022.

  16. Jumper J, Evans R, Pritzel A, Green T, Figurnov M, Ronneberger O, et al. Highly accurate protein structure prediction with AlphaFold. Nature. 2021;596(7873):583–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  17. Gorgulla C, Boeszoermenyi A, Wang ZF, Fischer PD, Coote PW, Padmanabha Das KM, et al. An open-source drug discovery platform enables ultra-large virtual screens. Nature. 2020;580(7805):663–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Wu SS, Joshi N, Sharrett J, Rao S, Shah A, Scharpf J, et al. Risk factors Associated with recurrence and death in patients with tall cell papillary thyroid Cancer: a Single-Institution Cohort Study with Predictive Nomogram. JAMA Otolaryngol Head Neck Surg. 2023;149(1):79–86.

    Article  PubMed  Google Scholar 

  19. Shin E, Koo JS. Cell Component and Function of Tumor Microenvironment in Thyroid Cancer. Int J Mol Sci. 2022;23(20).

  20. Kanehisa M, Goto S. KEGG: kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 2000;28(1):27–30.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  21. Kanehisa M. Toward understanding the origin and evolution of cellular organisms. Protein Sci. 2019;28(11):1947–51.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Kanehisa M, Furumichi M, Sato Y, Kawashima M, Ishiguro-Watanabe M. KEGG for taxonomy-based analysis of pathways and genomes. Nucleic Acids Res. 2023;51(D1):D587–d92.

    Article  PubMed  Google Scholar 

  23. Ramírez-Moya J, Santisteban P. miRNA-Directed regulation of the Main Signaling Pathways in thyroid Cancer. Front Endocrinol (Lausanne). 2019;10:430.

    Article  PubMed  Google Scholar 

  24. Schubert L, Mariko ML, Clerc J, Huillard O, Groussin L. MAPK Pathway Inhibitors in Thyroid Cancer: Preclinical and Clinical Data. Cancers (Basel). 2023;15(3).

  25. Boucai L, Saqcena M, Kuo F, Grewal RK, Socci N, Knauf JA et al. Genomic and transcriptomic characteristics of metastatic thyroid cancers with exceptional responses to radioactive iodine therapy. Clin Cancer Res. 2023.

  26. Cordioli MI, Moraes L, Alves MT, Delcelo R, Monte O, Longui CA, et al. Thyroid-specific genes expression uncovered Age-Related differences in Pediatric thyroid carcinomas. Int J Endocrinol. 2016;2016:1956740.

    Article  PubMed  PubMed Central  Google Scholar 

  27. Li J, Dong JN, Zhao Z, Lv Q, Yun B, Liu JQ, et al. Expression of sodium/iodide transporters and thyroid stimulating hormone receptors in thyroid cancer patients and its correlation with iodine nutrition status and pathology. Eur Rev Med Pharmacol Sci. 2018;22(14):4573–80.

    CAS  PubMed  Google Scholar 

  28. Ravi N, Yang M, Mylona N, Wennerberg J, Paulsson K. Global RNA Expression and DNA Methylation Patterns in Primary Anaplastic Thyroid Cancer. Cancers (Basel). 2020;12(3).

  29. So YK, Son YI, Baek CH, Jeong HS, Chung MK, Ko YH. Expression of sodium-iodide symporter and TSH receptor in subclinical metastatic lymph nodes of papillary thyroid microcarcinoma. Ann Surg Oncol. 2012;19(3):990–5.

    Article  PubMed  Google Scholar 

  30. Pu W, Shi X, Yu P, Zhang M, Liu Z, Tan L, et al. Single-cell transcriptomic analysis of the tumor ecosystems underlying initiation and progression of papillary thyroid carcinoma. Nat Commun. 2021;12(1):6058.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  31. Liu H, Huang J, Peng J, Wu X, Zhang Y, Zhu W, et al. Upregulation of the inwardly rectifying potassium channel Kir2.1 (KCNJ2) modulates multidrug resistance of small-cell lung cancer under the regulation of miR-7 and the Ras/MAPK pathway. Mol Cancer. 2015;14:59.

    Article  PubMed  PubMed Central  Google Scholar 

  32. Rezania S, Kammerer S, Li C, Steinecker-Frohnwieser B, Gorischek A, DeVaney TT, et al. Overexpression of KCNJ3 gene splice variants affects vital parameters of the malignant breast cancer cell line MCF-7 in an opposing manner. BMC Cancer. 2016;16:628.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Takanami I, Inoue Y, Gika M. G-protein inwardly rectifying potassium channel 1 (GIRK 1) gene expression correlates with tumor progression in non-small cell lung cancer. BMC Cancer. 2004;4:79.

    Article  PubMed  PubMed Central  Google Scholar 

  34. Jiang S, Zhu L, Yang J, Hu L, Gu J, Xing X, et al. Integrated expression profiling of potassium channels identifys KCNN4 as a prognostic biomarker of pancreatic cancer. Biochem Biophys Res Commun. 2017;494(1–2):113–9.

    Article  CAS  PubMed  Google Scholar 

  35. Zhu Z, Lei Z, Qian J, Zhang C, Gong Y, Yin G, et al. The Ion Channel-Related gene signatures correlated with diagnosis, prognosis, and Individualized Treatment in patients with Clear Cell Renal Cell Carcinoma. Front Pharmacol. 2022;13:889142.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Li H, Zhou X, Wang G, Hua D, Li S, Xu T, et al. CAR-T cells targeting TSHR demonstrate safety and potent preclinical activity against differentiated thyroid Cancer. J Clin Endocrinol Metab. 2022;107(4):1110–26.

    Article  PubMed  Google Scholar 

  37. van Staveren WC, Solís DW, Delys L, Duprez L, Andry G, Franc B, et al. Human thyroid tumor cell lines derived from different tumor types present a common dedifferentiated phenotype. Cancer Res. 2007;67(17):8113–20.

    Article  PubMed  Google Scholar 

  38. Pilli T, Prasad KV, Jayarama S, Pacini F, Prabhakar BS. Potential utility and limitations of thyroid cancer cell lines as models for studying thyroid cancer. Thyroid. 2009;19(12):1333–42.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  39. Meireles AM, Preto A, Rocha AS, Rebocho AP, Máximo V, Pereira-Castro I, et al. Molecular and genotypic characterization of human thyroid follicular cell carcinoma-derived cell lines. Thyroid. 2007;17(8):707–15.

    Article  CAS  PubMed  Google Scholar 

  40. Roepke TK, King EC, Reyna-Neyra A, Paroder M, Purtell K, Koba W, et al. Kcne2 deletion uncovers its crucial role in thyroid hormone biosynthesis. Nat Med. 2009;15(10):1186–94.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Purtell K, Paroder-Belenitsky M, Reyna-Neyra A, Nicola JP, Koba W, Fine E, et al. The KCNQ1-KCNE2 K+ channel is required for adequate thyroid I- uptake. Faseb j. 2012;26(8):3252–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Ramos HE, da Silva MR, Carré A, Silva JC Jr, Paninka RM, Oliveira TL, et al. Molecular insights into the possible role of Kir4.1 and Kir5.1 in thyroid hormone biosynthesis. Horm Res Paediatr. 2015;83(2):141–7.

    Article  CAS  PubMed  Google Scholar 

  43. McClenahan SJ, Kent CN, Kharade SV, Isaeva E, Williams JC, Han C, et al. VU6036720: the First Potent and Selective in Vitro inhibitor of Heteromeric Kir4.1/5.1 Inward Rectifier Potassium channels. Mol Pharmacol. 2022;101(5):357–70.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

We thank all members of the Thyroid and Breast Surgery Laboratory as well as various colleagues in the Department of General Surgery, Tongji Hospital, Tongji Medical College, Huazhong University of Science and Technology.

Funding

This work was supported by the National Natural Science Foundation of China (No.81802676) and the Wuhan Youth Cadre Project (2017zqnlxr01 and 2017zqnlxr02).

Author information

Authors and Affiliations

Authors

Contributions

YD and XL conceived of and designed the study. XY analyzed, interpreted the data and wrote the manuscript. HL, TX, MD, DKD, and CP provided the administrative support. GW, YW, SX, and XC provided the study materials and recruited the patients. All authors reviewed the manuscript.

Corresponding authors

Correspondence to Xingrui Li or Yaying Du.

Ethics declarations

Competing interests

The authors declare that they have no competing interests.

Ethics approval and consent to participate

All human samples were analyzed according to the guidelines of the Medical Ethics Committee of Tongji Hospital, Tongji Medical College, Huazhong University of Science and Technology (Institution Review Board Approval: TJ-C20180110). Each patient was requested to sign a written informed consent form, and the specimens were anonymized and handled according to the accepted ethical and legal standards.

Consent for publication

Not applicable.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Electronic supplementary material

Below is the link to the electronic supplementary material.

Supplementary Material 1

Supplementary Material 2

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Yang, X., Wu, Y., Xu, S. et al. Targeting the inward rectifier potassium channel 5.1 in thyroid cancer: artificial intelligence-facilitated molecular docking for drug discovery. BMC Endocr Disord 23, 113 (2023). https://doi.org/10.1186/s12902-023-01360-z

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12902-023-01360-z

Keywords