与第二个文件中的数组相比，解析文件到数组的方式更快
我目前有一个包含MS2光谱数据(QE_2706_229_sequest_high_conf.mgf)的MGF文件。文件模板在下面的链接中，以及示例的片段：
http://www.matrixscience.com/help/data_file_help.html
BEGIN IONS
TITLE=File3249 Spectrum10594 scans: 11084
PEPMASS=499.59366 927079.3
CHARGE=3+
RTINSECONDS=1710
SCANS=11084
104.053180 3866.360000
110.071530 178805.000000
111.068610 1869.210000
111.074780 10738.600000
112.087240 13117.900000
113.071150 7148.790000
114.102690 4146.490000
115.086840 11835.600000
116.070850 6230.980000
... ...
END IONS这个未加注释的光谱文件包含数千个这样的条目，总文件大小约为150MB。然后我就有了一系列需要解析的文本文件。每个文件都类似于上面的格式，第一列被读入一个numpy数组。然后为每个条目解析未注释的光谱文件，直到从注释的文本文件输入中找到匹配的数组。
(文件名GRPGPVAGHHQMPR)
m/z i matches
104.05318 3866.4
110.07153 178805.4
111.06861 1869.2
111.07478 10738.6
112.08724 13117.9
113.07115 7148.8
114.10269 4146.5
115.08684 11835.6
116.07085 6231.0一旦找到匹配项，就会写入一个MGF注释文件，该文件包含未注释文件中的完整条目信息，但其中一行指定了与该特定条目匹配的注释文本文件的文件名。输出如下：
BEGIN IONS
SEQ=GRPGPVAGHHQMPR
TITLE=File3249 Spectrum10594 scans: 11084
PEPMASS=499.59366 927079.3
... ...
END IONS可能有一种计算成本更低的解析方法。给定要搜索的2,000个带注释的文件，使用上面的大型未注释文件，在2.6 GHz四核英特尔Haswell cpu上解析当前需要大约12个小时。
下面是工作代码：
import numpy as np
import sys
from pyteomics import mgf
from glob import glob
def main():
"""
Usage: python mgf_parser
"""
pep_files = glob('*.txt')
mgf_file = 'QE_2706_229_sequest_high_conf.mgf'
process(mgf_file, pep_files)
def process(mgf_file, pep_files):
"""Parses spectra from annotated text file. Converts m/z values to numpy array.
If spectra array matches entry in MGF file, writes annotated MGF file.
"""
ann_arrays = {}
for ann_spectra in pep_files:
a = np.genfromtxt(ann_spectra, dtype=float, invalid_raise=False,
usemask=False, filling_values=0.0, usecols=(0))
b = np.delete(a, 0)
ann_arrays[ann_spectra] = b
with mgf.read(mgf_file) as reader:
for spectrum in reader:
for ann_spectra, array in ann_arrays.iteritems():
if np.array_equal(array, spectrum['m/z array']):
print '> Spectral match found for file {}.txt'.format(ann_spectra[:-4])
file_name = '{}.mgf'.format(ann_spectra[:-4])
spectrum['params']['seq'] = file_name[52:file_name.find('-') - 1]
mgf.write((spectrum,), file_name)
if __name__ == '__main__':
main()过去，它一次只能解析给定数量的文件。有没有更有效的解析方法的建议？
回答 1
Stack Overflow用户
发布于 2014-07-16 02:23:51
我看到了改进的空间，因为您正在为每个小文件重复解析整个MGF文件。如果您重构代码，使其只被解析一次，您可能会获得相当不错的加速。
下面是我如何调整代码，同时去掉bash循环，并使用mgf.write函数，该函数可能比np.savetxt慢一点，但更易于使用：
from pyteomics import mgf
import sys
import numpy as np
def process(mgf_file, pep_files):
ann_arrays = {}
for ann_spectra in pep_files:
a = np.genfromtxt(ann_spectra, invalid_raise=False,
filling_values=0.0, usecols=(0,))
b = np.delete(a, 0)
ann_arrays[ann_spectra] = b
with mgf.read(mgf_file) as reader:
for spectrum in reader:
for ann_spectra, array in ann_arrays.iteritems():
if np.allclose(array, spectrum['m/z array']):
# allclose may be better for floats than array_equal
file_name = 'DeNovo/good_training_seq/{}.mgf'.format(
ann_spectra[:-4])
spectrum['params']['seq'] = ann_spectra[
:ann_spectra.find('-') - 1]
mgf.write((spectrum,), file_name)
if __name__ == '__main__':
pep_files = sys.argv[1:]
mgf_file = '/DeNovo/QE_2706_229_sequest_high_conf.mgf'
process(mgf_file, pep_files)然后，为了实现与bash循环相同的效果，您可以将其调用为
python2.7 mgf_parser.py *.txt如果展开的参数列表太长，可以使用glob，而不是依靠bash来展开它：
from glob import iglob
pep_files = iglob(sys.argv[1])并这样命名，以防止bash的扩展：
python2.7 mgf_parser.py '*.txt'https://stackoverflow.com/questions/24751797
复制相似问题