setting Coldata levels in Deseq2

Question

[AP]

Hello everyone,

I am trying to see differential expression of a fungal genome with and without treatment of a chemical gossypol. I have one sample for control and 3 treated samples at 3 different timepoints. And I have 2 replications. So, alltogether I have 8 samples (4-R1 and 4-R2). I am going to see differetial expression of genes from the HTseq count data in DESeq2 and my script looks like this:

directory <- "C:\\Users\\hp\\Desktop\\RNA.seq.data"

sampleFiles <- c('count.1','count.2','count.3','count.4','count.5','count.6','count.7','count.8')

sampleCondition <- c('R1control','R1time1','R1time2','R1time4','R2control','R2time1','R2time2','R2time4')

sampleTable<-data.frame(sampleName=sampleFiles, fileName=sampleFiles, condition=sampleCondition)

ddsHTSeq<- DESeqDataSetFromHTSeqCount(sampleTable = sampleTable,
                                   directory = directory,
                                   design= ~ condition)
ddsHTSeq
colData(ddsHTSeq)$condition<-factor(colData(ddsHTSeq)$condition, levels=c('count'))

dds<-DESeq(ddsHTSeq)

res<-results(dds)

res<-res[order(res$padj),]

head(res)

I am am confused on how to assign levels in colData step. When I do colData(ddHTSeq) , it gives me a DataFrame with 8 rows and 1 column condition

  condition      factor
count.1           NA
count.2           NA
count.3           NA
count.4           NA
count.5           NA
count.6           NA
count.7           NA
count.8           NA

Error in designAndArgChecker(object, betaPrior) : variables in the design formula cannot have NA values.

[daniele.avancini]

May I ask what would you expect from colData(ddHTSeq)? I guess in the factor column you should have if the sample received or didn't receive the chemical gossypol. Right? Also, your code above is missing a ) at the end of colData.

[Devon Ryan]

I put the missing ) back, I may have accidentally deleted it when I reformatted this post.

Answer 1

[Devon Ryan]

You've thoroughly confused DESeq2 by telling it that you don't have replicates. Please replace your sampleTable line with the following:

sampleTable <- c("control", "t1", "t2", "t3", "control", "t1", "t2", "t3")

Your replicate number isn't a condition.

Also, you MUST delete the following line:

colData(ddsHTSeq)$condition<-factor(colData(ddsHTSeq)$condition, levels=c('count')

That is what's causing there to be no levels in your colData. You do not need to do anything like this, control will already default to being the base level for all comparisons.

[AP]

Devon Ryan,

Do you mean adding to Sample condition not sample table right?. I added in sample condition the following line and got the result as below:

sampleCondition <- c ("control", "t1", "t2", "t4", "control", "t1", "t2", "t4")

log2 fold change (MLE): condition t4 vs control 
Wald test p-value: condition t4 vs control 
DataFrame with 6 rows and 6 columns
               baseMean log2FoldChange     lfcSE      stat       pvalue
              <numeric>      <numeric> <numeric> <numeric>    <numeric>
XLOC_011966 298932.4717      11.167532 0.5369988  20.79620 4.685118e-96
XLOC_015099   4680.3330       8.129216 0.5742844  14.15538 1.729899e-45
XLOC_006378    897.0141       6.881496 0.5325741  12.92120 3.417696e-38
XLOC_003469    221.9459      -7.663667 0.6253951 -12.25412 1.596718e-34
XLOC_004664    808.9044       9.633222 0.8395317  11.47452 1.771611e-30
XLOC_015318   5509.5439       8.408940 0.7359340  11.42622 3.092968e-30
                    padj
               <numeric>
XLOC_011966 3.493692e-92
XLOC_015099 6.449928e-42
XLOC_006378 8.495253e-35
XLOC_003469 2.976681e-31
XLOC_004664 2.642180e-27
XLOC_015318 3.844043e-27

Here I got 6 rows and 6 columns. Is it keeping both control at same and other time points as separate columns? And does this looks ok?

[Devon Ryan]

Yes, sampleCondition rather than sampleTable. For the results of other comparisons, see help(results).

[AP]

Thank you so much Devon Ryan.