cucapra / cucapra/pollen

Bring flattening to BED files

Open
#206 2 comments 0 reactions 0 assignees View on GitHub
Dominant language
Rust
Stars
45
Forks
3
PR merge metrics
No merged PRs in 30d

Description

Let's take the same "data structure flattening" approach as FlatGFA and apply it to [BED files](https://bedtools.readthedocs.io/en/latest/content/general-usage.html). This would serve a few goals:

* Possibly an easy win that biologists would care about.
* Further prove our general principles by applying it to a second data format.
* Help develop general principles and library functionality for developing flat versions of _arbitrary_ file formats!
* Dovetail nicely with FlatGFA itself, because some of the input/output datasets involved in odgi-style analyses are in fact BED files.
* Get closer to demonstrating the advantages of avoiding serialization/deserialization to text by composing odgi-style and bed tools-style operators.

As a high-level vision, let's take this shell stanza from [this odgi tutorial](https://odgi.readthedocs.io/en/latest/rst/tutorials/detect_complex_regions.html):

```sh
odgi depth -i chr8.pan.og -r chm13#chr8 | \
bedtools makewindows -b /dev/stdin -w 5000 > chm13.chr8.w5kbps.bed
odgi depth -i chr8.pan.og -b chm13.chr8.w5kbps.bed --threads 2 | \
bedtools sort > chr8.pan.depth.w5kbps.bed
```

…and port it to FlatGFA-augmented Python:

```py
chr8 = flatgfa.load("chr8.flatgfa")
windows = chr8.depth(path="chm13#chr8").makewindows(5000)
res = chr8.depth(ranges=windows).sort()
res.top(5).print()
```

The type of `chr8` in this fantasy example is `FlatGFA`. The type of `windows` and `res` is something like `FlatBED`. We will need a binary format for this, some interval-oriented operations on this format, and some Python bindings.

Here is a step-by-step plan, broken down into phases.

## 1: Explore

To understand what’s going on here:

- Follow [that odgi tutorial](https://odgi.readthedocs.io/en/latest/rst/tutorials/detect_complex_regions.html), producing a few BED files.
- Using [the BED file format spec](https://bedtools.readthedocs.io/en/latest/content/general-usage.html) and the examples produced above, think carefully about how much of the BED format we need to support to get started with, how much flexibility we want to build in, etc.

## 2: Parser and Printer

Let’s start by *just* writing a BED parser/printer. Let’s not worry about the in-memory data structure too much at first; we’ll just use something flat-ish (i.e., avoid pointers) and worry about the details later.

1. Make up a simple Rust data structure for representing a BED file. Maybe it’s mostly just a `Vec<(u32, u32, u32)>` or similar, where each line is a name reference, a start, and an end. To keep the data structure flat, avoid putting strings directly in the data structure and refer to them in a separate pool of bytes. Put this in a file called, like, `bed.rs` in the FlatGFA crate.
2. Write a parser. For ideas about how to write a fast parser for these tab-delimited file formats, borrow from [`gfaline.rs`](https://github.com/cucapra/pollen/blob/main/flatgfa/src/gfaline.rs). In addition to ideas, consider also calling some of the actual functions there: for example, use `parse_num` to parse text integers.
3. Write a printer. Probably by implementing the `Display` trait. For inspiration, see [`print.rs`](https://github.com/cucapra/pollen/blob/main/flatgfa/src/print.rs).
4. Make a little command-line wrapper that just parses and re-prints BED files.
5. Add some tests to check that this round-trip reproduces BED files exactly, byte-for-byte.
6. Open a PR! Declare success!

## 3: Bytes and Files

Now let's invent a binary file format: call it FlatBED. The goal is round-trip conversion between this binary format and text BED files.

This phase is very similar to an equivalent phase in the #204 plan. Namely, break out [the zerocopy crate](https://docs.rs/zerocopy/latest/zerocopy/) and start reading and writing to files. Please see #204's Phase 2 for an outline of the strategy to pursue here; the steps are essentially identical.

## 4: Replicate Bedtools

Let's implement some of the basic operators from [bedtools][] so we can compare against it.

The fantasy here is that we "win" against bedtools and/or bedtk in simple situations _even if we implement the simplest possible algorithms_ for each operation. Even without _any_ asymptotic complexity advantage, we have a significant I/O advantage by skipping serialization/deserialization. At least for operators where the amount of computation to be done is small, we hope to run faster than the baselines by just cutting out the I/O time.

1. Implement equivalents of `bedtools makewindows` and `bedtools sort`, which are the commands used in our "target" script above. They should work in our command-line tool developed in the previous phase. Use the simplest possible O(n)-time approaches to doing the actual interval operations.
2. Set up some tests to compare the output with the real bedtools results.
3. Open a PR and let's get this merged!
4. Compare performance against bedtools, and possibly also [bedtk][] where appropriate.
5. Consider implementing a few additional subcommands from bedtools: just pick the low-hanging fruit. Expose each in our little command-line tool. Open PRs for each of these as we finish them.

[bedtools]: https://bedtools.readthedocs.io/en/latest/
[bedtk]: https://github.com/lh3/bedtk

## 5: Make FlatGFA Produce FlatBED

Now let's make `fgfa` itself interact with these new flattened BED-like files.

1. Add an option to `fgfa depth` that produces FlatBED files instead of ordinary tab-separated text output. (This is worth a PR of its own. Also needs some tests.)
2. Add a `-b` option to `fgfa depth` so it can consume these FlatBED files. (Also do some testing and open a PR.)
3. Finally, now that we have all of these components, we can reproduce the shell snippet above using invocations of our tools. Benchmark it against odgi/bedtools or odgi/bedtk.

## 6: Python Bindings

Let’s add Python bindings for all of the above, so we can write the fantasy Python program above.

The upshot should be that we can now avoid writing files to disk (we already avoided serialization/deserialization). Benchmark this little Python script against *both* ourselves (using separate invocations of the `fgfa` command-line tool) *and* the odgi/bedtools baseline.

Contributor guide

No contributing guide indexed for this repository

Research direction

Start with the FlatGFA crate and inspect the mentioned gfaline.rs and print.rs files before implementing the Phase 2 parser/printer work in bed.rs. The initial milestone is a command-line wrapper whose BED parsing and printing round-trip matches the input byte-for-byte, with tests covering that result; later phases expand this into FlatBED, bedtools operations, and Python bindings.

Written by the indexing model from the issue text.

Assessment

Tech stack
rust
Domain
cli, data
Issue type
Feature
Difficulty
5/5
Estimated time
Over a week
Activity status
Stale
Clarity
Mostly clear
Newbie friendliness
35/100

Get new issues in your inbox

A short digest of beginner-friendly GitHub issues.