Network-based bioinformatics analysis of HOX genes and their related genes in endometriosis

Background: Dysregulation of some HOX genes are important for the development of endometriosis. Endometriosis is a complex gynecologic disorder affecting as many as 10-15% of premenopausal women. The pathogenesis of endometriosis, as the presence of endometrium-like tissue outside the uterine cavity, is largely unknown. HOX genes, encoding homeodomain transcription factors, are dynamically expressed in endometrium, in which they are necessary for endometrial growth, differentiation, and implantation. Few investigations have been studied in this vast range of HOX genes with a network-based view. Network theory allows for a holistic understanding of the role of genes in diseases. The purpose of this study is to identify genes that cause this disease by focusing on the HOX family genes and their regulatory genes. Materials & Methods: Fifteen samples of eutopic and ectopic endometrium collected from patients in proliferation phase. PCR array was performed on sample groups for 84 HOX genes and their related genes. The statistical significance was determined by using t-test. PCA and Hierarchical clustering were carried out between different groups. The genes co-expression networks were constructed based on correlation between normalized gene expression data. Gene function was annotated based on Gene Ontology, the intervention in other diseases and expression in other tissues using the DAVID to clarify the function of important genes and the mechanism of endometriosis. Results: Among the analyzed genes, 30, 40 and 35 significantly differentially expressed genes were in eutopic, ectopic compared with normal and ectopic with eutopic samples, respectively. Twenty-eight genes had no significant alteration in none of the sample groups. Conclusion: The most significant genes and the genes in co-expression networks are effective genes in development, metabolic and neural processes. Another important point is that the some non-differently expressed genes in networks have many interactions with the significant genes, indicating that these non-significant genes also are likely contributed in endometriosis.

co-expression networks are effective genes in development, metabolic and neural processes. Another important point is that the some non-differently expressed genes in networks have many interactions with the significant genes, indicating that these non-significant genes also are likely contributed in endometriosis. Keywords: endometriosis, HOX, PCR-array, Co-expression network Background Endometriosis is known as the presence of endometrial tissues inside an abnormal anatomical location other than those of the uterus. It is a common, estrogendependent gynecological disorder that can lead to multiple manifestations, including dysmenorrhea, pelvic pain, pelvic mass, infertility and malignant behaviors [1][2][3][4], with prevalence rates of 10-15% in women of reproductive age [5,6]. Symptoms of endometriosis often affect social relationships, sexuality and mental health [7]. The understanding of the pathogenesis of endometriosis is poor.
The combination of retrograde menstruation of endometrial tissue fragments through the fallopian tubes into the peritoneal cavity and the immune dysfunctions are the most accepted causes of the pathogenesis of endometriosis. The live endometrial fragments in retrograde menstruation can hold in the ectopic site through adhesion, aggression and angiogenesis process and then form the ectopic lesions; therefore, the inherent biological properties of eutopic endometrium may play an important role in the pathogenesis of endometriosis [8].
Although the association between endometriosis and infertility is supported in many literature [9][10][11], the involved mechanisms have not been figured out. Several factors have been reported, including distortion of pelvic anatomy, hormone secretion abnormalities, alterations in peritoneal fluid and fertilization disorders. Accordingly, changes in endometrial development in patients with endometriosis may help to infertility associated with endometriosis. It is believed that patients with endometriosis have reduced the incidence of implantation. Ultrastructural defects in endometrium in women with endometriosis have been reported. Molecular markers of the endometrial receptors undergo change in patients with endometriosis; the expression patterns of the endometrium in the native endometrium of women with endometriosis are aberrant. Also, changes in biochemical or molecular markers in the tissue have been pointed out [12].
The Homeobox (HOX) gene family operates through a conserved homeodomain as transcriptional regulators for controlling fetal morphogenesis. In addition, the lack of proper expression of the HOX gene results in developmental anomalies. As transcription factors, the HOX genes regulate other downstream target genes leading to proper endometrial development and receptivity in implantation and function in different tissues of the fetus to provide spatial and temporal coordination for each cell [13]. Initially, HOX genes were expressed during embryonic development. However, the stable expression of the HOX genes has been observed in the reproductive system. HOX genes are essential for endometrial growth, differentiation and receptivity by mediating some functions of the sex steroids during each reproduction period. 39 HOX genes have been identified to divide into four HOX loci [14], each localized on a different chromosome and containing 9 to 11 orthologous genes. Genes in the four clusters can be categorized into 13 paralogous groups. Individual genes of the HOXA cluster are assigned a distinctive identity for each part of the paramesonephric duct, leading to the development of the fallopian tubes (hoxa9), uterus (hoxa10), lower part of the uterus and cervix (hoxa11) and upper vagina (hoxa13) [15]. HOXB has been detected to be closely related to the self-renewal of the hematopoietic stem cells and effective proliferation of hematopoietic progenitor cells [16]. Studies show a collinear activation of the HOXC genes in the limb ectoderm besides, possible implication of HOXC genes in the nail/claw/hoof transition [17]. While all four HOX gene clusters are initially activated during embryonic development, both HOXA and HOXD clusters are then reactivated during the development of appendicular skeletons, involved in the building of the limbs [18]. The identification of molecular markers of endometriosis is likely to help understanding its course and eventually, to uncover novel genes to be targeted by drugs.
In the present survey, we compared gene expression using PCR-array in a collection https://david.ncifcrf.gov) [19]. Co-expression gene network was constructed using correlation between genes expression data and visualized by Cytoscape software. In silico co-expression network construction can really provide useful information about pathogenesis of disease. Several parameters can be identified to help establish hierarchies between genes from a network or even belonging to preset remote networks. In general, the identification of new deregulated gene expression networks in endometriosis highlighted some of the interactions between the protein products of these genes and caused a significant reduction in infertility due to early detection of precursor lesions [20].  PCR-Array: PCR array was conducted for HOX family genes by using RT2 profiler PCR array human homeobox (HOX) genes kit (Qiagen, Cat.No: PAHS-083Z) and RT2 SYBR green ROX qPCR mastermix (Qiagen, Cat.No: 330502). There were primers for 84 tests and housekeeping genes on mentioned HOX PCR array kit according to manufacturer's instructions. Cycling condition on step one plus real time PCR system (ABI), was included an initial denaturation at 95°C for 10 minutes followed by 40 cycles at 95°C for 15 seconds, and 60°C for 1 minute. Results were expressed as values of cycle threshold (Ct) then normalized to GAPDH as housekeeping gene (control gene). To calculate fold changes in mRNA abundance 2 -∆∆Ct method was used.

Statistical analysis:
The statistical significance of difference in expression levels between groups was analyzed using Student's t-test with two tails. P-value <0.05 was considered statistically significant. PCA (Principal Component analysis) [21] was carried out using the Minitab 16 statistical software. PCA Was performed to simplify the large amount of data. PCA was applied algorithm on log2 of Differentially Expressed Genes (DEGs) and non-DEGs data to underlying cluster structures of endometrial samples.
Hierarchical clustering [22] was implemented by correlation coefficient and complete linkage of in Minitab 16 statistical software [23]. The methods were performed twice, using normalized gene expression data of DEGs and non-DEGs among sample groups in three clusters.

Construction of the gene co-expression network:
The genes co-expression networks were constructed of the basis of the normalized gene expression data.
Matlab was recruited to compute the Pearson correlation between each pair of genes , and then the significant correlation pairs (p<0.05) were imported into Cytoscape software [24] (version 3.6.1) for visualization, and the circle algorithms were used [25]. The co-expression networks were constructed of the normal tissues.
Coding gene functional analysis: Gene function was annotated based on the GO, the intervention in other diseases and expression in other tissues using the DAVID; https://david.ncifcrf.gov [19] to clarify the function of DEGs and the mechanism of endometriosis.

Results
PCR-array with primer sets targeting HOX genes and their associated genes used to evaluate the changes in their expression between ectopic and eutopic tissues with control tissues and ectopic with eutopic tissues. The genes expression which were significantly (P<0.05) altered (DEGs) in ectopic, eutopic with normal groups and ectopic with eutopic groups are illustrated in Table I. Among the analyzed genes, there are 33.33% (28) genes in non-DEG. The expression levels of 30, 41 and 35 genes were significantly altered in eutopic and ectopic compared with normal samples and ectopic rather than eutopic samples , respectively. Twelve genes were common between the three groups and 10 of 41 DEGs were exclusively related to ectopic versus normal tissue, 5 of 35 DEGs were solely related to eutopic versus normal tissue and and 3 of 35 DEGs were only related to ectopic versus eutopic tissue. One of the interesting points is that at least in one of the analysed groups, the expression out of five of the seven studied HOXB genes, were significantly upregulated. Three genes of the HOXA cluster were present. Only HOXA9 was significantly up-regulated in ectopic versus normal and eutopic samples. Seven genes of the HOXC cluster were studied and found that expression of five of them were significantly altered in three groups of analysis. Eight genes of the HOXD cluster were analyzed, and it was shown that expression of six of them underwent significant change in different studies ( Fig. 1 and Table I).
DEGs were categorized into six groups based on venn diagram (Fig. 1) The last group consists of the genes expression of which has siginificant ly altered in the three analysis ( Fig. 1 and Table I).
For further investigation, hierarchical clustering was performed on the nine studied samples using DEGs (Fig. 3A) and non-DEGs (Fig. 3B)  Seven co-expression networks were pointed out based on pairwise correlation of gene expression data and visualized using Cytoscape (Fig. 4). As pointed out above, expression about 64% (54/84) of the genes significantly changed between normal, eutopic and ectopic samples. All genes in the network B1 were significant genes in diseased tissues, and significant genes of ectopic tissues are abundant in the network A1 (71%) ( Table II, Fig. 4).
The signaling pathways present in the KEGG and Reactome databases searched for genes in the networks (Table III, Fig. 4), have identified two pathways of "Signaling pathways regulating pluripotency of stem cells" and "Transcriptional misregulation in cancer" from the KEGG database and pathway of "activation of anterior HOX genes during the early embryogenesis" from the Reactome database. Fifteen of which were in these pathways and the MEIS1 was present in all three pathways.

Functional analysis
The DEGs and the genes in seven networks (Fig. 4) were mapped to the DAVID database to investigate the GO, pathways, interference with other diseases and expression in other tissues. DAVID is the most popular tool in the field of functional annotation (Tables IV, V and S1).

Analysis of interference with other diseases:
According to analysis of interference with other diseases of genes, many DEGs in ectopic tissues play a role in developmental diseases class, such as cleft lip and clubfoot and metabolic diseases class like bone mineral density.
Many DEGs in eutopic tissues are also involved in neurological diseases class like parkinson's disease and autism in addition to developmental and metabolic diseases class. Many networks A1,2,3 and 4 and B2 genes affected developmental diseases.
Most genes in the network C cause neurological diseases. (Table IV).

GO annotation analysis:
On the basis of the GO annotation, the most of DEGs in the eutopic and ectopic tissues, and different networks were involved in the development or morphogenesis of various systems, especially skeletal system. Only most network C genes play roles in the development of the nervous system such as dopaminergic neuron and cerebellum (Table S1).

Discussion
The pathogenesis of endometriosis which was related to infertility remains unclear, but studies have shown that dysregulation HOX genes are important for the development of endometriosis [2,[26][27][28]. The investigation of the molecular mechanism of the disease is critically important for the diagnosis, treatment and prognosis. The expression analysis by PCR-array can provide information about the expression differences. In this study, the role of HOX genes and their associated genes in the disease was investigated. The data was screened obtained 16 upregulated and 25 downregulated DEGs of ectopic to normal samples, 7 upregulated and 23 downregulated DEGs of eutopic samples to normal samples and 21 upregulated and 14 downregulated DEGs of ectopic samples to eutopic samples from 9 samples of endometriosis patients and 9 samples of normal women (Fig. 1).
The GO and related disease analysis of DEGs revealed that most of them were significantly involved in development process especially skeletal system in ectopic and eutopic tissues as well as nervous system in eutopic tissues (Tables S1 and V).
Using PCA and hierarchical clustering, it was observed that these DEGs were accurately classified into three ectopic, eutopic, and normal groups. Accordingly, ectopic and normal samples were located in the same top node in the hierarchical clustering and, they are also closer together compered to eutopic sample in the PCA analysis (Figs. 2, 3). Evidence indicates the eutopic endometrium has greater ability in the formation of new blood vessels, cell metastasis and invasion, and the level of its gene expression is different compared to that of normal endometrium [8].
Therefore, abnormal alterations of gene expression of eutopic endometrium might be the source of the pathogenesis of endometriosis. Although the ectopic lesions were established from eutopic tissues, evidence indicated that the properties of each are very distinctive, that it suggests genetic regulation could be involved in the alteration of these phenotypes [29]. Therefore, despite the fact that eutopic tissue is similar in appearance to normal endometrial tissue, since eutopic cells can transform into other cells, the magnitude of gene expression changes was high in these cells. The genes such as EN1 that are in both CtoU and UtoN groups, the expression of which have significantly altered in eutopic tissue to normal and ectopic endometrium, are probably important in the initiation of the disease. Also, the genes like ARX that are in both CtoU and CtoN groups, whose expression has significantly altered in ectopic endometrium to normal and eutopic endometrium, are probably important in creating of the ectopic tissue and spreading the disease.
The importance of HOXB genes in endometriosis were also pointed out [26]. In the current survey, five genes out of 7 genes of HOXB cluster were of DEGs in the eutopic, ectopic or both of them and most of them upregulated (Fig. 1). It is possibly due to their role in sensory perception of pain, hemopoiesis, hematopoietic stem cell differentiation and angiogenesis [16].
In this study, our main concern was about the pathogenesis process and relationships between the DEGs. Thus co-expression network construction with expression data of all studied genes in normal samples identified seven networks (Fig. 4). The most genes of six networks were involved in developmental and metabolic processes, but genes in the network C involved in neural processes. Of the six genes forming this network, there are three DEGs dependent on the eutopic tissues. Probably the neural processes in the eutopic endometrium mentioned above (Table S1) were associated with the network C genes. Investigation of the tissues which express the genes of each network indicate that most genes in each network are expressed in a specific tissue, showing connections between the genes of each network together (Table V). Non-DEGs that were a minority in a network and had many connections with DEGs (Fig. 4) may also play a role in the development of endometriosis, but this role has not been shown a significant change in the expression of the genes. Hence, some non-DEGs genes like SIX1 may be effective in creating ectopic tissue and BARX1 genes which are likely to be effective in creating eutopic tissues.
In conclusion, our findings support the general notion, concerning the role of HOX family genes in development of endometriosis. In addition, the genes, such as VAX2, PHOX2B, LMX1B, DLX1and SIX4, have high significant alteration in gene expression of endometriosis tissue samples. Most of these genes have not previously been reported as endometriosis-related genes. In addition, since some non-DEGs have many connections to DEGs in co-expression networks, these non-DEGs also play important roles in endometriosis. Accordingly, the functional study showed that the significant genes engaged in developmental, neurological and metabolic pathways.

Ethics approval and consent to participate
The study protocol was approved by Ethical Committee of Royan Institute for Reproductive Biomedicine. Informed consent was obtained from all the participants for being included in the study.

Availability of data and material
Data sharing is not applicable to this article as no datasets were generated or analysed during the current study.

Competing interests
The authors declared no conflict of interest.    Figure 1 Venn diagram presentation of DEGs (P-value <0.05) in ectopic samples compared to normal