A bit-packed k-mer representation (and relevant utilities) for rust

Hi,

Actually library, store kmer in an array of u64 it could be nice to let user choose the type of array. I work with very short kmer, use and u64 to store it, would have a significant cost.

In this pull request I introduce a trait Encoder with 3 method encoding decoding and rev_comp, Kmer can use this Encoder method to perform encoding, decoding and reverse complement operation on array.

I also introduce 2 encoder, the first one Naive support any encoding type, the counter part is /slower/ operation. Xor10 operation is faster but support only one encoding (A -> 00, C -> 01, T -> 10, G -> 11).

We can add many other encoder, an encoder use smid instruction to convert 32 nucleotide in u64, @Daniel-Liu-c0deb0t create some algorithm to do this here nuc2bits. Usage of interface allow user to create Encoder fit perfectly to his needs.

In this pull request Kmer constructor accept u8 slice and an object implement Encoder, we can move this Encoder as Kmer type argument, but the type declaration became a little nightmare:

let kmer: Kmer<15, u16,  { word_for_k::<u16, 15>() }, Encoding::Naive::ACTG> = Kmer::new("ACTGAGAGAGACCAT");

This type complexity could by simplify by type aliasing and also when Rust get a more complete constant generic interface.

To simplify code base I also create trait Data to group all trait must be implemented by type use in array and add a method to_u8 required by Naive encoder, this trait and this method could be move or remove.

My implementation could probably be improve in many point, but I think general structure are already nice, intresting and ready to be discuss.

This pull request also containts PR #8

It's a try to solve #7.

I just had a P generic parameter to let user choose which type is use to store kmer.

I also add dependencies to bit_field crate, this crate provide a trait to manipulate bit in object, this trait is implemented for u8, u16, u32, u64, u128, usize, i8, i16, i32, i64, i128 and isize. I think this could be very usefull for us.

A nice side effect of use bit_field crate is P must implement bit_field::BitField so it's reduce number user choice of number.

I also replace macro word_for_k by a const function just for esthetic reason.

Run benchmark on sequence length between 2^i with i is between 8 and 16 to check behaviour of benchmark on different length (check linearity).

Add reverse complement benchmark.

Maybe benchmark code is too rustic, I can rewrite it in more readable way.

Mostly copied from the niffler configs.

There are a lot of jobs defined, but the main ones are

• tests, for running tests on stable, beta, windows (stable) and mac (stable)
• cross_testing, to make sure it works on some non-x86 targets (s390x is especially curious, since it is big endian)
• coverage using tarpaulin (and running on nightly)
• publish to make sure it is publishable on crates.io
• MSRV check (targeting 1.51 for now)
• wasm to check if it builds for webassembly
• test_all_feature_combinations checks that all feature combinations work together, since they should be additive (more info)

So after trying to write some code interfacing kmers with some k-mers parsed using needletail, I realized we are storing the k-mers in the opposite order. In kmers the leftmost nucleotide is stored in the highest order bit, in needletail and some of the other (C++) libraries I have used in the past it's the opposite.

Are there any strong arguments for one of these schemes over the other (cc @natir, @Daniel-Liu-c0deb0t, @luizirber)? Is this something we should also consider making part of the policy of how the k-mers are constructed?

This issue will provide a roadmap for the library, along with specific tasks (TODOs). Ideally we should break these tasks into short and long term tasks and, as the library becomes more mature, tie individual tasks to specific release candidates.

So, it turns out that currently (as of rust v1.51) with MVP const generics, we are not allowed to do simple computations with generic parameters. For example, it would be desirable to have something like this to automatically compute the size of the storage array we want to use based on the value of K provided to the class. However, such a capability is gated behind a feature flag on the nightly branch of rustc. We should determine how we want to handle this from a design perspective. I see a few options:

1. Go with > 1 const generic parameter, to avoid having to do const-generic arithmetic on the receiving end. Once that is available in stable, we can of course simplify the interface. This is a little bit onerous, since now we have to e.g. provide a helper (const) function or macro or some such so that the user doesn't have to think about the number of words that should be used for storage, which is, anyway, error prone.
2. Throw caution to the wind and require nightly rustc with the feature gate to allow ourselves to do the arithmetic we want with const generic parameters.
3. Some other and much more clever solution I've not yet considered.

I'd appreciate others' thoughts on this.

As a fundamental start point for any good software project, we need to decide on the name and generate a logo. I kept the repo name simple, but I was imagining we'd write it out as something like kme-rs or kme.rs or something cute like that ;P. I just wanted to open discussion on this very important topic here to invite opinions and thoughts.

I'm assuming we'd be using 2 bits per bp. There's a couple of common encodings, but I want to suggest

A -> 00
C -> 01
T -> 10
G -> 11

There's a couple of benefits:

1. These are the 2nd and 3rd bits of the ASCII encoding of the corresponding base pairs. Conversion from byte strings would be easy.
2. Complement by using XOR ...0101010

On a slightly unrelated note, I've worked on some sequence manipulation stuff that use SIMD (eg., here, here for a library that was abandoned). Many of these ideas could be applicable here as well. I'm assuming that we want scalar ops only here because SIMD registers are probably too wide (128 or 256 bits) for handling kmers that are relatively short.