Abstract
Suffixient arrays are recent structures that have attracted attention because they offer relevant pattern matching functionality in less asymptotic space than the Run-Length BWT, the de-facto standard to index highly repetitive string collections. Various algorithms exist for building them from the suffix array data structures. We present the first construction algorithm that is (i) linear-time, (ii) one-pass over the structures, and (iii) implemented and practical. This makes the construction particularly useful on large text collections, which we demonstrate empirically by showing that it dominates the space/time tradeoff map of the implemented constructions.
1 Introduction and Related Work↩︎
Many of the fastest-growing text collections are highly repetitive, meaning that most documents are near-copies of others: genome collections like the One Million Genome Initiative, software repositories like Software Heritage, and versioned text collections like Wikipedia are examples of this phenomenon. Effectively searching such huge collections requires using storage formats that exploit repetitiveness to reduce space while offering fast indexed searches. The Run-Length Burrows-Wheeler Transform (RLBWT) is at the root of many successful compressed indexes for this scenario [1]–[3], BW94?. The space of this family of indexes dominated by \(r\), the number of equal-letter runs in the BWT of the collection.
A very recent development [4] further pushes towards space reduction by defining a new measure \(\chi \le r\), the size of the smallest so-called suffixient set of the text. This is a subset of the text positions that, if sorted by the corresponding prefixes in colexicographic order (i.e., comparing right-to-left), yields the so-called suffixient array [5]. Provided with a method to access the compressed text, the suffixient array can efficiently search for one occurrence of a pattern in the text, or find all its maximal substrings that occur in the text. This is less than the classic pattern matching functionality of finding all the pattern occurrences, but sufficient for various applications while using less space.
Since smallest suffixient sets are meant to be used on large repetitive text collections, their efficient construction is key to their adoption. Cenzato et al. [5], [6] already proposed two efficient constructions, which work on top of the suffix array (SA), the Burrows-Wheeler Transform (BWT), and the longest common prefix array (LCP) of the reversed text. Their first algorithm [6], performs one single pass over those structures (i.e., it is “one-pass”) and takes time \(O(n+\overline{r}\log\sigma)\), where \(n\) is the length of the text, \(\overline{r}\) is the \(r\) value for the reverse text, and \(\sigma\) is the alphabet size. Their second algorithm [6], which we call LF, is linear-time but not one-pass, just as their third algorithm [5], which we call FM. FM is faster than LF but uses more space. Bonizzoni et al. [7] recently proposed the first algorithm that is both linear-time and one-pass, but no implementation is reported. Fujimaru et al. [8] present linear-time online algorithms that work directly on the text, processing them left-to-right or right-to-left. In another SPIRE submission, Fujimaru et al. [9] present another linear-time algorithm, which is however not one-pass nor implemented; it uses some structures built on the text and some on its reverse. They also present the first sublinear-time construction (also not implemented). Finally, Köppl and Kucherov [10] gave a near-real time construction algorithm, without reporting an implementation.
In this context, we present the first algorithm to build smallest suffixient sets that is (i) linear-time, (ii) one-pass, and (iii) implemented. See Table ¿tbl:tab:compar?. We empirically demonstrate its practicality by comparing it with the implemented alternatives: our algorithm outperforms previous algorithms, or matches the best ones, in both time and space.
| Algorithm | linear-time | one-pass | implemented |
|---|---|---|---|
| One-pass [6] | no | yes | yes |
| LF [6] | yes | no | yes |
| FM [5] | yes | no | yes |
| Bonizzoni et al. [7] | yes | yes | no |
| Online [8] | yes | yes*\(\!\!\!\) | yes |
| Linear [9] | yes | no | no |
| Near-linear [10] | no | yes*\(\!\!\!\) | no |
| Ours | yes | yes | yes |
2 Preliminaries↩︎
By \([i, j]\) we refer to all values \(i, i + 1, \dots, j\), if \(i \leq j\), and \([i, j] = \emptyset\), otherwise; and we define \([n] = [1, n]\). A string \(T[1, n] \in \Sigma^n\) is a concatenation of \(n\) elements of an alphabet \(\Sigma = \{1, \dots, \sigma\}\). The length of \(T\) is denoted by \(\vert T \vert = n\). The empty string \(\epsilon\) is the only string satisfying \(\vert \epsilon \vert = 0\). \(T[i]\) indicates the \(i\)-th element of \(T\), and \(T[i, j] = T[i] \dots T[j]\) is the substring of \(T\) starting at position \(i\) and ending at position \(j\), if \(1 \leq i \leq j \leq n\), and \(T[i, j] = \epsilon\) otherwise. A prefix of \(T\) is a substring of the form \(T[1, j]\) and a suffix is a substring of the form \(T[i, n]\). \(T^{\text{rev}}\) denotes the reversal of \(T\), that is, \(T^{\text{rev}} = T[n] \dots T[1]\).
Given two strings \(T, S \in \Sigma^+\), none a prefix of the other, and a total order over \(\Sigma\), the lexicographic order \(<_{lex}\) over \(\Sigma^+\) is defined as \(T <_{lex} S\) if and only if \(T[1,i] = S[1,i]\) and \(T[i + 1] < S[i + 1]\), for some \(i\). The suffix array of \(T[1, n]\) (\(\mathrm{SA}(T)\)) is a permutation of the elements of \(T\) such that \(T[\mathrm{SA}(T)[i], n] <_{lex} T[\mathrm{SA}(T)[j], n]\) for any \(1 \leq i < j \leq n\), namely \(\mathrm{SA}(T)[i]\) is the \(i\)-th suffix of \(T\) in lexicographic order. We assume \(T\) ends with a special character \(\$\), which is smaller than any other symbol in \(\Sigma\) and only occurs in \(T[n]\), so \(\mathrm{SA}\) is well defined. The longest common prefix array of \(T\) (\(\mathrm{LCP}(T)\)) is defined as \(\mathrm{LCP}(T)[i] = \ell\) if and only if \(2 \leq i\) and \(T[\mathrm{SA}(T)[i - 1], \mathrm{SA}(T)[i - 1] + \ell - 1] = T[\mathrm{SA}(T)[i], \mathrm{SA}(T)[i] + \ell - 1]\) and \(T[\mathrm{SA}(T)[i - 1] + \ell] \neq T[\mathrm{SA}(T)[i] + \ell]\), that is, \(\mathrm{LCP}(T)[i]\) is the length of the maximal prefix shared by the \((i - 1)\)-th and the \(i\)-th sufixes of \(T\) in the lexicographic order; we further define \(\mathrm{LCP}(T)[1] = 0\). The Burrows-Wheeler Transform array of \(T\) (\(\mathrm{BWT}(T)\)) is defined as \(\mathrm{BWT}(T)[i] = T[\mathrm{SA}(T)[i - 1]]\), if \(2 \leq i\), and \(\mathrm{BWT}(T)[1] = T[n]\), that is \(\mathrm{BWT}(T)[i]\) is the symbol preceding the \(i\)-th sufix of \(T\). By \(\mathrm{SA}_r\), \(\mathrm{LCP}_r\), and \(\mathrm{BWT}_r\) we refer to \(\mathrm{SA}(T^{\text{rev}})\), \(\mathrm{LCP}(T^{\text{rev}})\), and \(\mathrm{BWT}(T^{\text{rev}})\), respectively, where in this case we assume \(T^{\text{rev}} = (T[1,n-1])^{\text{rev}}\$\).
Given an integer array \(\mathcal{I}\), the previous smaller value array of \(\mathcal{I}\) (\(\mathrm{PSV}(\mathcal{I})\)) is a lenght-\(\vert \mathcal{I} \vert\) array defined as \(\mathrm{PSV}(\mathcal{I})[i] = j\) if and only if \(j = \max (\{k<i \mid \mathcal{I}[k] < \mathcal{I}[i] \} \cup \{0\})\), that is, \(j\) is the largest index in \([1, i - 1]\) for which \(\mathcal{I}\) is strictly smaller than \(\mathcal{I}[i]\), if such value exists, and \(j = 0\) otherwise. The next smaller value array, \(\mathrm{NSV}(\mathcal{I})\), is defined analogously for the smaller value with and index bigger than \(i\) (or \(|\mathcal{I}|+1\) if none exists). For brevity, by \(\mathrm{PSV}\) and \(\mathrm{NSV}\) we will refer to \(\mathrm{PSV}(\mathrm{LCP}_r)\) and \(\mathrm{NSV}(\mathrm{LCP}_r)\), respectively.
A run in \(\mathrm{BWT}_r\) is a maximal contiguous substring that repeats a single symbol. Given \(i \in [2, n]\) and \(a, c \in \Sigma\), with \(a \neq c\), we say that \(i\) is a \(c\)-run break if \(\mathrm{BWT}_r[i - 1, i] = ac\) or \(\mathrm{BWT}_r[i - 1, i] = ca\). By \(box(i) = [l_i, r_i]\) we refer to the maximal interval such that \(i \in [l_i, r_i]\) and \(\mathrm{LCP}_r[i] \leq \mathrm{LCP}_r[j]\), for any \(j \in [l_i, r_i]\), in other words, \(box(i) = [\mathrm{PSV}[i] + 1, \mathrm{NSV}[i] - 1]\). By \(text(i) = n - \mathrm{SA}_r[i] + 1\) we indicate the position in \(T\) corresponding to \(\mathrm{BWT}_r[i]\). Similarly, by \(bwt(j) = \mathrm{SA}_r^{-1}[n + 1 - j]\) we denote the position in \(\mathrm{SA}_r\) of the symbol \(T[j]\).
We say the substring \(T[i, j]\) (\(i - 1 \leq j\)) is right-maximal if \(j = n\) (so \(T[i, j]\) is a suffix of \(T\)) or there exists \(a, b \in \Sigma\) such that \(a \neq b\) and both \(T[i, j] \cdot a\) and \(T[i, j] \cdot b\) occur in \(T\). Note that, by definition, \(\epsilon\) is always right-maximal.
2.1 Suffixient sets↩︎
Suffixient sets were initially proposed by Depuydt et al. [4].
Definition 1 (Suffixient set [4]–[6]). A set \(S \subseteq [n]\) is suffixient* for a string \(T\) if, for every one-character right-extension \(T[i,j]\) (\(j\geq i\)) of every right-maximal string \(T[i,j-1]\), there exists \(x\in S\) such that \(T[i,j]\) is a suffix of \(T[1,x]\).*
They also proposed a compressed index based on a suffixient set of a text \(T\) that can efficiently find one occurrence of any given pattern \(P\) within \(T\), or say that \(P\) does not occur in \(T\). The index can also find all the maximal exact matches (MEMs) of \(P\). They show how to compute a suffixient set of size \(O(\overline{r})\), where \(\overline{r}\) is the number of runs in \(\mathrm{BWT}_r\).
Their index uses \(O(\vert S \vert + g)\) words of space, where \(g\) is the size of a grammar that gives random access to the input text and \(S\) is the suffixient set they build. It follows that the smaller the suffixient set built, the smaller their index is.
Cenzato et al. [6] followed up and studied methods to build suffixient sets of smallest cardinality. They introduced three algorithms to ahieve this goal. The first of them (Algorithm 1) runs in quadratic-time and is conceptually simple. The second one (Algorithm 2) runs in time \(O(n + \overline{r}\log\sigma)\). The algorithm is called one-pass in the paper because it needs just one sequential pass over \(\mathrm{SA}_r\), \(\mathrm{LCP}_r\), and \(\mathrm{BWT}_r\), which enables prefix free parsing, a technique which allows streaming those arrays in compressed space [11]. The third algorithm (Algorithm 4) runs in \(O(n)\) time and uses the same three arrays. Since it relies on the \(\mathrm{LF}\) mapping [12], it needs random access on \(\mathrm{LCP}_r\), so it does not run on top of PFP. Later [5], the same authors propose another linear-time algorithm (Algorithm 7), which uses \(\mathrm{PSV}\)/\(\mathrm{NSV}\) arrays and, in practice, is the fastest between the four algorithms, at the cost of using more working space. In the same paper, the authors proposed an index called suffixient array, a sampling of the positions of the prefix array (the mirroring version of \(\mathrm{SA}\), which uses lexicographic order over the reversed suffixes) which can be binary-searched and offers fast pattern matching, if random access to the text is still available. This index outperforms many well-known ones, such as the \(r\)-index [2], the fastest index in practice offering full pattern matching (i.e., finding all the pattern occurrences, not just one).
Several other construction algorithms were proposed after the work of Cenzato et al. [5], [6]. Bonizzoni et al. [7] proposed the first algorithm that is linear-time and one-pass. Fujimaru et al. [8] proposed an online algorithm that works directly on the text (so it avoids computing \(\mathrm{SA}_r\), \(\mathrm{LCP}_r\), and \(\mathrm{BWT}_r\)), and later [9] presented another linear-time algorithm that builds some structures on the text and on its reverse, but it is not linear-time nor one-pass. They also present the first sublinear-time construction, which relies on some complex data structures proposed by Kempa and Kociumaka [13]. Finally, Köppl and Kucherov [10] gave the first near-real time construction algorithm.
Suffixient sets have also become a field of study from other perspectives. Cenzato et al.[14] showed how to determine in linear time if an arbitrary set is suffixient of smallest cardinality. The measure \(\chi\) (the size of a smallest suffixient set) is also studied as a measure of repetitiveness, due to the combinatorial properties of suffixient sets [15]. It is known that \(\chi = O(\overline{r})\), and also that \(\chi = O(r)\) [15], so \(\chi = O(\min(\overline{r}, r))\).
3 A Simple Linear-Time Construction Algorithm↩︎
As pointed in Sec. 2.1, Cenzato et al. [5] give a simple algorithm to compute a suffixient set of smallest cardinality. They define a set \(B_{i, c}\) that contains the positions in \(box(i)\) that are \(c\)-run breaks in \(\mathrm{BWT}_r\). Then, for a \(c\)-run break \(i\), with \(\mathrm{BWT}_r[i'] = c\) for some \(i' \in \{i - 1, i\}\), the algorithm relies on adding to the resulting set the position \(text(i')\) if and only if it holds \[i = \max\{j\;:\;j\in B_{i, c} \wedge \mathrm{LCP}_r[j] = \max \mathrm{LCP}_r[B_{i, c}]\}.\label{eq:rightmost-maximum}\tag{1}\] That algorithm spends \(O(n)\) time to determine if Eq. (1 ) holds, and, since boxes may overlap, it spends quadratic time to run the whole algorithm. In this section we present a method to evaluate Eq. (1 ) in \(O(1)\) time, which leads to a linear-time version of the algorithm. As a first step, we introduce the notion of last c-candidate positions (Definition 2). Later, we prove that last c-candidate positions are exactly the positions that satisfy Eq. 1 (Proposition [prop:last-candidate]). Finally, we show a modified version of the quadratic algorithm using this equivalence (Algorithm 2).
In the rest of this section, we assume there are exactly \(k\) \(c\)-run breaks in \(\mathrm{BWT}_r\). Let \(1 \leq i_1 < i_2 < \dots < i_k \leq n\) be the positions of those run-breaks, and let \(box(i_h) = [l_h, r_h]\), for each \(h \in [k]\).
Definition 2 (last c-candidate). Let \(i_h\) be the position of a \(c\)-run break. We say \(i_h\) is a c-candidate* if one of the following conditions hold:*
\(h = 1\),
\(i_{h - 1} < l_h\), or
\(i_h \le r_{h - 1}\) and \(i_{h - 1}\) is a \(c\)-candidate.
Moreover, we say that \(i_h\) is a last c-candidate* if it is a c-candidate and one of the following conditions hold:*
\(h = k\), or
\(r_h < i_{h + 1}\).
Intuitively, a position is s c-candidate if it is a \(c\)-run break and its \(\mathrm{LCP}_r\) value is not smaller than the \(\mathrm{LCP}_r\) value of any previous \(c\)-run break (if any) within its box. On the other hand, a c-candidate position is a last c-candidate if there are no larger c-candidates positions within its box. The following proposition establishes the equivalence between finding \(c\)-run break positions that are the rightmost maximum within a box and last c-candidate positions.
Let \(i_h\) be a \(c\)-run break. We have \[i_h = \max\{i\;:\;i\in B_{i_h, c} \wedge \mathrm{LCP}_r[i] = \max \mathrm{LCP}_r[B_{i_h, c}]\}\] if and only if \(i_h\) is a last c-candidate.
Proof. Let \(B_{i_h, c} = \{i_{h - a}, \dots, i_h, \dots, i_{h + b}\}\), \(box(i_h) = [l_h, r_h]\), and define \(i_0 := 0\) and \(i_{k + 1} := n + 1\). By definition of \(B_{i_h, c}\), for every \(j \in [h - a, h + b]\) it holds \(\mathrm{LCP}_r[i_h] \leq \mathrm{LCP}_r[i_j]\). In addition, it holds \[\label{eq:lc-boxes}i_{h - a - 1} < l_h \leq i_{h - a} < \dots < i_h < \dots < i_{h + b} \leq r_h < i_{h + b + 1}\tag{2}\]
(\(\Rightarrow\)) Suppose \(i_h = \max\{i\;:\;i \in B_{i_h, c} \wedge \mathrm{LCP}_r[i] = \max \mathrm{LCP}_r[B_{i_h, c}]\}\). By the maximality of \(\mathrm{LCP}_r[i_h]\), for every \(i_j \in B_{i_h, c}\) it holds \(\mathrm{LCP}_r[i_h] = \mathrm{LCP}_r[i_j]\). This implies \(b = 0\) and \(box(i_{h - a}) = \dots = box(i_h)\). Since by Eq. (2 ) it holds \(i_{h - a - 1} < l_{h - a} \leq i_{h - a}\), \(i_{h - a}\) is a c-candidate. In addition, by Eq. (2 ), for \(d = 1, \dots, a\) we have \(i_{h - a + d} < r_{h - a + d - 1}\), then \(i_{h - a + d}\) is a c-candidate, so \(i_h\) is a c-candidate. Finally, by Eq. (2 ) it holds \(r_h < i_{h + 1}\), therefore, \(i_h\) is a last c-candidate.
(\(\Leftarrow\)) We proceed by contrapositive. Suppose \(i_h \neq i_m = \max\{i\;:\;i\in B_{i_h, c} \wedge \mathrm{LCP}_r[i] = \max \mathrm{LCP}_r[B_{i_h, c}]\}\). If \(m < h\), since \(i_m \in B_{i_h, c}\), there exists \(e \in [h - m]\) such that \(\mathrm{LCP}_r[i_{h - e}] > \mathrm{LCP}_r[i_{h - e + 1}] = \dots = \mathrm{LCP}_r[i_{h}]\). Because \(\mathrm{LCP}_r[i_{h - e }] > \mathrm{LCP}_r[i_{h - e + 1}]\), we have \(l_{h - e + 1} \le i_{h - e}\) and \(r_{h - e} < i_{h - e + 1}\), so \(i_{h - e + 1}\) is not a c-candidate. Now, for \(d = 1, \dots, e\) we have \(l_{h - e + d} \le i_{h - e + d - 1}\) and \(i_{h - e + d - 1}\) is not a c-candidate, then \(i_{h - e + d}\) is not a c-candidate, which implies \(i_h\) is not a c-candidate. On the other hand, if \(h < m \leq k\), then \(1 \leq b\) and \(i_{h + 1} \in B_{i_h, c}\), so we have \(\mathrm{LCP}_r[i_h] \leq \mathrm{LCP}_r[i_{h + 1}]\), which means \(i_{h + 1} \le r_h\). Therefore, \(i_h\) cannot be a last c-candidate.0◻ ◻
From Proposition [prop:last-candidate] we can develop a particular strategy to find local maxima inside boxes by finding last candidates within them. Actually, since the candidate condition propagates forward, a single sequential scan of the \(\mathrm{BWT}_r\) suffices to identify last candidate positions, given that we know boxes’ boundaries. As mentioned in Section 2, we can compute \(box(i) = [l_i, r_i] = [\mathrm{PSV}[i] + 1, \mathrm{NSV}[i] - 1]\) in \(O(1)\) time.
Following the above ideas, we show a linearizarion of the quadratic algorithm in Algorithm 1. We call this algorithm plain last-candidate (PLC) algorithm (the reason for “plain” will be made clear later). It uses a data structure \(R\) to save information about the candidate condition on the last \(c\)-run break seen up to now. The algorithm sequentially scans the \(\mathrm{BWT}_r\) looking for run-breaks. For a given \(c\)-run break \(i\), where \(\mathrm{BWT}_r[i'] = c\) for some \(i' \in \{i -1, i\}\), \(R[c]\) stores \((sa\_pos \gets i, text\_pos \gets n - \mathrm{SA}_r[i'] + 1, candidate \in \{true, false\})\), where \(candidate\) is set to \(true\) if and only if position \(i\) is a c-candidate.
Algorithm 1 shows the PLC algorithm. We prove next its correctness.
Lemma 1. Given \(T \in \Sigma^n\), Algorithm 1 computes a suffixient set of smallest cardinality \(\mathcal{S}\) in \(O(n)\) time and \(O(n)\) space.
Proof. For each \(c \in \Sigma\), \(R[c]\) is set to \((1, 0, true)\) at the beginning of the algorithm (line [line:lc46init]). Then, the first time we find a \(c\)-run break, the condition of line [line:lc46check32lc] is never satisfied, but the condition of line [line:lc46check32candidate] is true, so \(R[c]\) is updated and set as a c-candidate (because the first \(c\)-run break is always a c-candidate). Now, assume that we have scanned the \(\mathrm{BWT}_r\) until position \(i - 1\) and \(i\) is a \(c\)-run break, with \(\mathrm{BWT}_r[i'] = c\) for some \(i' \in \{i - 1, i\}\). Also, assume that \(j < i\) is the previous \(c\)-run break, \(j' \in \{j - 1, j\}\) and \(\mathrm{BWT}_r[j'] = c\). Then \(R[c].sa\_pos = j\), \(R[c].text\_pos = n - \mathrm{SA}_r[j'] + 1\). If \(R[c].candidate = true\) and \(\mathrm{NSV}[R[c].sa\_pos] \leq i\) (line [line:lc46check32lc]), then \(j\) is a last c-candidate and we add \(R[c].text\_pos\) to the suffixient set we have built up to now. On the other hand, if \(R[c].sa\_pos \leq \mathrm{PSV}[i]\) or \(j\) is a c-candidate and \(i < \mathrm{NSV}[R[c].sa\_pos]\) (line [line:lc46check32candidate]) we set \(i\) as a c-candidate (lines [line:lc46candidate32true] and [line:lc46update]), otherwise it is set as a non-\(\textit{c-candidate}\) (lines [line:lc32non32candidate] and [line:lc46update]). Now, we update the remaining attributes of \(R[c]\) according to \(i\)’s values: \(R[c].sa\_pos = i\), \(R[c].text\_pos = n - \mathrm{SA}_r[i'] + 1\). Finally, if \(i\) is the last \(c\)-run break, it is added to the output suffixient set in the final loop if it is a c-candidate (line [line:lc46eval]). Therefore, at the end of the algorithm execution, the output set \(\mathcal{S}\) consists of all (and only the) last candidate positions. Then, Algorithm 1 correctness follows from correctness of quadratic algorithm [4], [6], and Proposition [prop:last-candidate].
The arrays \(\mathrm{BWT}_r\), \(\mathrm{LCP}_r\), \(\mathrm{SA}_r\), \(\mathrm{NSV}\), and \(\mathrm{PSV}\) consume \(O(n)\) space. The data structure \(R\) uses \(O(\sigma)\) space and the set \(\mathcal{S}\) is also of size \(O(n)\). Thus, the space usage of Algorithm 1 is \(O(n)\). With regard to running time, computing \(\mathrm{BWT}_r\), \(\mathrm{LCP}_r\), \(\mathrm{SA}_r\) in line [line:lc46compute32arrays] takes \(O(n)\) time. \(\mathrm{PSV}\) and \(\mathrm{NSV}\) in line [line:lc46compute-sv] can also be computed in \(O(n)\) time [16]. Scanning \(\mathrm{BWT}_r\) takes \(O(n)\) time, and we perform \(O(1)\) operations for each \(\mathrm{BWT}_r\) position. The total running time is then bounded by \(O(n)\).0◻ ◻
Cenzato et al. [5] give two linear-time algorithms to compute a smallest suffixient set. One is based on the \(\mathrm{LF}\) array (so we refer to it as \(\mathrm{LF}\) algorithm) and the other is based on \(\mathrm{PSV}/\mathrm{NSV}\) arrays (and we call it FM algorithm). As pointed out in the paper, the FM algorithm is faster than the LF algorithm, at the cost of increasing space usage. As shown in Table 3, the \(\textrm{PLC}\) algorithm is neither faster than the FM algorithm nor uses less space. Interestingly, with a slight modification, we can turn it as fast as the FM algorithm using as little space as the LF algorithm, getting the best performance of both. We refer to the new algorithm as LC.
By the form in which Algorithm 1 works, we can see that it never asks for \(\mathrm{PSV}\)/\(\mathrm{NSV}\) values at positions that are not run-breaks. So, there is no need to compute \(\mathrm{PSV}\)/\(\mathrm{NSV}\) for other positions. Further, we can compute the needed values on the fly as we find the run-breaks, using a stack-based strategy similar to the one of Bonizzoni et al. [7]. In particular, we use two stacks, \(psv\) and \(nsv\). In \(psv\), we push each \(i\) as we iterate from \(i = 2\) to \(i = n\), that is, every \(i\) can be the \(\mathrm{PSV}\) value for some run-break. Then, for each run-break \(i\) we pop elements from \(psv\) until \(\mathrm{LCP}_r[psv.top()] < \mathrm{LCP}_r[i]\) or \(psv\) is empty, so \(\mathrm{PSV}[i] = psv.top()\) or \(\mathrm{PSV}[i] = 1\), respectively. Since we compute each needed \(\mathrm{PSV}\) value when it is requested, we store it in a variable \(psv_i\), so there is no need to create the array \(\mathrm{PSV}\), which reduces running time and space usage. On the other hand, in \(nsv\) we push a \(c\)-run break \(j\) to mean it is still pending to find its \(\mathrm{NSV}\) value, which will not be required until the next occurrence of a \(c\)-run break. At the beginning of iteration \(i\), while \(\mathrm{LCP}_r[i] < \mathrm{LCP}_r[nsv.top()]\) we set \(R[\mathrm{BWT}_r[nsv.top() - 1]].nsv = R[\mathrm{BWT}_r[nsv.top()]].nsv = i\) (because we had defined \(R[\mathrm{BWT}_r[j - 1]].sa\_pos = R[\mathrm{BWT}_r[j]].sa\_pos = j\) when we processed the run-break \(j\)) and call \(nsv.pop()\) (because we already found \(nsv.top()\) \(\mathrm{NSV}\) value). Since we store \(\mathrm{NSV}\) values only for the run-breaks, we can store those values in the data structure \(R\), adding a new field \(R[c].nsv\), which is set \(R[c].nsv = n + 1\) as a default value when we update \(R[c]\) on line [line:lc46update] until we compute its real value (which can be also \(n + 1\)). Then there is also no need to create an array \(\mathrm{NSV}\), which again reduces running time and space usage. Algorithm 2 applies these modifications to Algorithm 1.
4 Experimental Results↩︎
In this section we empirically compare our new algorithm with previous ones that are implemented.
The experiments were run using a workstation Intel(R) Pentium(R) IV, CPU @ 2.4 GHz with 8 cores, 10 MB cache, 256 gigabytes of RAM, 860 gigabytes of disk running Devuan GNU/LINUX 2.1. We compared the running time and the space consumption of algorithm
LC against PLC, FM, and LF, computing smallest suffixient set over DNA sequences from the Pizza&Chilli collection2. All the implementations of those algorithms are publicly available at https://github.com/regindex/suffixient. We used the elapsed time and the maximum resident set size metrics returned by the utility usr/bin/time3 as the running-time and the space consumption of the algorithms, respectively. We include in the construction cost the time to build the suffix array and other structures from the text. We also measure the performance of
the online algorithm [8] (called Onl, implementation is available upon request) as a way to compare a different paradigm to compute suffixient
sets.
The results on DNA are shown in Figure 3. We can see that LC used the same space than LF, 3.1–3.3 times less space than PLC and FM, and 17–19 times less space than Onl (except on dna, where the ratios are 2.3 and 11, respectively). On the other hand, LC was 14%–39% faster than LF, 6%–32% faster than PLC, and 6%–24% faster than FM. Compared to Onl, LC gave mixed results, being from 30% faster to 48% slower. From those results, we conclude that LC outperforms all the other algorithms in both time and space, with the exception of the running time of Onl, where the comparison is mixed. It must be noted that the good time performance of Onl comes at the price of using a large amount of space (nearly 20 times that of LC) which is prohibitive in many practical settings.
Figure 4 shows the case of some non-DNA representative texts: natural language and source code. The results stay in the same line as on DNA. We can see that Onl performs noticeably better when \(\chi\) is small, while the other constructions (including ours) are much more sensitive to \(n\).
5 Conclusions and Future Work↩︎
We have proposed the first algorithm that is (i) one-pass (over the suffix-array domain), (ii) linear-time, and (iii) implemented, using the same space as the linear-time implementation using least space [5] and faster than the fastest implementation in that paper. Our algorithm also shows to be more practical than an online algorithm [8], the only other implementation we know of, in the sense that uses much less space (by a factor around 20), being faster on some inputs and slower on others.
As a future line of research, we plan to evaluate if algorithm FM [5] admits the same transformation we made on PLC to obtain LC: computing \(\mathrm{PSV}\) and \(\mathrm{NSV}\) values on the fly and just for run-breaks. That would yield a one-pass version of FM as well. We believe that, with that transformation, FM would be even faster than LC while not using more space.
Being one-pass opens the door to running the construction in streaming mode, where the suffix-array related structures of the reverse text are visited sequentially and not maintained in main memory. A relevant measure for such a streamed construction is the amount of main memory needed. Our algorithm needs only \(O(\overline{r})\) main memory space for NSV because it is built only on run breaks, and \(O(\sigma)\) for other structures, but it stores an amount of PSV array cells that is proportional to the height of the suffix tree. While small in practice, this height can be up to \(\Theta(n)\) in pathologic cases. Reducing this extra space, say to \(O(\overline{r})\) as well, would yield a streamed construction whose main memory usage is sensitive to the compressibility of the text.
We also plan to code LC (and eventually the one-pass version of FM) on top of PFP [11]. Cenzato et al. [5] experiment with running the algorithms over the concatenation of many copies of large DNA sequences, producing huge highly-repetitive files. In such context, the PFP version of the one-pass algorithm performed much better than any of the linear-time algorithms. We believe the PFP-based version of LC can be even faster in that setting.