The Burrows-Wheeler Transform of an Elastic-Degenerate String

Abstract

Degenerate strings (DS) and elastic degenerate strings (EDS) are a way to represent, in a compact form, strings that contain a high degree of similarity. They can be particularly useful in some fields, such as text processing or the study of DNA mutations in computational biology, where it is necessary to eficiently manage several variations of a sequence. In practice, a degenerate string is a string whose symbols, called degenerate, can have several alternatives (hence a degenerate symbol is a set). In the literature difer ent constraints have been imposed on degenerate string symbols. For example, the symbol can only be i) a set of letters of the alphabet, ii) a set of strings of the same length, or iii) a set of strings of variable length (including the empty string). We consider the latter in its most general form, which is known as elastic degenerate strings. Our contribution is the introduction of the Burrows-Wheeler transform of an elastic-degenerate string (EDS-BWT). We show that EDS-BWT is reversible and that it can be used to solve the pattern matching problem, i.e., the problem of finding a standard string pattern within an EDS, by adapting the inner properties of the classical Burrows-Wheeler transform. Finally, we implemented the EDS-BWT encoding/decoding and the prototype edsBWTSearch to experimentally compare our pattern matching approach to other existing tools managing elastic degenerate strings.

Keywords

• eol>Burrows-Wheeler Transform
• elastic degenerate strings
• backward search
• pattern matching problem

1. Introduction

The Burrows–Wheeler transform (BWT) was introduced in [ 1 ] as a method for compressing a single string, and then, it was shown to be efe ctive in many other areas, with applications spanning beyond its original purpose [ 2 ]. For instance, it has been successfully used for compact text indexing [ 3, 4, 5, 6 ] and for bioinformatics applications, e.g., for sequence alignment [ 7 ], phylogenetic analysis [ 8 ], genome assembly [ 9 ] as well as for sequencing data compression [ 10 ].

Roughly, the BWT performs a permutation of the letters of an input string . First, the cyclic rotations of are sorted in lexicographic order, and then the bwt( ) is obtained by concatenating the last letters of the (sorted) cyclic rotations of . The BWT can also be define d by sorting the sufixes of $ [ 1 ], where $ is an end-marker symbol that does not appear in itself and it is lexicographically smaller than any of the symbols in .

Shifting the focus from a string to a collection of strings, the BWT of a string collection can be defined analogously. That is, either by sorting the cyclic rotations of all strings in the collection 1, as in [ 11 ], or by sorting their sufixes, as in [ 12 ]. In the latter case, a distinct end-marker symbol is appended to each string, giving an ordered collection of strings. The BWT and its extension to a collection of strings allow to define a map from symbols occurring in the transformed string to the string of lexicographically sorted symbols of . This mapping is known as LF-mapping, and it allows both the reversibility of the transform and the search for patterns (called backward search). The LF-mapping and the backward search are the key machinery in the FM-index [ 13 ].

Degenerate strings (DS) and elastic degenerate strings (EDS) have been introduced as a way to eficiently represent sequences that contain a high degree of similarity. For instance, in bioinformatics, they are used to represent pangenomes, which are collections of closely-related genomic sequences that one needs to analyze together [ 14 ].

A degenerate string (also known as indeterminate string) over an alphabet Σ is a sequence = 1 · · · where each is a subset of Σ. It represents any string that can be obtained by selecting one letter from each subset from left to right. In 2017, Iliopoulos et al. [ 15 ] defined a more general notion of degenerate strings: an elastic-degenerate string (EDS) over Σ is a of non-empty subsets of strings over Σ, where each is called sequence = 1 · · · degenerate string [ 16 ]. degenerate symbol. If instead each contains strings of the same length, is called general Here is an example of an EDS, which we use throughout this paper:

{︃ =

ATTGCT }︃{︃ }︃{︃

CTACGGACT }︃{︃ A }︃{︃

CTGT }︃ (1) We remark that Eq. (1) is a compact way to represent all the strings that are obtained by taking an element from each set and concatenating them in order.

Elastic-degenerate strings and their variants have been much studied in the literature in recent years, mainly for the pattern matching problem, which consists of finding the occurrences of a pattern in an ED string [ 17, 16, 18, 19, 20 ]. For instance, the pattern occurs in , across the first three degenerate symbols. A variety of methods and data structures have been used for the pattern matching problem on very similar strings. The following partial list gives a few examples [ 21, 22, 23, 24, 25, 26, 27 ] (see also references therein). 1.1. Our contribution In this paper, we introduce the Burrows-Wheeler transform of an elastic degenerate string (EDS-BWT), which applies the EBWT to a sequence of ordered collections of strings of any length, and show that this transform is reversible. The core of our method is a mix between the classical BWT [ 1 ] and the EBWT [ 12 ], together with a mapping between strings belonging to consecutive degenerate symbols. We also present an algorithm for solving the pattern matching problem on an elastic degenerate string. Specifically, we are able to return the number of starting positions of pattern occurrences and, for each pattern occurrence, the index of the degenerate symbol and the string position at which the occurrence starts. For instance, searching pattern

1In this case, one needs to use the -order defined in [11].

= in the string in Eq. (1), we return that there is at least one occurrence that starts in the first degenerate symbol at the last letter of the only string in it.

To the best of our knowledge, no other work does it, although some authors introduced the BWT for degenerate strings [ 27 ] and for closely-related sequences [ 28, 29, 5, 30, 31, 32 ] (see also references therein).

2. Background and Notation

Let Σ = {1, 2, . . . , } be a finite ordered alphabet Σ with 1 < 2 < . . . < , where < denotes the standard lexicographic order, and let be the empty string. Let be a string (or word) of length on Σ and [] its -th letter (or symbol). A substring [, ] of coincides with [] · [ + 1] · · · [], where · is the concatenation operator. For any 1 ≤ ≤ , the substring [1, ] is called a prefix of and [, ] a sufix of .

The Burrows-Wheeler Transform (BWT) [ 1 ]. The BWT is a well-known reversible transformation defined on a string that permutes its letters. By appending an end-marker symbol $ to and by sorting all the sufixes of $ in lexicographic order, the output of the BWT is a string bwt() of length + 1 obtained by concatenating the letters (circularly) preceding each sufix in the list of sorted sufixes. More precisely, for each , bwt()[] is the letter preceding the -th lexicographically smallest sufix of the string $, except for the sufix $, where the preceding letter is set to be $.

The BWT of a string collection [ 11, 12 ]. The BWT extended to a string collection , known as EBWT, is a reversible transformation that produces a string ebwt() that is a permutation of all the symbols of all strings in .

Let = {1, 2, . . . , } be a collection of strings on the alphabet Σ. We append to each string ∈ a diferent end-marker symbol $, not belonging to Σ and lexicographically smaller than any other symbol in Σ, by setting $ < $ for each < 2. That is, if has length , we define [ + 1] = $. In the following, we will always use to refer to a string of length + 1 terminating with the end-marker symbol $. Finally, let = ∑︀ =1 + be the number of letters of all strings in (including their end-marker symbol).

Note that, after appending a distinct end-marker symbol to each string in , the collection becomes an ordered collection. Exclusively for implementation purposes, a unique end-marker symbol $ is used for all strings in , even if each $ implicitly carries the index of the string to which it was appended.

LF-mapping. Let be the string bwt() and the string obtained by sorting all the symbols of lexicographically. The reversibility of the BWT (as well as of the EBWT) is based on the following two properties stated in [ 1 ]: • For all = 1, . . . , + 1, the symbol [] circularly follows the symbol [] in the string ; • For each symbol ∈ Σ, the -th occurrence of in corresponds to the -th occurrence of in . 2We note that in [ 11 ] the EBWT is defined without appending end-marker symbols to the strings: cyclic rotations of the strings in are sorted by means of an order relation, called -order, on infinite strings. From the second property, it follows that given a position in such that [] = , the position in corresponding to that occurrence of is given by [[]] + rank(, []), where is an array storing for any ∈ Σ the total number of symbols in that are smaller than , and rank(, []) is the number of occurrences of [] in the prefix [1, ]. This mapping gives a correspondence between symbol occurrences in and symbol occurrences in , and is known as LF-mapping [ 3, 13 ].

Bit vectors. A string of zeros and ones is called bitvector. Given a bitvector of length , the operators rank and select are defined as follows, for any 1 ≤ ≤ and ∈ {0, 1}: rank(, ) = |{ | 1 ≤ ≤ and [] = }| select(, ) = , with [] = and rank(, ) = , if such exists.

In other words, rank(, ) gives the number of occurrences of the bit in the prefix [1, ], while select(, ) returns the index of the -th occurrence of in (if it exists). A bitvector can be preprocessed in order to support rank and select queries in constant time [ 33 ].

Elastic Degenerate string. An elastic degenerate string (see [ 15, 17 ]) (or ED string) is a sequence of finite nonempty sets of strings (including ) 1 · · · of combined total length . Each is called degenerate symbol.

Given an elastic degenerate string = 1 · · · , for each = 1, . . . , , we denote the strings in by ,1, . . . , ,ℓ , where || = ℓ ≥ 1.

Essentially, an elastic degenerate string represents all possible strings that can be constructed by taking an element from each degenerate symbol and concatenating them. For example, {, }{, }{ } represent the strings , , , .

3. The BWT of an elastic degenerate string

In this section we define EDS-BWT, which is the BWT of an elastic degenerate string = 1 · · · , where = {,1, . . . , ,ℓ }, || = ℓ ≥ 1 and ℓ = ∑︀=1 ℓ. We need to append to each , an end-marker symbol $ ∈/ Σ, with = + ∑︀−=11 ℓ. In this way, denoting , $ by , we obtain a single ordered collection (of strings) = {1, 2, . . . , ℓ1 , ℓ1+1, . . . , ℓ1+ℓ2 , . . . , ℓ}. Roughly, we are appending an end-marker symbol at the end of each string of each degenerate symbol, proceeding left-to-right among the symbols and giving an arbitrary order to the strings within the degenerate symbols. Note that when , = , the resulting string is $, with as above.

Definition 3.1. Given an elastic degenerate string = 1 · · · the Burrows-Wheeler transform of (called EDS-BWT) is defined as the pair edsbwt() = (ebwt(), ℒ) where is the ordered collection = {1, . . . , ℓ} and ℒ is defined as follows: For = 1, . . . , ℓ, let be the position of $ in ebwt(), and let the associated string belong to the degenerate symbol , for some 1 ≤ ≤ .

Then, if > 1, we define ℒ() to be the interval of positions [, ] in ebwt() such that − 1 = {, . . . , }. Otherwise, ℒ() circularly gives the interval of positions [, ] such that = {, . . . , }.

The ED string of our running example Eq. (1) gives the ordered collection each > 1. Also, the next remark follows:

= { $1, $2, $3, $4, $5, $6, $7, $8}. together with its index can be returned.

Remark 3.2. As in [ 12 ] for the computation of the EBWT, we can use the same end-marker for all strings, so that the size of the alphabet increases by just one. However, each end-marker has a diferent (implicit) index determined by the order of the strings in the collection . In this way, during the construction of ebwt(), a list containing the position in ebwt() of each $

Note that, in the definition of

ℒ, = (ℓ1 + · · · + ℓ− 2) + 1 and = ℓ1 + · · · + ℓ− 1 for Remark 3.3. Let and ′ be two diferent positions in ebwt() such that ebwt[] = $ and ebwt[′] = $′ . If the associated strings and ′ belong to the same degenerate symbol , then ℒ() = ℒ(′) by definition. row F 1 $1 2 $2 3 $3 4 $4 5 $5 6 $6 7 $7 8 $8 9 A 10 A 11 A 12 A 13 A 14 A 15 A 16 C 17 C 18 C 19 C 20 C 21 C 22 G 23 G 24 G 25 G 26 T 27 T 28 T 29 T 30 T 31 T 32 T 33 T 34 T T A A A T A $7 T T T $4 $6 T G $1 A G A $2 $5 $8 G T C T C C G C $3 C T C A Reversibility To show that the EDS-BWT is reversible, we describe how LF-mapping works for this new transform. Note that LF-mapping for the EDS-BWT is diferent from that for the classical EBWT. Indeed, in the EBWT, the strings are independent from each other. Conversely, in our EDS-BWT, the strings in belonging to consecutive degenerate symbols need to be “linked” to reconstruct the elastic degenerate string .

Table 1 shows the transform of our running example. For instance, edsbwt()[ 20 ] = $5, and its associated string 5 belongs to 3, so ℒ(20) = [ 2, 4 ], which is an interval of three positions. Thus the previous degenerate symbol 2 is a set of three strings. Specifically, it is 2 = { , , }.

The following proposition states the properties of the LF-mapping for an ED string. The proof follows immediately from the original LF-mapping properties and the definition of ℒ. Proposition 3.1. Let be the collection of strings associated to an ED string and let edsbwt() = (ebwt(), ℒ). Let = ebwt() and be the sequence of the lexicographically sorted letters of . The following conditions hold: • For all = 1, . . . , , the letter [] is circularly followed by the letter [] in its associated string ; • For each letter ∈ Σ, its occurrences in appear in the same order as in , i.e. the -th occurrence of in corresponds to the -th occurrence of in . • For all = 1, . . . , such that [] = $ for some , the letter [] is preceded (in the ED string) by all the symbols in [, ], where ℒ() = [, ]. If any end-marker symbol appears in [, ], then [] is preceded (in the ED string) by the empty string.

Proposition 3.1 allows us to define the following permutation that maps each position of to and allows us to link the first symbol of a string in a degenerate symbol to the last symbol of any string within the previous degenerate symbol.

[] = {︃[[]] + (, []) if [] ̸= $ ℒ() if [] = $ (2) The of our running example is shown in Table 2.

By using this modified LF-mapping, we are able to reconstruct the ED string, as follows. We begin from the first end-marker symbol, which is $1, in position 15. Using ℒ we can ifnd the positions [ 8, 8 ] of the last letters in the final degenerate symbol, which in this case is [ 8 ] = . Note that, in this example, we have a single string in both the first and last degenerate symbols. In the general case, by Definition 3.3, ℒ applied to any end-marker symbol in 1 gives the interval containing the positions of the last letters of all strings in . Then, using the LF-mapping described in Eq. (2), we find that is preceded by [ 28 ] = , and continue in this way until we have reconstructed the string 8 = . This is correct since 5 = { }. This step stops at position 21, which is such that [ 21 ] = $8. Thus, we use ℒ to obtain the positions corresponding to the last letters in the previous degenerate symbol 4, ℒ(21) = [ 6, 7 ]. Since we have two positions now, we reconstruct two diferent strings using the same strategy as before; the first is , ending at position 12, while the second terminates immediately, since (7) = $7, indicating that 7 is the empty string. This is correct since 4 = {, }.

We do not need to compute both ℒ(7) and ℒ(12), since by Definition 3.3 they give the same result of [ 5, 5 ]. The reconstruction continues by alternating between the parallel reconstruction of all the strings in a degenerate symbol and the linking to the strings in the previous degenerate symbol. This process ends when it reaches the position of first end-marker, since that was the starting position.

Implementation of ℒ using bitvectors. Let = 1 · · · be an ED string. We define the bitvector () associated to as the concatenation of (1), . . . , (), where the bitvector () of length || has all zeroes except for its first bit, which is one. For instance, the bitvector of our running example is () = 1 100 1 10 1. An analogous definition of a bitvector to represent the underlying structure of a degenerate string appears in [ 34 ].

The following proposition shows how, using rank and select on (), we can compute ℒ() in constant time, provided that the index such that [] = $ is given. Proposition 3.2. Let = 1 · · · be an elastic degenerate string, and let = 1, . . . , ℓ the ordered collection of strings contained in the degenerate symbols and let 1 ≤ ≤ ℓ. For > 1, let , be the indexes such that, if ∈ , then − 1 = {, . . . , }. If = 1, let , such that = {, . . . , }. Then and can be computed in constant time using just and (). Specifically, for > 1, thus rank()(, 1) > 1, then = select()(rank()(, 1) − 1, 1), = select()(rank()(, 1), 1) − 1.

For = 1, thus rank()(, 1) = 1, then

= select()(rank()(ℓ, 1), 1), = ℓ.

Proof. Since () is the concatenation of the ()’s, and since each () contains exactly one 1, then rank()(, 1) gives the such that ∈ , for each . Now observe that, for > 1, and are the indexes of the first and last string of − 1, so ()[] = 1 and ()[ + 1] = 1. In particular, ()[] is the ( − 1)-th 1, while is the position preceding the -th 1. On the other hand, for = 1, and are the indexes of the first and last string of .

We can now compute and using select. Finally, rank and select can be performed in constant time on a bitvector, for example using the sdsl-lite library [ 33 ].

Note that we need the index for computing ℒ(). The mapping from to can be obtained during the construction of the ebwt() (see Definition 3.2). In our implementation, to compute in constant time, we use a bitvector such that [] = 1 if and only if ebwt()[] = $ for some 1 ≤ ≤ ℓ, and an associated array of length ℓ such that [] = , if rank(, 1) = , for each 1 ≤ ≤ ℓ. Hence = [rank(, 1)].

We conclude this section by encapsulating the previous observations in Algorithm 1. For simplicity, when describing the algorithms in the following, we just write $ ← [] to mean that we obtain the index from position , if such exists, and assign 0 to otherwise. Algorithm 1: link(, ()) ← ← 1 $ ← []; 2 if = 0 then 3 return ∅ 4 if $ ∈/ 1 then 5 ← select()(rank()(, 1) − 1, 1); 6 ← select()(rank()(, 1), 1) − 1; 7 else 8 select()(rank()(ℓ, 1), 1); 9 (()); 10 return [, ];

4. Exact pattern matching problem

In this section we show how the EDS-BWT can be used to resolve the exact pattern matching problem for an elastic degenerate string = 1 · · · .

In the classical problem, a string ∈ Σ* contains the pattern ∈ Σ* if and only if is a substring of . This concept is generalized to (elastic) degenerate strings in the natural way. Definition 4.1. Let ∈ Σ* , and let = 1 · · · be an ED string over Σ, with = {,1, . . . , ,ℓ } for each . We say that contains the pattern if there exists a collection of strings , ∈ , . . . , +,+ ∈ +, such that is a pattern of , · · · +,+ , for some , such that 0 ≤ ≤ − 1 and 1 ≤ ≤ − .

In this case, we say that contains an occurrence of at position (, , ) if begins at position in , .

Note that more than one occurrence of can start (and end) at the same starting (and ending) position, since they can select diferent strings in the degenerate symbols. For example the pattern in = {, , C}{, , }{, , }{T} can be obtained in two diferent ways starting at position (1, 3, 5) and ending at position (4, 1, 1) (in bold), taking either the red or the blue strings in the degenerate symbols 2 and 3. Conversely, there are no occurrences of , because and belong to the same degenerate symbol.

In [ 3 ], the authors showed how to search a pattern = [1, ] backwards by the output of the classical BWT. A backward search algorithm first searches for the []-interval (i.e. the interval in of the symbols associated to sufixes starting with []), then for the ( [ − 1] [])-interval (i.e. the interval in of the symbols associated to sufixes starting with [ − 1] []), and so on, until the whole pattern [1, ] is found, if there is one. Specifically, given an interval [, ] such that [, ] are all the letters followed by [, ] in the original string, then the letters followed by [ − 1, ] are in the interval [′, ′], with ′ = [] + ( − 1, ) + 1, ′ = [] + (, ).

We adapt the backward search algorithm defined in [ 3 ] in order to apply it to edsbwt() and solve the pattern matching problem. Algorithm 2: edsBWTSearch(edsbwt(), , (), [1, ]) 1 = [], = − 1; 2 = [] + 1, = [ + 1]; 3 = {[, ]} ; 4 while ̸= ∅ and > 0 do

// Compute ℒ on every interval 5 for [, ] ∈ do 6 for ∈ [, ] do 7 ← (, link(, ())) ; = − 1 ; ← ′; 17 if ̸= ∅ then 18 return "Found"; 19 else 20

return "Not found";

Algorithm 2 (edsBWTSearch) outlines the steps of the algorithm. edsBWTSearch loops over the symbols of , starting from the last, and for each iteration updates a list Intervals of the intervals of which have a positive match for the current symbol. The update happens in two steps: • the first step calls link (Algorithm 1). The algorithm finds all the end-marker symbols contained in the current intervals, and compute ℒ for each of those positions, adding new intervals to the list. We note the following: – the newly added intervals are checked again, since they could also contain an end-marker symbol, meaning that its corresponding string is the empty string; – the algorithm actually uses a non-circular version of link, because the pattern matching problem is not circular. Thus, when an end-marker symbol belonging to the first degenerate symbol is considered, the output of link is the empty interval (diferently from lines 8-9 in Algorithm 1); – the Merge called in line 7 (and again in line 14) is a modified version of the union operation. It adds the interval obtained from link to the list of intervals, but ifrst checks if is consecutive to or included in other intervals in Intervals. If this happens, a new interval is instead created by combining the two into one; • the second step applies the modified backward search to each interval in the list Intervals.

Note that Merge is applied also in this step.

Since link and the -mapping are called for each interval, then the Merge operation allows for a faster search by reducing the number of intervals (see Section 5 for an upper bound on the number of intervals at each iteration).

Table 3 shows the backward search of pattern on our running example. The search begins from interval [ 16, 21 ], which is the interval of positions corresponding to letter [ 3 ] = (marked in blue in the table). For readability, we will not color every interval in the table at every step, but we will just mark some “branchings” that visually show the link and update parts. Nonetheless, we fully describe the search in the following.

Since there are three end-marker symbols in [ 16, 21 ], link is recursively called on each of them, giving four additional intervals [ 1, 1 ], [ 2, 4 ] (in orange in the table), [ 5, 5 ], [ 6, 7 ]. Note that [ 5, 5 ] is obtained from the recursive call on [ 6, 7 ], since [ 7 ] = $7, which corresponds to the empty string in 4 of Eq. (1). The intervals are thus merged into one: [ 1, 7 ]. Thus = [[ 1, 7 ], [ 16, 21 ]].

The next step is to update each interval, searching for [ 2 ] = . Interval [ 1, 7 ] gives [ 9, 12 ] and interval [ 16, 21 ] gives [ 13, 14 ]. Thus = [[ 9, 14 ]].

Again, each end-marker symbol in the intervals is checked, giving [ 1, 1 ] (in green) and [ 5, 5 ]. Finally, [ 1 ] = is searched, obtaining five occurrences of , starting at positions 1, 5, 9, 10 and 13 of .

The positions in the original ED string can be recovered using the LF-mapping on each of the resulting positions until reaching an end-marker symbol. This gives the string index, the degenerate symbol index and the position in the string, which are (1, 1, 6), (2, 1, 2), (2, 2, 1), (3, 1, 2), (3, 1, 9).

Experiments. In order to evaluate our backward search applied to edsbwt(), we implemented3 a prototype in C++ that builds ℒ and then solves the pattern matching problem on an EDS. Since the EBWT is a well-known structure in string algorithms, we used existing tools to build the ebwt(). Note that the eficient construction of the EBWT has been the subject of extensive research (see for instance [ 12, 35, 36, 37, 38, 39 ]), that is beyond the goal of this paper. Therefore, the time needed to build the ebwt() depends on the tool used, and the resources available. Moreover, since our pattern matching strategy is of-line (it builds edsbwt()), the cost to pre-process , which however took a couple of seconds for the dataset with the longest EDS in our experiments, can be amortized for all the pattern searches.

In order to show the efectiveness of our pattern matching method, we have considered two existing tools, edsm [ 17 ] and eds_search [ 20 ], that solve the elastic-degenerate string matching problem in an on-line manner by taking as input an ED string on a biologic alphabet. In particular, edsm takes an ED string and a solid pattern , it builds the sufix tree of and scan the ED string left-to-right to return the indices of the degenerate symbols at which each occurrence of ends; while eds_search, given as input an ED string and a solid pattern , performs an on-line backward pattern matching by combining the traditional algorithms BNDM and SHIFT-AND, and returns the number of positions in which an occurrence of starts.

Note that there exist other tools for the pattern matching on closely-related sequences, but they encode the similar strings with data structures diferent from ED strings, see for instance [ 31, 29, 28, 40 ].

The experiments were conducted on a DELL PowerEdge R630 machine, 24-core, with Intel(R) Xeon(R) CPU E5-2620 v3 at 2.40 GHz, with 128 GB of shared memory. The system is Ubuntu 14.04.2 LTS.

We tested our tool on the synthetic elastic degenerate strings used in [ 17 ], which are 5 randomly generated strings of lengths 100000, 200000, 400000, 800000 and 1600000. We searched 40 patterns of lengths 8, 16, 32 and 644 (ten patterns each). The patterns were obtained by selecting random parts of the ED string and extracting a substring of the desired length, so that each pattern would occur at least once. This is important, since the absence of a pattern stops prematurely the execution of the algorithm, giving distorted timings.

Table 4 shows that the performance for all datasets of edsBWTSearch is comparable to the one of edsm. The tool eds_search is always the fastest, but it does not return the positions of the pattern occurrences. Excluding the pre-processing time needed to build edsbwt(), edsBWTSearch is always faster than edsm; however, the pre-processing time for the largest dataset is 4.72 seconds that, if added to the edsBWTSearch time of 2.69 seconds, gives 7.41 seconds. Finally, we note that edsBWTSearch is implemented in semi-external memory and stores on disk part of the index edsbwt(). In this way, edsBWTSearch is the tool that uses the least RAM on the largest dataset.

5. Conclusion, discussion and further work

In this work, we introduced a new transform EDS-BWT inspired by the BWT of a string and the EBWT of a string collection. The EDS-BWT permutes the letters of the strings in associated with an ED string by exploiting the () and also returns a function ℒ that allows to link the strings of a degenerate symbol to the strings of the previous one. The introduced transform works over any alphabet and it does not concatenate strings in , but keeps track of links between strings of consecutive degenerate symbols. Moreover, it allows, like the BWT on a string, to build an index on which to perform pattern searches.

We observe that the time and space usage for computing EDS-BWT depends mainly on the construction of ebwt(), because the construction of ℒ can be achieved by simply reading and producing the bit vector (). Since there are several tools for building ebwt(), one can choose the implementation that best suits their own resources.

The pattern matching algorithm involves iterations, where is the length of the pattern . It is easy to see that the first iteration may obtain at most + 1 intervals, because, as for the classical backward search, edsBWTSearch produces a single interval for the LF-mapping of letter [], while ℒ can return at most one interval for each degenerate symbol (see Definition 3.3). In the same fashion, at most intervals can be added at each subsequent

4eds_search does not support patterns of length greater than 32.

¯ 100000 200000 400000 800000 1600000 iteration, one for each degenerate symbol of the EDS. Therefore the worst case is to add intervals at each iteration, for a total of + 1. However, this is a theoretical upper bound that do not take into account the interval merging. Indeed, in the above case, the intervals of each iteration would be merged into one, giving just intervals, which is far from the worst case. Accounting for merging, the worst case occurs when all the intervals are non consecutive. Thus, the maximum amount of intervals added at each iteration is ⌈︀ 2 ⌉︀ , one every two degenerate symbols. Since we have iterations, this sums up to ⌈︀ ⌉︀ + 1 total maximum intervals. 2

During our experiments, we observed the number of intervals to be strictly decreasing after the first iteration. We believe this not to be a coincidence and plan to investigate the matter in subsequent work. Finally, we note that it is possible to reduce the number of intervals by building a partial index for patterns of small lengths.

Another future direction of work is to show that, like the EBWT, the EDS-BWT is dynamic, meaning that we can add/remove a string to a degenerate symbol without needing to rebuild the entire transform, but suitably adding/removing the symbols and updating ℒ.

Moreover, we aim to show that EDS-BWT also allows us to search for multiple patterns at the same time.

Acknowledgments

Work partially supported by the MUR PRIN 2022YRB97K PINC, by INdAM - GNCS Project, codice CUP_E53C23001670001, “Compressione, indicizzazione, analisi e confronto di dati biologici”, and by PNRR - M4C2 - Investimento 1.5, Ecosistema dell’Innovazione ECS00000017 - “THE Tuscany Health Ecosystem” - Spoke 6 “Precision medicine & personalized healthcare”, funded by the European Commission under the NextGeneration EU programme.

References

[1] M. Burrows , D. J. Wheeler , A Block Sorting Lossless Data Compression Algorithm , Technical Report 124 , Digital Equipment Corporation, 1994 .

[2] G. Rosone , M. Sciortino , The Burrows-Wheeler Transform between Data Compression and Combinatorics on Words , in: CiE, volume 7921 LNCS of LNCS , Springer, 2013 , pp. 353 - 364 .

[3] P. Ferragina , G. Manzini, Opportunistic data structures with applications , in: FOCS, IEEE Computer Society, 2000 , pp. 390 - 398 . doi: 10 .1109/SFCS. 2000 . 892127 .

[4] V. Mäkinen , G. Navarro, Succinct sufix arrays based on run-length encoding , Nordic J. of Computing 12 ( 2005 ) 40 - 66 .

[5] V. Mäkinen , G. Navarro, J. Sirén , N. Välimäki , Storage and retrieval of highly repetitive sequence collections , J. Comput. Biol . 17 ( 2010 ) 281 - 308 .

[6] T. Gagie , G. Navarro, N. Prezza , Fully functional sufix trees and optimal text searching in bwt-runs bounded space , J. ACM 67 ( 2020 ). doi: 10 .1145/3375890.

[7] H. Li , R. Durbin , Fast and accurate long-read alignment with Burrows-Wheeler transform , Bioinformatics 26 ( 2010 ) 589 - 595 . doi: 10 .1093/bioinformatics/btp698.

[8] V. Guerrini , A. Conte , R. Grossi , G. Liti, G. Rosone , L. Tattini, phyBWT2: phylogeny reconstruction via eBWT positional clustering , Algorithms Mol. Biol . 18 ( 2023 ) 11 . doi: 10 . 1186/S13015-023-00232-4.

[9] J. T. Simpson , R. Durbin , Eficient construction of an assembly string graph using the FM-index , Bioinform . 26 ( 2010 ) 367 - 373 .

[10] V. Guerrini , F. A. Louza , G. Rosone, Parallel lossy compression for large fastq files , in: Biomedical Engineering Systems and Technologies , Springer Nature Switzerland, Cham, 2023 , pp. 97 - 120 .

[11] S. Mantaci , A. Restivo , G. Rosone, M. Sciortino , An extension of the Burrows-Wheeler Transform , Theor. Comput. Sci . 387 ( 2007 ) 298 - 312 . doi: 10 .1016/j.tcs. 2007 . 07 .014.

[12] M. J. Bauer , A. J. Cox , G. Rosone, Lightweight algorithms for constructing and inverting the BWT of string collections , Theor. Comput. Sci . 483 ( 2013 ) 134 - 148 . Source code: https://github.com/BEETL/BEETL.

[13] P. Ferragina , G. Manzini, An experimental study of a compressed index , Information Sciences 135 ( 2001 ) 13 - 28 . doi: 10 .1016/S0020- 0255 ( 01 ) 00098 - 6 .

[14] T. C. P.-G. Consortium , Computational pan-genomics: status, promises and challenges , Briefings in Bioinformatics 19 ( 2016 ) 118 - 135 . doi: 10 .1093/bib/bbw089.

[15] C. S. Iliopoulos , R. Kundu , S. P. Pissis , Eficient pattern matching in elastic-degenerate texts , in: 11th International Conference on Language and Automata Theory and Applications (LATA) , volume 10168 of Springer LNCS , 2017 , pp. 131 - 142 .

[16] M. Alzamel , L. A. K. Ayad , G. Bernardini, R. Grossi , C. S. Iliopoulos , N. Pisanti , S. P. Pissis , G. Rosone, Comparing degenerate strings, Fundam. Informaticae 175 ( 2020 ) 41 - 58 . doi: 10 .3233/FI-2020- 1947 .

[17] R. Grossi , C. S. Iliopoulos , C. Liu, N. Pisanti , S. P. Pissis , A. Retha , G. Rosone , F. Vayani , L. Versari , et al., On-line pattern matching on similar texts , in: Proceedings of 28th Annual Symposium on Combinatorial Pattern Matching (CPM) , volume 78 , Schloss DagstuhlLeibniz -Zentrum für Informatik GmbH , 2017 , p. 1 .

[18] C. S. Iliopoulos , R. Kundu , S. P. Pissis , Eficient pattern matching in elastic-degenerate strings , Information and Computation 279 ( 2021 ) 104616 . doi: 10 .1016/j.ic. 2020 . 104616 .

[19] G. Bernardini , P. Gawrychowski , N. Pisanti , S. P. Pissis , G. Rosone, Elastic-degenerate string matching via fast matrix multiplication , SIAM Journal on Computing 51 ( 2022 ) 549 - 576 . doi: 10 .1137/20M1368033.

[20] P. Procházka , O. Cvacho , L. Krčál , J. Holub , Backward pattern matching on elasticdegenerate strings , SN Computer Science 4 ( 2023 ) 442 .

[21] H. Soldano , A. Viari , M. Champesme , Searching for flexible repeated patterns using a non-transitive similarity relation , Pattern Recognition Letters 16 ( 1995 ) 233 - 246 .

[22] N. Pisanti , H. Soldano , M. Carpentier, Incremental inference of relational motifs with a degenerate alphabet , in: CPM , volume 3537 of Springer LNCS , 2005 , pp. 229 - 240 .

[23] N. Pisanti , H. Soldano , M. Carpentier , J. Pothier , A relational extension of the notion of motifs: Application to the common 3d protein substructures searching problem , Journal of Computational Biology 16 ( 2009 ) 1635 - 1660 .

[24] K. Abrahamson , Generalized string matching , SIAM Journal of Computing 16 ( 1987 ) 1039 - 1051 .

[25] M. Crochemore , C. S. Iliopoulos , T. Kociumaka , J. Radoszewski , W. Rytter , T. Walen, Covering problems for partial words and for indeterminate strings , Theoretical Computer Science 698 ( 2017 ) 25 - 39 .

[26] C. S. Iliopoulos , J. Radoszewski , Truly Subquadratic-Time Extension Queries and Periodicity Detection in Strings with Uncertainties , in: 27th Annual Symposium on Combinatorial Pattern Matching (CPM) , volume 54 of LIPIcs, 2016 , pp. 8 : 1 - 8 : 12 .

[27] J. W. Daykin , R. Groult , Y. Guesnet , T. Lecroq , A. Lefebvre , M. Léonard , L. Mouchard, É. Prieur, B. W. Watson , Eficient pattern matching in degenerate strings with the Burrows-Wheeler transform , Information Processing Letters 147 ( 2019 ) 82 - 87 . doi: 10 .1016/ j.ipl. 2019 . 03 .003.

[28] J. C. Na , H. Kim , H. Park , T. Lecroq, M. Léonard , L. Mouchard , K. Park, FM-index of alignment: A compressed index for similar strings , Theoretical Computer Science 638 ( 2016 ) 159 - 170 .

[29] S. Maciuca , C. del Ojo Elias , G. McVean , Z. Iqbal , A Natural Encoding of Genetic Variation in a Burrows-Wheeler Transform to Enable Mapping and Genome Inference , in: Algorithms in Bioinformatics, Springer International Publishing, Cham, 2016 , pp. 222 - 233 .

[30] S. Huang , T. W. Lam , W. K. Sung , S. L. Tam , S. M. Yiu , Indexing similar DNA sequences , in: Algorithmic Aspects in Information and Management , Springer Berlin Heidelberg, Berlin, Heidelberg, 2010 , pp. 180 - 190 .

[31] L. Huang , V. Popic , S. Batzoglou , Short read alignment with populations of genomes , Bioinformatics 29 ( 2013 ) i361 - i370 . doi: 10 .1093/bioinformatics/btt215.

[32] T. Büchler , E. Ohlebusch, An improved encoding of genetic variation in a Burrows-Wheeler transform , Bioinformatics 36 ( 2019 ) 1413 - 1419 . doi: 10 .1093/bioinformatics/ btz782.

[33] S. Gog , T. Beller , A. Mofat , M. Petri, From theory to practice: Plug and play with succinct data structures , in: 13th International Symposium on Experimental Algorithms , (SEA 2014 ), 2014 , pp. 326 - 337 . doi: 10 .1007/978-3- 319 -07959-2_ 28 , source code: https: //github.com/simongog/sdsl-lite.

[34] P. Bille , I. L. Gørtz , T. Stordalen, Rank and select on degenerate strings , in: 2024 Data Compression Conference (DCC) , IEEE, 2024 , pp. 283 - 292 .

[35] L. Egidi , F. A. Louza , G. Manzini, G. P. Telles , External memory BWT and LCP computation for sequence collections with applications , Algorithms Mol. Biol . 14 ( 2019 ) 6: 1 - 6 : 15 . doi: 10 .1186/s13015-019-0140-0.

[36] P. Bonizzoni , G. Della Vedova , Y. Pirola , M. Previtali , R. Rizzi , Multithread Multistring Burrows-Wheeler Transform and Longest Common Prefix Array , Journal of computational biology 26 ( 2019 ) 948 - 961 . doi: 10 .1089/cmb. 2018 . 0230 .

[37] F. A. Louza , G. P. Telles , S. Gog , N. Prezza , G. Rosone, gsufsort: constructing sufix arrays, LCP arrays and BWTs for string collections , Algorithms Mol. Biol . 15 ( 2020 ) 18 . doi: 10 .1186/s13015-020-00177-y.

[38] N. Prezza , G. Rosone, Space-eficient construction of compressed sufix trees , Theoretical Computer Science 852 ( 2021 ) 138 - 156 . doi: 10 .1016/j.tcs. 2020 . 11 .024.

[39] D. Díaz-Domínguez , G. Navarro, Eficient Construction of the BWT for Repetitive Text Using String Compression , in: CPM 2022 , volume 223 of LIPIcs, 2022 , pp. 29 : 1 - 29 : 18 .

[40] J. Sirén , Indexing variation graphs , in: Proceedings of the Meeting on Algorithm Engineering and Experiments (ALENEX) , 2017 , pp. 13 - 27 . doi: 10 .1137/1.9781611974768.2.