python - 使用 python 更快地读取大型 fastq 文件

Question

我有几个平均有 500.000.000 行(125.000.000 个序列)的 fastq 文件。有没有一种快速的方法可以更快地读取这些fastq文件。

我想做的是读取每个序列并使用前 16 个序列作为条形码。然后统计每个文件中的条形码数量。

这是我的脚本，需要几个小时:

import os, errno
from Bio import SeqIO
import gzip
files = os.listdir(".")
for file in files[:]:
    if not file.endswith(".fastq.gz"):
        files.remove(file)

maps = {}
for file in files:
    print "Now Parsing file %s"%file
    maps[file] = {}
    with gzip.open(file,"r") as handle:
        recs = SeqIO.parse(handle,"fastq")
        for rec in recs:
            tag = str(rec.seq)[0:16]
            if tag not in map[file]:
                maps[file][tag] = 1
            else:
                maps[file][tag] += 1

我有 250 GB RAM 和 20 个 CPU，可用于多线程...

谢谢。

Answer 1

未经测试，但您可以通过以下方式以“令人尴尬的并行”方式做到这一点:

import multiprocessing as mp
import os, errno
from Bio import SeqIO
import gzip

def ImportFile(file):

    maps = {}
    with gzip.open(file,"r") as handle:
        recs = SeqIO.parse(handle,"fastq")
        for rec in recs:
            tag = str(rec.seq)[0:16]
            if tag not in maps.keys():
                maps[tag] = 1
            else:
                maps[tag] += 1

    return {file:maps}


files = os.listdir(".")
for file in files[:]:
    if not file.endswith(".fastq.gz"):
        files.remove(file)

# I'd test this with smaller numbers before using up all 20 cores
pool = mp.Pool(processes=10)
output = pool.map(ImportFile,files)