Maxime Crochemore
email: maxime.crochemore@kcl.ac.uk

Roberto Grossi
email: grossi@di.unipi.it

Juha Kärkkäinen
email: juha.karkkainen@cs.helsinki.fi

Gad M Landau
email: landau@cs.haifa.ac.il

Computing the Burrows-Wheeler Transform in Place and in Small Space

Keywords: Burrows-Wheeler transform, in-place algorithms, string algorithms, suffix sorting

des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

The Burrows-Wheeler Transform [START_REF] Burrows | A block-sorting lossless data compression algorithm[END_REF] (known as BWT) of a text string is at the heart of the bzip2 family of text compressors, and finds also applications in text indexing and sequence processing. Consider an input text string T ≡ T [0 . . n -1] and the set of its suffixes T i ≡ T [i . . n -1] (0 ≤ i < n) under the lexicographic order, where T [n -1] is an endmarker character $ smaller than any other character in T . The alphabet Σ from which the characters in T are drawn can be unbounded.

A classical way to define the BWT uses the n circular shifts of the text T = mississippi$ as shown in the first column of Table 1. We perform a lexicographic sort of these shifts, as shown in the second column: if we mark the last character from each of the circular shifts in this order, we obtain a sequence L of n characters that is called the BWT of T . Its relation with suffix sort is well known, as illustrated in the third column: the rth character in L is T [j -1] if and only if T j is the rth suffix in the sort (except the borderline case j = 0, for which we take T [n -1] as character).

As it can be seen in the example of Table 1, the BWT produces a text of the same length as the input text T . The transform is reversible since it is a one-to-one function when the input text is terminated by an endmarker $. Thus, not only we can recover T from L alone, but typically L is more compressible than T itself using 0th-order compressors [START_REF] Manzini | An analysis of the Burrows-Wheeler transform[END_REF]. There are now efficient methods that convert T to L and vice versa, taking O(n log n) time for unbounded alphabets in the worst case [START_REF] Adjeroh | The Burrows-Wheeler Transform: Data Compression, Suffix Arrays, and Pattern Matching[END_REF]. 1 The BWT is also a key element of some compressed text indexing implementations due to the small amount of space it requires: some examples are the solution by Ferragina and Manzini [START_REF] Ferragina | Indexing compressed text[END_REF] or that by Grossi et al. [START_REF] Grossi | High-order entropycompressed text indexes[END_REF], where the transform is associated with the techniques of wavelet trees and of succinct data structures using rank-select queries on binary sequences [START_REF] Munro | Tables[END_REF]. One of the prominent applications of the BWT is for software dealing with Next Generation Sequencing, where millions of short strings, called reads, are mapped onto a reference genome. Typical and popular software of this type are Bowtie [START_REF] Langmead | Ultrafast and memory-efficient alignment of short DNA sequences to the human genome[END_REF], BWA [START_REF] Li | Fast and accurate long-read alignment with Burrows-Wheeler transform[END_REF] and SOAP2 [START_REF] Lam | High Throughput Short Read Alignment via Bi-directional BWT[END_REF]. Here it is crucial that the genome is indexed in a compact manner to get reasonable running time. Space issues for computing the BWT are thus relevant: frequently the input data is so large that the input text T stays in main memory while any additional data structure of similar size cannot fit in the rest of the main memory [START_REF] Wing-Kai Hon | Breaking a time-and-space barrier in constructing full-text indices[END_REF].

All the previous work for computing the BWT of T relies on the fact that (a) we need first to store the suffix sorting of T (also known as suffix array [START_REF] Manber | Suffix arrays: A new method for on-line string searches[END_REF]), thus occupying n memory cells for storing integers, and (b) we need to output the BWT in another array storing n characters. Motivated by these observations, we want to study the case in which (a) and (b) are avoided, thus saving on the space occupied by them.

In this paper, our goal is to obtain the BWT by directly permuting T and using just O(1) memory cells, i.e., we aim at an in-place algorithm for computing the BWT. We consider the model in which the text T is stored as an array of n entries, where each entry stores exactly one character of T . Note that storing an integer usually takes more space than a character, so we assume that only the characters of T can be kept in the array T . Moreover, T is not read-only but it can be modified at any time, and just O(1) additional memory cells (besides T ) can be kept for storing auxiliary information. 2Note that our model represents some realistic situations in which one has to handle large text collections, or large genomic sequences, without relying on extra memory for (a) and (b). Hence it is crucial to maximize the amount of data that can fit into main memory: not storing explicitly (a) and (b) permits to save space, which is typically regarded as taking more than half of the total space required. For instance, DNA sequences are stored by using 2 bits per character and machine integers take 64 bits. Here we just need 2n bits to store the (genomic) text and save the 64n bits used for storing the intermediate suffix sorting in (a) and the 2n bits for storing the output of BWT in another array in (b): this means that during the BWT construction, we can fit almost 33 times more text using the same main memory size, thus eliminating the usage of the slower external memory for this time-consuming task in these cases.

From the combinatorial point of view, the in-place BWT is an interesting question to solve on strings. There are space saving approaches storing the suffix sorting in compressed form [START_REF] Wing-Kai Hon | A space and time efficient algorithm for constructing compressed suffix arrays[END_REF][START_REF] Chae | Alphabet-independent linear-time construction of compressed suffix arrays using o(nlogn)-bit working space[END_REF][START_REF] Wing-Kai Hon | Breaking a time-and-space barrier in constructing full-text indices[END_REF][START_REF] Okanohara | A linear-time Burrows-Wheeler transform using induced sorting[END_REF] or only partially at a time [START_REF] Kärkkäinen | Fast BWT in small space by blockwise suffix sorting[END_REF], but none of them provides an in-place algorithm. In-place selection and sorting does not seem to help either [START_REF] Dobkin | Optimal time minimal space selection algorithms[END_REF][START_REF] Franceschini | In-Place Suffix Sorting[END_REF][START_REF] Hoare | Algorithm 65: Find[END_REF][START_REF] Vinay | Selection from read-only memory and sorting with minimum data movement[END_REF][START_REF] Raman | Improved Upper Bounds for Time-Space Trade-offs for Selection[END_REF]. It is well known that in-place sorting requires the same comparison cost of Θ(n log n) as in standard sorting. But for the BWT, we only know its comparison cost of Θ(n log n) for the standard construction. As far as we know, no result is known for the in-place construction of BWT: a naive solution is not that simple, even if it results in exponential time. Indeed, any movement of a character T [j] to another position inside T at least changes the content of its suffixes T i for 0 ≤ i ≤ j, making the algorithmic flavor of this problem different from that of in-place sorting n elements.

The above discussion suggests that a careful orchestration of the move-ment of the characters inside T is needed to avoid losing the content of some suffixes before they contribute to the BWT. Our idea is to define a sequence of transformations B 0 , B 1 , . . . , B n , where B 0 is the input text T and B n is the final BWT of T . For 1 ≤ ℓ ≤ n, we have that B ℓ is the BWT of the last ℓ characters in T and is computed from B ℓ-1 (re)using just O(1) extra memory cells. We think that this sequence of transformations could be of independent interest for the community of string algorithms, and some of the combinatorial properties that we use can be found in [START_REF] Belazzougui | Linear time construction of compressed text indices in compact space[END_REF][START_REF] Wing-Kai Hon | Breaking a time-and-space barrier in constructing full-text indices[END_REF][START_REF] Lam | Compressed indexing and local alignment of DNA[END_REF][START_REF] Salson | A four-stage algorithm for updating a Burrows-Wheeler Transform[END_REF].

In this paper we propose an O(n 2 )-time approach that builds the above sequence of transformations using four integer variables and one character variable, taking O(n) time per transformation in the worst case. The resulting in-place algorithm is simple and can be easily encoded in few lines of C code or similar programming languages. However we do not claim any practicality of our solution due to its quadratic worst-case cost. Our contribution is that it could lay out the path towards faster methods for the space-efficient computation of the BWT: any method to compute B k from B k-1 in t(n) time (re)using s(n) space, would lead to a construction of the BWT in O(n • t(n)) time using O(1 + s(n)) space. To this end, it is worth noting that the inputs for BWT are typically large and a fast algorithm that is in-place or uses very low additional memory, would be relevant in practice.

Our theoretical study has also an impact on the practical algorithm design. A natural question is what we get if we allow for some extra space. We prove that using O(k • σ k ) additional space for any given parameter k ≤ n, where σ k ≤ min{|Σ|, k} is the maximum number of distinct characters found among k consecutive positions in T , we can compute the BWT (and its inverse) of a text of n characters in O((n 2 /k + n) log k) time in the comparison model. We observe that σ k is practically a constant in many applications. The practical implications of this trade-off in space versus time can be appreciated by observing that, for any arbitrarily small ǫ > 0, we can obtain the BWT of a text of n bytes in O(nǫ -1 log n) time using just ǫn extra bytes for a constant-size alphabet. This is useful when the text occupies a great part of the available main memory, and only ǫn free cells are available. This avoids using external-memory algorithms, which are clearly slower as I/O access takes several orders of magnitudes more time than main memory access.

The paper is organized as follows. We describe how to perform the inplace BWT in Section 2. We then discuss how to invert the BWT, so as to obtain the original text T , in Section 3. We illustrate the trade-off between space and time in Section 4. Finally, we draw some conclusions in Section 5.

In-Place BWT

Given the input text T = T [0 . . n -1] where T [n -1] = $, moving a single character inside T can change the content of many suffixes. The idea to circumvent this difficulty without using storage for the suffix sort is to proceed by induction from right to left in T , while maintaining the BWT of the current suffix T s , denoted by BWT(T s ). We assume 0 ≤ s ≤ n -3, since the last two suffixes of T are equal to their respective BWT.

To compute BWT(T s ), suppose that BWT(T s+1 ) has been already computed and stored in the last positions of T , i.e. T [s + 1 . . n -1]. Consider the current character c = T [s]: if we look at the content of T [s . . n -1], we no longer find T s , but the character c followed by the permutation BWT(T s+1 ) of T s+1 . Nevertheless, we still have enough information as we will show in the proof of Theorem 1 that the position of $ inside BWT(T s+1 ) is related to the rank of T s+1 among the suffixes T s+1 , . . . , T n-1 in lexicographic order. We exploit this fact in the following steps (see Figure 1). 1. Find the position p of the $ in T [s + 1 . . n -1]: note that p -s is the (local) rank of the suffix T s+1 that originally was starting at position s + 1.

a a α α b b b ′ b ′ β β γ γ c c $ $ BWT(T s ) BWT(T s+1 ) s s p p r r
2. Find the rank r of the suffix T s (originally in position s). Using character c, scan T [s + 1 . . n -1] and compute the sum of two counts: how many characters are strictly smaller than c, and how many occurrences of c appear in T [s + 1 . . p] (and add s as an offset to obtain r). use capital letters to denote the BWT partially built on the last positions of T . Suppose that we have already computed the BWT for the last 7 characters in T , namely, we have missiIPSPIS$. We then have p = 11 and, since there is one character ($) smaller than c = i, and two characters that are equal to c and occur before position p, we have r = s + 3 = 7. This means that we have to replace $ by c and shift IPS by one position left so as to insert $ in position r. The next configuration is missIPS$PISI, which ends with the BWT of T s . Proof : We prove first the correctness. Let T be the input text and T ′ be its modification at a generic iteration s, where 0

Store c into T [p] (thus replacing the $).

≤ s ≤ n -3. Note that T ′ [0 . . s] = T [0 . . s] while T ′ [s + 1 . . n -1] = BWT(T [s + 1 . . n -1]). By induction, the position p of $ in T ′ [s + 1 . . n -1]
indicates the rank p -s of T s+1 among the suffixes in {T s+1 , T s+2 , . . . , T n-1 } in lexicographic order. The base case for T n-2 and T n-1 is trivially satisfied. Hence, we show how to preserve this property for 0 ≤ s ≤ n -3.

First note that the character c = T [s] goes in position p, since it precedes T s+1 inside T . Next, we have to find the new position r for T s , so that r -(s-1) is its rank among the suffixes in S = {T s , T s+1 , . . . , T n-1 } in lexicographic order. First count how many characters smaller than c occur in T ′ [p . . n -1]: there are as many suffixes in S that are smaller than T s since their first character is smaller than c. To this quantity, add the number of occurrences of c in T ′ [s + 1 . . p -1]: these suffixes are also smaller since they start with c but have rank smaller than p, i.e. the rank of T s+1 . In this way, we discover how many suffixes are smaller than T s in S: inserting $ in the corresponding location r of T ′ , by shifting the characters in T ′ [s + 1 . . r] to the left, which thus occupy positions s, . . . , r -1 (see Figure 1), we maintain the induction. Hence,

T ′ [0 . . s -1] = T [0 . . s -1] and T ′ [s . . n -1] = BWT(T [s . . n -1]).
When s = 0, we obtain the BWT of T .

As for the complexity, note that each of the n -2 iterations requires O(n) time, since it can be implemented by O(1) scans of T ′ [s . . n -1]. This gives a total cost of O(n 2 ). We use four integer variables (i, p, r, s) and one character variable (c) in the C code shown in Fig. 2, and thus we need O(1) memory cells for the local variables.

Inverting the BWT

Reversing the permutation performed by the in-place BWT is called inverting the BWT. Initially we have the BWT of the original input text T , denoted BWT(T ). We want to invert the latter by permuting its characters. Thus we reverse the approach described in Section 2. We maintain the invariant that there is a pointer L to a certain position in the input buffer storing BWT(T ) so that, at any time, (a) the prefix of the buffer to the left of L stores the prefix of T obtained so far by the inverting process and (b) the remaining suffix of the buffer (pointed by L till the end of the input buffer) stores the portion of the BWT still to be inverted. For the sake of notation, we identify L with the entire suffix of the input buffer that still has to be inverted.

Under this invariant, which is initially true by setting L to the beginning of the input buffer for BWT(T ), we proceed as follows. We find the position p of $ in L, and then select the pth character in the alphabet order in the multiset given by the characters of L. Stability is needed, since equal characters should be considered in the order of their appearance in L, as detailed below.

1. Find the position p of the $ in L, and increment p (since array indexing starts from 0).

2. Let select be a selection algorithm that works on read-only input, i.e., it does not move elements around while finding the pth smallest element. Using select on L, select the pth character c in the multiset of the characters of L or, equivalently, the pth character in the sorted list of chracters of L.

3. Let q denote the position inside L of the f th occurrence of c, which we hit in a stable fashion when finding c in L. Here f is the difference between p and the number of characters c ′ of L such that c ′ < c.

4.

Replace the occurrence of c at position q by $, and remove the old occurrence of $ by shifting to the right the first p characters of L.

5. At this point, the first position in L is free: store the character c in it, and shorten L by one character at the beginning (i.e. advance the pointer L by one position towards the end of the input buffer).

The C code in Fig. 3 implements Steps 1-5 above, where END_MARKER denotes $. Note that it is a bit longer than the code for the in-place BWT in Fig. 2. As it can be seen below, the original text is reconstructed from left to right as a prefix of increasing length (indicated with small letters).

IPSSM$PISSII → mIPSS$PISSII → miIPSSPISSI$ → misIPSSPIS$I → missIPS$PISI → missiIPSPIS$ → missisIPSPI$ → mississIP$PI → mississiIPP$ → mississipIP$ → mississippI$ → mississippi$
The proof of correctness proceeds along the same lines as in the proof of Theorem 1, since we are reversing the procedure described there. As for the complexity, each of the n -2 iterations is dominated by the cost of select.

Let t s (n) be the time complexity in the comparison model and s s (n) be the space complexity required by select. Using the result in [START_REF] Vinay | Selection from read-only memory and sorting with minimum data movement[END_REF], we have t s (n) = O(n 1+ǫ ) in the worst case for any fixed small constant ǫ > 0 with s s (n) = O(1), and we have t s (n) = O(n 1+ǫ ) = O(n log log n) on the average (which meets the randomized lower bound in [START_REF] Chan | Comparison-based time-space lower bounds for selection[END_REF]), with s s (n) = O(1).

We can state the complexity in general terms.

Theorem 2. Let t s (n) be the time complexity in the comparison model and s s (n) be the space complexity required by select. Given the BWT of a text T of n characters, we can recover T by permuting the BWT (also known as inverse BWT) in O(n•t s (n)) time in the comparison model using O(1+s s (n)) additional memory cells.

We give some examples of the bounds that can be attained with Theorem 2.

Corollary 1. Given the BWT of a text T of n characters, we can recover T by inverting the BWT in O(n 2+ǫ ) time in the worst case, or O(n 2 log log n) time on the average, in the comparison model using O(1) additional memory cells.

Using slightly more additional space than a constant-literally speaking, the algorithm is no more in-place-and the result in [START_REF] Raman | Improved Upper Bounds for Time-Space Trade-offs for Selection[END_REF], where t s (n) = O(n(log n) 2 ) and s s (n) = O(log n), we derive the following.

Corollary 2. Given the BWT of a text T of n characters, we can recover T by inverting the BWT in O((n log n) 2 ) time in the comparison model using O(log n) additional memory cells.

Finally, for the special case in which the alphabet of the distinct characters in T is of constant size (as in DNA and ASCII texts), we obtain an improved bound since select can be immediately implemented by a simple scheme that employs O(|Σ|) = O(1) counters.

Practical Trade-Off between Space and Time

The inplace algorithms described so far have the drawback of requiring Ω(n 2 ) time, which make them unfeasible for long texts. A natural question is how much the latter bound can be improved using extra space. For example, using the dynamic wavelet tree data structure [START_REF] Navarro | Optimal Dynamic Sequence Representations[END_REF] in additional O(n + |Σ| log n) bits of space, we can maintain the BWT through insertion and deletion operations of individual symbols, supporting rank and select operations, with a cost of O(log n/ log log n) time per operation. Using the latter data structure, our algorithms in Sections 2 and 3 would give a bound of O(n log n) time with additional O(n + |Σ| log n) bits of space besides that needed for storing the n characters of the input text T . 3 However the resulting solution is not very practical as the data structure in [START_REF] Navarro | Optimal Dynamic Sequence Representations[END_REF] is quite sophisticated. We show next how to smoothly adapt our algorithms in Sections 2 and 3 to a situation where extra memory is allowed, producing some trade-off solutions that are amenable for implementation with a flexible parameter k for the additional space.

Theorem 3. Given a text T of n characters, we can compute its Burrows-Wheeler Transform (BWT) and its inverse in O((n 2 /k + n) log k) time in the comparison model using O(k • σ k ) additional space, where σ k ≤ min{|Σ|, k} is the maximum number of distinct characters found among k consecutive positions in T .

To appreciate the bound in Theorem 3 from a practical point of view, consider the situation in which the text T occupies a great part of the available memory, and the remaining free cells are a constant fraction of the text size. Our algorithm takes O(n log n) time by fixing k to be a suitable fraction of n. This avoids to use external-memory algorithms, which are clearly slower as I/O access takes several order of magnitudes with respect to main memory access. In general, if the available memory size is M , we obtain the following result by setting k = Θ(M -n).

Corollary 4. Let M be the number of available cells in main memory. Given a text T of n < M characters over a constant alphabet, we can compute its Burrows-Wheeler Transform (BWT) and its inverse in O(( n 2 M -n + n) log n) time in the comparison model using ≤ M total memory cells including those containing T . When M ≥ (1 + ǫ)n for a constant ǫ > 0, this gives O(n log n) time using just ǫn additional cells.

The idea to prove Theorem 3 is to have n/k batches. Each batch simulates k consecutive iterations in the external for loop on s in Figure 2, taking O((n + k) log k) time and using O(k • σ k ) space as follows.

Base case. Let s 1 be the largest multiple of k that is smaller than n. We can compute BWT(T [s 1 . . n -1]) and store it in

T [s 1 . . n -1] in O(k log k) time and O(k) additional space by observing that |T [s 1 . . n -1]| < k.
Inductive case. Suppose that BWT(T [s 1 . . n-1]) has been already stored in T [s 1 . . n -1], where s 1 > 0 is now a generic multiple of k. Letting s 0 = s 1 -k, we want to show how to store BWT(T [s 0 .

. n -1]) in T [s 0 . . n -1] in O((n + k) log k) time using O(k • σ k ) additional space.
Since we have one base case and ≤ n/k inductive cases, the final cost will be O(k

log k + (n/k) • (n + k) log k) = O((n 2 /k + n) log k) time using O(k • σ k )
additional space, as stated in Theorem 3.

Inductive case.

We can abstract the problem for a string X (i.e. T [s 0 . . n -1]) of length m, where the characters in X[k . . m -1] are already permuted according to their BWT, and the characters in X[0 . . k -1] are still in their original order. We want to compute BWT(X) in O((m+k) log k) time using O(k

•σ k ) additional space. Let Z denote X[k . . m -1]
where the $ character is virtually removed from its position, say j. Hence Z is of length m -k -1 and the pair $, j is a breakpoint for Z. In general, a breakpoint is a pair c, j such that c is virtually occupying position j of Z: if two or more breakpoints claim the same position j, there should be a relative order among them.

Our goal is to compute the k + 1 breakpoints for Z so that (a) their characters are those in X[0 . . k -1] plus $, and (b) flattening Z and these breakpoints correctly produces BWT(X) as follows. Given Z and an ordered list of k + 1 breakpoints B = c 0 , j 0 , . . . , c k , j k , where 0 ≤ j 0 ≤ • • • ≤ j k ≤ |Z|, flattening Z and B produces a string with the characters of Z suitably shifted to the left to make room for the characters in the breakpoints of B as follows. We scan Z starting with j = 0: if j = j r for the breakpoint c r , j r at the beginning of B, we output c r and remove c r , j r from B; else (j = j r ), we output Z[j] and increase j. (If c r = $, we do not really output it, but we retain its position j r for the next batch.) The computation ends when Z has been completely scanned and the list B has been emptied. The required time is O(m log k) and the computation can be performed using O(k • σ k ) additional space.

For this we need the following auxiliary data structures for string Z, which require O(k • σ k ) additional space. (Note that Z and X require just O(1) space as they originate from T .)

1. Static array C of σ k ≤ k entries, where α 1 < • • • < α σ k are the distinct characters in X[0 . . k -1]: entry C[i] is the number of positions j in Z such that Z[j] < α i 2.
Static rank data structure R 1 on Z supporting queries that, for an integer j ′ and a character α i , report how many positions j satisfy Z[j] = α i and j ≤ j ′ .

3. Dynamic list B of breakpoints, initially containing only the pair $, j .

4. Dynamic rank data structure R 2 on B supporting queries that, for a character α i , report how many breakpoints c, l to the left of $, j in B satisfy c = α i .

As it is clear, we want to populate the list B by simulating the algorithm described in Section 2. Namely, we want to find the breakpoint of character

X[s ′ ] for s ′ = k -1, k -2, . . . , 0 in Z.
Consider the breakpoint $, j , which exists in B by construction. Let α i = X[s ′ ], and r ′ = r 0 + r 1 + r 2 be sum of three quantities: the number r 0 of positions j ′ in Z such that Z[j ′ ] < α i , the number r 1 of positions j ′ in Z such that Z[j ′ ] = α i and j ′ ≤ j, and the number r 2 of breakpoints c, l that are to the left of $, j in B and have c = α i . Note that we can compute r 0 using entry C[i] in point 1, r 1 using the data structure R 1 in point 2, and r 2 using the data structures B and R 2 in points 3-4.

Fact 1. The value of r ′ is the rank of X[s ′ . . m -1] among X[s ′ + 1 . . m -1], X[s ′ + 2 . . m -1], . . . , X[m -1 . . m -1] in lexicographic order.
Proof : It follows from the fact that X ≡ T [s 0 . . n-1] and X[s ′ +d . . m-1] ≡ T s+d , where s = s 0 + s ′ and d ≥ 0.

After computing r ′ , we replace the breakpoint $, j by c, j , and create a new breakpoint $, r ′ to be inserted in B, updating the data structures in points 3-4 accordingly. Proof : We scan X[0 . . k -1] and find the distinct characters

α 1 < • • • < α σ k .
We then perform a scan of Z to store in C [START_REF] Adjeroh | The Burrows-Wheeler Transform: Data Compression, Suffix Arrays, and Pattern Matching[END_REF] the number of positions j such that Z[j] < α 1 and, for i > 1, to store in C[i] the number of positions j such that α i-1 ≤ Z[j] < α i . We then store in C its prefix sums, thus giving the wanted C. Time is O(m log σ k ) since during the scan we perform a binary search among the σ k characters. Space is O(σ k ) by definition of C.

Lemma 2. The static data structure R 1 can be stored in O(k • σ k ) space and built in O(m log σ k ) time, so that each query requires O(log σ k + m/k) time.

Proof : We store R 1 as a collection of σ k arrays R i 1 , for 1 ≤ i ≤ σ k . Entry R i 1 [h],
for 0 ≤ h ≤ k, stores the number of occurrences of character α i in the hth segment of m/k consecutive positions in Z: namely, R i 1 [0] = 0 and, for h ≥ 1, R i 1 [h] stores the number of positions j such that Z[j] = α i and (m/k) • (h -1) ≤ j ≤ max{m -1, (m/k)h -1}. After initializing the entries in all the arrays to zero, their correct value can be set with a single scan of Z in O(m log σ k ) time. 4 After that, we perform a post-processing and store in each R i 1 its prefix sums: in this way, R i 1 [h] stores the number of positions j such that Z[j] = α i and 0 ≤ j ≤ max{m-1, (m/k)h-1}. For a query with a character c ′ and an integer j ′ , it takes O(log σ k ) time to establish that c ′ = α i for a certain i, and O(1) time to find the largest h such that (m/k)h < j ′ . The correct answer for the query is then given by the sum of the content of R i 1 [h] and the number of α i 's found in Z[(m/k)h . . j ′ ] by its direct inspection. Scanning the latter takes O(m/k) time by definition of R i 1 .

Lemma 3. The dynamic list B and the dynamic data structure R 2 can be stored in O(k) space and built in O(k log k) time, so that each query and each update requires O(log k) time.

Proof : We handle the list B = c 0 , j 0 , . . . , c r , j r , where 0 ≤ r ≤ k, and the data structures R 2 together. In particular, R 2 is the wavelet tree of height O(log r) built on the sequence c 0 • • • c r of characters taken from B. The query for α i can be performed as a count query of characters α i in c 0 • • • c d , where d ≤ r and c d = $ (e.g. [START_REF] Navarro | Wavelet trees for all[END_REF]). Each update (insertion, deletion, replacement) can be handled in O(log r + log k) = O(log k) time [START_REF] Grossi | The wavelet trie: maintaining an indexed sequence of strings in compressed space[END_REF][START_REF] Navarro | Optimal Dynamic Sequence Representations[END_REF].

We now have all the ingredients to prove the time and space bounds for computing the BWT as stated in Theorem 3. By the inductive scheme and Fact 1, the computation is correctly performed. As for the time bounds, we use Lemma 1-3. Note that to invert the BWT, we can now implement the algorithm described in Section 3 in a simpler way, since R 1 and R 2 support also the selection of the pth symbol c = α i . Moreover, flattening Z and B removes the positions j stored in the pairs in B from Z, as this simulates the shift of some characters of Z to the right. The time and space analysis is similar to that of computing the BWT.

Conclusions

We presented an in-place BWT construction taking O(n 2 ) time in the comparison model. It would be interesting to improve this bound. Note that the while loop in our in-place BWT can be avoided using O(|Σ|) space, where Σ is the alphabet of characters occurring in T . Time can be further reduced to O(n 2 / log |Σ| n) by packing characters but this is still not useful for large text collections.

We do not know whether a lower bound better than Ω(n log n) holds for the problem in the comparison model since the space is very constrained. This is an interesting question to investigate.

On the practical side, our in-place algorithms can be adpated to provide a trade-off between space and time, when O(k • σ k ) extra memory cells are allowed, providing an O((n 2 /k + n) log k)-time algorithm to obtain (and invert) the BWT. For example in real systems where the alphabet size is a constant, for any arbitrarily small ǫ > 0, the BWT of a text of n bytes can be computed in O(nǫ -1 log n) time using just ǫn extra bytes.

APPENDIX: REPLY TO REVIEWERS' COMMENTS

We report the Reviewers' comments and our responses marked with "***".

REVIEWER #1

Comments on the revised version of the manuscript. While the contribution looks good, there are some issues with technical details that would warrant another revision.

In this revision of the manuscript the authors have made improvements, in particular added section 4 that describes time/space trade-offs for building and inverting the BWT with a small amount of extra space. It is essentially a variant of the algorithm of Hon et al. [START_REF] Grossi | The wavelet trie: maintaining an indexed sequence of strings in compressed space[END_REF], optimized for space usage. In particular, the improvements in space usage come from time/space trade-offs in the rank/select structure over BWT.

There are still some issues:

-The time bounds with \epsilon n bytes of extra space mentioned in the abstract and in sections 1, 4.1, and 5 assume a constant alphabet. This should be mentioned. There is also some confusion with s and k in section 5. *** done *** -Section 4.1: The initial discussion of breakpoints is a bit misleading, as it gives an impression that $ is going to be inserted back to position j. This is rectified only on the next page. The revised paper is stronger; it would still benefit from full proof-reading, suggestions: *** sorry for the typos *** + p3 line 9: "it is not rare the case when the input data is so large that" -> frequently the input data is so large that / not infrequently, ... / it is not rare for the input data to be so large that the input text T stays in main memory while any additional

Figure 1 :

 1 Figure 1: An illustration of Steps 1-4 of the in-place construction of the BWT.

4 .Figure 2 :

 42 Figure 2: In-place construction of BWT.

Theorem 1 .

 1 Given a text T of n characters, we can compute its Burrows-Wheeler Transform (BWT) in O(n 2 ) time in the comparison model using O(1) additional memory cells.

Corollary 3 .

 3 Given the BWT of a text T of n characters drawn from a constant-size alphabet, we can recover T by inverting the BWT in O(n 2 ) time in the comparison model using O(1) additional memory cells.

Lemma 1 .

 1 The static array C can be stored in O(k • σ k ) space and built in O(m log k) time.

* 20 -

 20 ** done *** -Section 4.1, paragraph 2: The length of Z is m-k-1. *** done *** -Section 4.1, paragraph 3: There may be a breakpoint at position |Z|, if the lexicographically largest suffix starts in X[0,k-1]. *** done *** Section 4.1, bullet 1: Is array C based on Z[j] \le \alpha_i or Z[j] < \alpha_i? *** done, it is Z[j] < \alpha_i *** -Section 4.1, lemmas 1 and 2: If the structure is built in O(n log k) time, the overall time complexity becomes O(n^2/k log k) instead of O(n^2/k + n log k). *** done, we wrote O((n^2/k + n) log k) *** -Section 4.1, lemma 2: For large alphabets and large values of k, the space usage may be larger than the time bound O(n log k). *** done, added a footnote sayng that it is a standard trick to initialize only what is REVIEWER #2

  *** done *** + p3 line -4: saving over the space -> saving on the space *** done *** + p4 line -6: even allowing for exponential time -> even if it results in exponential time *** done *** + p5 line -7: This avoids us to use external-memory algorithms, which are clearly slower as I/O access takes several order of magnitudes with respect to main memory access. -> This avoids using ........ several order_s_ of magnitude_ _more time_ / longer *** done *** + p16 line 8: The time and space analysis is similar that of computing the BWT. -> ... similar to that of ... *** done *** + p16 line 14: but still not useful -> but this is still not useful *** done *** REVIEWER #3 The paper looks acceptable now. Minor comments: Footnote 3 -is this a private communication? or give reference? *** done, it is a private communication *** The CPM 2013 paper should be included in the References. *** done ***

Table 1 :

 1

	cyclic shifts	sorted cyclic shifts	suffixes
		L	i T i
	mississippi$	$mississipp i	11 $
	$mississippi	i$mississip p	10 i$
	i$mississipp	ippi$missis s	7 ippi$
	pi$mississip	issippi$mis s	4 issippi$
	ppi$mississi	ississippi$ m	1 ississippi$
	ippi$mississ	mississippi $	0 mississippi$
	sippi$missis	pi$mississi p	9 pi$
	ssippi$missi	ppi$mississ i	8 ppi$
	issippi$miss	sippi$missi s	6 sippi$
	sissippi$mis	sissippi$mi s	3 sissippi$
	ssissippi$mi	ssippi$miss i	5 ssippi$
	ississippi$m	ssissippi$m i	2 ssissippi$

BWT L for the text T = mississippi$ and its relation with suffix sort.

As is standard in many string algorithms, we assume that any two characters in Σ can only be compared and this takes O(1) time. Hence, comparing characterwise any two suffixes may require O(n) time in the worst case.

In C code, we would declare T as unsigned char T[n] and use this storage plus O(1) local variables of constant size. A more formal model would say that each memory cells hosts a character from Σ and so an integer of log n bits requires log |Σ| n cells. We prefer to keep it simpler and say that an extra cell can contain an integer.

This theoretical solution has been suggested by Rossano Venturini (private communication).

If k • σ k = w(n log σ k ), it is a standard trick to allocate O(k • σ k ) memory and initialize only what is needed in in O(m log σ k ) time[START_REF] Aho | The Design and Analysis of Computer Algorithms[END_REF].

Acknowledgments

The second author is grateful to Gianni Franceschini for some preliminary discussions and to Venkatesh Raman for pointing out the results in [23,28] and Rossano Venturini for his comments and indicating the result in [3]. We also thank the anonymous reviewers for their careful reading of our manuscript.

✩ A preliminary version of the results in this paper appeared in [6]. The work of the second author has been supported in part by the Italian Ministry of Education, University, and Research (MIUR) under PRIN 2012C4E3KT project. The work of the third author has been supported by the Academy of Finland grant 118653 (ALGODAN). The work of the fourth author has been partially supported by the National Science Foundation Award 0904246, Israel Science Foundation grant 347/09, Yahoo, Grant No. 2008217 from the United States-Israel Binational Science Foundation (BSF) and DFG.

void inplaceIBWT( unsigned char L[ ], int n ){ int f, i, p, q, count; unsigned char c; /* step 1 */ p = 0; while( L[ p ] != END_MARKER ) p++; p++;

p = q+1; /* also the new END_MARKER has been shifted */ else p = q; } }