RNA-Seq Analysis with DEXSeq

Question

[justinkablan225]

Hello to all, I'm trying to follow the DEXSeq package workflow to do same basis analysis on RNA-Seq data. Then, at the step of constructing the DEXSeqDataSet object I've an error message. The gtf file given below is that what was use for sequence alignment and read count. All count reads are also OK. How can I fix it ?

countFiles = list.files(pattern=".txt$", full.names=TRUE)

basename(countFiles) 

[1] "C_M1.txt" "C_M2.txt" "C_M3.txt" "C_M4.txt" "C_S1.txt" "C_S2.txt" "C_S3.txt" "C_S4.txt"

flattenedFile = list.files(pattern="gtf$", full.names=TRUE)

basename(flattenedFile) 

[1] "Homo_sapiens.GRCh38.95.gtf"

dxd = DEXSeqDataSetFromHTSeq(
+   countFiles,
+   sampleData=sampleTable,
+   design= ~ sample + exon + condition:exon,
+   flattenedfile=flattenedFile) 

Error in FUN(X[[i]], ...) : subscript out of bounds

[WouterDeCoster]

I added code markup to your post for increased readability. You can do this by selecting the text and clicking the 101010 button. When you compose or edit a post that button is in your toolbar, see image below:

In addition, you posted this as a "Tutorial" while it is clearly a question. I have changed that as well.

Please keep these things in mind for future questions.

Answer 1

[shawn.w.foley]

DEXSeq is meant to be used to generate a new GFF file then run with that annotation, not the annotation that the mapping was performed with. The detailed workflow is here.

From the vignette, you need to run:

python /path/to/library/DEXSeq/python_scripts/dexseq_prepare_annotation.py Homo_sapiens.GRCh38.95.gtf Homo_sapiens.GRCh38.DEXSeq.gff

The dexseq_prepare_annotation.py file identifies exonic regions unique to different isoforms.

And rather than running HTSeq directly, you should run:

python /path/to/library/DEXSeq/python_scripts/dexseq_count.py Homo_sapiens.GRCh38.DEXSeq.gff input.sam output.txt

Then in the R session you should see:

basename(flattenedFile) 
[1] "Homo_sapiens.GRCh38.DEXSeq.gff"