treetime/clock_tree.py
def calc_rate_susceptibility(self, rate_std=None, params=None):
        """return the time tree estimation of evolutionary rates +/- one
        standard deviation form the ML estimate.

        Returns
        -------
        TreeTime.return_code : str
            success or failure
        """
        params = params or {}
        if rate_std is None:
            if not (self.clock_model['valid_confidence'] and 'cov' in self.clock_model):
                self.logger("ClockTree.calc_rate_susceptibility: need valid standard deviation of the clock rate to estimate dating error.", 1, warn=True)
                return ttconf.ERROR
            rate_std = np.sqrt(self.clock_model['cov'][0,0])

        current_rate = np.abs(self.clock_model['slope'])
        upper_rate = self.clock_model['slope'] + rate_std
        lower_rate = max(0.1*current_rate, self.clock_model['slope'] - rate_std)
        for n in self.tree.find_clades():
            if n.up:
                n._orig_gamma = n.branch_length_interpolator.gamma
                n.branch_length_interpolator.gamma = n._orig_gamma*upper_rate/current_rate

        self.logger("###ClockTree.calc_rate_susceptibility: run with upper bound of rate estimate", 1)
        self.make_time_tree(**params)
        self.logger("###ClockTree.calc_rate_susceptibility: rate: %f, LH:%f"%(upper_rate, self.tree.positional_joint_LH), 2)
        for n in self.tree.find_clades():
            n.numdate_rate_variation = [(upper_rate, n.numdate)]
            if n.up:
                n.branch_length_interpolator.gamma = n._orig_gamma*lower_rate/current_rate

        self.logger("###ClockTree.calc_rate_susceptibility: run with lower bound of rate estimate", 1)
        self.make_time_tree(**params)
        self.logger("###ClockTree.calc_rate_susceptibility: rate: %f, LH:%f"%(lower_rate, self.tree.positional_joint_LH), 2)
        for n in self.tree.find_clades():
            n.numdate_rate_variation.append((lower_rate, n.numdate))
            if n.up:
                n.branch_length_interpolator.gamma  = n._orig_gamma

        self.logger("###ClockTree.calc_rate_susceptibility: run with central rate estimate", 1)
        self.make_time_tree(**params)
        self.logger("###ClockTree.calc_rate_susceptibility: rate: %f, LH:%f"%(current_rate, self.tree.positional_joint_LH), 2)
        for n in self.tree.find_clades():
            n.numdate_rate_variation.append((current_rate, n.numdate))
            n.numdate_rate_variation.sort(key=lambda x:x[1]) # sort estimates for different rates by numdate

        return ttconf.SUCCESS
def calc_rate_susceptibility(self, rate_std=None, params=None):
        """return the time tree estimation of evolutionary rates +/- one
        standard deviation form the ML estimate.

        Returns
        -------
        TreeTime.return_code : str
            success or failure
        """
        params = params or {}
        if rate_std is None:
            if not (self.clock_model['valid_confidence'] and 'cov' in self.clock_model):
                self.logger("ClockTree.calc_rate_susceptibility: need valid standard deviation of the clock rate to estimate dating error.", 1, warn=True)
                return ttconf.ERROR
            rate_std = np.sqrt(self.clock_model['cov'][0,0])

        current_rate = np.abs(self.clock_model['slope'])
        upper_rate = self.clock_model['slope'] + rate_std
        lower_rate = max(0.1*current_rate, self.clock_model['slope'] - rate_std)
        for n in self.tree.find_clades():
            if n.up:
                n._orig_gamma = n.branch_length_interpolator.gamma
                n.branch_length_interpolator.gamma = n._orig_gamma*upper_rate/current_rate

        self.logger("###ClockTree.calc_rate_susceptibility: run with upper bound of rate estimate", 1)
        self.make_time_tree(**params)
        self.logger("###ClockTree.calc_rate_susceptibility: rate: %f, LH:%f"%(upper_rate, self.tree.positional_joint_LH), 2)
        for n in self.tree.find_clades():
            n.numdate_rate_variation = [(upper_rate, n.numdate)]
            if n.up:
                n.branch_length_interpolator.gamma = n._orig_gamma*lower_rate/current_rate

        self.logger("###ClockTree.calc_rate_susceptibility: run with lower bound of rate estimate", 1)
        self.make_time_tree(**params)
        self.logger("###ClockTree.calc_rate_susceptibility: rate: %f, LH:%f"%(lower_rate, self.tree.positional_joint_LH), 2)
        for n in self.tree.find_clades():
            n.numdate_rate_variation.append((lower_rate, n.numdate))
            if n.up:
                n.branch_length_interpolator.gamma  = n._orig_gamma

        self.logger("###ClockTree.calc_rate_susceptibility: run with central rate estimate", 1)
        self.make_time_tree(**params)
        self.logger("###ClockTree.calc_rate_susceptibility: rate: %f, LH:%f"%(current_rate, self.tree.positional_joint_LH), 2)
        for n in self.tree.find_clades():
            n.numdate_rate_variation.append((current_rate, n.numdate))
            n.numdate_rate_variation.sort(key=lambda x:x[1]) # sort estimates for different rates by numdate

        return ttconf.SUCCESS
def calc_rate_susceptibility ( self , rate_std = None , params = None ) : params = params or { } if rate_std is None : if not ( self . clock_model [ 'valid_confidence' ] and 'cov' in self . clock_model ) : self . logger ( "ClockTree.calc_rate_susceptibility: need valid standard deviation of the clock rate to estimate dating error." , 1 , warn = True ) return ttconf . ERROR rate_std = np . sqrt ( self . clock_model [ 'cov' ] [ 0 , 0 ] ) current_rate = np . abs ( self . clock_model [ 'slope' ] ) upper_rate = self . clock_model [ 'slope' ] + rate_std lower_rate = max ( 0.1 * current_rate , self . clock_model [ 'slope' ] - rate_std ) for n in self . tree . find_clades ( ) : if n . up : n . _orig_gamma = n . branch_length_interpolator . gamma n . branch_length_interpolator . gamma = n . _orig_gamma * upper_rate / current_rate self . logger ( "###ClockTree.calc_rate_susceptibility: run with upper bound of rate estimate" , 1 ) self . make_time_tree ( * * params ) self . logger ( "###ClockTree.calc_rate_susceptibility: rate: %f, LH:%f" % ( upper_rate , self . tree . positional_joint_LH ) , 2 ) for n in self . tree . find_clades ( ) : n . numdate_rate_variation = [ ( upper_rate , n . numdate ) ] if n . up : n . branch_length_interpolator . gamma = n . _orig_gamma * lower_rate / current_rate self . logger ( "###ClockTree.calc_rate_susceptibility: run with lower bound of rate estimate" , 1 ) self . make_time_tree ( * * params ) self . logger ( "###ClockTree.calc_rate_susceptibility: rate: %f, LH:%f" % ( lower_rate , self . tree . positional_joint_LH ) , 2 ) for n in self . tree . find_clades ( ) : n . numdate_rate_variation . append ( ( lower_rate , n . numdate ) ) if n . up : n . branch_length_interpolator . gamma = n . _orig_gamma self . logger ( "###ClockTree.calc_rate_susceptibility: run with central rate estimate" , 1 ) self . make_time_tree ( * * params ) self . logger ( "###ClockTree.calc_rate_susceptibility: rate: %f, LH:%f" % ( current_rate , self . tree . positional_joint_LH ) , 2 ) for n in self . tree . find_clades ( ) : n . numdate_rate_variation . append ( ( current_rate , n . numdate ) ) n . numdate_rate_variation . sort ( key = lambda x : x [ 1 ] ) # sort estimates for different rates by numdate return ttconf . SUCCESS
return the time tree estimation of evolutionary rates +/- one
        standard deviation form the ML estimate.

        Returns
        -------
        TreeTime.return_code : str
            success or failure
return the time tree estimation of evolutionary rates + / - one standard deviation form the ML estimate .