I am attempting to create a module-trait correlation heatmap from a WGCNA consensus network as described in this vignette: https://labs.genetics.ucla.edu/horvath/CoexpressionNetwork/Rpackages/WGCNA/Tutorials/Consensus-RelateModsToTraits.pdf
However when I calculate the correlation between the module eigengenes and my trait information, each value in the resulting correlation matrix comes back NA. I have two sets with four samples each (this is a test run before adding in all samples). My trait data is:
> Traits
[[1]]
[[1]]$data
CCF_200mM IMR_200mM
CCF_200_1_S35_L006_R1_001 1 0
CCF_200_2_S36_L006_R1_001 1 0
CCF_200_3_S37_L006_R1_001 1 0
CCF_200_4_S38_L006_R1_001 1 0
[[2]]
[[2]]$data
CCF_200mM IMR_200mM
IMR_200_1_S51_L006_R1_001 0 1
IMR_200_2_S52_L006_R1_001 0 1
IMR_200_3_S53_L006_R1_001 0 1
IMR_200_4_S54_L006_R1_001 0 1
Where 0 denotes NO and 1 denotes YES.
Additionally, all values within both consMEs (my module eigengene data variable) and my trait information are numeric:
> str(Traits)
List of 2
$ :List of 1
..$ data:'data.frame': 4 obs. of 2 variables:
.. ..$ CCF_200mM: num [1:4] 1 1 1 1
.. ..$ IMR_200mM: num [1:4] 0 0 0 0
$ :List of 1
..$ data:'data.frame': 4 obs. of 2 variables:
.. ..$ CCF_200mM: num [1:4] 0 0 0 0
.. ..$ IMR_200mM: num [1:4] 1 1 1 1
> str(consMEs)
List of 2
$ :List of 1
..$ data:'data.frame': 4 obs. of 76 variables:
.. ..$ ME10: num [1:4] -0.3744 0.8431 -0.0945 -0.3742
.. ..$ ME48: num [1:4] 0.448 0.542 -0.403 -0.586
.. ..$ ME68: num [1:4] 0.315 0.637 -0.332 -0.62
.. ..$ ME3 : num [1:4] -0.262 -0.5514 0.7918 0.0216
.. ..$ ME44: num [1:4] -0.509 -0.0171 0.8118 -0.2856
.. ..$ ME19: num [1:4] -0.662 0.108 0.723 -0.169
.. ..$ ME8 : num [1:4] -0.247 0.628 0.295 -0.677
.. ..$ ME64: num [1:4] -0.401 0.443 0.545 -0.587
.. ..$ ME7 : num [1:4] -0.536 -0.379 0.732 0.183
.. ..$ ME33: num [1:4] -0.74 0.31 0.578 -0.149
.. ..$ ME61: num [1:4] -0.331 0.704 0.216 -0.589
...
$ :List of 1
..$ data:'data.frame': 4 obs. of 76 variables:
.. ..$ ME10: num [1:4] 0.42 -0.487 -0.507 0.574
.. ..$ ME48: num [1:4] 0.576 -0.432 -0.557 0.413
.. ..$ ME68: num [1:4] 0.46 0.232 -0.844 0.152
.. ..$ ME3 : num [1:4] 0.818 -0.538 -0.105 -0.175
.. ..$ ME44: num [1:4] 0.551 -0.361 -0.618 0.429
.. ..$ ME19: num [1:4] 0.595 0.395 -0.494 -0.496
.. ..$ ME8 : num [1:4] 0.764 0.127 -0.41 -0.481
.. ..$ ME64: num [1:4] 0.707 -0.536 0.229 -0.401
.. ..$ ME7 : num [1:4] -0.247 -0.455 -0.142 0.844
.. ..$ ME33: num [1:4] -0.812 0.113 0.144 0.554
.. ..$ ME61: num [1:4] -0.766 0.427 0.465 -0.126
...
The code I am attempting is:
moduleTraitCor = list()
moduleTraitPvalue = list()
for (set in 1:nSets)
{
moduleTraitCor[[set]] = cor(consMEs[[set]]$data, Traits[[set]]$data, use = "p");
moduleTraitPvalue[[set]] = corPvalueFisher(moduleTraitCor[[set]], exprSize$nSamples[set]);
}
Which returns:
> moduleTraitCor
[[1]]
CCF_200mM IMR_200mM
ME10 NA NA
ME48 NA NA
ME68 NA NA
ME3 NA NA
ME44 NA NA
ME19 NA NA
ME8 NA NA
ME64 NA NA
ME7 NA NA
ME33 NA NA
ME61 NA NA
...
[[2]]
CCF_200mM IMR_200mM
ME10 NA NA
ME48 NA NA
ME68 NA NA
ME3 NA NA
ME44 NA NA
ME19 NA NA
ME8 NA NA
ME64 NA NA
ME7 NA NA
ME33 NA NA
ME61 NA NA
...
How can I get cor() to compute the correlation between Traits and my consensus module eigengenes?
> sessionInfo()
R version 3.4.3 (2017-11-30)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 14.04.5 LTS
Matrix products: default
BLAS: /usr/lib/libblas/libblas.so.3.0
LAPACK: /usr/lib/lapack/liblapack.so.3.0
locale:
[1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
[3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
[5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
[7] LC_PAPER=en_US.UTF-8 LC_NAME=C
[9] LC_ADDRESS=C LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] ggplot2_2.2.1 flashClust_1.01-2 WGCNA_1.62
[4] fastcluster_1.1.24 dynamicTreeCut_1.63-1
Thanks for your help.
You should make sure that you are using numbers in your argument. If one of the variables used in the equation is a character vector, the results would be NAs.