{'fraction_non_alphanumeric': 0.05011571568904947, 'fraction_numerical': 0.012620502489136045, 'mean_word_length': 3.7709905473032013, 'pattern_counts': {'":': 0, '<': 22, '<?xml version=': 0, '>': 13, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 17, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Lightweight LCP Construction for Very Large Collections of Strings

The longest common prefix array is a very advantageous data structure that, combined with the suffix array and the Burrows-Wheeler transform, allows to efficiently compute some combinatorial properties of a string useful in several applications, especially in biological contexts. Nowadays, the input data for many problems are big collections of strings, for instance the data coming from"next-generation"DNA sequencing (NGS) technologies. In this paper we present the first lightweight algorithm (called extLCP) for the simultaneous computation of the longest common prefix array and the Burrows-Wheeler transform of a very large collection of strings having any length. The computation is realized by performing disk data accesses only via sequential scans, and the total disk space usage never needs more than twice the output size, excluding the disk space required for the input. Moreover, extLCP allows to compute also the suffix array of the strings of the collection, without any other further data structure is needed. Finally, we test our algorithm on real data and compare our results with another tool capable to work in external memory on large collections of strings.


Introduction
The suffix array (SA), the longest common prefix array (LCP) and the Burrows-Wheeler transform (BWT) are data structures with many important applications in stringology [22].
The LCP array of a string contains the lengths of the longest common prefixes of the suffixes pointed to by adjacent elements of the suffix array of the string [19]. The most immediate utility of the LCP is to speed up suffix array algorithms and to simulate the more powerful, but more resource consuming, suffix tree. Indeed, the LCP array, combined with SA or BWT, simplifies the algorithms for some applications such as the rapid search for maximal exact matches, shortest unique substrings and shortest absent words [23,5,11,1].
Real world situations may present us with datasets that are not a single string but a large collection of strings, examples being corpora of web pages or the data coming from "next-generation" DNA sequencing (NGS) technologies. It is common for the latter in particular to generate aggregations of hundreds of millions of DNA strings (sometimes called "reads"), leading to a need for algorithms that deal not only with the collective nature of the data but also with its large overall size. In this paper we give an algorithm for the construction of LCP array that addresses both these requirements.
The importance of the LCP array has meant that its construction has been wellstudied in the literature. For instance, there are algorithms that work in linear time and in internal memory (cf. [15]). Other algorithms work in semi-external memory (see for instance [13]) or directly via BWT (see [5]). In [14], the first algorithm establishing that LCP array can be computed in O(sort(n)) complexity in the external memory model (i.e., the complexity of sorting n integers in external memory) is presented. Recently, in [6] an external memory algorithm to construct the suffix array of a string based on the induced sorting principle is provided. Such an algorithm can be augmented to also construct the LCP array. The overhead in time and I/O operations for this extended algorithm over plain suffix array construction is roughly two. In another recent paper [12], the authors introduced an external memory LCP array construction algorithm for a string. The main idea in [12] is to divide the string into blocks that are small enough to fit in RAM and then scan the rest of the string once for each block. Such a strategy needs 16n bytes of disk space, where n is the length of the string.
One's initial strategy for building the LCP array of a collection might therefore be to concatenate its members into a single string and then apply one of the above methods. However, some stumbling blocks become apparent, and it is not immediate how these single-string approaches can be adapted to circumvent them. First, many existing algorithms for computing the LCP array require data structures of size proportional to the input data to be held in RAM, which has made it impractical to compute the LCP array of a very large string. Moreover, by definition, the values in LCP array should not exceed the lengths of the strings, so one would ideally employ distinct end-marker symbols to act as separators between the strings of the collection. However, assigning a different end-marker to each string is not feasible when the number of strings in the collection is very large. On the other hand, the use of the same end-marker symbol throughout the collection could give rise to LCP values with the undesirable properties of sometimes exceeding the lengths of the strings and having dependencies on the order in which the strings are concatenated.
In the literature, the problem of the computation of LCP array for a collection of strings has been considered in [26] and in [18]. Moreover, a preliminary version of the results presented in this paper is in [4]. In particular, this paper includes full proofs of theorems from [4] and more detailed examples. Furthermore, additional experimental results are also described. Finally, in this paper we introduce a simple solution for dealing with strings having different lengths, allowing extLCP to work on any collection of strings.
Defining N and K as the sum of the lengths of all strings and the length of the longest string in the collection respectively, the approach in [26] requires O(N log K) time, but the O(N log N ) bits of internal memory needed to store the collection and its SA in internal memory prevents the method from scaling to massive data.
The external memory algorithm (called eGSA) in [18] builds both suffix and LCP arrays for a collection of strings. Such interesting strategy has an overhead in working disk space, so it seems that eGSA cannot be used for huge collections of strings.
The methodology presented in this paper attempts to overcome the above mentioned limitations. In fact, the aim of our paper is to design a lightweight algorithm for the computation, at the same time, of the longest common prefix array and the Burrows-Wheeler transform of a very large collection of strings having different or same length. The computation is realized by performing disk data accesses only via sequential scans. The total disk space usage never needs more than twice the output size, excluding the disk space required for the input.
In our approach, we directly compute both data structures, without needing to concatenate the strings of the collection and without requiring pre-computed auxiliary information such as the suffix array of the collection.
In order to do this, our algorithm is built upon the approach introduced in [3] related to an extension of the Burrows-Wheeler transform to a collection of strings defined in [20].
In particular, our algorithm (called extLCP) adds to the strategy in [3] (called BCR) some lightweight data structures, and allows the simultaneous computation of both the longest common prefix and the Burrows-Wheeler transform of a collection of m strings. Such a computation is performed in O((m+σ 2 ) log N ) bits of memory, with a worst-case time complexity of O(K(N + sort(m))), where sort(m) is the time taken to sort m integers, σ is the size of the alphabet, N is the sum of the lengths of all strings and K is the length of the longest string. Note that extLCP needs at most (2N − m)(log σ + log K) + N log σ bits of disk space and requires The low memory requirement enables our algorithm to scale to the size of dataset encountered in human whole genome sequencing datasets: in our experiments, we compute the BWT and LCP of collections as large as 800 million strings having length 100.
Moreover, extLCP allows to compute also the suffix array of the strings of a collection (called generalized suffix array). Such a further computation increases the number of I/O operations, but it does not need further data structures in internal memory.
Section 2 gives preliminaries that we will use throughout the paper, in Section 3 we define the main data structures for a collection of strings, Section 4 describes an incremental strategy for computing the BWT of a collection of strings of any length. Section 5 shows the sequential computation of the LCP array for the collection. We present details on the efficient implementation of the algorithm and its complexity in Sections 6 and 7, respectively. Computational results on real data are described in Section 8. Section 9 is devoted to some conclusions.

Preliminaries
Let Σ = {c 1 , c 2 , . . . , c σ } be a finite ordered alphabet with c 1 < c 2 < . . . < c σ , where < denotes the standard lexicographic order. We append to a finite string w ∈ Σ * an end-marker symbol $ that satisfies $ < c 1 . We denote its characters by w [1], w [2], . . . , w[k], where k is the length of w, denoted by |w|. Note that, , with a substring w [1, j] being called a prefix, while a substring w[i, k] is referred to as a suffix. A range is delimited by a square bracket if the correspondent endpoint is included, whereas the parenthesis means that the endpoint of the range is excluded.
We denote by S = {w 0 , w 1 , . . . , w m−1 } the collection of m strings of length at most K. We suppose that to each string w i is appended an end-marker symbol $ i smaller than c 1 , and $ i < $ j if i < j. Let us denote by N the sum of the lengths of all strings in S.
For j = 0, . . . , |w i | − 1, we refer to the suffix w i [|w i | − j, |w i |] of a string w i as its j-suffix ; the 0-suffix of w i contains $ i alone. The length of a j-suffix is equal to j (up to considering the end-marker). Let us denote by S j the collection of the j-suffixes of all the strings of S.
In our algorithm presented in this paper we use a unique end-marker $ = c 0 for all strings in S, because we set w s [|w s |] < w t [|w t |] if and only if s < t, so that if two strings w s and w t share the j-suffix, if and only if s < t. However, to ease the presentation distinct end-markers are shown.
We say that the symbol 3. SA, LCP and BWT of a collection of strings Suffix array, longest common prefix array and Burrows-Wheeler transform are all commonly defined with reference to a single string. This section describes the extension of such notions to a collection of strings.
The suffix array SA of a string w is an array containing the permutation of the integers 1, 2, . . . , |w| that arranges the starting positions of the suffixes of w into lexicographical order. There exist some natural extensions of the suffix array to a collection of strings (see [26]).
We define the generalized suffix array GSA of the collection S = {w 0 , w 1 , . . . , w m−1 } as the array of N pairs of integers (t, j), corresponding to the lexicographically sorted suffixes w j [t, |w j |], where 1 ≤ t ≤ |w j | and 0 ≤ j ≤ m − 1.
In particular, GSA[q] = (t, j) is the pair corresponding to the q-th smallest suffix of the strings in S, i.e. to the suffix w j [t, |w j |].
The longest common prefix array of a string contains the lengths of the longest common prefixes of the suffixes pointed to by adjacent elements of SA of the string [24]. The longest common prefix array LCP of a collection S of strings, denoted by lcp(S), is an array storing the length of the longest common prefixes between two consecutive suffixes of S in the lexicographic order. For every j = 1, . . . , N − 1, if GSA[j − 1] = (p 1 , p 2 ) and GSA[j] = (q 1 , q 2 ), LCP [j] is the length of the longest common prefix of suffixes starting at positions p 1 and q 1 of the words w p2 and w q2 , respectively. We set LCP [1] = 0.
For i < j, a range minimum query RM Q(i, j) on the interval [i, j] in the LCP array returns an index k such that LCP [k] = min{LCP [l] : i ≤ l ≤ j}. If GSA[r] = (p 1 , p 2 ) and GSA[s] = (q 1 , q 2 ), it is not difficult to show that the length of the longest common prefix between the suffixes starting at positions p 1 and q 1 of the words w p2 and w q2 corresponds to LCP [RM Q(r + 1, s)].
The suffix array of a string is related to the Burrows-Wheeler transform introduced in [7]. The original Burrows and Wheeler transform (BWT) on a string is described as follows: given a word w ∈ Σ * , the output of BWT is the pair (L, I) obtained by lexicographically sorting the list of the conjugates of w. In particular, L is the word obtained by concatenating the last symbol of each conjugate in the sorted list and I is the position of w in such a list. For instance, if w = mathematics then (L, I) = (mmihttsecaa, 7) 1 .
Actually, in several implementations of the BWT, in order to improve the efficiency of the computation, one can consider a variant of the BWT by the sorting the suffixes of w$ rather than the conjugates of w. To ensure the reversibility of the transform, one needs to append the symbol $ at the end of the input string . In this variant, the output bwt(w) is a permutation of w, obtained as concatenation of the letters that (circularly) precede the first symbol of the suffix in the lexicographically sorted list of its suffixes: for i = 1, . . . . In other words, the i-th symbol of the BWT is the symbol just before the i-th suffix. For instance, if w = mathematics$ then bwt(w) = smmihtt$ecaa. Note that, in this case, the second output of the BWT (the index I) is not useful, because one can use the position of the $-symbol for recovering the input string.
The Burrows-Wheeler transform can be extended to a collection of strings. In its original definition [20] (see also [21]), such a reversible transformation (called EBWT) produces a string that is a permutation of the characters of all strings in S and it does not make use of any end-marker. The EBWT of a collection S is a word (denoted by c-ebwt(S))) obtained by letter permutation of the words in S together a set of indexes (denoted by I) used to recover the original collection. In particular, c-ebwt(S) is obtained by concatenating the last symbol of each element in the sorted list of the conjugates of the words in S. The sorting exploits an order relation defined by using lexicographic order between infinite words. For instance, if S = {abac, cbab, bca, cba}, the output of EBWT of S is the couple (ccbbbcacaaabba, {1, 9, 13, 14}).
In this paper we use a more efficient variant of EBWT of a collection of strings, that needs to append a different end-marker to each string of S. In this case the word obtained as output, denoted by ebwt(S), is obtained by concatenating the symbols just preceding each suffix of the list of the suffixes of the words in S in according with the lexicographic order. The output ebwt(S) can be also defined in terms of the generalized suffix array of S. In particular, if GSA[i] = (t, j) then ebwt(S)[i] = w j [(t − 1)]; when t = 1, then ebwt(S)[i] = $ j . By using the above example, ebwt(S) = cbaacbb$ 0 bacca$ 2 ab$ 3 $ 1 . The external memory methods for computing the output ebwt(S) are given in [2]. Furthermore, in practice, such methods use a unique end-marker rather than m different end-markers. Note that ebwt(S) differs, for at least m symbols, from bwt applied to the string obtained by concatenating all strings in S.

Computing the EBWT of a collection of strings having any length
Our approach for computing the LCP array is built upon the BCR algorithm introduced in [3] to compute ebwt(S), where S is a collection of strings. Note that BCR in [3] is described for collections of strings of fixed length. In this section we focus on the description of BCR on collections of strings of any length.
The BCR algorithm [3] computes the ebwt(S) of the collection of strings S = {w 0 , w 1 , . . . , w m−1 } without concatenating the strings. In the sequel, we assume that 0 ≤ i ≤ m − 1. Note that we consider the symbol $ i appended to each string w i . We suppose that K is the maximal length (including the end-markers) of the strings in S and N the sum of their lengths (including the end-markers). We assume that K ≪ m.
In Subsection 4.1 we describe how the algorithm works on a collection of strings of any length, Subsection 4.2 is devoted to detail the involved data structures needed to allow a lightweight implementation and to reduce the I/O operations.
4.1. The algorithm. The basic idea of the BCR algorithm is to scan all the strings w 0 , w 1 , . . . , w m−1 in the collection S from right to left at the same time. This means that, at each iteration, it considers a "slice" of (at most) m characters from the strings in the collection. BCR builds incrementally, via K iterations, the Burrows-Wheeler transform of S by simulating, step by step, the insertion of all suffixes having the same length in the list of sorted suffixes.
At the end of each iteration j = 0, . . . , K − 1, BCR algorithm builds a partial ebwt(S) (denoted by ebwt j (S)). It is referred to as partial, because if one inserts the $-characters in their correct position (rather than the symbols that precede the j-suffixes), then one immediately obtains the ebwt(S j ) of all j-suffixes of S (if the length of some string is greater than j). For instance, if S = {AAT ACACT GT ACCAAC$ 0 , GAACAGAAAGCT C$ 1 } (including the distinct end-markers), ebwt 3 (S) corresponds to the computation of ebwt({AAC$ 0 , CT C$ 1 }) when the end-markers are inserted.
The key point of this strategy is to establish the positions where the new symbols associated with the j-suffixes must be inserted in ebwt j−1 (S) in order to obtain ebwt j (S). In other words, we have to find the position of each j-suffix in the list of sorted t-suffixes with t < j, without explicitly computing such a list.
In particular, at the step j = 0, we have to consider the symbols associated with the 0-suffixes and establish how they must be concatenated in order to obtain ebwt 0 (S). Since we use (implicit) distinct end-markers and $ i < $ j if i < j, then it is easy to verify that ebwt 0 (S) is obtained by considering the last symbol of each string w i and by concatenating them in the same order as the strings appear in the collection: At the iteration 1 ≤ j < K − 1, we need to retain ebwt j−1 (S), keep track of the positions within it of the symbols w i [|w i | − j] associated with the (j − 1)-suffixes of S and build ebwt j (S) by inserting at most m new symbols into ebwt j−1 (S), i.e. w i [|w i | − j − 1] for j < |w i | − 1 or $ for j = |w i | − 1. Such operations simulate the insertion of the j-suffixes into the list of sorted suffixes computed in the previous step. Note that if j ≥ |w i | then the string w i is not considered.
The process ends at the step K−1, when BCR considers the j-suffixes of the strings w i with |w i | = K and inserts the end-markers of such strings into ebwt K−2 (S), in order to obtain ebwt K−1 (S), i.e. the ebwt(S).
In order to find, at each step j ≥ 1, the positions where the new (at most) m symbols must be inserted into ebwt j−1 (S), we use the notions related to backward search, such as the table C and the rank function that have been extensively used in FM-index (see [8]). Formally, given a string w, C w [x] is a table of σ integers that, for each character x ∈ Σ, contains the number of occurrences of lexicographically smaller characters in the string w. Moreover, given a symbol x, an integer r and a string w, the function rank(x, r, w) returns the number of occurrences of character x in the prefix w [1, r].
In our context we suppose that, at the step j − 1 ≥ 0, we have inserted the symbol x = w i [|w i | − j] at the position r of ebwt j−1 (S). At the step j ≥ 1, for each i with j < |w i |, the position p, where we have to insert the new symbol associated with the j-suffix w i [|w i |− j, |w i |] into ebwt j−1 (S), is computed in the following way: Note that such a formula corresponds to the computation of the so-called LFmapping, a fundamental operation of the FM-index.

Data Structures.
In order to decrease the I/O operations, BCR does not use the table C and computes the rank function on segments of ebwt j−1 (S) rather than on the entire ebwt j−1 (S). More precisely, BCR considers ebwt j (S) as the concatenation of σ + 1 segments B j (0), B j (1), . . . , B j (σ), where the symbols in B j (0) are the characters preceding the lexicographically sorted 0-suffixes of S j (such suffixes consist in only the end-marker symbols) and the symbols in B j (h), with h ≥ 1, are the characters preceding the lexicographically sorted suffixes of S j starting with c h ∈ Σ. It is easy to see that, ebwt 0 (S) = B 0 (0) and the segments B 0 (h) (for h = 1, . . . , σ) does not contain any symbols. Now, we can omit the table C, because it is only useful for determining the segment where the suffix must be inserted. Indeed, the table C simply allows us to establish that if the j-suffix starts with the symbol c z = w i [|w i | − j] then it is larger than the suffixes starting with a symbol smaller than c z . This is equivalent to say that the j-suffix must be inserted in the block containing the suffixes starting with the symbol c z , i.e. in B j (z).
One can verify that the position r in ebwt j−1 (S) corresponds to a position s in (or the end-marker $ for the last step) must be inserted in the position r of B j (z), where r is obtained by computing the number of occurrences of Moreover, the computation of B j−1 (v) [1, s] corresponds to the computation of the function rank(c z , s, B j−1 (v)). Note that, in order to compute the occurrences in ) bits of memory can be used.
We remark that for each step j, our strategy computes, for each j-suffix τ , its definitive position in the sorted list of the elements of S j , regardless of the fact that some j-suffixes lexicographically smaller than τ have not been considered yet in the current step j. Actually, this means that the algorithm computes the definitive position (called absolute position) of the symbol associated to τ in the correspondent segment B j .
During each iteration j, we use the arrays U , N , P , Q containing at most m elements each. The values in such arrays are updated during each step. For ease of presentation we denote by U j , N j , P j , Q j the arrays at the step j, described as follows. For each q = 0, . . . , m − 1: • U j is an array that stores all the symbols, if they exist, located at the position j + 1 from the right. More formally, Note that # is a symbol that does not appear in any string of the collection and it is not involved in the computation of ebwt j (S). The array U j takes O(m log σ) bits of workspace.
Note that, at each step j ranging from 0 to K − 1, the algorithm considers only the first t values of arrays P j , Q j and N j where t is the number of the strings of the collection having length greater than or equal to j.
At the start of the iteration j, we compute the new values in Q j , P j , N j in the following way. We suppose that the The position r of the new symbol w i [|w i | − j − 1] (or the end-marker $ i for j = |w i |) that must be inserted in B j (z) is obtained by computing the number of Note that during the iteration j, we store the values in P j , Q j , N j and U j by replacing the corresponding values in P j−1 , Finally, we sort Q j , P j , N j where the first and the second keys of the sorting are the values in Q j and P j , respectively. So that we can insert, in sequential way, the symbols U j [q] (for q = 0, . . . , m − 1), into each segment B j (h) for h = 0, . . . , σ.
The algorithm starts, at the step j = 0, by computing the ebwt 0 (S). The 0-suffixes are $ 0 , $ 1 , and the new symbols that we have to insert in B 0 (0) (the segment associated with the suffixes starting with the end-marker) are C and C. We Then, we consider U 1 = [A, T ]. Both symbols in U 1 should be inserted into B 1 (2) because both the associated suffixes (C$ 0 and C$ 1 ) start with the symbol C. So, we set Q 1 = [2,2] and N 1 = [0, 1]. The position in B 1 (2) of the symbol A associated with the 1-suffix of w 0 is equal to 1 in the segment B 1 (2), indeed C$ 0 is the smallest suffix in S 1 starting with the letter C. The position in B 1 (2) of the symbol T associated with the 1-suffix of w 1 is equal to 2 in the segment B 1 (2), because the symbol C associated with $ 1 follows the symbol C associated with $ 0 . So we set P 1 = [1, 2] and obtain B 1 (2) by inserting A and T . Then ebwt 1 (S) = CCAT .
During the third step, the array U 2 = [A, C] is considered. Since the last inserted symbol of w 0 is A then the new symbol A associated with AC$ 0 must be inserted in the segment B 2 (1), whereas since the last inserted symbol of w 1 is T then the new symbol C associated with T C$ 1 must be inserted in the segment B 2 (4). So we set Since the number of occurrences of T in B 1 (0), B 1 (1) and in B 1 (2)[1, 2] is 1 then we have to insert the symbol C at the position 1 in B 2 (4). Since AC$ 0 and T C$ 1 are the smallest suffixes starting with A and T respectively, we set P 2 = [1, 1].
The first three iterations are depicted in Figure 1. The process continues via left extensions of suffixes until all symbols have been inserted.
Note that the list of sorted suffixes is inserted in the figure in order to make it easier to understand. In each table, the first column represents ebwt j (S) = B j (0)B j (1)B j (2)B j (3)B j (4) after each iteration. The new symbols inserted in each iteration are shown in bold. It is verify that at the end of iteration 0, we have

LCP computation of a collection of strings via EBWT
The main goal of this section consists in the description of the strategy for computing, by using the EBWT, the LCP array of a massive collection of strings via sequential scans of the disk data. In particular, the main theorem of the section enables the simultaneous computation of both LCP and EBWT of a string collection S = {w 0 , w 1 , . . . , w m−1 } of maximum length K. We recall that the last symbol of each string w i is the (implicit) end-marker $ i .
Our method follows the BCR algorithm in the sense that it scans all the strings from right-to-left in K steps and simulates, at the step j, the insertion of suffixes of length j into the sorted list of suffixes. This time, however, we wish to compute both the LCP and EBWT. So, at the step j, the longest common prefix array (denoted by lcp j (S)) of the collection S j of the suffixes having length at most j is computed alongside ebwt j (S). As well as computing an LCP value for the inserted suffix, we must also modify the LCP value for the suffix that comes after it to reflect the longest prefix common to it and the inserted suffix.
Our goal in this section is to frame these calculations in terms of RM Qs on sets of intervals within lcp j (S), serving as a preliminary to the next section, where we show how these computations may be arranged to proceed via sequential scans of the data.
It is easy to see that when j = K − 1, lcp j (S) coincides with lcp(S), the LCP array of S. Since all m end-markers are considered distinct, the longest common prefix of any pair of 0-suffixes is 0, so the first m positions into lcp j (S) are 0 for any j ≥ 0.
The general idea behind our method is described in the following and depicted in Figure 2. At each step j > 0, the value of a generic element lcp j (S)[r] (with r > m) have to be computed by taking into account the j-suffix τ that is placed at the position r of the sorted suffixes of S j . Such a value depends on the f -suffix τ 1 placed at the position r − 1 in the sorted list and, moreover, could lead an updating of the value lcp j (S)[r + 1] (if it exists) corresponding to the g-suffix τ 2 , as shown in Figure 2. By using our method, the computation of lcp j (S)[r] can be realized by using the values of the arrays lcp j−1 (S) and ebwt j−1 (S). In fact, if τ = cγ (where c is a symbol and γ is the (j − 1)-suffix placed at a certain position t in the sorted list of the suffixes of S j−1 ), then c = ebwt j−1 (S)[t]. Let us denote τ 1 = aγ 1 and τ 2 = bγ 2 , where a and b are symbols and γ 1 and γ 2 are, respectively, the (f −1)-suffix and the (g − 1)-suffix of S j−1 placed at the positions d 1 and d 2 . One can see that lcp j (S)[r] is equal to 0 if c = a, otherwise it is equal to 1 plus the longest common prefix between γ 1 and γ computed in the array lcp j−1 (S). Moreover lcp j (S)[r + 1] is equal to 0 if c = b, otherwise it is equal to 1 plus the longest common prefix between γ and γ 2 computed in the array lcp j−1 (S). In this section we show how such a computation can be sequentially performed.  Note that lcp j (S) can be considered the concatenation of σ + 1 arrays L j (0), L j (1), . . . , L j (σ) where, for h = 1, . . . , σ, the array L j (h) contains the values corresponding to the lengths of the longest common prefix of the suffixes of S j that start with c h ∈ Σ, while L j (0) (corresponding to the 0-suffixes) is an array of m zeroes. It is easy to see that lcp 0 (S) = L 0 (0) and that L 0 (h) is empty for h > 0. We note that, for each 1 ≤ h ≤ σ, L j (h) [1] = 0 and L j (h)[i] ≥ 1 for i > 1, because the suffixes associated with such values share at the least the symbol c h ∈ Σ.
As shown in Section 4, ebwt j (S) can be partitioned in an analogous way into segments B j (0), B j (1), . . . , B j (σ). Given the segments B j (h) and L j (h), h = 0, . . . , σ, for the symbol x occurring at position r of B j (h) we define the (j, h)-LCP Current Interval of x in r (denoted by LCI h j (x, r)) as the range (d 1 , r] in L j (h) (so we set LCI h j (x, r) = L j (h)(d 1 , r]), where d 1 is the greatest position smaller than r of the symbol x in B j (h), if such a position exists. If such a position does not exist, we de- . Analogously, we define for the symbol x the (j, h)-LCP Successive Interval of x in r (denoted by LSI h j (x, r)) as the range (r, In our context, if d 1 and d 2 exist, then B j (h)[t] = x for t = d 1 + 1, . . . , r − 1 and for t = r + 1, . . . , d 2 − 1, so it is easy to verify that d 1 = select(rank(x, r, B j (h)) − 1, x, B j (h)) and d 2 = select(rank(x, r, B j (h)) + 1, x, B j (h)).
We observe that the computation of the minimum value into LCI h j (x, r) and LSI h j (x, r) is equivalent to the computation of RM Q(d 1 , r) and RM Q(r, d 2 ), respectively. We cannot directly compute these values, because we build each B j (h) and L j (h) in sequential way, so we do not know the left extreme of LCI h j (x, r) and the right extreme of LSI h j (x, r) of all symbols that we have to insert in each B j (h), h = 1, . . . , σ.
The following theorem, related to [10, Lemma 4.1], shows how to compute the segments L j (h), with j > 0, by using L j−1 (h) and B j−1 (h) for any h > 0.
For each position (r b + 1) / ∈ I (where r b ∈ I and 0 ≤ b < l), then Proof. We consider a generic position r ∈ I corresponding to the position where the new symbol w i [|w i | − j − 1] (or w i [|w i |] = $ i ) must be inserted into B j (z) and the corresponding value must be inserted into L j (z). The new symbol w i [|w i | − j − 1] precedes the j-suffix w i [|w i |−j, |w i |]. Such a suffix is obtained by concatenating the symbol c z = w i [|w i | − j] with its (j − 1)-suffix. Let us suppose that the symbol c z associated with (j − 1)-suffix starting with the symbol c v is in position t in B j−1 (v).
Hence, t is also the position of the (j − 1)-suffix in the lexicographic order among the suffixes of S of length at most (j − 1) starting with c v .
In order to find L j (z)[r], we have to distinguish two cases.
In the first case, the j-suffix is the smallest suffix (in the lexicographic order) of length at most j starting with the symbol c z , i.e. r = 1. This means that c z does not exist in any segment B j−1 (h), h = 0, . . . , v − 1 and there does not exist in and L j (z) [1] = 0. In the second case, there exists some suffix starting with c z of length at most j that is lexicographically smaller than the j-suffix, i.e. r > 1. Recall that L j (z) [r] represents the length of the longest common prefix between the j-suffix and a fsuffix of a string w q (for some 0 ≤ q ≤ m− 1) , with f ≤ j, starting with the symbol c z , that immediately precedes the j-suffix in the lexicographic order.
If the longest common prefix between the f -suffix of w q and the j- and r > 1, so L j (z)[r] = 1. This means that the symbol c z associated to the (f − 1)-suffix of w q is not contained in the segment If the longest common prefix between the f -suffix of w q and the j-suffix of w i is longer than c z then both the (f − 1)-suffix of w q and the (j − 1)-suffix of w i start with the same symbol c v . So, we can suppose that the symbol associated with the (f − 1)-suffix of w q is at the position d 1 in the segment B j−1 (v). Remark that symbols in position d 1 and t are equal to c z . Then represents the length of the longest common prefix between the j-suffix and a g-suffix of a string w p (for some 0 ≤ p ≤ m − 1), with g ≤ j, starting with the symbol c z , that immediately follows the j-suffix in the lexicographic order.
If the longest common prefix between the j-suffix of w i and g-suffix of w p is c z then LSI v j−1 (c z , t) = L j−1 (v)[t], so L j (z)[r + 1] = 1. This means that the symbol c z associated to the (g − 1)-suffix of w p is not contained in the segment B j−1 (v) (but it is contained in some B j−1 (h), h = v + 1, . . . , σ).
If the longest common prefix between the j-suffix of w i and the g-suffix of w p is longer than c z then both the (j − 1)-suffix of w i and the (g − 1)-suffix of w p start with the same symbol c v . So, we can suppose that the symbol associated with the (g − 1)-suffix of w p is at the position d 2 in the segment B j−1 (v). Remark that symbols in position r and d 2 are equal to c z . Then LSI v j−1 (c z , t) = L j−1 (v)(t, d 2 ] and L j (z)[r + 1] = min(LSI v j−1 (c z , t)) + 1. Note that the position r + 1 in L j (z) does not exist when the the j-suffix is the greatest suffix (in the lexicographic order) of length at most j starting with the symbol c z . The suffix that, eventually, immediately follows the j-suffix in lexicographic order is involved in a segment B j (h), for some h = z + 1, . . . , σ, hence such suffix starts with a symbol greater than c z .
A consequence of the theorem is that the segments B j and L j can be constructed sequentially and stored in external files. This fact will be used in the next section.

Lightweight implementation via sequential scans
Based on the strategy described in the previous section, here we propose an algorithm (named extLCP) that simultaneously computes the EBWT and the LCP of a collection of strings S. Memory use is minimized by reading data sequentially from files held in external memory: only a small proportion of the symbols of S need to be held in internal memory. We could also add the computation of the generalized suffix array of S without adding further data structures.
Note that, at each iteration j, both the segments B j and L j , initially empty, are stored in different external files that replace the files used in the previous iteration. Consequently, both ebwt j (S) and lcp j (S) are updated accordingly.
6.1. Basic idea of the algorithm. The main part of the algorithm extLCP consists of K consecutive iterations. At iteration j, we consider all the j-suffixes of S and simulate their insertion in the sorted suffixes list. For each symbol associated with the j-suffix that we have to insert at the position r into B j (h), we also have to insert the new values related to the longest common prefix at position r and r + 1 into L j (h), where c h is the first symbol of the considered j-suffix.
In order to compute ebwt j (S) and lcp j (S), the algorithm extLCP needs to hold six arrays of m integers in internal memory. Four of these (P , Q, N and U ) are as employed by the algorithm BCR (see Section 4) and further two arrays (C and S) are needed to compute and update the values of the longest common prefixes at each iteration. As for P , Q, N and U arrays (see Section 4), for ease of presentation we denote by C j and S j the arrays at the j-th iteration. They contain exactly one integer for each string, i.e. they use O(m log K) bits of workspace. They are sequentially computed by using other auxiliary data structures described in Subsection 6.2.
More formally, if i = N j [q] then C j [q] stores the length of the longest common prefix between the j-suffix of w i and the previous suffix (if it exists) in the list of sorted suffixes with respect to the lexicographic order of all the suffixes of S of length at most j, whereas S j [q] contains the length of the longest common prefix between the j-suffix of w i and the successive suffix τ in the list of sorted suffixes (if it exists). Such values will be computed at the iteration j − 1 according to Theorem 5.1. Note that S j [q] is used when the suffix τ exists and P j [q] + 1 = P j [q + 1]. Example 6.1 provides an idea of the running of the algorithm extLCP and shows how the involved arrays are used. Example 6.1. Figure 3 illustrates an execution of our method on the collection S = {w 0 , w 1 } = {AAT ACACT GT ACCAAC$ 0 , GAACAGAAAGCT C$ 1 } at two consecutive iterations. Note that we have appended different end-marker to each string ($ 0 and $ 1 , respectively). In particular, we suppose that at the iteration j = 12, we have computed the arrays P 12 , Q 12 and N 12 . We recall that we have computed to terminate the iteration 12. We obtain that LCI 1 12 (G, 3) corresponds to the range (1, 3] in L 12 (1). So the minimum value is 2 and the value in lcp 13 (S) associated with the 13-suffix τ (i.e. τ = GAACAGAAAGCT C$ 1 ) of w 1 is 3 (this value is stored into C 13 ). Now we have to compute the value of the longest common prefix between the suffix τ and the suffix that immediately follows τ in the lexicographic order and to store this value into S 13 . Since the symbol G does not appear in the range L 12 (1)(3, 10], it means that there are not suffixes starting with GA lexicographically greater than τ , so such value will be less than 2. Because the symbol G does appear at least once in B 12 (h) with h > 1 (in this case B 12 (2)[7] = G), it means that there exists at least a suffix starting with G lexicographically greater than τ (in this case the suffix GCT C$ 1 ). So the value in lcp j+1 (S) of the suffix that follows τ must be updated to 1 (i.e. we store this value in S 13 ). Similarly, since LCI 2 12 (A, 4) = L 12 (2) (1,4], the minimum value is 1 and so the value in lcp 13 (S) for the 13-suffix ω (i.e. ω = ACACT GT ACCAAC$ 0 ) of w 0 is 2 (i.e. we insert this value in C 13 ). Moreover, LSI 2 12 (A, 4) = L 12 (2) (4,5], the minimum value is 2 and hence the value in lcp 13 (S) of the suffix that follows ω must be updated to 3 (i.e. we insert this value in S 13 ).
At the iteration j = 13, we compute the arrays P 13 , Q 13 and N 13 , whereas the arrays C 13 and S 13   Each iteration j > 0 can be divided into two consecutive phases. During the first phase we only read the segments B j−1 in order to compute the arrays P j , Q j N j and U j . Then we sort Q j , P j , N j , C j , S j , where the first and the second keys of the sorting are the values in Q j and P j respectively. We omit the description of the first phase, because it can be found in Section 4, so we focus on the second phase.
In the second phase, the segments B j−1 and L j−1 are read once sequentially both for the construction of new segments B j and L j and for the computation of the arrays C j+1 and S j+1 , as they will be used in the next iteration. Moreover, the computation of the segments L j is performed by using the arrays C j and S j constructed during the previous step.
Since the identical elements in Q j are consecutive, we open the pair files B j−1 (h) and L j−1 (h) (for h = 0, . . . , σ) at most once. Each of these files can be sequentially read, because the positions in P j are sorted in according with Q j .
In the sequel, we focus on a single segment B j (z), for z = 1, . . . , σ, by assuming that Q j [p] = z for each l ≤ p ≤ l ′ , with l ≥ 1 and l ′ ≤ m, i.e. we are considering the elements in B j (z) and L j (z) associated with the suffixes starting with c z . Since P j [l] < . . . < P j [l ′ ], we can sequentially build B j (z) and L j (z) by copying the old values from B j−1 (z) and L j−1 (z) respectively, by inserting each new symbol In the sequel, we say that the left extreme is the opening position and the right extreme is the closing position of some symbol x. A pseudo-code that realizes the second phase of a generic iteration j > 0 of the algorithm extLCP can be found in Figure 4 and it uses the following additional arrays of σ elements, defined as follows: • isMinLCIop and isMinLSIop are arrays of flags, where isMinLCIop [v] and isMinLSIop [v] indicate the opening of the LCI z j and LSI z j associated with c v , respectively.
• minLCI and minLSI are arrays of integers, where minLCI[v] and minLSI [v] store the minimum value among the values in L j (z) from the opening position of the intervals associated with c v and the current position s. • minLSInseq is an array of integers, where minLSInseq[v] contains the index g of the position in which minLSI[v] will be stored in S j+1 . This is useful, because when we close the LSI z j associated with the symbol c v at some position s, we have to store minLSI [v] in some position g of S j+1 , such that P j [g] < s. The routines insertNewSymbol, updateLCI and updateLSI are described in Figure 5 and Figure 6. When at any position s in B j (z) we insert the symbol c v (from U j or from B j−1 (z)), then the position s is assumed to be: • the opening position of LCI z j (c v , y), if another occurrence of c v will be inserted, as new symbol, at some next position y = P j [f ], for some f such Figure 6. Routines for computing the minimum values in LCI z j and LSI z j .
that P j [f ] > s. So, we set isMinLCIop[v] = 1 and the value minLCI[v] is updated (see updateLCI(x) routine in Figure 6). • the closing position of LSI z j (c v , P j [g]), if another occurrence of c v has been inserted, as new symbol, at some previous position P j [g], for some g such that P j [g] < s. In this case, LSI z j (c v , P j [g]) = L j (z)(P j [g], s] (because minLSInseq[v] = g) and we set S j+1 [g] = min(LSI z j (c v , P j [g])) + 1 according to Theorem 5.1. Such a minimum value is stored into minLSI [v]. We set isMinLSIop[v] = 0 (see updateLSI(x) routine in Figure 6). When B j (z) is entirely built, the closing position of some LSI z j (c v , P j [g]) for some l ≤ g ≤ l ′ could remain not found. So, we could have some value in isMinLSIop equal to 1. For each c v such that isMinLSIop[v] = 1, the last occurrence of c v appears at position P j [g] (we recall that we have set minLSInseq[v] to g). In this case, LSI z j (c v , P j [g]) = L j (z)[P j [g]] and we set S j+1 [g] = 1 according to Theorem 5.1 (see the for loop in Figure 4).
One can verify that these steps work in a sequential way. Moreover, one can deduce that, while the same segment is considered, for each symbol c v ∈ Σ at most one LCI z j (c v , P j [f ]) for some l ≤ f ≤ l ′ and at most one LSI z j (c v , P j [g]) for some l ≤ g ≤ l ′ will have not their closing position.

Complexity of the algorithm extLCP
The complexity of extLCP algorithm depends mainly on K, i.e. the length of the longest string in the collection, because the algorithm works in K passes and for each step it needs to build B j (z) and L j (z) from B j−1 (z) and L j−1 (z) for each z = 0, . . . , σ. At each step j, the used internal memory depends on the number of strings in S of length greater than j. Such a value is upper bounded by m for each step. Note that, at each step, the total size of σ + 1 files containing the partial ebwt(S) is increased by at most m symbols. Analogously the total size of σ + 1 files containing the partial lcp(S) is increased by at most m values. So, the used disk space mainly depends on m for each step.
Since our algorithm accesses disk data only by sequential scans, we analyze it counting the number of disk passes as in the standard external memory model (see [27]). We denote by B the disk block size and we assume that both the RAM size and B are measured in units of Θ(log N )-bit words.
From the size of the data structures and from the description of the phases of the extLCP algorithm given in previous sections, we can state the following theorem. Proof. At each iteration j, the main data structures used by extLCP are U j of size m log σ bits, Q j of size at most m log σ bits, N j of size at most m log m bits, P j of size at most m log N bits, C j and S j of size at most 2m log K bits. Moreover, we need σ 2 log N bits for computing the number of occurrences of each symbol in each segment B j (h), for h = 1, . . . , σ. The additional arrays take O(σ log m) (assuming that K ≪ m). So, the workspace is O((m + σ 2 ) log N ) bits. As we sort Q j , P j , N j , C j , S j where the first and the second keys of the sorting are the values in Q j and P j , we need O(mj + sort(m)) CPU time, where sort(m) is the time taken to sort m integers in internal memory. The total CPU time is O(K(N + sort(m)). We observe that the input is split into K files, so that each file contains m symbols, one for each string. At each step, such symbols will read and store into U j array. The amount of I/O operations for handling the input file is O( N B log σ N ). During the first phase of each iteration j > 0, we need to read at most m(j − 1)(log σ + log K) bits for computing Q j , P j , N j arrays. In the second phase we need to read at most m(j − 1)(log σ + log K) bits from B j−1 and L j−1 segments in order to obtain C j+1 , S j+1 arrays and B j and L j segments by writing mj(log σ + log K) bits of disk space. The disk I/O operations for obtaining the output of each iteration is O mj The internal memory of extLCP can be reduced significantly by observing that rather than storing P , Q, N and C to internal memory, they could reside on disk because these arrays are sequentially processed. In this way, the algorithm needs to store in internal memory the array S, isMinLCIop, isMinLSIop, minLCI, minLSI and minLSInseq that require random accesses. Clearly in this case one should use an external memory sorting algorithm.
Remark 7.2. It is easy to verify that the algorithm extLCP can be also adapted for computing the generalized suffix array of S together the EBWT. Such further computation increases the number of I/O operations, but one do not need of further data structures in internal memory. Note that, for each i, GSA[q] = (t, j) is the pair corresponding to the q-th smallest suffix of the strings in S, i.e. to the suffix w j [t, |w j |]. Since we know at each iteration of extLCP the values t and j, because they are store into P j and N j arrays respectively, it is enough to modify the extLCP code by adding the instructions for storing such values in external memory.

Computational experiments and discussion
Our algorithm extLCP represents the first lightweight method that simultaneously computes, via sequential scans, the LCP array and the EBWT of a vast collection of strings.
We developed an implementation of the algorithm described in Section 6, which is available upon request from the authors 2 .
Our primary goal has been to analyze the additional overhead in runtime and memory consumption of simultaneously computing both EBWT and LCP via extLCP compared with the cost of using BCR (cf. [3]) to compute only the EBWT. For this goal, we used a publicly available collection of human genome strings from the Sequence Read Archive [9] at ftp://ftp.sra.ebi.ac.uk/vol1/ERA015/ERA015743/srf/ and created subsets containing 100, 200 and 800 million reads, each read being 100 bases in length on the alphabet {A, C, G, T, N }.  Table 1. The input string collections were generated on an Illumina GAIIx sequencer, all strings are 100 bases long. Size is the input size in GiB, wall clock time-the amount of time that elapsed from the start to the completion of the instance-is given as microseconds per input base, and memory denotes the maximal amount of memory (in GiB) used during execution. The efficiency column states the CPU efficiency values, i.e. the proportion of time for which the CPU was occupied and not waiting for I/O operations to finish, as taken from the output of the /usr/bin/time command. Table 1 shows the results for the instances that we created. We show an increasing in runtime since extLCP writes the values of LCP after that the symbols in EBWT are written, so it effectively increases the I/O operations. So, a time optimization could be obtained if we read/write at the same time both the elements in EBWT and LCP by using two different disks. The tests on 100 and 200 million reads instances were done on the same machine, having 16GiB of memory and two quad-core Intel Xeon E5450 3.0GHz processors. The collection of 800 million reads was processed on a machine with 64GiB of RAM and four quad-core Intel Xeon E7330 2.4GHz processors. On both machines, only a single core was used for the computation.
For the experiments in the following, we use single identical nodes in a cluster of Dell PowerEdge M600 servers, each equipped with 8 Intel(R) Xeon(R) CPU X5460 3.16GHz processors 6144KiB cache and 32GiB of RAM, not used in exclusive mode.
As pointed out in Section 4, extLCP is not restricted to work on collections of strings of equal length. To examine the behaviour of our algorithm on datasets with reads of varying length and number, we have created five different datasets of DNA strings from datasets in Table 1: A 0 is composed of 4 million strings of 100bp, A 1 − A 4 have been obtained rearranging the characters in A 0 so that they form strings whose distribution of lengths is specified in Table 2.
strings strings length number  Table 2. Lengths distribution of the strings in the datasets A 0 , A 1 , A 2 , A 3 , A 4 containing DNA strings and in the datasets P 0 , P 1 , P 2 , P 3 , P 4 containing protein strings.
Datasets A 0 , A 3 and A 4 have the same number of characters and the same number of strings but different string lengths distribution and rising maximum string length. As can be observed from the first three rows of Table 3, a sensible increase in wall clock time is registered for greater maximum string length. This is consistent with the theoretical results on complexity described in Section 7 where we show that the total number of I/O operations depends on the maximal length of the strings. As expected, a corresponding decrease in CPU efficiency values is registered, while greater use of memory is due to an increment of the datatype size in C and S arrays. Datasets A 4 , A 1 and A 2 have approximately the same number of characters and exactly the same maximal string length but decreasing number of strings. The corresponding rows of Table 3 show a decrement in wall clock time and memory usage. Furthermore an increment of CPU efficiency values is registered, due to the fact that the number of strings in the collection affects I/O operations more than effective EBWT and LCP computation.
We have also tested the extLCP performances on datasets containing strings over alphabets with more than five symbols. For this purpose we have created, from  Table 3. The amount of time that elapsed from the start to the completion of the instance is given as microseconds per input base is the wall clock time. The efficiency column provides the CPU efficiency values, i.e. the proportion of time for which the CPU was occupied and not waiting for I/O operations to finish, as taken from the output of the /usr/bin/time command. The column mem is the maximum value returned by the RES field of top command representing the amount of used internal memory in MiB. The column aux disk space denotes the maximum amount of external memory in GiB used to accommodate all the working files, tot disk space includes the output files too.
UniProtKB/TrEMBL dataset 3 , five different datasets (denoted by P 0 − P 4 ) whose protein strings have the same lengths distribution of the strings in the datasets A 0 − A 4 (see Table 2). The results of Table 3 show that the same considerations made for the DNA string collections A 0 − A 4 also hold for the datasets P 0 − P 4 . Moreover, we notice a general increase in CPU efficiency and wall clock time, while the used internal and external space remains unchanged. Altogether, one can verify that both the varying string length and the size of the alphabet could have a significant influence on the extLCP performance.
Note that we have not compared our method with other lightweight approaches that compute the LCP array of a single string (see for instance [6]), because their implementations cannot be promptly tested on string collections, nor on a generic string produced as output by EBWT of a string collection. In general, an entirely like-for-like comparison would imply the concatenation of the strings of the collection by different end-markers. However, for our knowledge, the existing implementations do not support the many millions of distinct end-markers our test collections would require. An alternative is to concatenate each of strings with the same end-marker. The use of the same end-marker without additional information leads to values in the LCP exceeding the lengths of the strings and depending on the order in which the strings are concatenated.  Table 4. The column size is the input size in GiB; the columns strings length and mean length correspond to the range of lengths and the average length in the dataset respectively; the columns LCP mean and LCP max are the average and maximum values of the longest common prefix of the strings in the dataset.
At the best of our knowledge, the only lightweight tool capable of working on large string collections in lightweight memory is eGSA [18]. The algorithm eGSA 4 takes in input a string collection S and returns the generalized suffix and LCP arrays of S. Note that, as well as extLCP, eGSA uses (implicit) distinct end-markers.
For large collections eGSA works in two phases [17,16]. In the first phase, the collection S is partitioned into different subsets of equal dimension and each subset is treated as a single string, that is the concatenation of the strings of the subset. Then eGSA produces, in internal memory, the SA and the LCP for each concatenated string and writes them to external memory. In the second phase, eGSA merges the arrays previously computed to obtain suffix and LCP arrays of the entire collection S. The number of created subsets is related to the amount of internal memory available for the computation, such a value can be set through the parameter MEMLIMIT.
We point out that using subsets of equal dimension in the first phase implies that the performance of eGSA does not deteriorate when a large collection with strings of varying length is considered.
Moreover, while extLCP computes the LCP array via the EBWT and is parameter-free, eGSA builds the LCP array by producing the GSA for string collections and needs that MEMLIMIT (or the number of subsets) is set. The overall approach adopted by eGSA could be thought as "orthogonal" to that of extLCP, in the sense that extLCP builds the LCP array, incrementally, by inserting at each step a "slice" of (at most) m characters from the strings in the collection, while eGSA proceeds by building separately (in internal memory) the LCP array of each string (or subset) and then opportunely merging them. Actually, we were not able to compare the performance of extLCP and eGSA on very large datasets. Indeed, current implementation of eGSA shows an overhead in disk space usage that prevented us to make tests on very large string collections. For instance, for a collection of 400 millions of strings of length 100 eGSA needs 40 bytes per symbol, whereas extLCP uses 5 byte per symbol. So, the datasets in Table 1 are too demanding to be managed. However, experiments have been conducted using smaller real collections (see Table 4): the dataset D 1 (long jump library of Human Chromosome 14 5 ), the datasets D 2 and D 3 (Human Genome  Table 5. The wall clock time-the amount of time that elapsed from the start to the completion of the instance-is given as microseconds per input base. The efficiency column provides the CPU efficiency values, i.e. the proportion of time for which the CPU was occupied and not waiting for I/O operations to finish, as taken from the output of the /usr/bin/time command. The column % mem is the maximum value returned by %MEM field of top command representing the percentage of used internal memory. The column aux disk space denotes the maximum amount of external memory in GiB used to accommodate all the working files, tot disk space includes the output files too. sequences from Sequence Read Archive 6 ) are collections of DNA strings on the alphabet {A, C, G, T, N }. Moreover, in order to evaluate extLCP and eGSA when the dataset contains strings over alphabets with more than five symbols, we have created a new collection of proteins (called D 4 ), obtained from UniProtKB/TrEMBL dataset, by truncating the strings longer than 100 aminoacids, so that it has the same maximum string length as the datasets D 1 , D 2 , D 3 . For each dataset, extLCP is compared with eGSA using two different values of MEMLIMIT: an higher value is comparable with the size of memory required from the data structures used by extLCP; a lower value have been chosen to produce a number of subsets less than 1, 024. We have not set other parameters of eGSA code.
The results of our tests are described in Table 5. The first three experiments, in which we consider strings over the same alphabet, show that the bigger the string collection is, the better is the overall performance of extLCP with respect to eGSA. In fact, the wall clock time values of extLCP and eGSA become comparable, but the total disk usage for eGSA significantly increases. However, the fourth experiment in Table 5 shows that the alphabet size is a parameter that seems to have a significant impact on the wall clock time value for extLCP rather than for eGSA.
Finally, we observe that, as pointed in Remark 7.2, extLCP could produce as additional output the generalized suffix array. Although our algorithm is not optimized for this purpose, our tests have shown that the wall clock time is about twice because the I/O operations are doubled. However, the computation of the GSA does not produce an increase in the internal memory.

Conclusions
In this paper, we proposed extLCP which is a lightweight algorithm to construct, at the same time, the LCP array and the EBWT of a collection of strings.
Actually, the LCP and BWT are two of the three data structures needed to build a compressed suffix tree (CST) [25] of a string. The strategy proposed in this paper could enable the lightweight construction of CSTs of string collections for comparing, indexing and assembling vast datasets of strings when memory is the main bottleneck.
The algorithm extLCP is a tool parameter-free designed for very large collection, indeed Table 1 shows that it also works with collections of 74GiB and Tables 3 and  5 show that it can be also used on large collections with varying string lengths.
The experimental results show that our algorithm is a competitive tool for the lightweight simultaneous computation of LCP and EBWT on string collections.
Our current prototype can be further optimized in terms of memory by performing the sorting step in external memory. Further saving of the working space could be obtained if we embody our strategy in BCRext or BCRext++ (see [3]). These methods, although slower than BCR, need to store only a constant and (for the DNA alphabet) negligibly small number of integers in RAM regardless of the size of the input data.