annonars 0.44.2

Genome annotation based on Rust and RocksDB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70

//! Data structures for (de-)serialization as generated by `prost-build`.

use noodles::vcf::variant::record::AlternateBases;

pub use crate::pbs::helixmtdb::Record;
use noodles::vcf::variant::record_buf::info::field;

impl Record {
    /// Creates a new `Record` from a VCF record and allele number.
    pub fn from_vcf_allele(
        record: &noodles::vcf::variant::RecordBuf,
        allele_no: usize,
    ) -> Result<Self, anyhow::Error> {
        let chrom = record.reference_sequence_name().to_string();
        let pos: usize = record
            .variant_start()
            .expect("Telomeric breakends not supported")
            .get();
        let pos = i32::try_from(pos)?;
        let ref_allele = record.reference_bases().to_string();
        let alt_allele = record
            .alternate_bases()
            .iter()
            .nth(allele_no)
            .ok_or_else(|| anyhow::anyhow!("no such allele: {}", allele_no))??
            .to_string();
        let num_total =
            if let Some(Some(field::Value::Integer(num_total))) = record.info().get("AN") {
                *num_total
            } else {
                anyhow::bail!("missing INFO/AN in HelixMtDb record")
            };
        let num_het =
            if let Some(Some(field::Value::Integer(num_het))) = record.info().get("AC_het") {
                *num_het
            } else {
                anyhow::bail!("missing INFO/AC in HelixMtDb record")
            };
        let num_hom =
            if let Some(Some(field::Value::Integer(num_hom))) = record.info().get("AC_hom") {
                *num_hom
            } else {
                anyhow::bail!("missing INFO/AC_hom in HelixMtDb record")
            };
        let feature_type =
            if let Some(Some(field::Value::String(feature))) = record.info().get("FEATURE") {
                feature.to_string()
            } else {
                anyhow::bail!("missing INFO/FEATURE in HelixMtDb record")
            };
        let gene_name =
            if let Some(Some(field::Value::String(gene_name))) = record.info().get("GENE") {
                gene_name.to_string()
            } else {
                anyhow::bail!("missing INFO/GENE in HelixMtDb record")
            };

        Ok(Record {
            chrom,
            pos,
            ref_allele,
            alt_allele,
            num_total,
            num_het,
            num_hom,
            feature_type,
            gene_name,
        })
    }
}