San*_*dra 7 python performance bioinformatics fastq
我是Python的新手,我编写了一个(可能非常难看)脚本,它应该从fastq文件中随机选择一个序列子集.fastq文件以每行四行的块存储信息.每个块中的第一行以字符"@"开头.我用作输入文件的fastq文件是36 GB,包含大约14,000,000行.
我试图重写一个使用过多内存的现有脚本,并设法减少了很多内存使用量.但脚本需要永远运行,我不明白为什么.
parser = argparse.ArgumentParser()
parser.add_argument("infile", type = str, help = "The name of the fastq input file.", default = sys.stdin)
parser.add_argument("outputfile", type = str, help = "Name of the output file.")
parser.add_argument("-n", help="Number of sequences to sample", default=1)
args = parser.parse_args()
def sample():
linesamples = []
infile = open(args.infile, 'r')
outputfile = open(args.outputfile, 'w')
# count the number of fastq "chunks" in the input file:
seqs = subprocess.check_output(["grep", "-c", "@", str(args.infile)])
# randomly select n fastq "chunks":
seqsamples = random.sample(xrange(0,int(seqs)), int(args.n))
# make a list of the lines that are to be fetched from the fastq file:
for i in seqsamples:
linesamples.append(int(4*i+0))
linesamples.append(int(4*i+1))
linesamples.append(int(4*i+2))
linesamples.append(int(4*i+3))
# fetch lines from input file and write them to output file.
for i, line in enumerate(infile):
if i in linesamples:
outputfile.write(line)
Run Code Online (Sandbox Code Playgroud)
grep-step几乎没有时间,但是超过500分钟后,脚本仍然没有开始写入输出文件.所以我想这是grep和最后一个for循环之间需要很长时间的步骤之一.但我完全不明白哪一步,以及我能做些什么来加快它.
根据 的大小linesamples,if i in linesamples由于您要在列表中搜索每次迭代,因此将需要很长时间infile。您可以将其转换为 aset以缩短查找时间。而且,enumerate效率不是很高 - 我已经用line_num我们在每次迭代中递增的构造替换了它。
def sample():
linesamples = set()
infile = open(args.infile, 'r')
outputfile = open(args.outputfile, 'w')
# count the number of fastq "chunks" in the input file:
seqs = subprocess.check_output(["grep", "-c", "@", str(args.infile)])
# randomly select n fastq "chunks":
seqsamples = random.sample(xrange(0,int(seqs)), int(args.n))
for i in seqsamples:
linesamples.add(int(4*i+0))
linesamples.add(int(4*i+1))
linesamples.add(int(4*i+2))
linesamples.add(int(4*i+3))
# make a list of the lines that are to be fetched from the fastq file:
# fetch lines from input file and write them to output file.
line_num = 0
for line in infile:
if line_num in linesamples:
outputfile.write(line)
line_num += 1
outputfile.close()
Run Code Online (Sandbox Code Playgroud)