Optimizing the Performance of the FM-Index for Large-Scale Data

A string $T$ is a finite sequence of symbols $T=T[1..n]=T[1]\mathrm{\cdots }T[n]$ over an alphabet $\mathrm{\Sigma }=\{c_1,\mathrm{\dots },c_{\sigma }\}$ whose symbols can be unambiguously ordered. We denote $\epsilon$ as the empty string, $|T|$ as the length of $T$, $\sigma$ as the size of $\mathrm{\Sigma }$, and the rotation of $T$ is a transformation where the characters are circularly shifted so that the first character becomes the last character of the rotated string. A bitvector $B$ is defined as a finite sequence $B=B[1..n]=B[1]\mathrm{\cdots }B[n]$ over the alphabet $\mathrm{\Sigma }=\{0,1\}$. Unlike standard arrays, which typically require $\mathrm{\Theta }(n\log \sigma )$ bits to store $n$ elements on an alphabet of size $\sigma$, a bitvector uses only $n$ bits, making it a highly space-efficient representation when working over binary data.

We define the substring of $T[i]\mathrm{\cdots }T[j]$ of $T$ as $T[i..j]$ starting at position $i$ and ending in position $j$, with $T[i..j]=\epsilon$ if $i>j$. For a string $T$ and $1\le i\le n$, $T[1..i]$ is called the $i$-th prefix of $T$, and $T[i..n]$ is called the $i$-th suffix of $T$. We call a prefix $T[1..i]$ of $T$ a proper prefix if . Similarly, we call a suffix $T[i..n]$ of $T$ a proper suffix if . Given a set of strings $𝒮$, $𝒮$ is prefix-free if no string in $𝒮$ is a prefix of another string in $𝒮$.

We use $\prec$ to denote the lexicographic ordering of two strings, which is defined as follows. For two strings $T_2[1..m]$ and $T_1[1..n]$, $T_2\prec T_1$ if $T_2$ is a proper prefix of $T_1$, or there exists an index $1\le i\le min(n,m)$ such that $T_2[1..i-1]=T_1[1..i-1]$ and . We define the co-lexicographic order of two strings, T1 and T2, as the lexicographic order obtained by reading the strings from right to left. We assume that the string $T$ is terminated by a special symbol $\$$ that is lexicographically smaller than any other symbol in $\mathrm{\Sigma }$, unless otherwise specified. This provides the technical convenience of allowing all rotations of a string and all suffixes of the same string to have coinciding lexicographic ordering.

We define a count query and a locate query for some pattern $P$ as follows. count($P$) returns the number of occurrences of the pattern in $T$. locate($P$) returns all positions in $T$ where this pattern occurs. These queries can be answered without access to the text by using an index of the text.

Throughout this survey, the time and space complexities will be specified in the RAM model of computation, where programs are run on a random-access memory where words of $w=\log n$ bits are accessed and manipulated in constant time. Therefore, we will generally measure space in words as well. Moreover, in this model, we assume that all ordinary arithmetic and logical operations are performed in constant time.

For a string $T$ of length $n$, a symbol $c\in \mathrm{\Sigma }$, and a position $i$ in $T$, the rank function $T.{\mathrm{𝗋𝖺𝗇𝗄}}_c(i)$, is defined as the number of occurrences of $c$ in $T[1..i]$. Given an integer $i$, the select function $T.{\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}}_c(i)$, returns the position of the $i$-th occurrence of $c$ in $T$, if it exists.

For a bitvector $B$, we will instead define the rank function $B.\mathrm{𝗋𝖺𝗇𝗄}(i)$, as the number of occurrences of $1$ in $B[1..i]$, and the select function, $B.\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}(i)$, as the position in $B$ of the $i$-th occurrence of $1$. We note that for a bitvector we automatically assume that the character of interest is $1$.

The suffix tree [39, 24, 1] is a compressed trie storing the suffixes of a string $T$, where paths with a single child are compressed into a single node with the edge labeled by the concatenation of the compressed edge symbols. As there are $n$ suffixes and no suffix can be repeated twice, there are exactly $n$ leaves. Furthermore, since there are no 1-degree nodes and $n$ leaves, the suffix tree can be stored in the $𝖮(n)$ space. Assuming that $T$ is defined by a constant-size alphabet or an alphabet consisting of integers in a polynomial range, it can also be built in $𝖮(n)$ time  [39, 24, 38, 12].

The suffix array [23] (denoted as SA) was introduced to reduce the space consumption of the suffix tree. For a given string $T$ of length $n$, the SA is the array $\mathrm{𝖲𝖠}=[1..n]$ defined such that $\mathrm{𝖲𝖠}[i]$ is the text position of the $i$-th suffix of $T$ in lexicographical order. The theoretical space required by the SA is still $𝖮(n)$, but only a single integer is stored for each suffix, so the practical space is much smaller than that of the suffix tree. Another benefit of SA is that it is possible to build without having to build the suffix tree in $𝖮(n)$ time [19].

Both the suffix tree and the suffix array can be used to match a pattern $P=[1..m]$ to all occurrences in $T$. In the suffix tree, the process involves traversing the nodes from the root symbols matching $P$ to $T$, eventually ending at some node or on the edge between two nodes. All paths that follow from this point represent a suffix that starts with the pattern $P$ and as such the pattern has been matched. This is done in $𝖮(m)$ time assuming that the child node that extends the path can be located in constant time through the use of hashing or some $𝖮(1)$ lookup data structure. Given the suffix array of $T$, the query pattern can be found by performing a binary search over the suffix array. The binary search takes $𝖮(m\log n)$ time and returns an interval in the suffix array that represents all occurrences of the query pattern in $T$. A binary search can be used because the suffixes are sorted lexicographically, so suffixes that start with the same prefix are grouped together. This bound can be reduced to $𝖮(m+\log n)$ with the use of additional data structures  [23].

We now define the inverse of SA as inverse suffix array (ISA). Given a suffix array of $T$, $\mathrm{𝖲𝖠}$, the corresponding inverse suffix array, $\mathrm{𝖨𝖲𝖠}$, is defined as $\mathrm{𝖨𝖲𝖠}[i]=j$ if and only if $\mathrm{𝖲𝖠}[j]=i$. In other words, this structure returns the position of a suffix in the SA given its starting position in the text. Furthermore, many applications, such as FM-index, that make use of SA do not store the entire SA and instead store a sample of SA values, which is called the sampled SA.

The Burrows-Wheeler Transform ($\mathrm{𝖡𝖶𝖳}[1..n]$) [8] of a string $T[1..n]$ is a reversible permutation of the characters of the string. This means that it is possible to retrieve the original string, $T$, from its BWT. Let us denote the matrix of all sorted rotations of $T$ as ${\mathrm{M}}_{\mathrm{𝖡𝖶𝖳}}$, and F and L as the First and Last column of this matrix, respectively. By definition, L is the BWT of $T$ and F contains all the characters in $T$ in lexicographic order. If we were to only consider for each row of ${\mathrm{M}}_{\mathrm{𝖡𝖶𝖳}}$ from the first character of each rotation to the special character $ (the last character in $T$), we would have the suffixes of the original input in lexicographic order. That is exactly what is stored by the SA and this insight allows for the BWT to be constructed from the SA. More formally, this is defined as $\mathrm{𝖡𝖶𝖳}[i]=T[n]=\$$ if $\mathrm{𝖲𝖠}[i]=1$ and $\mathrm{𝖡𝖶𝖳}[i]=T[\mathrm{𝖲𝖠}[i]-1]$ otherwise. More precisely, the BWT is storing the characters that occur before each suffix in lexicographic order. The BWT can be built directly in $𝖮(n)$ time and $𝖮(n\log \sigma$) bits of space  [25], or built from the SA of $T$ using the formulas above.

A key property of the BWT is that the $i$-th occurrence of a character in L and the $i$-th occurrence of the same character in F are the same character in $T$. This property comes from the fact that the BWT is constructed by lexicographically ordered cyclic rotations. If we consider the row of the first occurrence of a character $c$ in L, the character, say $c^{\prime }$, in F in that same row is the character that follows $c$ in the original input. Since this is the first occurrence of $c$, this character, $c^{\prime }$, is also the lexicographically smallest character that follows any occurrence of $c$. During the sorting of F all $c$ occurrences are grouped together and their order is determined by what follows them and since $c^{\prime }$ is the lexicographically smallest, the first occurrence of $c$ in L is the first occurrence of $c$ in F. Next, we denote $𝖢[c]$ as the number of suffixes starting with a character smaller than $c$, which can be derived from F by counting the number of occurrences of each character. We now define the LF-mapping of the BWT as follows:

The LF-mapping function relates the row of the $i$-th occurrence of a character in L to the row of the $i$-th occurrence of the same character in $F$ and makes it possible to reconstruct the string $T$ from its BWT. To reconstruct $T$, we set an iterator $t=1$ and $T[n]=\text{BWT}[1]$ and for each $i=n-1,\mathrm{\dots },1$ do $t=\mathrm{𝖫𝖥}(t)$ and $T[i]=\mathrm{𝖡𝖶𝖳}[t]$. The way this works is that the LF-mapping is used to find the row of ${\mathrm{M}}_{\mathrm{𝖡𝖶𝖳}}$ where F contains the same character as $𝖫[t]=\mathrm{𝖡𝖶𝖳}[t]$ and since BWT stores the character that occurs before the character in F in $T$, then BWT in this row gives us the next character in reconstruction. The LF-mapping is essential in the FM-index and all related indexes, as will be shown.

The FM-index by Ferragina and Manzini  [14] is a compressed full-text index and the first self-index, an index that also replaces the text entirely. Formally, it makes use of the BWT and the SA to efficiently support count and locate queries, essential queries in pattern matching. The BWT of the input string supports count queries through the use of LF-mapping in a technique called backward search [14]. Backward search is formally defined as follows.

Using the result of the backward search, a count query can be answered simply by providing the size of the range, or $ep-sp+1$. Note that none of the steps of the backward search rely on the complete BWT matrix ${\mathrm{M}}_{\mathrm{𝖡𝖶𝖳}}$, it is sufficient to retain only the BWT, which is commonly stored in compressed form, and the array $𝖢[1..\sigma ]$. It is also easy to see that the backward search algorithm takes $𝖮(|P|)$ time, assuming that each $LF(i)$ computation can be performed in constant time. Consequently, a count query can be answered by the FM-index in $𝖮(|P|)$ time under this assumption [14, 13]. The SA of the input string enables all occurrences to be efficiently located. A range in BWT directly corresponds to a range in SA, which means that locating the occurrences would be trivial with the full SA. However, in FM-index only a sampled SA is stored. We formally define the process of locating occurrences with a sampled SA as follows.

The FM-index stores SA values at regular positions in the text $T$. Specifically, for some sampling parameter, $s$, each position $\mathrm{𝖲𝖠}[i]$ is kept if $\mathrm{𝖲𝖠}[i]mods\equiv 0$. As a result, locating each of the $occ$ occurrences requires traversing up to $s$ LF-mapping steps, yielding a total time complexity of $𝖮(occ\cdot s)$. In the original implementation, $s={\log }^{1+ϵ}(n)$ for some $ϵ>0$. The FM-index occupies $\frac{5n{H}_K}{\log n}+𝖮(\frac{n\log \sigma }{\log n})$ words of space.

The wavelet tree is a data structure by Grossi et al. [16] that can be used to store a string $T$ over an alphabet $\mathrm{\Sigma }$ of size $\sigma$ achieving $\frac{n\log \sigma }{\log n}+𝗈(\frac{n\log \sigma }{\log n})$ words of space and $𝖮(\log \sigma )$ time $\mathrm{𝗋𝖺𝗇𝗄}$ and $\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}$ operations over the input text. To achieve these bounds, we assume that the underlying bitvectors of the wavelet tree support constant-time $\mathrm{𝗋𝖺𝗇𝗄}$ and $\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}$. A wavelet tree is made up of $\log \sigma$ bitvectors organized in a balanced binary tree structure. The bitvector at the root, $B_{root}$, is the same size as the input $T$ and is defined by choosing a character $v\in \mathrm{\Sigma }$, usually the median character in lexicographic order, such that $B_{root}[i]=0$ if and $B_{root}[i]=1$ otherwise. The positions marked as $0$ form the subsequence that is used for the left child and the positions marked as $1$ form the subsequence for the right child. The bitvectors for each of these children are created in the same way as the root, keeping in mind that each child’s subsequence is on half the alphabet of the parent sequence. This partitioning is performed recursively until the bitvector of a node consists entirely of $0$s or $1$s, which denotes that this node is a leaf. Answering the query $T.{\mathrm{𝗋𝖺𝗇𝗄}}_c(i)$ is done by traversing the tree and computing $\mathrm{𝗋𝖺𝗇𝗄}$ on the bitvector at each level. We initialize $po{s}_0=i$ and continue down the left subtree with $po{s}_{j+1}=B.{\mathrm{𝗋𝖺𝗇𝗄}}_0(po{s}_j)$ if $c\prec v$ or down the right subtree with $po{s}_{j+1}=B.{\mathrm{𝗋𝖺𝗇𝗄}}_1(po{s}_j)$ otherwise. The query stops when a leaf is reached and the result is $po{s}_h$, where $h$ is the height of the tree. The query $T.{\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}}_c(i)$ is answered by first identifying the leaf node for $c$, which can be done using a rank query. Then, the tree is traversed from the leaf node to the root through $\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}$ queries on every node’s bitvector. We initialize $po{s}_h=i$ and if the current node is the left child of its parent, we traverse with $po{s}_j=B.{\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}}_0(po{s}_{j+1})$; otherwise, we traverse with $po{s}_j=B.{\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}}_1(po{s}_{j+1})$. The result of the query is then $po{s}_0$, which we get once we reach the root.

Prefix-Free Parsing (PFP) is a technique created by Boucher et al. [5] that compresses and replaces the input before constructing the index. This allows large data sets, such as pangenomes, to be indexed without storing the full input, reducing practical memory requirements. PFP takes as input a string $T$ of length $n$, and two parameters, which we denote $w$ and $p$. It produces a parse and a dictionary of $T$, denoted as P and D, respectively, which can be used to reconstruct the original input. The dictionary consists of unique overlapping phrases stored in their lexicographic ordering, and the parse gives the order in which these phrases appear in the input. The Prefix-Free in the name of the algorithm is because the output has the property that none of the suffixes of length greater than $w$ of the phrases in D is a prefix of any other. This concept is formalized in the following lemma. We refer to PFP of $T$ as $\text{PFP}(T)$, consisting of the dictionary D and the parse P.

In order to construct these data structures, PFP first appends $w$ copies of $\$$ – a special character that is lexicographically smaller than any element in $\mathrm{\Sigma }$ – to both the beginning and end of $T$. For the sake of explanation, we consider the string $T^{\prime }={\$}^{w}T{\$}^w$. Next, we define the set of trigger strings E, which determine how to parse $T^{\prime }$. Given a parameter $p$, we construct the set of trigger strings by computing the Karp-Rabin hash, $H_p(t)$, of substrings of length $w$ by sliding a window over $T^{\prime }$. E is then the set of substrings $t=T^{\prime }[i..i+w-1]$, where $H_p(t)\equiv 0$ or $t={\$}^w$. Then, using E, we can scan through $T^{\prime }$, store every substring that starts and ends with a trigger string in the dictionary, and then only keep the unique phrases. The parse is then a string of integers, where each integer is a position in the dictionary, and the string indicates the order in which these phrases appear in $T^{\prime }$. Note that by definition of the trigger strings and by appending $w$ copies of $\$$, that $T^{\prime }$ starts and ends with a trigger string. PFP can be used as a preprocessing step to build data structures such as the BWT, SA, and LCP.

Mäkinen and Navarro introduced a major improvement in the space requirements of the FM-index in 2004 [26, 22] with the RLFM index. The key idea behind the RLFM index is to exploit the structure of the run-length encoded Burrows–Wheeler Transform (RLBWT). The RLBWT is a data structure that encodes the BWT in $𝖮(r)$ words, where $r$ is the number of runs in the BWT, and supports $LF(i)$ queries in $𝖮(\log {\log }_w(\sigma +n/r))$ time. Here, a run in a string is defined as a maximal sequence of consecutive identical characters; thus, the number of runs corresponds to the number of such maximal contiguous blocks in the BWT. For example, if $\text{BWT}=\$aaacgggttacc$, then this BWT has seven runs. The construction of the BWT matrix groups similar prefixes together, indicating shared context in the input. This common context leads to the repeating of the same letter for several rows in BWT and is especially true in highly repetitive data where the number of runs, $r$, can be significantly less than the number of characters in the BWT. The RLBWT can be constructed indirectly or directly. The indirect construction of the RLBWT, which is the type of construction done in the RLFM index, is done by first constructing the BWT and then run-length encoding it. This can be done $𝖮(n)$ time and $𝖮(|\mathrm{𝖯𝖥𝖯}|)$ space [5] or $𝖮(n)$ time and $𝖮(\frac{n\log \sigma }{\log n})$ space[3, 25]. Direct construction of the RLBWT can be done in optimal time ($𝖮(n+r\log r)$ time) and BWT-runs bounded space ($𝖮(r)$ space)[29].

Taking advantage of the compressibility of BWT through run-length encoding, the RLFM index achieves substantial space savings while maintaining efficient count query performance. Mäkinen and Navarro run-length encode the BWT by replacing it with a bitvector, $B$, of length the size of the input that contains a $1$ at the start of every run in the input, and a string, $S$, which lists all the characters that make up each run in order. Consider the earlier example where $\text{BWT}=$ $aaacgggttacc. This sequence is replaced by two components: a bitvector $B=$1100110010110, where each 1 denotes the start of a new run, and a string $S=$ $acgtac, which stores the distinct characters corresponding to each run in the order of appearance in BWT. The authors demonstrate that the LF-mapping can still be computed using only these representations. To support this, a new bitvector $B^{\prime }$ is derived by performing a stable sort of the runs encoded in $B$, where the sorting is based on the lexicographic order of the associated character of each run in $S$, and stability ensures that the original relative order among runs of the same character is preserved. For the given example, this results in $B^{\prime }=$ 1100111010010. Additionally, instead of storing $S$ directly, a bitvector, $S_c$, is stored for each character $c$ storing a 1 at every index where $c$ occurs in $S$. In our example, we would then have the following bitmaps: $S_{\$}=$ 1000000, $S_{𝚊}=$ 0100010, $S_{𝚌}=$ 0010001, $S_{𝚐}=$ 0001000, and $S_{𝚝}=$ 0000100. Mäkinen and Navarro compute the LF-mapping in the RLFM index by presenting an equivalent equation to $𝖢[c]+\mathrm{𝖡𝖶𝖳}.{\mathrm{𝗋𝖺𝗇𝗄}}_c(i)$ that uses $B$, $B^{\prime }$, and the $S_c$ bitmaps. We denote that ${𝖢}_S$ is the same as the array $𝖢$, but over the string $S$ rather than the input. The authors present the following lemma.

To prove this lemma, the authors establish that since ${𝖢}_S[c]$ gives the number of runs of characters smaller than $c$ and $B^{\prime }$ contains the runs sorted lexicographically, then $B^{\prime }.\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}({𝖢}_S[c]+1)$ is the position in $B^{\prime }$ corresponding to the first run of $c$. It follows that $B^{\prime }.\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}({𝖢}_S[c]+1)-1$ is the sum of the run lengths of all characters smaller than $c$, which is equivalent to the number of characters smaller than $c$ in BWT, or $𝖢[c]$. The authors used this to then give the equation for the LF-mapping in RLBWT as follows.

The proof starts with the fact that $B.\mathrm{𝗋𝖺𝗇𝗄}(i)$ gives the position of the run that finishes at $i$ in $S$, which means that $S_c.\mathrm{𝗋𝖺𝗇𝗄}(B.\mathrm{𝗋𝖺𝗇𝗄}(i))$ is the number of runs of $c$ in the BWT from $1$ to $i$. Therefore, $C_S[c]+1+S_c.\mathrm{𝗋𝖺𝗇𝗄}(B.\mathrm{𝗋𝖺𝗇𝗄}(i))-1$ points to a position in $B^{\prime }$ within the runs of $c$ or to the character immediately after, which means that $B^{\prime }.\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}(C_S[c]+1+S_c.\mathrm{𝗋𝖺𝗇𝗄}(B.\mathrm{𝗋𝖺𝗇𝗄}(i)))-1$ is the starting position for the runs of $c$ plus the sum of the run lengths of $c$ in BWT$[1..i]$. This is the same as $C[c]+S_c.\mathrm{𝗋𝖺𝗇𝗄}(B.\mathrm{𝗋𝖺𝗇𝗄}(i))$. If BWT$[i]\ne c$, then this is the LF-mapping function. Otherwise, it is $B^{\prime }.\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}(C_S[c]+1+S_c.\mathrm{𝗋𝖺𝗇𝗄}(B.\mathrm{𝗋𝖺𝗇𝗄}(i)))+i-B.\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}(B.\mathrm{𝗋𝖺𝗇𝗄}(i))$, because the position we are looking for is within the run $c$, so the added term of $i-B.\mathrm{𝗌𝖾𝗅𝖾𝖼𝗍}(B.\mathrm{𝗋𝖺𝗇𝗄}(i))$ gives how far into the run $i$ is.

Mäkinen and Navarro presented an improved FM-index that when performing solely count queries required less space than the FM-index and achieved the same query time. The authors tested their implementation against the implementation of the FM-index, the original SA, the Compressed Suffix Array (CSA) [36], and the LZ self-index [27] using random patterns from the index text. Their index required 0.63 of the original input size, which is less than every competing method except for CSA with tuned parameters. In terms of count query time, RLFM was only slower than the implementation of the FM-index, requiring around 0.1 msec to count occurrences for the longest pattern length compared to around 0.05 msec. In order to answer locate queries, RLFM uses the same SA sampling techniques as the FM-index, which requires more than $𝖮(r)$ words of space. The RLFM with locate therefore requires more space than CSA and less than the FM-index, but had comparable locate query times across all pattern lengths.

Although count queries could be performed in $𝖮(r)$ space with the RLFM index, it was not until the development of the r-index index by Gagie et al. [15] that locate queries could also be supported within $𝖮(r)$ words of space. This is because the SA values are needed to answer a locate query, and storing the full SA requires $𝖮(n)$ space. Even when sampling the SA values at a constant rate, as was done in the original FM-index, the space required would be linearly proportional to the input size. Policriti and Prezza’s Toehold Lemma [35, 34] enabled SA values to be sampled once per run while still supporting locate queries. The lemma is as follows.

The Toehold Lemma states that by storing a predecessor data structure in $𝖮(r)$ words containing SA values at the start and end of each run in the BWT, it is possible to maintain a “toehold” – that is, an SA value for at least one position within the current BWT interval – at every step of the backward search. To prove this, assume that a pattern $P$, so far matched with $i+1..m$, the interval BWT$[s_p,e_p]$ has been found, and SA$p$ is known for some index $p$ in the BWT interval. If the next character that matches $P[i]$ is the same as BWT$[p]$, then the value of the SA for some index $p^{\prime }$ in the new interval $[s_p^{\prime },e_p^{\prime }]$ can be found by setting $p^{\prime }=\mathrm{𝖫𝖥}(p)$ and $\text{SA}[p^{\prime }]=\text{SA}[p]-1$. If instead BWT$[p]\ne P[i]$ and $P$ occurs in the text, then there must exist an index $u$ such that $s_p\le u\le e_p$ and BWT$[u]=P[i]$, as well as at least one index in $[s_p,e_p]$ for which BWT differs from $P[i]$. This implies that the next interval intersects a run boundary and, therefore, a sampled SA value is stored for this position. Finding this SA value requires a predecessor query to the earlier defined predecessor structure that can be completed in $𝖮(\log \log (n/r))$ time. Therefore, the lemma is proven using induction. Improvements to the Toehold lemma by Bannai et al. [2] show that the SA only needs to be sampled at run starts for the proof to remain valid, which helps improve practical performance but the $𝖮(r)$ words of space bound remains.

In their paper, Gagie et al. show that knowing a position in the BWT and its SA value (through Lemma 3.2), it is possible to compute the SA values of neighboring rows in $𝖮(\log \log (n/r))$ time. To do this, the authors used the definition of the $\varphi$-array of Kärkkäinen et al. [18].

The presented lemma states that by storing an $𝖮(r)$ space predecessor data structure it is possible to evaluate $\varphi$ and ${\varphi }^{-1}$ in $𝖮(\log \log (n/r))$ time. The predecessor data structure stores for every position $i$ in the text corresponding to a start of the run the positions in the text (the SA values) of the preceding and following position of the position $i$’s corresponding position in BWT, $q$. Essentially, the predecessor data structure, $N$, can be defined as follows for every $i$: for $i=\text{SA}[q]-1$. To then find the SA values neighboring for a given position, $k$, in the text (equivalently position $p$ in the BWT), first a predecessor query is done to find the stored position $i$ that starts the run that contains $k$ to get the entry . The neighboring values are defined as $\text{SA}[p-1]=\text{SA}[q-1]+k-i$ and $\text{SA}[p+1]=\text{SA}[q+1]+k-i$. The proof for this relies on the fact that contiguous positions in the BWT that are in the same run tend to remain contiguous after applying LF-mapping. This lemma shows that if a count query on a pattern returns a BWT interval and a SA value is known for a position in this interval, it is possible to identify the SA values for every position in the interval. Lemma 3.2 guarantees that a position in the BWT interval is always known. This allows for locate queries to be answered with queries to $𝖮(r)$ space predecessor data structures in $𝖮(\log \log (n/r))$ time. Gagie et al. also establish that it is possible to increase the size of index in order to achieve optimal time count and locate queries.

The r-index combines of the SA and BWT samples sampled at run boundaries requiring $𝖮(r)$ words of space, and supporting count queries in $𝖮(m\log \log (n/r))$ time and locate queries in $𝖮((m+occ)\log \log (n/r))$ time, where $m$ is the length of the query pattern and $occ$ is the number of occurrences of the pattern in the input. In their experiments, r-index outperforms all competing implementations of the FM-index with respect to space-time trade-offs, taking at most the same amount of space and being magnitudes faster in locate queries. Compared to indexes based on Lempel-Ziv compression or grammars, r-index uses more space but remains significantly faster for locate queries.

The subsampled r-index, or sr-index, was developed by Cobas et al.  [9] to reduce the space requirement of the r-index, which can be large compared to Lempel-Ziv-based indexes when $r$ is large. The main principle of sr-index is to reduce redundancy by removing run samples that are within a text distance of at most $s$ from another sample. Here, a run sample refers to a value stored for the start of a run in the BWT, whether this is a bit in a bitvector or a SA value. Doing this effectively reduces the space required by r-index from $𝖮(r)$ to $𝖮(min(r,n/s))$. Subsampling is performed by iterating over each run sample $t_i$ for $i=2,\mathrm{\dots },r-1$ and checking whether the distance between the run samples $t_{i+1}$ and $t_{i-1}$ is at most $s$ – that is, whether $t_{i+1}-t_{i-1}\le s$. If the condition is met, then the run sample $t_i$ is removed from the data structures in the r-index that stores information about it. In this way, all run samples are removed so that if the text were segmented into blocks of size $s+1$, there would be no more than two run samples in each block. The sr-index additionally introduces a bitvector of size $𝖮(r)$ bits, marking each removed run sample with a $1$. The space overhead of this bitvector is negligible compared to the savings gained by removing entries from the r-index data structures, as it is stored in bits rather than full words. Traditionally, the r-index relies on storing SA values at every run boundary to support locate queries. As a result, removing some of these run samples requires a modification to the way locate queries are performed. In particular, the process of maintaining a toehold SA value for every interval and the process of finding the SA values of all positions in a final BWT interval changes. For example, computing the SA value for the last position of a run is trivial if the run sample for the run was not removed by sub-sampling. However, in the case where the sample was removed, then in order to find the SA for the run boundary the LF-mapping must be applied recursively until a run sample that was not removed is found. The number of steps needed to find a sample is bounded by $s$, as this is guaranteed by the subsampling process.

The sr-index takes $𝖮(m+s)$ time to complete the backward search of a pattern of length $m$, rather than the expected $𝖮(m\cdot s)$ time. The authors managed to do this by taking advantage of the fact that the only SA values needed are the values of the final interval of the search. To do this, the authors keep track of the last time that the next character in the backward search did not match the character at the position in the BWT of the currently known SA value, i.e., BWT$[p]\ne P[i]$ as shown in Subsection 3.2. Then, the final SA value can be computed by subtracting the number of LF steps that have occurred since this mismatch from the SA value for this mismatched interval. The remaining step is to find the SA value for the mismatched interval, which can be located in at most $s$ LF steps, resulting in a query time of $𝖮(m+s)$. To locate all occurrences, the process involves breaking the final interval into maximal runs and using a single LF-mapping per run to calculate the LF-mapping for every position within the run, as the results are contiguous for positions in the same run. From the subsampling process, all runs will be completed in at most $s$ steps, at which point the SA values can be found by applying $\varphi$ to the LF-mapping results found earlier, given the location of all occurrences. The locate query therefore takes $𝖮((m+s\cdot occ)\log \log (n/r))$ time.

As the sampling distance increases, Cobas et al. showed that the space required by sr-index reduces but the query times jump to at least three times the r-index query times. The authors also compared the sr-index with Lempel-Ziv-based indexes. For most data sets, the Lempel-Ziv indexes required less space than sr-index except for a few cases where sr-index took less space when using the largest sampling distance. For every data set, sr-index was significantly faster in querying than the Lempel-Ziv indexes.

The r-index achieves significant space savings by leveraging $\varphi$ and ${\varphi }^{-1}$ to minimize the number of SA values that must be stored. However, in practice, this method of accessing SA is slow because it requires the SA value to be computed for every position between the desired position and the nearest run boundary position. This led to the development of the ${\varphi }^{-1}$ forest by Boucher et al. [6] that aims to improve the time required to access the SA by providing random access to the SA.

The ${\varphi }^{-1}$ forest, which works analogously for $\varphi$, exploits the iterative nature of ${\varphi }^{-1}$ to compute random access to SA. The primary idea behind this auxiliary data structure relates to the toehold lemma 3.2. In particular, $\mathrm{𝖲𝖠}[i+1]-\mathrm{𝖲𝖠}[i]$, or the difference in offsets of two consecutive lexicographical suffixes, does not change as we iteratively set $j=i$, when $\mathrm{𝖲𝖠}[j]=\mathrm{𝖲𝖠}[i]-1$ and $\mathrm{𝖡𝖶𝖳}[i]=\mathrm{𝖡𝖶𝖳}[i+1]$ for the current $j$. Consequently, given any text position, it is possible to compute successive values of ${\varphi }^{-1}$ by iterating over indices and applying the known differences between consecutive SA values, until reaching a run boundary where $\mathrm{𝖡𝖶𝖳}[i]\ne \mathrm{𝖡𝖶𝖳}[i+1]$. Furthermore, when we have a list of SA values at run ends, $ℰ$, and at run starts, $𝒮$, ${\varphi }^{-1}$ can be easily determined for each value in $ℰ$, since each occurs before some value in $𝒮$. Using these concepts, ${\varphi }^{-1}$ can be redefined as follows.

Using this idea, the ${\varphi }^{-1}$-graph is constructed with the nodes being the values of $ℰ$ and edges between nodes occurring when one node is the predecessor of the run following another node (i.e., $\overline{v_i{v}_j}$ exists when $pre{d}_{ℰ}(𝒮[i+1])=ℰ[j]$). Each edge is also labeled with a cost, $c$, and a limit, $l$, which represent the distances from the positions where a change in predecessor occurs. Each cost is defined by the number of steps the start of the next run takes from its predecessor, or $𝒮[i+1]-pre{d}_{ℰ}(i+1)$, and each limit is defined by the distance from the next run, or $suc{c}_{ℰ}(ℰ[i]+1)-ℰ[i]$. An important property of this graph is that every node has a single outgoing edge because the predecessor function is unique, which means that there is no ambiguity in the path to follow when traversing the graph. Boucher et al. show that by following this graph and utilizing the cost and limit values along the edges between nodes on the path, it is possible to determine ${\varphi }^{-1}$ for any given position in the text. Assume that we are given a SA value, $\mathrm{𝖲𝖠}[s]$, then with the ${\varphi }^{-1}$ graph, a total of $m$ iterations of ${\varphi }^{-1}$ can be calculated by traversing a $m$-length path starting from $v_i$, where $v_i=ℰ[i]=pre{d}_{ℰ}(\mathrm{𝖲𝖠}[s])$. The formula is as follows: ${\varphi }^{-m}(\mathrm{𝖲𝖠}[s])=𝒮[i_m+1]=ℰ[i_m]+{\sum }_{k=1}^{m-1}{c}_{i_k}+c_0$, where $c_0=\mathrm{𝖲𝖠}[s]-ℰ[i]$, as long as the requirement is met that from $j=2,\mathrm{\dots },m-1$, or the sum of costs is less than the limit at every step of the path. This condition determines the length $m$ of the path and also maintains the necessary condition that every node on the path is the predecessor of the SA value found by the previous node.

The ${\varphi }^{-1}$-graph is then used to construct the ${\varphi }^{-1}$-forest, which allows iterations of ${\varphi }^{-1}$ to be skipped during the process of finding a missing SA sample in a locate query. To do this, long non-overlapping paths from the ${\varphi }^{-1}$-graph are extracted and a balanced binary search tree is built on each path such that given $\mathrm{𝖲𝖠}[s]$ it is easy to find the largest index in the path for which ${\varphi }^{-m}(\mathrm{𝖲𝖠}[s])$ can be computed following the limit requirement laid out above. In the balanced binary search tree, each node stores a cost and limit value $c$ and $l$, respectively. The leaves of the tree are the edges in the path that the tree represents, and therefore their values are easily determined. The internal nodes then represent their cost as the sum of the cost of their two children and their limit as the minimum between the limit of the left child and the limit of the right child minus the cost of the left child, or $min(l_u,l_v-c_u)$. The tree can then be queried with a given SA sample, $\mathrm{𝖲𝖠}[s]$, by first finding the predecessor $ℰ[i]=pre{d}_{ℰ}(\mathrm{𝖲𝖠}[s])$ to compute $c_0=\mathrm{𝖲𝖠}[s]-ℰ[i]$ and determine the leaf $v$ that corresponds to the edge $e_{i_j}$, or the outgoing edge of $ℰ[i]$. The goal of the query is to find the largest $k$ such that ${\sum }_{x=j}^{m-1}{c}_{i_x}+c_0$ for all $m=j,\mathrm{\dots },k-1$. The algorithm proceeds in two phases: the first involves traversing the tree from $v_i$ and the second involves descending the tree to a leaf node that represents $k$. This fast-forward query can also be trivially extended to include a constraint, $d$, such that the result of the query is the largest $k\le d$. This constraint is important when considering realistic applications in which a specific value within a certain distance from $\mathrm{𝖲𝖠}[s]$ is desired.

The ${\varphi }^{-1}$-forest is composed of several of these binary trees, hence its name. To query the forest, with a given SA position, $i$, the first step is again to determine the predecessor at position, $j$, and then set $d=i-j$ and $c_0=\mathrm{𝖲𝖠}[j]-pre{d}_{ℰ}(\mathrm{𝖲𝖠}[j])$ and $v=e_j$. The trees in the ${\varphi }^{-1}$ forest are then searched to find one that contains $v$. If none is found, ${\varphi }^{-1}$ is computed normally, $v$ is updated to the new node, $d$ is decremented by one, and the trees are searched again for the new $v$. If $v$ is found in a tree, then the tree is queried as described above and from the result $v$ is updated to the new node, $d$ is subtracted from the number of edges traversed in the query, and the trees are searched again for the new $v$. This process continues until $d=0$. In this way, $\mathrm{𝖲𝖠}[i]$ can be computed in $𝖮((i-j)(\log \log (n/r)+\log r))$ time with $𝖮(r)$ additional space to the r-index.

Boucher et al. used their ${\varphi }^{-1}$-forest to augment the r-index and evaluated its performance against traditional r-index in both the Pizza&Chili repetitive corpus [33] and data sets comprising an increasing number of Chromosome 19 sequences. Their results showed comparable SA access times to the traditional r-index on the Chromosome 19 datasets, and up to a sixfold speedup on the Pizza&Chili datasets. Although this improvement came at the cost of slightly increased memory usage compared to r-index, it marked the first method to improve the calculation of $\mathrm{\Phi }$.

Deng et al. [10] introduced FM-Indexing Grammars Induced by Suffix Sorting (FIGISS), which aims to accelerate the FM-index by enhancing the backward search step using grammar compression. FIGISS uses a context-free grammar which we briefly define as follows. A context-free grammar is a formal system of production rules that recursively replace terminal and nonterminal symbols with a nonterminal symbol. This is formalized as grammar $G=(V,\mathrm{\Sigma },R,S)$, where $V$ is the set of nonterminal symbols, $\mathrm{\Sigma }$ is a set of terminal symbols, $R$ is the set of production rules that take the form $A\to \alpha$ with $A\in V$ and $\alpha \in {(\mathrm{\Sigma }\cup V)}^{\ast }$, and $S$ is the first rule $S\in V$. Their method uses Grammar Compression by Induced Suffix Sorting, known as GCIS [32], which creates a leftmost $S$-type suffixes (LMS) factorization of an input text $T$ using the linear-time SAIS algorithm [31]. This algorithm works by first padding the left and right of the input with characters that are lexicographically smaller than all the other characters in the string. Then, it labels each suffix of the input as an $S$-type suffix or an $L$-type suffix. This is done by creating a string the same size as the input, iterating through every position $i$ of the input, and adding an $S$ to the string if the suffix starting at that position is lexicographically smaller than the suffix starting at the following position, $i+1$, or by adding an $L$ otherwise. Suffixes that start at $T[0]$ and $T[n+1]$ are also labeled as $S$-type suffixes. Afterwards, all left-most $S$-type (LMS) suffixes are defined as every $S$-type suffix that occurs right after an $L$-type suffix, where the suffix starting at $T[0]$ is also considered an LMS suffix. The LMS factorization is formed by taking each substring of the created string that starts at one LMS suffix and ends right before the next LMS suffix. Note that the first and last positions of the created string are not included in any substring since they corresponded to the padding characters. Then, each LMS substring is replaced with a nonterminal symbol producing the grammar that is used by FIGGIS.

To perform pattern matching, FIGISS constructs the grammar, ${𝒢}_T$, on the input text $T$, then constructs the grammar, ${𝒢}_P$, on the query pattern $P$, and finally matches $P$ to $T$ by matching their nonterminals. When constructing ${𝒢}_P$, the same nonterminal as ${𝒢}_T$ is used when the right-hand side of an LMS substring in $T^{(1)}$, the LMS factorization of $T$, matches the right-hand side of an LMS substring in $P^{(1)}$, the LMS factorization of $P$. The idea of pattern matching in FIGISS is built on finding a core, which is a maximal length substring in $P^{(1)}$ that is contained in $T^{(1)}$ for every occurrence of $P$ in $T$. Note that the first and last nonterminals $P^{(1)}$ cannot be part of the core because the substrings those nonterminals represent do not appear in the same form in $T^{(1)}$. To find this core, FIGISS constructs a BWT on $T^{(1)}$, represented as a wavelet tree that is constructed using the colexicographic ranking of the nonterminals. Then, a backward search of $p-2$ steps is done, where $p$ is the number of nonterminals in $P^{(1)}$. The resulting interval represents all occurrences of $P^{(1)}[1..p-1]$ to be found in $T^{(1)}$, which correspondes to all occurrences of $P[1..|P|-1]$ in $T$. To then find all occurrences including the last nonterminal, $P^{(1)}[p]$, the backward search is started from the intervals defined by any nonterminal of $T^{(1)}$ that contains the LMS substring represented by $P^{(1)}[p]$ as a prefix rather than starting with the entire range of the BWT. In practice, Deng et al. limit the maximum possible size of any LMS substring with a parameter to improve performance and include a generalized suffix tree on the right-hand side of nonterminals to support pattern matching for queries shorter than the maximum length parameter.

Deng et al. evaluated their implementation with different aggregation parameters against RLFM using datasets from the Pizza&Chili repetitive corpus, a data set of 15 Chromosome 19 sequences and a data set of randomly generated strings. FIGISS requires less space than RLFM for every data set, and for every parameter tested, the space decreases as the parameter increases. When answering count queries, FIGISS is faster than RLFM for some aggregation parameters when the pattern length is particularly long, but otherwise it is slower than RLFM. Finally, FIGISS requires less total memory during index construction and less final index space but takes significantly more time to construct.

One of the main limitations of FIGISS was that a single round of LMS factorization led to short phrases that only helped improve the efficiency of count queries for particularly long patterns. This is because short phrases do not provide significant speedups from the reduction of backwards search steps compared to the increased overhead required unless the query pattern is long. PFP-FM, an accelerated FM-index, developed by Hong et al. [17] similarly accelerates the count query of the FM-index, but allows for acceleration across a wider range of pattern lengths. This is because PFP-FM uses Prefix-Free Parsing (PFP) which uses user-defined parameters to change the length of the phrases that are created. This allows for longer phrases, which provides a greater reduction in the number of backwards search steps than shorter phrases offer. A detailed description of how PFP works can be found in the preliminaries section 2.7.

The index constructed by PFP-FM consists of the FM-index of $T$, the FM-index of the parse $𝒫$, a bitvector $ℬ$, a hash function $h$, and a mapping $ℳ$. $ℬ$ has the same length as the BWT and stores a $1$ for every position for which a rotation of the input starting at that position is prefixed by a trigger string. The hash function $h$ hashes the phrases, and the map $ℳ$ is used to map the hash of each phrase to the lexicographic rank of the phrase. Hong et al. note that with a sufficiently large hash function, it is possible to store $ℳ$ in space proportional to the number of phrases while maintaining a negligible probability of a hash collision. This allows the hash function and the map to effectively replace the dictionary $𝒟$ of the input. Therefore, in order to check if a phrase occurs in the dictionary of the input, the phrase is first hashed and then the map $ℳ$ is queried with the phrase hash. A null return from the map means that the phrase does not occur in the input.

To perform a count query, PFP-FM first parses the pattern $P$, the same way as the input was parsed using the same rolling hash function. The main difference is that the pattern does not have characters appended to it and is not assumed to be circular, which means that the pattern will not necessarily start or end with a trigger string. Therefore, when $P$ is parsed, there will be up to three resulting components: a prefix that ends at the first trigger string in $P$, a sequence of complete phrases that are in $P$, and a suffix that starts at the last trigger string in $P$. First, if there is a suffix component, it is determined whether that suffix occurs in $T$ by performing a backward search for it in the FM-index of $T$. If no BWT range is identified, or if a BWT range is identified and the parsed $P$ consists solely of the suffix, then the query process is considered complete. Otherwise, the next step is to determine whether the sequence of phrases in $P$ occurs in $T$. To do this, the phrases are each hashed with $h$ and the hashes are placed on the map $ℳ$, which returns NULL if the phrase is not in the input, and their phrase ID if it is. If a NULL value is returned for any phrase, then that phrase and, by extension, $P$ do not occur in $T$, and the query is done. Otherwise, the sequence of phrases will be encoded in a sequence of phrase IDs. PFP-FM then performs a backward search for this sequence of phrase IDs in the FM-index of $𝒫$, which means that every backward search step matches a full phrase rather than a single character. The primary acceleration achieved by PFP-FM arises from this mechanism. Continuing with the backward search process, if a BWT range is found and the parsing of $P$ does not have a prefix, the query is complete. Otherwise, the final step is to find the prefix in the parsing of $P$ in $T$ through a backward search in the FM-index of $T$. In practice, Hong et al. demonstrate that the PFP parameters, $w$ and $p$, can be adjusted to systematically optimize the performance gains achieved by PFP-FM, as these parameters directly influence the average phrase length.

The authors performed experiments that compared PFP-FM with FIGISS [10], the standard RLFM[26, 22], RLCSA[37], Bowtie[21], and Bowtie2[20] on increasing size SARS-CoV-2 datasets and the human genome GRCH38. In the SARS-CoV-2 datasets, Hong et al. report that PFP-FM requires less cumulative time to complete 1,000 count queries on patterns of increasing length than all other competing methods. PFP-FM completes queries at least 451% faster than the second-fastest method, FIGISS. The results were the same on the GRCH38 data set, where PFP-FM was between 3.9 and 7.1 times faster than RLCSA, the second fastest method. This improvement comes at the cost of an increase in index size and memory required during construction.

Nishimoto and Tabei  [30] introduced the move data structure in 2021. The primary contribution of their paper was to introduce, OptBWTR, which is the first string index capable of performing locate, count, and decompression in $𝖮(r)$ words of space and optimal time. It also supports extract queries in optimal time with $𝖮(r+b)$ words of space, where $b$ is some sampling parameter, and prefix search in optimal time with $𝖮(r^{\prime }+d)$ words of space, where $r^{\prime }$ is the number of runs in the RLBWT of the input and $d$ is the size of the set of strings used in the prefix search. This outperforms the $r$-index, which can attain $𝖮(r)$ words of space or optimal query times, but not both simultaneously. Additionally, the r-index does not support anything other than locate and count queries. The index of Nishimoto and Tabei achieves this result through the use of the move data structure, which enables the $LF$ and ${\varphi }^{-1}$ functions – both essential for count and locate queries – to be computed in $𝖮(r)$ words of space and constant time. The previous best space-time tradeoff possible for these functions was in the r-index, which takes $𝖮(m\log \log (n/r))$ time for count queries and $𝖮((m+occ)\log \log (n/r))$ time for locate queries with $𝖮(r)$ words of space[15].

The primary idea behind the move data structure is that of disjoint interval sequences and move queries. A disjoint interval sequence is a way to turn a function’s domain and range into a set of non-overlapping input intervals and a set of non-overlapping output intervals. A move query takes a position $i$ from the function’s domain and the input interval that contains this position, $x$, and returns an $i^{\prime }$ and $x^{\prime }$ such that $i^{\prime }$ is a new position in the function’s domain and $x^{\prime }$ is the input interval that contains this new position. The way in which $i^{\prime }$ is determined is by computing $i^{\prime }=q_x+(i-p_x)$, where $p_x$ and $q_x$ are the start of the $x$-th input and output interval, respectively. This defines a bijective function that is related in its form to both the LF and the ${\varphi }^{-1}$ functions. This is significant, as Nishimoto and Tabei demonstrate that dividing the interval sequences to obtain a balanced interval sequence allows move queries to be answered in constant time when the move structure is built on such sequences. Since LF and ${\varphi }^{-1}$ have the same form as the bijective function defined by a disjoint interval sequence and a move query, then by building a move structure on the balanced interval sequence of these functions, it is possible to compute them in constant time.

A move data structure is composed of two arrays, pair and index. The pair stores all pairs from the original interval sequence. In the move data structure for the LF function, these pairs are defined as the BWT positions of all run heads and the LF mapping of these positions. The index array stores the index of the interval that contains the second value of the $i$-th pair. In terms of the LF move structure, this stores the index of the BWT run that contains the LF mapping of the $i$-th run head. For simplicity, we will now assume that the move structure is a table for which at every position $i$ all the information that would be contained in both arrays is stored. A move query on the LF move data structure of some position $i$ known to be in the input interval $x$ is as follows. First, the LF-mapping of $i$, let us say $i^{\prime }$, is determined using the information from the $x$-th entry of the table by computing $i^{\prime }=q_x+(i-p_x)$. Then, a linear search is conducted on the table to find the input interval containing $i^{\prime }$, let us say $x^{\prime }$, by starting at the input interval that contains the LF-mapping of the start of the $x$-th run, which is stored by the table. The final step is to return $(i^{\prime },x^{\prime })$ for use in future move queries that would be necessary in any other backward search. The second step seems like it could be linear in the number of input intervals that need to be searched, but Nishimoto and Tabei established that the number of intervals searched will be constant because the move structure is constructed on a balanced interval sequence.

In their paper, RLBWT Tricks, Brown et al. [7] implement the move structure defined by Nishimoto and Tabei, which had yet to be implemented, and refine their approach for practical applications.

First, Brown et al. modified the move table introduced by Nishimoto and Tabei by replacing absolute indices with run-relative start indices. As a result, the table representation requires one fewer column of integers. This change preserves the ability to compute the LF-mapping and $\varphi$ efficiently, but requires an additional sparse bitvector over BWT to mark run starts. This bitvector is used to convert run/offset pairs into absolute text positions at the end of the backward search. Each row in the table is composed of: the character of the run, the length of the run, the index of the run, the LF-mapping of the position at the start of the run, the index of the run that contains the LF-mapping, and the offset into that run it occurs. The main modification that Brown et al. do to the move data structure is to compress the table. Several different methods are presented, implemented, and compared to achieve this. Primarily, the table is partitioned into blocks of fixed size $B$, allowing each block to be easily loaded into the cache and taking advantage of the inherent locality of the tables for speed-ups. The fixed sizing specifically has the added benefit of positions within the blocks being mapped to easily. Additionally, by using uncompressed bitvectors as a wavelet tree to store the run characters in the block rank, and select queries can be done fast. The authors also augment each block with the position of the first run before and after the current block for each character to allow for quick transitions into other blocks while scanning.

To compress the run-length and offset columns, the authors use directly addressable codes (DACs), which provide a hierarchical encoding of integers [28]. Three different methods are presented to compress the LF-mapping columns. For two of these methods, the property that the LF-mapping values for multiple runs of the same character in a block form a non-decreasing subsequence is exploited. This enables the creation of sublists for each character, which can be represented using partial differences. One method stores these partial differences using DACs and samples the list of differences. The second method represents the partial differences as an uncompressed bitvector, where the number of 0’s before each 1 encodes a difference value. For example, to retrieve the fifth value in the difference list, a select query is performed on the bitvector to locate the fifth 1. The third method omits partial differences entirely and instead samples the list of differences while applying linear interpolation to recover missing values. Rather than using a single bitvector to represent run head indices, an approach more suitable for large alphabets, a separate uncompressed bitvector is maintained for each character. These per-character bitvectors support faster rank and select queries compared to wavelet trees, at the cost of increased space per character, which is acceptable for data such as DNA that has a small alphabet. Without having implemented the second half of the Nishimoto and Tabei paper, their implementations rely on sequential scans of the table, which takes $\mathrm{\Omega }(r)$ time in the worst case. To combat against large sequential scans, the authors additionally implemented a way to reduce the total number of runs by splitting large runs by some maximum run length during construction.

The implementation developed by Brown et al. is designed to efficiently support count queries. Multiple variations of their approach are evaluated, both internally and in comparison to established alternatives such as r-index and RLCSA. The evaluation is carried out on datasets of increasing size, specifically Chromosome 19 and SARS-CoV-2, using randomly selected substrings of varying lengths as query input. Among the tested methods, the variant that employs a bitvector to encode the partial differences list of the LF-mapping values – referred to as lookup-bv – demonstrates superior performance relative to all other approaches. For shorter query lengths, the standard lookup-bv implementation exhibits the best performance, while for longer queries, a variant that integrates lookup-bv with run-splitting techniques proves to be the most effective. In all data sets, the best performing variant of Brown et al. implementation achieves execution speeds comparable to or exceeding those of r-index, although with higher space consumption. Although it is slower than RLCSA, it compensates for it by requiring less space. Overall, the approach proposed by Brown et al. offers a balanced trade-off between time and space efficiency, positioning it between the two competing methods.

Movi, developed by Zakeri et al. [40], is the first full-text pangenome index to utilize Nishimoto and Tabei’s move data structure and builds on the work of Brown et al. as well. Movi takes advantage of the locality of reference of the move data structure to reduce the frequency of cache misses and employs a parallelization strategy to hide the latency of any cache misses that do occur. During LF-mapping steps, the move table stores the index of the run that contains the LF-mapping of the current run head. This other run is usually in a different block of the table, and, as such, accessing it results in a cache miss as the new block is retrieved. To hide this, Zakeri et al. use a memory prefetch command and work on multiple count queries at the same time. In this way, while one of the count queries is waiting for a block to be loaded into the cache, work is still being done on another count query, and the latency is not felt. Hiding this latency is crucial as the jump to a run in another block can occur every iteration of the backward search.

During the backward search, we would normally want to use the LF-mapping of the current BWT range indices to find the next BWT range indices. In the move table, we instead use the stored value that indicates the index of the run that contains the LF-mapping of the current run head. However, the BWT index in the current backward search step does not always correspond to a run head, so the run that contains the LF-mapping of the current run head may not contain the LF-mapping of the current BWT index. In order to find the appropriate run, a sequential search through the runs in the table is necessary. Movi provides two different methods for handling the search. The first is a sequential search through the table that checks runs of characters that match the character of the LF-mapping we are looking for to try and find the run that contains this LF-mapping. The second is a constant-time method in which pointers are stored for each character that point to the indices of runs of the character. The locality of reference allows for the sequential search to be practically inexpensive, especially when compared to a cache failure, but the constant-time method is the faster of the two, as the sequential search still has an unpredictable number of memory accesses. However, the constant-time method does have the added cost of being needed to store the pointers.

The authors compare the backward search capabilities of Movi’s two methods against the r-index using 10 million simulated 150 bp reads of the seven bacterial species in the Zymo community. The results show that both Movi methods were around 16 times faster than the r-index in completing the query, but required three times as much space in the default method (using sequential search) and 11 times as much space in the constant method (storing pointers).

Move-r is a static text index created by Bertram et al. [4] that incorporates the move data structure described in 5.1 with the r-index described in Subsection 3.2. This allows the Move-r index to be practically faster than other r-index implementations with a small increase in the space requirement. This improvement comes from the ability of the Move-r index to use the move data structure to compute LF and ${\varphi }^{-1}$ in $𝖮(1)$ through move queries.

One of the primary contributions is the inclusion of an algorithm to balance a disjoint interval sequence on which the move data structure is built. Having a balanced disjoint interval sequence is essential to the proof that allows move queries to be answered in constant time. The proposed algorithm works by using two balanced search trees (BST), one for the input disjoint intervals and one for output disjoint intervals. It iterates through every output interval in order and determines whether it is an $a$-heavy interval by checking if it contains $\ge 2a$ input intervals, where $a$ is some balancing parameter. If it does, then the output interval and related input interval are split, producing an $a$-balanced interval and another interval that is not necessarily balanced. The algorithm repeats this process until there are only $a$-balanced intervals in the binary search trees. Each step takes $𝖮(\log k)$ time for BST updates, where $k$ is the number of pairs of intervals in the disjoint interval sequence, and $𝖮(a)$ time to determine if an interval is $a$-heavy. The authors show that there are at most $\frac{k}{a-1}$ splits, so their algorithm takes $𝖮(k+\frac{k}{a}\log k)$ time and $𝖮(1)$ additional space. The algorithm of Nishimoto and Tabei for constructing the balanced disjoint interval sequence takes $𝖮(k\log k)$ time and $𝖮(k)$ space. The authors also provide a practically optimized locate algorithm with an asymptotically worse runtime compared to the locate algorithm provided by Nishimoto and Tabei [30]. Their algorithm saves up to $3m$ cache misses, where $m$ is the length of the query pattern, by precomputing and storing the nearest SA value for each run head, which allows constant time lookup during the query rather than accessing multiple different data structures in each iteration of the backward search as is done by OptBWTR.

Move-r is compared against the r-index[15], online-rlbwt [2], rcomp [29], lookup-bv [7], and block-rlbwt[11] using Pizza&Chili [33], concatenated samples of chr19, concatenated samples of the SARS-CoV-2 genome, and other repetitive texts. During construction of the index, block-rlbwt tended to take the least amount of time and the least amount of peak memory consumption, but move-r was 0.9-2 times as fast as the other methods, while requiring 0.33-1 times as much peak memory usage. For count queries, block-rlbwt offered the best trade-off between query throughput and index size, and move-r consistently had the highest query throughput but was the third largest index. For locate queries, move-r was 0.8-2.5 times as large and 2-35 times as fast.