diff options
| -rw-r--r-- | bh20seqanalyzer/main.py | 40 | ||||
| -rw-r--r-- | bh20sequploader/main.py | 9 | ||||
| -rw-r--r-- | bh20sequploader/qc_fasta.py | 64 | ||||
| -rw-r--r-- | scripts/docker/Dockerfile | 2 |
4 files changed, 73 insertions, 42 deletions
diff --git a/bh20seqanalyzer/main.py b/bh20seqanalyzer/main.py index 2fefa86..1746587 100644 --- a/bh20seqanalyzer/main.py +++ b/bh20seqanalyzer/main.py @@ -235,9 +235,10 @@ def upload_schema(api, workflow_def_project): return "keep:%s/schema.yml" % pdh -def print_status(api, uploader_project): +def print_status(api, uploader_project, fmt): pending = arvados.util.list_all(api.collections().list, filters=[["owner_uuid", "=", uploader_project]]) out = [] + status = {} for p in pending: prop = p["properties"] out.append(prop) @@ -245,7 +246,38 @@ def print_status(api, uploader_project): prop["status"] = "pending" prop["created_at"] = p["created_at"] prop["uuid"] = p["uuid"] - print(json.dumps(out, indent=2)) + status[prop["status"]] = status.get(prop["status"], 0) + 1 + if fmt == "html": + print( +""" +<html> +<body> +""") + print("<p>Total collections in upload project %s</p>" % len(out)) + print("<p>Status %s</p>" % status) + print( +""" +<table> +<tr><th>Collection</th> +<th>Sequence label</th> +<th>Status</th> +<th>Errors</th></tr> +""") + for r in out: + print("<tr valign='top'>") + print("<td><a href='https://workbench.lugli.arvadosapi.com/collections/%s'>%s</a></td>" % (r["uuid"], r["uuid"])) + print("<td>%s</td>" % r["sequence_label"]) + print("<td>%s</td>" % r["status"]) + print("<td><pre>%s</pre></td>" % "\n".join(r.get("errors", []))) + print("</tr>") + print( +""" +</table> +</body> +</html> +""") + else: + print(json.dumps(out, indent=2)) def main(): parser = argparse.ArgumentParser(description='Analyze collections uploaded to a project') @@ -264,7 +296,7 @@ def main(): parser.add_argument('--kickoff', action="store_true") parser.add_argument('--no-start-analysis', action="store_true") parser.add_argument('--once', action="store_true") - parser.add_argument('--print-status', action="store_true") + parser.add_argument('--print-status', type=str, default=None) args = parser.parse_args() api = arvados.api() @@ -284,7 +316,7 @@ def main(): return if args.print_status: - print_status(api, args.uploader_project) + print_status(api, args.uploader_project, args.print_status) exit(0) logging.info("Starting up, monitoring %s for uploads" % (args.uploader_project)) diff --git a/bh20sequploader/main.py b/bh20sequploader/main.py index 8555e2b..fd0278d 100644 --- a/bh20sequploader/main.py +++ b/bh20sequploader/main.py @@ -23,14 +23,16 @@ ARVADOS_API_TOKEN='2fbebpmbo3rw3x05ueu2i6nx70zhrsb1p22ycu3ry34m4x4462' UPLOAD_PROJECT='lugli-j7d0g-n5clictpuvwk8aa' def qc_stuff(metadata, sequence_p1, sequence_p2, do_qc=True): + failed = False try: log.debug("Checking metadata" if do_qc else "Skipping metadata check") if do_qc and not qc_metadata(metadata.name): - raise Exception("Failed metadata qc") + log.warning("Failed metadata qc") + failed = True except Exception as e: log.debug(e) print(e) - exit(1) + failed = True target = [] try: @@ -43,6 +45,9 @@ def qc_stuff(metadata, sequence_p1, sequence_p2, do_qc=True): except Exception as e: log.debug(e) print(e) + failed = True + + if failed: exit(1) return target diff --git a/bh20sequploader/qc_fasta.py b/bh20sequploader/qc_fasta.py index 8c6ebd3..37eb4e8 100644 --- a/bh20sequploader/qc_fasta.py +++ b/bh20sequploader/qc_fasta.py @@ -58,42 +58,36 @@ def qc_fasta(arg_sequence, check_with_clustalw=True): return ("sequence.fasta"+gz, seqlabel) with tempfile.NamedTemporaryFile() as tmp1: - refstring = pkg_resources.resource_string(__name__, "SARS-CoV-2-reference.fasta") - tmp1.write(refstring) - tmp1.write(submitlabel.encode("utf8")) - tmp1.write(("".join(submitseq)).encode("utf8")) - tmp1.flush() - subbp = 0 - refbp = 0 - similarity = 0 - try: - cmd = ["clustalw", "-infile="+tmp1.name, - "-quicktree", "-iteration=none", "-type=DNA"] - print("QC checking similarity to reference") - print(" ".join(cmd)) - result = subprocess.run(cmd, stdout=subprocess.PIPE) - res = result.stdout.decode("utf-8") - g1 = re.search(r"^Sequence 1: [^ ]+ +(\d+) bp$", res, flags=re.MULTILINE) - refbp = float(g1.group(1)) - g2 = re.search(r"^Sequence 2: [^ ]+ +(\d+) bp$", res, flags=re.MULTILINE) - subbp = float(g2.group(1)) - g3 = re.search(r"^Sequences \(1:2\) Aligned\. Score: (\d+(\.\d+)?)$", res, flags=re.MULTILINE) - similarity = float(g3.group(1)) + with tempfile.NamedTemporaryFile() as tmp2: + refstring = pkg_resources.resource_string(__name__, "SARS-CoV-2-reference.fasta") + tmp1.write(refstring) + tmp1.flush() + tmp2.write(submitlabel.encode("utf8")) + tmp2.write(("".join(submitseq)).encode("utf8")) + tmp2.flush() + subbp = 0 + refbp = 0 + similarity = 0 + try: + cmd = ["minimap2", "-c", tmp1.name, tmp2.name] + logging.info("QC checking similarity to reference") + logging.info(" ".join(cmd)) + result = subprocess.run(cmd, stdout=subprocess.PIPE) + result.check_returncode() + res = result.stdout.decode("utf-8") + mm = res.split("\t") + if len(mm) >= 10: + # divide Number of matching bases in the mapping / Target sequence length + similarity = (float(mm[9]) / float(mm[6])) * 100.0 + else: + similarity = 0 + except Exception as e: + logging.warn("QC against reference sequence using 'minimap2': %s", e, exc_info=e) - print(g1.group(0)) - print(g2.group(0)) - print(g3.group(0)) - except Exception as e: - logging.warn("QC against reference sequence using 'clustalw': %s", e) - - if refbp and (subbp/refbp) < .7: - raise ValueError("QC fail: submit sequence length is shorter than 70% reference") - if refbp and (subbp/refbp) > 1.3: - raise ValueError("QC fail: submit sequence length is greater than 130% reference") - if similarity and similarity < 70.0: - raise ValueError("QC fail: submit similarity is less than 70%") - if refbp == 0 or similarity == 0: - raise ValueError("QC fail") + if similarity and similarity < 70.0: + raise ValueError("QC fail: alignment to reference was less than 70%% (was %2.2f%%)" % (similarity)) + if similarity == 0: + raise ValueError("QC fail") return ("sequence.fasta"+gz, seqlabel) elif seq_type == "text/fastq": diff --git a/scripts/docker/Dockerfile b/scripts/docker/Dockerfile index 8811927..02829d4 100644 --- a/scripts/docker/Dockerfile +++ b/scripts/docker/Dockerfile @@ -3,7 +3,7 @@ FROM debian:10 RUN apt-get update && \ apt-get -yq --no-install-recommends -o Acquire::Retries=6 install \ python3 python3-pip python3-setuptools python3-dev python-pycurl \ - clustalw python3-biopython libcurl4-openssl-dev build-essential \ + minimap2 python3-biopython libcurl4-openssl-dev build-essential \ libssl-dev libmagic-dev python3-magic && \ apt-get clean |