Entering edit mode
7.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.show
orplt.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?