Analyzing the genes and pathways related to major depressive disorder via a systems biology approach

Abstract Introduction Major depressive disorder (MDD) is a mental disorder caused by the combination of genetic, environmental, and psychological factors. Over the years, a number of genes potentially associated with MDD have been identified. However, in many cases, the role of these genes and their relationship in the etiology and development of MDD remains unclear. Under such situation, a systems biology approach focusing on the function correlation and interaction of the candidate genes in the context of MDD will provide useful information on exploring the molecular mechanisms underlying the disease. Methods We collected genes potentially related to MDD by screening the human genetic studies deposited in PubMed (https://www.ncbi.nlm.nih.gov/pubmed). The main biological themes within the genes were explored by function and pathway enrichment analysis. Then, the interaction of genes was analyzed in the context of protein–protein interaction network and a MDD‐specific network was built by Steiner minimal tree algorithm. Results We collected 255 candidate genes reported to be associated with MDD from available publications. Functional analysis revealed that biological processes and biochemical pathways related to neuronal development, endocrine, cell growth and/or survivals, and immunology were enriched in these genes. The pathways could be largely grouped into three modules involved in biological procedures related to nervous system, the immune system, and the endocrine system, respectively. From the MDD‐specific network, 35 novel genes potentially associated with the disease were identified. Conclusion By means of network‐ and pathway‐based methods, we explored the molecular mechanism underlying the pathogenesis of MDD at a systems biology level. Results from our work could provide valuable clues for understanding the molecular features of MDD.

Till now, the cause of MDD is still poorly understood although much effort has been dedicated to explore the pathogenesis and molecular mechanisms of the disease via various approaches (CONVERGE consortium, 2015;Flint & Kendler, 2014;Kang et al., 2012;Mehta, Menke, Menke, & Binder, 2010). Physiologically, MDD is featured with symptom heterogeneity and changes in multiple biological systems are involved (Belmaker & Agam, 2008;Guo et al., 2012). Generally, MDD develops as a result of the combination of multiple factors, including the genetic factors, environmental, and psychological factors (Han, 2012). Actually, a large fraction of the risk of MDD can be attributed to genetics (American Psychiatric Association, 2013; Kendler et al., 2019;Ripke et al., 2013). For example, it is estimated that heritability for MDD is about 40% and the risk of developing depression for members from a family with depression history is 1.5-3 times higher than the normal population (Kendler, Gatz, Gatz, Gardner, & Pedersen, 2006;Pincus et al., 1999). As a polygenic disorder with divergent genetic architecture, many genetic factors, as well as gene-environment interactions, are believed to be among the risk factors of MDD (CONVERGE consortium, 2015;Ripke et al., 2013). A number of genes have been suggested to be associated with MDD, for example, the sodium-dependent serotonin transporter and solute carrier family 6 member 4 (SLC6A4), 5-hydroxytryptamine receptor 2A (5HT2A), apolipoprotein E (APOE), and brain-derived neurotrophic factor (BDNF; Bosker et al., 2011;Flint & Kendler, 2014;Lopez-Leon et al., 2008). Among them, SLC6A4 is one of the most extensively studied genes, which is responsible for transporting serotonin from the synaptic spaces into the presynaptic neurons and recycling it in a sodium-dependent manner. The 5-HTTLPR polymorphism of this gene is found to be associated with both depression and other mental disorders (Clarke, Flint, Flint, Attwood, & Munafò, 2010). As the main excitatory receptor of serotonin, the genetic variants of 5HT2A have been found to be related to several psychiatric disorders, including depression (Choi et al., 2004). The epsilon-4 type allele of APOE is found to be associated with depression in patients with Alzheimer's disease (Delano-Wood et al., 2008). BDNF is involved in activity-dependent neuronal plasticity, and evidence from clinical studies shows that decreased activity of BDNF occurs in the brain of patients with major depression (Lee & Kim, 2010). Similar to other complex mental disorders, genetic studies have suggested that for MDD, the individual differences may be caused by multiple genes and their variants. Genes with different functions may work cooperatively to increase the risk of MDD, with a relatively small effect exerted by each gene. In line with this view, more and more genes have been found to be potentially associated with MDD (Wray et al., 2018). For these genes, although a few plausible candidate genes have been partially replicated, some of them are considered to be problematic (Flint & Kendler, 2014). This is especially true as high-throughput methods like genome-wide association study (GWAS) are increasingly applied to genetic studies of the disease.
Under such circumstances, a comprehensive analysis of the potential causal genes of MDD within a pathway and/or a network framework may not only provide us important insights beyond the conventional single-gene analyses, but also offer consolidated validation for the individual candidate genes.
In the current study, we first collected the MDD-related genes from genetic association studies. Then, we conducted biological enrichment analyses to detect the significant biological themes within these genetic factors and investigated the interactions among the enriched biochemical pathways. In addition, a MDD-related subnetwork based on protein-protein interaction network was constructed and its topological characteristics were analyzed. This study could offer valuable hints for understanding the molecular mechanisms of MDD from a perspective of systems biology.

| Susceptibility gene set of MDD
As a polygenic disease, a number of genes potentially associated with the pathogenesis of MDD have been reported (Gatt, Burton, Burton, Williams, & Schofield, 2015;Manoharan, Shewade, Shewade, Rajkumar, & Adithan, 2016;Yin et al., 2016). In this study, the candidate genes for MDD were collected by searching the human genetic association studies deposited in PubMed (https ://www.ncbi.nlm.nih.gov/ pubme d/). Briefly, similar to previous studies (Wang & Li, 2010) NOT (Neoplasms [MeSH])." As of August 2017, we obtained a total of 1,514 publications related to MDD. Next, we reviewed the abstracts of these articles and kept only the association studies related to MDD with human subjects. From the selected publications, we narrowed our selection by focusing on those reporting a significant association of one or more genes with the disease. To reduce the number of potential false-positive findings, the studies reporting negative or insignificant associations were not included although some genes analyzed in these studies might be real pathogenic genes of MDD. Then, the full reports of the selected publications were examined to ensure the consistency of the conclusions and the contents. In the collected publications, several genome-wide association (GWA) studies on MDD were included, and genes reported to be significantly associated with MDD were selected. Via such a procedure, a list of 261 studies reporting the association of one or more candidate genes with MDD were obtained ( Figure 1). From these studies, genes reported to be associated with MDD were compiled for further analysis.

| Functional enrichment analysis
To reveal the major biological themes within the MDD-related genes, the function characteristics of these genes were explored. Briefly, gene ontology (GO; Fu et al., 2015) and pathway enrichment analysis were conducted on the MDD-related genes. Since in this study, we focused on the biological features underlying the candidate genes, only the GO category of biological process was analyzed. Biological pathways enriched in the MDD-related genes may be those with disturbed function in the pathogenesis of MDD. Both GO and pathway enrichment analysis were finished by the ToppFun module of ToppGene (http://toppg ene.cchmc.org; Chen, Bardes, Bardes, Aronow, & Jegga, 2009). For GO biological process analysis, items with 5 or more MDD-related genes and a false discovery rate (FDR) less than 0.05 were kept as significantly enriched ones. Then, the enriched items were subjected to REVIGO (Supek, Bosnjak, Bosnjak, Škunca, & Šmuc, 2011;http://revigo.irb.hr/) to remove the redundant GO terms and obtained a list of nonredundant GO biological process terms enriched in the candidate genes. For pathway analysis, the Kyoto Gene and Genome Encyclopedia (Du et al., 2016;KEGG) F I G U R E 1 PRISMA flow diagram illustrating search strategy and studies included in the analysis. PRISMA is Preferred Reporting Items for Systematic Reviews and Meta-Analyses (http://www.prisma-state ment.org/) PATHWAY was adopted as the pathway database, and a FDR threshold of 0.05 was used to define a significant pathway. False discovery rate was calculated via the method of Benjamini and Hochberg (Benjamini & Hochberg, 1995).

| Pathway cross talk analyses
The etiology and development of a complex disease are usually the result of simultaneous disturbance of multiple biological processes or pathways. Therefore, the relationship between the pathogenically abnormal pathways can provide useful clues to understand the molecular mechanisms of the disease. Through analyzing the network formed by correlated pathways, we are able to explore the biological pathways summarized from many different studies via a systematic approach, which may help us to understand the etiology and progression of a disease from a macro perspective. Here, we used the pathways enriched in the MDD-related genes to construct the pathway cross talk network, in which two pathways were defined as connected if they shared three or more overlapping MDD-related genes. The purpose of such definition was to reduce the false positives and ensure that the correlation between a pathway pair was biologically meaningful. To describe the overlap between a given pair of pathways, we adopted two measurements (Jia, Kao, Kao, Kuo, & Zhao, 2011;Liu, Fan, Fan, Liu, Cheng, & Wang, 2015), that is, the and the Overlap Coefficient= |A∩B| min(|A|,|B|) , with A and B being the lists of MDD-related genes included in the two tested pathways, and |A| and |B| representing the number of MDD-related genes contained in the two pathways. In addition, we used the arithmetic mean of these two coefficients to measure the significance of pathway correlation and arranged all pairs of pathway in descending order of the significance. Then, Cytoscape (Shannon. et al., 2003) was used to output a diagraphic representation of the cross talk relationship between the pathways.

| The construction of MDD subnetwork
Biomolecular network, especially the protein-protein interaction network, has become an effective tool to analyze the molecular relationship in complicated biomolecular systems (Li, Wang, Wang, Zhao, Wu, & Pan, 2016;Przulj, Wigle, Wigle, & Jurisica, 2004). In this study, we treated the genes/proteins and their interactions as nodes and edges, respectively; then, these nodes and edges were connected to form a molecular network. The protein-protein interaction network data used in this study were derived from direct physical interactions from six major common protein-protein interaction databases, that is, BioGM, Integrity, DIP, Peppermint, MIPS/ Mpact, and HPRD, with the self-interaction and redundant pairs excluded. Finally, a relatively complete human physical interaction group was obtained, which included 16,022 genes/protein and 228,122 interactions.

| MDD candidate gene sets
Based on the human genetic association studies, we compiled a list of 255 candidate genes reported to be associated with MDD (Table S1; referred to as MDDgene, hereafter). Among the candidate genes collected, there were some overlapping genes that were not only associated with MDD, but also involved in the occurrence and development of other neurological diseases. For example, some genes related to immune regulation and inflammation may be associated with Alzheimer's disease or depression (e.g., IL10 and IL1B), genes of the dopamine neurotransmitter system (e.g., DRD1 and DRD4), and members from the immunophilin protein family (e.g., FKBP4 and FKBP5) that may be associated with Alzheimer's disease or depressive disorders. In addition, there were also genes related to the serotonin neurotransmitter system and cell transport system, such as HTR2A, HTR6, TPH1, SLC1A2, SLC6A3, and SLC6A4. At the same time, the gene set included some specific genes related to MDD, such as ADCY9, ITPR1, and PCLO, which were involved in calcium signaling, binding, and salivary secretion biological pathways. Genes related to embryonic development (e.g., CHST11 and PTPRR), cellular stress response, and blood clotting (e.g., DNAJB2, EHD3) were also included. The diversity of MDDgene was consistent with the fact that MDD was a multigene and complex disease involving various physiological procedures.

| Functional enrichment analysis of MDDgene
Functional enrichment analysis revealed a more detailed biological function spectrum of these MDD-related genes (Table S2).
Among the GO terms overrepresented in MDDgene, those related to cell signaling, synaptic transmission, cell transport, endocrine system, or response to stimuli were included. GO terms associated with response to stimuli (e.g., multicellular organismal response to stress, response to wounding, response to light stimulus, and response to pain) were overrepresented. Such results were in line with previous findings that complicated correlations existed between the pathophysiological state of MDD and stress. Biological process terms related to synaptic transmission (e.g., trans-synaptic signaling; synaptic signaling; neuron-neuron synaptic transmission; positive regulation of synaptic transmission; synaptic transmission, glutamatergic; and synaptic transmission, GABAergic), dopamine signaling (dopamine transport, dopamine metabolic process, dopamine uptake involved in synaptic transmission, and dopamine uptake), and other neural functions (e.g., regulation of synaptic plasticity, long-term synaptic potentiation, neuron apoptotic process, and memory) were also enriched. Meanwhile, GO terms related to endocrine system (e.g., hormone secretion, insulin secretion, response to insulin, and response to hormone) were overrepresented. These results demonstrated that the members of MDDgene were diverse in molecular functions.

| Pathways enriched in MDD candidate genes
Pathway analysis identified 73 pathways with significant enrichment in MDDgene (Table 1). Several pathways related to neurotransmission or neural function modulation were identified, for example, neuroactive ligand-receptor interaction, glutamatergic synapse, serotonergic synapse, dopaminergic synapse, GABAergic synapse, cholinergic synapse, and retrograde endocannabinoid signaling. A number of pathways involved in cellular signaling cascade were enriched, for example, cAMP signaling pathway, MAPK signaling pathway, and calcium signaling pathway. In addition, pathways related to neurological disorders, such as morphine addiction, amphetamine addiction, Alzheimer's disease, and alcoholism, were significantly enriched. Moreover, immune responseassociated biological processes consisting of inflammatory bowel disease, inflammatory mediator regulation of TRP channels, interleukin-17 (IL-17) signaling pathway, and T-cell receptor signaling pathway were also significantly enriched, suggesting the immunological system was involved in the etiology and pathological process of MDD.
We further analyzed the cross talk between the enriched pathways that were significantly associated with MDD. Most of these pathways interacted with one or more other pathways, which resulted in a cross talk network with 68 nodes (i.e., pathways) and 325 edges (i.e., connection between two neighboring pathways; Figure 2). Based on the biological function and the relevance of these pathways, we could roughly divide the pathways into three modules. Pathways in the first module were mainly related to cellular signaling transduction (e.g., cAMP signaling pathway, calcium signaling pathway, cGMP-PKG signaling pathway, and phospholipase D signaling pathway) or the endocrine control (e.g., renin secretion, aldosterone synthesis and secretion, oxytocin signaling pathway, thyroid hormone synthesis, and estrogen signaling pathway). In the second module, many pathways were related to neuronal function like neurotransmission (e.g., cholinergic synapse, dopaminergic synapse, GABAergic synapse, glutamatergic synapse, and long-term depression), neurological disorders (e.g., amphetamine addiction, cocaine addiction, morphine addiction, nicotine addiction, alcoholism, amyotrophic lateral sclerosis, and Alzheimer's disease), endocrine, and metabolic diseases (e.g., type II diabetes mellitus and insulin resistance). The last module was largely concentrated in pathways related to the immune system, such as cytosolic DNA-sensing pathway, IL-17 signaling pathway, NOD-like receptor signaling pathway, T-cell receptor signaling pathway and Th17 cell differentiation, and Toll-like receptor signaling pathway. These three modules were not independent of each other; instead, they were interconnected by one or more pathways. In this cross talk network, a few other types of pathways related to biological processes such as aging, apoptosis, and environmental adaptation were also included. Thus, the etiology and development of MDD could be the consequence of the abnormality in multiple systems.

| MDD-specific network
To further explore the feature of genes associated with MDD, we constructed a subnetwork for the disease from the human protein-protein interaction network via the Steiner minimal tree algorithm (Li, Mao, Mao, & Wei, 2008;Sadeghi & Fröhlich, 2013), which tried to connect the largest number of input nodes (genes included in MDDgene in our case) via the least number of interlinking nodes At the same time, 35 genes outside of MDDgene were introduced into the MDD-specific molecular network (Table 2). Given these genes interacted closely with those known to be related to MDD, they might also be involved in the pathogenesis of the disease phenotype. Further functional enrichment analysis indicated that these genes were mainly involved in neuronal development, behavior, learning and memory, and glutamate receptor signaling.

| D ISCUSS I ON
Recent years, our understanding on the molecular mechanisms of MDD has been greatly improved. With the advancement and maturity of high-throughput technology, we are able to identify the elements related to this disease on much larger scales. Although more and more genes/proteins potentially involved in the disease have been reported, a thorough analysis of the biochemical processes associated with the pathogenesis of MDD from the molecular aspect is still missing. In such case, a systematic analysis of MDD-related genes via a pathway-and network-based analytical framework will provide us insight on the disease beyond the single candidate genebased analyses. In this study, we tried to pool and curate the genes related to MDD from human genetic studies, and systematically delineated the interconnection of these genes based on pathway and network analysis.   GABRB3, GABRD, GABRG2, AVPR1B, GHRHR, CNR1, VIPR2, DRD1,  HTR1A, DRD4, HTR1B, HTR2A, HTR2C, HTR4, HTR6, GRIA1, GRIA2,  GRIA4, GRIK1, GRIK4, GRIN2A, GRIN2B, NR3C1, GABBR2, GRM7,  GRM8, CRHR1, CRHR2, OPRM1,   can not only give us a more comprehensive view for the pathological mechanisms of MDD, but they also are more robust to the influence of false-positive genes.
As revealed by function enrichment analysis, genes related to MDD were diverse in function, mainly involved in cell signaling, immune system, metabolic process, drug response processes, and neurodevelopment. Gene ontology biological process terms such as reverse cholesterol transport, positive regulation of IL-6 production, response to ethanol, lipoprotein metabolic process, diol metabolic process, xenobiotic metabolic process, and regulation of neuronal synaptic plasticity were overrepresented among MDDgene, implying the important roles of these processes in the pathological processes of MDD. In addition, we noticed terms related to memory, visual learning, social behavior, sleep, axon regeneration, and axon guidance were also enriched in MDDgene, consistent with a priori biological findings on MDD.
Biological pathways enriched in MDDgene were involved in multiple biological systems, including the nervous system, immune system, endocrine systems, and signal transduction systems, or related to disorders like drug addiction and immune metabolism diseases. Actually, abnormality or dysregulation of many of these pathways has been known to be related to neurological diseases.
For example, calcium signaling pathway has been reported to be involved in diseases such as nicotine addiction (Wang & Li, 2010), Alzheimer's disease (Karttunen et al., 2011), bipolar disorder and schizophrenia (Berridge, 2014), and depression (Donev & Alawam, 2015;. Another example is the pathway of GABAergic synapse. As the most abundant inhibitory neurotransmitter in the mammalian central nervous system (Lloyd, Perrault, & Zivkovic, 2017;Zhang et al., 2018), the defect of GABAergic neurons in the frontal cortex may be responsible for the pathogenesis and development of MDD (Czéh et al., 2018). The identification of GABAergic synapse pathway in MDDgene provides

TA B L E 1 (Continued)
additional evidence that GABAergic dysfunction may lead to mood and cognitive symptoms of MDD. Interleukin-17 in the IL-17 signaling pathway plays a crucial role in acute and chronic inflammatory responses (Zhao, Li, Li, Wang, Manthari, & Wang, 2018). Neuroendocrine and immune system interactions play an important role in stress response (Ashley & Demas, 2017;Dantzer, 2018). Such results suggest that the immune system plays important roles in the onset of MDD. Both stress and inflammatory cytokine activation have been reported to have adverse effect on the neurogenesis and neural plasticity (Syed et al., 2018). By comparing our results with that of a meta-analysis on genes implicated in MDD (Gatt et al., 2015;Manoharan et al., 2016;Yin et al., 2016), we found that most of the pathways reported earlier were also identified in the current study. and CaMKII play an important role in long-term potentiation and long-term depression, and they connect multiple pathway genes, suggesting that CaM and CaMKII may play an important role in the development of synaptic plasticity. Perhaps it is the key factor that affects the development of MDD. In addition, the genes CLOCK and BMALL are essential in several pathways related to MDD (e.g., prolactin signaling and circadian rhythm), suggesting they may be involved in the development of MDD. Since these pathways are interconnected and they function cooperatively, dysfunction in one pathway may cause abnormality or dysregulation in others and eventually lead to the onset and development of MDD.
In the pathway cross talk network, there were several pathways related to other diseases, such as pathway of Alzheimer's disease, AGE-RAGE signaling pathway in diabetic complication, pathway of alcoholism, and pathway of dilated cardiomyopathy. Available evidence shows that each of these diseases has close correlation MDD. For example, it has been found that depression is associated with an increased risk of Alzheimer's disease, with MDD patients being 1.5 times more likely to develop Alzheimer's disease and 20% to 50% patients with Alzheimer's disease having depressive symptoms (Gibson et al., 2017;Saczynski. et al., 2010). Comparison of the molecules involved in the two diseases shows that they share a number of genes, regulatory elements like miRNAs, and quite several biological processes and pathways (Hu, Xin, Xin, Hu, Zhang, F I G U R E 2 Cross talk between pathways related to major depressive disorder (MDD). The circular nodes represent pathways significantly enriched in the genes in associated with MDD, and each edge represents the cross talk between the two connected pathways, with the width corresponding to strength of the cross talk (i.e., the average of the Jaccard Coefficient and the Overlap Coefficient). The nodes labeled with numbers represent the following pathways: 1, "vascular smooth muscle contraction"; 2, "dilated cardiomyopathy"; 3, "estrogen signaling pathway"; 4, "gap junction"; 5, "inflammatory mediator regulation of TRP channels"; 6, "long-term potentiation"; 7, "longevity regulated pathway-multiple species"; 8, "Rap1 signaling pathway"; 9, "neuroactive ligand-receptor interaction"; 10, "amyotrophic lateral sclerosis"; 11, "taste transduction"; 12, "insulin resistance"; 13, "apoptosis"; 14, "AGE-RAGE signaling pathway in diabetic complications"; and 15, "prolactin signaling pathway" & Wang, 2017;Mendes-Silva et al., 2016), which is consistent with the prior knowledge that depression may be a risk factor for Alzheimer's disease or part of the symptoms of dementia. We also detected pathways related to diabetes (i.e., insulin secretion, insulin resistance, and type II diabetes mellitus). Connection between diabetes and depression has been studied extensively, and there is clear symbiotic relationship between the two diseases (Han, 2012;Lloyd, Pambianco, Pambianco, & Orchard, 2010;Patterson, Khazall, Khazall, MacKay, Anisman, & Abizaid, 2013;Roy & Lloyd, 2012;Semenkovich, Brown, Brown, Svrakic, & Lustman, 2015). A possible explanation is that diabetes may affect the function of brain regions like hippocampus (Semenkovich et al., 2015), the abnormality in which may be involved in the pathogenesis of MDD (Colla et al., 2007;Ho, Sommers, Sommers, & Lucki, 2013).
In the subnetwork constructed by genes related to MDD, six genes outside of the MDDgene, that is, APP (amyloid beta precursor protein), HSP90AB1 (heat-shock protein HSP 90-beta), PRKACA (catalytic subunit α of protein kinase A), GRB2 (growth factor receptor-bound protein 2), PRKCA (protein kinase C alpha), and SP1 (transcription factor Sp1), were localized at the key positions in the subnetwork. Compared with other genes, they interacted with more genes in the network. We further extracted the genes interacting with these six genes to examine their connection with other genes ( Figure 5). In the genetic interaction network centered on these F I G U R E 3 Major depressive disorder specific network. The major depressive disorder (MDD)-specific subnetwork was constructed via the Steiner minimum algorithm, including 203 nodes and 415 edges. The circular nodes represent the known genes related to MDD, while the red triangular nodes represent the genes newly introduced to the subnetwork, which may be genes potentially related to MDD. The edge represents the interaction between genes genes, approximately 78% (54/69) of the genes were members of MDDgene. Functionally, pathways related to the immune system or the nervous system were enriched in these genes, implicating these genes may be involved in MDD through their connection with pathways related to the immune system or the nervous system.
Among these genes, APP encodes the precursor molecule of beta amyloid, the primary component of amyloid plaques found in the brains of patients with Alzheimer's disease. In our pathway enrichment analysis, pathway related to Alzheimer's disease was also enriched in MDDgene. Thus, even though the evidence on the role of APP in the pathogenesis of MDD is still limited, it may be closely related to MDD. Catalytic subunit α of protein kinase A and GRB2 have been reported to be related to MDD. Previous studies using human peripheral and postmortem brain tissue samples have shown that some depressed patient's exhibit reduced PRKACA activity (Kastenhuber. et al., 2017;Pandey et al., 2007). There are also some databases related to the genetic information of MDD, but no dataset specific for MDD. As MDD is a complex There are some limitations in the current study. First, there are some subjective factors in the procedure of MDD candidate gene collection.
For example, the collection of studies on MDD may be not comprehensive enough because of the specific screening conditions we used; in many candidate gene studies, the selection of genes could be biased as they are often chosen based on prior knowledge of the disease itself or related diseases. For such reason, the collected MDD-related genes may include a high fraction of genes also associated with other mental disorder, but we believe that with further improvement, the pathogenic genes for MDD will be supplemented and the dataset will become more and more reliable. Second, although multiple pathway databases are available, we only utilized the KEGG pathway database in pathway enrichment analysis, which might lead to bias in the result. But on the other hand, the definition of pathway may be different in various pathway databases, which means pathways with same or similar names may be not consistent in different databases. To avoid the potential confusion caused by merging multiple databases, we relied on KEGG pathway database for our analysis. Third, the current available human PPIN is still incomplete and may include false-positive data, which may have impact on our results. Although there are some shortcomings in the current study, we believe the results obtained by us should be reliable. Finally, several studies on MDD via GWA meta-analysis have been published recently (Howard et al., 2018(Howard et al., , 2019Wray et al., 2018).
Based on large sample sizes, a number of novel variants potentially associated with MDD have been identified in each study. These studies clearly demonstrated the power of GWA meta-analysis in detecting the genetic factors underlying complex disorders like MDD. However, due to the difficulties in data integration, we did not include the genes reported in these studies in our analysis.

F I G U R E 4
Diagram of the major pathways and genes related to major depressive disorder (MDD). MDD is a complex disease with a number of genes and pathways coordinated and interrelated by multiple systems. The nodes in the rectangle represent the genes involved in each pathway. Small elliptical nodes represent neurotransmitters such as GABA, serotonin, dopamine, and glutamate. The large ellipse represents the main pathway involved in MDD. The dashed line and the solid line represent the indirect and direct relationship between the parts; the line of the arrow or breakpoint indicates the activation and inhibition of the action, respectively

| CON CLUS ION
In this study, we conducted a systematic analysis on genes genetically associated with MDD. Based on the 255 disease-related genes collected, 73 significantly enriched pathways were identified. Pathway cross talk analysis indicated that three major modules were formed by these biological pathways, with each module including pathways related to cellular signaling transduction or the endocrine control, neuronal function or neurological disorders, and the immune system, respectively. Then, the disease-specific subnetwork was constructed and a number of novel genes potentially involved in MDD were identified. When more candidate genes associated with MDD are identified, the procedure outlined in this study should provide more detailed gene interaction and pathological molecular network on MDD. In addition, information on MDD from other sources can also be integrated into the framework used in this study; then, we will be able to obtain a more comprehensive and meaningful understanding on the molecular mechanisms on the pathogenesis of MDD.

ACK N OWLED G M ENTS
This study was supported in part by grants from National

CO N FLI C T O F I NTE R E S T
None declared.

DATA AVA I L A B I L I T Y S TAT E M E N T
Supporting data are provided as additional supporting files.

O RCI D
Ju Wang https://orcid.org/0000-0002-6514-0923 F I G U R E 5 Novel genes potentially related to major depressive disorder (MDD). The six essential genes were mapped to protein-protein interaction network (PPIN) to extract the genes interacting with them. The triangle nodes are non-MDD susceptibility gene, and the red triangle nodes are essential genes. The small circular nodes are susceptibility genes in the MDD gene set