Thank you... i'll give it a try.
• 0 views
•
link
written
by
lfkhajavi •
RSEM error: RSEM currently does not support gapped alignments, sorry!
hello,
i'm using STAR v2.5.2b and RSEM v1.2.31
STAR command:
for i in RNAseqFastq/*_R1_fastq.gz; do STAR --runThreadN 16 --genomeDir STAR/STAR_indices --sjdbGTFfile hg38_annotation.gtf --readFilesIn $i ${i%R1_fastq.gz}R2_fastq.gz --outFileNamePrefix
STAR/BAMfiles/$(basename $i fastq.gz) --readFilesCommand zcat --sjdbOverhang 50 --outFilterType BySJout --outFilterMultimapNmax 20 --outMultimapperOrder Random --alignSJoverhangMin 8 --alignSJDBoverhangMin 1
--outFilterMismatchNmax 999 --outFilterMismatchNoverLmax 0.04 --alignIntronMin 20 --alignIntronMax 1000000 --alignMatesGapMax 1000000 --quantMode TranscriptomeSAM --alignSoftClipAtReferenceEnds No
--outSAMstrandField intronMotif --outSAMmultNmax 1 --clip5pNbases 3 --outFilterIntronMotifs RemoveNoncanonical --outSAMtype BAM SortedByCoordinate; done
RSEM prepare reference:
rsem-prepare-reference -p 4 --gtf hg38_annotation.gtf hg38_genome.fa RSEM/hg38
RSEM calculate expression:
for i in STAR/*.toTranscriptome.out.bam; do rsem-calculate-expression --append-names --calc-ci --output-genome-bam --alignments --paired-end $i RSEM/hg38 RSEM/parALL/RNAseq_quals/$(basename $i
_Aligned.toTranscriptome.out.bam); done
ERROR: RSEM currently does not support gapped alignments, sorry!
I can't figure out where my problem lies... any help will be greatly appreciated. Thank you in advance, Leila
• 6,070 views
•
link
updated
by
Jake Warner
•
written
by
lfkhajavi •
•
written
by
lfkhajavi •
1 answer
The error is pretty verbose. RSEM does not support gapped alignment. So when you map to a genome (which has introns) you need to use an end to end aligner like bowtie or set up STAR appropriately. Try these options for STAR to prohibit gaps:
--alignEndsType EndToEnd --alignIntronMax 1 --alignIntronMin 2 --scoreDelOpen -10000 --scoreInsOpen -10000
• 0 views
•
link
written
by
Jake Warner
Before you do this consider @Alex's answer (which I had included above).
if you generated the genome indexes with annotations, you will splicing to annotated junctions, as --alignIntronMax only controls the annotated junctions. You can either (i) use genome with annotations and --alignSJDBoverhangMin 999 (any number > read length) while mapping, or (ii) re-generate the genome without annotations.
hello,
changing option --alignSJDBoverhangMin 1 to --alignSJDBoverhangMin 999 did not fix the problem... unfortunately. re-mapping using option --alignEndsType EndToEnd didn't remedy the problem
the commands listed were working well with RSEM until I re-ran star with the addition of --outMultimapperOrder Random & --clip5pNbases 3 It's with this new alignment that I'm having issues with RSEM...
Could the option --clip5pNbases 3 be introducing soft- or hard- clippings?
• 0 views
•
link
written
by
lfkhajavi •
Log in to answer this question.
STAR is producing spliced alignment which RSEM does not like. @Alex Dobin (author of STAR) recommends this to turn off spliced alignments.
From this page.