Hello, I'm very new to R and I'm trying to run DESeq2 analysis on a RNASeq dataset that has 6 different conditions.
The 6 conditions are Dat1_WT, Dat1_KO, Dat2_WT, Dat2_KO, Dat3_WT, Dat3_KO.
I want to generate volcano plots and extract the differentially expressed genes for the following pairwise comparisons: Dat1_WT vs Dat1_KO Dat2_WT vs Dat2_KO Dat3_WT vs Dat3_KO
The issue arises when I look at the calculated results sets using the resultsNames()
function and get different comparisons than what I'd like. I've looked at several other posts here but wasn't able to fix the issue (while keeping the apeglm instead of ashr.
Here is my entire workflow:
mydata<-read.table("data/july9_salmon.merged.gene_counts_length_scaled.tsv",sep="\t",header=T)
#convert to matrix
mymat<-as.matrix(mydata[,3:22])
rownames(mymat)<-mydata[,2]
# Making a vector describing the experimental conditions
condition<-factor(c("Dat1_WT", "Dat2_WT", "Dat2_KO", "Dat2_KO", "Dat2_KO", "Dat2_KO", "Dat3_WT", "Dat3_WT", "Dat3_WT", "Dat3_KO", "Dat3_KO", "Dat1_WT", "Dat3_KO", "Dat1_WT", "Dat1_KO", "Dat1_KO", "Dat1_KO", "Dat2_WT", "Dat2_WT", "Dat2_WT"))
#checking sample table
sampleTable<-data.frame(sample_name=colnames(mymat),condition=condition)
dds<-DESeqDataSetFromMatrix(countData=round(mymat),colData=sampleTable,design=~condition)
##Running DESeq2
dds<-DESeq(dds)
##Making PCA
# Make rlog (regularlized log) transform of counts
rld<-rlog(dds)
plotPCA(rld,intgroup="condition") + geom_text_repel(aes(label=colData(dds)$sample_name),size=3)
#Now here is the issue
resultsNames(dds)
#Not the comparisons I'd like
[1] "Intercept" "condition_Dat2_WT_vs_Dat2_KO" "condition_Dat3_KO_vs_Dat2_KO"
[4] "condition_Dat3_WT_vs_Dat2_KO" "condition_Dat1_KO_vs_Dat2_KO" "condition_Dat1_WT_vs_Dat2_KO"
### Calculating differential expression between different timepoints / conditions
Dat1<-results(dds,contrast=c("condition","Dat1_WT","Dat1_KO"),alpha=0.05)
Dat2<-results(dds,contrast=c("condition","Dat2_WT","Dat2_KO"),alpha=0.05)
Dat3<-results(dds,contrast=c("condition","Dat3_WT","Dat3_KO"),alpha=0.05)
## What do the results look like?
Dat1
Dat2
Dat3
res_my<-lfcShrink(dds,coef="condition_Dat1_WT_vs_Dat1_KO",res=res_my,type="apeglm")
Here at the end, I get the following error:
Error in lfcShrink(dds, contrast = contrast, res = dd1) : type='apeglm' shrinkage only for use with 'coef'
I need to use the apeglm, not ashr with contrast.
How do I fix this? Any help is highly appreciated.
When I run relevel(dds, ref = "ND_WT"), I get the following error:
I see in that post that you do the following:
What does that do?
You have to relevel the column in colData, not the dds object itself. The condition line just mocks up some factors. You need to do
dds$condition <- relevel(dds$condition, ref="ND_WT")