Skip to main content

ORIGINAL RESEARCH article

Front. Oncol., 15 February 2021
Sec. Cancer Genetics

Medication for Acromegaly Reduces Expression of MUC16, MACC1 and GRHL2 in Pituitary Neuroendocrine Tumour Tissue

Rihards Saksis&#x;Rihards Saksis1†Ivars Silamikelis&#x;Ivars Silamikelis1†Pola LaksaPola Laksa1Kaspars MegnisKaspars Megnis1Raitis PeculisRaitis Peculis1Ilona MandrikaIlona Mandrika1Olesja RogozaOlesja Rogoza1Ramona PetrovskaRamona Petrovska1Inga Balcere,Inga Balcere2,3Ilze Konrade,Ilze Konrade2,3Liva SteinaLiva Steina4Janis StukensJanis Stukens4Austra BreiksaAustra Breiksa4Jurijs NazarovsJurijs Nazarovs4Jelizaveta SokolovskaJelizaveta Sokolovska5Valdis Pirags,Valdis Pirags4,5Janis KlovinsJanis Klovins1Vita Rovite*Vita Rovite1*
  • 1Latvian Biomedical Research and Study Centre, Riga, Latvia
  • 2Riga East Clinical University Hospital, Riga, Latvia
  • 3Riga Stradins University, Riga, Latvia
  • 4Pauls Stradins Clinical University Hospital, Riga, Latvia
  • 5University of Latvia Faculty of Medicine, Riga, Latvia

Acromegaly is a disease mainly caused by pituitary neuroendocrine tumor (PitNET) overproducing growth hormone. First-line medication for this condition is the use of somatostatin analogs (SSAs), that decrease tumor mass and induce antiproliferative effects on PitNET cells. Dopamine agonists (DAs) can also be used if SSA treatment is not effective. This study aimed to determine differences in transcriptome signatures induced by SSA/DA therapy in PitNET tissue. We selected tumor tissue from twelve patients with somatotropinomas, with half of the patients receiving SSA/DA treatment before surgery and the other half treatment naive. Transcriptome sequencing was then carried out to identify differentially expressed genes (DEGs) and their protein–protein interactions, using pathway analyses. We found 34 upregulated and six downregulated DEGs in patients with SSA/DA treatment. Three tumor development promoting factors MUC16, MACC1, and GRHL2, were significantly downregulated in therapy administered PitNET tissue; this finding was supported by functional studies in GH3 cells. Protein–protein interactions and pathway analyses revealed extracellular matrix involvement in the antiproliferative effects of this type of the drug treatment, with pronounced alterations in collagen regulation. Here, we have demonstrated that somatotropinomas can be distinguished based on their transcriptional profiles following SSA/DA therapy, and SSA/DA treatment does indeed cause changes in gene expression. Treatment with SSA/DA significantly downregulated several factors involved in tumorigenesis, including MUC16, MACC1, and GRHL2. Genes that were upregulated, however, did not have a direct influence on antiproliferative function in the PitNET cells. These findings suggested that SSA/DA treatment acted in a tumor suppressive manner and furthermore, collagen related interactions and pathways were enriched, implicating extracellular matrix involvement in this anti-tumor effect of drug treatment.

Introduction

Pituitary neuroendocrine tumors (PitNETs) are benign intracranial endocrine tumors with potentially high prevalence in the population. Clinically significant PitNET affects approximately 0.1% of the general population, but clinically insignificant or undiagnosed PitNET can be found in approximately 16.7% (1, 2).

Somatotropinoma, which develops from anterior pituitary somatotroph cells, are characterized by increased synthesis and secretion of growth hormone (GH). They constitute 10–15% of all clinically significant somatotropinomas and usually cause acromegaly in adults or gigantism in children with additional comorbidities (3). Acromegaly is a rare and chronic endocrine disorder and is characterized by abnormal growth of the extremities, cardiovascular diseases, and metabolic disorders, such as diabetes mellitus that is caused by increased levels of insulin growth factor 1 (IGF-1), activated by high GH levels (4).

The goal of therapeutic strategies for this condition is to normalize GH and IGF-1 levels, remove tumor mass and/or stabilize tumor growth while maintaining normal pituitary function (5). Somatotropinoma resection is recommended as the primary therapy, but often medical treatment is chosen as the first-line therapy to reduce tumor mass (4, 5). After PitNET resection, often the drug treatment is then used to better control the disease (5). Somatostatin receptors (SSTRs) are the main targets for this therapy, and 90% of somatotropinomas express predominantly SSTR2 and SSTR5 subtypes. If somatostatin analogs (SSAs) fail to control IGF-1 levels, dopamine agonists (DAs) that target dopamine receptor 2 (D2R) are used as complementary management options (6). Although SSA and DA have demonstrated an inhibitory effect on the secretion of both hormones and cell proliferation, approximately one-third of acromegalic patients are resistant to cabergoline and octreotide treatment (6).

The molecular mechanisms involved in the downstream signaling and pharmacological resistance after SSA/DA therapy are still poorly understood (6). Their principle inhibitory effect on hormone secretion is mediated through the inhibition of calcium channels and activation of potassium channels (through SSTR1, SSTR2, SSTR5) and also by inhibition of the adenylyl cyclase-cAMP signaling pathway (through SSTR3, SSTR4) (7, 8). However, stimulation of other second messenger signaling molecules, such as protein tyrosine phosphatases, plays an important role in controlling growth caused by somatostatin/somatostatin receptor ligands that prevent cell proliferation by inhibition of the phosphoinositide 3-kinase (PI3K)/AKT signaling pathway, mainly through SSTR2 or via the inhibition/activation of MAPK signal pathways by SSTR1, SSTR2, and SSTR5 (8).

Studies have found that long-term SSA therapy is associated with lower beta-arrestin expression levels (9, 10), by potentially working as natural desensitization factors for the SSTR, and expression profiles of PitNETs can change due to SSA treatment. It has been demonstrated that several factors that act as tumor suppressors and influence gene expression in PitNET development, such as aryl hydrocarbon receptor-interacting protein (AIP) and tumor suppressor gene PLAG1 like zinc finger 1 (ZAC1), modulate the antiproliferative effects of SSA (1113). SSAs can induce an increase in the expression of AIP and ZAC1 (14, 15), suggesting that SSA treatment affects factors related to PitNET development and clinical parameters, that could affect not only molecular signaling but also the transcriptional profiles of PitNET cells.

In this study, we compared differences in the total transcriptome of PitNET patients with and without SSA/DA treatment before surgery. Our objective was to determine whether the transcriptome signature could distinguish those acromegaly patients that had been on therapy with those that had not, and also to identify candidate genes whose expression is affected by SSA/DA treatment.

Materials and Methods

Sample and Data Collection

All somatotropinoma tissue samples were obtained from patients who underwent planned resections at Pauls Stradins Clinical University Hospital, Latvia. All patients were recruited from the Genome Database of the Latvian population (LGDB), a government-funded national biobank (16), and all biobank activities and research in this article comply with the Declaration of Helsinki. Broad informed consent for LGDB and project-specific consent for research involving the pituitary tumors were obtained from all patients. Both the biobank and PitNET research was approved by the Central Medical Ethics Committee of Latvia (protocol No. 22.03.07/A7 and 01.29.1/28, respectively). Sample collection and the research process both complied to the Declaration of Helsinki. After resection, PitNET tissue samples were stored in RNAlater® Solution (Thermo Fisher Scientific, USA) and stored for up to 24 h at +4°C and then at −20°C for up to 2 months, until DNA/RNA extraction.

Total RNA Extraction From Pituitary Neuroendocrine Tumor Tissue

Twenty to 30 mg of PitNET tissues was homogenized using a FastPrep-24™ homogenizer in RLT Plus buffer in Lysing Matrix D 1.5 ml tubes (MP Biomedicals, USA), and total RNA was extracted with AllPrep DNA/RNA Mini Kit (Qiagen, Germany) and stored at −80°C. Total RNA quality control was carried out using Agilent 2100 bioanalyzer and an Agilent RNA 6000 Pico Kit (Agilent Technologies, USA), and concentration was measured using Qubit 2.0 fluorometer and Qubit™ RNA HS Assay Kit (Thermo Scientific, USA).

Transcriptome Library Preparation and Sequencing

Transcriptome libraries were prepared with the Low Input RiboMinus™ Eukaryote System v2 and Ion Total RNA-Seq Kit v2 (Thermo Fisher Scientific, USA) following the manufacturer’s instructions. The average concentration of total RNA was 88.15 ng/µl (range: 35.2–138 ng/µl). Input amount for the transcriptome library preparation was 500, ng and the average library concentration was 23.55 ng/µl (8.91–31.8 ng/µl) with an average fragment length of 214 bp (189–236 bp). The preparation of transcriptome libraries was carried out in relation to 1:1 based on patient samples who had or had no medical therapy before PitNET resection, to reduce any batch effects. Libraries were sequenced at the Latvian Biomedical Research and Study centre, Genome Centre Core facility using the Ion Proton™ System for Next-Generation Sequencing (Ion OneTouch and Ion PI™ Hi-Q™ OT2 Reagents 200 Kit for emulsion PCR and Ion PI™ Hi-Q™ Sequencing 200 Solutions kit and Ion PI™ Chip V3 chips for sequencing, all from Thermo Fisher Scientific, USA). Libraries were sequenced three to four times to acquire at least 20 M reads per sample, and at least 10 M would represent uniquely mapped reads, so that the rRNA content would not reach 50% of the sample. For sample PA05, 20 M of total reads was not achieved, and uniquely mapped reads were close to 7 M. Nevertheless, the relevant sample was included in the data analysis based on a 2014 study by Liu et al., analyzing the effects of read number and biological repetitions on differentially expressed genes (DEGs). They concluded that by increasing the amount of biological repetitions independent of the number of reads from the sample library, they were able to increase the number of DEGs and logFC resulting in increased accuracy of expression levels and greater efficiency in the resulting analysis (17).

Data Analysis

Transcriptome data has been deposited in an online repository (GEO Submission GSE160195). Differentially expressed genes were obtained using the DESeq2 (v1.26.0) package (18) from the Bioconductor (v3.10) project (19) and R (v3.6.1) software (www.r-project.org). Read counts were first filtered based on read count frequency in all samples. Genes with at least ten reads in at least three samples were included in the analysis. To detect sample outliers, read count data were transformed with variant stabilizing transformations (VSTs), considering the design of the experiment and visualized with sample distance heat mapping, PCA and MDS methods. Read count density visualization was used to check for problematic samples. Counts were replaced with trimmed mean values for genes which were marked as outliers based on their dispersion Cook’s distance values, which were calculated as the.99 percentile of the F (p, m − p) distribution for each gene. To account for batch effects, surrogate variables (SVAs) were calculated using the sva (v3.34.0) package (Bioconductor project) and added to the design matrix. DESeq function call was then used on raw data but filtered by read count frequency counts to detect DEGs. “Independent Hypothesis Weighing” (IHW) was specified as the independent filtering method, which is not the default in DESeq2, to gain statistical power for testing, using the IHW (v1.14.0) package (20), which uses the base mean value of counts for each gene, in this case, as a weight for the adjusted p-value calculation. Next, we performed log fold shrinkage with the apeglm (v1.8.0) package (21), which uses the heavy tailed Cauchy prior distribution, to account for extra variability that comes from genes with low counts and high dispersion values. All parameters were left at their default values, except for “coef” which lets us denote which design coefficient we want to perform shrinkage on. Genes with FDR <0.05 and transformed absolute logFC value >1.5 were deemed as significant and selected. To visualize differential expression results, a heatmap was constructed with the pheatmap (v1.0.12) package using the differences from means across all samples of the VST normalized counts for each gene. Heatmap data was clustered both by genes and samples using the kmeans algorithm. A volcano plot was also plotted with the Enhanced Volcano (v1.4.0) package using the same log fold shrinkage values obtained from apeglm. Kendall’s Tau correlation analysis between read count matrix and Knosp classification index and between read count matrix and Ki-67 proliferation index was performed to test whether these factors affect the expression levels of detected DEG’s and whether they should be included in the model (22). To further assess potential relationships between the significant differentially expressed genes and their link to a common signaling pathway or functional protein–protein interaction, enrichment analysis was carried out between statistically significant DEGs, using the STRING (v11.0) database (23), with the significance threshold set at 0.4. Experiments, databases, co-expression, neighborhood, gene fusion and co-occurrences were all set as sources for active interactions. The enrichment was performed by calculating the probabilities for each of the selected evidence sources, a prior was set to account for random interactions between two proteins and then the probabilities are added. A homology correction was added to the co-occurrence score. Finally, we annotated the DEGs using the Gene Ontology (GO) (release: 2020-01-01) (24) and GO enrichment analysis tools (25) and functional annotation tool from the DAVID Bioinformatics Resources 6.8 database (26), as well as gene set enrichment analysis (GSEA) with fgsea (v1.14.0) (minGSsize = 10, maxGSsize = 1,000) and clusterProfiler (v3.16.1) packages, from the Bioconductor project. All scripts used in the study are in Supplementary Material 1 and deposited in GitHub portal (https://github.com/rsak-384/PitNET-after-therapy).

Selection of the Validation Cohort

We searched online data repositories and the literature for PitNET transcriptome data, with SSA and/or SSA/DA treatment status before surgery. From the studies closely matching our own, we found transcriptomes of ten somatotropinomas with and without SSA/DA therapy, from the “Pangenomic Classification of Pituitary Neuroendocrine Tumors” study (27) which were carefully chosen using the “ArrayExpress” archive of functional genomics data. Both therapy and nontherapy group samples were screened using the mentioned metadata, to ensure that the samples were as similar as possible between both groups, in terms of clinical data, in order to maximize the likelihood of detecting a true signal.

Cell Line Culture and Stimulation

The rat pituitary GH3 cell line was obtained from ATTC. GH3 cells were maintained in F-12K medium supplemented with 2.5% fetal bovine serum, 15% horse serum, penicillin (100 U/ml) and streptomycin (100 µg/ml).Cells were cultured at 37°C in a humidified incubator containing 5% CO2. For the gene and protein expression experiments, cells were grown on a 6-well plate at a density of 1.5 × 106 cells/well and treated with 0.1 µM and 1 µM octreotide for 4, 8, and 24 h.

RNA Extraction and Quantitative Real-Time Polymerase Chain Reaction

Total RNA was isolated from GH3 cells usingTRIzol reagent, according to the manufacturer’s protocol. First-strand cDNA was synthesized from 1.6 µg RNA with RevertAidH Minus kit (ThermoFisherScientific). qRT-PCR was performed using Absolute Blue SYBR green MasterMix reagent (ThermoFisherScientific) with a ViiA7 real-time PCR detection system (AppliedBiosystems). The primer pairs that were employed for MUC16 were forward 5′-GCCTAGGAAGAACCAAAACTCA-3′ and reverse 5′-TCCAATGTGTAGTTCCCCAGT-3′; for GRHL2 were forward 5′-CCTCTGCCTGAGTCAAGACC-3′ and reverse 5′-TAGGAGCTGTGGCTGGCTAT-3′; for MACC1 were forward 5′- CCTGGATGCCTTAGGTGGTA-3′ and reverse 5′-CCCACCCAGGACTCTGATTA-3′; for GUSB were forward 5′-GACTGATCCTTCCATGTATCCCA-3′ and reverse 5′-CCCGCATAGTTGAAGAAGTCG-3′. mRNA levels were quantified and normalized to levels of reference gene GUSB using the 2−ΔΔCt method and presented as relative expression compared with values of untreated cells.

Western Blot

Total protein samples from control and octreotide treated GH3 cells were extracted by RIPA buffer (50 mM Tris-HCl; pH 8.0, 150 mM NaCl, 5 mM EDTA, 1% Igepal CA-630, 0.5% sodium deoxycholate, 0.1% SDS, and supplemented with Halt™ protease inhibitor cocktail; Thermo Fisher Scientific). Protein concentration was quantified with BCA Protein Assay reagent (Thermo Fisher Scientific), according to the manufacturer’s instructions. Lysates (45 µg) were electrophoresed through a 10% SDS-PAGE gel and transferred to Hybond-C-extra nitrocellulose membranes (Amersham Biosciences) for 45 min in a semi-dry transfer system. Membranes were blocked with 5% non-fat milk in TBS-T buffer (20 mM Tris, pH 7.6, 137 mM NaCl, and 0.05% Tween 20) for 1 h at room temperature. Membranes were incubated with 1 µg/ml MACC1(PA5-20758, Thermo Fisher Scientific) and 1 µg/ml β-actin (ab8224, Abcam) primary antibodies for 1 h at room temperature and incubated with horseradish peroxidase-conjugated anti-mouse IgG (sc-2005, Santa Cruz Biotechnology) and anti-rabbit (sc-2004, Santa Cruz Biotechnology) secondary antibodies for 1 h at room temperature. Membranes were then washed in TBS-T buffer three times and developed with Pierce ECL western blotting substrate (Thermo Fisher Scientific). UVP software VisionWorks LS (Thermo Fisher Scientific) was used to capture signals.

Results

Characterization of the Study Group

For all the patients within this study (Table 1), it was their first time diagnosis of somatotropinoma, displaying clinical symptoms of acromegaly (ICD-10, E22.0). One patient (PA11) had tumor with complementary prolactin expression (Supplementary Material Figure 1). All patients had primary PitNET having a tumor over 10 mm, at least in one dimension. Patient’s mean age at the time of diagnosis of PitNET was 43 years (age distribution varied from 22 to 65 years). Ten were females, and two were males. Three patients had SSA/DA therapy and three SSA treatment before PitNET resection, and six patients did not have SSA/DA treatment. We are fully aware that the use of DAs in our study group was more frequent than usual. This could be explained by the fact that not all patients agreed to surgical therapy immediately, and therefore, the preoperative treatment period in most cases in our cohort was longer than 3–6 months. In turn, if IGF-1 remains modestly elevated during SSA administration, the proposed algorithm recommends adding cabergoline to the treatment plan.

TABLE 1
www.frontiersin.org

Table 1 Clinical characteristics of PitNET patients.

Differentially Expressed Genes

A total of 18,266 genes remained after sequencing read count filtering, based on the presence of gene expression in at least three samples. Raw counts were replaced for 691 genes as they were flagged as outliers. Three statistically significant surrogate variables were detected with the SVA package and were included in the analysis model. Differential expression analysis between the non-therapy group and the therapy group with DESeq2 and log fold shrinkage with the apeglm algorithm detected 40 DEGs after applying thresholds of logFC >1.5 and P-adj <0.05. In total 34 (85%) DEGs were upregulated with a median logFC of 3.02 (IQR = 1.73), and six (15%) genes were downregulated with a median logFC of −2 (IQR = 0.64). Even though two distinct sample outliers were apparent in MDS, PCA and sample distance heatmap plots (Supplementary Material Figures 2, 3 and 4), they were not removed because both samples passed sample quality evaluation using FastQC. Two distinct sample groups were observed in the heatmap (Figure 1A) that were derived from clustering. These visual results are concordant with the two patient groups, with and without therapy, that were used in the study design.

FIGURE 1
www.frontiersin.org

Figure 1 Heatmap (A) of the statistically significant differentially expressed genes. Heatmap intensities were obtained by filtering apeglm transformed log fold change values with the following thresholds: abs. lfc >1.5 and padj <0.05. The color scale represents the approximate difference in each individual sample from the average, vst normalized read counts for that gene. Data is clustered both by genes and samples using the kmeans algorithm included in the pheatmap package. Volcano plot (B) of differential expression results with apeglm transformed log fold change values. Dashed vertical lines represent absolute log fold change threshold of 1.5 and horizontal dashed line represents p values threshold of 1.31e-4 or FDR ~0.05. Red dots represent significantly DEGs (FDR < 0.05 and absolute L2FC > 1.5).

Volcano plot using apeglm transformed log-fold change values presents the distribution of differentially expressed genes based on their log fold changes while being controlled for genes with low counts or abnormally high dispersion values. These manipulations confirmed that most of the DEGs had increased expression when comparing samples receiving therapy to those with no therapy (Figure 1B).

The three most upregulated genes in the treated group were metallophosphoesterase domain containing 1 (MPPED1), transmembrane protein 132C (TMEM132C), and septin 14 (SEPTIN14); in contrast, the three most downregulated genes in the SSA/DA and SSA therapy groups were grainyhead like transcription factor 2 (GRHL2), mucin 16 (MUC16), and MET transcriptional regulator (MACC1—Metastasis-Associated in Colon Cancer) (Table 2, Figure 2). Furthermore, OLFM2 displayed the most consistent statistically significant differential upregulation in all the SSA/DA and SSA therapy group samples as indicated by low variance, high logFC difference, and P-value. logFC, P values, and box plots for all the significant DEGs are presented in Supplementary Material Tables 1 and 2.

TABLE 2
www.frontiersin.org

Table 2 Top three differentially expressed genes based on their logFC values.

FIGURE 2
www.frontiersin.org

Figure 2 Box plot diagrams displaying raw read count distribution at the log10 scale for both the SSA/DA treated and untreated groups for the seven DEGs (OLFM2, MACC1, GRHL2, MUC16, TMEM132C, SEPTIN14, MPPED1), with the top L2FC positive and negative changes and for DEGs with the lowest FDR value (OLFM2).

By performing Kendall’s Tau correlation analysis between the read count matrix and Knosp classification index and Ki-67 proliferation index, no coherent correlation trend could be observed across all DEG’s for either of the potential factors; therefore there is no reason for their inclusion in the differential expression model Supplementary Material Table 3).

Protein–Protein Interactions and Pathway Analysis

By visualizing protein–protein interactions between all of the significant DEGs using STRING database, we were able to identify five potential interactions between the proteins of the following genes: COL8A2, COL16A1, SLC8A2, COL6A1, SLC6A1, ADAMTSL2, and ADAMTS10. According to the database protein interaction analysis results, the number of expected interactions (edges) in the network was one, but the achieved number was five, with a PPI enrichment p-value of 0.00022, meaning that our network holds more reciprocal interactions than would be expected if the same amount of random proteins where to be drawn from the genome (Figure 3).

FIGURE 3
www.frontiersin.org

Figure 3 Visual representation of protein–protein interaction enrichment analysis between the statistically significant 39 out of 40 DEGs using STRING (v11.0) database. The line thickness between the connected nodes indicates the strength of the supporting data. Experiments, databases, co-expression, neighborhood, gene fusion, and co-occurrence.

The protein for the gene LINC01529 could not be identified for Homo Sapiens; therefore, it was not included in the resulting network hence the final network consisted of protein interactions for 39 DEGs. Specific sources and their contribution for each pair of interactions are presented in Supplementary Material Table 4.

To determine which molecular mechanisms could be affected by the DEGs, GO enrichment analysis was performed along with functional annotation using DAVID. Several genes from the submitted list of DEGs (FDR < 0.05) were found to be significantly enriched in multiple pathways categorized as a part of the “Cellular component” annotation data set. The majority of these results are associated with either the plasma membrane, extracellular matrix (ECM), or synapses (Table 3).

TABLE 3
www.frontiersin.org

Table 3 GO enrichment “cellular component” data set top results enriched by SSA/DA therapy in PitNET patients, ranked by statistical significance.

By performing KEGG pathway analysis as a part of DAVID, we obtained one significantly (FDR < 0.05) enriched pathway for protein digestion and absorption. However, by performing Reactome pathway analysis in the same run, we obtained multiple significant hits, which were associated with collagen processes, glycosylation, and cell surface interactions (Table 4).

TABLE 4
www.frontiersin.org

Table 4 Pathway enrichment analysis results (FDR < 0.05).

When applying GSEA analysis six signaling pathways were detected as significantly enriched (P-adjusted < 0.05) in the discovery cohort with the Molecular Signature Database (v7.2) (MSigDB) hallmark gene sets, and one signaling pathway was detected when using the KEGG gene set from the cluster Profiler package. Protein digestion and absorption were the only signaling pathway to overlap when using results from DAVID and those of GSEA. We also assessed the involvement in signaling pathways related to cell proliferation, growth, and apoptosis of individual DEGs identified in our study (full pathway list in Supplementary Material Table 5) and identified that QPRT, SEPTIN14, CLU, PTGS2 are involved in following pathways ALCALA_APOPTOSIS, GO_CELL_CYCLE, HALLMARK_APOPTOSIS, GO_CELL_CYCLE and GO_RESPONSE_TO_TUMOR_NECROSIS_FACTOR (Supplementary Material Table 5).

Validation by Independent Data Set

To test for biological reproducibility of the experiment, the aforementioned 10 samples were run through the exact same workflow to determine differential gene expression detection as before, resulting in 88 significant DEGs (FDR < 0.05 and logFC> ± 1.5), 70 (79.55%) of which were upregulated, with a median logFC of 2.71 (IQR = 1.67) and 18 (20.45%) downregulated, with a median logFC of −2.1 (IQR = 0.59). From these 85 DEGs, amyloid beta precursor protein binding Family A member 2 (APBA2) was the only one differentially expressed in both discovery and validation cohorts, with logFC of 2.83 in the discovery cohort and 3.31 in the validation cohort.

Gene set enrichment analysis was performed on the validation cohort to test whether there may be any overlapping signaling pathways between both sets of samples. No overlapping signaling pathways, however, were detected although 478 signaling pathways were found in the publication repository data when using the MSigDB hallmark set and five when using the KEGG pathways analysis (27).

Discovery set DEGs were also cross validated with leading edge genes from the significant validation cohort pathways. Three genes were confirmed to be involved in multiple pathways, which were THBD, BCAM, and APBA2 (Supplementary Material Table 6).

Functional Characterization in GH3 Cells

To functionally characterize potential effects of SSA on MUC16, GRHL2, and MACC1, we treated GH3 cells with octreotide and observed changes in expression by real-time PCR for all three candidates and Western blot for MACC1. The results demonstrate that octreotide inhibits RNA expression of all three factors with most inhibition after 4 h of treatment. The MACC1 protein expression was also inhibited after octreotide stimulation and the most decrease in MACC1 level was observed after 24 h (Figure 4).

FIGURE 4
www.frontiersin.org

Figure 4 Effect of octreotide on MUC16, GRHL2, and MACC1 gene expression in GH3 cells. (A) Cells were treated with 0.1 µM and 1 µM octreotide for 4, 8, and 24 h. Gene expression was analyzed by real-time PCR and expressed relative to value from untreated control cells, which was defined as 1. Data were normalized to GUSB and values were calculated using a comparative (2-ΔCt) method. Results are mean ± SD from two independent experiments. (B) Western blot analysis of MACC1 protein expression in octreotide treated GH3 cells. β-actin was used as a loading control. The experiments were performed twice with similar results.

Discussion

The first-line medical management of somatotropinomas is the use of SSAs, causing reduced tumor cell proliferation and hormone secretion. If the control of IGF-1 is not achieved, DAs are added to the treatment regimen. The targets of this therapy are SSTRs and D2Rs that alter intracellular signaling to promote an antiproliferative effect, and there is indirect evidence that this may lead to altered gene transcription. However, the effects of SSA/DA medication on the overall cellular transcriptional profile and consequent functional changes have not been elucidated. This is one of the first reports that have evaluated transcriptome profiles from PitNET tumor tissue in two groups of patients with and without medical therapy, which was applied before tumor resection. As a result, we were able to distinguish tumor preoperative treatment status based on transcriptome profiles. We found 40 DEGs potentially affected by SSA/DA treatment, with significant downregulation seen in several tumorigenesis promoting factors. We extended these findings by applying pathway analysis, which revealed an alteration in the expression of collagen related genes.

Many studies have analyzed the influence of SSA on cell signaling components, and these involve major effects on SSTR1, SSTR2, and SSTR5 binding and alteration of the PI3K/AKT and MAPK pathways (8). It has also been demonstrated the SSA can reduce the Ki-67 proliferation index in tumor tissue from acromegaly patients (28), but the underlying causes were not clear. The only study to date that has assessed the impact of SSA as a pre-operative treatment on transcriptome profiles from PitNETs has been carried out on three somatotropinomas and was compared to untreated PitNETs (29). The authors detected decreased Ki-67 expression in the pre-treated group and increased MUC1 and CD40 expression (29). In our data, we did not observe a significant correlation of Ki-67 with treatment status; however, we observed that overall read counts in PitNET patients without treatment were higher (Supplementary Material Table 7); to identify more pronounced association the larger sample group would be required but the tendency could be observed.

Our data indicated that three genes having pronounced effects on tumorigenesis were strongly downregulated by SSA/DA in somatotropinoma tissues (MUC16, MACC1, and GRHL2) and confirmed these findings in functional GH3 cell experiments. MUC16 is a high molecular weight O-glycosylated protein. Elevated levels of MUC16 in the blood serves as a prognostic biomarker for ovarian cancer (30). MUC16 has also been implicated in the development of other neoplasms and as a potential marker in pancreatic, breast and lung cancers (3133).

The second candidate, MACC1 was involved in epithelial–mesenchymal transition and was potent regulator of cancer metastasis and cell invasion. Increased expression of MACC1 leads to greater proliferation, induction of angiogenesis, and is antiapoptotic in various cancers (34, 35). Therefore, over the twofold decrease of MACC1expression observed in our data due to SSA/DA treatment is in line with previous results in the literature.

The third factor identified in this study GRHL2 was DNA-binding nuclear protein, containing a highly conserved DNA-binding domain (36). It was involved in the regulation of various developmental processes such as closure of the neural tube and the regulation of progenitor cell functions during the development of the pituitary gland. GRHL2 is expressed in pituitary progenitor cells and in mice with decreased pituitary progenitor cell numbers (37). GRHL2 has been implicated in many cancers such as breast, lung, colorectal, gastric, ovarian, and prostate. However, its role as tumor promoter or as an anti-tumor factor appears to be tissue dependent (36).

The three genes that were most strongly upregulated by SSA/DA treatment in somatotropinoma tissues were MPPED1, TMEM132C, and SEPTIN14, potentially having more indirect evidence for their involvement in anti-tumor effect.

The only study looking at the role for MPPED1 has found that it is involved in the development of the central nervous system (38). No other functional studies have been performed on this candidate. However, as a member of metallophosphoesterase proteins, the related MPPED2 has been shown to have pronounced tumor suppressor activities. Downregulation of MPPED2 has been demonstrated in neuroblastoma and breast cancer (39, 40). MPPED1 and MPPED2 show over 71% similarity based on Ensembl.org data, but to determine whether they have related functions needs further investigation.

Information concerning a potential function for the second upregulated candidate, TMEM132C, is modest. TMEM132 family proteins are thought to be involved in cell adhesion (41), while other members of this protein family have been associated with insomnia, hearing loss, panic disorder, and anxiety behavior (41). Somatic variants in the related factor TMEM132D have been found in small-cell lung and pancreatic cancer (42, 43). One report has implicated the cg04475027 methylation site on TMEM132C, as a marker for breast cancer (44).

The third factor, SEPTIN14 belongs to the septin protein family of GTP-binding membrane-interacting proteins and has a function in cytoskeleton organization, cytokinesis, apoptosis, cell polarity, and cell cycle regulation (45). Studies have found that aberrant expression of septin, may induce antiproliferative and tumor suppressive effects (46, 47), and somatic variants of SEPTIN14 have been demonstrated in skin and gastrointestinal cancers (48).

As SSA treatment effect induces shrinkage of the tumor, it would be expected that SSA leads to alterations in cell proliferation, growth, apoptosis, and related pathways. Indeed, four of the individual DEGs (QPRT, SEPTIN14, CLU, PTGS2) in our study are involved in these pathways (Supplementary Material Table 5); it shows that there are perturbed parts of the relevant cell proliferation pathways, but what effect these DEGs constitute on whole pathway in the cells of PitNETs remains to be discovered.

The altered expression of CLU gene that encodes for custerin involved in apoptosis and clearance of cell debris has been also identified related to invasiveness traits in two PitNET studies (49, 50). We, however, did not observe this on the pathway level in the results of enrichment analysis and protein–protein interactions, but our data supported ECM related pathways (Table 4). Further investigation is needed to elucidate whether the period of SSA/DA treatment could influence the intrinsic transcriptomic regulation of PitNETs over the course of the therapy. As at the very start of the therapy, the signaling pathways could be active that slow the tumor growth, but afterwards, other pathways involved in the reduction of tumor mass that could involve the ECM reorganization could be more pronounced. It has been shown that the most reduction of tumor mass occurred during the first year of SSA treatment (5). These assumptions should be further studied as there is only one study to date considering the impact of SSA preoperative treatment of PitNET transcriptome (29). This study has also found altered expression of one of the mucin family factors MUC1, while in our study we found MUC16 as one of the most significant DEGs.

Other transcriptomic studies that could be compared to our results are reports that have assessed the invasiveness of the PitNETs. The invasiveness could also be related to increased cell growth, and in PitNETs, it has been demonstrated that resistance to SSA is accompanied with reduced tumor shrinkage (51). In the transcriptome studies comparing invasive and non-invasive PitNETs, there are various candidates identified that are related to cell proliferation, growth, apoptosis pathways, but still rarely these candidates overlap between the studies (49, 50, 5258) (Supplementary Material Table 5). These numerous studies demonstrate various intrinsic mechanisms that in large scale studies lead to heterogeneity in the observed DEGs. Mostly this difference arises from variation in design and diagnoses, but also notable is the usage of small groups to obtain conclusions. Small group comparison (which is one of the weaknesses of this study) is more likely to lead to spurious findings, although it’s clear that obtaining large homogeneous groups of the rare heterogeneous condition is a challenging task in itself. The lack of reproducibility in these literature data could also be attributed to sample size variation and use of the different methodologies.

Although the majority of DEGs detected here were upregulated (34 DEGs) after treatment, three of six genes that were downregulated appear to have functions more related to tumor pathogenesis, according to current literature. Thus, we propose that the influence of SSA/DA treatment at the transcriptome level was directed toward suppression of tumor promoting genes (GRHL2, MUC16, MACC1), and the protective effects of the upregulated DEGs were small or indirect (MPPED1, TMEM132C, SEPTIN14); however, further studies are needed to confirm this hypothesis.

Interestingly, olfactomedin 2 (OLFM2) was found to be significantly upregulated in all tumor tissues from the SSA/DA treated group (Figure 1B). OLFM2 is involved in vascular remodeling, scar tissue formation in blood vessels, and smooth muscle differentiation. Elevated levels of OLFM2 can be found in blood plasma interventional therapy with postoperative restenosis (5961). Vascular plasticity is an important hallmark of malignant tumors (62), but whether SSA administration can lead to blood vessel remodeling events in PitNET tissue, promoting protective anti-tumoral effects needs to be functionally assessed. Furthermore, OLFM2 is involved in interactions with the transcription factor Runx2 (59), which has been reported to be involved in the regulation of pituitary tumor growth (63, 64).

Protein–protein interactions and pathway analysis revealed potential interactions between COL8A2, COL16A1, SLC8A2, COL6A1, SLC6A1, ADAMTSL2, ADAMTS10 (Figure 3) and several significantly enriched pathways belonging to the “Cellular component” category (Table 3), which again demonstrated the lowest FDR scores for pathways including COL16A1, COL6A1, and COL8A2 suggesting collagen related effects upon SSA/DA treatment (Table 4).

Collagens are key components of the ECM and can affect the behavior of tumor cells, for example proliferation, differentiation, motility, and metastasis. Tumor mass usually consists of tumor cells and stroma, which is composed of the ECM. Tumor cells and stroma interact with each other and affect tumorigenesis and cellular characteristics. It has been reported that PitNETs can have elevated fibrotic scar tissue primarily composed of stromal collagens (6567). Furthermore, research on collagen-producing cells in PitNETs and normal pituitary suggests that properties of collagen production during tumorigenesis can change, leading to the formation of fibrosis in PitNETs (68). Interestingly, collagen type VI alpha 6 (COL6A6) is downregulated in PitNETs, and overexpression of this collagen caused anti-tumoral effects and decreased metastatic capacity (69).

Some studies investigating preoperative SSA administration have demonstrated that this treatment helps to reduce tumor mass and soften the tumor tissue, facilitating tumor resection (70, 71). However, other reports have not observed this (72), or whether the tissue softening could be related to ECM alterations. Overall, our data suggested that ECM remodelling, in particular collagen regulation, might be involved in SSA/DA treatment, but further investigations are still required to determine the underlying interactions.

The ECM is an essential part of the tumor microenvironment, that along with other participating factors, may affect cell signaling, molecular transport, metabolic pathways, and immune resistance mechanisms largely contributing to tumor growth and invasiveness (73). Some tumors are composed of even up to 60% of ECM that is due to fibroblast infiltration that produces excessive amounts of ECM and is associated with worth prognostics (74, 75). Various collagen dysregulation has been linked with tumor properties of many malignant cancers (76, 77). In our study the SSA/DA treatment could impact the expression of these proteins due to remodeling and restructuration of the ECM; both ADAMTS proteins and collagens could be responding to treatment induced reduction of the tumor mass, but further studies are needed to elaborate this hypothesis.

We also discovered “protein digestion and absorption” with both pathway enrichment methods that might imply molecular mechanisms’ underlying effects of SSA/DA treatment. Besides tumor microenvironment, also metabolic reprogramming could affect properties of tumor cells. It has been recognized that alterations in metabolic regulation of tumor cells occur that drives accelerated proliferation, growth, and survival via an increase in glycolysis, amino acid, and lipid metabolism, and mitochondrial biogenesis (78, 79). There are many studies considering the role of conventionally administered medications (metformin, sulfasalazine, proton pump inhibitors and other) for non-malignant disorders that target metabolic reprogramming of tumors and have beneficial effects on cancer treatments (78, 79); whether SSA could have a similar effect is an interesting topic for investigation. It is known that SSA has effects on PI3K/Akt signaling pathway (8) that is one of the important factors inducing cancer metabolic reprogramming (78). Therefore, the SSA could exert its activity on beneficial anti-tumoral effects via shifts in tumor cell metabolism, but this needs to be investigated further.

From our computational validation of the data from online available data sets of PitNET tissue transcriptomes, from patients with and without therapy (27), we discovered one candidate factor and (APBA2) several potentially involved pathways (related to THBD, BCAM and APBA2). THBD is expressed on endothelial cells that were found to be involved in blood coagulation (80). Further investigations are still needed, however, in order to determine the role of these factors during SSA/DA treatment. BCAM is a plasma membrane glycoprotein that can bind ECM proteins and has been reported to be expressed in both fetal and adult rat pituitary (81), however, the precise role of this protein in the pathogenesis and therapy of PitNET needs further study. The only candidate that was found to be differentially expressed in our results and in validation data is APBA2; this protein interacts with amyloid precursor proteins which are involved in the development of Alzheimer’s disease. It has been reported that APBA2 is involved in gonadotropin-releasing hormone-1 neuronal migration to the pituitary and neurogenesis (82). The decreased expression of APBA2 has been identified in superior temporal gyrus in schizophrenia patients (83) and gender-specific alterations in the expression of APBA2 have been found in dopamine neurons (84). According to Protein Atlas information, this protein is involved in signal transduction and vesicular trafficking in the central nervous system, and these actions could be related to different levels of PitNET functioning, development, or progression as well as response to SSA/DA treatment. But the precise functionality of APBA2 within the nervous system signaling or trafficking as well as the pathogenesis of PitNETs is not clear and needs to be further studied.

The potential candidates involved in SSA/DA effects in our data are not directly linked to previous reports relating to resistance mechanisms of SSA. We did not find any significant profile changes in SSTR, beta-arrestin, or cytoskeleton protein filamin A that was associated with drug resistance mechanisms (85). It has also been reported that SSTR2 and D2DR agonists can reduce migration and invasiveness of PitNET cells that are mediated via cofilin and filamin A mechanisms (8688). We did not observe a significant alteration in expression of these factors, however, in our transcriptome data. Nevertheless, our results indicated novel factors targeted by SSA/DA (GRHL2, MUC16, MACC1) that merit further exploration, so as to characterize their role in SSA responder and non-responder groups that could give more insight into the relation of these factors to drug resistance.

The limitation of this study was the small number of samples used; however, our main goal was to demonstrate that SSA/DA could cause distinct alterations in transcriptome profiles from PitNET tumors, and this has been convincingly demonstrated. Although the results indicate proof of principle and highlight some of the novel factors involved in the antiproliferative effects of SSA/DA treatment, it is also clear that an increased sample size would strengthen our observation.

In conclusion, we have detected changes in transcriptional profiles induced by SSA/DA therapy in PitNET tissue. The tumor promoting factors MUC16, MACC1 and GRHL2 were downregulated in PitNETs after SSA/DA therapy and in GH3 cell following octreotide treatment. Collagen related interactions were detected after analyses of pathways and enrichment, implicating ECM involvement in the anti-tumoral effects of drug treatment. Further functional analyses are needed to determine the impact of these molecules and their potential role in response to SSA/DA treatment in patients with PitNET.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: GEO, GSE160195.

Ethics Statement

The studies involving human participants were reviewed and approved by the Central Medical Ethics Committee, Ministry of Health of the Republic of Latvia. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

RS performed the data analysis and wrote the manuscript writing. IS performed the data analysis. PL and KM prepared the sequencing library and wrote the manuscript. RPec performed the data analysis and data interpretation. IM, OR, and RPet carried out the functional experiments. IB and LS obtained the clinical data. IK obtained the clinical data and wrote the manuscript. JSt performed PitNET resection and obtained the tumor tissue. AB and JN performed the ICH characterization. JSo participated in the study design and manuscript writing. VP and JK participated in the study design and manuscript editing. VR participated in the study design, data interpretation, and manuscript writing. All authors contributed to the article and approved the submitted version.

Funding

This research was supported by the European Regional Development Fund within the project “Molecular markers of pituitary tumor development, progression and therapy response” (1.1.1.1/16/A/066).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

The authors acknowledge the Latvian Biomedical Research and Study Centre and the Genome Database of the Latvian Population for providing the infrastructure, biological material, and data.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2020.593760/full#supplementary-material

References

1. Ezzat S, Asa SL, Couldwell WT, Barr CE, Dodge WE, Vance ML, et al. The prevalence of pituitary adenomas. Cancer (2004) 101:613–9. doi: 10.1002/cncr.20412

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Aflorei ED, Korbonits M. Epidemiology and etiopathogenesis of pituitary adenomas. J Neurooncol (2014) 117:379–94. doi: 10.1007/s11060-013-1354-5

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Syro LV, Rotondo F, Serna CA, Ortiz LD, Kovacs K. Pathology of GH-producing pituitary adenomas and GH cell hyperplasia of the pituitary. Pituitary (2017) 20:84–92. doi: 10.1007/s11102-016-0748-8

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Molitch ME. Diagnosis and Treatment of Pituitary Adenomas: A Review. JAMA (2017) 317:516–24. doi: 10.1001/jama.2016.19699

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Colao A, Auriemma RS, Pivonello R. The effects of somatostatin analogue therapy on pituitary tumor volume in patients with acromegaly. Pituitary (2016) 19:210–21. doi: 10.1007/s11102-015-0677-y

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Peverelli E, Treppiedi D, Giardino E, Vitali E, Lania AG, Mantovani G. Dopamine and Somatostatin Analogues Resistance of Pituitary Tumors: Focus on Cytoskeleton Involvement. Front Endocrinol (Lausanne) (2015) 6:187. doi: 10.3389/fendo.2015.00187

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Ben-Shlomo A, Melmed S. Pituitary somatostatin receptor signaling. Trends Endocrinol Metab (2010) 21:123–33. doi: 10.1016/j.tem.2009.12.003

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Gadelha MR, Wildemberg LE, Bronstein MD, Gatto F, Ferone D. Somatostatin receptor ligands in the treatment of acromegaly. Pituitary (2017) 20:100–8. doi: 10.1007/s11102-017-0791-0

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Gatto F, Feelders R, van der Pas R, Kros JM, Dogan F, van Koetsveld PM, et al. β-Arrestin 1 and 2 and G protein-coupled receptor kinase 2 expression in pituitary adenomas: role in the regulation of response to somatostatin analogue treatment in patients with acromegaly. Endocrinology (2013) 154:4715–25. doi: 10.1210/en.2013-1672

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Gatto F, Biermasz NR, Feelders RA, Kros JM, Dogan F, van der Lely A-J, et al. Low beta-arrestin expression correlates with the responsiveness to long-term somatostatin analog treatment in acromegaly. Eur J Endocrinol (2016) 174:651–62. doi: 10.1530/EJE-15-0391

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Ibáñez-Costa A, Korbonits M. AIP and the somatostatin system in pituitary tumours. J Endocrinol (2017) 235:R101–16. doi: 10.1530/JOE-17-0254

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Fukuda T, Tanaka T, Hamaguchi Y, Kawanami T, Nomiyama T, Yanase T. Augmented Growth Hormone Secretion and Stat3 Phosphorylation in an Aryl Hydrocarbon Receptor Interacting Protein (AIP)-Disrupted Somatotroph Cell Line. PLoS One (2016) 11:e0164131. doi: 10.1371/journal.pone.0164131

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Varrault A, Dantec C, Le Digarcher A, Chotard L, Bilanges B, Parrinello H, et al. Identification of Plagl1/Zac1 binding sites and target genes establishes its role in the regulation of extracellular matrix genes and the imprinted gene network. Nucleic Acids Res (2017) 45:10466–80. doi: 10.1093/nar/gkx672

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Chahal HS, Trivellin G, Leontiou CA, Alband N, Fowkes RC, Tahir A, et al. Somatostatin analogs modulate AIP in somatotroph adenomas: the role of the ZAC1 pathway. J Clin Endocrinol Metab (2012) 97:E1411–20. doi: 10.1210/jc.2012-1111

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Hollingshead BD, Patel RD, Perdew GH. Endogenous hepatic expression of the hepatitis B virus X-associated protein 2 is adequate for maximal association with aryl hydrocarbon receptor-90-kDa heat shock protein complexes. Mol Pharmacol (2006) 70:2096–107. doi: 10.1124/mol.106.029215

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Rovite V, Wolff-Sagi Y, Zaharenko L, Nikitina-Zake L, Grens E, Klovins J. Genome Database of the Latvian Population (LGDB): Design, Goals, and Primary Results. J Epidemiol (2018) 28:353–60. doi: 10.2188/jea.JE20170079

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Liu Y, Zhou J, White KP. RNA-seq differential expression studies: more sequence or more replication? Bioinformatics (2014) 30:301–4. doi: 10.1093/bioinformatics/btt688

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol (2014) 15:550. doi: 10.1186/s13059-014-0550-8

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Huber W, Carey VJ, Gentleman R, Anders S, Carlson M, Carvalho BS, et al. Orchestrating high-throughput genomic analysis with Bioconductor. Nat Methods (2015) 12:115–21. doi: 10.1038/nmeth.3252

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Ignatiadis N, Klaus B, Zaugg JB, Huber W. Data-driven hypothesis weighting increases detection power in genome-scale multiple testing. Nat Methods (2016) 13:577–80. doi: 10.1038/nmeth.3885

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Zhu A, Ibrahim JG, Love MI. Heavy-tailed prior distributions for sequence count data: removing the noise and preserving large differences. Bioinformatics (2019) 35:2084–92. doi: 10.1093/bioinformatics/bty895

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Puka L. Kendall’s Tau. In: Lovric M, editor. International Encyclopedia of Statistical Science. Berlin, Heidelberg: Springer Berlin Heidelberg (2014). p. 713–5. doi: 10.1007/978-3-642-04898-2_324

CrossRef Full Text | Google Scholar

23. Szklarczyk D, Gable AL, Lyon D, Junge A, Wyder S, Huerta-Cepas J, et al. STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets. Nucleic Acids Res (2019) 47:D607–13. doi: 10.1093/nar/gky1131

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Carbon S, Thomas PD, Albou LP, Hill DP, Gaudet P, Van Auken K, et al. The Gene Ontology Consortium. The Gene Ontology Resource: 20 years and still GOing strong. Nucleic Acids Res (2019) 47:D330–8. doi: 10.1093/nar/gky1055

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Mi H, Muruganujan A, Ebert D, Huang X, Thomas PD. PANTHER version 14: more genomes, a new PANTHER GO-slim and improvements in enrichment analysis tools. Nucleic Acids Res (2019) 47:D419–26. doi: 10.1093/nar/gky1038

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Huang DW, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat Protoc (2009) 4:44–57. doi: 10.1038/nprot.2008.211

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Neou M, Villa C, Armignacco R, Jouinot A, Raffin-Sanson M-L, Septier A, et al. Pangenomic Classification of Pituitary Neuroendocrine Tumors. Cancer Cell (2020) 37:123–34.e5. doi: 10.1016/j.ccell.2019.11.002

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Selek A, Cetinarslan B, Canturk Z, Tarkun I, Hanazay Y, Vural C, et al. The effect of somatostatin analogues on Ki-67 levels in GH-secreting adenomas. Growth Horm IGF Res (2019) 45:1–5. doi: 10.1016/j.ghir.2019.01.001

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Salomon MP, Wang X, Marzese DM, Hsu SC, Nelson N, Zhang X, et al. The Epigenomic Landscape of Pituitary Adenomas Reveals Specific Alterations and Differentiates Among Acromegaly, Cushing’s Disease and Endocrine-Inactive Subtypes. Clin Cancer Res (2018) 24:4126–36. doi: 10.1158/1078-0432.CCR-17-2206

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Felder M, Kapur A, Gonzalez-Bosquet J, Horibata S, Heintz J, Albrecht R, et al. MUC16 (CA125): tumor biomarker to cancer therapy, a work in progress. Mol Cancer (2014) 13:129. doi: 10.1186/1476-4598-13-129

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Haridas D, Chakraborty S, Ponnusamy MP, Lakshmanan I, Rachagani S, Cruz E, et al. Pathobiological implications of MUC16 expression in pancreatic cancer. PLoS One (2011) 6:e26839. doi: 10.1371/journal.pone.0026839

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Avula LR, Rudloff M, El-Behaedi S, Arons D, Albalawy R, Chen X, et al. Mesothelin Enhances Tumor Vascularity in Newly Forming Pancreatic Peritoneal Metastases. Mol Cancer Res (2020) 18:229–39. doi: 10.1158/1541-7786.MCR-19-0688

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Aithal A, Rauth S, Kshirsagar P, Shah A, Lakshmanan I, Junker WM, et al. MUC16 as a novel target for cancer therapy. Expert Opin Ther Targets (2018) 22:675–86. doi: 10.1080/14728222.2018.1498845

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Radhakrishnan H, Walther W, Zincke F, Kobelt D, Imbastari F, Erdem M, et al. MACC1-the first decade of a key metastasis molecule from gene discovery to clinical translation. Cancer Metastasis Rev (2018) 37:805–20. doi: 10.1007/s10555-018-9771-8

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Wu Z-Z, Chen L-S, Zhou R, Bin J-P, Liao Y-L, Liao W-J. Metastasis-associated in colon cancer-1 in gastric cancer: Beyond metastasis. World J Gastroenterol (2016) 22:6629–37. doi: 10.3748/wjg.v22.i29.6629

PubMed Abstract | CrossRef Full Text | Google Scholar

36. He J, Feng C, Zhu H, Wu S, Jin P, Xu T. Grainyhead-like 2 as a double-edged sword in development and cancer. Am J Transl Res (2020) 12:310–31. doi: 10.1038/onc.2017.178

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Edwards W, Nantie LB, Raetzman LT. Identification of a novel progenitor cell marker, grainyhead-like 2 in the developing pituitary. Dev Dyn (2016) 245:1097–106. doi: 10.1002/dvdy.24439

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Chen C-M, Wang H-Y, You L-R, Shang R-L, Liu F-C. Expression analysis of an evolutionarily conserved metallophosphodiesterase gene, Mpped1, in the normal and beta-catenin-deficient malformed dorsal telencephalon. Dev Dyn (2010) 239:1797–806. doi: 10.1002/dvdy.22293

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Pellecchia S, Sepe R, Federico A, Cuomo M, Credendino SC, Pisapia P, et al. The Metallophosphoesterase-Domain-Containing Protein 2 (MPPED2) Gene Acts as Tumor Suppressor in Breast Cancer. Cancers (Basel) (2019) 11:797. doi: 10.3390/cancers11060797

CrossRef Full Text | Google Scholar

40. Liguori L, Andolfo I, de Antonellis P, Aglio V, di Dato V, Marino N, et al. The metallophosphodiesterase Mpped2 impairs tumorigenesis in neuroblastoma. Cell Cycle (2012) 11:569–81. doi: 10.4161/cc.11.3.19063

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Sanchez-Pulido L, Ponting CP. TMEM132: an ancient architecture of cohesin and immunoglobulin domains define a new family of neural adhesion molecules. Bioinformatics (2018) 34:721–4. doi: 10.1093/bioinformatics/btx689

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Iwakawa R, Kohno T, Totoki Y, Shibata T, Tsuchihara K, Mimaki S, et al. Expression and clinical significance of genes frequently mutated in small cell lung cancers defined by whole exome/RNA sequencing. Carcinogenesis (2015) 36:616–21. doi: 10.1093/carcin/bgv026

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Forbes SA, Beare D, Gunasekaran P, Leung K, Bindal N, Boutselakis H, et al. COSMIC: exploring the world’s knowledge of somatic mutations in human cancer. Nucleic Acids Res (2015) 43:D805–11. doi: 10.1093/nar/gku1075

PubMed Abstract | CrossRef Full Text | Google Scholar

44. de Almeida BP, Apolónio JD, Binnie A, Castelo-Branco P. Roadmap of DNA methylation in breast cancer identifies novel prognostic biomarkers. BMC Cancer (2019) 19:219. doi: 10.1186/s12885-019-5403-0

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Neubauer K, Zieger B. The Mammalian Septin Interactome. Front Cell Dev Biol (2017) 5:3. doi: 10.3389/fcell.2017.00003

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Gonzalez ME, Makarova O, Peterson EA, Privette LM, Petty EM. Up-regulation of SEPT9_v1 stabilizes c-Jun-N-terminal kinase and contributes to its pro-proliferative activity in mammary epithelial cells. Cell Signal (2009) 21:477–87. doi: 10.1016/j.cellsig.2008.11.007

PubMed Abstract | CrossRef Full Text | Google Scholar

47. García-Fernández M, Kissel H, Brown S, Gorenc T, Schile AJ, Rafii S, et al. Sept4/ARTS is required for stem cell apoptosis and tumor suppression. Genes Dev (2010) 24:2282–93. doi: 10.1101/gad.1970110

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Angelis D, Spiliotis ET. Septin Mutations in Human Cancers. Front Cell Dev Biol (2016) 4:122. doi: 10.3389/fcell.2016.00122

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Lekva T, Berg JP, Fougner SL, Olstad OK, Ueland T, Bollerslev J. Gene expression profiling identifies ESRP1 as a potential regulator of epithelial mesenchymal transition in somatotroph adenomas from a large cohort of patients with acromegaly. J Clin Endocrinol Metab (2012) 97:E1506–14. doi: 10.1210/jc.2012-1760

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Chen Y, Chuan H-L, Yu S-Y, Li C-Z, Wu Z-B, Li G-L. Zhang Y-Z. A Novel Invasive-Related Biomarker in Three Subtypes of Nonfunctioning Pituitary Adenomas. World Neurosurg (2017) 100:514–21. doi: 10.1016/j.wneu.2017.01.010

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Dai C, Liu X, Ma W, Wang R. The Treatment of Refractory Pituitary Adenomas. Front Endocrinol (Lausanne) (2019) 10:334. doi: 10.3389/fendo.2019.00334

PubMed Abstract | CrossRef Full Text | Google Scholar

52. Kim YH, Kim JH. Transcriptome Analysis Identifies an Attenuated Local Immune Response in Invasive Nonfunctioning Pituitary Adenomas. Endocrinol Metab (Seoul Korea) (2019) 34:314–22. doi: 10.3803/EnM.2019.34.3.314

CrossRef Full Text | Google Scholar

53. Falch CM, Sundaram AYM, Øystese KA, Normann KR, Lekva T, Silamikelis I, et al. Gene expression profiling of fast- and slow-growing non-functioning gonadotroph pituitary adenomas. Eur J Endocrinol (2018) 178:295–307. doi: 10.1530/EJE-17-0702

PubMed Abstract | CrossRef Full Text | Google Scholar

54. Cao C, Wang W, Ma C, Jiang P. Computational analysis identifies invasion-associated genes in pituitary adenomas. Mol Med Rep (2015) 12:1977–82. doi: 10.3892/mmr.2015.3564

PubMed Abstract | CrossRef Full Text | Google Scholar

55. Wierinckx A, Auger C, Devauchelle P, Reynaud A, Chevallier P, Jan M, et al. A diagnostic marker set for invasion, proliferation, and aggressiveness of prolactin pituitary tumors. Endocr Relat Cancer (2007) 14:887–900. doi: 10.1677/ERC-07-0062

PubMed Abstract | CrossRef Full Text | Google Scholar

56. Yu S-Y, Hong L-C, Feng J, Wu Y-T, Zhang Y-Z. Integrative proteomics and transcriptomics identify novel invasive-related biomarkers of non-functioning pituitary adenomas. Tumour Biol (2016) 37:8923–30. doi: 10.1007/s13277-015-4767-2

PubMed Abstract | CrossRef Full Text | Google Scholar

57. de Araújo LJT, Lerario AM, de Castro M, Martins CS, Bronstein MD, Machado MC, et al. Transcriptome Analysis Showed a Differential Signature between Invasive and Non-invasive Corticotrophinomas. Front Endocrinol (Lausanne) (2017) 8:55. doi: 10.3389/fendo.2017.00055

PubMed Abstract | CrossRef Full Text | Google Scholar

58. Galland F, Lacroix L, Saulnier P, Dessen P, Meduri G, Bernier M, et al. Differential gene expression profiles of invasive and non-invasive non-functioning pituitary adenomas based on microarray analysis. Endocr Relat Cancer (2010) 17:361–71. doi: 10.1677/ERC-10-0018

PubMed Abstract | CrossRef Full Text | Google Scholar

59. Shi N, Li C-X, Cui X-B, Tomarev SI, Chen S-Y. Olfactomedin 2 Regulates Smooth Muscle Phenotypic Modulation and Vascular Remodeling Through Mediating Runt-Related Transcription Factor 2 Binding to Serum Response Factor. Arterioscler Thromb Vasc Biol (2017) 37:446–54. doi: 10.1161/ATVBAHA.116.308606

PubMed Abstract | CrossRef Full Text | Google Scholar

60. Shi N, Guo X, Chen S-Y. Olfactomedin 2, a novel regulator for transforming growth factor-β-induced smooth muscle differentiation of human embryonic stem cell-derived mesenchymal cells. Mol Biol Cell (2014) 25:4106–14. doi: 10.1091/mbc.E14-08-1255

PubMed Abstract | CrossRef Full Text | Google Scholar

61. Li H, Gui H, Yuan G, Zheng X, Gao C, Yuan H. Increased plasma olfactomedin 2 after interventional therapy is a predictor for restenosis in lower extremity arteriosclerosis obliterans patients. Scand J Clin Lab Invest (2018) 78:269–74. doi: 10.1080/00365513.2018.1452287

PubMed Abstract | CrossRef Full Text | Google Scholar

62. Farnsworth RH, Lackmann M, Achen MG, Stacker SA. Vascular remodeling in cancer. Oncogene (2014) 33:3496–505. doi: 10.1038/onc.2013.304

PubMed Abstract | CrossRef Full Text | Google Scholar

63. Righi A, Morandi L, Leonardi E, Farnedi A, Marucci G, Sisto A, et al. Galectin-3 expression in pituitary adenomas as a marker of aggressive behavior. Hum Pathol (2013) 44:2400–9. doi: 10.1016/j.humpath.2013.05.020

PubMed Abstract | CrossRef Full Text | Google Scholar

64. Zhang H-Y, Jin L, Stilling GA, Ruebel KH, Coonse K, Tanizaki Y, et al. RUNX1 and RUNX2 upregulate Galectin-3 expression in human pituitary tumors. Endocrine (2009) 35:101–11. doi: 10.1007/s12020-008-9129-z

PubMed Abstract | CrossRef Full Text | Google Scholar

65. Paez-Pereda M, Kuchenbauer F, Arzt E, Stalla GK. Regulation of pituitary hormones and cell proliferation by components of the extracellular matrix. Braz J Med Biol Res Rev Bras Pesqui Med E Biol (2005) 38:1487–94. doi: 10.1590/s0100-879x2005001000005

CrossRef Full Text | Google Scholar

66. Li J, Wang L, Mamon H, Kulke MH, Berbeco R, Makrigiorgos GM. Replacing PCR with COLD-PCR enriches variant DNA sequences and redefines the sensitivity of genetic testing. Nat Med (2008) 14:579–84. doi: 10.1038/nm1708

PubMed Abstract | CrossRef Full Text | Google Scholar

67. Tofrizal A, Fujiwara K, Yashiro T, Yamada S. Alterations of collagen-producing cells in human pituitary adenomas. Med Mol Morphol (2016) 49:224–32. doi: 10.1007/s00795-016-0140-9

PubMed Abstract | CrossRef Full Text | Google Scholar

68. Tofrizal A, Fujiwara K, Azuma M, Kikuchi M, Jindatip D, Yashiro T, et al. Tissue inhibitors of metalloproteinase-expressing cells in human anterior pituitary and pituitary adenoma. Med Mol Morphol (2017) 50:145–54. doi: 10.1007/s00795-017-0155-x

PubMed Abstract | CrossRef Full Text | Google Scholar

69. Long R, Liu Z, Li J, Yu H. COL6A6 interacted with P4HA3 to suppress the growth and metastasis of pituitary adenoma via blocking PI3K-Akt pathway. Aging (Albany NY) (2019) 11:8845–59. doi: 10.18632/aging.102300

PubMed Abstract | CrossRef Full Text | Google Scholar

70. Abe T, Lüdecke DK. Effects of preoperative octreotide treatment on different subtypes of 90 GH-secreting pituitary adenomas and outcome in one surgical centre. Eur J Endocrinol (2001) 145:137–45. doi: 10.1530/eje.0.1450137

PubMed Abstract | CrossRef Full Text | Google Scholar

71. Colao A, Ferone D, Cappabianca P, del Basso De Caro ML, Marzullo P, Monticelli A, et al. Effect of octreotide pretreatment on surgical outcome in acromegaly. J Clin Endocrinol Metab (1997) 82:3308–14. doi: 10.1210/jcem.82.10.4283

PubMed Abstract | CrossRef Full Text | Google Scholar

72. Carlsen SM, Lund-Johansen M, Schreiner T, Aanderud S, Johannesen O, Svartberg J, et al. Preoperative octreotide treatment in newly diagnosed acromegalic patients with macroadenomas increases cure short-term postoperative rates: a prospective, randomized trial. J Clin Endocrinol Metab (2008) 93:2984–90. doi: 10.1210/jc.2008-0315

PubMed Abstract | CrossRef Full Text | Google Scholar

73. Henke E, Nandigama R, Ergün S. Extracellular Matrix in the Tumor Microenvironment and Its Impact on Cancer Therapy. Front Mol Biosci (2019) 6:160. doi: 10.3389/fmolb.2019.00160

PubMed Abstract | CrossRef Full Text | Google Scholar

74. Naba A, Clauser KR, Hoersch S, Liu H, Carr SA, Hynes RO. The matrisome: in silico definition and in vivo characterization by proteomics of normal and tumor extracellular matrices. Mol Cell Proteomics (2012) 11:M111.014647. doi: 10.1074/mcp.M111.014647

PubMed Abstract | CrossRef Full Text | Google Scholar

75. Schober M, Jesenofsky R, Faissner R, Weidenauer C, Hagmann W, Michl P, et al. Desmoplasia and chemoresistance in pancreatic cancer. Cancers (Basel) (2014) 6:2137–54. doi: 10.3390/cancers6042137

PubMed Abstract | CrossRef Full Text | Google Scholar

76. Fang S, Dai Y, Mei Y, Yang M, Hu L, Yang H, et al. Clinical significance and biological role of cancer-derived Type I collagen in lung and esophageal cancers. Thorac Cancer (2019) 10:277–88. doi: 10.1111/1759-7714.12947

PubMed Abstract | CrossRef Full Text | Google Scholar

77. Miskolczi Z, Smith MP, Rowling EJ, Ferguson J, Barriuso J, Wellbrock C. Collagen abundance controls melanoma phenotypes through lineage-specific microenvironment sensing. Oncogene (2018) 37:3166–82. doi: 10.1038/s41388-018-0209-0

PubMed Abstract | CrossRef Full Text | Google Scholar

78. Phan LM, Yeung S-CJ, Lee M-H. Cancer metabolic reprogramming: importance, main features, and potentials for precise targeted anti-cancer therapies. Cancer Biol Med (2014) 11:1–19. doi: 10.7497/j.issn.2095-3941.2014.01.001

PubMed Abstract | CrossRef Full Text | Google Scholar

79. Yoshida GJ. Metabolic reprogramming: the emerging concept and associated therapeutic strategies. J Exp Clin Cancer Res (2015) 34:111. doi: 10.1186/s13046-015-0221-y

PubMed Abstract | CrossRef Full Text | Google Scholar

80. Demeulenaere M, Devreese K, Vanbelleghem H, De Zaeytijd J, Vande Walle J, Van Biesen W, et al. Thrombomodulin and Endothelial Dysfunction: A Disease-Modifier Shared between Malignant Hypertension and Atypical Hemolytic Uremic Syndrome. Nephron (2018) 140:63–73. doi: 10.1159/000490201

PubMed Abstract | CrossRef Full Text | Google Scholar

81. Azuma M, Tsukada T, Inagaki T, Casmad F, Jindatip D, Tofrizal A, et al. Immunohistochemical Study of the Laminin α5 Chain and Its Specific Receptor, Basal Cell Adhesion Molecule (BCAM), in both Fetal and Adult Rat Pituitary Glands. Acta Histochem Cytochem (2018) 51:145–52. doi: 10.1267/ahc.18014

PubMed Abstract | CrossRef Full Text | Google Scholar

82. Forni PE, Fornaro M, Guénette S, Wray S. A role for FE65 in controlling GnRH-1 neurogenesis. J Neurosci (2011) 31:480–91. doi: 10.1523/JNEUROSCI.4698-10.2011

PubMed Abstract | CrossRef Full Text | Google Scholar

83. Bowden NA, Scott RJ, Tooney PA. Altered gene expression in the superior temporal gyrus in schizophrenia. BMC Genomics (2008) 9:199. doi: 10.1186/1471-2164-9-199

PubMed Abstract | CrossRef Full Text | Google Scholar

84. Simunovic F, Yi M, Wang Y, Stephens R, Sonntag KC. Evidence for gender-specific transcriptional profiles of nigral dopamine neurons in Parkinson disease. PLoS One (2010) 5:e8856. doi: 10.1371/journal.pone.0008856

PubMed Abstract | CrossRef Full Text | Google Scholar

85. Paragliola RM, Corsello SM, Salvatori R. Somatostatin receptor ligands in acromegaly: clinical response and factors predicting resistance. Pituitary (2017) 20:109–15. doi: 10.1007/s11102-016-0768-4

PubMed Abstract | CrossRef Full Text | Google Scholar

86. Peverelli E, Giardino E, Treppiedi D, Catalano R, Mangili F, Locatelli M, et al. A novel pathway activated by somatostatin receptor type 2 (SST2): Inhibition of pituitary tumor cell migration and invasion through cytoskeleton protein recruitment. Int J Cancer (2018) 142:1842–52. doi: 10.1002/ijc.31205

PubMed Abstract | CrossRef Full Text | Google Scholar

87. Mantovani G, Treppiedi D, Giardino E, Catalano R, Mangili F, Vercesi P, et al. Cytoskeleton actin-binding proteins in clinical behavior of pituitary tumors. Endocr Relat Cancer (2019) 26:R95–108. doi: 10.1530/ERC-18-0442

PubMed Abstract | CrossRef Full Text | Google Scholar

88. Peverelli E, Giardino E, Treppiedi D, Locatelli M, Vaira V, Ferrero S, et al. Dopamine receptor type 2 (DRD2) inhibits migration and invasion of human tumorous pituitary cells through ROCK-mediated cofilin inactivation. Cancer Lett (2016) 381:279–86. doi: 10.1016/j.canlet.2016.08.005

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: somatostatin/dopamine (SSA/DA) therapy, acromegaly, transcriptome, next generation sequencing (NGS), somatotropinoma

Citation: Saksis R, Silamikelis I, Laksa P, Megnis K, Peculis R, Mandrika I, Rogoza O, Petrovska R, Balcere I, Konrade I, Steina L, Stukens J, Breiksa A, Nazarovs J, Sokolovska J, Pirags V, Klovins J and Rovite V (2021) Medication for Acromegaly Reduces Expression of MUC16, MACC1 and GRHL2 in Pituitary Neuroendocrine Tumour Tissue. Front. Oncol. 10:593760. doi: 10.3389/fonc.2020.593760

Received: 18 August 2020; Accepted: 11 December 2020;
Published: 15 February 2021.

Edited by:

Shicheng Guo, University of Wisconsin-Madison, United States

Reviewed by:

Pedro Marques, Queen Mary University of London, United Kingdom
Junwei Han, Harbin Medical University, China

Copyright © 2021 Saksis, Silamikelis, Laksa, Megnis, Peculis, Mandrika, Rogoza, Petrovska, Balcere, Konrade, Steina, Stukens, Breiksa, Nazarovs, Sokolovska, Pirags, Klovins and Rovite. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Vita Rovite, vita.rovite@biomed.lu.lv

These authors have contributed equally to this work

Disclaimer: All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.