This is very good question, but as far as I know there really isn't a good answer. I will attempt to provide some comments and tips.
However, I'm trying to do similar with a dataset which may consist of repeats (I have previously aligned the dataset using bowtie, and a lot of reads aligned to repeats).
If this dataset is also PE and the alignment is outputting either uniquely mapped reads or a single (random) position for multimappers, it is pretty much as good as it gets. A lot of repeat regions are not exact copies, so using PE there is a good change that one read of the pair mapping to a single position is helping with mapping of the other read.
I need to know reasonably confidently where precisely those repeats are located
It is probably a matter of phrasing, but what one needs is to know the precise genomic location of the reads. The location of the repeat elements can be obtained from the Repeat Masker track of UCSC, but you probably know this already.
With multimapped reads this step could be done by seeing which of the pair is uniquely mapped and choosing the mate which maps closest to it.
bowtie already does this. Did you map the read pairs together or separately?
- There will always be a number of multimapping reads (fragments) that will be discarded. The information of how many are being discarded for your samples would be helpful.
- Mapping can be improved using both reads and tunning
bowtie settings (see bellow)
Finding a valid paired-end alignment where both mates align to
repetitive regions of the reference can be very time-consuming. By
default, Bowtie avoids much of this cost by imposing a limit on the
number of "tries" it makes to match an alignment for one mate with a
nearby alignment for the other. The default limit is 100. This causes
bowtie to miss some valid paired-end alignments where both mates lie
in repetitive regions, but the user may use the --pairtries or -y
options to increase Bowtie's sensitivity as desired. source
The immediate step after aligning the pairs is to flatten the pair into the DNA fragment that has been sequenced.
As indicated above, the aligner will take care of that sort of. It will gave the mapping information of each read pair, which is then used by
ShortStack which uses the coverage of uniquely mapped reads to place multimapping reads for ChIP-seq, and the results were mixed. That said, it was SE data and the code has been through some changes, so it be worth a try.
When requested, the authors of
RepEnrich provided me with the scripts for an utility called
RepConsensus which maps reads to a consensus reference of repeat regions. The caveats are:
- mappings are better for elements whose sequences have not diverged a lot yet
- it will tell the position of the protein in the repeat sequence, but no information about the surrounding genomic regions because there isn't any.
Whatever you do, do not use all alignments, that is double counting of reads, and think very carefully before using custom methods. A blog post detailing some of the potential issues can be found here (I wrote it).