Weighted Gene Co-expression Network Analysis of Endometriosis and Identification of Functional Modules Associated With Its Main Hallmarks

Abstract

Although many genes have been identified using high throughput technologies in endometriosis (ES), only a small number of individual genes have been analyzed functionally. This is due to the complexity of the disease that has different stages and is affected by various genetic and environmental factors. Many genes are upregulated or downregulated at each stage of the disease, thus making it difficult to identify key genes. In addition, little is known about the differences between the different stages of the disease. We assumed that the study of the identified genes in ES at a system-level can help to better understand the molecular mechanism of the disease at different stages of the development. We used publicly available microarray data containing archived endometrial samples from women with minimal/mild endometriosis (MMES), mild/severe endometriosis (MSES) and without endometriosis. Using weighted gene co-expression analysis (WGCNA), functional modules were derived from normal endometrium (NEM) as the reference sample. Subsequently, we tested whether the topology or connectivity pattern of the modules was preserved in MMES and/or MSES. Common and specific hub genes were identified in non-preserved modules. Accordingly, hub genes were detected in the non-preserved modules at each stage. We identified sixteen co-expression modules. Of the 16 modules, nine were non-preserved in both MMES and MSES whereas five were preserved in NEM, MMES, and MSES. Importantly, two non-preserved modules were found in either MMES or MSES, highlighting differences between the two stages of the disease. Analyzing the hub genes in the non-preserved modules showed that they mostly lost or gained their centrality in NEM after developing the disease into MMES and MSES. The same scenario was observed, when the severeness of the disease switched from MMES to MSES. Interestingly, the expression analysis of the new selected gene candidates including CC2D2A, AEBP1, HOXB6, IER3, and STX18 as well as IGF-1, CYP11A1 and MMP-2 could validate such shifts between different stages. The overrepresented gene ontology (GO) terms were enriched in specific modules, such as genetic disposition, estrogen dependence, progesterone resistance and inflammation, which are known as endometriosis hallmarks. Some modules uncovered novel co-expressed gene clusters that were not previously discovered.

Keywords: endometriosis; expression; hub genes; module; weighted gene co-expression network.

FIGURE 1

Biweight midcorrelation (bicor) between gene expression levels in normal endometrium and ES disease stages. Each dot corresponds to a gene and the red line represents to the best-fit regression line. Correlation coefficients are reported in respective graph. Gene expression patterns were strongly correlated between normal endometrium and ES stages (A,B) as well as between ES stages (C).

FIGURE 2

Biweight midcorrelation (bicor) between connectivity patterns in normal endometrium and ES disease stages; x- and y- axes denote the average number of connectivity of a gene and the red line represents to the best-fit regression line. Connectivity pattern was less preserved between NEM and ES stages (A,B) compared to MMES vs. MSES (C).

FIGURE 3

Identification of gene co-expression modules in normal endometrium data (reference data) using average hierarchical linkage clustering; the y-axis denotes the co-expression distance and the x-axis corresponds to genes. Dynamic tree cutting was applied to identify modules by dividing the dendrogram at significant branch points. Modules are displayed with different colors in the horizontal bar immediately below the dendrogram, with gray representing unassigned genes.

FIGURE 4

Composite preservation statistics of reference modules after assignment in minimal/mild (MMES) and mild/severe (MSES) ES data (as test data): (A,B) the composite statistic medianRank (y-axis) as a function of the module size; each point represents a module, labeled by color. The dashed red line indicates the threshold medianRank = 10; (C,D) the summary statistic Z (y-axis), as a function of the module size; each point represents a module, labeled by color as in (A). The dashed red line indicates the threshold Z = 5. A module was considered as non-preserved if it had medianRank ≥ 8 or Zsummary < 5.

FIGURE 5

The percentage of hub genes in normal endometrium (NEM), minimal/mild endometriosis (MMES), and mild/severe endometriosis (MSES) in each non-preserved module. The percentages of specific hub genes for each condition, and common hub genes between them were calculated. The total number of hub genes and other details for each module can be found in Supplementary File S5.

FIGURE 6

The results of RT-PCR analysis of eight gene candidates in NEM, MMES, and MSES; ^∗ and ^∗∗ means are significantly different at 0.05 and 0.01 probability level, respectively.