Obsessive-Compulsive Disorder Interactome Profile Analysis: A Perspective From Molecular Mechanism.

INTRODUCTION: Obsessive-Compulsive Disorder (OCD) is one of the complex neuropsychiatric conditions. This disorder disables individuals in many different aspects of their personal and social life. Interactome analysis may provide a better understanding of this disorder's molecular origin and its underlying mechanisms.
METHODS: In this study, the OCD-associated genes were extracted from the literature. The criterion for gene selection was to choose genes with at least one significant report. Furthermore, by applying Cytoscape and its plugins, protein-protein interaction network, and gene ontology of the 31 candidate genes related to OCD from genetic association studies is examined. The cross-validation method was used for network centrality assessment.
RESULTS: A scale-free network, including 1940 nodes and 3269 edges for 31 genes, was constructed. According to the network centrality evaluation, ESR1, TNFα, DRD2, DRD4, HTR1B, HTR2A, and CDH2 showed the highest values and can be considered hub-bottlenecks elements. It is also confirmed by the number of 123 cross-validation tests that the frequency of these essential genes remains unaltered against the initial seed genes' changes with the accuracy of 0.962. Besides, enrichment analysis identified four highlighted biological processes related to the 31 candidate genes. The top biological processes are determined as dopamine transport, learning, memory, and monoamine transport.
CONCLUSION: Among 31 initial genes, 7 were introduced as crucial elements for onset and development in OCD and can be suggested for further investigations. Furthermore, the complex molecular origin of OCD requires high-throughput screening for diagnosis and treatment goals. The findings are a possible valuable source to establish molecular-based diagnostic tools for OCD. Copyright

Seven genes were suggested to have topological features in the Obsessive-Compulsive Disorder (OCD) interaction network.

Dopamine transport, learning, memory, and monoamine transport were the critical contributing biological processes in the OCD network.

ESR1 could play an indenciple role in a specific subtype of OCD for women.

The presence of OCD in different subtypes could be linked to the combination of malfunction of gene groups.

Plain Language Summary

The results of this study suggest that topological information of interactome profile of OCD reported genes from genetic studies; however, it is in its initial phase. Our PPI network analysis can serve as a paradigm for later experiments such as genomic, proteomic, and metabolomic investigations.

Introduction

Obsessive-Compulsive Disorder (OCD) is known to be the fourth most frequent heterogeneous mental condition in the world. The lifetime prevalence of OCD is about 1% to 3% (Tang et al., 2014), and it is typically featured with unwanted thoughts and compulsive behavior (Zamanian-Azodi et al., 2015). The complex nature of OCD is correlated to the interactions of many genetic and environmental factors (Taylor, Asmundson, & Jang, 2016). These interactions can lead to complex clinical demonstration and symptom overlap of OCD with other related disorders, including Tourette syndrome, chronic hair pulling, trichotillomania, and anxiety, hindering proper therapeutic approaches (Mahboubi et al., 2018). Many genes with their polymorphisms have been introduced for OCD etiology (Taylor, 2013). However, the contribution of each gene is known to be partial (Jaffe et al., 2014). This fact implies the importance of examining genetic elements in a whole interaction system rather than singly (Rezaei-Tavirani, Tavirani, & Azodi, 2019).

Any minor changes in these elements may influence the molecular interaction pattern and biological processes necessary for everyday functions. These alterations can promote the onset and development of a disease (Chavez, Weisbrod, Zheng, Eng, & Bruce, 2013). In this regard, protein-protein interaction network analysis can help map the interactome profile of a disease (Zamanian-Azodi et al., 2018). Identifying central genes in a predicted network is valuable for better understanding potential elements of disease underlying mechanisms and accordingly introducing possible therapeutic targets (Mahboubi et al., 2018). As mentioned earlier, many genes are reported to be associated with OCD risk. However, there are inconsistencies in their relations among studies. This might be due to various factors, including the difference in sample size, subtypes, Studied single Nucleotide Polymorphisms (SNPs), and the applied statistical criteria. In this regard, a comprehensive meta-analysis may provide more information. Based on one meta-analysis, there are significant relationships between OCD and some genes from genetic association studies (Taylor, 2013).

Nevertheless, this approach may not be adequate for network analysis since some of the other genes that were clinically showed association with OCD might be missed. For instance, dopaminergic system genes showed a lack of significant associations with OCD in the mentioned meta-analysis, while treatment with dopamine receptor antagonism proved to be effective (Koo, Kim, Roh, & Kim, 2010; Lochner et al., 2016). The heterogeneous basis of the genetic association studies made it difficult for the PPI network dataset design. At least 31 genes have been introduced for OCD onset and development, but there is not an explicit priority between them. Therefore, determination of the crucial genes in OCD by interactome analysis is the aim of this investigation.

Methods

Data visualization in the form of complex networks can be obtained by applying Cytoscape (an open-source network analyzer tool; http://cytoscape.org) (Saito et al., 2012). For analyzing genes from genetic association studies, 31 genes for OCD patients were searched through popular literature motor engines (Google Scholar) and specialized biomedical databank (PubMed). The first step for network analysis was to retrieve UniProt accession codes (Homo sapiens) of the 31 genes from http://www.uniprot.org. Then the codes were used for PPI mapping by Cytoscape v: 3.2.1. For PPI network construction, PSICQUIC (Proteomics Standard Initiative Common QUery InterfaCe) was applied. This algorithm is a universal web service client for importing interactions from public databases (Aranda et al., 2011). The applied public databases in this study were Mentha, Reactome-Fls, and Reactome (Lotia, Montojo, Dong, Bader, & Pico, 2013).

Further analysis was performed based on functional topological parameters (centrality), including node degree and betweenness centrality, using Network Analyzer Application. This software is a versatile and user-friendly tool for the analysis of biological networks. This plugin is well-integrated into Cytoscape and computes a comprehensive list of complex and straightforward topology parameters using resourceful graphic algorithms (Assenov, Ramírez, Schelhorn, Lengauer, & Albrecht, 2008). In the PPI network, the nodes represent genes, while the edges represent physical interactions (Luo et al., 2014). Degree (K) implies the number of connected edges to a consider vertices. Betweenness Centrality (BC) equals the number of shortest paths that pass through each node (Assenov et al., 2008). Genes with large K and high BC are called hub-bottlenecks. For hub-bottleneck validation, test networks were obtained based on omitting (1 to 4) top nodes of the main network.

It should be noted that for deleting more than one node, a combination of one hub-bottleneck gene and random selection from the remained top genes with five times repetition was performed. In the end, the robustness of the backbone network was calculated (Nair, Ghatge, Kakkar, & Shanker, 2014). ClueGO plugin is an up-to-date functional analyzer that enriches a list of genes based on gene ontology (Bindea et al., 2009). Here, GlueGO2.1.7 annotates genes for biological process. ClueGO presents gene ontology ranging from general to specific ones. Biological terms were assessed for the designated 31 OCD genes. Genes were grouped based on a significant functional association between terms and gene-sets. For redundancy reduction, the fusion of the relevant terms with the related genes was applied. The association strength between the terms is calculated by Kappa scoring. The Kappa score level was adjusted to 0.5 (medium) for Biological Process (BP) analysis. In this analysis, the minimum number and percentage of genes per term were adjusted to 3 and 4 as the default options, respectively. Likewise, minimum and maximum levels of ontology were set to 3 and 8 as the default option, respectively. The P value was also set to ≤0.05, which was corrected by the Bonferroni step-down method. The enrichment/depletion test for terms was set to two-sided (enrichment/depletion) based on hypergeometric (Bindea, Galon, & Mlecnik, 2013).

Results

Thirty-one OCD seed genes with their properties were retrieved from the genetic association studies as tabulated in Table 1.

Table 1

A list of thirty-one candidate genes for OCD derived from previously reported genetic studies

UA Code	Gene Name	Protein Name
P21964	COMT (Taylor, 2013)	Catechol-O-methyltransferase
P31645	SLC6A4 (Taylor, 2013)	Solute carrier family 6 (neurotransmitter transporter, serotonin), member 4
P46098	HTR3A (Lennertz et al., 2014)	5-hydroxytryptamine (serotonin) receptor 3A
O95264	HTR3B (Lennertz et al., 2014)	5-hydroxytryptamine (serotonin) receptor 3B
Q8WXA8	HTR3C (Lennertz et al., 2014)	5-hydroxytryptamine receptor 3 subunit C
Q70Z44	HTR3D (Lennertz et al., 2014)	5-hydroxytryptamine (serotonin) receptor 3 family member D
A5X5Y0	HTR3E (Lennertz et al., 2014)	5-hydroxytryptamine receptor 3E
O14490	DLGAP1 (Li et al., 2015)	Discs, large homolog-associated protein 1
O95886	DLGAP3 (SAPAP3) (Boardman et al., 2011)	Disks large-associated protein 3
P35462	DRD3 (Nicolini et al., 1996)	Dopamine receptor D3
P14416	DRD2 (Nicolini et al., 1996)	Dopamine receptor D2
P23560	BDNF (Zai et al., 2015)	Brain-derived neurotrophic factor
P21917	DRD4 (Millet et al., 2003)	Dopamine receptor D4
P03372	ESR1 (Alonso et al., 2011)	Estrogen receptor 1
P01375	TNFA (Taylor, 2013)	Tumor necrosis factor
P43005	SLC1A1 (Stewart et al., 2013)	Solute carrier family 1 (neuronal/epithelial high-affinity glutamate transporter, system Xag), member 1
Q13516	OLIG2 (Taylor, 2013)	Oligodendrocyte transcription factor 2
Q13224	GRIN2B (Arnold et al., 2004)	Glutamate receptor, ionotropic, N-methyl D-aspartate 2B
Q05329	GAD2 (Mas et al., 2014)	Glutamate decarboxylase 2
P28222	HTR1B (Mas et al., 2014)	5-hydroxytryptamine (serotonin) receptor 1B, G protein-coupled
P21397	MAOA (Taylor, 2013)	Monoamine oxidase A
P28223	HTR2A (Taylor, 2013)	5-hydroxytryptamine (serotonin) receptor 2A
Q13002	GRIK2 (Taylor, 2013)	Glutamate receptor ionotropic
Q16653	MOG (Zai et al., 2004)	Myelin-oligodendrocyte glycoprotein
Q96PX8	SLITRK1 (Ozomaro et al., 2013)	SLIT and NTRK-like family, member 1
Q16620	NTRK2 (Taylor, 2013)	BDNF/NT-3 growth factors receptor
Q16288	NTRK3 (Alonso et al., 2008)	Neurotrophic tyrosine kinase, receptor, type 3
P28221	HTR1D (Mundo et al., 2000)	5-hydroxytryptamine receptor 1D
P19022	CDH2 (McGregor et al., 2016)	Neuronal cadherin gene
O75751	EMT (SLC22A3) (Lazar et al., 2008)	Solute carrier family 22 member 3

Asterisk sign indicates 8 significant genes obtained from a meta-analysis (Taylor, 2013); UA refers to the UniProt accession number.

The PPI network of this study is formed by 1940 nodes and connected via 3269 edges. Graph topological properties help gain insight into the global network structure (Assenov et al., 2008). Node distribution based on centrality parameters such as K and BC shows the presence of genes with high centrality values (Figure 1). The top ten central genes are presented in Table 2. The higher the value, the more central the nodes are in the network (Scardoni, Petterlini, & Laudanna, 2009).

Figure 1

Graphical show of topological parameters of the scale-free network with the 31 related genes

A: The K (the number of connected edges) distribution; B: The BC (betweenness centrality) distribution.

Table 2

The top ten central genes in the network are listed

Row	Vertices	K	BC
1	ESR1	899	0.7
2	DRD2	258	0.06
3	HTR1B	258	0.06
4	DRD3	254	0.04
5	DRD4	254	0.05
6	HTR1D	245	0.03
7	TNFA	219	0.18
8	HTR2A	201	0.11
9	CDH2	199	0.16
10	HTR3C	84	0.08

The seven genes with asterisk signs indicate key genes (hub-bottlenecks). The genes’ prioritization is based on a cut-off of node K (the number of connected edges) >100 and BC (betweenness centrality) >0.05. The nodes are descended based on degree.

In Figure 1, K values are significantly heterogeneous. Just a few nodes show well-linked, whereas others possess a small number of connections. This distribution implies the presents of genes with high centrality values. The line indicates the power law. The R2 value is computed on logarithmic values, which is equal to 0.53, and the correlation is 0.84. Genes with a high degree are in the right down region of the plot (out of linear range). B: R=0.002 and correlation=0.38, The BT range is from 0-1. Nodes in the top and middle regions of the plot (out of linear range) have high BC values.

The amounts of degree and betweenness centralities for the introduced key genes are represented in Table 2. Based on Table 2, degree values for the hub-bottleneck nodes (except ESR1) vary between 200 and 250.

A backbone network of the introduced 7 hub-bottleneck nodes with the related topological features was obtained (the data are not shown). Furthermore, to confirm the backbone network, we used the cross-validation method. Accordingly, 123 test networks (18+3×35) were constructed by deleting 1 to 4 nodes from 18 top nodes (largest degree). At first, the deletion of one out of 18 top nodes for each test resulted in the construction of 18 networks. Then, as deleting more than 2 genes would result in many different combinations, the test was done five times by selecting one hub-bottleneck gene and randomly selecting from other top genes (Table 3, Figure 2).

Table 3

The frequency of hub-bottleneck genes and their accuracy in 123 tests

Number of Omitted Genes	Frequency of Hub-bottleneck Genes							Accuracy	Number of Test Networks
	ESR1	DRD2	HTR1B	DRD4	TNF	HTR2A	CDH2
1	18	18	18	18	17	18	17	0.984	18
2	35	35	35	35	30	35	30	0.959	35
3	35	35	35	35	30	35	30	0.959	35
4	34	35	35	35	30	35	30	0.955	35
Overall	122	123	123	123	107	123	107	0.962	123

Figure 2

Accuracy of hub-bottleneck genes in the test network

The x-axis shows the number of omitted genes in the test network, while the y-axis represents the accuracy of the corresponding test network.

A backbone network was built from Mentha, Reactome, and Reactome-Fls sources for 7 hub-bottleneck elements. This network has 1609 nodes, and 2279 edges denote that 83% of the nodes relative to the main network are presented in the backbone network.

Different subtypes have been reported for OCD that may appear due to the interaction of different genes and environmental factors. Central genes in this study are reported with different subtypes (Table 4).

Table 4

The list of central genes with the corresponding subtypes

Gene Ontology Term	GO Level	Number of Genes	% Associated Genes	Term P-value Corrected	Associated Genes Found in Term
Dopamine transport	[6, 7]	6	16.2	6.7E-10	DRD2, DRD3, DRD4, HTR1B, HTR2A, SLC22A3
Learning	[4, 6]	5	3.5	1.8E-06	DRD2, DRD3, DRD4, NTRK2, SLC6A4
Memory	[4, 6]	5	4.46	2.1E-05	DRD2, DRD4, HTR2A, NTRK2, SLC6A4
Monoamine transport	[4, 5]	7	9.09	7.4E-10	DRD2, DRD3, DRD4, HTR1B, HTR2A, SLC22A3, SLC6A4

The P-values s are obtained from the previously published reports.

Application of ClueGO plugin for functional enrichment analysis of OCD seed genes based on BP terms is shown in Figure 3 and Table 5.

Figure 3

Overview chart presentation of the main biological processes

Table 5

Detailed information of the most significant terms of highlighted groups resulted by ClueGO analysis

Gene Name	Polymorphism In OCD	Gender	Disease Subtype
ESR1 (Estrogen receptor)	rs34535804, five SNPs haplotype, located at the 5’ end of intron 1 (Alonso et al., 2011)	Female/Male especially in females	Contaminate on obsessions and cleaning compulsions
HTR1B	rs2000292, rs6296,rs6298,rs4140535,rs9359271 (Mas et al., 2014)	Male	Symmetry/Ordering (Ickowicz et al., 2007)
HTR2A	-1438A (Enoch et al., 2001)	Female	Sexually dimorphic
TNFA	rs1800629rs361525A–A haplotype (Hounie et al., 2008)	Female/Male	Tic, Rheumatic Fever (RF)
DRD2	rs1800497 (Nicolini et al., 1996)	Male	Tic
DRD4	Seven-repeat allele (R7) (Viswanath et al., 2013)	Female/Male	Tic
CDH2	rs1120154 and rs12605662 (McGregor et al., 2016)	Female/Male	Repetitive behaviors

The P-value correction method is Bonferroni’s step down; P-value ≤ 0.01 is reported for all data.

Groups are labeled based on the leading term (the most significant). The most significant term is defined as the term with the highest number and percentage of initial genes. Kappa Score ≥0.5. Figure 3 showed that dopamine transport, learning, memory, and monoamine transport are the main BPs.

Discussion

Molecular studies of Obsessive-Compulsive Disorder (OCD) have just begun (Zamanian-Azodi et al., 2015). Most of the recent investigations are focused on the genetic basis of OCD. Multiple genetic variants with small effect size are linked to OCD risk (Gruenblatt, Hauser, & Walitza, 2014). While individual examining of genes is informative, interactome evaluation is essential as well. Genes are in a complex interaction structure that any genetic alterations may lead to functional phenotype changes of the whole system (Chavez et al., 2013). PPI network analysis can identify specific interactome changes in a disease state. In this regard, examining the OCD interactome pattern can help gain further insight into the underlying mechanisms through topological analysis (Mahboubi et al., 2018). In this study, a total number of 31 human candidate genes from genetic association studies are reviewed. Among them, 8 genes were previously reported significant for OCD by a meta-analysis (Taylor, 2013).

On the other hand, none of the single candidate genes has been proved to be significantly associated with OCD by genome-wide association studies (Gruenblatt et al., 2014). Moreover, mRNA expression profiling of postmortem brain of OCD patients showed no significant associations of these genes and mRNA expression profile (Jaffe et al., 2014). Therefore, a network from reported genes may be inaccurate since they are incomplete and may show some bias towards neurotransmitter systems genes (Gruenblatt et al., 2014). To somehow overcome this limitation, it is tried to select genes from genetic-based studies that were with at least one significant association report with OCD.

There are about 1.5 edges per node; however, this value’s variation is considerable for different nodes. For example, ESR1, the top node, is characterized by 899 links. This gene is discriminated by about 650 links from the second top score node. In a PPI network, central genes (hub-bottlenecks) are essential for network integrity and may be considered for more evaluations (Mahboubi et al., 2018). Based on hub and bottleneck definition, when the function of one gene impairs, it may dramatically affect the function of other linked genes. Network analysis in Figure 1 (A & B) reveals that the location of the hub and bottleneck elements are well-defined in the right down region of degree and top region of BC distribution plots, respectively. Centrality analysis implies the significant contribution of ESR1, DRD2, HTR1B, DRD4, TNFA, HTR2A, and CDH2 in the OCD network. In the next step, a backbone network of central genes revealed 7 essential genes whose pivotal role in network integrity requires further analysis. About 83% of the nodes of the main network are presented in the backbone network. However, about 70% of the edges are constructed. As it is tabulated in Table 3, cross-validation confirms the robustness of the network since the frequency of hub-bottleneck genes was not affected by the changes in the initial genes. Therefore, the backbone network shows excellent accuracy for four deletions, and the pivotal role of hub-bottleneck genes is supported. According to Table 3 and Figure 2, the accuracy values vary from 95% to 98%. Besides, detailed characteristic of each vital elements as collected in Table 4, indicates the importance of these genes in OCD subtypes. In this light, ESR1 as the top hub-bottleneck gene in the network may play a major role in OCD. On the other hand, studies showed that ESR1 is associated with other mental disorders as well. ESR1 is seen in the several regions of the brain that are connected to cognitive function and emotional behavior. Steroid hormones regulate compulsive behavior, mood, and cognition (Sundermann, Maki, & Bishop, 2010). Estrogen hormone has a critical regulatory activity on serototonergic, glutamate, and dopaminergic pathways (Alonso et al., 2011). Some kinds of variants in ESR1 can cause changes in estrogen responsivity in the brain (Rettberg, Yao, & Brinton, 2014), resulting in the manifestation of many kinds of neuropsychiatric conditions (Sundermann et al., 2010).

On the other hand, the critical associations of ESR1 polymorphisms in women’s neuropsychiatric disorders are suggested by many studies (Sundermann et al., 2010; Millet et al., 2003). Besides, gender plays a vital role in OCD phenotype (Mathis et al., 2011). Women are more susceptible to certain types of OCD than men. For instance, cleaning obsessions and contamination compulsions are more evident in females that may be connected to ESR1 malfunction (Nissen et al., 2016). It is also suggested by many clinical assessments that there is a vital connection between OCD onset and severity and reproductive cycle activity, in which the level of hormones varies greatly (Vulink, Denys, Bus, & Westenberg, 2006; Labad, Menchón, & Alonso, 2005; Forray, Focseneanu, Pittman, McDougle, & Epperson, 2010; Guglielmi et al., 2014). It seems that ESR1, as part of reproductive cycle events, may play a key role in the pathogenesis of OCD.

Moreover, its central role in the OCD network reflects the possible crucial role of female reproductive cycle activity in OCD onset and development. TNF-α, the second essential gene in the network, has a noticeable function in immunologic reactions that suggests the vital contribution of immunologic reactions in OCD risk. It is also reported in previous investigations that TNF-α has an essential role in some subtypes of OCD (Hounie et al., 2008). DRD2 and DRD4 also implied key roles in the OCD network. It is also confirmed by previous studies that the dopamine system can be an appropriate target for OCD treatment. HTR1B and HTR2A are serotonin receptors, and their relation to OCD is gender-specific. HTR1B has a remarkable contribution in male patients (Mas et al., 2014), while HTR2A is mostly involved in females (Enoch et al., 2001). Another central gene, CDH2, has a primary role in early brain morphogenesis, long-term potentiation, synaptic plasticity, and synaptogenesis. This gene, as a producer of cell-cell adhesion protein, belongs to the cadherin gene superfamily. In OCD, dysregulation of the glutamate receptor availability may attribute to CDH2 polymorphisms (McGregor et al., 2016).

According to the ClueGO findings (Figure 3 and Table 5), dopamine transport, learning, memory, and monoamine transport processes may play an essential role in OCD etiology. These BPs may be influenced by the different contributing polymorphisms. Most genes in these processes are among the ranked genes. The highest-ranked process in this analysis is dopamine transport. On the other hand, most of the prescribed medicines for OCD belong to the Selective Serotonin Reuptake Inhibitors (SSRIs) family. However, not all patients’ responses to this drug family and other treatment options are essential to be evaluated. One of the potential therapeutic targets is the dopamine system; dopaminergic antagonists are designed and showed promising treatment-resistant cases (Koo et al., 2010). Here, the contribution of four central genes in the dopamine transport process suggests the important role of this process in OCD pathogenesis. Another nominated critical process in OCD pathogenesis is learning, suggesting that OCD patients may have impairments in learning. It has been reported that patients with OCD are slower in learning (Gross-Isseroff et al., 1996). Memory is another involved process in OCD seed genes analysis. This fact can be justified by previous clinical findings of memory impairments in OCD (Penadés, Catalán, Andrés, Salamero, & Gastó, 2005). In other words, our analysis provides more supports for memory process dysregulation in OCD. Overall, the investigation introduced some feasible candidates showing substantial participation in our network integrity and functions. Variation in OCD genes may lead to the formation of different combinations of genes’ functions. As a result, each of these combinations can correspond to a specific OCD subtype model. Activation of a particular combination is attributed to a specific internal or external stimulation. Hence, the onset of a specific subgroup for each individual depends on the origin of the trigger and the influenced pathway.

Our study results suggest topological information of interactome profile of OCD reported genes from genetic studies; however, it is in its initial phase. Our PPI network analysis can serve as a paradigm for later experiments such as genomic, proteomic, and metabolomic investigations. Thus, by enriching interaction sets, a more in-depth understanding of OCD etiology may be reachable. Finally, screening 31 genes associated with OCD leads to introducing a panel, including 7 crucial genes. This panel may potentially be used for diagnostic and therapeutic proposes of OCD after extensive validation tests.

Ethical Considerations

Compliance with ethical guidelines

All ethical principles are considered in this article.