Commit 505abbd0 authored by Lukas Jelonek's avatar Lukas Jelonek
Browse files

Added a module that runs hmmscan vs PFAM-A

parent 867e3500
#!/usr/bin/python3
import sys
import json
import argparse
parser = argparse.ArgumentParser(description='Convert hmmscan results to json documents')
parser.add_argument('--result', '-r', required=True, help='The run_hmmer result directory')
parser.add_argument('--output', '-o', required=True, help='The converted results json file')
args = parser.parse_args()
filename = args.result + "/domtblout.tsv"
documents = {}
tool = {'name': None, 'version': None, 'database': None}
with open(filename) as f:
# first scan for tool information
for line in f:
if line.startswith("# Program:"):
tool['name'] = line.split()[2]
if line.startswith("# Version:"):
tool['version'] = line.split()[2]
if line.startswith("# Target file:"):
tool['database'] = line.split()[3]
with open(filename) as f:
# second scan for data
for line in f:
if not line.startswith("#"):
split = line.split(None, 22)
query_id = split[3]
if not query_id in documents:
documents[query_id] = {"id": query_id, "computations": [{'tool': tool, 'results':[]}]}
results = documents[query_id]['computations'][0]["results"]
results.append({
'target': {
'name': split[0],
'accession': split[1],
'length': split[2],
'score': float(split[13]),
'start': int(split[15]),
'end': int(split[16]),
'description': split[22].rstrip()
},
'query': {
'start': int(split[17]),
'end': int(split[19])
}
})
output_filename = args.output
with open(output_filename, 'w') as o:
json.dump(documents, o)
#!/usr/bin/env python3
import argparse
import config
from os import system, makedirs
hmmscan_tool = config.load_config()['tools'].get('hmmscan', 'hmmscan')
parser = argparse.ArgumentParser(description='Search sequences against a profile database')
parser.add_argument('--fasta', '-f', required=True, help='A fasta file with aminoacid sequences')
parser.add_argument('--database', '-d', required=True, help='Database to search in')
parser.add_argument('--output', '-o', required=True, help='The result directory')
parser.add_argument('--evalue', '-e', default='0.0001', help='Evalue cutoff')
args = parser.parse_args()
makedirs(args.output, exist_ok=True)
system(hmmscan_tool +
" -E " + args.evalue +
" -o " +args.output+ "/hmmscan.out " +
" --tblout " + args.output + "/tblout.tsv " +
" --domtblout " + args.output + "/domtblout.tsv " +
" --pfamtblout " + args.output + "/pfamtblout.tsv " +
args.database + " " + args.fasta)
......@@ -14,6 +14,7 @@ profiles:
modules:
signalp:
blastp_swissprot:
hmmer_pfam_a:
- name: 'complete'
info: 'Profile that uses all available tools'
modules:
......
# Module manifest for the blastp against swissprot analysis
# The name of the module. Is needed for the list-analyses option, for custom
# configurations and custom profiles
name: 'hmmer_pfam_a'
# Short description of the analysis
info: 'hmmscan analysis against PFAM-A'
# The name of the script for the analysis step. Must take a --fasta parameter
analysis_script: 'run_hmmer.py'
# The name of the result to json converter script. Must take one parameter, the
# result file from the analysis_script
converter_script: 'convert_hmmer.py'
parameters:
database: '/vol/biodb/pfam/Pfam-A.hmm'
evalue: 1e-10
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment