I am trying to understand how DESeq2 calculates Log2foldchange, because it does not match my calculations of fold change from my normalized counts.
My script:
CountTable_POS = read.table("DIF_pos_counts", header=TRUE, sep="")
Design = data.frame( row.names = colnames( CountTable_POS ),
condition = c( "untreated", "untreated", "untreated",
"treated", "treated", "treated" ),
type = c( "single-end", "single-end", "single-end",
"single-end", "single-end","single-end" ) )
dds_pos<-DESeqDataSetFromMatrix(countData = CountTable_POS,
colData = Design,
design = ~ condition)
dds_pos<-DESeq(dds_pos, fitType="local")
norm<-counts(dds_pos, normalized=T)
write.csv(as.data.frame(norm), file="norm_counts.csv")
res_pos<-results(dds_pos, contrast=c("condition", "treated", "untreated"))
write.csv(as.data.frame(res_pos), file="DESeq2_out.csv")
RESULTS
from norm_counts.csv:
Locus Control1 Control2 Control3 Treated1 Treated2 Treated3 Mean_Control Mean_Treated Basemean Fold change
2R_188654 3.8127 5.7220 5.7191 1.9351 4.9898 2.9518 5.0846 3.2922 4.1884 0.6475
2R_188736 3.8127 3.8147 6.6723 1.9351 4.9898 5.9036 4.7666 4.2762 4.5214 0.8971
From DESeq2_out:
Locus Basemean Log2foldchange P value Fold change
2R_188654 4.1884 -0.3108 0.4093 0.8062
2R_188736 4.5214 -0.0819 0.8279 0.9448
As you can see, the basemeans (mean of control and treated counts) match- which tells me there is no internal error- but the fold changes do not. Why is this???
My script:
CountTable_POS = read.table("DIF_pos_counts", header=TRUE, sep="")
Design = data.frame( row.names = colnames( CountTable_POS ),
condition = c( "untreated", "untreated", "untreated",
"treated", "treated", "treated" ),
type = c( "single-end", "single-end", "single-end",
"single-end", "single-end","single-end" ) )
dds_pos<-DESeqDataSetFromMatrix(countData = CountTable_POS,
colData = Design,
design = ~ condition)
dds_pos<-DESeq(dds_pos, fitType="local")
norm<-counts(dds_pos, normalized=T)
write.csv(as.data.frame(norm), file="norm_counts.csv")
res_pos<-results(dds_pos, contrast=c("condition", "treated", "untreated"))
write.csv(as.data.frame(res_pos), file="DESeq2_out.csv")
RESULTS
from norm_counts.csv:
Locus Control1 Control2 Control3 Treated1 Treated2 Treated3 Mean_Control Mean_Treated Basemean Fold change
2R_188654 3.8127 5.7220 5.7191 1.9351 4.9898 2.9518 5.0846 3.2922 4.1884 0.6475
2R_188736 3.8127 3.8147 6.6723 1.9351 4.9898 5.9036 4.7666 4.2762 4.5214 0.8971
From DESeq2_out:
Locus Basemean Log2foldchange P value Fold change
2R_188654 4.1884 -0.3108 0.4093 0.8062
2R_188736 4.5214 -0.0819 0.8279 0.9448
As you can see, the basemeans (mean of control and treated counts) match- which tells me there is no internal error- but the fold changes do not. Why is this???
Comment