Question: RSEM error: RSEM currently does not support gapped alignments, sorry!
1
gravatar for lfkhajavi
23 months 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 23 months ago by Jake Warner770 • written 23 months 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 23 months ago • written 23 months ago by genomax73k
2
gravatar for Jake Warner
23 months ago by
Jake Warner770
Jake Warner770 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 23 months ago • written 23 months ago by Jake Warner770

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

ADD REPLYlink written 23 months 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 23 months ago by genomax73k

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 23 months ago • written 23 months 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: 2125 users visited in the last hour