From d29dfd593233541b85c1cefb239650279d57d59f Mon Sep 17 00:00:00 2001
From: Peter Amstutz
Date: Mon, 20 Apr 2020 13:41:56 -0400
Subject: Relabel sequences to match metadata subjects.

Arvados-DCO-1.1-Signed-off-by: Peter Amstutz <peter.amstutz@curii.com>
---
 workflows/pangenome-generate/relabel-seqs.cwl | 19 +++++++++++++++++++
 1 file changed, 19 insertions(+)
 create mode 100644 workflows/pangenome-generate/relabel-seqs.cwl

(limited to 'workflows/pangenome-generate/relabel-seqs.cwl')

diff --git a/workflows/pangenome-generate/relabel-seqs.cwl b/workflows/pangenome-generate/relabel-seqs.cwl
new file mode 100644
index 0000000..b5b7231
--- /dev/null
+++ b/workflows/pangenome-generate/relabel-seqs.cwl
@@ -0,0 +1,19 @@
+cwlVersion: v1.1
+class: CommandLineTool
+inputs:
+  readsFA: File[]
+  subjects: string[]
+outputs:
+  relabeledSeqs:
+    type: stdout
+requirements:
+  InlineJavascriptRequirement: {}
+  InitialWorkDirRequirement:
+    listing:
+      - entry: {$include: relabel-seqs.py}
+        entryname: relabel-seqs.py
+hints:
+  DockerRequirement:
+    dockerPull: commonworkflowlanguage/cwltool_module
+stdout: relabeledSeqs.fasta
+baseCommand: [python, relabel-seqs.py]
-- 
cgit v1.2.3


From 9ddcfeacb3191638f42b08af999889d867f0f81c Mon Sep 17 00:00:00 2001
From: Peter Amstutz
Date: Mon, 20 Apr 2020 14:57:25 -0400
Subject: Better handling of duplicate sequences

Also save original fasta label in metadata
---
 bh20sequploader/bh20seq-schema.yml                  | 11 ++++++++++-
 workflows/pangenome-generate/merge-metadata.cwl     |  2 ++
 workflows/pangenome-generate/merge-metadata.py      | 21 ++++++++++++++++++++-
 workflows/pangenome-generate/pangenome-generate.cwl | 10 ++++------
 workflows/pangenome-generate/relabel-seqs.cwl       | 10 ++++++++--
 workflows/pangenome-generate/relabel-seqs.py        | 12 +++++++-----
 workflows/pangenome-generate/seqkit-rmdup.cwl       |  4 ++--
 workflows/pangenome-generate/testjob.yml            | 16 ++++++++++++++++
 8 files changed, 69 insertions(+), 17 deletions(-)
 create mode 100644 workflows/pangenome-generate/testjob.yml

(limited to 'workflows/pangenome-generate/relabel-seqs.cwl')

diff --git a/bh20sequploader/bh20seq-schema.yml b/bh20sequploader/bh20seq-schema.yml
index 64008f2..982447c 100644
--- a/bh20sequploader/bh20seq-schema.yml
+++ b/bh20sequploader/bh20seq-schema.yml
@@ -18,6 +18,7 @@ $graph:
         jsonldPredicate:
           _id: http://www.ebi.ac.uk/efo/EFO_0000532
           _type: "@id"
+          identity: true
     host_id:
         doc: Identifer for the host. If you submit multiple samples from the same host, use the same host_id for those samples
         type: string
@@ -29,6 +30,7 @@ $graph:
         jsonldPredicate:
           _id: http://purl.obolibrary.org/obo/PATO_0000047
           _type: "@id"
+          identity: true
     host_age:
         doc: Age of the host as number (e.g. 50)
         type: int?
@@ -40,6 +42,7 @@ $graph:
         jsonldPredicate:
           _id: http://purl.obolibrary.org/obo/NCIT_C42574
           _type: "@id"
+          identity: true
     host_health_status:
         doc: A condition or state at a particular time
         type: string?
@@ -79,12 +82,14 @@ $graph:
       jsonldPredicate:
           _id: http://purl.obolibrary.org/obo/OBI_0001479
           _type: "@id"
+          identity: true
     specimen_source2:
       doc: Method how the specimen was derived as NCIT IRI, e.g. http://purl.obolibrary.org/obo/NCIT_C155835 (=throat swabb)
       type: string?
       jsonldPredicate:
           _id: http://purl.obolibrary.org/obo/OBI_0001479
           _type: "@id"
+          identity: true
     collection_date:
       doc: Date when the sample was taken
       type: string
@@ -96,6 +101,7 @@ $graph:
       jsonldPredicate:
         _id: http://purl.obolibrary.org/obo/GAZ_00000448
         _type: "@id"
+        identity: true
     sample_storage_conditions:
       doc: Information about storage of a specified type, e.g.  frozen specimen, paraffin, fresh ....
       type: string?
@@ -126,6 +132,7 @@ $graph:
       jsonldPredicate:
           _id: http://edamontology.org/data_1875
           _type: "@id"
+          identity: true
     virus_strain:
       doc: Name of the virus strain
       type: string?
@@ -141,12 +148,14 @@ $graph:
       jsonldPredicate:
         _id: http://purl.obolibrary.org/obo/OBI_0600047
         _type: "@id"
+        identity: true
     sample_sequencing_technology2:
       doc: Technology that was used to sequence this sample (e.g Sanger, Nanopor MiniION)
       type: string?
       jsonldPredicate:
         _id: http://purl.obolibrary.org/obo/OBI_0600047
         _type: "@id"
+        identity: true
     sequence_assembly_method:
       doc: Protocol which provides instructions on the alignment of sequencing reads to reference genome
       type: string?
@@ -215,7 +224,7 @@ $graph:
       jsonldPredicate:
           _id: http://semanticscience.org/resource/SIO_000115
           _type: "@id"
-          noLinkCheck: true
+          identity: true
 
 - name: MainSchema
   type: record
diff --git a/workflows/pangenome-generate/merge-metadata.cwl b/workflows/pangenome-generate/merge-metadata.cwl
index 9164c09..fcefe32 100644
--- a/workflows/pangenome-generate/merge-metadata.cwl
+++ b/workflows/pangenome-generate/merge-metadata.cwl
@@ -7,6 +7,8 @@ inputs:
   metadata: File[]
   metadataSchema: File
   subjects: string[]
+  dups: File?
+  originalLabels: File
 outputs:
   merged: stdout
 stdout: mergedmetadata.ttl
diff --git a/workflows/pangenome-generate/merge-metadata.py b/workflows/pangenome-generate/merge-metadata.py
index 64275b1..bfec781 100644
--- a/workflows/pangenome-generate/merge-metadata.py
+++ b/workflows/pangenome-generate/merge-metadata.py
@@ -1,9 +1,13 @@
+import re
 import schema_salad.schema
 import schema_salad.jsonld_context
+import json
 
 metadataSchema = '$(inputs.metadataSchema.path)'
 metadata = $(inputs.metadata)
 subjects = $(inputs.subjects)
+dups = json.loads('''$(inputs.dups)''')
+originalLabels = $(inputs.originalLabels)
 
 (document_loader,
  avsc_names,
@@ -11,7 +15,22 @@ subjects = $(inputs.subjects)
  metaschema_loader) = schema_salad.schema.load_schema(metadataSchema)
 
 for i, m in enumerate(metadata):
-    doc, metadata = schema_salad.schema.load_and_validate(document_loader, avsc_names, m["path"], True)
+    doc, metadata = schema_salad.schema.load_and_validate(document_loader, avsc_names, m["path"], False, False)
     doc["id"] = subjects[i]
     g = schema_salad.jsonld_context.makerdf(subjects[i], doc, document_loader.ctx)
     print(g.serialize(format="ntriples").decode("utf-8"))
+
+import logging
+
+if dups:
+    sameseqs = open(dups["path"], "rt")
+    for d in sameseqs:
+        logging.warn(d)
+        g = re.match(r"\\d+\\t(.*)", d)
+        logging.warn("%s", g.group(1))
+        sp = g.group(1).split(",")
+        for n in sp[1:]:
+            print("<%s> <http://biohackathon.org/bh20-seq-schema/has_duplicate_sequence> <%s> ." % (n.strip(), sp[0].strip()))
+
+orig = open(originalLabels["path"], "rt")
+print(orig.read())
diff --git a/workflows/pangenome-generate/pangenome-generate.cwl b/workflows/pangenome-generate/pangenome-generate.cwl
index 896f936..0cb1368 100644
--- a/workflows/pangenome-generate/pangenome-generate.cwl
+++ b/workflows/pangenome-generate/pangenome-generate.cwl
@@ -26,15 +26,11 @@ steps:
     in:
       readsFA: inputReads
       subjects: subjects
-    out: [relabeledSeqs]
+    out: [relabeledSeqs, originalLabels]
     run: relabel-seqs.cwl
-  common:
-    in: {readsFA: relabel/relabeledSeqs}
-    out: [duplicatedReads]
-    run: seqkit-common.cwl
   dedup:
     in: {readsFA: relabel/relabeledSeqs}
-    out: [readsMergeDedup]
+    out: [readsMergeDedup, dups]
     run: seqkit-rmdup.cwl
   overlapReads:
     in: {readsFA: dedup/readsMergeDedup}
@@ -63,5 +59,7 @@ steps:
       metadata: metadata
       metadataSchema: metadataSchema
       subjects: subjects
+      dups: dedup/dups
+      originalLabels: relabel/originalLabels
     out: [merged]
     run: merge-metadata.cwl
diff --git a/workflows/pangenome-generate/relabel-seqs.cwl b/workflows/pangenome-generate/relabel-seqs.cwl
index b5b7231..2b780d4 100644
--- a/workflows/pangenome-generate/relabel-seqs.cwl
+++ b/workflows/pangenome-generate/relabel-seqs.cwl
@@ -5,7 +5,13 @@ inputs:
   subjects: string[]
 outputs:
   relabeledSeqs:
-    type: stdout
+    type: File
+    outputBinding:
+      glob: relabeledSeqs.fasta
+  originalLabels:
+    type: File
+    outputBinding:
+      glob: originalLabels.ttl
 requirements:
   InlineJavascriptRequirement: {}
   InitialWorkDirRequirement:
@@ -15,5 +21,5 @@ requirements:
 hints:
   DockerRequirement:
     dockerPull: commonworkflowlanguage/cwltool_module
-stdout: relabeledSeqs.fasta
+stdout:
 baseCommand: [python, relabel-seqs.py]
diff --git a/workflows/pangenome-generate/relabel-seqs.py b/workflows/pangenome-generate/relabel-seqs.py
index 32f2386..b558fe2 100644
--- a/workflows/pangenome-generate/relabel-seqs.py
+++ b/workflows/pangenome-generate/relabel-seqs.py
@@ -1,13 +1,15 @@
-import sys
-
 reads = $(inputs.readsFA)
 subjects = $(inputs.subjects)
 
+relabeled_fasta = open("relabeledSeqs.fasta", "wt")
+original_labels = open("originalLabels.ttl", "wt")
+
 for i, r in enumerate(reads):
     with open(r["path"], "rt") as fa:
-        fa.readline()
-        print(">"+subjects[i])
+        label = fa.readline()
+        original_labels.write("<%s> <http://biohackathon.org/bh20-seq-schema/original_fasta_label> \\"%s\\" .\\n" % (subjects[i], label[1:].strip().replace('"', '\\\\"')))
+        relabeled_fasta.write(">"+subjects[i]+"\\n")
         data = fa.read(8096)
         while data:
-            sys.stdout.write(data)
+            relabeled_fasta.write(data)
             data = fa.read(8096)
diff --git a/workflows/pangenome-generate/seqkit-rmdup.cwl b/workflows/pangenome-generate/seqkit-rmdup.cwl
index 07184c3..071fa66 100644
--- a/workflows/pangenome-generate/seqkit-rmdup.cwl
+++ b/workflows/pangenome-generate/seqkit-rmdup.cwl
@@ -1,14 +1,14 @@
 cwlVersion: v1.1
 class: CommandLineTool
 inputs:
-  readsFA: File[]
+  readsFA: File
 outputs:
   readsMergeDedup:
     type: File
     outputBinding:
       glob: readsMergeDedup.fasta
   dups:
-    type: File
+    type: File?
     outputBinding:
       glob: dups.txt
 requirements:
diff --git a/workflows/pangenome-generate/testjob.yml b/workflows/pangenome-generate/testjob.yml
new file mode 100644
index 0000000..a48aff8
--- /dev/null
+++ b/workflows/pangenome-generate/testjob.yml
@@ -0,0 +1,16 @@
+inputReads:
+  - class: File
+    location: ../../example/sequence.fasta
+  - class: File
+    location: ../../example/sequence.fasta
+metadata:
+  - class: File
+    location: ../../example/metadata.yaml
+  - class: File
+    location: ../../example/metadata.yaml
+metadataSchema:
+  class: File
+  location: ../../bh20sequploader/bh20seq-schema.yml
+subjects:
+  - http://arvados.org/keep/seq1
+  - http://arvados.org/keep/seq2
-- 
cgit v1.2.3


From f4c3da88c1233802fea46cc972a81dc3b5b51185 Mon Sep 17 00:00:00 2001
From: Peter Amstutz
Date: Tue, 21 Apr 2020 15:37:58 -0400
Subject: Work around CWL content size limit by chunking

Arvados-DCO-1.1-Signed-off-by: Peter Amstutz <peter.amstutz@curii.com>
---
 bh20sequploader/main.py                       |  1 +
 workflows/pangenome-generate/relabel-seqs.cwl | 31 +++++++++++++++++++++++----
 workflows/pangenome-generate/relabel-seqs.py  | 22 +++++++++++++------
 3 files changed, 44 insertions(+), 10 deletions(-)

(limited to 'workflows/pangenome-generate/relabel-seqs.cwl')

diff --git a/bh20sequploader/main.py b/bh20sequploader/main.py
index 2fda347..4c4711d 100644
--- a/bh20sequploader/main.py
+++ b/bh20sequploader/main.py
@@ -63,6 +63,7 @@ def main():
     external_ip = urllib.request.urlopen('https://ident.me').read().decode('utf8')
 
     properties = {
+        "sequence_label": seqlabel,
         "upload_app": "bh20-seq-uploader",
         "upload_ip": external_ip,
         "upload_user": "%s@%s" % (getpass.getuser(), socket.gethostname())
diff --git a/workflows/pangenome-generate/relabel-seqs.cwl b/workflows/pangenome-generate/relabel-seqs.cwl
index 2b780d4..01196f6 100644
--- a/workflows/pangenome-generate/relabel-seqs.cwl
+++ b/workflows/pangenome-generate/relabel-seqs.cwl
@@ -3,6 +3,10 @@ class: CommandLineTool
 inputs:
   readsFA: File[]
   subjects: string[]
+  script:
+    type: File
+    default: {class: File, location: relabel-seqs.py}
+    inputBinding: {}
 outputs:
   relabeledSeqs:
     type: File
@@ -15,11 +19,30 @@ outputs:
 requirements:
   InlineJavascriptRequirement: {}
   InitialWorkDirRequirement:
-    listing:
-      - entry: {$include: relabel-seqs.py}
-        entryname: relabel-seqs.py
+    listing: |
+          ${
+          var i = 0;
+          var b = 1;
+          var out = [];
+          for (; i < inputs.readsFA.length; i++) {
+            var block = [];
+            for (; i < (b*100) && i < inputs.readsFA.length; i++) {
+              block.push(inputs.readsFA[i]);
+            }
+            out.push({
+              entryname: "block"+b,
+              entry: JSON.stringify(block)
+            });
+            b++;
+          }
+          out.push({
+            entry: JSON.stringify(inputs.subjects),
+            entryname: "subjects"
+          });
+          return out;
+          }
 hints:
   DockerRequirement:
     dockerPull: commonworkflowlanguage/cwltool_module
 stdout:
-baseCommand: [python, relabel-seqs.py]
+baseCommand: [python]
diff --git a/workflows/pangenome-generate/relabel-seqs.py b/workflows/pangenome-generate/relabel-seqs.py
index 1188ceb..970540f 100644
--- a/workflows/pangenome-generate/relabel-seqs.py
+++ b/workflows/pangenome-generate/relabel-seqs.py
@@ -1,5 +1,15 @@
-reads = $(inputs.readsFA)
-subjects = $(inputs.subjects)
+import os
+import json
+
+reads = []
+b = 1
+while os.path.exists("block%i" % b):
+    with open("block%i" % b) as f:
+        reads.extend(json.load(f))
+    b += 1
+
+with open("subjects") as f:
+    subjects = json.load(f)
 
 relabeled_fasta = open("relabeledSeqs.fasta", "wt")
 original_labels = open("originalLabels.ttl", "wt")
@@ -7,12 +17,12 @@ original_labels = open("originalLabels.ttl", "wt")
 for i, r in enumerate(reads):
     with open(r["path"], "rt") as fa:
         label = fa.readline()
-        original_labels.write("<%s> <http://biohackathon.org/bh20-seq-schema/original_fasta_label> \\"%s\\" .\\n" % (subjects[i], label[1:].strip().replace('"', '\\\\"')))
-        relabeled_fasta.write(">"+subjects[i]+"\\n")
+        original_labels.write("<%s> <http://biohackathon.org/bh20-seq-schema/original_fasta_label> \"%s\" .\n" % (subjects[i], label[1:].strip().replace('"', '\\"')))
+        relabeled_fasta.write(">"+subjects[i]+"\n")
         data = fa.read(8096)
         while data:
             relabeled_fasta.write(data)
-            endswithnewline = data.endswith("\\n")
+            endswithnewline = data.endswith("\n")
             data = fa.read(8096)
         if not endswithnewline:
-            relabeled_fasta.write("\\n")
+            relabeled_fasta.write("\n")
-- 
cgit v1.2.3


From 61726edb9293fe529e6efbe5bb6f1cc953bb3c4e Mon Sep 17 00:00:00 2001
From: Peter Amstutz
Date: Tue, 21 Apr 2020 16:20:10 -0400
Subject: Workaround CWL limit by chunking file list

Arvados-DCO-1.1-Signed-off-by: Peter Amstutz <peter.amstutz@curii.com>
---
 workflows/pangenome-generate/merge-metadata.cwl | 46 +++++++++++++++++++++----
 workflows/pangenome-generate/merge-metadata.py  | 33 ++++++++++++------
 workflows/pangenome-generate/relabel-seqs.cwl   | 12 ++++---
 workflows/pangenome-generate/relabel-seqs.py    | 18 +++++-----
 4 files changed, 79 insertions(+), 30 deletions(-)

(limited to 'workflows/pangenome-generate/relabel-seqs.cwl')

diff --git a/workflows/pangenome-generate/merge-metadata.cwl b/workflows/pangenome-generate/merge-metadata.cwl
index fcefe32..4d9c808 100644
--- a/workflows/pangenome-generate/merge-metadata.cwl
+++ b/workflows/pangenome-generate/merge-metadata.cwl
@@ -5,16 +5,48 @@ hints:
     dockerPull: commonworkflowlanguage/cwltool_module
 inputs:
   metadata: File[]
-  metadataSchema: File
   subjects: string[]
-  dups: File?
-  originalLabels: File
+  metadataSchema:
+    type: File
+    inputBinding: {position: 2}
+  originalLabels:
+    type: File
+    inputBinding: {position: 3}
+  dups:
+    type: File?
+    inputBinding: {position: 4}
+  script:
+    type: File
+    inputBinding: {position: 1}
+    default: {class: File, location: merge-metadata.py}
 outputs:
   merged: stdout
 stdout: mergedmetadata.ttl
 requirements:
+  InlineJavascriptRequirement: {}
   InitialWorkDirRequirement:
-    listing:
-      - entry: {$include: merge-metadata.py}
-        entryname: merge-metadata.py
-baseCommand: [python3, merge-metadata.py]
+    listing: |
+          ${
+          var i = 0;
+          var b = 1;
+          var out = [];
+          for (; i < inputs.metadata.length; i++) {
+            var block = [];
+            var sub = [];
+            for (; i < (b*150) && i < inputs.metadata.length; i++) {
+              block.push(inputs.metadata[i]);
+              sub.push(inputs.subjects[i]);
+            }
+            out.push({
+              entryname: "block"+b,
+              entry: JSON.stringify(block)
+            });
+            out.push({
+              entryname: "subs"+b,
+              entry: JSON.stringify(sub)
+            });
+            b++;
+          }
+          return out;
+          }
+baseCommand: python
diff --git a/workflows/pangenome-generate/merge-metadata.py b/workflows/pangenome-generate/merge-metadata.py
index bfec781..65d08a6 100644
--- a/workflows/pangenome-generate/merge-metadata.py
+++ b/workflows/pangenome-generate/merge-metadata.py
@@ -2,12 +2,27 @@ import re
 import schema_salad.schema
 import schema_salad.jsonld_context
 import json
+import sys
+import os
+import logging
+
+metadataSchema = sys.argv[1]
+originalLabels = sys.argv[2]
+dups = None
+if len(sys.argv) == 4:
+    dups = sys.argv[3]
+
+def readitems(stem):
+    items = []
+    b = 1
+    while os.path.exists("%s%i" % (stem, b)):
+        with open("%s%i" % (stem, b)) as f:
+            items.extend(json.load(f))
+        b += 1
+    return items
 
-metadataSchema = '$(inputs.metadataSchema.path)'
-metadata = $(inputs.metadata)
-subjects = $(inputs.subjects)
-dups = json.loads('''$(inputs.dups)''')
-originalLabels = $(inputs.originalLabels)
+metadata = readitems("block")
+subjects = readitems("subs")
 
 (document_loader,
  avsc_names,
@@ -20,17 +35,15 @@ for i, m in enumerate(metadata):
     g = schema_salad.jsonld_context.makerdf(subjects[i], doc, document_loader.ctx)
     print(g.serialize(format="ntriples").decode("utf-8"))
 
-import logging
-
 if dups:
-    sameseqs = open(dups["path"], "rt")
+    sameseqs = open(dups, "rt")
     for d in sameseqs:
         logging.warn(d)
-        g = re.match(r"\\d+\\t(.*)", d)
+        g = re.match(r"\d+\t(.*)", d)
         logging.warn("%s", g.group(1))
         sp = g.group(1).split(",")
         for n in sp[1:]:
             print("<%s> <http://biohackathon.org/bh20-seq-schema/has_duplicate_sequence> <%s> ." % (n.strip(), sp[0].strip()))
 
-orig = open(originalLabels["path"], "rt")
+orig = open(originalLabels, "rt")
 print(orig.read())
diff --git a/workflows/pangenome-generate/relabel-seqs.cwl b/workflows/pangenome-generate/relabel-seqs.cwl
index 01196f6..c1f17a4 100644
--- a/workflows/pangenome-generate/relabel-seqs.cwl
+++ b/workflows/pangenome-generate/relabel-seqs.cwl
@@ -26,19 +26,21 @@ requirements:
           var out = [];
           for (; i < inputs.readsFA.length; i++) {
             var block = [];
-            for (; i < (b*100) && i < inputs.readsFA.length; i++) {
+            var sub = [];
+            for (; i < (b*150) && i < inputs.readsFA.length; i++) {
               block.push(inputs.readsFA[i]);
+              sub.push(inputs.subjects[i]);
             }
             out.push({
               entryname: "block"+b,
               entry: JSON.stringify(block)
             });
+            out.push({
+              entryname: "subs"+b,
+              entry: JSON.stringify(sub)
+            });
             b++;
           }
-          out.push({
-            entry: JSON.stringify(inputs.subjects),
-            entryname: "subjects"
-          });
           return out;
           }
 hints:
diff --git a/workflows/pangenome-generate/relabel-seqs.py b/workflows/pangenome-generate/relabel-seqs.py
index 970540f..6b022a0 100644
--- a/workflows/pangenome-generate/relabel-seqs.py
+++ b/workflows/pangenome-generate/relabel-seqs.py
@@ -1,15 +1,17 @@
 import os
 import json
 
-reads = []
-b = 1
-while os.path.exists("block%i" % b):
-    with open("block%i" % b) as f:
-        reads.extend(json.load(f))
-    b += 1
+def readitems(stem):
+    items = []
+    b = 1
+    while os.path.exists("%s%i" % (stem, b)):
+        with open("%s%i" % (stem, b)) as f:
+            items.extend(json.load(f))
+        b += 1
+    return items
 
-with open("subjects") as f:
-    subjects = json.load(f)
+reads = readitems("block")
+subjects = readitems("subs")
 
 relabeled_fasta = open("relabeledSeqs.fasta", "wt")
 original_labels = open("originalLabels.ttl", "wt")
-- 
cgit v1.2.3