A recent paradigm shift in bioinformatics from a single reference genome to a pangenome brought with it several graph structures. These graph structures must implement operations, such as eficient construction from multiple genomes and read mapping. Read mapping is a well-studied problem in sequential data, and, together with data structures such as sufix array and Burrows-Wheeler transform, allows for eficient computation. Attempts to achieve comparatively high performance on graphs bring many complications since the common data structures on strings are not easily obtainable for graphs. In this work, we introduce prefix-free graphs, a novel pangenomic data structure; we show how to construct them and how to use them to obtain well-known data structures from stringology in sublinear space, allowing for many eficient operations on pangenomes.
which makes them computationally expensive. Here
we present a new class of graphs, prefix-free graphs,
The term pangenome was first used by Tettelin which are orders of magnitude faster to construct.
et al. [
Several approaches to pangenomic graph construc- sive knowledge. Using the sufix array from a prefixtion exist, such as variation graphs [7, 8, 9], cactus free graph, we can obtain several stringological data graphs [10, 11], and Wheeler graphs [12, 13]. Most structures which are not easily obtainable for graphs. of these approaches require an initial local alignment This feature of prefix-free graphs was implicitly of similar regions or a multiple sequence alignment, demonstrated in several articles [16, 17, 13, 18], where the authors used similar techniques as presented here. We think that explicitly defining and framing the prefix-free graph as a standalone pangenomic data structure can bring several benefits: • reduction in the complexity of the presentation of several space-eficient algorithms • support of theoretical research by clearly de
limiting the relevant terms • improved focus on the optimization of algo
rithms related to prefix-free graphs • enabling bringing prefix-free graphs closer to the biological data
In this work, we define prefix-free graphs and show how they can be constructed from the pangenome in its textual representation. Furthermore, we show in detail how prefix-free graphs can be used to generate the sufix array of a pangenome in sublinear space T using the automaton. Each time we encounter and linear time. Finally, we implement the pre- a trigger word, we recognize a new segment from sented algorithms as two binaries for easy construc- the start of the previous trigger word to the end of tion of prefix-free graphs from a set of sequences the current trigger word. If the segment’s sequence in FASTA format and from a pangenomic graph in was not yet observed during the scan, we add it to GFA format. Furthermore, we implement the rust the set of segments and assign a unique ID. Each library for working with prefix-free graphs. This time we append the corresponding ID to the path library contains an iterator, which can be directly representing the original sequence. used to generate the sufix array in sublinear space. Two special cases happen during the sequence scan, one at the beginning, when no previous trigger word was encountered, and another at the end, 2. Prefix-free graphs when the last k characters are sentinels. These are addressed by simply starting the first segment at The idea of prefix-free graphs is inspired by a the start of a sequence and ending the last at the technique used in the tool rsync named Context- sequence end.
Triggered Piecewise Hashing (CTPH) [19]. CTPH Notably, the adjacent segments overlap by exactly uses a rolling hash to partition a string into sub- k characters. Furthermore, trigger words occur strings such that long repeated substrings are par- only at the beginnings or ends of segments because titioned the same way. These substrings are then any occurrence of a trigger word in the middle of hashed with a traditional hash function and stored a segment would break it into two. This feature as a string signature. The signatures of several files and the choice of sentinels outside the sequence’s are then compared to determine changes. alphabet guarantee that the set of segments is prefix
In prefix-free graphs, we partition each sequence free. A sketch of the proof is shown in Figure 1. of a given pangenome into segments. These seg- After the previous steps, the set of segments and ments form nodes of the prefix-free graph, and the list of paths already represent a prefix-free graph. their adjacencies in the original sequences constitute However, to simplify the usage of prefix-free graphs, edges. The sequences are represented as paths in we recommend normalizing them. During the northe graph. malization, we sort segments lexicographically and
The segments have two essential characteristics change their IDs to correspond to the lexicographimaking them a good choice for nodes of a pange- cal ranks. This relabeling is then also propagated nomic graph. Similarly to CTPH, long repeated to paths accordingly. All construction steps can be sequences will be partitioned the same way. Fur- performed in space proportional to the sum of segthermore, no segment is a prefix of another, making ment lengths and the sum of path lengths, which is a set of segments prefix-free. The second charac- expected to be significantly smaller than the length teristic is crucial for connecting prefix-free graphs of the pangenome. and sufix arrays, as will be presented in the next To illustrate the entire procedure, consider a set of section. sequences {CACGTACT, CACACT, CACGACT} and a set of
To create a prefix-free graph from a given trigger words T = {AC, CG}. After the partitioning,
pangenome, we define a set of trigger words T , we obtain a set of segments with IDs {0:CAC, 1:ACG,
where each trigger word is a string of length k. For 2:CGTAC, 3:ACT.., 4:ACAC, 5:CGAC} and a list of
this set T , we build an Aho-Corasick automaton paths [[
pangenome from a prefix-free graph, we need to create several data structures. First, we concate
From this representation, original sequences of a nate all the segments into a single string using a pangenome can be reconstructed by expanding the separator # and append a sentinel $. We will call segment IDs in a particular path, ignoring the last this concatenation segment join. An example of a k characters of each segment. segment join is in Figure 3.
Next, we calculate the segment join’s sufix array and the longest common prefix array [ 26, 27]. For 3. Sufix array construction both of these arrays, there exist algorithms with linear time and space complexity which we can use.
A sufix array is a permutation of string positions We note that these linear complexities are proporwhich lexicographically sorts the sufixes of the tional to the length of the segment join, which is string starting at that position. For a set of strings, usually much smaller than the original pangenome. we consider a concatenation of the strings as shown Next, for each sufix of the segment join, we need in Figure 4. In the following text, we will refer to to calculate the corresponding segment ID and pothis concatenation as the pangenome. sition values. The value segment ID represents in
The sufix array is an influential data structure what segment the current sufix starts, and the value with many applications in eficient string algorithms segment position represents at what position in that solving problems such as exact pattern matching, particular segment the current sufix starts. These repeat finding, maximum exact match (MEM) find- arrays can be computed using an inverse permutaing, document retrieval and many more. There exist tion of a sufix array ISA (Equation 1) of a segment several algorithms for sufix array construction in join in linear time. linear time [21, 22, 23] with several practical implementations [24, 25]. Despite their linear time ISA[SA[i]] = i (1) complexity, these algorithms become bottlenecks in some applications because of their linear space complexity. This observation is especially relevant in pangenomics, where the datasets often do not fit in the computer memory.
Here, we show another crucial advantage of prefixfree graphs. Although they do not ofer any improvement of theoretical guarantees in the worst case, in practice, they often represent the pangenome in a substantially smaller space and allow us to generate the sufix array values one by one, possibly using the values directly in subsequential computation or storing them in compressed form. This iteration can be done without ever expanding the pangenome to its full textual representation in space proportional to the sum of segment lengths and the sum of path lengths.
join of our running example from Figure 3. Each position of the join can be assigned a segment ID and a position in the current segment by linearly scanning the segment join and incrementing the ID and position accordingly. Then, applying the ISA to these arrays changes the order of computed values in correspondence to the sorted sufixes. The resulting sufix array (SA), longest common prefix array (LCP), segment ID array (ID) and segment positions array (pos) are stored in a sufix table as shown in Table 1.
In the sufix table, one row can represent multiple positions of the pangenome. To identify these positions, we store some additional information in a segment table. For each segment of the prefixfree graph, we store its length, starting positions in the pangenome and ranks of the right contexts of example of a path join for our running example is these positions. To calculate the starting positions in Figure 4. Then, starting positions can be calcuand the ranks of the right contexts, we use a path lated by cumulatively summing the lengths of the join. Similarly to segment join, a path join is a segments in path join and subtracting the overlaps. concatenation with delimiters # and sentinel $, but The computation of ranks is more involved. It now constructed by concatenating the paths. An uses the normalized form of prefix-free graphs since it relies on a lexicographically smaller ID in a path representing a lexicographically smaller segment.
We construct the sufix array of the path join and ifnd its inverse permutation ISA. ISA gives us ranks for each position in the path join. To determine the rank of the right context for position i, we take the value of ISA[i + 1]. Finally, we store the starts and ranks sorted by the rank values in the segment table as shown in Table 2. • a singleton block with segment sufix occur
ring several times in the pangenome • a non-singleton block
sufix array value. Given the row index i, this value can be calculated with Equation 2.
SA value = starts[ID[i]] + pos[i] (2)
3.2. Iteration the first row yielding a SA value. Its segment ID is 0, and from Table 2, we see only one occurrence of With the previous tables stored in memory, we have segment 0 with starting position 9 in the pangenome. all the necessary ingredients to generate the sufix The ofset from the start of a segment pos[13] is 0. array value by value. Summing these two values, we get the first value
Each row in the sufix table represents a sufix of a of a sufix array 9 + 0 = 9 corresponding to the particular segment. There are four cases of what the lexicographically smallest sufix of a pangenome ifrst position of these sufixes can represent within P [9..] = ACACT. the segments: The second, slightly more complex case is a single• the sentinel $ ton block representing a segment sufix with several • a separator # occurrences in the pangenome. In this case, we must • a position within the last k characters of a report as many sufix array values as the number of segment occurrences. Because the starting segment positions • a position outside the last k characters, sep- in the segment table are sorted based on their right arator and sentinel context rank, we can iterate through these starting Since the pangenome has no corresponding posi- positions and apply Equation 2 to each of them. tion for the sentinel or separator characters, we can As an example, consider the row 14 in Table skip the first rows representing them. 1. This sufix occurs twice in the pangenome in
In the third case, the position is inside the trigger segments starting at positions 15 and 1. Since the word or the sentinels appended during the graph ofset from the start of a segment pos[14] is 0, we construction. The positions inside the trigger words report a sufix array values 15+0 = 15 and 1+0 = 1, are represented twice in the sufix table, once at the corresponding to the sufixes P [15..] = ACGACT and end of a segment and a second time at the beginning P [1..] = ACGTACT. of the following segment in the pangenome. These In the last case, we have a non-singleton block ending positions can violate the prefix-free prop- representing sufixes of several segments, possibly erty of the segment sufixes and, therefore, can be with multiple occurrences. These sufixes represent sorted incorrectly. Skipping through these positions identical substrings in the pangenome. Here, we ensures the prefix-free property for the rest of the report a sufix array value for each of the substrings. sufixes and also avoids double reporting. Therefore, To identify the first value, we must find the starting if the length of a current segment sufix is smaller or position with the smallest right context rank. Beequal to the size of the trigger words k, we skip the cause the ranks are sorted, this procedure is similar row as in the previous cases. This choice also plays to the merging phase of a merge sort. Therefore, to nicely with the previous choice of appending k sen- iterate through all sufix array values in the block, tinels during the construction of a prefix-free graph, we always identify the segment start with the next as these positions will not get reported either. smallest right context rank and apply Equation 2
Finally, in the last case, we report the sufix array to this segment start. values. The sufix table can be partitioned into As an example, consider the block of rows [20..21] blocks of the same segment sufixes. For example, in Table 1. The relevant segment IDs are 0 and consider rows 20 and 21, which form a single block. 3, with segment starts at positions 9, 8, 14 and All other blocks in the running example consist of 0. The right context ranks from smallest to highsingle rows; therefore, we call them singletons. est are 4, 5, 6, 9 with the corresponding segment This partitioning leads to three cases: starts 8, 14, 0, 9. Applying Equation 2 to these segment starts yields a sufix array values 8 + 0 = 8, 14 + 0 = 14, 0 + 0 = 0 and 9 + 1 = 10, represent• a singleton block with segment sufix occurring only once in the whole pangenome
We implemented the prefix-free graphs as a Rust package. This package contains two binary crates and one library crate. Binary crates are executables that serve the purpose of creating prefix-free graphs from FASTA and GFA formats. The binary crates are named fasta2pfg and gfa2pfg, and their usage is as follows: fasta2pfg -t triggers.txt < pangenome.fna > , → pfg.gfa gfa2pfg -t triggers.txt < pangenome.gfa > , → pfg.gfa
In addition to the sufix array value, the iterator provides the segment ID and position values. These are useful for the computation of several related data structures. Similarly to Boucher et al. [16], we can use them to compute the Burrows-Wheeler transform BWT[i] by storing the preceding characters of each segment and then reporting the character at position pos[i] − 1 of a segment id[i]. This computation allows for space-eficient construction of r-index [28] similar to MONI implementation [17], spaceeficient construction of Wheeler graphs similar in nature to the implementation in Goga and Baláž [13], and with the use of predecessor queries, spaceeficient construction of a tag array introduced in Goga et al. [18].
We compared the running time of the construc- leads to the need for new data structures to utilize tion algorithm with several pangenome construction all the contained information to its full potential. tools, namely with the PanGenome Graph Builder Pangenome graphs have been proposed as a solution [9], VG [8], and Minigraph-Cactus [11]. For com- capable of representing these datasets, exploiting parison, we used up to 256 SARS-CoV-2 sequences. the repetitiveness of a collection of related genomes These sequences have lengths of around 30 kbp and for eficient storage while preserving and sometimes high nucleotide similarity. We used stop codons even highlighting the underlying variations. How(TAA, TAG, and TGA) as trigger words to construct ever, this shift in the bioinformatics field from linear prefix-free graphs. The resulting running times are sequences to pangenomic graphs brings challenges. shown in Table 3. The enormous sizes of pangenomes make the con
The library crate provides an interface for working struction of pangenomic graphs nontrivial and often with prefix-free graphs, mainly an iterator of a sufix one of the major bottlenecks in the analysis. Curarray in sublinear space and linear time. The library rent tools rely on well-known computational steps, can be used from within the Rust programming such as all-against-all alignment [9], multiple selanguage as follows: quence alignment [11], or variant calling [8]. All these steps are computationally expensive and do let pfg = PFG::load("pfg.gfa"); not scale well with rapidly growing datasets. Furthermore, variant calling relies on a linear reference, for (i, (sa_i, id_i, pos_i)) in which can introduce a reference bias. , → pfg.iter().enumerate() { Here, we introduced the prefix-free graph as a println!("{}\t{}\t{}\t{}", i, sa_i, standalone data structure and showed how to build , → id_i, pos_i); it from currently available pangenomic datasets. } In comparison to the other tools, it ofers several crucial advantages. The construction of prefix-free graphs does not require any alignment or variant calling. Instead, it avoids these expensive steps by employing a set of trigger words, which split the sequences according to their contexts. The time complexity of this construction is linear with respect to the input size, and the space complexity is sublinear, proportional only to the size of the resulting data structure.
The novelty of prefix-free graphs brings several possible directions for future research, with the two main directions being the connection of prefixfree graphs to stringology and the choice of trigger words. In Chapter 3, we showed how to use prefixfree graphs to build a sufix array of a pangenome in sublinear space. Sufix arrays are at the core of many eficient and well-established string algorithms in bioinformatics, such as read mapping and pattern matching. We expect this connection will facilitate the development of similar algorithms on pangenomic datasets, supporting the paradigm shift.
This direction will require exploring which additional string algorithms can use prefix-free graphs to improve their applicability to these vast and repetitive datasets. An exciting attempt may be to map reads to prefix-free graphs in a similar fashion as the popular tool BWA [29] maps reads to a linear reference.
The second direction, the choice of trigger words, is, so far, mostly undiscovered territory. Trigger words ofer great flexibility in the construction process, so a better understanding of their selection is crucial for prefix-free graphs. From the computational perspective, a set of trigger words minimizing the size of the resulting graph is desirable since it would allow the analysis of larger datasets. On the other hand, the construction flexibility could be used to create prefix-free graphs based on biologically significant strings. In our experiments, we used stop codon sequences as the trigger words.
However, experiments with other motives, such as diferent binding sites, recombination hotspots, or repetitive elements, could illuminate the possibility of capturing particular biological phenomena with prefix-free graphs. Moreover, integrating the strand information in the graph construction may be beneficial for further reducing the graph size and capturing biological phenomena such as inversion.
This integration could be achieved by considering both the forward and reverse complement of the trigger words.
We believe the characteristics of prefix-free graphs are highly valueable and will have a significant impact when dealing with massive datasets, moving us closer to the ultimate goal of computational pangenomics.
ropean Union’s Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No 956229 (ALPACA) and by a grant from Slovak Research Grant Agency VEGA 1/0538/22. [6] Computational pan-genomics: status, promises (2022).
and challenges, Briefings in bioinformatics 19 [19] J. Kornblum, Identifying almost identical files (2018) 118–135. using context triggered piecewise hashing, Dig[7] D. M. Church, V. A. Schneider, K. M. Stein- ital investigation 3 (2006) 91–97. berg, M. C. Schatz, A. R. Quinlan, C.-S. Chin, [20] A. V. Aho, M. J. Corasick, Eficient string P. A. Kitts, B. Aken, G. T. Marth, M. M. matching: an aid to bibliographic search, ComHofman, et al., Extending reference assembly munications of the ACM 18 (1975) 333–340. models, Genome biology 16 (2015) 1–5. [21] J. Kärkkäinen, P. Sanders, Simple linear [8] E. Garrison, J. Sirén, A. M. Novak, G. Hickey, work sufix array construction, in: Automata, J. M. Eizenga, E. T. Dawson, W. Jones, Languages and Programming: 30th InternaS. Garg, C. Markello, M. F. Lin, et al., Varia- tional Colloquium, ICALP 2003 Eindhoven, tion graph toolkit improves read mapping by The Netherlands, June 30–July 4, 2003 Prorepresenting genetic variation in the reference, ceedings 30, Springer, 2003, pp. 943–955.
Nature biotechnology 36 (2018) 875–879. [22] G. Nong, S. Zhang, W. H. Chan, Linear sufix [9] E. Garrison, A. Guarracino, S. Heumos, F. Vil- array construction by almost pure inducedlani, Z. Bao, L. Tattini, J. Hagmann, S. Vor- sorting, in: 2009 data compression conference, brugg, S. Marco-Sola, C. Kubica, et al., Build- IEEE, 2009, pp. 193–202. ing pangenome graphs, bioRxiv (2023) 2023– [23] S. J. Puglisi, W. F. Smyth, A. H. Turpin, A tax04. onomy of sufix array construction algorithms, [10] B. Paten, M. Diekhans, D. Earl, J. S. John, acm Computing Surveys (CSUR) 39 (2007) J. Ma, B. Suh, D. Haussler, Cactus graphs 4–es. for genome comparisons, Journal of Computa- [24] F. A. Louza, S. Gog, G. P. Telles, Inducing tional Biology 18 (2011) 469–481. enhanced sufix arrays for string collections, [11] G. Hickey, J. Monlong, J. Ebler, A. M. Novak, Theoretical Computer Science 678 (2017) 22– J. M. Eizenga, Y. Gao, T. Marschall, H. Li, 39.
B. Paten, Pangenome graph construction from [25] Y. Mori, libdivsufsort, 2008. URL: https:// genome alignments with minigraph-cactus, Na- github.com/y-256/libdivsufsort, accessed on ture Biotechnology (2023) 1–11. 2023-06-23. [12] T. Gagie, G. Manzini, J. Sirén, Wheeler graphs: [26] G. M. Landau, T. Kasai, G. Lee, H. Arimura, A framework for bwt-based data structures, S. Arikawa, K. Park, Linear-time longestTheoretical computer science 698 (2017) 67– common-prefix computation in sufix arrays 78. and its applications, in: Combinatorial Pat[13] A. Goga, A. Baláž, Prefix-free parsing for tern Matching: 12th Annual Symposium, CPM building large tunnelled wheeler graphs, in: 2001 Jerusalem, Israel, July 1–4, 2001 Proceed22nd International Workshop on Algorithms in ings 12, Springer, 2001, pp. 181–192. Bioinformatics (WABI 2022), Schloss Dagstuhl- [27] G. Manzini, Two space saving tricks for linear Leibniz-Zentrum für Informatik, 2022. time lcp array computation, in: Scandina[14] M. Burrows, A block-sorting lossless data com- vian Workshop on Algorithm Theory, Springer, pression algorithm, SRC Research Report, 124 2004, pp. 372–383.
(1994). [28] T. Gagie, G. Navarro, N. Prezza, Fully func[15] P. Ferragina, G. Manzini, Opportunistic data tional sufix trees and optimal text searching in structures with applications, in: Proceedings bwt-runs bounded space, Journal of the ACM 41st annual symposium on foundations of com- (JACM) 67 (2020) 1–54.
puter science, IEEE, 2000, pp. 390–398. [29] H. Li, R. Durbin, Fast and accurate short [16] C. Boucher, T. Gagie, A. Kuhnle, B. Lang- read alignment with burrows–wheeler transmead, G. Manzini, T. Mun, Prefix-free parsing form, bioinformatics 25 (2009) 1754–1760. for building big bwts, Algorithms for Molecular
Biology 14 (2019) 1–15. [17] M. Rossi, M. Oliva, B. Langmead, T. Gagie,
C. Boucher, Moni: A pangenomic index for ifnding maximal exact matches, Journal of
Computational Biology 29 (2022) 169–187. [18] A. Goga, A. Baláž, A. Petescia, T. Gagie,
Maria: Multiple-alignment r-index with aggregation, arXiv preprint arXiv:2209.09218