Thuật toán Knuth-Morris-Pratt -- KMP string matching

Trong loạt bài này mình giới thiệu các thuật toán thao tác với các xâu kí tự. Chúng ta bắt đầu bằng bài toán cơ bản nhất: tìm sự xuất hiện của một xâu kí tự trong một văn bản.

Problem 1: Cho một văn bản T[1,2,\ldots,n] có chiều dài n với mỗi T[i] là một kí tự trong bảng chữ cái (alphabet) \Sigma và một xâu mẫu P[1,2,\ldots,m] có chiều dài m (m \leq n). Tìm (các) vị trí xuất hiện của xâu P (như là một xâu liên con liên tục) trong văn bản T.

 
Ví dụ: đoạn văn bản T[1,2,\ldots,15] = \mathrm{abcbabababacaab} và xâu mẫu P[1,2,\ldots, 7] = \mathrm{ababaca}, P xuất hiện ở vị trí 7 của văn bản P[1,\ldots, 7] = T[7,\ldots, 13]. Trong ví dụ này bảng chữ cái \Sigma gồm các kí tự latin viết thường.

Chúng ta bắt đầu với thuật toán trâu bò (brute force), đôi khi còn gọi là vét cạn. Các thuật toán trâu bò thường chậm tuy nhiên dễ thiết kế, cài đặt và thường là bước đầu tiên trước khi chúng ta phát triển các thuật toán cải tiến cao cấp hơn.

Thuật toán trâu bò

Ý tưởng của thuật toán trâu bò khá đơn giản như sau: bắt đầu từ kí tự đầu tiên của văn bản T, dịch dần xâu P về bên phải từng kí tự một. Nếu ta dịch đến kí tự thứ i, so sánh T[i,i+1,\ldots, i+{m-1}] với P[1,2,\ldots, m]. Nếu ta tìm thấy xâu mẫu, i.e, T[i,i+1,\ldots, i+{m-1}] = P[1,2,\ldots, m], ta trả về vị trí i. Nếu không ta tiếp tục dịch sang vị trí i+1. Giả mã như sau:

BFMatcher(T[1,2,\ldots, n], P[1,2,\ldots,m]):
    for i \leftarrow 1 to n-m + 1
        \mathrm{matched} \leftarrow \mathrm{TRUE}
        s \leftarrow 0
        while \mathrm{matched} and s \leq m-1
            if T[i + s] \not= P[s+1]
                \mathrm{matched} \leftarrow \mathrm{FALSE}         [[found a mismatched]]
            else
                s \leftarrow s+1
        if \mathrm{matched}
            return i
    return \mathrm{NONE}

 
Code của thuật toán bằng C:

// assume both string T and P are indexed from 1
// you can force this condition by adding a sentinel character at the beginning of both strings.
int BFmatcher(int n, int m){
	int i = 0, s = 0;
	int matched = 1;
	for(i = 1;  i <= n-m+1; i++){
		s = 0;
		matched = 1;
		while((matched) && (s <= m-1)){
			if (T[i+s] != P[s+1]){
					matched = 0;
			}else {
				s++;
			}
		}
		if(matched){
			return i;
		}
	}
	return -1; // -1 is NONE
}

Có thể thấy với mỗi vòng lặp, chúng ta mất tối đa m phép so sánh để kiểm tra xâu P[1,2,\ldots,m] có xuất hiện trong văn bản hay không. Do đó:

Lemma 1: Thời gian tìm xâu mẫu P[1,2,\ldots,m] trong văn bản T[1,2,\ldots,n] của giải thuật trâu bò là O((n-m)m) = O(mn).

 

Thuật toán Knuth-Morris-Pratt

Thuật toán KMP được Knuth, Morris và Pratt phát triển vào năm 1977 là thuật toán tìm kiếm có thời gian tuyến tính O(n+m). Thuật toán này chạy khá tốt trong thực tế, tuy nhiên cũng khá khó hiểu, vì thế khó gỡ lỗi. Trước hết mình trình bày một biến thể của thuật toán (biến thể này được trình bày trong [1]) với thời gian tính cỡ O(n + |\Sigma|m). Thuật toán này có hai điểm mạnh:

  1. Thuật toán này ,theo mình, đơn giản hơn thuật toán gốc và bao gồm ý tưởng chính của thuật toán gốc.
  2. Nếu kích thước của bảng chữ cái không quá lớn (256 với kí tự ASCII), thuật toán này có thể coi là tuyến tính O(n+m) (mình đã test thuật toán này với bảng chữ cái kích thước 26 trên spoj và đã được accept).

Note: Mình khuyến khích bạn đọc thêm bài về máy trạng thái hữu hạn để hiểu hơn về thuật toán KMP.

Thuật toán O(n+|\Sigma|m)

Ý tưởng của thuật toán KMP dựa trên quan sát rằng trong thuật toán trâu bò có nhiều bước so sánh dư thừa. Cụ thể, khi so sánh phần tử T[i] với P[j], nếu T[i] \not = P[j], ta đã có thông tin:

T[i-1] = P[j-1], T[i-2] = P[j-1], \ldots, T[i-j+1] = P[1]  \qquad (1)

Thuật toán trâu bò dịch sang P lên 1 đơn vị bằng cách và thiết lập biến chạy trên xâu mẫu j = 0không tận dụng thông tin này. Trong phần này, thay vì nói "dịch xâu mẫu lên x đơn vị" ta sẽ nói "lùi biến chạy j về vị trí s". Hai cách nói này về cơ bản là tương đương.

Như vậy tại sao thông tin này lại có ích? Xét văn bản T và xâu mẫu P trong ví dụ đầu tiên.
kmp-mis
Khi so sánh T[5,6,\ldots, 11] với P[1,2,\ldots, 7], chúng ta thấy hai xâu này không khớp tại vị trí T[10] \not= P[6] (i = 10, j = 6). Thay vì thiết lập biến chạy  j = 0 như trong thuật toán trâu bò, ta thiết lập j = 4 và tiếp tục so sánh T[i+1] = T[11] với T[j+1] = T[5]. Cứ làm như vậy, thì mỗi vòng lặp biến chạy trên văn bản i luôn tăng và do đó thời gian tính toán sẽ là O(n) nếu như ta có thể tính được vị trí cần lùi lại của biến chạy j trong thời gian O(1).

kmp-patter-shift

Tại sao trong ví dụ trên ta lại lùi biến chạy j về vị trí 4. Quan sát kĩ ta thấy chúng ta lùi j về vị trí 4 là do T[10,9,8,7] = P[4,3,2,1], tức là tiền tố (prefix) kết thúc tại j=4 của P khớp với hậu tố (suffix) bắt đầu từ 7 của T[1,2,\ldots, 10]. Tổng kết lại, quy luật lùi biến chạy j như sau:

Quy luật 1: Lùi j về vị trí s sao cho tiền tố kết thúc tại s của P[1,2,\ldots,j] (chính là P[1,2,\ldots, s]) khớp với hậu tố bắt đầu tại i-s của T[1,2,\ldots,i] (chính là T[i-s,i-s+1,\ldots, i]) và độ dài của đoạn khớp này (chính là s) lớn nhất.

Ta gọi kí tự T[i] là kí tự chốt (pivot). Trong ví dụ trên, nếu ta lùi lại j=5 thì P[1,2,\ldots, 5] không khớp khới với T[6,7,8,9,10]. Nếu ta lùi lại j = 4 thì P[1,2,3, 4] khớp khới với T[7,8,9,10] và độ dài đoạn khớp này là 4. Tương tự, nếu ta lùi j = 2 thì P[1,2] khớp khới với T[9,10] và độ dài đoạn khớp là 2. Theo luật trên, ta lùi j về vị trí sao cho độ dài đoạn khớp là lớn nhất, do đó, ta lùi j = 4.

Nhận xét thấy vị trí lùi lại j = s phụ thuộc vào hai yếu tố:

  1. Kí tự chốt T[i].
  2. Vị trí xảy ra không khớp trong P (chính là j).

Trường hợp kí tự chốt T[i] không xuất hiện trong P, ta lùi lại j = 0 (hay nói cách khác ta dịch P lên j đơn vị). Ví dụ nếu T[10] = \mathrm{f} thì ta có thể dịch P lên 6 đơn vị. Dựa trên nhận xét trên, ta sẽ tính trước mảng S[1,2,\ldots, |\Sigma|][1,2,\ldots, m] trong đó S[p,q] là giá trị ta sẽ lùi lại nếu T[i] là kí tự thứ p trong bảng chữ cái và vị trí xảy ra không khớp là q trong xâu mẫu P. Ví dụ trong ví dụ trên T[i]= \mathrm{b} (tương đương với p = 2 nếu \Sigma = \{a,b,c\}) và q = j = 6, thì S[2,6] = 4. Giả sử ta đã tính được bảng S bằng thủ tục ComputeBackup (ở dưới), thuật toán tìm xâu sẽ như sau:

KMPMatcher(T[1,2,\ldots, n], P[1,2,\ldots,m]):
    ComputeBackup(\Sigma, P[1,2,\ldots, m])
    j \leftarrow 1
    for i  \leftarrow 1 to n
        if T[i] \not= P[j]         [[found a mismatched]]
            p \leftarrow index of T[i] in \Sigma
            j \leftarrow S[p,j]            [[the restarting value]]
        j \leftarrow j+1
        if j = m+1
            return i-m+1.
    return \mathrm{NONE}

 
Code của giả mã bằng C:

// assume both string T and P are indexed from 1
// you can force this condition by adding a sentinel character at the beginning of both strings.
char Alph[28] = " abcdefghijklmnopqrstuvwxyz"; // the alphabet

int KMPmatcher(int n, int m){
	compute_backup(m);
	int i = 1, j = 1, p = 0;
	for(i = 1; i <=n; i++){
		if(T[i] != P[j]){
			p = T[i]-96; the index of T[i] in the alphabet
			j = S[p][j];
		}
		else j++;
		if(j == m+1){
			return i-m+1;
		}
	}
	return -1;
}
 

Với mỗi vòng lặp ta mất thời gian O(1), do đó, tổng thời gian là O(n + T(|\Sigma|,m)) trong đó T(|\Sigma|,m) là thời gian tính bảng S[1,2,\ldots, |\Sigma|][1,2,\ldots,m]. Để tính bảng S[1,2,\ldots, |\Sigma|][1,2,\ldots,m], ta sử dụng chính ý tưởng lùi biến chạy ở trên kết hợp với quy hoạch động.

Gọi X[1,2,\ldots,m] là mảng trong đó X[j] là chỉ số s < j-1 lớn nhất sao cho P[1,2,\ldots, s] khớp với P[j-s+1,j-s+2,\ldots, j]. Gọi p là chỉ số của P[j] trong bảng chữ cái \Sigma. S[k,j]X[j] có thể được tính thông qua X[j-1] dựa trên công thức quy hoạch động sau:

\begin{array}{lcl} S[k,j] &= & S[k,X[j-1]+1] \\ S[p,j] &= &j+1 \\ X[j] &=& S[p,X[j-1]+1] - 1 \qquad \mbox{why -1?}  \end{array}

Khởi tạo S[k,j] = 0 với mọi 1 \leq k \leq |\Sigma|X[0] = -1 (why -1?). Giả mã như sau:

ComputeBackup(\Sigma, P[1,2,\ldots, m]):
    X[0]\leftarrow -1
    S[1,2,\ldots,|\Sigma|][0] \leftarrow 1
    for j\leftarrow 1 to m
        for k\leftarrow 1 to |\Sigma|
            S[k,j] \leftarrow  S[k,X[j-1]+1]
        p \leftarrow index of P[j] in \Sigma
        S[p,j] = j+1
        X[j] = S[p,X[j-1]+1]-1

 
Code của giả mã bằng C:

void compute_backup(int m){
	int X = -1;
	int j = 0 ,k = 0, p = 0;
	for(k = 0; k < 27; k++){ // 26 is the size of the alphabet
		memset(S[k], 0, sizeof(S[k])); // set everything to 0
	}
	for(k = 0; k < 27 ; k ++){
		S[k][0] = 1;
	}
	for(j = 1; j <=m ; j++){
		for(k = 1; k < 27; k++){
			S[k][j] = S[k][X+1];
		}
		p = P[j]-96; // the index of P[j] in the alphabet
		S[p][j] = j+1;
		X = S[p][X+1]-1;
	}
}

Do X[j]S[k][j] chỉ phụ thuộc vào X[j-1], ta có thể tiết kiệm một chút bộ nhớ bằng cách chỉ dùng một biến X thay vì dùng cả mảng X[1,2,\ldots,m]. Dễ thấy thuật toán trên mất thời gian O(|\Sigma|m) để cập nhật bảng S, do đó:

Lemma 2: Thời gian tìm xâu mẫu P[1,2,\ldots,m] trong văn bản T[1,2,\ldots,n] của giải thuật KMPMatcherO(n + |\Sigma|m).

 

Thuật toán O(n+m)

Trong thuật toán trước, sự phụ thuộc thời gian vào kích thước bảng chữ cái dường như dễ hiểu vì kí tự chốt T[i] có thể là bất kì kí tự nào trong bảng chữ cái. Để xóa bỏ sự phụ thuộc này, ta chọn kí tự chốt là T[i-1]. Do T[i-1] = P[j-1], kí tự này luôn xuất hiện trong P. Quy luật lùi biến chạy có thể phát biểu lại như sau:

Quy luật 3/2: Lùi j về vị trí s sao cho tiền tố kết thúc tại s-1 của P[1,2,\ldots, j-1] (chính là P[1,2,\ldots, s-1]) khớp với hậu tố bắt đầu tại i-s của T[1,2,\ldots,i-1] (chính là T[i-s,i-s+1,\ldots, i-1]) và độ dài của đoạn khớp này (chính là s) lớn nhất.

Sự khác biệt cơ bản giữa quy luật 1 và quy luật 3/2 đó là: sau mỗi lần lùi biến chạy j, quy luật 1 luôn đảm bảo T[i] = P[j] và do đó ta tiếp tục xét kí tự tiếp theo i+1. Còn quy luật 3/2 không đảm bảo điều này, do đó, sau mỗi lần lùi biến chạy j, nếu T[i] \not= P[j], ta phải tiếp tục lùi biến chạy j một lần nữa (nghe hơi giống đệ quy). Dường như thuật toán này mất nhiều phép so sa nhơn thuật toán trước. Tuy nhiên, phân tích dưới đây chỉ ra rằng số lần lùi lại như vậy không nhiều.

Do T[i-j+1,i-j+2, \ldots, i-1] = P[1,2,\ldots,j-1], quy luật trên có thể được phát biểu lại như sau:

Quy luật 2: Lùi j về vị trí s sao cho tiền tố kết thúc tại s-1 của P[1,2,\ldots, j-1] (chính là P[1,2,\ldots, s-1]) khớp với hậu tố bắt đầu tại j-s của P[1,2,\ldots,j-1] (chính là P[j-s,j-s+1,\ldots, j-1]) và độ dài của đoạn khớp này (chính là s) lớn nhất.

Trong ví dụ trên, nếu ta lùi lại j=5 thì tiền tố P[1,2,\ldots, 4] = \mathrm{abab} của P[1,2,\ldots, 6] không khớp khới với hậu tố P[2,3,4,5] = \mathrm{baba}. Nếu ta lùi lại j = 4 thì tiền tố P[1,2,3] = \mathrm{aba} của P[1,2,\ldots, 6] khớp với hậu tố P[3,4,5] = \mathrm{aba} của P[1,2,\ldots, 6 và độ dài đoạn khớp này là 3. Tương tự, nếu ta lùi j = 2 thì P[1] khớp với P[5] và độ dài đoạn khớp là 1. Theo luật trên, ta lùi j về vị trí sao cho độ dài đoạn khớp là lớn nhất, do đó, ta lùi j = 4.

kmp-shift-2

Theo quy luật 2, ta có thể thấy vị trí lùi lại j = s chỉ phụ thuộc vào vị trí xảy ra không khớp trong P (chính là j), không phụ thuộc vào T. Gọi F[j] là vị trí j phải lùi về khi T[j] \not=  P[i]. Giả sử ta đã tính được F[1,2,\ldots,m] bằng thủ tục ComputeFailure (ở dưới), thuật toán tìm xâu sẽ như sau:

FastKMPMatcher(T[1,2,\ldots, n], P[1,2,\ldots,m]):
    ComputeFailure(P[1,2,\ldots, m])
    j \leftarrow 1
    for i  \leftarrow 1 to n
        while T[i] \not= P[j] and j  data-recalc-dims= 0" />
            j \leftarrow F[j]         [[backing up the pattern]]
        if j = m
            return i-m+1
        j \leftarrow j+1
    return \mathrm{NONE}

 
Code bằng C của giả mã:

int fastKMPmatcher(int n, int m){
	compute_failure(m);
	int i = 1, j = 1;
	for(i = 1; i <=n ; i++){
		while(T[i] != P[j] && j > 0){ // found a mismatch
			j = F[j];   // backup the pattern
		}
		if(j == m){
			return i-m+1;
		}
		j++;
	}
	return -1;

}

Giả sử thủ tục ComputeFailure(P[1,2,\ldots, m]) sử dụng T(m) phép so sánh, ta sẽ chứng minh rằng số phép so sánh của thuật toán FastKMPMatcher(T[1,2,\ldots, n], P[1,2,\ldots,m]) là O(n + T(m)). Không dễ để có thể thấy được điều này vì sự xuất hiện vòng lặp while bên trong vòng lặp for. Trước hết ta thấy sau mỗi phép so sánh bằng (T[i] = P[j]), cả ij đều tăng lên 1. Do đó, số phép so sánh bằng tối đa là n-1. Vòng lặp while giảm giá trị của j mỗi khi ta bắt gặp một phép so sánh không bằng (T[i] \not= T[j]). Tuy nhiên, số lượng giảm của j không thể vượt quá giá trị mà j đã tăng lên (sau mỗi phép so sánh bằng) qua các vòng lặp trước đó. Điều đó có nghĩa là tổng số phép so sánh không bằng không thể vượt quá tổng số phép so sánh bằng, do đó, tối đa là n-1. Như vậy, tổng số phép so sánh là 2(n-1) + T(m) = O(n + T(m)).

F[1,2,\ldots,m] có thể được tính dựa vào định nghĩa của F[j] và quy luật 2 trong thời gian O(m^3). Cụ thể, để tính F[j], ta sẽ kiểm tra tất cả các chỉ số s từ 1 đến j-2 sao cho P[1,2,\ldots, s-1] khớp với P[j-s,j-s+1,\ldots, j-1] và chọn ra chỉ số s lớn nhất thỏa mãn điều kiện này. Tuy nhiên, Knuth-Morris-Pratt chỉ ra rằng bảng F[1,2,\ldots,m] có thể được tính trong thời gian O(m) bằng cách sửa đổi thuật toán FastKMPMatcher đối sánh P với chính nó.

ComputeFailure(P[1,2,\ldots, m]):
    k \leftarrow 0
    F[1] \leftarrow 0
    for j  \leftarrow 1 to m-1
        while P[j] \not= P[k] and k  data-recalc-dims= 0" />
            k \leftarrow F[k] \qquad(*)
        k \leftarrow k+1
        F[j+1] \leftarrow k

 

Code bằng C của giả mã:

void compute_failure(int m){
	int k = 0, j = 1;
	F[1] = 0;
	for(j=1 ; j < m; j++){
		while(P[j] != P[k] && k > 0){
			k =  F[k];
		}
		k++;
		F[j+1] = k;
	}
}

Tại sao đối sánh P với chính nó lại cho ra kết quả ta mong muốn. Về mặt trực quan, để ý quy luật 3/2, nếu thay T trong quy luật này bởi P, và j bởi k, giá trị s trong quy luật 3/2 chính là F[j]. Thủ tục FastKMPMatcher về cơ bản thực thi quy luật này, do đó, ta có thể sử dụng chính thủ tục này để tính F[j].

Một cách hiểu khác của thủ tục này là như sau: giả sử F[j] = s, ta biết rằng P[1,2,\ldots, s-1] =  P[j-s, \ldots, j-1]. Do đó, ta tìm các ứng cử viên k, bắt đầu từ ứng viên lớn nhất (chính là F[j-1]), sao cho:

P[1,2,\ldots, k-1] =  P[j-k, \ldots, j-2]

Nếu P[k] = P[j-1] thì k+1 chính là F[j]. Nếu P[k] \not= P[j-1], ứng viên tiếp theo sẽ là F[k] (theo định nghĩa của bảng F), do đó ở dòng (*), ta gán k \leftarrow F[k]. Do ta xét các ứng cử viên từ lớn đến nhỏ, giá trị tìm được sẽ tương ứng với đoạn khớp lớn nhất. Phân tích tương tự như trên, số phép so sánh trong thủ tục ComputeFailure2(m-1) = O(m). Do đó,

Lemma 3: Thời gian tìm xâu mẫu P[1,2,\ldots,m] trong văn bản T[1,2,\ldots,n] của giải thuật FastKMPMatcherO(n + m).

 

Ví dụ: bảng F[1,2,\ldots, 7] ứng với mẫu P = \mathrm{ababaca} như sau:

P[i] a b a b a c a
F[i] 0 1 1 2 3 4 1

 

Code: kmp-all-in-one, kmp-accpted-code,

Thao khảo
[1] Sedgewick, R., Wayne, K. (2011). Algorithms, 4th Edition, Chapter 5. Addison-Wesley. ISBN: 978-0-321-57351-3
[2] Knuth, Donald E., James H. Morris, Jr, and Vaughan R. Pratt. Fast pattern matching in strings. SIAM Journal on Computing 6.2 (1977): 323-350.
[3] Jeff Erickson, String Matching lecture note, UIUC.
[4] Cormen, Thomas H.; Leiserson, Charles E., Rivest, Ronald L., Stein, Clifford. Introduction to Algorithms (2nd ed.), Chapter 32. MIT Press and McGraw-Hill (2001). ISBN 0-262-03293-7

Facebook Comments

Tags: , ,

Reply

Your email address will not be published. Required fields are marked *