various changes
This commit is contained in:
parent
80e00b0dfa
commit
93f6707b8f
12 changed files with 159 additions and 231 deletions
|
|
@ -1,6 +1,7 @@
|
|||
#!/usr/bin/env python3
|
||||
import subprocess
|
||||
import os
|
||||
import sys
|
||||
import numpy as np
|
||||
import pandas as pd
|
||||
from contextlib import suppress
|
||||
|
|
@ -8,6 +9,8 @@ from pathlib import Path
|
|||
import re
|
||||
import csv
|
||||
import argparse
|
||||
import shutil
|
||||
import time
|
||||
#https://realpython.com/python-pathlib/
|
||||
|
||||
# FIXME
|
||||
|
|
@ -22,8 +25,8 @@ homedir = os.path.expanduser('~')
|
|||
|
||||
# set working dir
|
||||
os.getcwd()
|
||||
os.chdir(homedir + '/git/LSHTM_analysis/foldx/')
|
||||
os.getcwd()
|
||||
#os.chdir(homedir + '/git/LSHTM_analysis/foldx/')
|
||||
#os.getcwd()
|
||||
|
||||
#=======================================================================
|
||||
#%% command line args
|
||||
|
|
@ -37,7 +40,7 @@ arg_parser.add_argument('-i', '--input_dir', help = 'Input dir containing pdb fi
|
|||
arg_parser.add_argument('-o', '--output_dir', help = 'Output dir for results. By default, it assmes homedir + <drug> + output')
|
||||
arg_parser.add_argument('-p', '--process_dir', help = 'Temp processing dir for running foldX. By default, it assmes homedir + <drug> + processing. Make sure it is somewhere with LOTS of storage as it writes all output!') #FIXME
|
||||
|
||||
arg_parser.add_argument('-pdb', '--pdb_file', help = 'PDB File to process. By default, it assmumes a file called <gene>_complex.pdb in input_dir')
|
||||
arg_parser.add_argument('-pdb', '--pdb_file', help = 'PDB File to process. By default, it assmumes a file called <gene>_complex.pdb in input_dir') # DO NOT specify an absolute path
|
||||
arg_parser.add_argument('-m', '--mutation_file', help = 'Mutation list. By default, assumes a file called <gene>_mcsm_snps.csv exists')
|
||||
|
||||
# FIXME: Doesn't work with 2 chains yet!
|
||||
|
|
@ -77,34 +80,45 @@ if not datadir:
|
|||
datadir = homedir + '/' + 'git/Data'
|
||||
|
||||
if not indir:
|
||||
indir = datadir + '/' + drug + '/input'
|
||||
indir =
|
||||
datadir + '/' + drug + '/input'
|
||||
|
||||
if not outdir:
|
||||
outdir = datadir + '/' + drug + '/output'
|
||||
|
||||
#TODO: perhaps better handled by refactoring code to prevent generating lots of output files!
|
||||
if not process_dir:
|
||||
process_dir = datadir + '/' + drug +'/' + 'processing'
|
||||
process_dir = datadir + '/' + drug + '/processing'
|
||||
|
||||
# Make all paths absolute in case the user forgot
|
||||
indir = os.path.abspath(indir)
|
||||
process_dir = os.path.abspath(process_dir)
|
||||
outdir = os.path.abspath(outdir)
|
||||
datadir = os.path.abspath(datadir)
|
||||
|
||||
#=======
|
||||
# input
|
||||
#=======
|
||||
# FIXME
|
||||
if pdb_filename:
|
||||
pdb_filename = os.path.abspath(pdb_filename)
|
||||
pdb_name = Path(pdb_filename).stem
|
||||
infile_pdb = pdb_filename
|
||||
else:
|
||||
pdb_filename = gene.lower() + '_complex.pdb'
|
||||
pdb_name = Path(pdb_filename).stem
|
||||
|
||||
infile_pdb = indir + '/' + pdb_filename
|
||||
infile_pdb = indir + '/' + pdb_filename
|
||||
|
||||
actual_pdb_filename = Path(infile_pdb).name
|
||||
|
||||
if mut_filename:
|
||||
mutation_file = mut_filename
|
||||
mutation_file = os.path.abspath(mut_filename)
|
||||
infile_muts = mutation_file
|
||||
print('User-provided mutation file in use:', infile_muts)
|
||||
else:
|
||||
mutation_file = gene.lower() + '_mcsm_formatted_snps.csv'
|
||||
|
||||
infile_muts = outdir + '/' + mutation_file
|
||||
infile_muts = outdir + '/' + mutation_file
|
||||
print('WARNING: Assuming default mutation file:', infile_muts)
|
||||
|
||||
#=======
|
||||
# output
|
||||
|
|
@ -116,6 +130,7 @@ print('Arguments being passed:'
|
|||
, '\nDrug:', args.drug
|
||||
, '\ngene:', args.gene
|
||||
, '\ninput dir:', indir
|
||||
, '\nprocess dir:', process_dir
|
||||
, '\noutput dir:', outdir
|
||||
, '\npdb file:', infile_pdb
|
||||
, '\npdb name:', pdb_name
|
||||
|
|
@ -124,6 +139,9 @@ print('Arguments being passed:'
|
|||
, '\nchain1:', args.chain1
|
||||
, '\noutput file:', outfile_foldx
|
||||
, '\n=============================================================')
|
||||
|
||||
#### Delay for 10 seconds to check the params ####
|
||||
time.sleep(10)
|
||||
#=======================================================================
|
||||
|
||||
def getInteractionEnergy(filename):
|
||||
|
|
@ -184,6 +202,19 @@ def loadFiles(df):
|
|||
f.close()
|
||||
return np.asarray(resultList, dtype=np.float32)
|
||||
|
||||
# TODO: put the subprocess call in a 'def'
|
||||
#def repairPDB():
|
||||
# subprocess.call(['foldx'
|
||||
# , '--command=RepairPDB'
|
||||
# , '--pdb-dir=' + indir
|
||||
# , '--pdb=' + actual_pdb_filename
|
||||
# , '--ionStrength=0.05'#
|
||||
# , '--pH=7'
|
||||
# , '--water=PREDICT'
|
||||
# , '--vdwDesign=1'
|
||||
# , 'outPDB=true'
|
||||
# , '--output-dir=' + process_dir])
|
||||
|
||||
#=======================================================================
|
||||
def main():
|
||||
pdbname = pdb_name
|
||||
|
|
@ -195,37 +226,128 @@ def main():
|
|||
nmuts = len(mutlist)
|
||||
print(nmuts)
|
||||
print(mutlist)
|
||||
print('start')
|
||||
#subprocess.check_output(['bash','repairPDB.sh', pdbname, process_dir])
|
||||
subprocess.check_output(['bash','repairPDB.sh', indir, actual_pdb_filename, process_dir])
|
||||
|
||||
print('end')
|
||||
output = subprocess.check_output(['bash', 'runfoldx.sh', pdbname, process_dir])
|
||||
print('start')
|
||||
# some common parameters for foldX
|
||||
foldx_common=' --ionStrength=0.05 --pH=7 --water=PREDICT --vdwDesign=1 '
|
||||
|
||||
print('\033[95mSTAGE: repair PDB (foldx subprocess) \033[0m')
|
||||
print('Running foldx RepairPDB for WT')
|
||||
subprocess.call(['foldx'
|
||||
, '--command=RepairPDB'
|
||||
, foldx_common
|
||||
, '--pdb-dir=' + indir
|
||||
, '--pdb=' + actual_pdb_filename
|
||||
, 'outPDB=true'
|
||||
, '--output-dir=' + process_dir])
|
||||
print('\033[95mCOMPLETED STAGE: repair PDB\033[0m')
|
||||
print('\n==========================================================')
|
||||
|
||||
|
||||
print('\033[95mSTAGE: Foldx commands BM, PN and SD (foldx subprocess) for WT\033[0m')
|
||||
print('Running foldx BuildModel for WT')
|
||||
subprocess.call(['foldx'
|
||||
, '--command=BuildModel'
|
||||
, foldx_common
|
||||
, '--pdb-dir=' + process_dir
|
||||
, '--pdb=' + pdbname + '_Repair.pdb'
|
||||
, '--mutant-file="individual_list_' + pdbname +'.txt"'
|
||||
, 'outPDB=true'
|
||||
, '--numberOfRuns=1'
|
||||
, '--output-dir=' + process_dir], cwd=process_dir)
|
||||
|
||||
print('Running foldx PrintNetworks for WT')
|
||||
subprocess.call(['foldx'
|
||||
, '--command=PrintNetworks'
|
||||
, '--pdb-dir=' + process_dir
|
||||
, '--pdb=' + pdbname + '_Repair.pdb'
|
||||
, '--water=PREDICT'
|
||||
, '--vdwDesign=1'
|
||||
, '--output-dir=' + process_dir], cwd=process_dir)
|
||||
|
||||
print('Running foldx SequenceDetail for WT')
|
||||
subprocess.call(['foldx'
|
||||
, '--command=SequenceDetail'
|
||||
, '--pdb-dir=' + process_dir
|
||||
, '--pdb=' + pdbname + '_Repair.pdb'
|
||||
, '--water=PREDICT'
|
||||
, '--vdwDesign=1'
|
||||
, '--output-dir=' + process_dir], cwd=process_dir)
|
||||
print('\033[95mCOMPLETED STAGE: Foldx commands BM, PN and SD\033[0m')
|
||||
print('\n==========================================================')
|
||||
|
||||
|
||||
print('\033[95mSTAGE: Print Networks (foldx subprocess) for MT\033[0m')
|
||||
for n in range(1,nmuts+1):
|
||||
print(n)
|
||||
with suppress(Exception):
|
||||
subprocess.check_output(['bash', 'runPrintNetworks.sh', pdbname, str(n), process_dir])
|
||||
|
||||
print('\033[95mNETWORK:\033[0m', n)
|
||||
print('Running foldx PrintNetworks for mutation', n)
|
||||
subprocess.call(['foldx'
|
||||
, '--command=PrintNetworks'
|
||||
, '--pdb-dir=' + process_dir
|
||||
, '--pdb=' + pdbname + '_Repair_' + str(n) + '.pdb'
|
||||
, '--water=PREDICT'
|
||||
, '--vdwDesign=1'
|
||||
, '--output-dir=' + process_dir], cwd=process_dir)
|
||||
print('\033[95mCOMPLETED STAGE: Print Networks (foldx subprocess) for MT\033[0m')
|
||||
print('\n==========================================================')
|
||||
|
||||
|
||||
print('\033[95mSTAGE: Rename Mutation Files (shell)\033[0m')
|
||||
for n in range(1,nmuts+1):
|
||||
print(n)
|
||||
print('\033[95mMUTATION:\033[0m', n)
|
||||
print('\033[96mCommand:\033[0m mutrenamefiles.sh %s %s %s' % (pdbname, str(n), process_dir ))
|
||||
#FIXME: bad design and needs to be done in a pythonic way
|
||||
with suppress(Exception):
|
||||
subprocess.check_output(['bash', 'mutrenamefiles.sh', pdbname, str(n), process_dir])
|
||||
|
||||
print('\033[95mCOMPLETED STAGE: Rename Mutation Files (shell)\033[0m')
|
||||
print('\n==========================================================')
|
||||
|
||||
|
||||
print('\033[95mSTAGE: Rename Files (shell) for WT\033[0m')
|
||||
# FIXME: this is bad design and needs to be done in a pythonic way
|
||||
out = subprocess.check_output(['bash','renamefiles.sh', pdbname, process_dir])
|
||||
|
||||
print('\033[95mCOMPLETED STAGE: Rename Files (shell) for WT\033[0m')
|
||||
print('\n==========================================================')
|
||||
|
||||
|
||||
if comp=='y':
|
||||
print('\033[95mSTAGE: Running foldx AnalyseComplex (foldx subprocess) for WT\033[0m')
|
||||
chain1=chainA
|
||||
chain2=chainB
|
||||
with suppress(Exception):
|
||||
subprocess.check_output(['bash','runcomplex.sh', pdbname, chain1, chain2, process_dir])
|
||||
for n in range(1,nmuts+1):
|
||||
with suppress(Exception):
|
||||
subprocess.check_output(['bash','mutruncomplex.sh', pdbname, chain1, chain2, str(n), process_dir])
|
||||
subprocess.call(['foldx'
|
||||
, '--command=AnalyseComplex'
|
||||
, '--pdb-dir=' + process_dir
|
||||
, '--pdb=' + pdbname + '_Repair.pdb'
|
||||
, '--analyseComplexChains=' + chain1 + ',' + chain2
|
||||
, '--water=PREDICT'
|
||||
, '--vdwDesign=1'
|
||||
, '--output-dir=' + process_dir], cwd=process_dir)
|
||||
|
||||
interactions = ['Distances','Electro_RR','Electro_MM','Electro_SM','Electro_SS','Disulfide_RR','Disulfide_MM','Disulfide_SM','Disulfide_SS',
|
||||
'Hbonds_RR','Hbonds_MM','Hbonds_SM','Hbonds_SS','Partcov_RR','Partcov_MM','Partcov_SM','Partcov_SS','VdWClashes_RR','VdWClashes_MM',
|
||||
'VdWClashes_SM','VdWClashes_SS','Volumetric_RR','Volumetric_MM','Volumetric_SM','Volumetric_SS']
|
||||
# FIXME why would we ever need to do this?!? Cargo-culted from runcomplex.sh
|
||||
ac_source = process_dir + '/Summary_' + pdbname + '_Repair_AC.fxout'
|
||||
ac_dest = process_dir + '/Summary_' + pdbname + '_Repair_AC.txt'
|
||||
shutil.copyfile(ac_source, ac_dest)
|
||||
print('\033[95mCOMPLETED STAGE: foldx AnalyseComplex (subprocess) for WT:\033[0m', n)
|
||||
|
||||
for n in range(1,nmuts+1):
|
||||
print('\033[95mSTAGE: Running foldx AnalyseComplex (foldx subprocess) for mutation:\033[0m', n)
|
||||
subprocess.call(['foldx'
|
||||
, '--command=AnalyseComplex'
|
||||
, '--pdb-dir=' + process_dir
|
||||
, '--pdb=' + pdbname + '_Repair_' + str(n) + '.pdb'
|
||||
, '--analyseComplexChains=' + chain1 + ',' + chain2
|
||||
, '--water=PREDICT'
|
||||
, '--vdwDesign=1'
|
||||
, '--output-dir=' + process_dir], cwd=process_dir)
|
||||
|
||||
# FIXME why would we ever need to do this?!? Cargo-culted from runcomplex.sh
|
||||
ac_mut_source = process_dir + '/Summary_' + pdbname + '_Repair_' + str(n) +'_AC.fxout'
|
||||
ac_mut_dest = process_dir + '/Summary_' + pdbname + '_Repair)' + str(n) +'_AC.txt'
|
||||
shutil.copyfile(ac_mut_source, ac_mut_dest)
|
||||
print('\033[95mCOMPLETED STAGE: foldx AnalyseComplex (subprocess) for mutation:\033[0m', n)
|
||||
print('\n==========================================================')
|
||||
|
||||
interactions = ['Distances','Electro_RR','Electro_MM','Electro_SM','Electro_SS','Disulfide_RR','Disulfide_MM','Disulfide_SM','Disulfide_SS', 'Hbonds_RR','Hbonds_MM','Hbonds_SM','Hbonds_SS','Partcov_RR','Partcov_MM','Partcov_SM','Partcov_SS','VdWClashes_RR','VdWClashes_MM','VdWClashes_SM','VdWClashes_SS','Volumetric_RR','Volumetric_MM','Volumetric_SM','Volumetric_SS']
|
||||
|
||||
dGdatafile = process_dir + '/Dif_' + pdbname + '_Repair.txt'
|
||||
dGdata = pd.read_csv(dGdatafile, sep = '\t')
|
||||
|
|
@ -267,8 +389,7 @@ def main():
|
|||
print(d)
|
||||
data[i+1] = d
|
||||
|
||||
interactions = ['ddG', 'Distances','Electro_RR','Electro_MM','Electro_SM','Electro_SS','Disulfide_RR','Disulfide_MM','Disulfide_SM','Disulfide_SS', 'Hbonds_RR','Hbonds_MM','Hbonds_SM','Hbonds_SS','Partcov_RR','Partcov_MM','Partcov_SM','Partcov_SS','VdWClashes_RR','VdWClashes_MM',
|
||||
'VdWClashes_SM','VdWClashes_SS','Volumetric_RR','Volumetric_MM','Volumetric_SM','Volumetric_SS']
|
||||
interactions = ['ddG', 'Distances','Electro_RR','Electro_MM','Electro_SM','Electro_SS','Disulfide_RR','Disulfide_MM','Disulfide_SM','Disulfide_SS', 'Hbonds_RR','Hbonds_MM','Hbonds_SM','Hbonds_SS','Partcov_RR','Partcov_MM','Partcov_SM','Partcov_SS','VdWClashes_RR','VdWClashes_MM','VdWClashes_SM','VdWClashes_SS','Volumetric_RR','Volumetric_MM','Volumetric_SM','Volumetric_SS']
|
||||
|
||||
print(interactions)
|
||||
|
||||
|
|
@ -292,8 +413,7 @@ def main():
|
|||
print(len(IE))
|
||||
data = np.append(data,[IE], axis = 0)
|
||||
print(data)
|
||||
interactions = ['ddG','Distances','Electro_RR','Electro_MM','Electro_SM','Electro_SS','Disulfide_RR','Disulfide_MM','Disulfide_SM','Disulfide_SS', 'Hbonds_RR','Hbonds_MM','Hbonds_SM','Hbonds_SS','Partcov_RR','Partcov_MM','Partcov_SM','Partcov_SS','VdWClashes_RR','VdWClashes_MM',
|
||||
'VdWClashes_SM','VdWClashes_SS','Volumetric_RR','Volumetric_MM','Volumetric_SM','Volumetric_SS','Interaction Energy']
|
||||
interactions = ['ddG','Distances','Electro_RR','Electro_MM','Electro_SM','Electro_SS','Disulfide_RR','Disulfide_MM','Disulfide_SM','Disulfide_SS','Hbonds_RR','Hbonds_MM','Hbonds_SM','Hbonds_SS','Partcov_RR','Partcov_MM','Partcov_SM','Partcov_SS','VdWClashes_RR','VdWClashes_MM','VdWClashes_SM','VdWClashes_SS','Volumetric_RR','Volumetric_MM','Volumetric_SM','Volumetric_SS','Interaction Energy']
|
||||
|
||||
mut_file = process_dir + '/individual_list_' + pdbname + '.txt'
|
||||
with open(mut_file) as csvfile:
|
||||
|
|
@ -331,6 +451,7 @@ def main():
|
|||
#outputfilename = 'foldx_results_' + pdbname + '.csv'
|
||||
#results.to_csv(outputfilename)
|
||||
results2.to_csv(outputfilename, index = False)
|
||||
print ('end')
|
||||
|
||||
if __name__ == '__main__':
|
||||
main()
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue