Handling of N characters

[RagnarGrootKoerkamp]

Would be good to have better policies around this, and ideally the user can choose:

• panic on N
• split sequence on N
• or return one allocated object, and a list of N-free slices?
• implicit cast from N to ACTG, at least for packed-seq.

[rob-p]

These are all reasonable policies. What do you think about one or two different functions depending on what the user wants? I think panic on N is a harsh default policy, which I'd prefer to avoid. Implicit cast is probably the easiest, though the user then loses some control. I like either splitting on N or returning N-free segments with a list of N occurrences. In practice, for the latter, you might want to provide an NRun struct or some such, as I suspect, in some cases, you'll just have long stretches of N (e.g. at the telomere or centromere of chromosomes).

[RagnarGrootKoerkamp]

Ok so some thoughts:

• By default, I would prefer panic on N, so we never silently do things the user didn't expect. Unfortunately though, it's kind of annoying to detect while bitpacking that we are only packing ACTG characters. So in practice, I think we now panic only on the very first and/or last characters for the bit that's not a multiple of 32 or so. This is also very inconsistent of course, and we should fix it.
• Yes, implicit cast is for sure easy to provide, but then I'd like it as PackedSeqVec::from_ascii_implicit_N(seq: &[u8]) or so, so that this is very explicitly a user choice.
• I'm not fully sure yet on 'split on N' vs 'return a list of N'. The drawback of split-on-N is that short slices don't work nearly as well with all the SIMD since you pay a 'startup cost' per slice. So indeed, returning a list of (ranges of) N positions sounds nicer, so that the results corresponding to those slices can be filtered away later.
• Another alternative would be to produce an additional bitstream of position is non-ACTG, but since N should be rare, space-wise it's probably better to store the indices instead.

@imartayan got anything to add?

Anyway, this is mostly blocked on finding some time to work on the actual implementation for this.

Also, some folks at INRIA are working on a fast fasta parser and it may make sense to go directly from fasta to PackedSeqVec, without first checking for \n and collecting into a single contiguous vec.

[RagnarGrootKoerkamp]

With the release of 4.0.0, there is now PackedNSeq::from_ascii that uses a 2+1 bit encoding: 2 bits for the classic ACTG values, and one additional indicator bit for anything ambiguous (non-ACTG).

I'll leave this open though, since the story around detecting N characters for plain PackedSeq is not quite there yet. Currently it selently accepts everything, but we should do better.

Specifically now that we have PackedNSeq, I could imagine something like:

• PackedNSeq only stores the indicator vector if needed at all. (When omitted, the ambiguous_kmers iterator trivially returns false always).
• PackedSeq panics on non-ACTG characters.