Hey!

I am going through the very helpful WGCNA tutorials and constructed a co-expression network, but now I am stuck with further steps and interpretation of my modules.

At first some information regarding my data and progress:

•  I have 60 samples from humans (assigned to 2 groups a 30 people)
•  Microarray data, pre-processed (batch effects, normalization, log2 transformation)
•  as input I used ~ 5000 genes, which I previously filtered variance-based
•  I decided to use a signed & weighted network (followed step-wise tutorial)
•  I ended up with 16 modules

Then, I tried to correlate the eigengene values of the modules with the external trait information. I was insecure if I can use correlation for binary information (my only "trait" information is if the sample is from human with/ without disease). However, I´ve read all posts I was able to find about this and there it was discussed that this is possible and basically a t-test between the groups.

So my code was this: moduleTraitCor = cor(MEs, Anno_180_T, use = "p"); Where I have my eigengenes in MEs and my sample class info in Anno_180_T.

So do you think this approach is ok so far?

My problem is: No single correlation between module eigengenes and the binary trait is significant.

I don´t know what I could do know, even functional enrichment does not really make sense when there is no relationship at all...

My question of interest was if I could find differentially co-expressed modules between disease/healthy persons. And I have additional data from 2 other time points which I thought to maybe analyze with consensus analysis/ eigengene networks...

Thank you for any help/ suggestions!!

PS: I also had a look at CEMiTool, which is somewhat based on WGCNA I think... There they seem to do gene set enrichment analysis for group differences. As far as I understood it mathematically, the difference is that there the expression of the whole module goes in, not only the eigengene. I tried this for my sample and every module was significant... I am really confused.

It is not uncommon to find no relationship between disease and co-expression modules. The first thing I would do is run a individual gene DE analysis (using say limma, but you could also use standardScreeningBinaryTrait in WGCNA) to see how much signal there is on the individual gene level. I would also perhaps not filter by variance on log-transformed data since, after log-transformation, the low-expressed probes tend to have more variance than the higher-expressed ones. I'd filter by mean expression and try to keep at least 10k genes in the analysis.

If you have additional time points, you could run the same analysis in the additional time points, especially the individual gene association analysis to see whether the individual gene associations among the three sets are at least somewhat concordant. Then you can try a consensus network analysis.

Hi Lin, thanks for using CEMiTool!

CEMiTool's automatic filter is, as you may have noticed, quite stringent, and it tends to leave in the analysis mostly genes with significant activity in their given pathways. This is probably why you found that all your modules ended up being differentially expressed between your sample classes. You could try making the filter less stringent by playing around with the filter_pval parameter in the cemitool function, set by default to 0.1. 

Also, since you applied a log-transformation to your data, it is possible that you introduced a dependency between the mean and the variance, which affects correlation-based analyses. You can check this by observing a straight line when using the plot_mean_var function on your CEMiTool object, which if present, you can try to correct by setting the apply_vst parameter in the cemitool function to TRUE.