Multi mapping reads

[Evi-050]

Hello,

I would really use some help with the following:
I am trying some different parameters in STAR (version 2.5.3a) to see how the percent of multi-mapping reads affects my final goal. I am trying the following combinations:

1. --outSAMmultNmax -1 \ #default --outFilterMultimapNmax 1000 \

2.--outSAMmultNmax 1 \ --outFilterMultimapNmax 1000 \

3. --outSAMmultNmax 1 \

The rest of the flags are the same for all 3 approaches:
--sjdbOverhang 48
--quantMode GeneCounts
--twopassMode Basic
--runThreadN 8
--genomeDir ${starIndexesDir}
--sjdbGTFfile ${gtfFile}
--alignSJDBoverhangMin 1
--readFilesIn TCGA-97-8174_R1.fastq TCGA-97-8174_R2.fastq
--outFileNamePrefix $dir/TCGA-97-8174-01A
--outSAMtype BAM SortedByCoordinate
--outTmpDir $dir/STAR_tmp_TCGA-97-8174-01A
--outSAMunmapped Within
--bamRemoveDuplicatesType UniqueIdentical
--outMultimapperOrder Random
--outFilterMismatchNoverLmax 0.04
--outFilterMatchNminOverLread 0.33
--outFilterScoreMinOverLread 0.33
--limitBAMsortRAM 55000000000
--limitGenomeGenerateRAM 35000000000
--chimOutType WithinBAM
--chimSegmentMin 10
--chimJunctionOverhangMin 10;

From the manual:
" The --outSAMmultNmax parameter limits the number of output alignments (SAM lines) for multimappers. For instance, --outSAMmultNmax 1 will output exactly one SAM line for each mapped read.

default: -1 int: max number of multiple alignments for a read that will be output to the SAM/BAM files. "
-1 all alignments (up to –outFilterMultimapNmax) will be output.

I would guess that --outFilterMultimapNmax is only needed and works when the --outSAMmultNmax is on default. But instead of seeing a difference between approach 1 and 2, I see between 1/2 and 3. The absence of --outFilterMultimapNmax 1000 flag increases the % of reads mapped to too many loci from 0% (approach 1,2 in presence of --outFilterMultimapNmax) to 17.2% (approach 3).

I f you have any idea why this is happening, it would be really enlightening!

Thank you in advance,
Evi

[alexdobin]

Hi Evi,

--outFilterMultimapNmax is a filter - it determines whether read alignments are output at all. If a read maps to more loci than --outFilterMultimapNmax, none of the read alignments will be output.
If a read maps to fewer loci than --outFilterMultimapNmax, then up to --outSAMmultNmax alignments will be written to the SAM/BAM file.

So you approaches 1 and 2 yield the same mapping stats because the filter is the same - however, you will see only one SAM alignment for each read in the 2nd approach.
The 3 approach uses default --outFilterMultimapNmax 10 and thus you will see a lot of reads considered "mapped to too many loci" and not output at all.

Note that for large --outFilterMultimapNmax you also need to increase the --winAnchorMultimapNmax to at least the same value (from default 50), which may slow down the mapping.

Cheers
Alex

Cheers
Alex

[Evi-050]

Hello,

Thank you very much for the detailed response!!

Evi