Central and Peripheral Changes in FOS Expression in Schizophrenia Based on Genome-Wide Gene Expression

Schizophrenia is a chronic, debilitating neuropsychiatric disorder. Multiple transcriptomic gene expression profiling analysis has been used to identify schizophrenia-associated genes, unravel disease-associated biomarkers, and predict clinical outcomes. We aimed to identify gene expression regulation, underlying pathways, and their roles in schizophrenia pathogenesis. We searched the Gene Expression Omnibus (GEO) database for microarray studies of fibroblasts, lymphoblasts, and post-mortem brains of schizophrenia patients. Our analysis demonstrated high FOS expression in non-neural peripheral samples and low FOS expression in brain tissues of schizophrenia patients compared with healthy controls. FOS exhibited predictive value for schizophrenia patients in these datasets. Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis revealed that “amphetamine addiction” was among the top 10 significantly enriched KEGG pathways. FOS and FOSB, which are implicated in the amphetamine addiction pathway, were up-regulated in schizophrenia fibroblast samples. Protein–protein interaction (PPI) network analysis revealed that proteins closely interacting with FOS-encoded protein were also involved in the amphetamine addiction pathway. Pearson correlation test indicated that FOS showed positive correlation with genes in the amphetamine pathway. The results revealed that FOS was acceptable as a biomarker for schizophrenia and may be involved in schizophrenia pathogenesis.

Introduction

Schizophrenia is a chronic and debilitating neuropsychiatric disorder affecting 1% of the population, posing a severe social and economic burden on societies worldwide (McGrath et al., 2008; Howes and Murray, 2014). Common symptoms include positive and negative symptoms and cognitive deficits (Lewis et al., 2012). Besides its complex symptomatology, schizophrenia is considered a neurodevelopmental disorder with heterogeneous, polygenic, and highly heritable etiology (Quadrato et al., 2016). Schizophrenia affects gross architectural structures, specific cell types, and ion channels across different brain regions, including the prefrontal cortex, thalamus, thalamic reticular nucleus, and basal ganglia (Heyes et al., 2015; Mouchlianitis et al., 2016; Goff et al., 2017). However, there is heterogeneity in the molecular and genetic phenotypes of patients. The precise etiology and pathogenesis underlying schizophrenia are not fully known. Identifications of gene changes as biomarkers for schizophrenia may be helpful for diagnostic assessment of patients.

Strong evidence suggests that dysfunction of multiple neurotransmitter systems may contribute to the pathophysiology of schizophrenia, including dopamine (DA), glutamate, and serotonin neurotransmission (de Bartolomeis et al., 2014). For example, N-methyl-D-aspartate (NMDA) receptor antagonists such as phencyclidine (PCP), ketamine, and MK-801 have psychotomimetic effects and have been used to generate pharmacological animal models of schizophrenia (Zuo et al., 2009). Another well-established pre-clinical schizophrenia rodent model is based on amphetamine (AMPH)-induced dopaminergic dysregulation (Renard et al., 2016). AMPH-induced schizophrenia-like sensorimotor cognitive deficits are severely disrupted in schizophrenia (Pedrazzi et al., 2015; Renard et al., 2016). Genes involved in these pathways have been examined to better understand the pathogenesis of schizophrenia and points to new targets for therapeutic investigation (Karam et al., 2010; Hall et al., 2015; Network et al., 2015). Environmental risk factors and their interactions with gene expression also play important roles in schizophrenia pathophysiology (Van Winkel et al., 2010). Understanding the affected pathways and specific gene expression profiles in the pathogenesis of schizophrenia may help to uncover disease-associated biomarkers for risk assessment, regulatory mechanisms, and personalization of treatment.

The recently developed technique of gene expression profiling analysis of the whole transcriptome has been widely used to identify schizophrenia-associated genes, unravel disease-associated biomarkers, and predict clinical outcomes. Several whole-genome expression studies have utilized gene expression data from human fibroblasts, blood, and post-mortem brain samples to identify gene alterations in schizophrenia patients compared with healthy controls (Rollins et al., 2010; Cattane et al., 2015). Expression studies on post-mortem brains and peripheral cells demonstrated overlapping gene expression results; however, contrasting results have also been observed in different studies. The heterogeneity between gene expression profiles in peripheral cells and in post-mortem brains can be attributed to intrinsic differences in expression levels between the central nervous system (CNS) and peripheral tissues; and confounding factors in patients themselves, including the course of disease, age, living habits, environmental events, and clinical medications. These problems can be partly solved by integrated analysis of gene expression data from multiple studies and multiple tissues.

On this basis, the aim of our study was to elucidate gene expression changes in the pathogenesis of schizophrenia and to acquire new potential biomarkers for diagnostic prediction. We searched for microarray data for schizophrenia from the Gene Expression Omnibus (GEO) database. Different gene expression microarray studies with schizophrenia samples were selected and compared to perform a reliable genome-wide gene expression profiling analysis. First, we analyzed different microarray gene expression studies with a sample of fibroblasts from schizophrenia patients and controls to identify differentially expressed genes. Second, we compared the gene alteration results observed in different schizophrenia post-mortem brain expression studies. Third, we compared the results by analyzing gene expression changes in samples from schizophrenia mice models. Our analysis enables identification of gene expression regulation in different body areas of schizophrenia patients, disease progression, and confirmation of pre-clinical studies. This study may help us identify risk genes for schizophrenia and the underlying pathways. Our study may be applicable in elucidating the pathogenesis of other neuropsychiatric disorders.

Materials and Methods

Inclusion Criteria

To identify specific gene expression changes in schizophrenia, we performed a systematic search in the GEO database for schizophrenia studies. The keywords used for the GEO database search were: “schizophrenia” AND “expression profiling by array” for “study type” AND “tissue” for “attribute name” AND “homo sapiens” for “organism.” Details of the search strategy and analysis process are outlined in Supplementary Figure S1. The final studies were selected based on the criteria that (1) transcriptomic profiles were analyzed in skin fibroblasts of schizophrenia patients; (2) original data in CEL format could be downloaded; and (3) matched samples of healthy controls were used.

After screening 26 studies in GEO, the GEO dataset GSE62333¹ was selected, which studied 20 schizophrenia patients matched with 20 healthy controls with human skin fibroblast samples based on the GPL11532 [HuGene-1_1-st] Affymetrix Human Gene 1.1 ST Array (Cattane et al., 2015). Array GPL11532 has 33,297 probes. All schizophrenia patients fulfilled the DSM-IV criteria for schizophrenia. The control samples consisted of healthy volunteers without drug or alcohol abuse and without family history of psychiatric diseases. Other diseases such as hypothyroidism or hyperthyroidism, metabolic disorders, and serious illnesses were also excluded in both groups.

To maintain consistency, we also selected GEO datasets with lymphoblast samples from schizophrenia patients. Finally, dataset GSE73129², which is based on platform GPL570 [HG-U133_Plus_2] Affymetrix Human Genome U133 Plus 2.0 Array was selected for validation. Array GPL570 contained 54,675 probes. Human lymphoblast samples from 36 schizophrenia patients and 41 healthy controls were included. We compared gene expression levels in human prefrontal cortex samples in GSE92538³, comprising 45 schizophrenia patients and 46 healthy controls. We also validated the results with a genome-wide transcriptomics analysis of mouse prefrontal cortex and hippocampus samples. GSE10784⁴ was selected, which contained 10 Df(16)A/+ mice and 10 control mice.

Data Analysis

Gene Expression Omnibus⁵ provides raw data (CEL file) and normalized data downloads of the: (1) citation information, (2) geo-annotation, and (3) statistical matrix. Four GEO datasets were included in the analysis, GSE62333, GSE73129, GSE92538, and GSE10784. All of the original data came from publicly available datasets.

R language was mainly used for statistical analysis. The original CEL data was subjected to background correction, normalization, and expression calculation before further analysis. We completed the data preprocessing by RMA (R package “affy”) (Gautier et al., 2004). Probes that did not match any known genes were removed. Only probes with the highest interquartile range were included for further analysis, if more than one probe matched with a gene. The selected probes were subjected to annotation using the R annotation the package “annotate.” The gene expression profiling data was transformed by log₂ and extracted for DEG identification and analysis. The DEGs were identified using the R package “limma.” Spearman correlation analysis was performed using the “circlize” package (Gu et al., 2014), ROC curves were derived using the “pROC” package (Robin et al., 2011). Principal component analysis (PCA) was performed using the prcomp function of R “stats” package, and visualization was done using the ggplot2 package. The heatmap of DEGs was generated by clustering (using the R package “pheatmap”) with p-value <0.05 between schizophrenia and control samples. Associated KEGG (Kyoto Encyclopedia of Genes and Genomes) enrichment analysis of the 250 DEGs was generated by Omics Bean and illustrates the results of KEGG pathway analysis. STRING database (version 9.1) was used to analyze protein–protein interaction (PPI) networks of the candidate genes. Dot plots of expression levels of identified genes were drawn by GraphPad Prism. All statistical tests were two-sided. A p-value less than 0.05 is considered to be statistically significant.

Results

Identification and Functional Analysis of Differentially Expressed Genes (DEGs) in Schizophrenia and Healthy Control

The GEO dataset GSE62333 was selected to identify DEGs in human skin fibroblast samples from schizophrenia patients and matched healthy controls. After the samples were normalized by the robust multichip averaging (RMA) process, the expression levels (transformed by log₂) of all included samples were extracted. The DEGs between two groups were calculated by the Classical Bayesian algorithm with 33,297 probes included. After screening using the criteria of | FC| > 1.2 and p-value < 0.05, 1,022 probes with significantly different expression levels between the two groups (517 up-regulated probes and 505 down-regulated probes) were identified. Probes that did not match any known genes were removed in further analysis. We selected the 250 probes with the highest inter-quartile range (125 up-regulated probes and 125 down-regulated probes, respectively). A heatmap that shows the 250 top-regulated genes matched with 250 probes between schizophrenia and control groups is presented in Figure 1A. The hierarchical clustering of these 250 genes demonstrated significant differences between schizophrenia and control groups.

FIGURE 1

Figure 1. Identification and functional analysis of 250 differentially expressed genes (DEGs) between schizophrenia and healthy control samples in GEO dataset GSE62333. (A) Heatmap showing 250 DEGs (| FC| > 1.3, p-value < 0.05) between schizophrenia and control samples. Each column represents the expression level of a gene, and each row represents a sample. The color scale below the heatmap represents the raw Z-score ranging from green to red (low to high expression level). Dendrograms above correspond to the Pearson correlation-based hierarchical clustering of the 250 genes. (B) Associated KEGG signaling pathways of the 250 DEGs. The network model generated by Omics Bean illustrates the results of KEGG pathway and biological process enrichment analysis. (C) A histogram according to the KEGG functional analysis displaying the percentage of genes affected in these pathways.

We performed KEGG enrichment analysis of the 250 DEGs for functional analysis to identify the biological processes, pathways, and networks shared by these genes. The top 10 significantly enriched KEGG pathways (p-value < 0.05) were “AMPH addiction,” “oxytocin signaling,” “pathways in cancer,” “estrogen signaling,” “vascular smooth muscle contraction,” “Ras signaling,” “AGE-RAGE signaling pathway in diabetic complications,” “choline metabolism in cancer,” “Chagas disease,” and “Rap1 signaling.” The network model generated by Omics Bean illustrates the results of KEGG pathway enrichment analysis (Figure 1B). A histogram was plotted based on the functional analysis to display the percentage of genes affected in these pathways (Figure 1C).

FOS Up-Regulation in Schizophrenia Peripheral Samples

After confirming the pathways of 250 DEGs by KEGG enrichment analysis, we ranked the top 10 DEGs between the schizophrenia and control samples to determine which pathways these genes were involved in. The top five up-regulated genes (FOSB, MMP1, FOS, RANBP3L, and SNORA38B) and top five down-regulated genes (SNORA68, SNORA23, SNORA20, ACTG2, and SNORA3A) were arranged by degree of change. Detailed information on these genes is listed in Table 1. Dot plots illustrated that all 10 DEGs were differentially expressed between the two groups (p-value < 0.05), suggesting that these genes may be useful as biomarkers of schizophrenia (Figure 2). We also performed PCA and while there was a separation of some samples into diagnostic groups, overall there was no clear separation between control and SCZ samples (Supplementary Figure S2).

TABLE 1

Table 1. Top 10 differentially expressed genes, either up-regulated or down-regulated.

FIGURE 2

Figure 2. Scatter plot of expression levels of the identified top 10 differentially expressed genes (DEGs) in GEO dataset GSE62333. The top five up-regulated genes (FOSB, MMP1, FOS, RANBP3L, and SNORA38B) and top five down-regulated genes (SNORA68, SNORA23, SNORA20, ACTG2, and SNORA3A) were arranged by degree of change. Detailed information on these genes is listed in Table 1.

Among the top 10 DEGs, five of them (SNORA38B, SNORA68, SNORA23, SNORA20, and SNORA3A) were non-coding genes. Five of the top 10 DEGs were protein-coding genes (RANBP3L, ACTG2, MMP1, FOSB, and FOS). Two genes in the FOS family (FOS and FOSB) were up-regulated in schizophrenia samples. Among the top 10 DEGs, FOS was most closely related to schizophrenia pathogenesis. Previous research has shown that NMDA receptor antagonist-induced and AMPH-induced pharmacological animal models of schizophrenia can up-regulate FOS gene expression (Zuo et al., 2009; Renard et al., 2016). And we have observed that “AMPH addiction” was one of the most significantly enriched KEGG pathways in schizophrenia samples. This suggested that FOS may have some connection with AMPH addiction pathway in schizophrenia samples.

Based on these findings, we constructed PPI networks of FOS to explore the known and predicted proteins interacted with FOS, as well as the underlying pathways using STRING. The STRING database is a biological database and web resource which can be used for exploring the known and predicted interaction networks of a particular protein. Ten gene-encoded proteins (JUN, CREB1, ATF2, JUNB, JUND, MAPK1, MAPK3, MAPK8, MAPK9, and IL2) were in the interaction network of FOS (Figure 3A). The edges between FOS and its 10 predicted functional partners represent protein–protein associations. The edges between the proteins do not necessarily imply binding interactions; edges in different colors point to different methods of identifying interactions. Already-known interactions can be found from curated databases or experimental results. Predicted interactions arise from text mining, co-expression, and protein homology. All proteins interacting with FOS were generated from known results from curated databases and experiments. Several interactions were predicted by text mining or co-expression. Further functional enrichment analysis of the protein network in KEGG pathways revealed that FOS protein is involved in the AMPH pathway. Three proteins that interact with FOS protein (JUN, CREB1, and ATF2) were also implicated in the AMPH pathway (Figure 3B). Proteins closely interacting with FOS protein were also found to participate in other psychoactive drug pathways, such as cocaine addiction and alcoholism. KEGG functional enrichment analysis of all genes in the network is shown in Figure 3C.

FIGURE 3

Figure 3. Functional protein–protein network analysis of FOS gene using STRING database. (A) The protein–protein interaction network of FOS and 10 other proteins. Proteins are presented by network nodes of different colors. The edges between FOS and its predicted functional partners represent protein–protein associations. Categorized identification of protein–protein interactions using different methods are labeled by circles. Already known interactions can be found from curated databases or experimental results. Predicted interactions can come from text mining, co-expression, and protein homology. (B) Further functional enrichment analysis of the protein network in KEGG pathways. FOS, JUN, CREB1, and ATF2 are involved in the amphetamine pathway. Proteins closely interacting with FOS protein also participate in other psychoactive drug pathways, such as cocaine addiction and alcoholism. (C) KEGG functional enrichment analysis of all genes in the protein–protein interaction network.

Our finding from PPI networks of FOS suggested FOS is tightly related to genes involved in FOS. Therefore, we further explicate the roles of the FOS gene and the AMPH pathway in the pathogenesis of schizophrenia. We calculated genes involved in the AMPH pathway from the 1,022 probes with significantly different gene expression levels (using the R package “pheatmap”). Seven genes in the AMPH pathway are shown in the heatmap (Figure 4A). To examine the relationship between FOS and other genes in the AMPH pathway, we calculated a Spearman correlation of gene expression levels of fold changes using the GSE62333 dataset. The results indicated that FOS was positively correlated with other genes in the AMPH pathway, especially FOSB, JUN, and GRIA3 (Figure 4B). Corrgram is a visual display technique to represent the pattern of correlations. To get further understanding the relationship between FOS gene and the AMPH pathway, Corrgram was derived according to r-value between FOS and six genes in the AMPH pathway (Figure 4C). There was a positive correlation between FOS and the following genes: FOSB (r = 0.73, p = 1.02E-07), JUN (r = 0.66, p = 0.000003) GRIA3 (r = 0.38, p = 0.014374). There was a negative correlation between FOS and CALM3 (r = −0.49, p = 0.001441). There was no correlation between FOS and the following genes: PRKCA (r = 0.28, p = 0.076489), GRIN3A (r = −0.006, p = 0.966528).

FIGURE 4

Figure 4. FOS was positively correlated with FOSB, JUN, and GRIA3 in the amphetamine pathway in GEO dataset GSE62333. Spearman correlation of related gene expression levels of fold changes were calculated for further analysis. (A) Heatmap showing genes involved in the amphetamine pathway from 1,022 probes with significantly different gene expression levels (using the R package “pheatmap”). (B) Correlation between FOS and other genes in the amphetamine pathway (FOSB, JUN, PRKCA, GRIA3, GRIN3A and CALM3). (C) Corrgram shows the correlation of FOS with other six genes in the amphetamine pathway. The cells of the matrix were shaded to show the correlation value. The correlation r-value between FOS and six genes in the amphetamine pathway were FOSB (r = 0.73, p = 1.02E-07), JUN (r = 0.66, p = 0.000003) GRIA3 (r = 0.38, p = 0.014374), PRKCA (r = 0.28, p = 0.001441), GRIN3A (r = –0.006, p = 0.966528), CALM3 (r = –0.49, p = 0.001441), respectively.

To overcome the limitations of individual studies, we validated our findings in dataset GSE73129. The results showed a significant difference of FOS expression between schizophrenia and control lymphoblast samples. No differences were found in other genes involved in the AMPH pathway (Figure 5A). The results from two GEO datasets GSE62333 and GSE73129 demonstrated up-regulated expression of FOS in fibroblast and lymphoblast samples from schizophrenia patients, which were involved in the AMPH pathway. These non-neural peripheral cells may be useful for studying molecular signatures in psychiatric disorders (Cattane et al., 2015; Horiuchi et al., 2016). Post-mortem brains are also valuable tools to identify molecular alterations in these diseases. To further address the role of FOS in schizophrenia pathogenesis, we analyzed human postmortem dorsolateral prefrontal cortex samples and brain tissues of a schizophrenia 22q11-deletion mouse model (Stark et al., 2008; Hagenauer et al., 2016). FOS gene expression was significantly down-regulated in CNS tissues, regardless of the source of postmortem tissue (UCDavis or UMichigan) or the region of mouse brain (prefrontal cortex or hippocampus) (Figures 5B,C). To measure the predictive value of FOS in the datasets we used, ROC curves for FOS expression in schizophrenia samples and control samples were performed. The area under the curve (AUC) is 78.6% and 71.8% in GSE62333 and GSE73129 dataset, respectively. The AUC in GSE92538 is 77.6% (UCDavis) and 79.4% (UMichigan), respectively. The AUC in GSE10784 is 86.0% (prefrontal cortex) and 71.0% (hippocampus), respectively. The highest AUC was observed in the prefrontal cortex samples of GSE10784 dataset. These results suggested that FOS is acceptable as a biomarker of schizophrenia. Details were presented in Figure 6.

FIGURE 5

Figure 5. Validation of expression levels of FOS and other genes involved in the amphetamine pathway in GEO datasets GSE73129, GSE92538, and GSE10784. The expression of FOS, FOSB, JUN, PRKCA, GRIA3, and CALM3 between schizophrenia and control lymphoblast samples have been presented in the histograms. MRNA expression level means the intensity of the detected fluorescence intensity transformed by log₂. NS indicates not significant, ^∗ indicates p-value < 0.05, ^∗∗ indicates p-value < 0.01, ^∗∗∗indicates p-value < 0.001. (A) Expression levels of FOS and other genes involved in the amphetamine pathway between schizophrenia and control lymphoblast samples in GEO dataset GSE73129. (B) Expression levels of FOS and other genes involved in the amphetamine pathway in human post-mortem dorsolateral prefrontal cortex samples from UCDavis or UMichigan. Samples were collected from schizophrenia patients and healthy controls. (C) FOS and other genes involved in the amphetamine pathway expression changes in the prefrontal cortex (PFC) and hippocampus (HPC) of the Df(16)A+/– schizophrenia mouse model.

FIGURE 6

Figure 6. ROC curve analysis reveals the predictive value of FOS for schizophrenia patients in different GEO datasets. Area under the curve (AUC) is measured in order to determine the predictive value of the used datasets. (A) ROC curve of FOS in GEO dataset GSE62333 which analyzed lymphoblast samples from schizophrenia patients and controls (P < 0.05). (B) ROC curve of FOS in GEO dataset GSE73129 which analyzed lymphoblast samples from schizophrenia patients and controls (P < 0.05). (C) ROC curve of FOS in GEO dataset GSE92538 which analyzed human post-mortem dorsolateral prefrontal cortex samples from UCDavis (UD) or UMichigan (UM). Samples were collected from schizophrenia patients and healthy controls (P < 0.05). (D) ROC curve of FOS in GEO dataset GSE10784 which analyzed the prefrontal cortex (PFC) and hippocampus (HPC) samples from the Df(16)A+/– schizophrenia mouse model (P < 0.05).

Collectively, our findings revealed that FOS, which is involved in the AMPH pathway, is significantly up-regulated in human fibroblast samples and can be validated in peripheral lymphoblast samples. Further analysis showed contrasting expression of FOS in CNS tissues. FOS was down-regulated in human postmortem dorsolateral prefrontal cortex samples and brain tissues of a schizophrenia 22q11-deletion mouse model.

Discussion

In this study, we performed a systematic search on the GEO database for microarray studies of schizophrenia with samples from patients’ fibroblasts, lymphoblasts, post-mortem human brains, and brains from a mouse model of schizophrenia. Unlike previous studies, which analyzed peripheral or CNS samples only from schizophrenia patients, we identified FOS expression changes in different tissues. Skin fibroblasts, lymphocytes, and lymphoblasts are frequently used in schizophrenia gene expression analysis because of the difficulty in obtaining CNS tissues or cells from schizophrenia patients. Postmortem brains are another good source that help to address this problem. Our study, which analyzed FOS gene expression in peripheral and CNS tissues, provides a new window into the molecular changes underlying schizophrenia pathogenesis.

In our analysis, we re-analyzed publicly available microarray gene expression data and obtained different conclusions. A previous study of GEO datasets GSE62333 also reported the up-regulation of FOSB and FOS in the fibroblasts of schizophrenia patients and mainly focused on EGR1 and other genes as biomarkers for disease diagnosis (Cattane et al., 2015). Our analysis focused on gene expression changes and related pathways contributing to the pathogenesis of schizophrenia. The analysis of the top 10 DEGs in schizophrenia and healthy samples revealed that FOS and FOSB (which are involved in the AMPH addiction pathway), were significantly up-regulated in schizophrenia fibroblast samples. Protein functional network analysis revealed proteins closely interacting with FOS that were also involved in the AMPH addiction pathway. These findings were consistent with our KEGG enrichment analysis results, which listed “AMPH addiction” as one of the top 10 significantly enriched KEGG pathways. Additionally, Spearman correlation analyzing samples in dataset GSE62333 indicated that FOS was positively correlated with FOSB, JUN, and GRIA3 in the AMPH pathway, and negatively correlated with CALM3 in the AMPH pathway. Thus, our results support previous findings and reveal a previously unrecognized pathway connected with FOS in the pathogenesis of schizophrenia.

As one of the most studied immediate early genes in the brain, FOS is highly expressed in addiction to psychoactive drugs such as AMPH, alcohol, and cocaine (Gallo et al., 2018). It also has been reported to be differentially expressed in schizophrenia samples, either through experimental or expression profiling studies (Zuo et al., 2009; Cattane et al., 2015; Renard et al., 2016; Monfil et al., 2018). Besides, FOS is thought to play an important role in the pathophysiology of schizophrenia (Zuo et al., 2009; Renard et al., 2016; Monfil et al., 2018). We observed that FOS expression was highly up-regulated in schizophrenia group of dataset GSE62333 which was then validated in another dataset GSE73129. ROC curves indicated fair predictive values of the datasets we analyzed. These combined results suggested FOS is acceptable as a biomarker of schizophrenia. Detection of FOS in blood samples may be helpful for schizophrenia diagnosis. Further exploration using postmortem human brains and 22q11-deletion mouse brain samples suggested that FOS is up-regulated in non-neural peripheral samples and down-regulated in brain tissues of schizophrenia patients compared with those of healthy controls. The alteration of FOS expression in peripheral and CNS tissues of schizophrenia indicates that this gene is sensitive for schizophrenia.

Among the top 10 DEGs we detected, five of them (SNORA38B, SNORA68, SNORA23, SNORA20, and SNORA3A) were non-coding genes. The ratio of non-coding to protein-coding genes is considered a function of developmental complexity. Prokaryotes have less than 25% non-coding DNA, while humans have approximately 98.5% non-coding DNA. RNA-based regulation is devoted to the majority of advances in human genomic programming (Mattick, 2004). Increasing evidence suggests that non-coding RNAs play an important role in neural development and function, and neurological disease progression (Mehler and Mattick, 2006). Although non-coding RNAs such as small nucleolar RNAs (snoRNAs), microRNAs and long non-coding RNAs (LncRNAs), have been studied in disease etiology (Barry, 2014; Lai et al., 2016), the involvement of SNORA38B, SNORA68, SNORA23, SNORA20, and SNORA3A in the pathogenesis of schizophrenia has yet to be reported. However, snoRNAs may play important roles in brain development and neurological disease progression (Pang et al., 2006). Recent studies have indicated that snoRNAs are highly expressed in tumor cells and are involved in cell behavior and oncogenesis (Williams and Farzaneh, 2012; Gong et al., 2017; Zhou et al., 2017).

Five of the top 10 DEGs were protein-coding genes (RANBP3L, ACTG2, MMP1, FOSB, and FOS). RANBP3L is a gene involved in detoxification which mediates bone morphogenetic protein (BMP)-specific nuclear export of SMAD to terminate BMP signaling (Chen et al., 2015). The involvement of RANBP3L in schizophrenia has yet to be reported, but exon microarray analysis of human dorsolateral prefrontal cortex revealed up-regulation of this gene in alcoholism (Manzardo et al., 2014). ACTG2 encodes gamma 2 enteric actin, which is a smooth muscle actin found in enteric tissues and is crucial for enteric muscle contraction (Thorson et al., 2014; Wangler et al., 2014; Matera et al., 2016). ACTG2 has yet to be associated with schizophrenia. Another actin-coding gene, ACTG1, is associated with brain development and hearing loss (Perrin et al., 2010; Park et al., 2013). Matrix metalloproteases (MMPs) are involved in the degradation of basement membrane and extracellular matrix (ECM) components in physiological and pathological processes, such as tissue development, wound healing, and tumor invasiveness (Poola et al., 2005). As a member of the MMP family, MMP1 encodes an enzyme that breaks down the ECM, and promotes tumor cell division and metastasis. Studies have shown that activator protein-1 (AP-1) genes FOS and JUN can regulate MMP1 gene expression in invasive tumors (Poola et al., 2005; Uhlirova and Bohmann, 2006; Baker et al., 2018). MMP-mediated ECM disruption is also involved in the pathogenesis of schizophrenia. MMP modulators may therefore be a potential therapeutic target for the treatment of schizophrenia (Chopra et al., 2015). Besides, FOS proteins can dimerize with JUN proteins to form the transcription factor complex AP-1, which regulates cell proliferation, differentiation, and transformation (Milde-Langosch, 2005). Few studies have reported functions of FOSB in the CNS. However, a truncated splice variant of FOSB named Delta-FOSB is involved in behavior and addiction to drugs (Nestler, 2008; Ruffle, 2014). FOS was up-regulated in schizophrenia samples and is a recognized marker of neural activation (Gallo et al., 2018).

There are several limitations in our study. First, the gene expression profiles could be affected by many factors, such as differences in methodology and/or sample preparation, patient characteristics, platform, and data analysis (Sullivan et al., 2006; Mistry and Pavlidis, 2010; Kumarasinghe et al., 2012; Mistry et al., 2013). Studies on schizophrenia are often associated with small sample sizes which may result in low statistical power. Second, several studies reported that peripheral cells may not reflect the gene expression profile of the human brain. Furthermore, using postmortem brains to understand dynamic changes in disease progression and development of complications is challenging (Rollins et al., 2010; Tylee et al., 2013; Cattane et al., 2015; Horiuchi et al., 2016). Finally, confirmatory experiments and comparison of results with microarray gene expression modifications from re-analysis will enable validation of our results.

Conclusion

To conclude, our results indicate that FOS was significantly up-regulated in schizophrenia fibroblast and lymphoblast samples. Exploration of FOS gene expression in CNS tissues revealed that this gene was largely down-regulated in datasets GSE92538 and GSE10784. The signatures we identified are consistent with current hypotheses of molecular dysfunction in schizophrenia, including alteration of the FOS gene and involvement of the AMPH addiction pathway. FOS and AMPH-related genes may thus represent novel biomarkers for diagnosis of schizophrenia in clinical practice.

Author Contributions

JH and FL performed data extraction and prepared the manuscript. BW, HT, and ZT reviewed the literature and performed the quality assessment. YQ and LL reviewed the revised manuscript and checked the grammars. HW and JC supervised all the work.

Funding

This work was supported by the National Natural Science Foundation of China (81270019) and National Science Foundation for Young Scientists of China (81501163). The funders had no involvement in any aspect of this study or manuscript preparation.

Conflict of Interest Statement

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 would like to thank the reviewers for their valuable comments and suggestions to improve the quality of the paper.

Supplementary Material

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

FIGURE S1 | Flowchart showing the identification and selection of microarray studies from Gene Expression Omnibus (GEO) and the data analysis procedure.

FIGURE S2 | Transcriptomic analysis via principle component analysis including control and SCZ samples in GEO dataset GSE62333 using top 10 differentially expressed genes.

Footnotes

1. ^ https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE62333
2. ^ https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE73129
3. ^ https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE92538
4. ^ https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE10784
5. ^ https://www.ncbi.nlm.nih.gov/geo/

References

Baker, J., Falconer, A. M. D., Wilkinson, D. J., Europe-Finner, G. N., Litherland, G. J., and Rowan, A. D. (2018). Protein kinase D3 modulates MMP1 and MMP13 expression in human chondrocytes. PLoS One 13:e0195864. doi: 10.1371/journal.pone.0195864

PubMed Abstract | CrossRef Full Text | Google Scholar

Barry, G. (2014). Integrating the roles of long and small non-coding RNA in brain function and disease. Mol. Psychiatry 19, 410–416. doi: 10.1038/mp.2013.196

PubMed Abstract | CrossRef Full Text | Google Scholar

Cattane, N., Minelli, A., Milanesi, E., Maj, C., Bignotti, S., Bortolomasi, M., et al. (2015). Altered gene expression in schizophrenia: findings from transcriptional signatures in fibroblasts and blood. PLoS One 10:e0116686. doi: 10.1371/journal.pone.0116686

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, F., Lin, X., Xu, P., Zhang, Z., Chen, Y., Wang, C., et al. (2015). Nuclear export of smads by RanBP3L regulates bone morphogenetic protein signaling and mesenchymal stem cell differentiation. Mol. Cell Biol. 35, 1700–1711. doi: 10.1128/MCB.00121-15

PubMed Abstract | CrossRef Full Text | Google Scholar

Chopra, K., Baveja, A., and Kuhad, A. (2015). MMPs: a novel drug target for schizophrenia. Expert. Opin. Ther. Targets 19, 77–85. doi: 10.1517/14728222.2014.957672

PubMed Abstract | CrossRef Full Text | Google Scholar

de Bartolomeis, A., Latte, G., Tomasetti, C., and Iasevoli, F. (2014). Glutamatergic postsynaptic density protein dysfunctions in synaptic plasticity and dendritic spines morphology: relevance to schizophrenia and other behavioral disorders pathophysiology, and implications for novel therapeutic approaches. Mol. Neurobiol. 49, 484–511. doi: 10.1007/s12035-013-8534-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Gallo, F. T., Katche, C., Morici, J. F., Medina, J. H., and Weisstaub, N. V. (2018). Immediate early genes, memory and psychiatric disorders: focus on c-Fos, Egr1 and Arc. Front. Behav. Neurosci. 12:79. doi: 10.3389/fnbeh.2018.00079

PubMed Abstract | CrossRef Full Text | Google Scholar

Gautier, L., Cope, L., Bolstad, B. M., and Irizarry, R. A. (2004). Affy analysis of Affymetrix GeneChip data at the probe level. Bioinformatics 20, 307–315. doi: 10.1093/bioinformatics/btg405

PubMed Abstract | CrossRef Full Text | Google Scholar

Goff, D. C., Falkai, P., Fleischhacker, W. W., Girgis, R. R., Kahn, R. M., Uchida, H., et al. (2017). The long-term effects of antipsychotic medication on clinical course in schizophrenia. Am. J. Psychiatry 174, 840–849. doi: 10.1176/appi.ajp.2017.16091016

PubMed Abstract | CrossRef Full Text | Google Scholar

Gong, J., Li, Y., Liu, C. J., Xiang, Y., Li, C., Ye, Y., et al. (2017). A pan-cancer analysis of the expression and clinical relevance of small nucleolar RNAs in human cancer. Cell Rep. 21, 1968–1981. doi: 10.1016/j.celrep.2017.10.070

PubMed Abstract | CrossRef Full Text | Google Scholar

Gu, Z., Gu, L., Eils, R., Schlesner, M., and Brors, B. (2014). circlize implements and enhances circular visualization in R. Bioinformatics 30, 2811–2812. doi: 10.1093/bioinformatics/btu393

PubMed Abstract | CrossRef Full Text | Google Scholar

Hagenauer, M. H., Schulmann, A., Li, J. Z., Vawter, M. P., Walsh, D. M., Thompson, R. C., et al. (2018). Inference of cell type content from human brain transcriptomic datasets illuminates the effects of age, manner of death, dissection, and psychiatric diagnosis. PLoS ONE 13:e0200003. doi: 10.1371/journal.pone.0200003

PubMed Abstract | CrossRef Full Text | Google Scholar

Hall, J., Trent, S., Thomas, K. L., O’Donovan, M. C., and Owen, M. J. (2015). Genetic risk for schizophrenia: convergence on synaptic pathways involved in plasticity. Biol. Psychiatry 77, 52–58. doi: 10.1016/j.biopsych.2014.07.011

PubMed Abstract | CrossRef Full Text | Google Scholar

Heyes, S., Pratt, W. S., Rees, E., Dahimene, S., Ferron, L., Owen, M. J., et al. (2015). Genetic disruption of voltage-gated calcium channels in psychiatric and neurological disorders. Prog. Neurobiol. 134, 36–54. doi: 10.1016/j.pneurobio.2015.09.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Horiuchi, Y., Kondo, M. A., Okada, K., Takayanagi, Y., Tanaka, T., Ho, T., et al. (2016). Molecular signatures associated with cognitive deficits in schizophrenia: a study of biopsied olfactory neural epithelium. Transl. Psychiatry 6:e915. doi: 10.1038/tp.2016.154

PubMed Abstract | CrossRef Full Text | Google Scholar

Howes, O. D., and Murray, R. M. (2014). Schizophrenia: an integrated sociodevelopmental-cognitive model. Lancet 383, 1677–1687. doi: 10.1016/S0140-6736(13)62036-X

PubMed Abstract | CrossRef Full Text | Google Scholar

Karam, C. S., Ballon, J. S., Bivens, N. M., Freyberg, Z., Girgis, R. R., Lizardi-Ortiz, J. E., et al. (2010). Signaling pathways in schizophrenia: emerging targets and therapeutic strategies. Trends Pharmacol. Sci. 31, 381–390. doi: 10.1016/j.tips.2010.05.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Kumarasinghe, N., Tooney, P. A., and Schall, U. (2012). Finding the needle in the haystack: a review of microarray gene expression research into schizophrenia. Aust. N. Z. J. Psychiatry 46, 598–610. doi: 10.1177/0004867412442405

PubMed Abstract | CrossRef Full Text | Google Scholar

Lai, C. Y., Lee, S. Y., Scarr, E., Yu, Y. H., Lin, Y. T., Liu, C. M., et al. (2016). Aberrant expression of microRNAs as biomarker for schizophrenia: from acute state to partial remission, and from peripheral blood to cortical tissue. Transl. Psychiatry 6:e717. doi: 10.1038/tp.2015.213

PubMed Abstract | CrossRef Full Text | Google Scholar

Lewis, D. A., Curley, A. A., Glausier, J. R., and Volk, D. W. (2012). Cortical parvalbumin interneurons and cognitive dysfunction in schizophrenia. Trends Neurosci. 35, 57–67. doi: 10.1016/j.tins.2011.10.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Manzardo, A. M., Gunewardena, S., Wang, K., and Butler, M. G. (2014). Exon microarray analysis of human dorsolateral prefrontal cortex in alcoholism. Alcohol. Clin. Exp. Res. 38, 1594–1601. doi: 10.1111/acer.12429

PubMed Abstract | CrossRef Full Text | Google Scholar

Matera, I., Rusmini, M., Guo, Y., Lerone, M., Li, J., Zhang, J., et al. (2016). Variants of the ACTG2 gene correlate with degree of severity and presence of megacystis in chronic intestinal pseudo-obstruction. Eur. J. Hum. Genet. 24, 1211–1215. doi: 10.1038/ejhg.2015.275

PubMed Abstract | CrossRef Full Text | Google Scholar

Mattick, J. S. (2004). RNA regulation: a new genetics? Nat. Rev. Genet. 5, 316–323. doi: 10.1038/nrg1321

PubMed Abstract | CrossRef Full Text | Google Scholar

McGrath, J., Saha, S., Chant, D., and Welham, J. (2008). Schizophrenia: a concise overview of incidence, prevalence, and mortality. Epidemiol. Rev. 30, 67–76. doi: 10.1093/epirev/mxn001

PubMed Abstract | CrossRef Full Text | Google Scholar

Mehler, M. F., and Mattick, J. S. (2006). Non-coding RNAs in the nervous system. J. Physiol. 575, 333–341. doi: 10.1113/jphysiol.2006.113191

CrossRef Full Text | Google Scholar

Milde-Langosch, K. (2005). The Fos family of transcription factors and their role in tumourigenesis. Eur. J. Cancer 41, 2449–2461. doi: 10.1016/j.ejca.2005.08.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Mistry, M., Gillis, J., and Pavlidis, P. (2013). Genome-wide expression profiling of schizophrenia using a large combined cohort. Mol. Psychiatry 18, 215–225. doi: 10.1038/mp.2011.172

PubMed Abstract | CrossRef Full Text | Google Scholar

Mistry, M., and Pavlidis, P. (2010). A cross-laboratory comparison of expression profiling data from normal human postmortem brain. Neuroscience 167, 384–395. doi: 10.1016/j.neuroscience.2010.01.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Monfil, T., Vazquez Roque, R. A., Camacho-Abrego, I., Tendilla-Beltran, H., Iannitti, T., Meneses-Morales, I., et al. (2018). Hyper-response to novelty increases c-Fos expression in the hippocampus and prefrontal cortex in a rat model of Schizophrenia. Neurochem. Res. 43, 441–448. doi: 10.1007/s11064-017-2439-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Mouchlianitis, E., McCutcheon, R., and Howes, O. D. (2016). Brain-imaging studies of treatment-resistant schizophrenia: a systematic review. Lancet Psychiatry 3, 451–463. doi: 10.1016/S2215-0366(15)00540-4

CrossRef Full Text | Google Scholar

Nestler, E. J. (2008). Review. Transcriptional mechanisms of addiction: role of DeltaFosB. Philos. Trans. R. Soc. Lond. B Biol. Sci. 363, 3245–3255. doi: 10.1098/rstb.2008.0067

PubMed Abstract | CrossRef Full Text | Google Scholar

Network, T., O’Dushlaine, C., Rossin, L., Lee, P. H., Duncan, L., Parikshak, N. N., et al. (2015). Psychiatric genome-wide association study analyses implicate neuronal, immune and histone pathways. Nat. Neurosci. 18:199. doi: 10.1038/nn.3922

PubMed Abstract | CrossRef Full Text | Google Scholar

Pang, K. C., Frith, M. C., and Mattick, J. S. (2006). Rapid evolution of noncoding RNAs: lack of conservation does not mean lack of function. Trends Genet. 22, 1–5. doi: 10.1016/j.tig.2005.10.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Park, G., Gim, J., Kim, A. R., Han, K. H., Kim, H. S., Oh, S. H., et al. (2013). Multiphasic analysis of whole exome sequencing data identifies a novel mutation of ACTG1 in a nonsyndromic hearing loss family. BMC Genomics 14:191. doi: 10.1186/1471-2164-14-191

PubMed Abstract | CrossRef Full Text | Google Scholar

Pedrazzi, J. F., Issy, A. C., Gomes, F. V., Guimaraes, F. S., and Del-Bel, E. A. (2015). Cannabidiol effects in the prepulse inhibition disruption induced by amphetamine. Psychopharmacology 232, 3057–3065. doi: 10.1007/s00213-015-3945-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Perrin, B. J., Sonnemann, K. J., and Ervasti, J. M. (2010). beta-actin and gamma-actin are each dispensable for auditory hair cell development but required for Stereocilia maintenance. PLoS Genet. 6:e1001158. doi: 10.1371/journal.pgen.1001158

PubMed Abstract | CrossRef Full Text | Google Scholar

Poola, I., DeWitty, R. L., Marshalleck, J. J., Bhatnagar, R., Abraham, J., and Leffall, L. D. (2005). Identification of MMP-1 as a putative breast cancer predictive marker by global gene expression analysis. Nat. Med. 11, 481–483. doi: 10.1038/nm1243

PubMed Abstract | CrossRef Full Text | Google Scholar

Quadrato, G., Brown, J., and Arlotta, P. (2016). The promises and challenges of human brain organoids as models of neuropsychiatric disease. Nat. Med. 22, 1220–1228. doi: 10.1038/nm.4214

PubMed Abstract | CrossRef Full Text | Google Scholar

Renard, J., Loureiro, M., Rosen, L. G., Zunder, J., de Oliveira, C., Schmid, S., et al. (2016). Cannabidiol counteracts amphetamine-induced neuronal and behavioral sensitization of the mesolimbic dopamine pathway through a Novel mTOR/p70S6 kinase signaling pathway. J. Neurosci. 36, 5160–5169. doi: 10.1523/JNEUROSCI.3387-15.2016

PubMed Abstract | CrossRef Full Text | Google Scholar

Robin, X., Turck, N., Hainard, A., Tiberti, N., Lisacek, F., Sanchez, J.-C., et al. (2011). pROC: an open-source package for R and S+ to analyze and compare ROC curves. BMC Bioinformatics 12:77. doi: 10.1186/1471-2105-12-77

PubMed Abstract | CrossRef Full Text | Google Scholar

Rollins, B., Martin, M. V., Morgan, L., and Vawter, M. P. (2010). Analysis of whole genome biomarker expression in blood and brain. Am. J. Med. Genet. B. Neuropsychiatr. Genet. 153B, 919–936. doi: 10.1002/ajmg.b.31062

PubMed Abstract | CrossRef Full Text | Google Scholar

Ruffle, J. K. (2014). Molecular neurobiology of addiction: what’s all the (Delta)FosB about? Am. J. Drug Alcohol. Abuse 40, 428–437. doi: 10.3109/00952990.2014.933840

PubMed Abstract | CrossRef Full Text | Google Scholar

Stark, K. L., Xu, B., Bagchi, A., Lai, W. S., Liu, H., Hsu, R., et al. (2008). Altered brain microRNA biogenesis contributes to phenotypic deficits in a 22q11-deletion mouse model. Nat. Genet. 40, 751–760. doi: 10.1038/ng.138

PubMed Abstract | CrossRef Full Text | Google Scholar

Sullivan, P. F., Fan, C., and Perou, C. M. (2006). Evaluating the comparability of gene expression in blood and brain. Am. J. Med. Genet. B Neuropsychiatr. Genet. 141B, 261–268. doi: 10.1002/ajmg.b.30272

PubMed Abstract | CrossRef Full Text | Google Scholar

Thorson, W., Diaz-Horta, O., Foster, J. II, Spiliopoulos, M., Quintero, R., Farooq, A., et al. (2014). De novo ACTG2 mutations cause congenital distended bladder, microcolon, and intestinal hypoperistalsis. Hum. Genet. 133, 737–742. doi: 10.1007/s00439-013-1406-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Tylee, D. S., Kawaguchi, D. M., and Glatt, S. J. (2013). On the outside, looking in: a review and evaluation of the comparability of blood and brain “-omes”. Am. J. Med. Genet. B Neuropsychiatr. Genet. 162B, 595–603. doi: 10.1002/ajmg.b.32150

PubMed Abstract | CrossRef Full Text | Google Scholar

Uhlirova, M., and Bohmann, D. (2006). JNK- and Fos-regulated Mmp1 expression cooperates with Ras to induce invasive tumors in Drosophila. EMBO J. 25, 5294–5304. doi: 10.1038/sj.emboj.7601401

PubMed Abstract | CrossRef Full Text | Google Scholar

Van Winkel, R., Esquivel, G., Kenis, G., Wichers, M., Collip, D., Peerbooms, O., et al. (2010). REVIEW: Genome-wide findings in schizophrenia and the role of gene-environment interplay. CNS Neurosci. Ther. 16, e185–e192. doi: 10.1111/j.1755-5949.2010.00155.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Wangler, M. F., Gonzaga-Jauregui, C., Gambin, T., Penney, S., Moss, T., Chopra, A., et al. (2014). Heterozygous de novo and inherited mutations in the smooth muscle actin (ACTG2) gene underlie megacystis-microcolon-intestinal hypoperistalsis syndrome. PLoS Genet. 10:e1004258. doi: 10.1371/journal.pgen.1004258

PubMed Abstract | CrossRef Full Text | Google Scholar

Williams, G. T., and Farzaneh, F. (2012). Are snoRNAs and snoRNA host genes new players in cancer? Nat. Rev. Cancer 12, 84–88. doi: 10.1038/nrc3195

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, F., Liu, Y., Rohde, C., Pauli, C., Gerloff, D., Kohn, M., et al. (2017). AML1-ETO requires enhanced C/D box snoRNA/RNP formation to induce self-renewal and leukaemia. Nat. Cell Biol. 19, 844–855. doi: 10.1038/ncb3563

PubMed Abstract | CrossRef Full Text | Google Scholar

Zuo, D. Y., Cao, Y., Zhang, L., Wang, H. F., and Wu, Y. L. (2009). Effects of acute and chronic administration of MK-801 on c-Fos protein expression in mice brain regions implicated in schizophrenia with or without clozapine. Prog. Neuropsychopharmacol. Biol. Psychiatry 33, 290–295. doi: 10.1016/j.pnpbp.2008.12.002

PubMed Abstract | CrossRef Full Text | Google Scholar