I have an optimization problem where I have a degenerate nucleotide sequence I want to align against subsets of a reference genome (exons, specifically, to make the problem more tractable).
The sequence itself contains subregions of variable length.
The sequence looks like this, containing degenerate bases:
L C R 5' - SSS......SSS|NNN.....NNN|WWW......WWW - 3'
Per IUPAC notation,
N denotes any base, and
I am marking each block
R to broadly denote "left", "center", and "right" regions, oriented from 5' to 3' ("forward"-strand oriented).
The widths of the
R blocks can be from 20 to 25 bases long. The width of the
C block can be 14 to 16 bases long.
One programmatic way I could think to do this would be to brute-force a solution, generating all combinations of sequences of varied lengths,
blastn-aligning each of such combinations to the sequence for the exon, and then looking for the best alignments among such sequences.
The output would be the scores for regions from sequences that align, to be able to say, from point X to point Y of the genome, a sequence here aligns better than other candidate sequences.
To avoid reinventing tools, and avoid brute-force methods, I'd like to know if there is a tool or even a wrapper around tools, which does this, in case it is already written, or another approach that is more elegant/efficient.