library(beeswarm)
library(gplots)
## Warning: package 'gplots' was built under R version 4.3.2
##
## Attaching package: 'gplots'
## The following object is masked from 'package:stats':
##
## lowess
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(RColorBrewer)
library(tidyverse)
## Warning: package 'ggplot2' was built under R version 4.3.2
## Warning: package 'tidyr' was built under R version 4.3.2
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ forcats 1.0.0.9000 ✔ readr 2.1.5
## ✔ ggplot2 3.5.1 ✔ stringr 1.5.1
## ✔ lubridate 1.9.3 ✔ tibble 3.2.1
## ✔ purrr 1.0.2 ✔ tidyr 1.3.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(WGCNA)
## Loading required package: dynamicTreeCut
## Loading required package: fastcluster
##
## Attaching package: 'fastcluster'
##
## The following object is masked from 'package:stats':
##
## hclust
##
##
##
## Attaching package: 'WGCNA'
##
## The following object is masked from 'package:stats':
##
## cor
library(flashClust)
##
## Attaching package: 'flashClust'
##
## The following object is masked from 'package:fastcluster':
##
## hclust
##
## The following object is masked from 'package:stats':
##
## hclust
#import WGCNA outputs
rootdir <- '/Users/usri/Desktop/Camk2a.WGCNA.8june/july/'
outputfigs <- '/Users/usri/Desktop/Camk2a.WGCNA.8june/july/'
FileBaseName <- 'camk2a.Wgcna.module'
cleanDat <- readRDS("~/Desktop/Camk2a.WGCNA.8june/july/cleanDat.after.outlier.removal.15samples.rds")
dim(cleanDat)
## [1] 1544 15
head(cleanDat, 3)
## X.106P2fractionIP X.120P2fractionIP X.117P2fractionIP
## AAK1|Q2M2I8 23.48784 23.92167 23.18945
## AAK1|Q2M2I8.1 25.21936 25.97150 25.12306
## ABAT|P80404 27.55898 27.40838 27.95464
## X.123P2fractionIP X.137P2fractionIP X.131P2fractionIP
## AAK1|Q2M2I8 23.31690 23.23492 22.92868
## AAK1|Q2M2I8.1 25.51681 25.18611 25.24855
## ABAT|P80404 27.51978 27.19804 27.35241
## X.143P2fractionIP X.133P2fractionIP X.106HomogenateIP
## AAK1|Q2M2I8 23.61881 22.99170 22.93829
## AAK1|Q2M2I8.1 25.07017 25.34287 25.46261
## ABAT|P80404 27.27474 27.16163 27.14351
## X.120HomogenateIP X.117HomogenateIP X.123HomogenateIP
## AAK1|Q2M2I8 23.87091 23.85305 23.11254
## AAK1|Q2M2I8.1 25.80004 25.22436 25.56972
## ABAT|P80404 26.87173 26.93971 26.87173
## X.131HomogenateIP X.143HomogenateIP X.133HomogenateIP
## AAK1|Q2M2I8 23.40135 19.97449 23.80321
## AAK1|Q2M2I8.1 25.29442 25.52288 25.54158
## ABAT|P80404 27.12597 27.27474 26.92225
net <- readRDS("~/Desktop/Camk2a.WGCNA.8june/july/camk2a.final.net.rds")
MEs <- net$MEs
MEs<-tmpMEs<-data.frame()
MEList = moduleEigengenes(t(cleanDat), colors = net$colors)
MEs = orderMEs(MEList$eigengenes)
colnames(MEs)<-gsub("ME","",colnames(MEs)) #let's be consistent in case prefix was added, remove it.
tmpMEs <- MEs #net$MEs
colnames(tmpMEs) <- paste("ME",colnames(MEs),sep="")
MEs[,"grey"] <- NULL
tmpMEs[,"MEgrey"] <- NULL
MEs <- tmpMEs
colnames(MEs)<-gsub("ME","",colnames(MEs))
MEs
## tan turquoise yellow greenyellow
## X.106P2fractionIP 0.002286292 -0.140061765 -0.27955616 -0.002570635
## X.120P2fractionIP -0.120859965 -0.031801206 -0.32547276 0.408460405
## X.117P2fractionIP -0.183685117 -0.108613605 -0.44123605 0.423884673
## X.123P2fractionIP -0.201417732 -0.009744013 -0.34286651 0.412058843
## X.137P2fractionIP -0.412131752 -0.383416336 -0.08059125 -0.265240843
## X.131P2fractionIP 0.359446246 -0.398738082 -0.10381111 -0.360464667
## X.143P2fractionIP -0.020596037 -0.301034794 -0.06961930 -0.153300089
## X.133P2fractionIP -0.396108527 -0.368437479 -0.13239815 -0.317628307
## X.106HomogenateIP 0.231776392 0.236850352 0.23919658 0.043715728
## X.120HomogenateIP -0.076339854 0.235162511 0.25794282 -0.138742698
## X.117HomogenateIP -0.195773544 0.239193849 0.17775195 -0.153816139
## X.123HomogenateIP 0.043415824 0.238581003 0.27887421 0.092151133
## X.131HomogenateIP 0.316620387 0.274321956 0.23482013 -0.116749780
## X.143HomogenateIP 0.470288615 0.271085934 0.33531843 0.268866448
## X.133HomogenateIP 0.183078771 0.246651674 0.25164717 -0.140624069
## pink magenta purple green blue
## X.106P2fractionIP 0.01887541 -0.19343566 0.261704219 0.14537444 0.2737437
## X.120P2fractionIP 0.17798158 0.43125082 0.502950915 -0.34025051 0.1520824
## X.117P2fractionIP 0.21990967 -0.37439010 -0.092544078 0.38555799 0.4707747
## X.123P2fractionIP 0.32838544 0.09049241 0.374687993 0.08071588 0.2301341
## X.137P2fractionIP -0.42986549 -0.12314501 -0.037425105 0.27410872 0.1484854
## X.131P2fractionIP -0.23091230 -0.10470751 -0.401823883 0.33794708 0.1676142
## X.143P2fractionIP -0.38813302 -0.26765805 0.038611506 0.34054666 0.2480337
## X.133P2fractionIP -0.45988970 -0.34163156 -0.245486895 0.23140313 0.1265096
## X.106HomogenateIP 0.22760904 0.07101496 -0.170744596 -0.26169039 -0.2183401
## X.120HomogenateIP 0.27838637 0.39525708 0.318379295 -0.17152034 -0.1980493
## X.117HomogenateIP 0.25784141 -0.29850892 -0.398736958 -0.06224406 -0.2018206
## X.123HomogenateIP 0.09274704 0.27694708 -0.058753328 -0.16971057 -0.3364152
## X.131HomogenateIP -0.08044073 0.06756535 -0.107710398 -0.22059390 -0.2554036
## X.143HomogenateIP 0.01079054 0.26708948 0.020100045 -0.34637129 -0.2358256
## X.133HomogenateIP -0.02328525 0.10385964 -0.003208732 -0.22327284 -0.3715233
## red brown black salmon
## X.106P2fractionIP 0.26595078 -0.35269404 0.25149258 -0.16269702
## X.120P2fractionIP 0.27800499 -0.38742119 0.36770638 0.39502611
## X.117P2fractionIP 0.50297334 -0.19016970 -0.27106026 -0.33982731
## X.123P2fractionIP 0.41260075 -0.36160886 0.01582102 -0.26720920
## X.137P2fractionIP -0.02007925 0.29636985 0.33917169 0.36540251
## X.131P2fractionIP -0.28062443 0.45243526 0.27711135 0.16767152
## X.143P2fractionIP 0.14761594 0.29257079 -0.14160357 -0.07940971
## X.133P2fractionIP 0.04719025 0.38342670 0.32639851 0.20657914
## X.106HomogenateIP -0.21758651 -0.10447807 0.01611855 -0.12455437
## X.120HomogenateIP 0.01676919 -0.11221520 -0.48714456 0.25612715
## X.117HomogenateIP -0.23341681 -0.03910131 -0.07488413 -0.54553163
## X.123HomogenateIP -0.19301170 0.05463966 -0.22127167 -0.02874865
## X.131HomogenateIP -0.32487724 -0.02710729 -0.07498642 -0.09345293
## X.143HomogenateIP -0.20704422 0.03941463 -0.33394578 0.13758446
## X.133HomogenateIP -0.19446510 0.05593878 0.01107632 0.11303993
## create orderedModules object from MEs
# Replace this with your actual data or method of creating orderedModules
orderedModules <- data.frame(
ModuleName = colnames(MEs),
Color = colnames(MEs))
colnames(orderedModules) <- c("ModuleName", "Color")
orderedModules
## ModuleName Color
## 1 tan tan
## 2 turquoise turquoise
## 3 yellow yellow
## 4 greenyellow greenyellow
## 5 pink pink
## 6 magenta magenta
## 7 purple purple
## 8 green green
## 9 blue blue
## 10 red red
## 11 brown brown
## 12 black black
## 13 salmon salmon
##create toplot object from MEs (no gray module included)
toplot <- MEs
colnames(toplot) <- colnames(MEs)
rownames(toplot) <- rownames(MEs)
toplot <- t(toplot) #may have to transpose sometimes: t(toplot)
#rownames(toplot) <- paste(orderedModules[match(colnames(MEs),orderedModules[,2]),1]," ",rownames(toplot),sep="")
toplot
## X.106P2fractionIP X.120P2fractionIP X.117P2fractionIP
## tan 0.002286292 -0.12085997 -0.18368512
## turquoise -0.140061765 -0.03180121 -0.10861361
## yellow -0.279556163 -0.32547276 -0.44123605
## greenyellow -0.002570635 0.40846040 0.42388467
## pink 0.018875406 0.17798158 0.21990967
## magenta -0.193435658 0.43125082 -0.37439010
## purple 0.261704219 0.50295091 -0.09254408
## green 0.145374437 -0.34025051 0.38555799
## blue 0.273743679 0.15208244 0.47077470
## red 0.265950779 0.27800499 0.50297334
## brown -0.352694041 -0.38742119 -0.19016970
## black 0.251492576 0.36770638 -0.27106026
## salmon -0.162697024 0.39502611 -0.33982731
## X.123P2fractionIP X.137P2fractionIP X.131P2fractionIP
## tan -0.201417732 -0.41213175 0.3594462
## turquoise -0.009744013 -0.38341634 -0.3987381
## yellow -0.342866509 -0.08059125 -0.1038111
## greenyellow 0.412058843 -0.26524084 -0.3604647
## pink 0.328385442 -0.42986549 -0.2309123
## magenta 0.090492407 -0.12314501 -0.1047075
## purple 0.374687993 -0.03742510 -0.4018239
## green 0.080715884 0.27410872 0.3379471
## blue 0.230134075 0.14848538 0.1676142
## red 0.412600748 -0.02007925 -0.2806244
## brown -0.361608862 0.29636985 0.4524353
## black 0.015821024 0.33917169 0.2771113
## salmon -0.267209201 0.36540251 0.1676715
## X.143P2fractionIP X.133P2fractionIP X.106HomogenateIP
## tan -0.02059604 -0.39610853 0.23177639
## turquoise -0.30103479 -0.36843748 0.23685035
## yellow -0.06961930 -0.13239815 0.23919658
## greenyellow -0.15330009 -0.31762831 0.04371573
## pink -0.38813302 -0.45988970 0.22760904
## magenta -0.26765805 -0.34163156 0.07101496
## purple 0.03861151 -0.24548689 -0.17074460
## green 0.34054666 0.23140313 -0.26169039
## blue 0.24803367 0.12650964 -0.21834013
## red 0.14761594 0.04719025 -0.21758651
## brown 0.29257079 0.38342670 -0.10447807
## black -0.14160357 0.32639851 0.01611855
## salmon -0.07940971 0.20657914 -0.12455437
## X.120HomogenateIP X.117HomogenateIP X.123HomogenateIP
## tan -0.07633985 -0.19577354 0.04341582
## turquoise 0.23516251 0.23919385 0.23858100
## yellow 0.25794282 0.17775195 0.27887421
## greenyellow -0.13874270 -0.15381614 0.09215113
## pink 0.27838637 0.25784141 0.09274704
## magenta 0.39525708 -0.29850892 0.27694708
## purple 0.31837929 -0.39873696 -0.05875333
## green -0.17152034 -0.06224406 -0.16971057
## blue -0.19804932 -0.20182060 -0.33641520
## red 0.01676919 -0.23341681 -0.19301170
## brown -0.11221520 -0.03910131 0.05463966
## black -0.48714456 -0.07488413 -0.22127167
## salmon 0.25612715 -0.54553163 -0.02874865
## X.131HomogenateIP X.143HomogenateIP X.133HomogenateIP
## tan 0.31662039 0.47028861 0.183078771
## turquoise 0.27432196 0.27108593 0.246651674
## yellow 0.23482013 0.33531843 0.251647170
## greenyellow -0.11674978 0.26886645 -0.140624069
## pink -0.08044073 0.01079054 -0.023285251
## magenta 0.06756535 0.26708948 0.103859637
## purple -0.10771040 0.02010005 -0.003208732
## green -0.22059390 -0.34637129 -0.223272839
## blue -0.25540361 -0.23582560 -0.371523340
## red -0.32487724 -0.20704422 -0.194465096
## brown -0.02710729 0.03941463 0.055938779
## black -0.07498642 -0.33394578 0.011076316
## salmon -0.09345293 0.13758446 0.113039934
numericMeta <- readRDS("~/Desktop/Camk2a.WGCNA.8june/JUNE30/numericMeta.rds")
numericMeta
## SampleID Treatment Group
## X.106P2fractionIP 106 P2fractionIP CamK2a + LPS P2fractionIP
## X.120P2fractionIP 120 P2fractionIP CamK2a + saline P2fractionIP
## X.117P2fractionIP 117 P2fractionIP CamK2a + LPS P2fractionIP
## X.123P2fractionIP 123 P2fractionIP CamK2a + saline P2fractionIP
## X.137P2fractionIP 137 P2fractionIP CamK2a + LPS P2fractionIP
## X.131P2fractionIP 131 P2fractionIP CamK2a + saline P2fractionIP
## X.143P2fractionIP 143 P2fractionIP CamK2a + LPS P2fractionIP
## X.133P2fractionIP 133 P2fractionIP CamK2a + saline P2fractionIP
## X.106HomogenateIP 106 HomogenateIP CamK2a + LPS HomogenateIP
## X.120HomogenateIP 120 HomogenateIP CamK2a + saline HomogenateIP
## X.117HomogenateIP 117 HomogenateIP CamK2a + LPS HomogenateIP
## X.123HomogenateIP 123 HomogenateIP CamK2a + saline HomogenateIP
## X.131HomogenateIP 131 HomogenateIP CamK2a + saline HomogenateIP
## X.143HomogenateIP 143 HomogenateIP CamK2a + LPS HomogenateIP
## X.133HomogenateIP 133 HomogenateIP CamK2a + saline HomogenateIP
## Condition Condition.Test Sample
## X.106P2fractionIP P2fractionIP.CamK2a + LPS 1 X.106P2fractionIP
## X.120P2fractionIP P2fractionIP.CamK2a + saline 0 X.120P2fractionIP
## X.117P2fractionIP P2fractionIP.CamK2a + LPS 1 X.117P2fractionIP
## X.123P2fractionIP P2fractionIP.CamK2a + saline 0 X.123P2fractionIP
## X.137P2fractionIP P2fractionIP.CamK2a + LPS 1 X.137P2fractionIP
## X.131P2fractionIP P2fractionIP.CamK2a + saline 0 X.131P2fractionIP
## X.143P2fractionIP P2fractionIP.CamK2a + LPS 1 X.143P2fractionIP
## X.133P2fractionIP P2fractionIP.CamK2a + saline 0 X.133P2fractionIP
## X.106HomogenateIP HomogenateIP.CamK2a + LPS 1 X.106HomogenateIP
## X.120HomogenateIP HomogenateIP.CamK2a + saline 0 X.120HomogenateIP
## X.117HomogenateIP HomogenateIP.CamK2a + LPS 1 X.117HomogenateIP
## X.123HomogenateIP HomogenateIP.CamK2a + saline 0 X.123HomogenateIP
## X.131HomogenateIP HomogenateIP.CamK2a + saline 0 X.131HomogenateIP
## X.143HomogenateIP HomogenateIP.CamK2a + LPS 1 X.143HomogenateIP
## X.133HomogenateIP HomogenateIP.CamK2a + saline 0 X.133HomogenateIP
## X Genotype
## X.106P2fractionIP P2fractionIP.CamK2a.LPS P2fractionIP.CamK2a + LPS
## X.120P2fractionIP P2fractionIP.CamK2a.saline P2fractionIP.CamK2a + saline
## X.117P2fractionIP P2fractionIP.CamK2a.LPS P2fractionIP.CamK2a + LPS
## X.123P2fractionIP P2fractionIP.CamK2a.saline P2fractionIP.CamK2a + saline
## X.137P2fractionIP P2fractionIP.CamK2a.LPS P2fractionIP.CamK2a + LPS
## X.131P2fractionIP P2fractionIP.CamK2a.saline P2fractionIP.CamK2a + saline
## X.143P2fractionIP P2fractionIP.CamK2a.LPS P2fractionIP.CamK2a + LPS
## X.133P2fractionIP P2fractionIP.CamK2a.saline P2fractionIP.CamK2a + saline
## X.106HomogenateIP P2fractionIP.CamK2a.LPS HomogenateIP.CamK2a + LPS
## X.120HomogenateIP P2fractionIP.CamK2a.saline HomogenateIP.CamK2a + saline
## X.117HomogenateIP P2fractionIP.CamK2a.LPS HomogenateIP.CamK2a + LPS
## X.123HomogenateIP P2fractionIP.CamK2a.saline HomogenateIP.CamK2a + saline
## X.131HomogenateIP P2fractionIP.CamK2a.saline HomogenateIP.CamK2a + saline
## X.143HomogenateIP P2fractionIP.CamK2a.LPS HomogenateIP.CamK2a + LPS
## X.133HomogenateIP P2fractionIP.CamK2a.saline HomogenateIP.CamK2a + saline
regvars <- data.frame(as.factor(numericMeta[,"Genotype"]))
colnames(regvars) <- c("Genotype")
lm <- lm(data.matrix(MEs)~Genotype,data=regvars)
pvec <- rep(NA,ncol(MEs))
for (i in 1:ncol(MEs)) {
f <- summary(lm)[[i]]$fstatistic ## Get F statistics
pvec[i] <- pf(f[1],f[2],f[3],lower.tail=F) ## Get the p-value corresponding to the whole model
}
names(pvec) <- colnames(MEs)
Groupingx <- numericMeta$Genotype
metdat <- data.frame(Status=Groupingx)
pvec
## tan turquoise yellow greenyellow pink magenta
## 3.139295e-01 2.442273e-04 9.386552e-05 9.344510e-01 4.842534e-01 1.097167e-01
## purple green blue red brown black
## 6.698333e-01 6.714808e-03 2.440795e-06 3.897617e-02 9.949859e-01 7.543990e-02
## salmon
## 5.048735e-01
metdat
## Status
## 1 P2fractionIP.CamK2a + LPS
## 2 P2fractionIP.CamK2a + saline
## 3 P2fractionIP.CamK2a + LPS
## 4 P2fractionIP.CamK2a + saline
## 5 P2fractionIP.CamK2a + LPS
## 6 P2fractionIP.CamK2a + saline
## 7 P2fractionIP.CamK2a + LPS
## 8 P2fractionIP.CamK2a + saline
## 9 HomogenateIP.CamK2a + LPS
## 10 HomogenateIP.CamK2a + saline
## 11 HomogenateIP.CamK2a + LPS
## 12 HomogenateIP.CamK2a + saline
## 13 HomogenateIP.CamK2a + saline
## 14 HomogenateIP.CamK2a + LPS
## 15 HomogenateIP.CamK2a + saline
par(mfrow=c(3,3))
par(mar=c(4.5,6,4.5,1.5))
#pdf(file=paste0(outputfigs,"/BoxPlots_",FileBaseName,".pdf"),width=18,height=11.25)
# Replace this with your actual data or method of creating orderedModules
orderedModules <- data.frame(
ModuleName = colnames(MEs),
Color = colnames(MEs))
colnames(orderedModules) <- c("ModuleName", "Color")
orderedModules
## ModuleName Color
## 1 tan tan
## 2 turquoise turquoise
## 3 yellow yellow
## 4 greenyellow greenyellow
## 5 pink pink
## 6 magenta magenta
## 7 purple purple
## 8 green green
## 9 blue blue
## 10 red red
## 11 brown brown
## 12 black black
## 13 salmon salmon
library(beeswarm)
library(gplots)
for (i in 1:(nrow(toplot))) { # grey already excluded, no -1
titlecolor<-if(signif(pvec,2)[i] <0.05) { "red" } else { "black" }
boxplot(toplot[i,]~factor(numericMeta$Genotype),col=colnames(MEs)[i],ylab="Eigenprotein Value",main=paste0(orderedModules[match(colnames(MEs)[i],orderedModules[,2]),1]," ",colnames(MEs)[i],"\nANOVA p = ",signif(pvec,2)[i]),xlab="Genotype",las=1,col.main=titlecolor) #no outliers: ,outline=FALSE)
transcol=paste0(col2hex(colnames(MEs)[i]),"99")
beeswarm(toplot[i,]~factor(numericMeta$Genotype),method="swarm",add=TRUE,corralWidth=0.5,vertical=TRUE,pch=21,bg=transcol,col="black",cex=0.8,corral="gutter") }
