X-Git-Url: http://woldlab.caltech.edu/gitweb/?p=erange.git;a=blobdiff_plain;f=cistematic%2Fgenomes%2Fmmusculus.py;h=72f1350070ff7027118e5a9c569e1da23d5e8e41;hp=4e90cfe7e271f5b704d6cb152bb83da2b0264ad9;hb=HEAD;hpb=4ad5495359e4322da39868020a7398676261679e diff --git a/cistematic/genomes/mmusculus.py b/cistematic/genomes/mmusculus.py index 4e90cfe..72f1350 100644 --- a/cistematic/genomes/mmusculus.py +++ b/cistematic/genomes/mmusculus.py @@ -1,7 +1,7 @@ ########################################################################### # # # C O P Y R I G H T N O T I C E # -# Copyright (c) 2003-10 by: # +# Copyright (c) 2003-13 by: # # * California Institute of Technology # # # # All Rights Reserved. # @@ -34,36 +34,58 @@ from cistematic.core.geneinfo import geneinfoDB from os import environ if environ.get("CISTEMATIC_ROOT"): - cisRoot = environ.get("CISTEMATIC_ROOT") + cisRoot = environ.get("CISTEMATIC_ROOT") else: cisRoot = "/proj/genome" geneDB = "%s/M_musculus/mmusculus.genedb" % cisRoot -def loadChromosome(db, chromID, chromPath, chromOut): - seqArray = [] - mmGenome = Genome("mmusculus", dbFile=db) - inFile = open(chromPath, "r") - line = inFile.readline() - for line in inFile: - seqArray.append(line.strip()) +def buildMmusculusDB(db=geneDB, downloadDir="%s/download" % cisRoot): + genePath = "%s/seq_gene.md" % downloadDir # ftp://ftp.ncbi.nih.gov/genomes/M_musculus/mapview/seq_gene.md + goDefPath = "%s/GO.terms_and_ids" % downloadDir # ftp://ftp.geneontology.org/pub/go/doc/GO.terms_and_ids + goPath = "%s/gene2go" % downloadDir # ftp://ftp.ncbi.nih.gov/gene/DATA/gene2go.gz + # chromosomes are from ftp://hgdownload.cse.ucsc.edu/goldenPath/mm9/chromosomes + # but ignoring all random chromosomes + chromList = ["1", "2", "3", "4", "5", "6", "7", "8", "9", "10", + "11", "12", "13", "14", "15", "16", "17", "18", "19", + "X", "Y", "M" + ] - seq = string.join(seqArray, "") - seqLen = len(seq) - if seqLen < 1: - print "Problems reading sequence from file" + print "Creating database %s" % db + createDBFile(db) - print "writing to file %s" % chromOut - outFile = open("%s%s" % (cisRoot, chromOut), "w") - outFile.write(seq) - outFile.close() - mmGenome.addChromosomeEntry(chromID, chromOut, "file") + print "Adding gene entries" + loadGeneEntries(db, genePath, chromList) + + print "Adding gene features" + loadGeneFeatures(db, genePath, chromList) + + print "Adding gene annotations" + loadGeneAnnotations(db) + + print "Adding gene ontology" + loadGeneOntology(db, goPath, goDefPath) + + for chromID in chromList: + print "Loading chromosome %s" % chromID + chromPath = "%s/chr%s.fa" % (downloadDir, chromID) + loadChromosome(db, chromID, chromPath, "/M_musculus/chromo%s.bin" % chromID) + + print "Creating Indices" + createDBindices(db) + + print "Finished creating database %s" % db + + +def createDBFile(db): + mmGenome = Genome("mmusculus", dbFile=db) + mmGenome.createGeneDB(db) def loadGeneEntries(db, gFile, cDict): - """ FIXME - NEED TO DEAL WITH ALTERNATIVE SPLICING ENTRIES - """ + #TODO: - NEED TO DEAL WITH ALTERNATIVE SPLICING ENTRIES + geneEntries = [] alreadySeen = [] mmGenome = Genome("mmusculus", dbFile=db) @@ -184,7 +206,7 @@ def loadGeneOntology(db, goPath, goDefPath): synonyms = idb.geneIDSynonyms(gID) if len(synonyms) >0: for entry in synonyms: - gene_name += "," + gene_name += "," gene_name += entry else: gene_name = " " @@ -198,66 +220,26 @@ def loadGeneOntology(db, goPath, goDefPath): mmGenome.addGoInfoBatch(goArray) -def createDBFile(db): - mmGenome = Genome("mmusculus", dbFile=db) - mmGenome.createGeneDB(db) - - -def createDBindices(db): +def loadChromosome(db, chromID, chromPath, chromOut): + seqArray = [] mmGenome = Genome("mmusculus", dbFile=db) - mmGenome.createIndices() - - -def buildMmusculusDB(db=geneDB, downloadDir="%s/download" % cisRoot): - genePath = "%s/seq_gene.md" % downloadDir # ftp://ftp.ncbi.nih.gov/genomes/M_musculus/mapview/seq_gene.md - goDefPath = "%s/GO.terms_and_ids" % downloadDir # ftp://ftp.geneontology.org/pub/go/doc/GO.terms_and_ids - goPath = "%s/gene2go" % downloadDir # ftp://ftp.ncbi.nih.gov/gene/DATA/gene2go.gz - # chromosomes are from ftp://hgdownload.cse.ucsc.edu/goldenPath/mm9/chromosomes - # but ignoring all random chromosomes - chromDict = {"1": "%s/chr1.fa" % downloadDir, - "2": "%s/chr2.fa" % downloadDir, - "3": "%s/chr3.fa" % downloadDir, - "4": "%s/chr4.fa" % downloadDir, - "5": "%s/chr5.fa" % downloadDir, - "6": "%s/chr6.fa" % downloadDir, - "7": "%s/chr7.fa" % downloadDir, - "8": "%s/chr8.fa" % downloadDir, - "9": "%s/chr9.fa" % downloadDir, - "10": "%s/chr10.fa" % downloadDir, - "11": "%s/chr11.fa" % downloadDir, - "12": "%s/chr12.fa" % downloadDir, - "13": "%s/chr13.fa" % downloadDir, - "14": "%s/chr14.fa" % downloadDir, - "15": "%s/chr15.fa" % downloadDir, - "16": "%s/chr16.fa" % downloadDir, - "17": "%s/chr17.fa" % downloadDir, - "18": "%s/chr18.fa" % downloadDir, - "19": "%s/chr19.fa" % downloadDir, - "X": "%s/chrX.fa" % downloadDir, - "Y": "%s/chrY.fa" % downloadDir, - "M": "%s/chrM.fa" % downloadDir - } - - print "Creating database %s" % db - createDBFile(db) - - print "Adding gene entries" - loadGeneEntries(db, genePath, chromDict) - - print "Adding gene features" - loadGeneFeatures(db, genePath, chromDict) - - print "Adding gene annotations" - loadGeneAnnotations(db) + inFile = open(chromPath, "r") + line = inFile.readline() + for line in inFile: + seqArray.append(line.strip()) - print "Adding gene ontology" - loadGeneOntology(db, goPath, goDefPath) + seq = string.join(seqArray, "") + seqLen = len(seq) + if seqLen < 1: + print "Problems reading sequence from file" - for chromID in chromDict.keys(): - print "Loading chromosome %s" % chromID - loadChromosome(db, chromID, chromDict[chromID], "/M_musculus/chromo%s.bin" % chromID) + print "writing to file %s" % chromOut + outFile = open("%s%s" % (cisRoot, chromOut), "w") + outFile.write(seq) + outFile.close() + mmGenome.addChromosomeEntry(chromID, chromOut, "file") - print "Creating Indices" - createDBindices(db) - print "Finished creating database %s" % db \ No newline at end of file +def createDBindices(db): + mmGenome = Genome("mmusculus", dbFile=db) + mmGenome.createIndices()