Merge pull request #41 from jcha40/master

MS Timing normalization & remove numba

Author: iglc

BuildTree/CellPopulationEngine.py

@@ -4,7 +4,6 @@
 from collections import defaultdict
 import collections
 import numpy as np
-import numba as nb
 import itertools
 import operator
 import logging
@@ -46,7 +45,6 @@ def sample_ccf(xk, pk):
             return custm.rvs(size=1)[0]
 
     @staticmethod
-    @nb.njit()
     def logSumExp(ns):
         max_ = np.max(ns)
         ds = ns - max_

BuildTree/ClusterObject.py

@@ -1,6 +1,5 @@
 import logging
 import numpy as np
-import numba as nb
 
 
 class Cluster:
@@ -112,7 +111,6 @@ def remove_mutation(self, mutation, update_cluster_hist=True):
                                                                                                self._identifier))
 
     @staticmethod
-    @nb.njit()
     def logSumExp(ns):
         max_ = np.max(ns)
         ds = ns - max_

BuildTree/Tree.py

@@ -1,6 +1,5 @@
 import logging
 import numpy as np
-import numba as nb
 import itertools
 import functools
 
@@ -265,7 +264,6 @@ def get_possible_configurations(potential_children):
         return new_children
 
     @staticmethod
-    @nb.njit()
     def logSumExp(ns):
         max_ = np.max(ns)
         ds = ns - max_

Dockerfile

@@ -1,7 +1,7 @@
 from bitnami/minideb
 RUN install_packages python-pip build-essential python-dev r-base r-base-dev git graphviz python-tk
 RUN pip install setuptools wheel
-RUN pip install numpy scipy matplotlib pandas numba
+RUN pip install numpy scipy matplotlib pandas
 COPY req /tmp/req
 RUN apt-get -y upgrade
 RUN apt-get -y update

SinglePatientTiming/TimingEngine.py

@@ -53,6 +53,7 @@ def get_concordant_cn_states(self):
             states_across_samples = {}  # hold union of all copy number states in all samples
             for sample in self.sample_list:
                 if region in sample.cn_states:
+                    # Coerce copy number states to integer states
                     if 0. <= sample.cn_states[region].cn_a1 < 1.:
                         cn_a1 = 0.
                     elif 1. < sample.cn_states[region].cn_a1 <= 2.:
@@ -182,7 +183,7 @@ def _get_cluster_ccfs(self):
                 cluster_ccfs.setdefault(c, np.zeros((n_samples, 101)))
                 cluster_ccfs[c] += np.log(mut.ccf_dist + 1e-10)
         for c in cluster_ccfs:
-            cluster_ccfs[c] = np.exp(cluster_ccfs[c] - logsumexp(cluster_ccfs[c]))
+            cluster_ccfs[c] = np.exp(cluster_ccfs[c] - logsumexp(cluster_ccfs[c], axis=1, keepdims=True))
         return cluster_ccfs
 
     def time_events(self):
@@ -375,7 +376,7 @@ def call_wgd(self, use_concordant_states=False):
                     regions_supporting_WGD.append(cn_state)
                 if cn_state.cn_a1 >= 2 and cn_state.cn_a2 >= 2:
                     regions_both_arms_gained.append(cn_state)
-            if len(regions_both_arms_gained) >= 5 or len(regions_supporting_WGD) * 2 >= \
+            if len(regions_both_arms_gained) >= 5 and len(regions_supporting_WGD) * 2 >= \
                     len(self.arm_regions) - len(self.missing_arms):
                 supporting_arm_states = [TimingCNState([self], s.chrN, s.arm, (s.cn_a1, s.cn_a2), s.purity, supporting_muts=s.supporting_muts) for
                                          s in supporting_arm_states]

req

@@ -1,6 +1,5 @@
 intervaltree==2.1.0
 scikit-learn==0.18.1
 networkx==1.11
-numba==0.45.1
 seaborn