Comparing RNA-Seq and Ribo-Seq

Question

I'm attempting to compare total RNA-seq with Ribo-Seq, to determine if changes in Ribo-Seq are due to changes in transcriptional expression (akin to the analysis performed by Anota2Seq). I am, however, struggling to understand how to accomplish this in DESeq2.

My design is as follows:

cts <- read.delim("Interaction.txt", row.names=1)
coldata <- read.delim("Sample_Sheet.txt", row.names=1)
cts <- round(cts)
dds <- DESeqDataSetFromMatrix(countData = cts, colData = coldata, design = ~Genotype + Type + Genotype:Type)
dds <- DESeq(dds)
resultsNames(dds)

> coldata
                Genotype Type
C_24h_Ribo_1     Control Ribo
C_24h_Ribo_2     Control Ribo
C_24h_Ribo_3     Control Ribo
CHIK_24h_Ribo_1    Viral Ribo
CHIK_24h_Ribo_2    Viral Ribo
CHIK_24h_Ribo_3    Viral Ribo
C_24h_RNA_1      Control  RNA
C_24h_RNA_2      Control  RNA
C_24h_RNA_3      Control  RNA
CHIK_24h_RNA_1     Viral  RNA
CHIK_24h_RNA_2     Viral  RNA
CHIK_24h_RNA_3     Viral  RNA

> resultsNames(dds)
[1] "Intercept"                 "Genotype_Viral_vs_Control" "Type_RNA_vs_Ribo"          "GenotypeViral.TypeRNA"

I do not understand which contrast is the one I need, or how to access the results i'm looking for (or interpret them). Any help?

Answer 1

Adding to what benn pointed out, at this line:

dds <- DESeqDataSetFromMatrix(countData = cts, colData = coldata, design = ~Genotype + Type + Genotype:Type)

You specified that the counts, can be a explained by Genotype (being Control or Viral) or Type (RNAseq or Ribo) or an interaction effect

What you need is this coefficient, Genotype:Type. To use it, go

res <- results(dds,name="GenotypeViral.TypeRNA")

Based on your setup, you are always comparing RNA to Ribo (see "Type_RNA_vs_Ribo" in resultsNames). If changes to the abundance of a mRNA is similar in RNAseq or Ribo, then the coefficient of the interaction term will be zero and it will be zero in the log2foldchange column in results(res). A positive fold change indicates the log2foldchange between RNA and Ribo is higher in Viral.

If you would like to compare ribo versus and rna-seq, then set RNA as the reference by:

coldata$Type = factor(coldata$Type,levels=c("RNA","Ribo"))
# run deseq again
dds <- DESeqDataSetFromMatrix(countData = cts, colData = coldata, design = ~Genotype + Type + Genotype:Type)
dds <- DESeq(dds)

With this you will compare by default Ribo versus RNA, and the interaction term will tell you how much more enriched or depleted the ribo fraction is, in the Viral samples.

Would recommend checking the vignette out here:

http://bioconductor.org/packages/devel/bioc/vignettes/DESeq2/inst/doc/DESeq2.html#contrasts