Hi, I am new in R and I have performed DESEQ2 analysis and I am confused at a particular step
I try to plot the normalized counts into box plot but it does not match with the fold change results generated by Deseq2
For instance, using CADM2 gene as an example and I compared OAC wBE with normal mucosa. The log2 fold change shows a positive 1.6 fold change higher in OAC wBE than normal mucosa.
But when I plot them into a box plot, it shows that the normal mucosa has higher expression than OAC wBE. I have tried to extract the normalized counts of CADM2 and calculate the ratio of mean between the two groups, the result is similar with the log2foldchange result generated by DEseq2. However, when I plot it, the mean (shown in hexagonal white box) is higher in normal mucosa than OAC wBE. May I know, what is the problem here?
Here is the R code that I used
countsTable <- read.csv("merge_normal_cancer_5.csv",row.names = 1,sep=",") summary(colSums(countsTable)) condition <- rep("Normal_Muc",121) condition[which(grepl("woBE",names(countsTable)))] <- "OAC_woBE" condition[which(grepl("wBE",names(countsTable)))] <- "OAC_wBE" colours = c("darksalmon","lightblue2", "olivedrab2") colData<-data.frame(condition,colours=colours[c(as.numeric(as.factor(condition)))]) row.names(colData)<-colnames(countsTable) colData colData$condition <- relevel(factor(colData$condition), ref = "Normal_Muc") library(DESeq2) dds<-DESeqDataSetFromMatrix(countData= countsTable,colData=colData,design=~condition) dds<-DESeq(dds) res <- results(dds) head(results(dds, tidy=TRUE)) summary(res) head(res) head(dds) View(dds) normalized.counts <- as.data.frame(counts(dds, normalized=TRUE )) write.csv(normalized.counts, "merge_normalized_counts.csv") res_woBE_muc <- results(dds,contrast=c("condition","OAC_woBE","Normal_Muc")) res_wBE_muc <- results(dds, contrast=c("condition", "OAC_wBE","Normal_Muc")) res_woBE_wBE <- results(dds, contrast=c("condition", "OAC_woBE","OAC_wBE")) write.csv(res_woBE_muc, "res_woBE_muc.csv") write.csv(res_wBE_muc, "res_wBE_muc.csv") write.csv(res_woBE_wBE, "res_woBE_wBE.csv") res_woBE_muc <- res_woBE_muc[order(res_woBE_muc$padj),] res_wBE_muc <- res_wBE_muc[order(res_wBE_muc$padj),] rstudioapi::writeRStudioPreference("data_viewer_max_columns", 1000L) res_woBE_muc["CRIP3",] res_wBE_muc["CADM2",] res_woBE_wBE["SPP1",] ab <- plotCounts(dds, gene="SPP1", intgroup="condition", returnData=TRUE) pq <- ggplot(ab,aes(x=condition, y=count, color=condition)) +geom_boxplot() pq + geom_jitter(shape=16, position=position_jitter(0.2)) + labs(title="plot",x="group", y="count SPP1")+scale_y_continuous(trans='log10')+stat_summary(fun.y=mean,geom="point",shape=23, size=4) # include your problematic code here with any corresponding output # please also include the results of running the following in an R session sessionInfo( )