bio_seq/
kmer.rs

1// Copyright 2021-2024 Jeff Knaggs
2// Licensed under the MIT license (http://opensource.org/licenses/MIT)
3// This file may not be copied, modified, or distributed
4// except according to those terms.
5
6//! Encoded sequences of static length
7//!
8//! Generally, the underlying storage type of `Kmer` should lend itself to optimisation. The default `Kmer` instance is packed into a `usize`, which can be efficiently `Copy`ed on the stack.
9//!
10//! `k * A::BITS` must fit in the storage type, e.g. `usize` (64 bits).
11//!
12//! ```
13//! use bio_seq::prelude::*;
14//!
15//! for (amino_kmer, amino_string) in Seq::<Amino>::try_from("SSLMNHKKL").unwrap()
16//!         .kmers::<3>()
17//!         .zip(["SSL", "SLM", "LMN", "MNH", "NHK", "HKK", "KKL"])
18//!     {
19//!         assert_eq!(amino_kmer, amino_string);
20//!     }
21//! ```
22//!
23//! Kmers can be copied from other sequence types:
24//!
25//! ```
26//! # use bio_seq::prelude::*;
27//! let kmer: Kmer<Dna, 8> = dna!("AGTTGGCA").try_into().unwrap();
28//! ```
29
30// storage types wider than usize are deliberately narrowed
31#![allow(clippy::cast_possible_truncation)]
32
33use crate::Bs;
34use crate::codec::{self, Codec};
35use crate::prelude::ParseBioError;
36use crate::seq::{Seq, SeqArray, SeqSlice};
37use crate::{
38    Complement, ComplementMut, Reverse, ReverseComplement, ReverseComplementMut, ReverseMut,
39};
40use bitvec::field::BitField;
41use bitvec::view::BitView;
42use core::fmt;
43use core::hash::{Hash, Hasher};
44use core::marker::PhantomData;
45use core::ops::Deref;
46use core::ptr;
47use core::str::FromStr;
48
49pub(crate) mod integral;
50
51#[cfg(feature = "serde")]
52use serde_derive::{Deserialize, Serialize};
53
54const fn make_2bit_table() -> [u8; 256] {
55    let mut table = [0u8; 256];
56    let mut i: usize = 0;
57    while i < 256 {
58        let b0: u8 = (i as u8 & 0b11_00_00_00) >> 6;
59        let b1: u8 = (i as u8 & 0b00_11_00_00) >> 2;
60        let b2: u8 = (i as u8 & 0b00_00_11_00) << 2;
61        let b3: u8 = (i as u8 & 0b00_00_00_11) << 6;
62
63        table[i] = b3 | b2 | b1 | b0;
64        i += 1;
65    }
66    table
67}
68
69const REV_2BIT: [u8; 256] = make_2bit_table();
70
71pub(crate) mod sealed {
72    use crate::Bs;
73
74    pub trait KmerStorage: Copy + Clone + PartialEq + std::fmt::Debug {
75        const BITS: usize;
76        type BaN: AsRef<Bs> + AsMut<Bs>;
77
78        fn to_bitarray(self) -> Self::BaN;
79        fn from_bitslice(bs: &Bs) -> Self;
80
81        fn shiftr(&mut self, n: u32);
82
83        fn complement(&mut self, mask: usize);
84        fn rev_blocks_2(&mut self);
85    }
86}
87
88pub trait KmerStorage: sealed::KmerStorage {}
89
90impl KmerStorage for usize {}
91
92impl KmerStorage for u64 {}
93
94impl KmerStorage for u128 {}
95
96/// By default k-mers are backed by `usize` and `Codec::BITS` * `K` must be <= 64 on 64-bit platforms
97///
98/// ```compile_fail
99/// # use bio_seq::prelude::*;
100/// // 40 * 2 bits does not fit in a usize
101/// let kmer: Kmer<Dna, 40> = "ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGT".parse().unwrap();
102/// ```
103///
104/// ```
105/// # use bio_seq::prelude::*;
106/// // ..but 80 bits does fit into a u128
107/// let kmer: Kmer<Dna, 40, u128> = "ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGT".parse().unwrap();
108/// assert_eq!(kmer.len(), 40);
109/// ```
110#[derive(Debug, PartialEq, Eq, PartialOrd, Ord, Copy, Clone)]
111#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))]
112#[repr(transparent)]
113pub struct Kmer<C: Codec, const K: usize, S: KmerStorage = usize> {
114    pub(crate) _p: PhantomData<C>,
115    pub(crate) bs: S,
116}
117
118impl<A: Codec, const K: usize, S: KmerStorage> Kmer<A, K, S> {
119    // This error message can be formatted with constants in nightly (const_format)
120    const ASSERT_K: () = assert!(
121        K * A::BITS as usize <= S::BITS,
122        "`KmerStorage` not large enough for `Kmer`",
123    );
124
125    const ASSERT_K_NONZERO: () = assert!(K > 0, "`K` must be greater than 0");
126
127    const fn assert_k() {
128        let () = Self::ASSERT_K;
129        let () = Self::ASSERT_K_NONZERO;
130    }
131
132    const BITS: usize = K * A::BITS as usize;
133
134    pub fn len(&self) -> usize {
135        K
136    }
137
138    pub fn is_empty(&self) -> bool {
139        // This is recommended by clippy since we have `len`
140        // Kmers are never empty if K > 0
141        false
142    }
143
144    pub fn rotated_left(&self, n: u32) -> Self {
145        let n: usize = (n as usize % K) * A::BITS as usize;
146        let mut ba = self.bs.to_bitarray();
147        let bs: &mut Bs = ba.as_mut();
148        bs[..Self::BITS].rotate_left(n);
149
150        Kmer {
151            _p: PhantomData,
152            bs: S::from_bitslice(&bs[..Self::BITS]),
153        }
154    }
155
156    pub fn rotated_right(&self, n: u32) -> Self {
157        let n: usize = (n as usize % K) * A::BITS as usize;
158        let mut ba = self.bs.to_bitarray();
159        let bs: &mut Bs = ba.as_mut();
160        bs[..Self::BITS].rotate_right(n);
161
162        Kmer {
163            _p: PhantomData,
164            bs: S::from_bitslice(&bs[..Self::BITS]),
165        }
166    }
167
168    /// Shift bases to the right and push a base onto the end.
169    ///
170    /// ```
171    /// use bio_seq::prelude::*;
172    /// use bio_seq::codec::dna::Dna;
173    ///
174    /// let k = kmer!("ACGAT");
175    /// assert_eq!(k.pushr(Dna::T).to_string(), "CGATT");
176    /// ```
177    pub fn pushr(self, base: A) -> Self {
178        let mut ba = self.rotated_left(1).bs.to_bitarray();
179        let bs: &mut Bs = ba.as_mut();
180
181        let start = Self::BITS - A::BITS as usize;
182        let end = start + A::BITS as usize;
183
184        bs[start..end].store(base.to_bits());
185
186        Kmer {
187            _p: PhantomData,
188            bs: S::from_bitslice(bs),
189        }
190    }
191
192    /// Push a base from the left
193    pub fn pushl(self, base: A) -> Self {
194        let mut ba = self.rotated_right(1).bs.to_bitarray();
195        let bs: &mut Bs = ba.as_mut();
196
197        bs[..A::BITS as usize].store(base.to_bits());
198
199        Kmer {
200            _p: PhantomData,
201            bs: S::from_bitslice(bs),
202        }
203    }
204
205    /// Create Kmer from sequence without checking length
206    pub fn unsafe_from_seqslice(seq: &SeqSlice<A>) -> Self {
207        Self::assert_k();
208        debug_assert!(K == seq.len(), "K != seq.len()");
209        Kmer {
210            _p: PhantomData,
211            bs: S::from_bitslice(&seq.bs),
212        }
213    }
214
215    fn complement(&mut self) {
216        self.bs.complement(K * A::BITS as usize);
217    }
218
219    fn rev_blocks_2(&mut self) {
220        // TODO: assert K == 2
221        self.bs.rev_blocks_2();
222        self.bs.shiftr((S::BITS - (A::BITS as usize * K)) as u32);
223    }
224}
225
226impl<A: Codec, const K: usize> From<usize> for Kmer<A, K, usize> {
227    fn from(i: usize) -> Kmer<A, K, usize> {
228        Self::assert_k();
229        Kmer {
230            _p: PhantomData,
231            bs: i,
232        }
233    }
234}
235
236impl<A: Codec, const K: usize> From<u64> for Kmer<A, K, u64> {
237    fn from(i: u64) -> Kmer<A, K, u64> {
238        Self::assert_k();
239        Kmer {
240            _p: PhantomData,
241            bs: i,
242        }
243    }
244}
245
246impl<A: Codec, const K: usize> From<usize> for Kmer<A, K, u64> {
247    fn from(i: usize) -> Kmer<A, K, u64> {
248        Self::assert_k();
249        Kmer {
250            _p: PhantomData,
251            bs: i as u64,
252        }
253    }
254}
255
256/*
257impl<A: Codec, const K: usize, S: KmerStorage> From<&SeqSlice<A>> for Kmer<A, K, S> {
258    fn from(seq: &SeqSlice<A>) -> Self {
259        Kmer {
260            _p: PhantomData,
261            bs: S::from_bitslice(&seq.bs),
262        }
263    }
264}
265*/
266
267impl<S: KmerStorage + Into<usize>, A: Codec, const K: usize> From<&Kmer<A, K, S>> for usize {
268    fn from(kmer: &Kmer<A, K, S>) -> usize {
269        kmer.bs.into()
270    }
271}
272
273impl<A: Codec, const K: usize> Deref for Kmer<A, K, usize> {
274    type Target = SeqSlice<A>;
275
276    fn deref(&self) -> &Self::Target {
277        let bs: &Bs = &self.bs.view_bits()[0..(K * A::BITS as usize)];
278        let bs: *const Bs = ptr::from_ref::<Bs>(bs);
279        unsafe { &*(bs as *const SeqSlice<A>) }
280    }
281}
282
283impl<A: Codec, const K: usize> AsRef<SeqSlice<A>> for Kmer<A, K, usize> {
284    fn as_ref(&self) -> &SeqSlice<A> {
285        self
286    }
287}
288
289impl<A: Codec, const K: usize, S: KmerStorage> fmt::Display for Kmer<A, K, S> {
290    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
291        let mut s = String::new();
292        let ba = self.bs.to_bitarray();
293        let bs: &Bs = &ba.as_ref()[0..(K * A::BITS as usize)];
294
295        bs.chunks(A::BITS as usize).for_each(|chunk| {
296            s.push(A::unsafe_from_bits(chunk.load_le::<u8>()).to_char());
297        });
298        write!(f, "{s}")
299    }
300}
301
302/// An iterator over all kmers of a sequence with a specified length
303pub struct KmerIter<'a, A: Codec, const K: usize> {
304    pub(crate) slice: &'a SeqSlice<A>,
305    pub(crate) index: usize,
306    pub(crate) len: usize,
307    pub(crate) _p: PhantomData<A>,
308}
309
310impl<A: Codec, const K: usize, S: KmerStorage> Kmer<A, K, S> {
311    fn unsafe_from(seq: &SeqSlice<A>) -> Self {
312        Self::assert_k();
313        debug_assert!(K == seq.len(), "K != seq.len()");
314        Kmer {
315            _p: PhantomData,
316            bs: S::from_bitslice(&seq.bs),
317        }
318    }
319}
320
321impl<A: Codec, const K: usize> Iterator for KmerIter<'_, A, K> {
322    type Item = Kmer<A, K>;
323    fn next(&mut self) -> Option<Kmer<A, K>> {
324        let i = self.index;
325        if self.index + K > self.len {
326            return None;
327        }
328        self.index += 1;
329        Some(Kmer::<A, K>::unsafe_from(&self.slice[i..i + K]))
330    }
331
332    fn size_hint(&self) -> (usize, Option<usize>) {
333        let n = (self.len.saturating_sub(self.index) + 1).saturating_sub(K);
334        (n, Some(n))
335    }
336}
337
338impl<A: Codec, const K: usize> ExactSizeIterator for KmerIter<'_, A, K> {}
339
340/// ```
341/// use bio_seq::prelude::*;
342/// use std::hash::{Hash, Hasher, DefaultHasher};
343///
344/// let mut hasher1 = DefaultHasher::new();
345/// kmer!("AAA").hash(&mut hasher1);
346/// let hash1 = hasher1.finish();
347///
348/// let mut hasher2 = DefaultHasher::new();
349/// kmer!("AAAA").hash(&mut hasher2);
350/// let hash2 = hasher2.finish();
351///
352/// assert_ne!(hash1, hash2);
353/// ```
354impl<A: Codec, const K: usize, S: KmerStorage> Hash for Kmer<A, K, S> {
355    fn hash<H: Hasher>(&self, state: &mut H) {
356        // Length prefix as 8 little-endian bytes (see `SeqSlice` for rationale).
357        state.write(&(K as u64).to_le_bytes());
358        let ba = self.bs.to_bitarray();
359        let bs: &Bs = ba.as_ref();
360        // Only the meaningful bits are hashed; the unused high bits of the
361        // storage word/array are padding and must not influence the digest.
362        let n_bits = K * A::BITS as usize;
363        crate::hash::hash_bits(&bs[..n_bits], state);
364    }
365}
366
367impl<A: Codec, const K: usize, S: KmerStorage> TryFrom<&SeqSlice<A>> for Kmer<A, K, S> {
368    type Error = ParseBioError;
369
370    fn try_from(seq: &SeqSlice<A>) -> Result<Self, Self::Error> {
371        if seq.len() == K {
372            Ok(Kmer::<A, K, S>::unsafe_from(&seq[0..K]))
373        } else {
374            Err(ParseBioError::MismatchedLength(K, seq.len()))
375        }
376    }
377}
378
379impl<A: Codec, const K: usize> TryFrom<Seq<A>> for Kmer<A, K> {
380    type Error = ParseBioError;
381
382    fn try_from(seq: Seq<A>) -> Result<Self, Self::Error> {
383        Self::try_from(seq.as_ref())
384    }
385}
386
387impl<A: Codec, const K: usize, S: KmerStorage> PartialEq<SeqArray<A, K, 1>> for Kmer<A, K, S> {
388    fn eq(&self, seq: &SeqArray<A, K, 1>) -> bool {
389        if seq.len() != K {
390            return false;
391        }
392        &Kmer::<A, K, S>::unsafe_from(seq.as_ref()) == self
393    }
394}
395
396impl<A: Codec, const K: usize, S: KmerStorage> PartialEq<&SeqArray<A, K, 1>> for Kmer<A, K, S> {
397    fn eq(&self, seq: &&SeqArray<A, K, 1>) -> bool {
398        if seq.len() != K {
399            return false;
400        }
401        &Kmer::<A, K, S>::unsafe_from(seq.as_ref()) == self
402    }
403}
404
405impl<A: Codec, const K: usize> PartialEq<Seq<A>> for Kmer<A, K> {
406    fn eq(&self, seq: &Seq<A>) -> bool {
407        if seq.len() != K {
408            return false;
409        }
410        &Kmer::<A, K>::unsafe_from(seq.as_ref()) == self
411    }
412}
413
414impl<A: Codec, const K: usize, S: KmerStorage> PartialEq<SeqSlice<A>> for Kmer<A, K, S> {
415    fn eq(&self, seq: &SeqSlice<A>) -> bool {
416        if seq.len() != K {
417            return false;
418        }
419        &Kmer::<A, K, S>::unsafe_from(seq) == self
420    }
421}
422
423impl<A: Codec, const K: usize, S: KmerStorage> PartialEq<&SeqSlice<A>> for Kmer<A, K, S> {
424    fn eq(&self, seq: &&SeqSlice<A>) -> bool {
425        if seq.len() != K {
426            return false;
427        }
428        &Kmer::<A, K, S>::unsafe_from(seq) == self
429    }
430}
431
432impl<A: Codec, const K: usize> PartialEq<&str> for Kmer<A, K> {
433    fn eq(&self, seq: &&str) -> bool {
434        &self.to_string() == seq
435    }
436}
437
438impl<A: Codec, const K: usize, S: KmerStorage> FromStr for Kmer<A, K, S> {
439    type Err = ParseBioError;
440
441    fn from_str(s: &str) -> Result<Self, Self::Err> {
442        if s.len() != K {
443            return Err(ParseBioError::MismatchedLength(K, s.len()));
444        }
445        let seq: Seq<A> = Seq::from_str(s)?;
446        Kmer::<A, K, S>::try_from(seq.as_ref())
447    }
448}
449
450impl<A: Codec, const K: usize> From<Kmer<A, K, usize>> for Seq<A> {
451    fn from(kmer: Kmer<A, K, usize>) -> Self {
452        let mut seq: Seq<A> = Seq::with_capacity(K);
453        seq.extend(kmer.iter());
454        seq
455    }
456}
457
458impl<const K: usize> ComplementMut for Kmer<codec::dna::Dna, K, usize> {
459    fn comp(&mut self) {
460        self.complement();
461    }
462}
463
464impl<const K: usize> Complement for Kmer<codec::dna::Dna, K, usize> {}
465
466impl<A: Codec, const K: usize> ReverseMut for Kmer<A, K, usize> {
467    fn rev(&mut self) {
468        self.rev_blocks_2();
469    }
470}
471
472impl<A: Codec, const K: usize> Reverse for Kmer<A, K, usize> {}
473
474impl<const K: usize> ReverseComplementMut for Kmer<codec::dna::Dna, K, usize> {}
475
476impl<const K: usize> ReverseComplement for Kmer<codec::dna::Dna, K, usize> {}
477
478/// Convenient compile time kmer constructor
479///
480/// This is a wrapper for the `dna!` macro that returns a `Kmer`:
481/// ```
482/// # use bio_seq::prelude::*;
483/// let kmer: Kmer<Dna, 8> = kmer!("ACGTACGT");
484/// ```
485#[macro_export]
486macro_rules! kmer {
487    ($seq:expr) => {
488        $crate::kmer::Kmer::<$crate::codec::dna::Dna, { $seq.len() }>::unsafe_from_seqslice(dna!($seq))
489    };
490    ($seq:expr, $storage:ty) => {
491        $crate::kmer::Kmer::<$crate::codec::dna::Dna, { $seq.len() }, $storage>::unsafe_from_seqslice(
492            dna!($seq),
493        )
494    };
495}
496
497#[cfg(test)]
498mod tests {
499    use crate::prelude::*;
500    use crate::seq::SeqArray;
501
502    #[test]
503    #[allow(clippy::cast_sign_loss)]
504    fn kmer_to_usize() {
505        let s: &'static SeqSlice<Dna> = dna!("AACTT");
506        println!("{s}");
507
508        for (kmer, index) in s.kmers::<2>().zip([0b00_00, 0b01_00, 0b11_01, 0b11_11]) {
509            println!("{kmer}");
510            assert_eq!(index as usize, (&kmer).into());
511        }
512    }
513    #[test]
514    fn pushl_test() {
515        let k = kmer!("ACGT");
516        let k1 = k.pushl(Dna::G);
517        let k2 = k1.pushl(Dna::A);
518        let k3 = k2.pushl(Dna::T);
519        let k4 = k3.pushl(Dna::C);
520        let k5 = k4.pushl(Dna::C);
521
522        assert_eq!(k1, kmer!("GACG"));
523        assert_eq!(k2, kmer!("AGAC"));
524        assert_eq!(k3, kmer!("TAGA"));
525        assert_eq!(k4, kmer!("CTAG"));
526        assert_eq!(k5, kmer!("CCTA"));
527    }
528    #[test]
529    fn pushr_test() {
530        let k = kmer!("ACGT");
531        let k1 = k.pushr(Dna::G);
532        let k2 = k1.pushr(Dna::A);
533        let k3 = k2.pushr(Dna::T);
534        let k4 = k3.pushr(Dna::C);
535        let k5 = k4.pushr(Dna::C);
536
537        println!("{k1}");
538        assert_eq!(k1, kmer!("CGTG"));
539        assert_eq!(k2, kmer!("GTGA"));
540        assert_eq!(k3, kmer!("TGAT"));
541        assert_eq!(k4, kmer!("GATC"));
542        assert_eq!(k5, kmer!("ATCC"));
543    }
544
545    #[test]
546    #[allow(clippy::cast_sign_loss)]
547    fn amino_kmer_to_usize() {
548        for (kmer, index) in Seq::<Amino>::try_from("SRY")
549            .unwrap()
550            .kmers::<2>()
551            .zip([0b0010_0001_1000, 0b0100_1100_1000])
552        {
553            assert_eq!(index as usize, usize::from(&kmer));
554        }
555    }
556    #[test]
557    fn big_kmer_shiftr() {
558        let mut kmer: Kmer<Dna, 32, u64> = kmer!("AATTTGTGGGTTCGTCTGCGGCTCCGCCCTTA", u64);
559        for base in dna!("TACTATGAGGACGATCAGCACCATAAGAACAAA") {
560            kmer = kmer.pushr(base);
561        }
562        assert_eq!(kmer!("ACTATGAGGACGATCAGCACCATAAGAACAAA", u64), kmer);
563    }
564
565    #[test]
566    fn big_kmer_shiftl() {
567        let mut kmer: Kmer<Dna, 32, u64> = kmer!("AATTTGTGGGTTCGTCTGCGGCTCCGCCCTTA", u64);
568        for base in dna!("GTACTATGAGGACGATCAGCACCATAAGAACAAA") {
569            kmer = kmer.pushl(base);
570        }
571        assert_eq!(kmer!("AAACAAGAATACCACGACTAGCAGGAGTATCA", u64), kmer);
572    }
573
574    #[test]
575    fn amino_kmer_iter() {
576        for (kmer, target) in Seq::<Amino>::try_from("SSLMNHKKL")
577            .unwrap()
578            .kmers::<3>()
579            .zip(["SSL", "SLM", "LMN", "MNH", "NHK", "HKK", "KKL"])
580        {
581            assert_eq!(kmer, target);
582        }
583    }
584
585    #[test]
586    fn test_rotations() {
587        let kmer: Kmer<Dna, 9> = Kmer::try_from(dna!("ACTGCGATG")).unwrap();
588
589        for (shift, rotation) in vec![
590            "ACTGCGATG",
591            "CTGCGATGA",
592            "TGCGATGAC",
593            "GCGATGACT",
594            "CGATGACTG",
595            "GATGACTGC",
596            "ATGACTGCG",
597            "TGACTGCGA",
598            "GACTGCGAT",
599            "ACTGCGATG",
600            "CTGCGATGA",
601            "TGCGATGAC",
602        ]
603        .into_iter()
604        .enumerate()
605        {
606            //            println!("{} {} {}", shift, kmer.rotated_left(shift as u32), rotation);
607            assert_eq!(kmer.rotated_left(shift as u32), rotation);
608        }
609
610        for (shift, rotation) in vec![
611            "ACTGCGATG",
612            "GACTGCGAT",
613            "TGACTGCGA",
614            "ATGACTGCG",
615            "GATGACTGC",
616            "CGATGACTG",
617            "GCGATGACT",
618            "TGCGATGAC",
619            "CTGCGATGA",
620            "ACTGCGATG",
621            "GACTGCGAT",
622        ]
623        .into_iter()
624        .enumerate()
625        {
626            //            println!("{} {} {}", shift, kmer.rotated_right(shift as u32), rotation);
627            assert_eq!(kmer.rotated_right(shift as u32), rotation);
628        }
629
630        let kmer: Kmer<Dna, 8> = Kmer::try_from(dna!("ACTGCGAT")).unwrap().rotated_left(1);
631
632        assert_ne!(kmer.to_string(), "ACTGCGAT");
633        assert_eq!(kmer.to_string(), "CTGCGATA");
634
635        let kmer: Kmer<Dna, 8> = kmer!("ACTGCGAT").rotated_right(1);
636
637        assert_ne!(kmer.to_string(), "ACTGCGAT");
638        assert_eq!(kmer.to_string(), "TACTGCGA");
639
640        let kmer: Kmer<Dna, 9> = Kmer::from_str("ACTGCGATG").unwrap().rotated_left(0);
641
642        assert_eq!(kmer.to_string(), "ACTGCGATG");
643        assert_ne!(kmer.to_string(), "ACTGCGATGA");
644
645        let kmer: Kmer<Dna, 9> = Kmer::try_from(dna!("ACTGCGATG")).unwrap().rotated_right(0);
646
647        assert_eq!(kmer.to_string(), "ACTGCGATG");
648        assert_ne!(kmer.to_string(), "ACTGCGATGA");
649
650        let kmer: Kmer<Dna, 9> = Kmer::from_str("ACTGCGATG").unwrap().rotated_left(9 * 3307);
651
652        assert_eq!(kmer.to_string(), "ACTGCGATG");
653        assert_ne!(kmer.to_string(), "ACTGCGATGA");
654
655        let kmer: Kmer<Dna, 9> = kmer!("ACTGCGATG").rotated_right(9 * 3307);
656
657        assert_eq!(kmer.to_string(), "ACTGCGATG");
658        assert_ne!(kmer.to_string(), "ACTGCGATGA");
659    }
660
661    #[test]
662    fn eq_functions() {
663        assert_eq!(kmer!("ACGT"), dna!("ACGT"));
664
665        // this should be a compiler error:
666        // assert_ne!(kmer!("ACGT"), dna!("ACGTA"));
667
668        let kmer: Kmer<Iupac, 4> = Kmer::from_str("ACGT").unwrap();
669        assert_eq!(kmer, iupac!("ACGT"));
670        assert_ne!(kmer, iupac!("NCGT"));
671    }
672
673    #[test]
674    fn kmer_iter() {
675        //let seq = dna!("ACTGA");
676        let cs: Vec<Kmer<Dna, 3>> = dna!("ACTGA").kmers().collect();
677        assert_eq!(cs[0], "ACT");
678        assert_eq!(cs[1], "CTG");
679        assert_eq!(cs[2], "TGA");
680        assert_eq!(cs.len(), 3);
681    }
682
683    #[test]
684    fn k_check() {
685        let _kmer = Kmer::<Dna, 32>::from(0);
686        let _kmer = Kmer::<Amino, 10>::from(0);
687        let _kmer = Kmer::<Iupac, 14>::from(0);
688    }
689
690    #[test]
691    fn kmer_revcomp() {
692        assert_eq!(kmer!("ACGT"), kmer!("ACGT").to_revcomp());
693        assert_ne!(kmer!("GTCGTA"), kmer!("TACGAC"));
694
695        let rc = kmer!("GTCGTA").to_revcomp();
696
697        assert_eq!(rc, kmer!("TACGAC"));
698
699        assert_eq!(
700            kmer!("GCTATCGATCTGATCG"),
701            kmer!("CGATCAGATCGATAGC").to_revcomp()
702        );
703    }
704
705    #[test]
706    fn kmer_deref() {
707        let kmer: Kmer<Dna, 3> = kmer!("ACG");
708        let seq: &SeqSlice<Dna> = &kmer;
709
710        assert_eq!(*seq, *kmer);
711        assert_eq!(seq.to_string(), "ACG");
712
713        let kmer: Kmer<Dna, 8> = kmer!("AAAAAAAA");
714        let seq: &SeqSlice<Dna> = &kmer;
715
716        assert_eq!(seq.to_string(), "AAAAAAAA");
717
718        let kmer: Kmer<Dna, 8> = kmer!("TTTTTTTT");
719        let seq: &SeqSlice<Dna> = &kmer;
720
721        assert_eq!(seq.to_string(), "TTTTTTTT");
722
723        let kmer: Kmer<Dna, 16> = kmer!("AGCTAGCTAGCTAGCT");
724        let seq: &SeqSlice<Dna> = &kmer;
725
726        assert_eq!(seq.to_string(), "AGCTAGCTAGCTAGCT");
727    }
728
729    #[test]
730    fn kmer_as_ref() {
731        let kmer: Kmer<Dna, 4> = kmer!("ACGT");
732        let seq: &SeqSlice<Dna> = kmer.as_ref();
733
734        assert_eq!(seq.to_string(), "ACGT");
735
736        let kmer: Kmer<Dna, 16> = kmer!("AGCTAGCTAGCTAGCT");
737        let seq: &SeqSlice<Dna> = kmer.as_ref();
738
739        assert_eq!(seq.to_string(), "AGCTAGCTAGCTAGCT");
740    }
741
742    #[test]
743    fn kmer_rev() {
744        let mut kmer: Kmer<Dna, 4> = kmer!("ACGT");
745
746        kmer.rev();
747
748        assert_eq!(kmer.to_string(), "TGCA");
749
750        kmer.rev();
751
752        assert_eq!(kmer.to_string(), "ACGT");
753
754        kmer.rev();
755
756        assert_eq!(kmer.to_string(), "TGCA");
757    }
758
759    #[test]
760    #[allow(clippy::needless_borrow)]
761    fn kmer_storage_types() {
762        let s1 = "AACGTAGCCGCGAACTTACGTAGCCGCGAAAA";
763        let s2 = "AACGTAGCCGCGAACTTACGTAGCCGCGAAA";
764        let s3 = "ACGTAGCCGCGAACTTACGTAGCCGCGAAAA";
765
766        let s4 = "AACGTAGCCGCGAACTTACGTAGCCGCGAAAAAACGTAGCCGCGAACTTACGTAGCCGCGAAAA";
767        let s5 = "AACGTAGCCGCGAACTTACGTAGCCGCGAAAAAACGTAGCCGCGAACTTACGTAGCCGCGAAAAA";
768
769        assert_eq!(s1.len(), 32);
770        assert_eq!(s2.len(), 31);
771        assert_eq!(s3.len(), 31);
772        assert_eq!(s4.len(), 64);
773        assert_eq!(s5.len(), 65);
774
775        let kmer1_64 = Kmer::<Dna, 32, u64>::from_str(&s1).unwrap();
776        let kmer2_64 = Kmer::<Dna, 31, u64>::from_str(&s2).unwrap();
777        let kmer3_64 = Kmer::<Dna, 31, u64>::from_str(&s3).unwrap();
778
779        let kmer1 = Kmer::<Dna, 32, u64>::from_str(&s1).unwrap();
780        let kmer2 = Kmer::<Dna, 31, u64>::from_str(&s2).unwrap();
781        let kmer3 = Kmer::<Dna, 31, u64>::from_str(&s3).unwrap();
782
783        let kmer4_128 = Kmer::<Dna, 64, u128>::from_str(&s4).unwrap();
784
785        let seq5: Seq<Dna> = s5.try_into().unwrap();
786
787        assert_eq!(kmer4_128, &seq5[..64]);
788        assert_ne!(kmer4_128, &seq5[1..]);
789
790        assert_eq!(kmer1, &seq5[..32]);
791        assert_eq!(kmer1, &seq5[32..64]);
792
793        assert_eq!(kmer1_64, &seq5[..32]);
794        assert_eq!(kmer1_64, &seq5[32..64]);
795
796        assert_ne!(kmer1, &seq5[..31]);
797        assert_ne!(kmer1, &seq5[32..]);
798
799        assert_ne!(kmer1_64, &seq5[1..33]);
800        assert_ne!(kmer1_64, &seq5[33..]);
801
802        assert_eq!(kmer4_128, kmer4_128);
803        assert_eq!(kmer1_64, kmer1_64);
804        assert_eq!(kmer2, kmer2);
805
806        assert_ne!(kmer2, kmer3);
807        assert_ne!(kmer2_64, kmer3_64);
808        // PartialEq is not implemented for different storgage types
809        /*
810                assert_ne!(kmer2, kmer3_64);
811                assert_eq!(kmer2, kmer2_64);
812                assert_eq!(kmer1, kmer1_64);
813        */
814    }
815
816    #[test]
817    fn try_from_seq() {
818        let seq: Seq<Dna> = Seq::try_from("ACACACACACACGT").unwrap();
819        assert_eq!(
820            Kmer::<Dna, 8>::try_from(&seq[..8]).unwrap().to_string(),
821            "ACACACAC"
822        );
823        assert_eq!(
824            Kmer::<Dna, 8>::try_from(&seq[1..9]).unwrap().to_string(),
825            "CACACACA"
826        );
827
828        let err: Result<Kmer<Dna, 8>, ParseBioError> = Kmer::try_from(&seq[2..9]);
829        assert_eq!(err, Err(ParseBioError::MismatchedLength(8, 7)));
830
831        let err: Result<Kmer<Dna, 8>, ParseBioError> = Kmer::try_from(seq);
832        assert_eq!(err, Err(ParseBioError::MismatchedLength(8, 14)));
833
834        let seq: Seq<Dna> = Seq::try_from("ACACACACACACGT").unwrap();
835
836        assert_eq!(
837            Kmer::<Dna, 14>::try_from(seq).unwrap().to_string(),
838            "ACACACACACACGT"
839        );
840    }
841}
bio_seq/kmer.rs

bio_seq/
kmer.rs
