Question: RSEM error: RSEM currently does not support gapped alignments, sorry!
1
gravatar for lfkhajavi
3.3 years ago by
lfkhajavi30
lfkhajavi30 wrote:

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

ADD COMMENTlink modified 3.3 years ago by Jake Warner810 • written 3.3 years ago by lfkhajavi30

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.

However, please note that RSEM does * not * support gapped alignments. So make sure that your aligner does not produce alignments with intersions/deletions.

ADD REPLYlink modified 3.3 years ago • written 3.3 years ago by GenoMax96k
2
gravatar for Jake Warner
3.3 years ago by
Jake Warner810
Jake Warner810 wrote:

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
ADD COMMENTlink modified 3.3 years ago • written 3.3 years ago by Jake Warner810

Thank you... i'll give it a try.

ADD REPLYlink written 3.3 years ago by lfkhajavi30
1

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.

ADD REPLYlink written 3.3 years ago by GenoMax96k

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?

ADD REPLYlink modified 3.3 years ago • written 3.3 years ago by lfkhajavi30
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 2432 users visited in the last hour
_