I have an experiment which I did in blocks (week 1 vs week 2) and 4 groups (A, B, C and D).
My understanding is that I can do the analysis controlling for block by doing this.
dds_results <- DESeqDataSetFromMatrix(countData = genes.matrix, colData = sample_info, design = ~ block + group)
I can then extract results of the contrast of interest:
results.AB <- results(dds_results, contrast=c("group", "A", "B"), alpha=0.05, pAdjustMethod = "BH", lfcThreshold = lfc)
Which gives me some number of up- and down-regulated DE genes (lets say 50)
Problem 1 I used lfcShrink to get shrunken log2foldChanges but when I do the summary() it seems to indicate far fewer genes are significant now.
results.AB <- lfcShrink(dds_results, contrast=c("group", "A", "B"), alpha=0.05, pAdjustMethod = "BH", lfcThreshold = lfc)
The p-values between results() and lfcShrink() are not the same -- shouldn't the p-values stay the same?
Does this mean that I need to adjust type=? but I get an error that type can only be adjusted for coeff.
Problem 2 I wanted to see how consistent blocks were but I reran the above splitting the data by block and essentially do the same as above and the 2 blocks yields drastically different numbers, though there is a high degree of gene overlap.
block 1: 60 DE genes block 2: 25 DE genes
I did some descriptives on the data and both look the same, maybe overall levels of expression are lower in block 2.
Why would this occur? And what is the best way to control for block?