BACKGROUND
the main technological advances that accompanied the genomic and post-genomic eras are high-throughput sequencing and hybridization microarrays.
sequencing technology enabled scientists to obtain the full genomic sequence for many species, including the human and many model organisms.
sequencing technology is also being used to selectively re-sequence the human genome to detect genome variations such as single nucleotide polymorphisms or large-scale structural variations.
because understanding these variations can immediately impact medical sciences, making sequencing more efficient and accessible is imperative.
however, traditional methodologies used to sequence the first mammalian genomes remain expensive, time consuming, and labor intensive.
oligonucleotide microarray technology, used to interrogate the rna or dna content of a sample, has emerged as a widely accessible and effective tool for studying gene expression or detecting protein-dna interactions .
microarrays have also been used to detect genome variations, like snps or structural variations.
array-based hybridization also has limitations.
for example, probes can behave in highly non-uniform ways, and the effects of cross-hybridization and resolution limits are poorly understood.
although significant research efforts are focused on these problems, they remain inherent in all hybridization-based methods.
the new sequencing technology referred to as "second-generation" shows promise to eliminate many of the problems associated with traditional sequencing technology and also those with oligonucleotide microarray technology.
second-generation sequencers are able to sequence more quickly and at lower cost in terms of both money and labor.
new sequencing technologies are developed to sequence at greater depth, meaning that a clone can be sequenced from a sample even when that clone exists at very low abundance .
several second-generation technologies have been developed using diverse methods.
two recent second-generation sequencers are from  <dig> life sciences  and solexa .
the  <dig> sequencers use an emulsion method for dna amplification and a pyrosequencing protocol for sequencing by synthesis  at picolitre-scale volumes.
current  <dig> sequencers can produce 25– <dig> m nt of sequence in a single run, in the form of reads with length up to  <dig> nt , enabling this technology to be used for de-novo sequencing in addition to re-sequencing  <cit> .
solexa/illumina 1g sequencers also use sequencing by synthesis, with dna amplified on the surface of a flow cell, resulting in a random array of dense clusters  <cit> .
the solexa technology is faster and cheaper than that used in  <dig> sequencers, producing 1g nucleotides of sequence in one run but producing much shorter reads.
each individual read is roughly  <dig> to  <dig> bases in length .
the solexa sequencing technology has recently started producing breakthrough results.
research described in  <cit>  and  <cit>  employed solexa sequencers to obtain high-resolution genomic maps of several histone modifications, as well as localization data for the dna-binding proteins.
effectiveness of chip-sequencing has also been demonstrated by  <cit> , who used solexa sequencing to obtain locations of stat <dig> binding sites in hela s <dig> cells before and following ifn-γ stimulation.
mapping reads from the solexa sequencer presents an obvious algorithmic challenge: tens of millions of reads must be mapped to a large  genome in a reasonable amount of time.
strong efforts to design short-read mapping algorithms have resulted in methods that are effective in particular contexts.
the mapping algorithm implemented as part of the solexa analysis pipeline is named eland .
eland is optimized to map very short reads, with length at most  <dig> nt, and ignores the additional bases when the sequenced reads are longer.
eland also only allows at most two mismatches between the read and the genomic sequence to which it maps, which will clearly be too few for longer reads.
despite these restrictions, eland remains very useful for many mapping tasks because it is extremely efficient.
the sxoligosearch algorithm  can quickly map reads of varying length, using different criteria, with performance depending on both the read length and mapping criteria.
the performance of sxoligosearch depends on use of the proprietary synabase data structure.
this data structure is a heavily compressed and annotated index for the reference genome that retains all non-redundant information.
the gains in mapping speed by using such a data structure come at a cost in terms of the memory required for the synabase data structure.
extreme memory requirements of the data structure makes sxoligosearch unsuitable for use on hardware available in most labs.
in most re-sequencing applications accuracy of mapping is a primary concern.
we must know with great accuracy what part of the genome was actually sequenced.
there are several reasons why it might be difficult to determine the location in the reference genome from which a read was derived, or even if a read was derived from the reference genome.
these include problems with the experiments, such as sequencing errors or sample contamination.
mapping is also made more difficult by repeats in the genome, and by polymorphisms.
while mapping algorithms cannot be expected to be robust to all such problems, effort should be made to make the algorithms as robust as possible.
two sources of information with the potential to improve mapping accuracy are the 3' ends of longer reads, which are often ignored because they contain a higher frequency of errors, and the base-call quality scores.
the quality scores describe the confidence of bases in each read.
sequencing quality scores, introduced in the phred algorithm  <cit> , assign a probability to the four possible nucleotides for each sequenced base.
the solexa analysis pipeline, for example, includes a program called bustard to calculate quality scores.
because the bases with lower quality scores are more likely to be sequencing errors, any potential mapping for a read should be penalized less for mismatching at positions with lower scores.
the quality scores are especially important for mapping longer reads, since the 3' ends of longer reads are known to have a higher frequency of sequencing errors.
to investigate whether these sources of information can be used to improve accuracy when mapping reads, we developed the rmap tool, which can map reads having a wide range of lengths and allows base-call quality scores to determine which positions in each read are more important when mapping.
the only requirement on the base-call quality scores for use in rmap is that they increase monotonically with the inverse of the error probability for a particular base call.
our results indicate that significant gains in mapping performance can be achieved by considering the information in 3' ends of longer reads, and appropriately using the quality scores.
RESULTS
the mapping criteria
we designed rmap to use two different mapping criteria, both based on approximate matching of the read and the reference genome.
the first criterion is a simple count of mismatches between a read and the aligned genomic segment.
under this criterion, any unknown nucleotides in the reference genome  will induce a mismatch with any nucleotide.
uncalled positions, where the sequencing was unable to determine the nucleotide, also induce a mismatch.
for a fixed read length, by allowing a greater number of mismatches, more reads can be mapped to reference genome.
we refer to this simple mismatch criterion as rmapm .
the second criterion, also based on mismatch-counts, makes use of the base-call quality scores.
a cutoff for the base-call quality score is used to designate positions as either high-quality  or low-quality , depending on whether the quality score of the highest-scoring base at that position exceeds the cutoff.
low-quality positions always induce a match .
to prevent the possibility of trivial matches, a quality control step eliminates reads with too many low-quality positions.
as with the first criterion, mapping accuracy can be controlled by manipulating the number of allowed mismatches when mapping.
but for the second criterion, manipulating the quality-score cutoff provides another means of adjusting sensitivity and specificity, and allows positions to contribute when they are of high-quality, but not be penalized if they are low-quality.
we refer to this criterion as rmapq .
evaluating mapping accuracy
measuring mapping accuracy
in measuring mapping accuracy, we want to quantify both sensitivity and specificity by using reads sequenced from dna samples from selected genomic regions instead of the entire genome.
by mapping those reads to the genome, we can evaluate how accurately they are mapped to the target region.
however, there are theoretical and practical limits to how well these can be measured.
inability to map a read correctly can be attributed to sequencing errors , to variation between the sampled genome and the reference genome, or can result from ambiguities caused by repeats in the reference genome.
these diverse sources of error make it difficult to measure accuracy in terms of traditional sensitivity and specificity.
ambiguous reads, under a given mapping criterion, are reads that map to more than one location in the reference genome.
reads that map to a single location are called uniquely mappable  reads.
all reads that are not mapped uniquely to some location in the reference genome are said to be unmappable .
we define target region coverage  as the number of bases in the target region covered by at least one mappable read, divided by the total number of bases in the target region.
in order to compare the coverage values for different read lengths, we use only the first base of each read to represent that read.
by counting bases covered, rather than number of reads that map to the target region, greater target region coverage is achieved when the reads map uniformly in the target region.
we define mapping selectivity as the number of mappable reads that map inside the target region, divided by the total number of mappable reads.
a read is said to map inside the target region if any part of the read overlaps the target region.
the selectivity shows how well the mapping criterion places mappable reads inside the target region.
in this study, when we refer to mapping accuracy, we are referring to both coverage and mapping selectivity .
evaluation data
we used the data from samples of two bacs provided for sequencer validation by solexa, which covers  <dig> kb of the chromosome  <dig> mhc region in an a1-b8-dr <dig> alternate haplotype assembly based on sequence data from the cox library  <cit> .
we chose one lane of reads sequenced by the 1g sequencers at the cshl genome center.
the total number of raw  <dig> nt reads in this data set is  <dig>  million, with the quality score of each base ranging from - <dig> to  <dig> .
as reference genome we used hg <dig>  all chromosomes except chr <dig>  which we replaced entirely with chr6_cox_hap <dig>  the a1-b8-dr <dig> alternate haplotype assembly.
we chose to include this alternate haplotype in the reference because it is the origin of the bac that was sequenced.
we excluded the ordinary chr <dig> because it has high similarity with chr6_cox_hap <dig>  and including both of these would have resulted in a high proportion of reads mapping ambiguously to chr <dig> and chr6_cox_hap <dig> .
evaluation procedure
to investigate how information is distributed within the reads, we ran rmap on all reads with lengths ranging from 25– <dig> nt, allowing mismatches in the range of  <dig> to  <dig>  and using both the rmapm and rmapq criteria.
for the rmapq criterion, we chose { <dig>   <dig>   <dig>   <dig>   <dig>  24} as the set of quality-score cutoffs to evaluate.
reads with fewer than  <dig> contiguous hq bases  were considered unmappable and removed from consideration, as the algorithm requires a minimum number of high-quality bases for efficiency .
any read lacking  <dig> consecutive high-quality bases would likely have a very high overall amount of error.
mapping longer reads with more mismatches increases accuracy
the solexa sequencer can produce reads of more than  <dig> bases, and longer reads contain more sequence information.
although it is known that the quality of sequenced bases in reads decreases toward the 3' end of the read, especially as read length increases, it remains to be shown how much useful information may still exist in bases at the 3' ends of longer reads.
making use of any additional bases is only expected to improve mapping accuracy if the additional bases contain information of sufficient quality.
when the bac reads were mapped to the human genome using the rmapm criterion, with length from 25– <dig> nt, and different number of allowed mismatches, we found that there is generally a great deal of information in 3' end bases up to  <dig> nt.
these results are presented in figure  <dig> and supplementary table  <dig>  .
the bac coverage always increased with length of mapped reads, except when only one or zero mismatches are allowed.
the mapping selectivity decreases monotonically with read length when zero or at most one mismatch is allowed.
when multiple mismatches are allowed, the mapping selectivity first increases with read length, then decreases slightly.
taking the mean of these two measures as overall mapping accuracy, we see that the combined target region coverage and mapping selectivity is maximized when read length is  <dig> nt and up to  <dig> mismatches are allowed.
comparing read lengths between  <dig> nt and  <dig> nt, the mapping selectivity increases  <dig> % while the bac region coverage increases  <dig> %.
by extrapolating our results, even greater improvements are expected as read lengths increase beyond  <dig> nt.
using quality score information increases accuracy
we tested different cutoffs for defining high quality bases to estimate the ability of the rmapq criterion to most effectively use the quality information in the reads.
we achieve better mapping performance in both bac coverage and mapping selectivity when using longer reads with high quality score cutoff of  <dig> to  <dig> than without using quality score filtering .
this is due to a larger number of reads being mapped unambiguously to the genome by the rmapq criterion, and a larger number of those being mapped to the target  regions.
the best mapping accuracy was achieved when reads were of length  <dig> nt, a quality-score cutoff of  <dig> was used, and when at most one mismatch was allowed.
for these settings, the mapping selectivity further improved almost 3% with the same bac coverage, compared with the best accuracy of rmap using the rmapm criterion.
these results demonstrate that the way in which the quality scores are incorporated into rmap results in improved accuracy.
the reads derived from the bacs provided by solexa for the purpose of validation are of very high quality.
for a given  <dig> kb bac region, most of these  <dig>  million reads can be mapped almost perfectly.
this introduces a saturation effect, leaving limited space for improvement, as highly accurate mapping is achieved easily.
rmap still makes significant improvements within this narrow range.
in many applications, without this ceiling effect, the improvement is expected to be even more pronounced.
we also compared the mapping accuracy of rmap with another available, but presently unpublished, method called maq  <cit> .
using default parameter setting, we ran maq to map the set of  <dig>  m reads from the bac.
although maq had speed and memory usage similar to rmap, we found maq to have lower mapping accuracy  when allowing at most  <dig> mismatches.
discussion
widely-accessible second-generation sequencing technologies promise to revolutionize many areas of bio-medical research.
in addition to de novo sequencing of new species, these technologies make targeted re-sequencing a reality.
re-sequencing will provide more accurate means of interrogating the oligo-nucleotide content of samples, and of identifying important genome variations, such as disease-related snps.
mapping reads to genomes is a critical step in re-sequencing data analysis, and both the algorithmic and software technology must keep pace with surging advances in the throughput of sequencing instruments.
in order to maximize the use of available information in mapping solexa reads, we developed the rmap tool, which incorporates base-call quality scores to improve accuracy.
rmap responds to an urgent need for such an algorithm by providing both the accuracy to handle emerging mapping tasks.
our results in applying rmap have shown that more reads can be mapped into the target regions when using the rmapm criterion to map longer reads and allow more mismatches.
we have also shown that the way in which quality scores are used in rmap  significantly increases both coverage and mapping selectivity.
although second-generation sequencing technology is currently producing many important results, there is still little understanding of how this technology should behave with respect to sequencing errors and what are the general properties of typical re-sequencing data sets.
as more knowledge accumulates about typical results from this new sequencing technology, more information can be incorporated into algorithms for mapping reads and other associated analysis tasks.
in theory we could move toward an ideal mapping criterion by predictive modeling, where a model would be trained to identify the location from which each read was derived.
although the best mapping criteria may not be amenable to high-throughput computation, some approximation of those criteria could be developed.
cross-validation and the use of a wide range of data sets could be used to ensure that the trained criteria are sufficiently general.
in practice such a procedure would require high-quality training data, and an extreme amount of computing time to train and evaluate such models.
rmap does not consider insertions or deletions , which are potentially important in certain sequencing applications .
the straight-forward strategy for handling indels is to extend initial seed matches using a smith-waterman-style alignment, as is commonly done in database search programs like blast  <cit> .
for short reads, with length ≤  <dig> bases, providing this greater flexibility will require careful investigation into scoring the alignments, because using simple scoring schemes for indels may result in higher rate of false-positive mappings and ambiguities.
in addition, because indels have nothing analogous to the base-call quality scores, it will be more difficult to distinguish errors from real genotypic variations during the mapping stage of analysis.
we have observed that data sets containing a high proportion of low-complexity reads, such as single nucleotide repeats, cause a decrease in execution speed of rmap.
this is explained by such low complexity seeds resulting in a large increase in the number of full comparisons , that are induced by these highly-common seeds.
the exclusion strategy used in rmap is also heavily used by popular programs for searching sequence databases.
research effort toward improving the structure of seeds used in these algorithms has also led to improvements in the algorithms themselves  <cit> , and similar improvements might be observed by developing seeds specifically for the short read mapping problem.
the filtration strategy used in rmap can also be extended to "multiple filtration", as described by  <cit> .
this uses multiple criteria for excluding possible mappings, and would result in the algorithm performing fewer full comparisons between reads and genomic regions.
unlike many other applications of approximate matching, the extreme volume of data that must be mapped means that algorithms for mapping reads must be conscious of the memory required.
in the framework of rmap, the most straight-forward implementation of multiple filtration would use larger  or additional hash-tables, which generally requires a great deal more memory.
for short reads, as produced by the solexa sequencer, the full comparisons can be computed very fast, so another concern is keeping the work required to implement more powerful filtration less than that required to do the extra comparisons.
aside from further processing on the set of reads, efficiency improvements can be gained by indexing the reference genome.
this is the strategy used by sxoligosearch, through the proprietary synabase data structure.
because the reference genome will likely remain static for many mapping tasks , this indexing can be done off-line, so the time required to build such an index is not important.
the usual problem with indexing the genome is that such index structures are often several times larger than the genome itself.
the more information built into the index to facilitate searches, the larger the index must be.
as stated above, hardware currently available for mapping usually has a few gb of memory, and this restricts the kinds of indexing that can be done on entire genomes.
to make genome index structures practical, they must either be small, or support on-line processing so that only a small portion of the structure must reside in memory at once.
it is likely that special purpose data structures for indexing the genome can be developed to fit these requirements, and greatly improve efficiency of mapping reads.
specifically with respect to rmap, adding any means of eliminating exact genomic repeats during mapping will greatly improve efficiency.
CONCLUSIONS
our results indicate that significant gains in solexa read mapping performance can be achieved by considering the information in 3' ends of longer reads, and appropriately using the base-call quality scores.
the rmap tool we have developed will enable researchers to effiectively exploit this information in targeted re-sequencing projects.
availability and requirements
the rmap tool is freely available for downloading at .

 the main technological advances that accompanied the genomic and post-genomic eras are high-throughput sequencing and hybridization microarrays. sequencing technology enabled scientists to obtain the full genomic sequence for many species, including the human and many model organisms. sequencing technology is also being used to selectively re-sequence the human genome to detect genome variations such as single nucleotide polymorphisms or large-scale structural variations. because understanding these variations can immediately impact medical sciences, making sequencing more efficient and accessible is imperative. however, traditional methodologies used to sequence the first mammalian genomes remain expensive, time consuming, and labor intensive. oligonucleotide microarray technology, used to interrogate the rna or dna content of a sample, has emerged as a widely accessible and effective tool for studying gene expression or detecting protein-dna interactions. microarrays have also been used to detect genome variations, like snps or structural variations. array-based hybridization also has limitations. for example, probes can behave in highly non-uniform ways, and the effects of cross-hybridization and resolution limits are poorly understood. although significant research efforts are focused on these problems, they remain inherent in all hybridization-based methods. the new sequencing technology referred to as "second-generation" shows promise to eliminate many of the problems associated with traditional sequencing technology and also those with oligonucleotide microarray technology. second-generation sequencers are able to sequence more quickly and at lower cost in terms of both money and labor. new sequencing technologies are developed to sequence at greater depth, meaning that a clone can be sequenced from a sample even when that clone exists at very low abundance. several second-generation technologies have been developed using diverse methods. two recent second-generation sequencers are from life sciences and solexa. the sequencers use an emulsion method for dna amplification and a pyrosequencing protocol for sequencing by synthesis at picolitre-scale volumes. current sequencers can produce 25– m nt of sequence in a single run, in the form of reads with length up to nt, enabling this technology to be used for de-novo sequencing in addition to re-sequencing. solexa/illumina 1g sequencers also use sequencing by synthesis, with dna amplified on the surface of a flow cell, resulting in a random array of dense clusters. the solexa technology is faster and cheaper than that used in sequencers, producing 1g nucleotides of sequence in one run but producing much shorter reads. each individual read is roughly to bases in length. the solexa sequencing technology has recently started producing breakthrough results. research described in and employed solexa sequencers to obtain high-resolution genomic maps of several histone modifications, as well as localization data for the dna-binding proteins. effectiveness of chip-sequencing has also been demonstrated by, who used solexa sequencing to obtain locations of stat binding sites in hela s cells before and following ifn-γ stimulation. mapping reads from the solexa sequencer presents an obvious algorithmic challenge: tens of millions of reads must be mapped to a large genome in a reasonable amount of time. strong efforts to design short-read mapping algorithms have resulted in methods that are effective in particular contexts. the mapping algorithm implemented as part of the solexa analysis pipeline is named eland. eland is optimized to map very short reads, with length at most nt, and ignores the additional bases when the sequenced reads are longer. eland also only allows at most two mismatches between the read and the genomic sequence to which it maps, which will clearly be too few for longer reads. despite these restrictions, eland remains very useful for many mapping tasks because it is extremely efficient. the sxoligosearch algorithm can quickly map reads of varying length, using different criteria, with performance depending on both the read length and mapping criteria. the performance of sxoligosearch depends on use of the proprietary synabase data structure. this data structure is a heavily compressed and annotated index for the reference genome that retains all non-redundant information. the gains in mapping speed by using such a data structure come at a cost in terms of the memory required for the synabase data structure. extreme memory requirements of the data structure makes sxoligosearch unsuitable for use on hardware available in most labs. in most re-sequencing applications accuracy of mapping is a primary concern. we must know with great accuracy what part of the genome was actually sequenced. there are several reasons why it might be difficult to determine the location in the reference genome from which a read was derived, or even if a read was derived from the reference genome. these include problems with the experiments, such as sequencing errors or sample contamination. mapping is also made more difficult by repeats in the genome, and by polymorphisms. while mapping algorithms cannot be expected to be robust to all such problems, effort should be made to make the algorithms as robust as possible. two sources of information with the potential to improve mapping accuracy are the 3' ends of longer reads, which are often ignored because they contain a higher frequency of errors, and the base-call quality scores. the quality scores describe the confidence of bases in each read. sequencing quality scores, introduced in the phred algorithm, assign a probability to the four possible nucleotides for each sequenced base. the solexa analysis pipeline, for example, includes a program called bustard to calculate quality scores. because the bases with lower quality scores are more likely to be sequencing errors, any potential mapping for a read should be penalized less for mismatching at positions with lower scores. the quality scores are especially important for mapping longer reads, since the 3' ends of longer reads are known to have a higher frequency of sequencing errors. to investigate whether these sources of information can be used to improve accuracy when mapping reads, we developed the rmap tool, which can map reads having a wide range of lengths and allows base-call quality scores to determine which positions in each read are more important when mapping. the only requirement on the base-call quality scores for use in rmap is that they increase monotonically with the inverse of the error probability for a particular base call. our results indicate that significant gains in mapping performance can be achieved by considering the information in 3' ends of longer reads, and appropriately using the quality scores.  the mapping criteria we designed rmap to use two different mapping criteria, both based on approximate matching of the read and the reference genome. the first criterion is a simple count of mismatches between a read and the aligned genomic segment. under this criterion, any unknown nucleotides in the reference genome will induce a mismatch with any nucleotide. uncalled positions, where the sequencing was unable to determine the nucleotide, also induce a mismatch. for a fixed read length, by allowing a greater number of mismatches, more reads can be mapped to reference genome. we refer to this simple mismatch criterion as rmapm. the second criterion, also based on mismatch-counts, makes use of the base-call quality scores. a cutoff for the base-call quality score is used to designate positions as either high-quality or low-quality, depending on whether the quality score of the highest-scoring base at that position exceeds the cutoff. low-quality positions always induce a match. to prevent the possibility of trivial matches, a quality control step eliminates reads with too many low-quality positions. as with the first criterion, mapping accuracy can be controlled by manipulating the number of allowed mismatches when mapping. but for the second criterion, manipulating the quality-score cutoff provides another means of adjusting sensitivity and specificity, and allows positions to contribute when they are of high-quality, but not be penalized if they are low-quality. we refer to this criterion as rmapq. evaluating mapping accuracy measuring mapping accuracy in measuring mapping accuracy, we want to quantify both sensitivity and specificity by using reads sequenced from dna samples from selected genomic regions instead of the entire genome. by mapping those reads to the genome, we can evaluate how accurately they are mapped to the target region. however, there are theoretical and practical limits to how well these can be measured. inability to map a read correctly can be attributed to sequencing errors, to variation between the sampled genome and the reference genome, or can result from ambiguities caused by repeats in the reference genome. these diverse sources of error make it difficult to measure accuracy in terms of traditional sensitivity and specificity. ambiguous reads, under a given mapping criterion, are reads that map to more than one location in the reference genome. reads that map to a single location are called uniquely mappable reads. all reads that are not mapped uniquely to some location in the reference genome are said to be unmappable. we define target region coverage as the number of bases in the target region covered by at least one mappable read, divided by the total number of bases in the target region. in order to compare the coverage values for different read lengths, we use only the first base of each read to represent that read. by counting bases covered, rather than number of reads that map to the target region, greater target region coverage is achieved when the reads map uniformly in the target region. we define mapping selectivity as the number of mappable reads that map inside the target region, divided by the total number of mappable reads. a read is said to map inside the target region if any part of the read overlaps the target region. the selectivity shows how well the mapping criterion places mappable reads inside the target region. in this study, when we refer to mapping accuracy, we are referring to both coverage and mapping selectivity. evaluation data we used the data from samples of two bacs provided for sequencer validation by solexa, which covers kb of the chromosome mhc region in an a1-b8-dr alternate haplotype assembly based on sequence data from the cox library. we chose one lane of reads sequenced by the 1g sequencers at the cshl genome center. the total number of raw nt reads in this data set is million, with the quality score of each base ranging from - to. as reference genome we used hg all chromosomes except chr which we replaced entirely with chr6_cox_hap the a1-b8-dr alternate haplotype assembly. we chose to include this alternate haplotype in the reference because it is the origin of the bac that was sequenced. we excluded the ordinary chr because it has high similarity with chr6_cox_hap and including both of these would have resulted in a high proportion of reads mapping ambiguously to chr and chr6_cox_hap. evaluation procedure to investigate how information is distributed within the reads, we ran rmap on all reads with lengths ranging from 25– nt, allowing mismatches in the range of to and using both the rmapm and rmapq criteria. for the rmapq criterion, we chose { 24} as the set of quality-score cutoffs to evaluate. reads with fewer than contiguous hq bases were considered unmappable and removed from consideration, as the algorithm requires a minimum number of high-quality bases for efficiency. any read lacking consecutive high-quality bases would likely have a very high overall amount of error. mapping longer reads with more mismatches increases accuracy the solexa sequencer can produce reads of more than bases, and longer reads contain more sequence information. although it is known that the quality of sequenced bases in reads decreases toward the 3' end of the read, especially as read length increases, it remains to be shown how much useful information may still exist in bases at the 3' ends of longer reads. making use of any additional bases is only expected to improve mapping accuracy if the additional bases contain information of sufficient quality. when the bac reads were mapped to the human genome using the rmapm criterion, with length from 25– nt, and different number of allowed mismatches, we found that there is generally a great deal of information in 3' end bases up to nt. these results are presented in figure and supplementary table. the bac coverage always increased with length of mapped reads, except when only one or zero mismatches are allowed. the mapping selectivity decreases monotonically with read length when zero or at most one mismatch is allowed. when multiple mismatches are allowed, the mapping selectivity first increases with read length, then decreases slightly. taking the mean of these two measures as overall mapping accuracy, we see that the combined target region coverage and mapping selectivity is maximized when read length is nt and up to mismatches are allowed. comparing read lengths between nt and nt, the mapping selectivity increases % while the bac region coverage increases %. by extrapolating our results, even greater improvements are expected as read lengths increase beyond nt. using quality score information increases accuracy we tested different cutoffs for defining high quality bases to estimate the ability of the rmapq criterion to most effectively use the quality information in the reads. we achieve better mapping performance in both bac coverage and mapping selectivity when using longer reads with high quality score cutoff of to than without using quality score filtering. this is due to a larger number of reads being mapped unambiguously to the genome by the rmapq criterion, and a larger number of those being mapped to the target regions. the best mapping accuracy was achieved when reads were of length nt, a quality-score cutoff of was used, and when at most one mismatch was allowed. for these settings, the mapping selectivity further improved almost 3% with the same bac coverage, compared with the best accuracy of rmap using the rmapm criterion. these results demonstrate that the way in which the quality scores are incorporated into rmap results in improved accuracy. the reads derived from the bacs provided by solexa for the purpose of validation are of very high quality. for a given kb bac region, most of these million reads can be mapped almost perfectly. this introduces a saturation effect, leaving limited space for improvement, as highly accurate mapping is achieved easily. rmap still makes significant improvements within this narrow range. in many applications, without this ceiling effect, the improvement is expected to be even more pronounced. we also compared the mapping accuracy of rmap with another available, but presently unpublished, method called maq. using default parameter setting, we ran maq to map the set of m reads from the bac. although maq had speed and memory usage similar to rmap, we found maq to have lower mapping accuracy when allowing at most mismatches. discussion widely-accessible second-generation sequencing technologies promise to revolutionize many areas of bio-medical research. in addition to de novo sequencing of new species, these technologies make targeted re-sequencing a reality. re-sequencing will provide more accurate means of interrogating the oligo-nucleotide content of samples, and of identifying important genome variations, such as disease-related snps. mapping reads to genomes is a critical step in re-sequencing data analysis, and both the algorithmic and software technology must keep pace with surging advances in the throughput of sequencing instruments. in order to maximize the use of available information in mapping solexa reads, we developed the rmap tool, which incorporates base-call quality scores to improve accuracy. rmap responds to an urgent need for such an algorithm by providing both the accuracy to handle emerging mapping tasks. our results in applying rmap have shown that more reads can be mapped into the target regions when using the rmapm criterion to map longer reads and allow more mismatches. we have also shown that the way in which quality scores are used in rmap significantly increases both coverage and mapping selectivity. although second-generation sequencing technology is currently producing many important results, there is still little understanding of how this technology should behave with respect to sequencing errors and what are the general properties of typical re-sequencing data sets. as more knowledge accumulates about typical results from this new sequencing technology, more information can be incorporated into algorithms for mapping reads and other associated analysis tasks. in theory we could move toward an ideal mapping criterion by predictive modeling, where a model would be trained to identify the location from which each read was derived. although the best mapping criteria may not be amenable to high-throughput computation, some approximation of those criteria could be developed. cross-validation and the use of a wide range of data sets could be used to ensure that the trained criteria are sufficiently general. in practice such a procedure would require high-quality training data, and an extreme amount of computing time to train and evaluate such models. rmap does not consider insertions or deletions, which are potentially important in certain sequencing applications. the straight-forward strategy for handling indels is to extend initial seed matches using a smith-waterman-style alignment, as is commonly done in database search programs like blast. for short reads, with length ≤ bases, providing this greater flexibility will require careful investigation into scoring the alignments, because using simple scoring schemes for indels may result in higher rate of false-positive mappings and ambiguities. in addition, because indels have nothing analogous to the base-call quality scores, it will be more difficult to distinguish errors from real genotypic variations during the mapping stage of analysis. we have observed that data sets containing a high proportion of low-complexity reads, such as single nucleotide repeats, cause a decrease in execution speed of rmap. this is explained by such low complexity seeds resulting in a large increase in the number of full comparisons, that are induced by these highly-common seeds. the exclusion strategy used in rmap is also heavily used by popular programs for searching sequence databases. research effort toward improving the structure of seeds used in these algorithms has also led to improvements in the algorithms themselves, and similar improvements might be observed by developing seeds specifically for the short read mapping problem. the filtration strategy used in rmap can also be extended to "multiple filtration", as described by. this uses multiple criteria for excluding possible mappings, and would result in the algorithm performing fewer full comparisons between reads and genomic regions. unlike many other applications of approximate matching, the extreme volume of data that must be mapped means that algorithms for mapping reads must be conscious of the memory required. in the framework of rmap, the most straight-forward implementation of multiple filtration would use larger or additional hash-tables, which generally requires a great deal more memory. for short reads, as produced by the solexa sequencer, the full comparisons can be computed very fast, so another concern is keeping the work required to implement more powerful filtration less than that required to do the extra comparisons. aside from further processing on the set of reads, efficiency improvements can be gained by indexing the reference genome. this is the strategy used by sxoligosearch, through the proprietary synabase data structure. because the reference genome will likely remain static for many mapping tasks, this indexing can be done off-line, so the time required to build such an index is not important. the usual problem with indexing the genome is that such index structures are often several times larger than the genome itself. the more information built into the index to facilitate searches, the larger the index must be. as stated above, hardware currently available for mapping usually has a few gb of memory, and this restricts the kinds of indexing that can be done on entire genomes. to make genome index structures practical, they must either be small, or support on-line processing so that only a small portion of the structure must reside in memory at once. it is likely that special purpose data structures for indexing the genome can be developed to fit these requirements, and greatly improve efficiency of mapping reads. specifically with respect to rmap, adding any means of eliminating exact genomic repeats during mapping will greatly improve efficiency.  our results indicate that significant gains in solexa read mapping performance can be achieved by considering the information in 3' ends of longer reads, and appropriately using the base-call quality scores. the rmap tool we have developed will enable researchers to effiectively exploit this information in targeted re-sequencing projects. availability and requirements the rmap tool is freely available for downloading at. 
to investigate whether these sources of information can be used to improve accuracy when mapping reads, we developed the rmap tool, which can map reads having a wide range of lengths and allows base-call quality scores to determine which positions in each read are more important when mapping. the rmap tool we have developed will enable researchers to effectively exploit this information in targeted re-sequencing projects. our results indicate that significant gains in solexa read mapping performance can be achieved by considering the information in 3' ends of longer reads, and appropriately using the base-call quality scores. accurately mapping the reads back to a reference genome is a critical task in almost all applications. two sources of information that are often ignored when mapping reads from the solexa technology are the 3' ends of longer reads, which contain a much higher frequency of sequencing errors, and the base-call quality scores. we applied rmap to analyze data re-sequenced from two human bac regions for varying read lengths, and varying criteria for use of quality scores. the increase in speed, sensitivity and availability of sequencing technology brings demand for advances in computational technology to perform associated analysis tasks.
 second-generation sequencing has the potential to revolutionize genomics and impact all areas of biomedical science. new technologies will make re-sequencing widely available for such applications as identifying genome variations or interrogating the oligonucleotide content of a large sample. the increase in speed, sensitivity and availability of sequencing technology brings demand for advances in computational technology to perform associated analysis tasks. the solexa/illumina 1g sequencer can produce tens of millions of reads, ranging in length from ~25– nt, in a single experiment. accurately mapping the reads back to a reference genome is a critical task in almost all applications. two sources of information that are often ignored when mapping reads from the solexa technology are the 3' ends of longer reads, which contain a much higher frequency of sequencing errors, and the base-call quality scores.  to investigate whether these sources of information can be used to improve accuracy when mapping reads, we developed the rmap tool, which can map reads having a wide range of lengths and allows base-call quality scores to determine which positions in each read are more important when mapping. we applied rmap to analyze data re-sequenced from two human bac regions for varying read lengths, and varying criteria for use of quality scores. rmap is freely available for downloading at.  our results indicate that significant gains in solexa read mapping performance can be achieved by considering the information in 3' ends of longer reads, and appropriately using the base-call quality scores. the rmap tool we have developed will enable researchers to effectively exploit this information in targeted re-sequencing projects. 