我正在寻找一种在基因组序列中发现短串联重复序列的算法。
基本上,给定一个只能由 4 个字符“ATCG”组成的非常长的字符串,我需要找到 2-5 个字符长且彼此相邻的短重复。
例如: TACATGAGATCATGATGATGATGATGGAGCTGTGAGATC公司 会给 ATGATGATG 或 ATG 重复 3 次
该算法需要扩展到 100 万个字符的字符串,因此我试图尽可能接近线性运行时间。
我目前的算法: 由于重复可以是 2-5 个字符长,我逐个字符地检查字符串,看看第 N 个字符是否与第 N+X 个字符相同,X 是 2 到 5。每个 X 都有一个计数器来计算顺序匹配并在不匹配时重置,我们知道当 X = 计数器时是否有重复。然后可以手动检查后续重复。
最佳答案
你正在查看给你 O(n)
的每个字符,因为你在每个字符上比较下一个(最大)五个字符,这给你一个常量 c
:
var data = get_input();
var compare = { `A`, `T`, `G`, `A`, `T` } // or whatever
var MAX_LOOKAHEAD = compare.length
var n
var c
for(n = data_array.length; n < size; i++) { // Has runtime O(n)
for(c = 0; c < MAX_LOOKAHEAD; c++) { // Maximum O(c)
if( compare[c] != data[i+c] ) {
break;
} else {
report( "found match at position " + i )
}
}
}
很容易看出这运行了 O(n*c)
次。由于 c
非常小,因此可以将其忽略 - 我认为人们无法摆脱该常量 - 这导致总运行时间为 O(n)
。
好消息:
您可以通过并行化来加快速度。例如。您可以将其拆分为 k
间隔,并通过为它们提供适当的开始和结束索引让多个线程为您完成这项工作。这可以为您提供线性加速。
如果您这样做,请确保将交叉点视为特殊情况,因为如果您的间隔将匹配分成两部分,您可能会错过匹配。
例如n = 50000
:
4 个线程的分区:(n/10000) - 1 = 4
。第 5 个线程不会有太多工作要做,因为它只处理交叉点,这就是为什么我们不需要考虑它的(在我们的例子中很小)开销。
1 10000 20000 40000 50000
|-------------------|-------------------|-------------------|-------------------|
| <- thread 1 -> | <- thread 2 -> | <- thread 3 -> | <- thread 4 -> |
|---| |---| |---|
|___________________|___________________|
|
thread 5
这就是它的样子:
var data;
var compare = { `A`, `T`, `G`, `A`, `T` };
var MAX_LOOKAHEAD = compare.length;
thread_function(args[]) {
var from = args[0];
var to = args[1];
for(n = from ; n < to ; i++) {
for(c = 0; c < MAX_LOOKAHEAD; c++) {
if( compare[c] != data[i+c] ) {
break;
} else {
report( "found match at position " + i )
}
}
}
}
main() {
var data_size = 50000;
var thread_count = 4;
var interval_size = data_size / ( thread_count + 1) ;
var tid[]
// This loop starts the threads for us:
for( var i = 0; i < thread_count; i++ ) {
var args = { interval_size * i, (interval_size * i) + interval_size };
tid.add( create_thread( thread_function, args ) );
}
// And this handles the intersections:
for( var i = 1; i < thread_count - 1; i++ ) {
var args = { interval_size * i, (interval_size * i) + interval_size };
from = (interval_size * i) - compare.length + 1;
to = (interval_size * i) + compare.length - 1;
for(j = from; j < to ; j++) {
for(k = 0; k < MAX_LOOKAHEAD; k++) {
if( compare[k] != data[j+k] ) {
break;
} else {
report( "found match at position " + j )
}
}
}
}
wait_for_multiple_threads( tid );
}
关于string - 寻找连续重复序列的算法,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/23574738/