Transcriptome profiling of cervical cancer cells acquired resistance to cisplatin by deep sequencing

Abstract Cervical cancer is one of the most fatal malignancies in females. Acquired resistance to chemotherapeutic agent is one reason behind this lethality. In this study, we developed cisplatin resistance cell line, subsequently examined the molecular mechanisms linked. Transcriptome sequencing technology was utilized to compare the various expression models between the cisplatin-resistant cell line (Hela/DDP) and its parental cell line human cervical adenocarcinoma Hela. The present study has identified 2,312 differentially expressed genes (DEGs). Results showed there were 1,437 up-regulated genes and 875 down-regulated ones. Databases analysis including Gene ontology (GO), Cluster of Orthologous Groups of proteins (COG) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway were performed to reveal potential molecular mechanisms. We studied AKT3, a crucial gene in the PI3K/AKT pathway which clustered the most DEGs. Silencing AKT3 in Hela/DDP could enhance its sensibility to cisplatin. Quantitative real-time reverse transcription PCR (qRT-PCR) and western blot experiments were showed that expression of AKT3 was decreased after siRNA interference and inhibitor treatment. CCK-8 experiments showed that low expression of Akt3/pAkt enhanced the sensitivity of drug-resistant cells to cisplatin. Apoptotic analysis demonstrated that inhibition of AKT3 increased the rate of Hela/DDP apoptosis. Our results suggest a novel mechanism by which upregulated expression of AKT3 in cervical cancer may lead to resistance to cisplatin.


Introduction
The incidence of cervical cancer ranks third among female malignant tumours, and its mortality rate ranks second. Therefore, it seriously harms women's health. For 2017, the estimated new cancer cases are 12,820, and deaths of uterine are 4,210 in America. Cervical cancer mortality increased by nearly 3% annually between 2010 and 2014 [1]. Cisplatin is a first-line chemotherapy for the treatment of various solid tumours. This agent remains as the cornerstone in the systemic antineoplastic management of the malignancy. Combing cisplatin with paclitaxel is one of most commonly utilized regimens in the recurrent or metastatic disease setting [2] although cisplatin acts with tumours to induce apoptosis. However, in actual clinical treatment, tumour cells will develop resistance. Studies have shown that patients with recurrent or metastatic swelling have a relatively poor prognosis, and their 1-year survival rate is only 15-20% [3]. Resistance mechanism of cisplatin in tumours are multifactorial according to research studies, including drug metabolism: abnormal expression of transporters, cell detoxification mechanism: enhanced intracellular detoxification, gene level: increased the capacity of DNA repair and increase the ability to inhibit apoptosis. The above are the most significant recognized mechanisms [4]. However, the molecular mechanism of how cisplatin resistance develops involved in cervical cancer need further study. Therefore, it is tremendously significant to clarify the potential molecular mechanism underlying resistance of cisplatin in cervical cancer cells.
PI3K/AKT signalling pathway involves several biological processes with cellular survival and apoptosis, tumour development as well as progression [5,6]. Activation of PI3K/AKT signalling pathway has been found in many various kinds of cancer. Moreover, it confers resistance to various cancer therapies [7]. Akt, a kind of serine/threonine kinase which acts downstream of PI3K, which is well-known as a critical factor in PI3K/AKT signalling [8], also recognized as another name: protein kinase B (PKB). It comprises three highly homologous isoforms: AKT1 (PKBa), AKT2 (PKBb) and AKT3 (PKBc). It acts pivotal parts in cellular survival and apoptosis [9]. Acquired anti-tumour drug resistance cells of tumour tend to increase Akt activation. Past studies have shown that AKT overexpression regulates drug resistance in cancer cells during chemotherapy [10]. AKT activity has been recognized as a factor associated with drug resistance closely in other type of cancer. For example, in breast cancer cells there was an activated AKT1 expression in drug resistance cells [11]. AKT3 is one of the most pivotal members of AKT family; it is contributed to the process of proliferation, tumorigenesis, cellular differentiation, and apoptosis [12][13][14][15][16][17]. AKT3 is also associated with chemoresistance of cancer cells [18,19]. Research study suggested that knockdown of AKT3 increases migration along with sensibility of breast cancer cells to cisplatin [20]. Besides, AKT3 knockdown also induces mitochondrial dysfunction in A549 lung cancer cells [21] However, the relationship between up-regulation of AKT3 expression in cervical cancer cells and its resistance of cisplatin need further research.
Next-generation sequencing (NGS) or high throughput sequencing technologies and parallel genomic sequencing facilities are utilized in a number of parts of cancer research [22]. Previously, NGS has been reported to study the genetic regulation of cancer development, drug resistance and prognosis of cervical cancer research [23]. Since RNA-Seq offers virtually unlimited possibilities in modern bioanalysis [24], we used this technique to study the difference in gene expression levels between Hela/DDP Hela. The aim of our study was to figure out the mechanism underlying cisplatin resistance in cervical cancer. Results showed that AKT3 fromPI3K/ AKT signalling pathway was up-regulated in the cisplatin resistant cell line Hela/DDP. With this found, we take AKT3 as a target gene involves the process of how resistance of cisplatin occurs. Furthermore, our work suggested that both inhibition of AKT3 and knockdown it contributes to Hela/DDP cells increase the sensitivity to cisplatin by inducing cellular apoptosis and reducing prefoliation. Tumour resistance to cisplatin is one of the causes of chemotherapy failure. Exploring the molecular mechanism of drug resistance is an important part of overcoming tumour resistance and will provide a theoretical basis for reversing tumour resistance and improving the efficacy of cancer chemotherapy. Research findings have confirmed that the PI3K/AKT signalling pathway is aberrantly activated in many solid tumours, especially in chemotherapyresistant cells. Therefore, it is necessary to further explore the relationship between PI3K/AKT signalling pathway and chemoresistance from a new perspective.

Cell culture
Human cervical squamous carcinoma cell line Hela was obtained from the American Type Culture Collection (ATCC, Rockville, MD, USA) and cultured in Dulbecco's Modified Eagle Medium (DMEM) (Hyclone) supplemented with 10% FBS (Gibco) and antibiotics (Gibco). Hela cells were cultured in a couple of increasing concentrations of cisplatin (Sigma) for over 6 months to establish cisplatin-resistant cervical cancer cell lines, Hela/DDP. All cells were maintained at 37 C in a humidified 5% CO 2 atmosphere.

Cytotoxicity assay
Assessment of cell viability was carried out using CCK-8 assay. An adherence of Hela and Hela/DDP cells (1 Â 10 5 cells/well) was seeded into 96-well plate followed by incubating for 24 h. The cells were treated with 100 lL of media containing DDP of different concentrations, which were 0, 0.2, 0.39, 0.78, 1.56, 3.13, 6.25, 12.5, 25, 50, 100 lg/mL and further incubated for 48 h. Then, wells were replaced with 100 lL media which containing 0.5 mg/mL CCK-8 and incubated at 37 C for additional 4 h. Absorbance at 450 nm was measured using a microplate reader. The inhibition rate for cell growth was calculated according to the following formula: inhibition rate (%) ¼ (1Àmean absorbance of treatment group/mean absorbance of untreated group) Â 100%.

Transcriptome sequencing
RNA degradation and contamination were checked on 1% agarose gels. RNA purity was tested using the NanoPhotometer spectrophotometer (IMPLEN, CA, USA). RNA concentration was measured using Qubit RNA Assay Kit in Qubit 2.0 Flurometer (Life Technologies, CA, USA). RNA integrity was measured utilizing the RNA Nano 6000 Assay Kit of the Agilent Bioanalyzer 2100 system (Agilent Technologies, CA, USA). For RNA sample preparation, a total quantity of 1 lg RNA per sample was used as input material. Sequencing libraries were produced using NEBNext UltraTM RNA Library Prep Kit for Illumina (NEB, USA) based on the instructions from manufacturer. In brief, mRNA was purified from total RNA using poly-T oligo-attached magnetic beads. Fragmentation was carried out using divalent cations under elevated temperature in NEBNext First Strand Synthesis Reaction Buffer (5X). We used AMPure XP system (Beckman Coulter, Beverly, USA) to select cDNA fragments of especially 200-250 bp in length. Then, 3 ll USER Enzyme (NEB, USA) was utilized with size-selected, adaptor-ligated cDNA at 37 C for 15 min followed by 5 min at 95 C before PCR. PCR was performed with Phusion High-Fidelity DNA polymerase, Universal PCR primers and Index (X) Primer. Finally, PCR products were purified using AMPure XP system and library quality was assessed on the Agilent Bioanalyzer 2100 system.
The clustering of the index-coded samples was performed on a cBot Cluster Generation System using TruSeq PE Cluster Kit v4-cBot-HS (Illumia) according to the manufacturer's recommendations. After cluster generation, the library preparations were sequenced on an Illumina Hiseq 2500 platform and paired-end reads were generated.

Differential expression analysis
Differential expression analysis of two conditions/groups was performed using the DESeq R package (1.10.1). p Values were adjusted using the Benjamini and Hochberg's approach for controlling the false discovery rate.
p Values were adjusted using p Value < .005 and |log 2 (fold change) | !1 was set as the threshold for significantly differential expression.

Gene functional annotation
Gene function was annotated based on the following databases: Nr (NCBI non-redundant protein sequences)；Nt (NCBI non-redundant nucleotide sequences)； Pfam (Protein family)； KOG/COG (Clusters of Orthologous Groups of proteins)；Swiss-Prot (A manually annotated and reviewed protein sequence database)；KO (KEGG Ortholog database)；GO (Gene Ontology).

RNA isolation, reverse transcription and qRT-PCR
Total RNA was extracted from cells using Trizol Reagent (Life Technologies). To detect different the relatively mRNA expression of 14 randomly chosen different expressive genes which were found after RNA sequencing, cDNA was synthesized using the Prime-Script RT reagent Kit (TaKaRa, Dalian, China). Quantitative reverse transcription-polymerase chain reaction (qRT-PCR) was carried out by utilizing SYBR Premix ExTaq (TaKaRa) with the ABI ViiA 7 Real Time PCR System (Thermo Lifetech). The relative expression of these 14 different expressive genes was quantified with the 2 ÀDDCT method. The mRNA level of tested genes was normalized to the level of ACTB, and was calculated by the 2 ÀDDCT method. PCR reactions for each gene were performed three times. Primers were used in qPCR reaction are presented in Table 1.

siRNA synthesis and cell transfection assays
Three human siRNAs (GenePharma, Shanghai, China) designed against AKT3 (GenBank NM_181690). There was one negative random siRNA, one negative control FAM siRNA and one GAPDH positive control siRNA (GenePharma, Shanghai, China). The sequences for AKT3-siRNAs and control siRNA are presented in Table 2.
The cells were first divided into seven groups, the mock group with only the transfection reagent, the negative control group added with a non-targeting control siRNA or the transfection reagent, and the other three groups supplemented with different AKT3 siRNA and the transfection reagents. Cells were seeded on 6-well plates (5 Â 10 5 cells for Hela/DDP) and were allowed to attach for 24 h, and then treated with 100 pmol siRNA per well. Using equimolar amounts of siRNAs were and Lipofectamine 2000 Transfection Reagent from Invitrogen (Madison, WI, United States). Transfected cells were grown at 37 C for 6 h, followed by incubation with complete medium. Cells were maintained for 24 h for further experiments

Western blot
Total protein extracted from cell lines were parted by 12% sodium dodecyl sulfate gel electrophoresis (20 mg per lane) and transferred onto a nitrocellulose membrane by electroblotting in transfer buffer for 150 min at 120 V. To prevent non-specific binding of antibodies by blocking with 5% skim milk the membrane was incubated with primary antibodies that specific for Akt3 (Santa Cruz Biotechnology Inc, USA), and GAPDH (Santa Cruz Biotechnology Inc, USA) at 1:1,000 dilution in TBS solution at 4 C overnight. The secondary antibody (CoWin Biosciences, China) at 1:15,000 dilution in skim milk was infrared at room temperature for 2 h. Result bands were visualized using an enhanced chemiluminescent substrate (Thermo Fisher Scientific, Inc., USA).

Cell apoptosis assays
Cells were seeded at a density of 5 Â 10 5 cells/60-mm Petri dish in complete medium for 24 h. Followed by treating with the designated concentrations of cisplatin for 24 h. To assess apoptotic cells, Apoptotic Cells Detection kits were used. Cells were accumulated and resuspended in 500 ll of binding buffer, 5 ll of Annexin V-FITC and 5 ll of PI. Then, they were cultured for 15 min at room temperature. CellQuest software (BD Biosciences) was used to analyze samples by flow cytometer.

Statistical analysis
p Values and false discovery rate (FDR) filtering were used to identify the mRNAs expression was significantly different between Hela and Hela/DDP cell lines. These differentially expressed mRNAs were identified through fold change filtering. Results from qPCR were analyzed using GraphPad Prism 5.

Cisplatin exposure induces drug-resistance in cervical cancer
To study mechanism of cisplatin resistance in cervical cancer, we established in vitro cisplatin-resistance cervical model (Hela/DDP) through stepwise exposure Hela cells to cisplatin. As shown in Figure 1. Hela/DDP displayed significant Table 1. Primer sequences used for qRT-PCR.

Analysis of transcriptome sequencing
After removing reads with adapters, unknown nucleotides, and low-quality reads of sequencing, we acquired a total of 14.54 Gb clean data. As shown in Table 3, The Q30 percentage (accuracy rate of base identification > 99.9% and error rate of sequencing < 1‰) was more than 85.07%. Clean reads were mapped to reference through using TopHat2 [25] which based on Bowtie2 [26]. The average mapping ratio with reference genes was 71.51%. Table 4 lists respective mapping rates of two samples. We normalized the number of Mapped Reads in the samples and the length of the transcript. FPKM (Fragments Per Kilobase of transcript) was an indicator of gene expression. Take Fold Change >2 and FDR < 0.01 as the screening criterion in the process of gene detection, the transcriptome sequencing analysis identified 2,312 differently expressed genes between Hela and Hela/ DDP cells, including 1,437 up-regulated genes and 875 down-regulated ones. The volcano plot which illustrated different expression is shown in Figure 2.

Functional enrichment analysis of DEGs
In Gene ontology analysis of DEGs, all of them were categorized into three major categories: biological processes, cellular components and molecular functions. The present study identified that the 2,137 DEGs were classified into 62 functional categories, 22 in biological processes, 21 in molecular functions and 19 in cellular components including reproduction, cellular processes and metabolic processes (p < .05; Figure 3). COG (Cluster of Orthologous Groups of proteins) database is a classification system based on orthologous relationships Table 2. Sequences of AKT3-siRNAs.

AKT3-siRNA1
Sense5 0 -GGAUGCCUCUACAACCCAUTT-3 0 antisense5 0 -AUGGGUUGUAGAGGCAUCCTT-3 0 AKT3-siRNA12 Sense5 0 -GCACACACUCUAACUGAAATT-3 0 antisense 5 0 -UUUCAGUUAGAGUGUGUGCTT-3 0 AKT3-siRNA13  Clean reads were paired-end reads of clean data.  Figure 2. Volcano graphs of all DEG enrichment. Each point in the volcano diagram represents one differentially expressed gene, and the x-axis represents the log 2 (fold change) and the y-axis shows the -log 10 (1 À P). The greater the absolute value of the horizontal coordinate, the greater the difference between the two samples. The larger the ordinate value, the differential expression is more significant, and the differentially expressed gene is more reliable. In the green dots, the differentially expressed gene is down-regulated, and the red dots represent the differentially expressed genes. The black dots represent the non-differentially expressed genes.   between genes of prokaryotes and eukaryotes. The results of the differential expression genes classified on COG system are shown in Figure 4.

KEGG pathway analysis of DEGs
To investigate high-level functions of DEGs between Hela and Hela/DDP cells after transcriptome sequencing, signalling pathway analysis was performed. Totally, 1,344 DEGs were annotated into 277 pathways in the KEGG (Kyoto Encyclopedia of Genes and Genomes) database. The annotation results of gene KEGG were classified according to the pathway types in KEGG, and the top 50 enriched pathways are shown in Figure 5(A). Except for pathways in cancer, the top 5 signalling pathways annotated with high numbers of DEGs were: PI3K/AKT signalling pathway, Rap1 signalling pathway, MAPK signalling pathway, Focal adhesion and Cytokinecytokine receptor interaction. Pathway significance enrichment analysis used pathway in the KEGG database as a unit and applied hypergeometric test to find pathways of significant enrichment in the differential expression genes compared with the whole genome. The result of the statistics of pathway enrichment is shown in Figure 5(B), showing 20 pathways with the least significant Q value. The pathway which involves the most DEGs annotation is displayed in Figure 5(C).

Verification of significant DEGs
We  Table 5. A evaluation was also performed between results of NGS and qRT-PCR, as shown in Figure 6.
AKT3-siRNA and MK-2206 attenuated AKT3 expression AKT3 has been suggested to be potential target of drug design as it plays a role in tumour migration [27]. To investigate the potential role of AKT3 on drug resistance, and the genes responsible for the significant enrichment in PI3K/AKT signalling pathway, we performed qRT-PCR to assessed its mRNA expression levels in Hela/DDP cells and parental Hela cells. The relative expression of mRNA between transcriptome sequencing and qPCR were highly consistent, as shown in Figure 7(A). In order to study whether AKT3 is related with cisplatin resistance in cervical cancer, we down-regulated AKT3 in Hela/DDP through siRNA or inhibitor. The efficiency of AKT3-siRNA was examined using real-time PCR and western blot. We demonstrated a decrease in AKT3 mRNA in Hela/DDP cells treated with AKT3-siRNAs at 24 h after transfection, respectively (28.38%, 68.22% and 4.17% decrease) (Figure 7(B)). Furthermore, AKT3-siRNA2 was used for high transfection rate.
To test variation in the protein level, we also measured the activation of Akt using western blot. The result showed that Hela/DDP cells had a high level of phosphorylated Akt (p-Akt) while additional administration of AKT3-siRNA, MK-2206 significantly suppressed pAkt of Hela/DDP cells compared with control (Figure 7(C)). The data suggested that PI3K/AKT signalling pathway were activated in Hela/DDP cells.

MK-2206 treatment increased sensitivity of Hela/DDP to cisplatin and AKT3-siRNA treatment induced apoptosis
The CCK-8 assay revealed that Hela/DDP cells transfected with AKT3-siRNA exhibited slightly increased sensitivity to

Discussion
In the present study, our aim was to obtain potential therapeutic strategies of resistance to cisplatin in cervical cancer by comparing transcriptome differences associated with cisplatin resistance cells (Hela/DDP) and Hela cell lines. Firstly, our data from next-generation sequencing pointed out 2,312 genes were significantly differentially expressed between the cisplatin-resistant cell line (Hela/DDP) and its parental cell line (Hela), including 1,437 up-regulated genes and 873 downregulated genes in Hela/DDP compared with Hela. Second, the function of these DEGs was annotated and enriched by reference to the database like GO, COG and KEGG. The top-five DEGs-enriched relative signalling pathways in the KEGG database were investigated in order to find the molecular mechanism of cisplatin resistance. Based on the RNA-seq and statistical results, AKT3 was identified as the target gene which was possible involved in a potential resistance mechanism. Lastly, qRT-PCR and western Blot were performed to certify the AKT3 mRNA and protein expression. In addition, apoptotic rate was measured using performing flow cytometry. Cisplatin is a cell cycle non-specific cytotoxic and a firstline drug for cervical cancer chemotherapy [28]. It plays a role of anticancer via DNA damage response and apoptosis by generating DNA lesions [29]. However, resistance to cisplatin results in limitations during treatment [4]. Resistance of cisplatin has been proved to be associated with cellular contents up-taking, drug efflux, increasing the ability of detoxification, inhibition of apoptosis and increased DNA repair [30].
In this work, we revealed a connection between the PI3K/ AKT pathway and cisplatin resistance in cervical cancer. This specific signalling pathway enriched most DEGs. It shows that the activation state of the PI3K/AKT pathway results in tumours developing resistance to targeted therapies as well as to chemical/radiotherapy [31]. Inhibition of certain components of the PI3K pathway not only synergizes with chemotherapy, radiation therapy, hormone therapy and targeted drugs in cancer therapy, but also overcomes its resistance [32,33]. In our study, there were three key genes in the pathway which were upregulated in cisplatin resistant cells: HSP90B1, AKT3 and CREB5.
Therefore, the connection between AKT3 and cisplatin resistance in cervical cells is further studied in this work. We first reported that AKT3 takes a crucial part in cisplatin resistance in cervical cells. It has been observed that AKT3 is highly expressed in cisplatin resistance Hela/DDP cells than cisplatin sensitive Hela cells both in mRNA and protein level. qPCR analysis also showed that mRNA measured matched their sequencing results. In addition, knockdown of AKT3 increased sensitivity to cisplatin and apoptosis. After treating with MK-2206 2HCl, Hela/DDP cells showed more sensitive to cisplatin as well. Western blot analysis further indicated that MK-2206 2HCl significantly decreased expression of p-AKT. Results suggested that AKT3 takes an important part in many aspects of cancer biology which involves proliferation, survival, apoptosis and drug resistance of cancer cells.
The activity of AKT3 may depend on upstream or downstream molecules that have network with it in particular kind of cancer cell. For instance, it has been reported that AKT3 contributed endocrine therapy resistance of ErbB2 (þ) breast tumour cells by regulating the expression of ErbB2 [34]. In human ovarian cancer, miR-489 has been demonstrated to inhibit cisplatin resistance and cell growth, and promote apoptosis by suppressing AKT3 expression [35]. Additionally, it has been shown that expression of AKT3 is negatively correlated with miR-505, which, in turn, regulates drug sensitivity in MCF7-ADR [36]. Therefore, microRNA-mediated AKT3 expression may play important roles in cisplatin resistance in cervical cancer cells.
In conclusion, our work analyzed transcriptome-related changes in cisplatin-resistant Hela cells and their progenitors. We found three new promising biomarkers in the PI3K/AKT signalling pathway (HSP90B1, AKT3 and CREB5). In addition, inhibition of AKT3 may have effects on the development of new therapeutic clinic strategies to overcome acquired cisplatin resistance in cervical cancer.