Merge branch 'master' of github.com:rust-bio/rust-bio-tools

Author: johanneskoester

Cargo.toml

@@ -1,6 +1,6 @@
 [package]
 name = "rust-bio-tools"
-version = "0.10.3-alpha.0"
+version = "0.10.4-alpha.0"
 authors = ["Johannes Köster <johannes.koester@tu-dortmund.de>", "Erik Clarke <ecl@pennmedicine.upenn.edu>"]
 description = "A set of fast and robust command line utilities for bioinformatics tasks based on Rust-Bio."
 license-file = "LICENSE.md"

src/bcf/annotate_dgidb.rs

@@ -33,18 +33,28 @@ pub fn annotate_dgidb(
     vcf_path: &str,
     api_path: String,
     field_name: &str,
+    datasources: Option<Vec<&str>>,
     genes_per_request: usize,
 ) -> Result<(), Box<dyn Error>> {
-    let gene_drug_interactions = request_interaction_drugs(vcf_path, api_path, genes_per_request)?;
+    let gene_drug_interactions =
+        request_interaction_drugs(vcf_path, api_path, datasources, genes_per_request)?;
     modify_vcf_entries(vcf_path, gene_drug_interactions, field_name)
 }
 
 fn request_interaction_drugs(
     vcf_path: &str,
     api_path: String,
+    datasources_opt: Option<Vec<&str>>,
     genes_per_request: usize,
 ) -> Result<Option<HashMap<String, Vec<(String, Vec<String>)>>>, Box<dyn Error>> {
     let mut genes = collect_genes(vcf_path)?;
+    let datasources = if let Some(entries) = datasources_opt {
+        let mut b = String::from("&interaction_sources=");
+        b.push_str(entries.join(",").as_str());
+        b
+    } else {
+        String::new()
+    };
     if genes.is_empty() {
         return Ok(None);
     }
@@ -57,6 +67,7 @@ fn request_interaction_drugs(
     {
         let mut slice_api_path = api_path.clone();
         slice_api_path.push_str(gene_slice.join(",").as_str());
+        slice_api_path.push_str(datasources.as_str());
         let res: Dgidb = reqwest::get(&slice_api_path)?.json()?;
 
         for term in res.matched_terms {

src/bcf/oncoprint.rs

@@ -12,7 +12,10 @@ use tera::{self, Context, Tera};
 
 use rust_htslib::bcf::{self, Read};
 
-pub fn oncoprint(sample_calls: &HashMap<String, String>) -> Result<(), Box<dyn Error>> {
+pub fn oncoprint(
+    sample_calls: &HashMap<String, String>,
+    vep_annotation: bool,
+) -> Result<(), Box<dyn Error>> {
     let mut data = HashMap::new();
     for (sample, path) in sample_calls.iter().sorted() {
         let mut genes = HashMap::new();
@@ -55,8 +58,16 @@ pub fn oncoprint(sample_calls: &HashMap<String, String>) -> Result<(), Box<dyn E
                         }
 
                         let gene = str::from_utf8(fields[3])?;
-                        let dna_alteration = str::from_utf8(fields[9])?;
-                        let protein_alteration = str::from_utf8(fields[10])?;
+                        let dna_alteration = if vep_annotation {
+                            str::from_utf8(fields[10])?
+                        } else {
+                            str::from_utf8(fields[9])?
+                        };
+                        let protein_alteration = if vep_annotation {
+                            str::from_utf8(fields[11])?
+                        } else {
+                            str::from_utf8(fields[10])?
+                        };
 
                         let rec = genes
                             .entry(gene.to_owned())

src/cli.yaml

@@ -187,6 +187,14 @@ subcommands:
               short: f
               default_value: dgiDB_drugs
               help: Info field name to be used for annotation.
+          - datasources:
+              long: sources
+              short: s
+              value_name: str
+              multiple: true
+              help: |
+                A list of data sources included in query. If omitted all sources are considered.
+                A list of all sources can be found at http://dgidb.org/api/v2/interaction_sources.json
           - genes-per-request:
               long: genes-per-request
               short: g
@@ -208,6 +216,9 @@ subcommands:
             multiple: true
             value_name: NAME=FILE
             help: VCF/BCF files to include (single-sample). Name is the sample name that will be used in the oncoprint.
+        - vep-annotation:
+            long: vep-annotation
+            help: Annotation field gets parsed as definded by VEP.
 
   - collapse-reads-to-fragments:
       about: |

src/main.rs

@@ -65,6 +65,10 @@ fn main() -> Result<(), Box<dyn Error>> {
             &matches.value_of("vcf").unwrap(),
             matches.value_of("api-path").unwrap().to_string(),
             &matches.value_of("field").unwrap(),
+            match matches.is_present("datasources") {
+                true => Some(matches.values_of("datasources").unwrap().collect()),
+                false => None,
+            },
             value_t!(matches, "genes-per-request", usize).unwrap(),
         ),
         ("oncoprint", Some(matches)) => {
@@ -74,7 +78,7 @@ fn main() -> Result<(), Box<dyn Error>> {
                 sample_calls.insert(e[0].to_owned(), e[1].to_owned());
             }
 
-            bcf::oncoprint::oncoprint(&sample_calls)
+            bcf::oncoprint::oncoprint(&sample_calls, matches.is_present("vep-annotation"))
         }
         ("collapse-reads-to-fragments", Some(matches)) => match matches.subcommand() {
             ("fastq", Some(matches)) => {

tests/expected/annotate_dgidb_drugbank_test.bcf

@@ -249,6 +249,19 @@ fn test_vcf_annotate_dgidb() {
     );
 }
 
+#[test]
+fn test_vcf_annotate_dgidb_drugbank() {
+    assert!(
+        Command::new("bash")
+            .arg("-c")
+            .arg("target/debug/rbt vcf-annotate-dgidb tests/annotate_dgidb_test.vcf -s DrugBank> /tmp/annotate_dgidb_drugbank_test.bcf")
+            .spawn().unwrap().wait().unwrap().success());
+    test_output(
+        "/tmp/annotate_dgidb_drugbank_test.bcf",
+        "tests/expected/annotate_dgidb_drugbank_test.bcf",
+    );
+}
+
 #[test]
 fn test_stats_fasta_file() {
     assert!(Command::new("bash")