TSS like old code

Use get_id() not name

Author: davmlaw

bin/pyreference_gtf_to_json.py

@@ -7,9 +7,10 @@
 from collections import defaultdict, Counter
 import gzip
 import json
+import os
 from pyreference.settings import CHROM, START, END, STRAND, IS_CODING, CODING_FEATURES, \
     PYREFERENCE_JSON_VERSION_KEY, PYREFERENCE_JSON_VERSION
-from pyreference.utils.file_utils import name_from_file_name
+from pyreference.utils.file_utils import name_from_file_name, file_md5sum
 
 
 class SetEncoder(json.JSONEncoder):
@@ -158,6 +159,8 @@ def main():
     add_UTR_features(transcripts_by_id, transcript_cds_by_id)
 
     data = {PYREFERENCE_JSON_VERSION_KEY : PYREFERENCE_JSON_VERSION,
+            "reference_gtf" : {"path" : os.path.abspath(args.gtf),
+                               "md5sum" : file_md5sum(args.gtf)},
             "genes_by_id" : genes_by_id,
             "transcripts_by_id" : transcripts_by_id,
             "gene_id_by_name" : gene_id_by_name,

pyreference/gene.py

@@ -66,7 +66,7 @@ def get_longest_transcript(self, coding_only=False):
             # We also want NM_007041 (len 2209) over NM_001001976 (len 2209)
             # Which is annoyingly zero padded - so use smallest ID length, then only if equal do alpha sort 
             def min_transcript_key(t):
-                (sys.maxint - t.length, len(t.get_id()), t.get_id())
+                return (sys.maxint - t.length, len(t.get_id()), t.get_id())
 
             longest_transcript = min(transcripts, key=min_transcript_key)
         return longest_transcript

pyreference/genomic_region.py

@@ -30,7 +30,10 @@ def iv(self):
 
     @lazy    
     def tss(self):
-        return self.iv.start_d_as_pos
+        ''' (Representative) Transcript Start Site
+            This is NOT the most 5' position (use iv.start_d_as_pos for that) '''
+        transcript_iv = self.get_representative_transcript().iv
+        return transcript_iv.start_d_as_pos
 
     def get_promoter_iv(self, promoter_range=1000):
         return iv_from_pos_range(self.tss, promoter_range)

pyreference/pyreference_config.py

@@ -42,8 +42,8 @@ def load_params_from_config(build=None, config=None):
             raise_from(configparser.NoOptionError(msg), noe)
 
     if not cfg.has_section(build):
-        params = {"build" : build, "config" : config}
-        msg = "Build='%(build)s', no section [%(build)s] in config file '%(config)s" % params 
+        msg_params = {"build" : build, "config" : config}
+        msg = "Build='%(build)s', no section [%(build)s] in config file '%(config)s" % msg_params 
         raise ValueError(msg)
 
     for f in GLOBAL_FLAGS:
@@ -52,6 +52,9 @@ def load_params_from_config(build=None, config=None):
         except configparser.NoOptionError as noe:
             pass
 
+    params["build"] = build
+    params["config"] = config
+    
     for k in defaults.keys():
         params[k] = cfg.get(build, k)
 

pyreference/reference.py

@@ -104,6 +104,11 @@ def __init__(self, build=None, config=None, **kwargs):
             if kwargs:
                 six.raise_from(ValueError("No 'genes_json' in passed kwargs"), config_exception)
             raise config_exception
+        
+        # Store this so we can ask about config later
+        self.build = params["build"]
+        self._args = {"build" : build, "config" : config}
+        self._build_params = params
 
 
     @lazy
@@ -281,7 +286,7 @@ def get_transcripts_in_iv(self, iv):
         return list(transcripts)
 
     def get_transcript_ids(self, iv):
-        return [feature.name for feature in self.get_transcripts_in_iv(iv)]
+        return [feature.get_id() for feature in self.get_transcripts_in_iv(iv)]
 
     def get_gene_names_array(self, iv):
         return list(set([t.get_gene_id() for t in self.get_transcripts_in_iv(iv)]))
@@ -355,4 +360,15 @@ def has_chr(self):
         chrom = some_transcript["chr"]
         return chrom.startswith("chr")
 
+    def __repr__(self):
+        return "PyReference (%s)" % self.build
+
+    @lazy
+    def config(self):
+        params = {"build" : self.build,
+                  "args" : self._args,
+                  "reference_gtf" : self._genes_dict["reference_gtf"],
+                  "build_params" : self._build_params.copy()}
+        return params
+
 

pyreference/settings.py

@@ -6,7 +6,7 @@
 
 
 # Change this when you introduce breaking changes
-PYREFERENCE_JSON_VERSION = 4
+PYREFERENCE_JSON_VERSION = 5
 PYREFERENCE_JSON_VERSION_KEY = "pyreference_json_version"
 
 CODING_FEATURES = {"CDS", "start_codon", "stop_codon"}

pyreference/transcript.py

@@ -176,3 +176,6 @@ def get_genomic_position(self, pos_on_transcript):
 
         raise NotOnTranscriptException("%s didn't contain %s" % (self.get_id(), pos_on_transcript))
 
+    def __repr__(self):
+        coding_str = " (coding)" if self.is_coding else ""
+        return "Transcript %s%s: length %d" % (self.get_id(), coding_str, self.length)

pyreference/utils/file_utils.py

@@ -3,6 +3,7 @@
 
 @author: dlawrence
 '''
+from md5 import md5
 import os
 
 
@@ -34,3 +35,9 @@ def file_or_file_name(f, mode='r'):
         return f # Already a File object
     else:
         raise ValueError("'%s' (%s) not a file or string" % (f, type(f)))
+
+def file_md5sum(filename):
+    m = md5()
+    with open(filename, "rb") as f:
+        m.update(f.read())
+    return m.hexdigest()

setup.py

@@ -4,7 +4,7 @@
 
 setup(name = 'pyreference',
       packages = find_packages(),
-      version = '0.2',
+      version = '0.3',
       description = 'Library for working with reference genomes',
       author = 'David Lawrence',
       author_email = 'davmlaw@gmail.com',