Entering edit mode
8.0 years ago
Burgenix
•
0
Hi,
Running MIRA for a contig threshold here. I am trying to plot the contig length and the number of iterations, but somehow, the plotter stops the program after a iteration. Any idea why?
import os
import sys
import shutil
import subprocess
import matplotlib.pyplot as plt
import numpy as np
name = sys.argv[1] + '-0'
ref_loc = sys.argv[2]
read_loc = sys.argv[3]
manifest = open('manifest-0.conf', 'w')
manifest.write('project = ' + name + '\njob = genome,mapping,accurate\nparameters = -GE:not=6 -NW:mrnl=0 -AS:nop=1 SOLEXA_SETTINGS -CO:msr=no COMMON_SETTINGS -SB:tor=no\nreadgroup\nis_reference\ndata = ' + ref_loc + '\nstrain = ' + name + '\nreadgroup = reads\ndata = ' + read_loc + '\ntechnology=solexa')
manifest.close()
avg_cov = 0
count = 0
new_avg_cov = 0
contig_length = 0
while new_avg_cov <= avg_cov * 2 or avg_cov == 0: #set threshold for average coverage points
command1 = 'mira manifest-' + str(count) + '.conf'
print(command1)
avg_cov = new_avg_cov
subprocess.call(command1, shell=True)
print('MIRA is really done!')
with open(name + '_assembly/' + name +'_d_info/'+ name +'_info_contigstats.txt') as f:
print('Fishing for coverage:')
for line in f:
if not line.startswith('#'):
new_avg_cov = float(line.split('\t')[5])
with open(name + '_assembly/' + name +'_d_info/'+ name +'_info_contigstats.txt') as f:
if not line.startswith('#'):
contig_length = int(line.split('\t')[1])
plt.scatter(count, contig_length)
plt.xlabel("Number of iterations", fontsize=14)
plt.ylabel("Length of assembled contig", fontsize=14)
plt.axis([0,100, 0,200000])
plt.pause(0.001)
plt.show()
new_bait = name + '_assembly/' + name +'_d_results/'+ name +'_out_AllStrains.unpadded.fasta'
print('new bait is ' + new_bait)
count += 1
name = name.split('-')[0] + '-' + str(count)
manifest = open('manifest-' + str(count) +'.conf', 'w')
print('new manifest is ' + 'manifest-' + str(count) +'.conf')
manifest.write('project = ' + name +'\njob = genome,mapping,accurate\nparameters = -GE:not=6 -NW:mrnl=0 -AS:nop=1 SOLEXA_SETTINGS -CO:msr=no COMMON_SETTINGS -SB:tor=no\nreadgroup\nis_reference\ndata = /home/radu/Desktop/new_test' + new_bait + '\nstrain = ' + name + '\nreadgroup = reads\ndata = ' + read_loc + '\ntechnology=solexa')
manifest.close()
if os.path.exists('./manifest-' + str(count-2) + '.conf'):
os.remove('./manifest-' + str(count-2) + '.conf')
shutil.rmtree('./' + name.split('-')[0] + '-' + str(count-2) + '_assembly/' )
Tried to use a line instead of a scatterplot, same result. I think It might be a problem with the manifest file, but i cant figure it out.
Have you tried running the program without
plt.showorplt.pause? Both of these functions may cause the program to block until some GUI action is taken. Similarly, have you checked that the conditions in the while loop are not satisfied?