Drug repositioning in non-small cell lung cancer (NSCLC) using gene co-expression and drug-gene interaction networks analysis.

Lung cancer is the most common cancer in men and women. This cancer is divided into two main types, namely non-small cell lung cancer (NSCLC) and small cell lung cancer (SCLC). Around 85 to 90 percent of lung cancers are NSCLC. Repositioning potent candidate drugs in NSCLC treatment is one of the important topics in cancer studies. Drug repositioning (DR) or drug repurposing is a method for identifying new therapeutic uses of existing drugs. The current study applies a computational drug repositioning method to identify candidate drugs to treat NSCLC patients. To this end, at first, the transcriptomics profile of NSCLC and healthy (control) samples was obtained from the GEO database with the accession number GSE21933. Then, the gene co-expression network was reconstructed for NSCLC samples using the WGCNA, and two significant purple and magenta gene modules were extracted. Next, a list of transcription factor genes that regulate purple and magenta modules' genes was extracted from the TRRUST V2.0 online database, and the TF-TG (transcription factors-target genes) network was drawn. Afterward, a list of drugs targeting TF-TG genes was obtained from the DGIdb V4.0 database, and two drug-gene interaction networks, including drug-TG and drug-TF, were drawn. After analyzing gene co-expression TF-TG, and drug-gene interaction networks, 16 drugs were selected as potent candidates for NSCLC treatment. Out of 16 selected drugs, nine drugs, namely Methotrexate, Olanzapine, Haloperidol, Fluorouracil, Nifedipine, Paclitaxel, Verapamil, Dexamethasone, and Docetaxel, were chosen from the drug-TG sub-network. In addition, nine drugs, including Cisplatin, Daunorubicin, Dexamethasone, Methotrexate, Hydrocortisone, Doxorubicin, Azacitidine, Vorinostat, and Doxorubicin Hydrochloride, were selected from the drug-TF sub-network. Methotrexate and Dexamethasone are common in drug-TG and drug-TF sub-networks. In conclusion, this study proposed 16 drugs as potent candidates for NSCLC treatment through analyzing gene co-expression, TF-TG, and drug-gene interaction networks.

Introduction

Lung cancer is one of the leading cancer death causes[1] worldwide. This type of cancer occurs when a cancerous tumor grows inside the lungs. Lung cancer contains two main types: non-small cell lung cancer (NSCLC) and small cell lung cancer (SCLC). NSCLC is the most common lung cancer[2]. Histopathological grading has identified about 85% to 90% of lung cancers as NSCLC and 15% to 20% as SCLC[3]. This cancer includes three different types of Adenocarcinoma, Squamous cell carcinoma, and large cell carcinoma.

Different studies based on computational approaches and network analysis have been undertaken to find biomarker genes for early NSCLC detection. Moreover, scientists have evaluated and discussed the effect of current drugs on this cancer. Based on a co-expression network analysis, Ling Kui et al.4 proposed several important genes as biomarkers for NSCLC treatment. Xiujuan Gao et al.[4] applied the gene expression profile of NSCLC samples and, based on a systems biology approach, reported Estrogen receptors (ERs) as promoters of NSCLC progression. In another study, Mei Zhao et al.[5] introduced five genes, including FGF2, GOLM1, GPC3, IL6, and SPP1, which deregulated in NSCLC tissues. They introduced these 5 genes for NSCLC prognosis in patients. A computational approach based on protein–protein interaction (PPI) network analysis was used in a similar study, and Stratifin had an important role in NSCLC[6] development. Furthermore, Yun-Qiang Zhang et al.[7] proposed HIST1H2BH and PLK1 as prognostic biomarkers for NSCLC patients.

Drug repositioning (DR) is utilized as a time- and cost-effective method to discover new drugs[8-12]. Drug repositioning is also referred to as drug repurposing, drug therapeutic, drug recycling, and drug reprofiling[13,14]. There are usually three kinds of methods for drug repurposing, including experimental biological methods, computational methods, and mixed methods[10,15,16]. Computational methods can be referred to as molecular docking, network mapping, signature matching, genetic association, and retrospective clinical analysis[13,17,18]. In the current study, a computational drug repositioning method is applied to identify candidate drugs to treat NSCLC.

Lately, various studies based on the network approach for drug repurposing have been carried out. Network-based strategy is one of the important computational methods in drug repurposing[19,20]. SAveRUNNER[21] is a network-based algorithm in this field. This algorithm predicts drug-disease relations based on a similarity measure. This method was provided as an R programing language package[22]. Xing Li and colleagues [23] proposed a network-based approach to discover lncRNA biomarkers in human lung adenocarcinoma. Furthermore, a computational approach for drug repurposing based on the system biology approach was proposed by Azam Peyvandipour and collages [24] in 2018. In another study, Wei-Feng Guo et al.[25] proposed a network controllability-based algorithm called combinatorial drug identification algorithm (CPGD). Besides, Albert Li and collages[26] proposed a network-based method, namely LncTx, to repurpose drugs in lung cancer. In a recent study by Zahra and her colleagues[27], they proposed a novel network-based method to discover candidate drugs for bladder cancer.

Anisha et al. [28] presented an overview of drug repositioning for anti-cancer applications, and they proposed a novel drug repurposing technique to target the MAPK signaling pathway in NSCLC. In a similar study, Muthu Kumar and colleagues[29] introduced another drug repurposing method for NSCLC, and they hypothesized that Nebivolol is an excellent candidate for inhibiting MEK in NSCLC patients. In another study, Joelle C. Boulos and colleagues[30] repurposed ALK Inhibitor Crizotinib for NSCLC, Acute Leukemia, and Multiple Myeloma Cells. Compared to the mentioned methods, this study applies a novel computational model based on gene co-expression and TF–TG interaction networks. Moreover, two drug–gene interaction networks, including drug-TF and drug-TG, were studied that have not been studied in previous studies.

Gene co-expression network analysis is one of the important network-based approaches in systems biology[11,31,32]. Different studies based on gene co-expression network analysis were done on different transcriptomic datasets. In this study, a gene co-expression network analysis was applied on the NSCLC transcriptomics dataset to repurpose some potent candidate drugs for NSCLC. Xue-Tao Li [33] and colleagues applied gene co-expression modules analysis in order to predict non-small cell lung cancer survivals. In a similar project, Guanghui Wang et al.[34] applied gene co-expression modules analysis on NSCLC metastases.

Weighted gene co-expression network analysis (WGCNA[35]) is a bioinformatics and systems biology tool that is employed to construct and analyze co-expression networks. This tool is an R programming language package and contains different functions for network construction, visualization, data simulation and gene selection and can be applied for detecting modules (clusters) of highly correlated genes[36]. In the present study, WGCNA was utilized to reconstruct and analyze the gene co-expression network for NSCLC transcriptomic dataset. Xuting Xu and colleagues [37] applied WGCNA to identify hub genes as biomarkers in lung cancer and introduced CCNB1, CCNE2, MCM7, and PCNA as hub biomarker genes. In a similar study, Binglin Chen et al.[38] applied WGCNA on the NSCLC transcriptomics dataset and identified four hub genes (AURKB, CDC20, TPX2 , and KIF2C) as NSCLC prognostic biomarkers based on co-expression network analysis. Moreover, WGCNA was utilized to discover prognostic markers in lung cancer by Bo Ling colleague [39].

The current study aimed to discover potent candidate drugs for NSCLC treatment by analyzing gene co-expression, TF–TG, and drug–gene interaction networks. To this end, at first, a gene co-expression network was reconstructed based on the WGCNA for the NSCLC transcriptome dataset. Then, two significant gene modules, named purple and magenta, were discovered from the reconstructed gene co-expression network. Next, a list of transcription factor (TF) genes regulating purple and magenta modules' genes was gathered from the TRRUST V2.0[40] online database. Afterward, a TF-gene interaction network was reconstructed for the gathered TFs and their target genes. This network is named the TF–TG network. Simultaneously, Gene Ontology (GO) and pathway enrichment analysis were done using the David bioinformatics Resources 6.8[41] for purple and magenta modules' genes, and the results were reported. Subsequently, in order to identify the existing drugs targeting TF–TG network genes, the DGIdb V4.0[42] online database was utilized. After obtaining a list of drugs targeting the TF–TG genes, we reconstructed two drug–gene interaction networks, including drug-TF and drug-TG. Consequently, for each of the drug-TF and Drug-TG networks, nine high-degree drugs (hub drugs) were selected and reported as potent candidate drugs. METHOTREXATE is a hub node in the drug-TG interaction network and regulates 6 genes of the purple and magenta modules. The highest degree drug node in the drug-TF interaction network is CISPLATIN, which regulates 11 TF genes.

In summary, the current study consists of the following three main steps: (1) Gene co-expression network reconstruction, (2) TF_TG interaction network analysis, and (c) Drug- Target interaction network analysis. Compared to the other studies, steps (2) and (3) are novel in our project and have not been applied for drug repurposing for NSCLC before. Moreover, the current study analyses interactions between drugs and both of TFs and non-TF genes (Drug-TF and Drug-TG), which have not been studied before for NSCLC treatment. Figure 1 shows the workflow diagram of the proposed approach.

Figure 1

The workflow diagram of the proposed method. This study applies a gene co-expression network and a drug–gene regulatory network analysis to reposition candidate drugs for NSCLC treatment. (a,b) At first, a transcriptome profile for normal and NSCLC samples was downloaded from the GEO database with the accession number GSE21933. (c,d) Then, a gene co-expression network was reconstructed for the differentially expressed genes (p_value < 0.01) of normal and NSCLC groups using the WGCNA package in the R programming environment, and two significant gene modules (purple and magenta) were extracted from the NSCLC co-expression network. (e) Next, a list of transcription factor genes, which regulate purple and magenta modules' genes, were obtained from the Trrust V2.0 [40] online database. (f,g) Subsequently, two drug–gene interaction networks, named drug-TG (target gene) and drug-TF (transcription factor gene), were drawn using the DGIdb V4.0[42] online database. (e) Finally, 18 candidate drugs are proposed for NSCLC treatment.

Result

Module analysis

For 4218 differentially expressed genes between normal and NSCLC, the gene co-expression network was reconstructed for NSCLC transcriptomics data using the WGCNA. Accordingly, 21 gene modules were discovered from this network (Supplementary Fig. C). The darkorange module is the smallest module with 47 genes, while the largest module is blue with 326 genes. The grey module shows genes that are not assigned to any other detected modules. This module is not considered for further analysis.

Comparing the modules between NSCLC and normal groups

Those modules that have changed significantly between NSCLC and normal groups could deregulate some biological processes and cause disease. Therefore, no-preserve modules between NSCLC and normal groups may cause the NSCLC. As described in the method section, the modules with are considered as no preservation modules for additional analysis. In this regard, the purple and magenta modules have and , respectively and are considered as no preservation modules between NSCLC and normal groups (see Fig. 2). These modules can represent cancer progression from normal to NSCLC stage. Table 1 shows all extracted modules along with their .

Figure 2

Magenta (a) and Purple (b) modules. The circle nodes represent genes (this figure was drawn in the Cytoscape[43] v.3.8.2 software).

Table 1

The of NSCLC co-expression modules compared to the normal gene expression data.

Module name	Size	Z_summary
Purple	167	0.93
Magenta	183	1.3
Orange	52	2.2
Darkgreen	81	2.7
Red	222	3.5
Grey60	101	4
Midnightblue	112	4.3
Greenyellow	160	4.6
Lightgreen	92	4.6
Cyan	133	5.5
Darkred	84	5.5
Lightcyan	108	5.8
Lightyellow	90	5.8
Darkturquoise	70	6.7
Brown	317	6.9
Blue	326	7.4
Royalblue	87	8.1
Darkorange	47	10
Salmon	140	14
Gold		17
Black	192	27
Grey	42	0.46

Enrichment analysis of the gene modules

In order to study the biological functions of the genes in purple and magenta modules, functional enrichment analysis was performed using the DAVID[5] (Database for Annotation Visualization and Integrated Discovery) database. Gene Ontology (GO) enrichment analysis shows that the genes of purple and magenta modules are enriched in 55 and 72 significant (p_value < 0.05) terms, respectively. The results show that the genes in the purple module are significantly enriched in some biological processes related to respiration and lung including: lung epithelial cell differentiation (p_value < 0.001), lung cell differentiation(p_value < 0.001), lung epithelium development (p_value < 0.004), respiratory system development (p_value < 0.006), and lung development (p_value < 0.01). As well as, the results show that the genes in the magenta module are not significantly enriched in biological processes related to respiration and lung. Therefore, the purple module genes are closer to NSCLC than the magenta module. More details for the GO results are reported in supplementary file S2.

Moreover, to investigate biological pathways related to the purple and magenta modules, the pathway enrichment analysis was done based on the REACTOME[44] database. The results revealed that the purple module is significantly enriched in the regulation of the insulin secretion pathway. In addition, the magenta module is significantly enriched in five biological pathways, including Gap junction assembly, TP53 Regulates Metabolic Genes, Tandem of pore domain in a weak inwardly rectifying K + channels (TWIK), Tight junction interactions, and Synthesis of 12-eicosatetraenoic acid derivatives (see supplementary file S2).

TF–TG regulatory network

In order to identify a list of transcription factor (TF) genes that regulate magenta and purple modules' genes, the TF–TG regulatory network was reconstructed. Regulatory information of TFs and TGs was retrieved from the TRRUST[40] online database. After reconstructing the TF–TG regulatory network for magenta and purple modules, we obtained a network with 178 nodes and 182 regulatory interactions. This network contains 107 TFs and 71 TGs. In this network, MUC1 with 11 input degrees and SP1 with 16 output degrees are high TG and TF nodes, respectively. Figure 3 shows the TF–TG regulatory network. A list of TF–TG regulatory interactions is reported in supplementary file S3.

Figure 3

The TF–TG interaction network. This network contains 178 nodes and 182 regulatory interactions. Out of 178 nodes, 107 and 71 nodes are TFs and TGs, respectively. All of the TG nodes are from the magenta and purple modules. The red circles and green triangles represent TGs and TFs, respectively (this figure was drawn in the Cytoscape[43] v.3.8.2 software).

Drug-TG and Drug-TF Interaction networks

The Drug Gene Interaction Database (DGIdb[42]) was used to detect potential drugs for NSCLC treatment. This database is comprehensive and contains drug–gene interaction information. Using DGIdb, we found 277 candidate drugs that target purple and magenta modules' genes. These drugs could have a regulatory effect on NSCLC progression. The drug–gene interaction network was reconstructed based on the obtained drugs and the purple and magenta modules' genes. The Cytoscape[43] v.3.8.2 software was used to reconstruct and visualize this network. This network is shown in Fig. 4, and further details are reported in Supplementary file S4. This network shows nine drugs, including Methotrexate, Olanzapine, Haloperidol, Fluorouracil, Nifedipine, Paclitaxel Verapamil, Dexamethasone, and Docetaxel, are high-degree nodes. These nine drugs, along with target genes, are selected from the network, and then a sub-network is drawn for these drugs and genes (see Fig. 5). In this sub-network, expression levels of genes in NSCLS samples compared to normal samples are shown with blue (Dow-Regulation) to red (Up-Regulation) colors. Among these genes, UGT1A9 has the highest up-regulation, and ATP1A2 has the highest down-regulation expression level in NSCLC group compared to the normal group. METHOTREXATE is a hub node in this sub-network and regulates 6 genes of the purple and magenta modules. High-degree drugs in the network regulate more genes and can have important regulatory effects. The details of target genes' expression level in 9 drugs of NSCLC group compared to the normal group are reported in Supplementary file S6.

Figure 4

The drug–gene interaction network. Totally, 277 candidate drugs were identified as regulators of the purple and magenta modules of the NSCLC network. The red circle shapes and blue hexagon shapes represent genes and drugs, respectively (this figure was drawn in the Cytoscape[43] v.3.8.2 software).

Figure 5

The expression level of hub drugs' target genes in the NSCLC group compared to the normal group. The circle and hexagon shapes represent genes and drugs, respectively. The size of a node indicates its degree (this figure was drawn in the Cytoscape[43] v.3.8.2 software).

Moreover, a list of drugs that target those TFs regulating magenta and purple genes was retrieved from the DGIdb database. A list of TFs with regulatory relationships with magenta and purple modules' genes is available in supplementary files S3. Supplementary Fig. D shows the Drug-TF interaction network, and the details of this network are reported in Supplementary file S5. This network contains 723 nodes, including 675 drugs and 48 TFs. The highest degree drug node is Cisplatin which regulates 11 TF genes, including NFE2L2, TP53, ESR1, BRCA1, ATM, MYC, E2F1, SMAD4, MYCN, TP73, and STAT1. Daunorubicin is the second highest degree drug node that regulates ten TF genes. Among all drugs, those with degree 7 or above along with target TFs were selected from the network, and then a sub-network was drawn for these drugs and TFs. Figure 6 shows this drug-TF sub-network. In this sub-network, the expression level of TF genes in NSCLC samples compared to normal samples is demonstrated with blue (Dow-Regulation) to red (Up-Regulation) colors. Furthermore, the TFs' expression level in 9 drugs of the NSCLC group compared to the normal group is reported in Supplementary file S6.

Figure 6

The expression level of hub drugs' target TFs in NSCLC group compared to the normal group. The triangle and hexagon shapes represent TF genes and drugs, respectively. The size of a node indicates its degree (this figure was drawn in the Cytoscape[43] v.3.8.2 software).

In order to investigate and confirm interactions of the candidate drugs and candidate target genes, the DrugBank[45] database was used. Information for some candidate drugs and candidate target genes is obtained from this database and reported in Table 2. For some other drugs there were no interaction information.

Table 2

Confirmation of the candidate drugs and candidate target genes thanks to the DrugBank database.

Drug name	Type	Target gene
Methotrexate	Transporter	Folate receptor alpha (FOLR1)
Methotrexate	Transporter	Solute carrier organic anion transporter family member 1B3 (SLCO1B3)
Olanzapine	Target	5-Hydroxytryptamine receptor 3A (HTR3A)
Paclitaxel	Transporter	Solute carrier organic anion transporter family member 1B3 (SLCO1B3)
Docetaxel	Transporter	Solute carrier organic anion transporter family member 1B3 (SLCO1B3)
Vorinostat	Target	Histone deacetylase 1 (HDAC1)
Vorinostat	Target	Histone deacetylase 1 (HDAC2)

The literature review of the recent articles shows that most of the proposed candidate drugs have significant effects on NSCLC. Methotrexate and Curcumin are introduced as novel therapeutic strategies to treat NSCLC[46]. The Methotrexate component of MTX-Gd is reported as a chemotherapeutic drug in cancer therapies[47]. Li-Qing Du and colleague noticed that this drug inhibit the expression of RAD51 in cancer cells[48]. Daye Zhang et al. reported that Lenvatinib and Dexamethasone inhibit the invasion and migration of NSCLC[49]. According to Haiyan Ge and colleagues, pemetrexed-induced senescence alleviates in NSCLC by Dexamethasone[50]. In another study, Tatjana Sarcev and colleagues concluded that Dexamethasone significantly decreases weight and appetite in lung cancer patients[51]. Furthermore, Juan P Cata et al. demonstrated that intraoperative Dexamethasone administration to NSCLC patients is not related to its impact on recurrence-free survival (RFS) and overall survival (OS)[52].

Xin Wang and colleagues revealed that the combination of Ondansetron and Olanzapine has better efficacy in preventing vomiting and chemotherapy-induced nausea in NSCLC patients[53]. According to Thierry André and colleagues, combining Oxaliplatin, Fluorouracil, and Leucovorin could improve colon cancer treatment[54]. In a similar study, Herbert Hurwitz et al. reported that Bevacizumab and Fluorouracil composition significantly improved the survival among patients with metastatic colorectal cancer[55]. Furthermore, the combination of Fluorouracil and Curcumin was studied in cancer treatment by Yumeng Wei and colleagues[56]. Barbora Chovancova and colleagues reported that calcium channel blocker Nifedipine inhibits immune escape and colorectal cancer progression[57]. Moreover, in several studies, it has been proved that Nifedipine can promote breast cancer[58,59].

According to Alan Sandler and colleagues, the combination of Paclitaxel, Bevacizumab, and Carboplatin has a significant survival benefit with the risk of increased treatment-related deaths for NSCLC patients[60]. Moreover, Atsuto Mouri and colleagues reported that the combination of Carboplatin and Paclitaxel could be effective and feasible in patients with SCLC, especially those with interstitial lung disease[61]. In another study, Dongjie Ma et al. showed that Paclitaxel increases the sensitivity of lung cancer cells to lobaplatin[62].

Chundi Zhang and colleagues reported that Verapamil might change the expression level of NW23 and EGFR in lung cancer by post-transcriptional and transcriptional levels, respectively[63]. In addition, S Merry et al. studied the role of Verapamil in overcoming cytotoxic drug resistance in human lung cancer[64]. Zhiyuan Shen and colleagues expressed that circular RNA Foxo3 reduction promotes chemoresistance and prostates cancer progression to Docetaxel[65]. In another study, Hai-Hong Zhou and colleagues recounted that combining Docetaxel and erastin may offer an effective administration for chemo-resistant ovarian cancer patients[66]. Furthermore, Marta Prieto-Vila et al. reported that Quercetin and Docetaxel combination could be a promising therapeutic approach in breast cancer treatment[67]. Juan Valle and colleagues introduced Cisplatin plus Gemcitabine as an effective option for treating advanced biliary cancer[68]. In another study, Deborah K Armstrong and colleagues noted that Cisplatin and Paclitaxel combination improved survival in patients with ovarian cancer[69]. Moreover, Kazumasa Noda et al. reported that Cisplatin plus Irinotecan could effectively treat small-cell lung cancer[70]. According to Ana Catarina Alves and colleagues, Daunorubicin coactions with membranes of cancer cells[71]. Furthermore, Yuanyuan Wang and colleagues found that Daunorubicin can be an effective strategy in NSCLC[71] treatment. In a study by Jia Guo and colleagues, the Daunorubicin and Tamoxifen combination was reported as an option to eliminate both cancer stem cells and breast cancer cells[72]. Lilia Antonova and colleagues reported that the expression of the breast cancer susceptibility gene BRCA1 was down-regulated by stress hormone Hydrocortisone in mouse cell line[73]. Yuan Hong and colleagues reported that Doxorubicin and Curcumin combination could be a method for Lung cancer therapy[74]. In a similar study, Abolfazl Akbarzadeh et al. reported the combination of Doxorubicin β-elemene co-loaded as a way to treat lung cancer[75]. Moreover, Vanesa Gregorc and colleagues showed that the NGR-hTNF plus Doxorubicin could be a way for SCLC[76] treatment. According to Yang Yang and colleagues' report, Trichostatin and Azacitidinecan amalgamation decreased tumorigenic of lung cancer cells[77]. Taofeek K Owonikoko and colleagues found that Vorinostat increased Carboplatin and Paclitaxel activity in NSCLC cells[78]. In Sang Eun Park and colleagues' study, Vorinostat and EGFR‑TKI combination was evaluated in NSCLC to reverse EGFR‑TKI resistance[79]. Furthermore, Chun-Hao Pan and colleagues reported that Vorinostat increased the cisplatin-mediated anticancer effects in SCLC cells[80]. Moreover, Doxorubicin Hydrochloride and Haloperidol were tested on different cancer treatments in humans and other organisms[81-87].

Gene set enrichment analysis and candidate drugs validation

In order to validate the proposed drugs for NSCLC treatment, the GSEA was performed based on the Enrichr[88] database. We considered high-degree drug nodes from the drug-TG sub-network, including Methotrexate, Olanzapine, Haloperidol, Fluorouracil, Nifedipin, Paclitaxel, Verapamil, Dexamethasone, and Docetaxel. In addition, high-degree drug nodes from the drug-TF sub-networks, including Cisplatin, Daunorubicin, Dexamethasone, Methotrexate, Hydrocortisone, Doxorubicin, Azacitidine, Vorinostat, and Doxorubicin Hydrochloride, were considered. Out of these 18 drugs, 2 drugs, including Dexamethasone and Methotrexate, are common between drug-TG and drug-TF sub-networks. Therefore, 16 drugs are assumed as potent candidate drugs for NSCLC treatment, and the CMAP analysis was performed for these drugs.

The results show that Methotrexate and Paclitaxel downregulate GLS2 and NFE2, respectively, and Haloperidol and Dexamethasone up-regulate HTR3A and GLS, respectively. Moreover, the Azacitidine up-regulates the DNMT 3A TF gene. Moreover, there was no information regarding other drugs. Romero-Benitez and colleagues studied the impact of paclitaxel on NFE2 in vivo, and they revealed that the expression of NFE2 was up-regulated on day 3[89]. In another study, Anna Schuhmacher et al. [90] assessed functional and coding variants of the HTR3A subunits in response to haloperidol. Moreover, Takuma Kusabe and colleagues[91] reported that the expression level of GLS is reduced by treatment dexamethasone.

Method

Dataset and preprocessing

In this study, the transcriptomics data with accession number GSE21933 was downloaded from the Gene Expression Omnibus (GEO) database. This data contains 42 male samples, including 21 normal and 21 primary non-small cell lung cancer (NSCLC) samples. The mean and standard deviation of the age for all samples in healthy and NSCLC are about 70 and 7.8, respectively. The annotation file with accession number GPL6254 was used to assign probes to gene IDs.

Hierarchical clustering was done for the normal and NSCLC samples independently to check outlier samples. Results show no outlier samples among the normal and NSCLC samples (see supplementary Fig. A). Therefore, all 42 samples, including normal and NSCLC, are considered for further analysis.

Gene co-expression network and gene modules

First of all, differentially expressed genes (GEGs) were calculated between the normal and NSCLC groups applying the adjusted p-value and Benjamini & Hochberg's method based on the GEO2R tool. Overall, 4218 genes with adjusted p_value less than 0.01 were considered as the initial gene list (see supplementary file S1). This gene list was used in gene co-expression network reconstruction.

Then, the gene co-expression network from NSCLC expression data was reconstructed through the Weighted Gene Co-expression Network Analysis (WGCNA[35]) package. This package can reconstruct the gene co-expression network in three different ways: "signed", "unsigned", and "signed hybrid". In this project, the type of gene co-expression network is signed hybrid. To adjust the scale-free property of the network, the β (soft thresholding power beta) parameter is applied in this package.

The soft threshold power beta is determined according to the standard scale-free network[92]. This parameter was set to 7 in NSCLC network (see supplementary Fig. B) to gain the scale independency of the network, where the scale-free index was 0.9. To extract modules for the gene co-expression network, the hierarchical clustering algorithm was applied in WGCNA (see supplementary Fig. C).

Module preservation analysis

To analyze module preservation, the score was used. The modules with is considered as no preservation[93]. The calculation of the is shown in Eq. (1) [94]. In this equation, and are the connectivity and density of the subnetwork, respectively[95]. The score in NSCLC network compared to normal data expression is calculated for all extracted modules (see Table 1). It should be noted that the grey module shows the genes which are not assigned to other detected modules. Two modules including purple and magenta (see Fig. 2) have and considered as no preservation modules. The gene list of these modules is reported in supplementary file S1. These genes could have crucial rules in NSCLC.

Enrichment analysis

To identify the biological mechanisms of the genes in purple and magenta modules, we used functional enrichment analysis based on the DAVID[96,97] (The Database for Annotation, Visualization, and Integrated Discovery) database. Moreover, a pathway enrichment analysis was done for these modules' genes using the Reactome[98] pathway database.

TF–TG regulatory relationships

In order to obtain regulatory information of Transcription factor(TF) genes and target genes(TG), the TRRUST[40] V2.0 online database was utilized. TRRUST is a manually curated database containing transcriptional regulatory information for mice and humans. This version of TRRUST contains 8444 TF–TG regulatory information of 800 human TFs. After obtaining TF–TG regulatory relationships, a TF–TG network, which contained TFs regulating magenta and purple modules' genes, was reconstructed.

Drug–gene interaction network

To identify the candidate drugs that target purple and magenta genes, the DGIdb[42] (Drug Gene Interaction Database) was used. This database is connected to 22 other related databases. This database brings back the target genes based on 24 related databases. In the current project, to identify drug–gene interactions information, only experimentally validated interactions were considered.

Gene set enrichment analysis (GSEA)

The gene set enrichment analysis (GSEA) was performed as a validation method to test whether the proposed candidate drugs can counteract the gene expression perturbations caused by NSCLC. To this end, the Connectivity Map (CMAP[99]) analysis was performed using the Enrichr[88] database. To perform the CMAP analysis, the genes of purple and magenta modules were submitted to the Enrichr[88] database to retrieve up-regulated or down-regulated genes in the cells treated with different drugs. Two datasets of CMAP-up and CMAP-down, which contained the genes up-regulated or down-regulated by different drugs, were extracted. We quested for our proposed candidate drugs in CMAP-up and CMAP-down datasets. Totally, 16 drugs were proposed as potent candidate drugs for NSCLC treatment, which was evaluated using the Enrichr database.

Guangda Li and colleagues[100] identified some hub genes by combining WGCNA, DEG analysis, and functional enrichment analysis in NSCLC. Moreover, in vitro experiments along with the CMAP database were applied to predict and verify small molecule drugs in NSCLC. These researchers reported cephaeline and Emetine with the potential to overcome resistance using CMAP database. In another study, Ying Zheng et al. [101] applied the CMAP database to predict the anesthetic drugs that regulate the differential expression of RNA binding proteins in cervical squamous cell carcinoma.

They reported 65 differentially expressed RNA binding proteins in cervical squamous cell carcinoma. Moreover, they obtained four anesthetics containing procaine, tetracaine, benzocaine, and pentoxyverine. Mengnan Zhao and colleagues[102] have done a study in order to identify a prognostic ferroptosis and iron-metabolism signature for esophageal squamous cell carcinoma and they identified 20 potential compounds using CMAP database.

Moreover, Hang Yang et al.[103] have done multi-omics-based research and they used CMAP for chemotherapy drug analysis and screening for drugs which reduce the expression of high-risk genes. In other study, Zetian Gong and colleague[104], explored several potential small molecule drugs using CMAP based on the mRNAs co-expressed with autophagy-related lncRNAs.

Discussion

This study applied a gene co-expression network analysis to identify potent candidate drugs for NSCLC treatment. To this end, at first, transcriptomics profiles of normal and NSCLC samples were collected, and 4218 genes with a significantly different expression between normal and NSCLC samples were selected for future analysis. Then, a gene co-expression network analysis was reconstructed based on the WGCNA package. Then, two significant gene modules named purple and magenta were identified. Next, a list of transcription factor genes regulating these two modules' genes was gathered from the TRUST V2.0 online database, and a TF–TG regulatory network was drawn. Subsequently, a list of existing drugs that target TF–TG network genes was collected from the DGIdb V4.0 database, and then two drug–gene interaction networks, including drug-TF and drug-TG, were drawn. In data collection, 675 and 278 drugs were identified for the drug-TF and drug-TG networks, respectively. Consequently, nine high-degree drugs from the drug-TF and drug-TG networks were selected separately and introduced as potent candidate drugs for NSCLC treatment. Eventually, 16 drugs were introduced as potent candidate drugs to treat NSCLC. Out of 16 selected drugs, nine drugs (Methotrexate, Olanzapine, Haloperidol, Fluorouracil, Nifedipine, Paclitaxel, Verapamil, Dexamethasone, and Docetaxel) were selected from the drug-TG network, and nine drugs (Cisplatin, Daunorubicin, Dexamethasone, Methotrexate, Hydrocortisone, Doxorubicin, Azacitidine, Vorinostat and Doxorubicin Hydrochloride) were selected from the drug-TF sub-network. Out of these 18 hub drugs, Methotrexate and Dexamethasone are common in drug-TF and drug-TG networks.

In order to evaluate the gene ontology and biological pathways for purple and magenta modules' genes, the DAVID online tool was used. Magenta and purple modules were enriched in 72 and 55 Go terms with a p_value < 0.05, respectively. The results showed that purple and magenta modules were more significantly enriched in phospholipid translocation biological process with a p_value 0.0007 and skin development biological process with a p_value < 1.015e−9, respectively. Moreover, five significant biological process terms of purple module are related to lung and respiratory. These five significant terms are: lung epithelial cell differentiation (p_value < 0.001), lung cell differentiation(p_value < 0.001), lung epithelium development (p_value < 0.004), respiratory system development (p_value < 0.006), and lung development (p_value < 0.01). Whereas, none of the significant biological process terms of magenta module are related to lung and respiratory. In conclusion, the purple module genes can be important compared to the magenta module in NSCLC studies.

In addition, a pathway enrichment analysis was done for these two modules based on the REACTOME database. The results show that the purple module was significantly enriched in the "regulation of the insulin secretion" pathway. Three genes of the purple module, including CACNA1C, RAPGEF3, and GNAI2, are involved in the regulation of the insulin secretion pathway. Talip Zengin et al.[105] introduced the RAPGEF3 for prognostic risk prediction according to overall survival time for lung adenocarcinoma patients. Xiao Wang and colleagues[106] have done genome sequencing analysis for lung adenocarcinoma and introduced CACNA1C as a cancer-related gene. Moreover, they reported that this gene was mutated in lung adenocarcinoma tumor tissue. Furthermore, the magenta module was significantly enriched in five biological pathways, including: "Gap junction assembly", "TP53 Regulates Metabolic Genes","Tandem of pore domain in a weak inwardly rectifying K + channels (TWIK)", "Tight junction interactions", and "Synthesis of 12-eicosatetraenoic acid derivatives".

The "Gap junction assembly" pathway involves four magenta module genes (GJB2, GJB4, GJB5, GJB6). Deng Yun Li et al.[107] and Seon-Sook Han et al. [108] reported that GJB2 expression is aberrantly higher in Lung adenocarcinoma than in control tissue. In a study that Yi-Pei Lin and colleagues[109] have done, GJB4 was reported as a novel biomarker for lung cancer. "TP53 Regulates Metabolic Genes" pathway involves five genes of the magenta module containing GPX2, SESN3, GLS2, SFN, and TP63. In their research, Kui Liu et al.[110] revealed that up-regulation of GPx2 is correlated with worse overall survival for NSCLC patients. Besides, Shuhao Li and colleagues[111] reported that SESN3 has high expression in lung cancer patients compared to healthy patients.

Moreover, this gene was reported as an oncogene in esophageal squamous cell carcinoma cells[112]. Rakibul Islam et al.[107,113] have done a survival analysis, and their results show a worse overall survival value for SFN, and Outcomes show that SFN may play a crucial role in the development of NSCLC. "Tandem of pore domain in a weak inwardly rectifying K+ channels (TWIK)" pathway involves two genes of the magenta module, including KCNK7 and KCNK1. Wen Wang and colleagues[114] constructed a ceRNA network, and they concluded that KCNK1 is specific to LINC00467 in Lung adenocarcinoma. The "Tight junction interactions" pathway involves three genes of the magenta module containing PRKCI, CLDN20, and PARD6G. Yongfeng Wu et al. [115] demonstrated that mutation of PRKCI and some other genes are identified to be correlated with NSCLC metastasis.

Similarly, Fei Yuan and colleagues[116] reported that PARD6G is differentially expressed between Lung adenocarcinoma and lung squamous cell cancer. Finally, the last significant pathway for the magenta module is "Synthesis of 12-eicosatetraenoic acid derivatives". This pathway contains two genes of the magenta module containing GPX2 and ALOX12B. Szymon Zmorzyński et al.[117] showed that the changes in the activity of the GPX2 isoform might be associated with other cancers development. In another study, Chao Ma et al.[118] reported that ALOX12B could predict lung adenocarcinoma accurately.

Conclusion

In conclusion, we used a gene co-expression network analysis to identify potent candidate drugs for the NSCLC treatment in this study. To this end, at first, a gene co-expression network was reconstructed for the transcriptomics data of the NSCLC patients. Then, two significant gene modules, namely magenta and purple, were discovered from the constructed co-expression network. After that, a TF-TG regulatory network was drawn for magenta and purple modules' genes and the TFs targeting these modules' genes. Next, two drug–gene interaction networks, namely drug-TG and drug-TF, were constructed. Subsequently, from each drug-TG and drug-TF network, nine high-degree drugs were selected and reported as potent candidates for NSCLC treatment. Consequently, 16 drugs, including Methotrexate, Olanzapine, Haloperidol, Fluorouracil, Nifedipine, Paclitaxel, Verapamil, Dexamethasone, Docetaxel, Cisplatin, Daunorubicin, Hydrocortisone, Doxorubicin, Azacitidine, Vorinostat, and Doxorubicin Hydrochloride, were introduced as potent candidate drugs to treat NSCLC. Moreover, gene ontology and pathway enrichment analyses were run for the magenta and purple modules.

Supplementary Information 1.

Supplementary Information 2.

Supplementary Information 3.

Supplementary Information 4.

Supplementary Information 5.

Supplementary Information 6.

Supplementary Figures.