Алгоритм Ахо-Корасик

Вступление

В посте я постарался избежать сложных дефиниций и строгих матетематических доказательств, а некоторые вещи вообще понятны интуитивно. Алгоритм удобно разбивается взаимосвязные части, поэтому и уловить принцип его работы не должно составлять труда.

Начальное описание

Алгоритм Ахо-Корасик реализует эффективный поиск всех вхождений всех строк-образцов в заданную строку. Был разработан в 1975 году Альфредом Ахо и Маргарет Корасик.
Опишем формально условие задачи. На вход поступают несколько строк pattern[i] и строка s. Наша задача — найти все возможные вхождения строк pattern[i] в s.

Суть алгоритма заключена в использование структуры данных — бора и построения по нему конечного детерминированного автомата. Важно помнить, что задача поиска подстроки в строки тривиально реализуется за квадратичное время, поэтому для эффективной работы важно, чтоб все части Ахо-Корасика ассимптотически не превосходили линию относительно длинны строк. Мы вернемся к оценке сложности в конце, а пока поближе посмотрим на составляющие алгоритма.

Построение бора по набору строк-образцов

Структура бора

Что же такое бор? Строго говоря, бор — это дерево, в котором каждая вершина обозначает какую-то строку (корень обозначает нулевую строку — ε). На ребрах между вершинами написана 1 буква (в этом его принципиальное различие с суффиксными деревьями и др.), таким образом, добираясь по ребрам из корня в какую-нибудь вершину и контангенируя буквы из ребер в порядке обхода, мы получим строку, соответствующую этой вершине. Из определения бора как дерева вытекает также единственность пути между корнем и любой вершиной, следовательно — каждой вершине соответствует ровно одна строка (в дальнейшем будем отождествлять вершину и строку, которую она обозначает).
Строить бор будем последовательным добавление исходных строк. Изначально у нас есть 1 вершина, корень (root) — пустая строка. Добавление строки происходит так: начиная в корне, двигаемся по нашему дереве, выбирая каждый раз ребро, соответствующее очередной букве строки. Если такого ребра нет, то мы создаем его вместе с вершиной. Вот пример построенного бора для строк: 1)acab, 2)accc, 3)acac, 4)baca, 5)abb, 6)z, 7)ac.

^{(Рис. 1)}

Обратите внимание на добавление строки 7. Она не создает новых вершин и ребер, а процесс ее добавления останавливается во внутренней вершине. Отсюда видно, что для каждой строки необходимо дополнительно хранить признак того является она строкой из условия или нет (красные круги).

^{(Рис. 2)}

Отметим также что, две строки в боре имеют общие ребра при условии наличия у них общего префикса. Крайний случай — все строки образцы попарно не имеют одинаковой начальной части. Значит верхняя оценка для числа вершин в боре — сумма длин всех строк + 1(корень).

Реализация

Будем хранить бор как массив вершин, где каждая вершина имеет свой уникальный номер, а корень имеет нулевое значение (root = 0). Возможное описание структуры вершины:

//k - размер алфавита
struct bohr_vrtx{
   int next_vrtx[k],pat_num;
   bool flag;
};

next_vrtx[i] — номер вершины, в которую мы придем по символу с номером i в алфавите, flag — бит, указывающий на то, является ли наша вершина исходной строкой, pat_num — номер строки-образца, обозначаемого этой вершиной.
Предподсчет длин всех добавляемых строк — лишние затраты по памяти. Будем использовать структуру данных из STL — vector. В нем память выделяется динамически, следовательно дополнительные затраты будут нулевыми. Явным образом вытекает процедура добавление строки (используем 26-буквенный строчный латинский алфавит => k=26).

vector<bohr_vrtx> bohr;

bohr_vrtx make_bohr_vrtx(){
   bohr_vrtx v;
   //(255)=(2^8-1)=(все единицы в каждом байте памяти)=(-1 в дополнительном коде целого 4-байтного числа int)
   memset(v.next_vrtx, 255, sizeof(v.next_vrtx));
   v.flag=false;
   return v;
}

void bohr_ini(){
   //добавляем единственную вершину - корень
   bohr.push_back(make_bohr_vrtx());
}   

void add_string_to_bohr(const string& s){
   int num=0; //начинаем с корня   
   for (int i=0; i<s.length(); i++){
      char ch=s[i]-'a'; //получаем номер в алфавите
      if (bohr[num].next_vrtx[ch]==-1){ //-1 - признак отсутствия ребра
         bohr.push_back(make_bohr_vrtx());
         bohr[num].next_vrtx[ch]=bohr.size()-1;         
         }
      num=bohr[num].next_vrtx[ch];
   }
   bohr[num].flag=true;
   pattern.push_back(s);
   bohr[num].pat_num=pattern.size()-1;
}

bool is_string_in_bohr(const string& s){
   int num=0;   
   for (int i=0; i<s.length(); i++){
      char ch=s[i]-'a';
      if (bohr[num].next_vrtx[ch]==-1){
         return false;         
         }
      num=bohr[num].next_vrtx[ch];
   }
   return true;
}

Построение автомата по бору

Описание принципа работы

Наша задача — построить конечный детерминированный автомат. Что за штука такая можно посмотреть здесь. Вкратце, состояние автомата — это какая-то вершина бора. Переход из состояний осуществляется по 2 параметрам — текущей вершине v и символу ch. по которорому нам надо сдвинуться из этой вершины. Поконкретнее, необходимо найти вершину u, которая обозначает наидлиннейшую строку, состоящую из суффикса строки v (возможно нулевого) + символа ch. Если такого в боре нет, то идем в корень.

Зачем это нам надо? Предположим, что мы можем вычислить такую вершину быстро, за константное время. Пусть, мы стоим в некой вершине бора, соответствующей подстроке [i..j] строки s, вхождения в которую мы ищем. Теперь найдем все строки бора, суффиксы s[i..j]. Утверждается, что их можно искать быстро (описано далее). После этого, просто перейдем из состояния автомата v в состояние u по символу s[j+1] и продолжим поиск.

^{(Рис. 3)}
Для реализации автомата нам понадобится понятие суффиксной ссылки из вершины.

Суффиксные ссылки

Назовем суффиксной ссылкой вершины v указатель на вершину u, такую что строка u — наибольший cобственный суффикс строки v, или, если такой вершины нет в боре, то указатель на корень. В частности, ссылка из корня ведет в него же. Нам понадобятся суффиксные ссылки для каждой вершины в боре, поэтому немного изменим структуру вершины и процедуру создание вершины, введя дополнительную переменную suff_link.

struct bohr_vrtx{
   //...
   int suff_link; //suff_link - суффиксная ссылка
};

bohr_vrtx make_bohr_vrtx(){
   bohr_vrtx v;
   //...
   v.suff_link=-1; //изначально - суф. ссылки нет
   return v;
}

Вот пример расстановки суф. ссылок для бора на рис. 1:

^{(Рис. 4)}

Реализация автомата

Вернемся к задача быстрого перехода между состояниями автомата. Очевидно, что всего возможных переходов существует bohr.size()*k, так как для для каждой возможной вершиной и каждым возможным символом в алфавите нужно знать переход. Предподсчет существенно снизит среднее время работы алгоритма, поэтому воспользуемся идеями ленивой динамики — будем считать по необходимости и запоминать уже сосчитанные значения.
Введем вычисляемую функцию для перехода (v,ch). Идея тут вот в чем: если из текущей вершины есть ребро c символом ch, то пройдем по нему, в обратном случаем пройдем по суффиксной ссылке и запустимся рекурсивно от новой вершины. Почему это работает, догадаться не трудно.

^{(Рис. 5)}

Вопрос лишь в корректном получение суф. ссылки от вершины. В этой задаче тоже можно использовать ленивую динамику. Эвристика заключена в следующем: для получения суф. ссылки вершины v (строки s[i..j]) спустимся до ее предка par, пройдем по суф. ссылке par и запустим переход от текущей вершины t по символу symb, который написан на ребре от par до v. Очевидно, что сначала мы попадем в наибольший суфикс s[i..j-1] такой что, он имеет ребро с символом symb, потом пройдем по этому ребру. По определению, получившаяся вершина и есть суффикная ссылка из вершин v.

^{(Рис. 6)}

Итак, видно, что функции получения суффиксной ссылки и перехода из состояния автомата взаимосвязаны. Их удобная реализация представляет 2 функции, каждая из которых рекурсивно вызывает другую. База обоих рекурсий — суф. ссылка из корня или из сына корня ведет в корень.

struct bohr_vrtx{
   //...
   int auto_move[k],par; //auto_move - запоминание перехода автомата, par - вершина-отец в дереве
   char symb; //символ на ребре от par к этой вершине 
};

bohr_vrtx make_bohr_vrtx(int p,char c){ //передаем номер отца и символ на ребре в боре
   bohr_vrtx v;
   //...
   memset(v.auto_move, 255, sizeof(v.auto_move));
   v.par=p;
   v.symb=c;
   return v;
}

int get_auto_move(int v, char ch);
 
int get_suff_link(int v){
   if (bohr[v].suff_link==-1) //если еще не считали
      if (v==0||bohr[v].par==0) //если v - корень или предок v - корень
         bohr[v].suff_link=0;
      else
         bohr[v].suff_link=get_auto_move(get_suff_link(bohr[v].par), bohr[v].symb);
   return bohr[v].suff_link;
}
 
int get_auto_move(int v, char ch){
   if (bohr[v].auto_move[ch]==-1)
      if (bohr[v].next_vrtx[ch]!=-1)
         bohr[v].auto_move[ch]=bohr[v].next_vrtx[ch];
      else
         if (v==0)
            bohr[v].auto_move[ch]=0;
         else
            bohr[v].auto_move[ch]=get_auto_move(get_suff_link(v), ch);
   return bohr[v].auto_move[ch];
}

Выявление «хороших» суффиксных ссылок

С автоматом несложно определить сам алгоритм: считываем строку, двигаемся из состояния в состояние по символам строки, в каждом из состояний двигаемся по суф. ссылками, то есть по суффиксам строки в позиции автомата, проверяя при этом наличие их в боре.
Все бы ничего, но оказывается, что этот вариант Ахо-Корасика имеет квадратичную ассимптотику относительно N — длинны считываемой строки s. Действительно, для строки из состояния v можно найти v.length() суффиксов, а переход из состояний может просто увеличивать на 1 длину этой строки. Цимес в том, чтобы двигаясь по суф. ссылкам попадать только в заведомо имеющиеся среди строк-образцов. Введем понятие «хороших» суф. ссылок suff_flink. Так, bohr[v].suff_flink — это ближайший суффикс, имеющийся в боре, для которого flag=true. Число «скачков» при использование таких ссылок уменьшится и станет пропорционально количеству искомых вхождений, оканчивающихся в этой позиции.

^{(Рис. 7)}

struct bohr_vrtx{
   //...
   int suff_flink; //suff_flink - "хорошая" суф. ссылка
};

bohr_vrtx make_bohr_vrtx(int p,char c){
   bohr_vrtx v;
   //...
   v.suff_flink=-1;
   return v;
}

Вычислять их довольно просто, все той же ленивой динамикой. Введем функцию подсчета «хорошой» суф. ссылки. Если для вершине по суф. ссылке flag=true, то это и есть искомая вершина, в ином случае рекурсивно запускаемся от этой же вершины.

int get_suff_flink(int v){ 
   if (bohr[v].suff_flink==-1){
      int u=get_suff_link(v);
      if (u==0) //либо v - корень, либо суф. ссылка v указывает на корень 
         bohr[v].suff_flink=0;
      else
         bohr[v].suff_flink=(bohr[u].flag)?u:get_suff_flink(u);
   }
   return bohr[v].suff_flink;
}

Реализация поиска по автомату

Поиск реализуется тривиально. Нам потребуется процедура хождения по «хорошим» суф. ссылкам check(v,i) из текущей позиции автомата v учитывая, что эта позиция оканчивается на i-ую букву в слове s.

void check(int v,int i){
    for(int u=v;u!=0;u=get_suff_flink(u)){
        if (bohr[u].flag) cout<<i-pattern[bohr[u].pat_num].length()+1<<" "<<pattern[bohr[u].pat_num]<<endl;
    }
}

void find_all_pos(const string& s){
    int u=0;
    for(int i=0;i<s.length();i++){
        u=get_auto_move(u,s[i]-'a');
        check(u,i+1);
    }
}

   bohr_ini();
   add_str_to_bohr("abc");
   add_str_to_bohr("bcdc");
   add_str_to_bohr("cccb");
   add_str_to_bohr("bcdd");
   add_str_to_bohr("bbbc");
   find_all_pos("abcdcbcddbbbcccbbbcccbb");   

Оценка сложности и способы хранения

Существующий вариант алгоритм проходит циклом по длине s (N=s.length()), откуда его уже можно оценить как O(N*O(check)), но так как check прыгает только по заведомо помеченным вершинам, для которых flag=true, то общую ассимптотику можно оценить как O(N+t), где t — количество всех возможных вхождений всех строк-образцов в s. Если быть точным и учитывать вычисления автомата и суф. ссылок, то алгоритм работает O(M*k+N+t), где M=bohr.size(). Память — константные массивы размера k для каждой вершины бора, откуда и выливается оценка O(M*k).

Оказывается, другой способ хранение, а конкретно, обращения к алфавиту, способен изменить эту оценку. Будем использовать отображение map<char,int> вместо массива. Читаем здесь и здесь, видим, что структура данных map из STL реализована красно-черным деревом, а время обращения к его элементам пропорционально логарифму числа элементов. В нашем случае — двоичному логарифму размера алфавита k (что практически константа). Общее время — O((M+N)*log k+t). На практике, это значительнее быстрее массива. Map вовсе не хранит лишних ячеек памяти для элементов, поэтому память пропорциональна количеству ребер в боре (а следовательно и количеству вершин в боре, т.к. в дереве с M вершин — M-1 ребер). Количество вычислений переходов автомата, очевидно, пропорционально длине строки. Получившаяся оценка — O((M+N)*log k).

_{PS: Вот ссылка на полностью реализованный рабочий алгоритм: Aho-Corasick alg}.

Автор: rmq

Источник