TechOverflow Logo
TechOverflow

在 Python 中解析 NCBI GeneInfo

English Deutsch 简体中文

问题：

你需要用 Python 解析 NCBI GeneInfo 格式的文件，如可从 NCBI FTP GENE_INFO 目录下载的文件。你想避免任何依赖。

解决方案

这是一个仅使用 Python 标准库的解析器。它已使用 2013/11/30 的 All_Data.gene_info.gz 测试，并将每行解析为不可变的 GeneInfo 对象。

同义词列表或数据库 XRefs 等字段分别转换为列表或字典。NCBI 表示中由 “-“ 表示的缺失值替换为 None 以便更容易过滤。

代码包含一个示例，不依赖于 Python 2.6 之外的任何软件或库。建议你修改解析器以适应不同但相似的格式或满足你的需求。解析器基于生成器，如果正确使用需要恒定内存。

注意在 All_Data.gene_info.gz 上的快速测试中，PyPy 比 CPython 2.7.6 快约 60%

parse_ncbi_geneinfo.py
#!/usr/bin/env python3
"""
NCBI gene info 格式的简单解析器。

版本 1.1：Python3 就绪
"""
from __future__ import with_statement
from collections import namedtuple
import gzip
from datetime import datetime

__author__  = "Uli Köhler"
__license__ = "Apache License v2.0"
__version__ = "1.1"

#初始化 GeneInfo 命名元组。注意：namedtuple 是不可变的
geneInfoFields = ["tax_id", "gene_id", "symbol", "locus_tag", "synonyms", "db_xrefs", "chromosome", "map_location", "description", "type_of_gene", "symbol_from_nomenclature_authority", "full_name_from_nomenclature_authority", "nomenclature_status", "other_designations", "modification_date"]
GeneInfo = namedtuple("GeneInfo", geneInfoFields)

def parseDBXrefs(xrefs):
    """将 HGNC:5|MIM:138670 这样的 DB xref 字符串解析为字典"""
    #按 | 分割，按 : 分割结果。创建字典（python 2.6 兼容方式）。
    if xrefs == "-": return {}
    return dict([(xrefParts[0], xrefParts[2])
                for xrefParts in (xref.partition(":")
                  for xref in xrefs.split("|"))])

def parseNCBIGeneInfo(filename):
    """
    NCBI gene info 格式解析器。
    生成包含单个基因信息的对象。

    支持透明 gzip 解压缩
    """
    #以透明解压缩解析
    openFunc = gzip.open if filename.endswith(".gz") else open
    with openFunc(filename) as infile:
        for line in infile:
            if line.startswith("#"): continue
            parts = line.strip().split("\t")
            #如果失败，格式不兼容标准
            assert len(parts) == len(geneInfoFields)
            #规范化数据
            normalizedInfo = {
                "tax_id": int(parts[0]),
                "gene_id": int(parts[1]),
                "symbol": parts[2],
                "locus_tag": None if parts[3] == "-" else parts[3],
                "synonyms": [] if parts[4] == "-" else parts[4].split("|"),
                "db_xrefs": parseDBXrefs(parts[5]),
                "chromosome": parts[6],
                "map_location": parts[7],
                "description": parts[8],
                "type_of_gene": parts[9],
                "symbol_from_nomenclature_authority": None if parts[10] == "-" else parts[10],
                "full_name_from_nomenclature_authority": None if parts[11] == "-" else parts[11],
                "nomenclature_status": None if parts[12] == "-" else parts[12],
                "other_designations": None if parts[13] == "-" else parts[13],
                "modification_date": datetime.strptime(parts[14], "%Y%m%d")
            }
            #或者，你可以在此处发出字典，如果你需要可变性：
            #    yield normalizedInfo
            yield GeneInfo(**normalizedInfo)

if __name__ == "__main__":
    import argparse
    parser = argparse.ArgumentParser()
    parser.add_argument("geneinfo_file", help="NCBI GeneInfo 输入文件（允许 .gz）")
    parser.add_argument("--print-records", action="store_true", help="打印所有 GeneInfo 对象，不仅仅是")
    args = parser.parse_args()
    #执行解析器
    recordCount = 0
    for geneInfo in parseNCBIGeneInfo(args.geneinfo_file):
        if args.print_records:
            print (geneInfo)
        #像这样访问记录：my_gene_id = geneInfo.gene_id
        recordCount += 1
    print ("Total records: %d" % recordCount)
Check out similar posts by category: Bioinformatics, Python