HtSEq count multimappers

[tellsparck]

Hi,
I am new to RNASeq analysis and have some sequences aligned by STAR.
I am planning to use HTSeq count to generate an input file for DeSeq, but learned now that HtSeq will discard multi-mappers. Is there a way to make it accept multimappers? Or is there another program like easyRNASeq that does not discard multimappers? I do have 10-35% multimappers in my data set.
Thanks!

[chadn737]

Why do you want to keep them. You cannot know where a multimapped read comes from and retaining it increases the potential error.

[tellsparck]

Thanks Chadn.Can't we set an upper limit for the multimappers (I guess cufflink uses 10) so that many genuine multigene family members are counted and repetititive sequences are not?

[chadn737]

But that makes the assumption that your reads are indeed coming from annotated genes and not the repetitive sequences. If you try doing what cufflinks does, and disperse the value of the read amongst all of its possible mappings, then this will affect the variation in all the genes, because that read still comes from only one of them. I would especially be worried about this if you are dealing with samples that have large differences in the number of multi-mapped genes. Because that will affect the variation of one sample more than the other. Yes, you potentially loose data in discarding these reads, but you can also be more confident in the results.

Sometimes when you really go back and look at these multimapped reads, they are a mess and you realize that your results are the better for discarding them.

[tellsparck]

Yes it is not straightforward. But the junk reads among multimappers will end up not being annotated while meaningful ones will. Wonder whether anyone has tried this approach to get extract the useful reads out of multimappers.But the problem of how to assign reads among various positions still remain...