Самый быстрый способ найти позиции несоответствия между двумя строками одинаковой длины
У меня есть миллионы пар строк одинаковой длины, которые я хочу сравнить и найти позицию, в которой они не совпадают.
Например для каждого $str1 а также $str2 мы хотим найти несоответствующую позицию с $str_source:
$str_source = "ATTCCGGG";
$str1 = "ATTGCGGG"; # 1 mismatch with Str1 at position 3 (0-based)
$str2 = "ATACCGGC"; # 2 mismatches with source at position 2 and 7
Есть ли быстрый способ сделать это. В настоящее время у меня есть метод стиля C, который я зацикливаю каждую позицию в обеих строках, используя функцию 'substr'. Но этот подход ужасно медленный.
my @mism_pos;
for $i (0 .. length($str_source)) {
$source_base = substr($str_source,$i,1);
$str_base = substr($str2,$i,$1);
if ($source_base ne $str_base) {
push @mism_pos,$i;
}
}
9 ответов
Инлайн:: C
Вычисления просты, сделайте это с Inline:: C (для документации http://search.cpan.org/perldoc/Inline::C-Cookbook и perldoc Inline::C):
use Inline C => << '...';
void find_diffs(char* x, char* y) {
int i;
Inline_Stack_Vars;
Inline_Stack_Reset;
for(i=0; x[i] && y[i]; ++i) {
if(x[i] != y[i]) {
Inline_Stack_Push(sv_2mortal(newSViv(i)));
}
}
Inline_Stack_Done;
}
...
@diffs= find_diffs("ATTCCGGG","ATTGCGGG"); print "@diffs\n";
@diffs= find_diffs("ATTCCGGG","ATACCGGC"); print "@diffs\n";
Вот вывод этого скрипта:
> script.pl
3
2 7
PDL
Если вы хотите быстро обрабатывать много данных в Perl, изучите PDL ( Документация):
use PDL;
use PDL::Char;
$PDL::SHARE=$PDL::SHARE; # keep stray warning quiet
my $source=PDL::Char->new("ATTCCGGG");
for my $str ( "ATTGCGGG", "ATACCGGC") {
my $match =PDL::Char->new($str);
my @diff=which($match!=$source)->list;
print "@diff\n";
}
(Тот же вывод, что и в первом скрипте.)
Примечания: Я очень счастливо использовал PDL при обработке геномных данных. Вместе с отображением в памяти доступа к данным, хранящимся на диске, можно быстро обрабатывать огромные объемы данных: вся обработка выполняется в высоко оптимизированных циклах C. Кроме того, вы можете легко получить доступ к тем же данным через Inline:: C для любых функций, отсутствующих в PDL.
Однако обратите внимание, что создание одного вектора PDL довольно медленное (постоянное время, это приемлемо для больших структур данных). Итак, вы скорее хотите создать один большой объект PDL со всеми вашими входными данными за один раз, чем зацикливаться на отдельных элементах данных.
Они похожи на генные последовательности. Если все строки состоят из 8 символов, а область возможных кодов - ( A, C, G, T), вы можете как-то преобразовать данные перед их обработкой. Это даст вам только 65536 возможных строк, так что вы можете специализировать свою реализацию.
Например, вы пишете метод, который принимает 8-символьную строку и отображает ее в целое число. Запомните это, чтобы операция была быстрой. Затем напишите функцию сравнения, которая дает два целых числа, рассказывает, как они различаются. Вы бы назвали это в подходящей циклической конструкции с числовым тестом на равенство unless ( $a != $b ) перед вызовом сравнения - короткое замыкание на идентичные коды, если хотите.
Вот сценарий бенчмаркинга, чтобы выяснить, есть ли различия в скорости различных подходов. Просто имейте в виду, что при первом запуске сценария, использующего Inline:: C, при запуске компилятора C будет происходить задержка и т. Д. Итак, запустите сценарий один раз, а затем выполните тестирование.
#!/usr/bin/perl
use strict;
use warnings;
use Benchmark qw( cmpthese );
my ($copies) = @ARGV;
$copies ||= 1;
my $x = 'ATTCCGGG' x $copies;
my $y = 'ATTGCGGG' x $copies;
my $z = 'ATACCGGC' x $copies;
sub wrapper {
my ($func, @args) = @_;
for my $s (@args) {
my $differences = $func->($x, $s);
# just trying to ensure results are not discarded
if ( @$differences == 0 ) {
print "There is no difference\n";
}
}
return;
}
cmpthese -5, {
explode => sub { wrapper(\&where_do_they_differ, $y, $z) },
mism_pos => sub { wrapper(\&mism_pos, $y, $z) },
inline_c => sub {
wrapper(\&i_dont_know_how_to_do_stuff_with_inline_c, $y, $z) },
};
sub where_do_they_differ {
my ($str1, $str2) = @_;
my @str1 = split //, $str1;
my @str2 = split //, $str2;
[ map {$str1[$_] eq $str2[$_] ? () : $_} 0 .. length($str1) - 1 ];
}
sub mism_pos {
my ($str1, $str2) = @_;
my @mism_pos;
for my $i (0 .. length($str1) - 1) {
if (substr($str1, $i, 1) ne substr($str2, $i, 1) ) {
push @mism_pos, $i;
}
}
return \@mism_pos;
}
sub i_dont_know_how_to_do_stuff_with_inline_c {
[ find_diffs(@_) ];
}
use Inline C => << 'EOC';
void find_diffs(char* x, char* y) {
int i;
Inline_Stack_Vars;
Inline_Stack_Reset;
for(i=0; x[i] && y[i]; ++i) {
if(x[i] != y[i]) {
Inline_Stack_Push(sv_2mortal(newSViv(i)));
}
}
Inline_Stack_Done;
}
EOC
Результаты (с использованием VC++ 9 в Windows XP с AS Perl 5.10.1) с $copies = 1:
Оцените взорваться mism_pos inline_c взорваться 15475/ с - -64% -84% mism_pos 43196/s 179% --56% inline_c 98378/ с 536% 128% -
Результаты с $copies = 100:
Оцените взорваться mism_pos inline_c взорваться 160 / с - -86% -99% mism_pos 1106/s 593% --90% inline_c 10808 / с 6667% 877% -
Похоже, это может быть важной частью вашего приложения. В этом случае вы можете рассмотреть возможность написания метода расширения C для сравнения.
Perl предоставляет механизм расширения XS, который делает это достаточно простым.
Самый быстрый способ сравнить строки, чтобы найти различия, - это XOR каждого байта из них вместе, а затем проверить на ноль. Если бы мне пришлось это сделать, я бы просто написал программу на C, чтобы выполнить разную работу, а не написал расширение C для Perl, тогда я бы запустил свою программу на C как подпроцесс Perl. Точный алгоритм будет зависеть от длины строк и количества данных. Однако это не займет более 100 строк C. Фактически, если вы хотите максимизировать скорость, программа на XOR байтах строк фиксированной длины и проверка на ноль могут быть написаны на ассемблере.
Вы делаете 2 вызова к substr для каждого сравнения символов, что, вероятно, замедляет вас.
Несколько оптимизаций я бы сделал
@source = split //,$str_source #split first rather than substr
@base = split //, $str_base
for $i (0 .. length($str_source)) {
$mism_pos{$1} = 1 if ($source[$i] ne $base); #hashing is faster than array push
}
return keys $mism_pos
Я не знаю, насколько он эффективен, но вы всегда можете сделать xor для двух подходящих строк и найти индекс первого несоответствия.
#! /usr/bin/env perl
use strict;
use warnings;
use 5.10.1;
my $str_source = "ATTCCGGG";
my $str1 = "ATTGCGGG";
my $str2 = "ATACCGGC";
my $str3 = "GTTCCGGG";
# this returns the index of all of the mismatches (zero based)
# it returns an empty list if the two strings match.
sub diff_index{
my($a,$b) = @_;
my $cmp = $a^$b;
my @cmp;
while( $cmp =~ /[^\0]/g ){ # match non-zero byte
push @cmp, pos($cmp) - 1;
}
return @cmp;
}
for my $str ( $str_source, $str1, $str2, $str3 ){
say '# "', $str, '"';
my @ret = diff_index $str_source, $str;
if( @ret ){
say '[ ', join( ', ', @ret), ' ]';
}else{
say '# match';
}
}
# "ATTCCGGG"
# match
# "ATTGCGGG"
[ 3 ]
# "ATACCGGC"
[ 2, 7 ]
# "GTTCCGGG"
[ 0 ]
Запуск его через B:: Concise показывает, что дорогостоящие операции процессора выполняются как одиночные коды операций. Это означает, что эти операции выполняются в C.
perl -MO=Concise,-exec,-compact,-src,diff_index test.pl |
perl -pE's/^[^#].*? \K([^\s]+)$/# $1/' # To fix highlighting bugs
main::diff_index:
# 15: my($a,$b) = @_;
1 <;> nextstate(main 53 test.pl:15) # v:%,*,&,$
2 <0> pushmark # s
3 <$> gv(*_) # s
4 <1> rv2av[t3] # lK/3
5 <0> pushmark # sRM*/128
6 <0> padsv[$a:53,58] # lRM*/LVINTRO
7 <0> padsv[$b:53,58] # lRM*/LVINTRO
8 <2> aassign[t4] # vKS
# 16: my $cmp = $a^$b;
9 <;> nextstate(main 54 test.pl:16) # v:%,*,&,$
a <0> padsv[$a:53,58] # s
b <0> padsv[$b:53,58] # s
c <2> bit_xor[t6] # sK <----- Single OP -----
d <0> padsv[$cmp:54,58] # sRM*/LVINTRO
e <2> sassign # vKS/2
# 18: my @cmp;
f <;> nextstate(main 55 test.pl:18) # v:%,*,&,{,$
g <0> padav[@cmp:55,58] # vM/LVINTRO
# 20: while( $cmp =~ /[^\0]/g ){ # match non-zero byte
h <;> nextstate(main 57 test.pl:20) # v:%,*,&,{,$
i <{> enterloop(next->r last->v redo->j) # v
s <0> padsv[$cmp:54,58] # s
t </> match(/"[^\\0]"/) # sKS/RTIME <----- Single OP -----
u <|> and(other->j) # vK/1
# 21: push @cmp, pos($cmp) - 1;
j <;> nextstate(main 56 test.pl:21) # v:%,*,&,$
k <0> pushmark # s
l <0> padav[@cmp:55,58] # lRM
m <0> padsv[$cmp:54,58] # sRM
n <1> pos[t8] # sK/1
o <$> const(IV 1) # s
p <2> subtract[t9] # sK/2
q <@> push[t10] # vK/2
r <0> unstack # v
goto # s
v <2> leaveloop # vK/2
# 24: return @cmp;
w <;> nextstate(main 58 test.pl:24) # v:%,*,&,{,$
x <0> pushmark # s
y <0> padav[@cmp:55,58]
z <@> return # K
10 <1> leavesub[1 ref] # K/REFC,1
Некоторые классические строки сравнения оптимизации:
оптимальное несоответствие - начать сравнение в конце строки поиска. например, поиск ABC в ABDABEABF, если вы сравниваете в начале, вы будете перемещаться по шаблону по одному символу за раз. Если вы будете искать с конца, вы сможете прыгнуть вдоль трех символов
плохая эвристика персонажа - выберите наименее часто встречающийся символ и выполните поиск по нему первым. например, на английском языке символ "z" встречается редко, и хорошие функции поиска по строке будут искать "лабиринт" и начинать сравнение с 3-го символа
Я собирался сказать, "напишите это на C" тоже.
Оказавшись там, вы можете использовать оптимизацию, например, сравнивая 4 символа одновременно (как 32-разрядные целые числа).
Или измените свое представление (4-х буквенное, верно?), Чтобы использовать 2-битное представление для базы (?), Чтобы вы могли сравнивать 16 символов одновременно.