TechOverflow Logo
TechOverflow

NCBI GeneInfo in Python parsen

English Deutsch

Problem:

Du musst Dateien im NCBI GeneInfo-Format parsen, wie sie beispielsweise aus dem NCBI FTP GENE_INFO-Verzeichnis heruntergeladen werden können, in Python. Du möchtest alle Abhängigkeiten vermeiden.

Lösung

Hier ist ein Parser, der nur die Python-Standardbibliothek verwendet. Er wurde mit der All_Data.gene_info.gz vom 30.11.2013 getestet und parst jede Zeile in unveränderliche GeneInfo-Objekte.

Felder wie die Synonymliste oder die Datenbank-XRefs werden in Listen bzw. Dictionaries umgewandelt. Fehlende Werte, die in der NCBI-Darstellung durch "-" gekennzeichnet sind, werden durch None ersetzt, um das Filtern zu erleichtern.

Der Code enthält ein Beispiel und hängt von keiner Software oder Bibliothek über Python 2.6 hinaus ab. Es wird empfohlen, den Parser an unterschiedliche aber ähnliche Formate anzupassen oder an die eigenen Bedürfnisse anzupassen. Der Parser ist generatorbasiert und benötigt bei korrekter Verwendung konstanten Speicher.

Beachte, dass PyPy in einem kurzen Test mit All_Data.gene_info.gz etwa 60 % schneller war als CPython 2.7.6.

parse_ncbi_geneinfo.py
#!/usr/bin/env python3
"""
Ein einfacher Parser für das NCBI Gene Info-Format.

Version 1.1: Python3-ready
"""
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"

#Initialisiert GeneInfo Named Tuple. Hinweis: namedtuple ist unveränderlich
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):
    """Parst einen DB-Xref-String wie HGNC:5|MIM:138670 in ein Dictionary"""
    #Nach | aufteilen, Ergebnisse nach : aufteilen. Ein Dict erstellen (Python 2.6-kompatibler Weg).
    if xrefs == "-": return {}
    return dict([(xrefParts[0], xrefParts[2])
                for xrefParts in (xref.partition(":")
                  for xref in xrefs.split("|"))])

def parseNCBIGeneInfo(filename):
    """
    Ein Parser für das NCBI Gene Info-Format.
    Yielded Objekte, die Infos über ein einzelnes Gen enthalten.

    Unterstützt transparente Gzip-Dekomprimierung
    """
    #Mit transparenter Dekomprimierung parsen
    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")
            #Wenn dies fehlschlägt, ist das Format nicht standardkompatibel
            assert len(parts) == len(geneInfoFields)
            #Daten normalisieren
            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")
            }
            #Alternativ kann das Dictionary hier ausgegeben werden, wenn Mutabilität benötigt wird:
            #    yield normalizedInfo
            yield GeneInfo(**normalizedInfo)

if __name__ == "__main__":
    import argparse
    parser = argparse.ArgumentParser()
    parser.add_argument("geneinfo_file", help="Die NCBI GeneInfo-Eingabedatei (.gz erlaubt)")
    parser.add_argument("--print-records", action="store_true", help="Alle GeneInfo-Objekte ausgeben, nicht nur")
    args = parser.parse_args()
    #Parser ausführen
    recordCount = 0
    for geneInfo in parseNCBIGeneInfo(args.geneinfo_file):
        if args.print_records:
            print (geneInfo)
        #Auf Datensätze wie folgt zugreifen: my_gene_id = geneInfo.gene_id
        recordCount += 1
    print ("Total records: %d" % recordCount)
Check out similar posts by category: Bioinformatics, Python