I am trying to look at post-transcriptional regulation using exon and intron reads as discussed in this paper (http://www.nature.com/nbt/journal/vaop/ncurrent/full/nbt.3269.html). Essentially they look for changes in exon reads between two samples and changes in intron reads between two samples and classify transcripts as post-transcriptionally regulated when there is a discrepancy between changes in exons and changes in introns.
They use a linear model ~ region + condition + region:condition. Where region is either exon(ex) or intron(in) and condition is treatment or control. I can get post-transcriptionally regulated genes from the default results(), but I also want to plot ∆exons vs ∆introns for treatment vs control and I wanted to check that the contrasts I'm using are correct.
layout <- data.frame(row.names = colnames(countMatrix), condition = c(rep('control',3), rep('treatment',3)), region = rep(c("ex","in"),each=ncol(cntEx))) dds <- DESeqDataSetFromMatrix(countData = countMatrix, colData = layout, design = ~ region*condition) dds <- DESeq(dds, betaPrior = FALSE) results <- results(dds, alpha=0.1) results_exons <- results(dds, contrast=c('condition','treatment','control')) results_introns <- results(dds, contrast=list(c('condition_treatment_vs_control','regionin.conditiontreatment'))) plot(results_exons$log2FoldChange, results_introns$log2FoldChange)