I’ve been running DEXSeq and I’m having trouble understanding the results that I’m getting. I’ve used the suggested model
~ sample + exon + condition:exon
Which should give me differential exon usage between my two conditions.
When plotting the exon usage of my two conditions and the exon expression of my two conditions, I’ve generally found that they’re the inverse of each other. For example, if I have an exon in a gene that looks to be differentially expressed between the two conditions, it doesn’t look to have differential exon usage. Vice versa, if an exon in a gene appears to have no differential expression between two conditions, it’s shown to have differential usage.
I’m struggling to understand the difference between these two metrics, does anyone have a good explanation?
Secondly, If I wanted to test for differential exon expression, what would I have to change my model to, to achieve this?
condition libType countName F1_S1 A paired-end raw_counts//sorted_Flowcell_A_1.bam.dexseq.count F1_S2 B paired-end raw_counts//sorted_Flowcell_A_2.bam.dexseq.count F1_S3 A paired-end raw_counts//sorted_Flowcell_A_3.bam.dexseq.count F1_S4 A paired-end raw_counts//sorted_Flowcell_A_4.bam.dexseq.count F1_S5 B paired-end raw_counts//sorted_Flowcell_A_5.bam.dexseq.count F1_S6 A paired-end raw_counts//sorted_Flowcell_A_6.bam.dexseq.count F1_S7 B paired-end raw_counts//sorted_Flowcell_A_7.bam.dexseq.count F1_S8 A paired-end raw_counts//sorted_Flowcell_A_8.bam.dexseq.count F2_S1 A paired-end raw_counts//sorted_Flowcell_B_1.bam.dexseq.count F2_S2 B paired-end raw_counts//sorted_Flowcell_B_2.bam.dexseq.count F2_S3 A paired-end raw_counts//sorted_Flowcell_B_3.bam.dexseq.count F2_S4 A paired-end raw_counts//sorted_Flowcell_B_4.bam.dexseq.count F2_S5 B paired-end raw_counts//sorted_Flowcell_B_5.bam.dexseq.count F2_S6 B paired-end raw_counts//sorted_Flowcell_B_6.bam.dexseq.count F2_S7 A paired-end raw_counts//sorted_Flowcell_B_7.bam.dexseq.count F2_S8 A paired-end raw_counts//sorted_Flowcell_B_8.bam.dexseq.count
dxd <- DEXSeqDataSetFromHTSeq(list.files("raw_counts/", full.names=T), sampleData = sampleTable, design = ~ sample + exon + condition:exon, flattenedfile = "../scrips/genome.gff") dxd <- estimateSizeFactors(dxd) dxd <- estimateDispersions(dxd) dxd <- testForDEU(dxd) dxd <- estimateExonFoldChanges(dxd, fitExpToVar="condition") dxr1 <- DEXSeqResults(dxd)