treetime/treeanc.py
def optimize_branch_length(self, mode='joint', **kwargs):
        """
        Perform optimization for the branch lengths of the entire tree.
        This method only does a single path and needs to be iterated.

        **Note** this method assumes that each node stores information
        about its sequence as numpy.array object (node.sequence attribute).
        Therefore, before calling this method, sequence reconstruction with
        either of the available models must be performed.

        Parameters
        ----------

         mode : str
            Optimize branch length assuming the joint ML sequence assignment
            of both ends of the branch (:code:`joint`), or trace over all possible sequence
            assignments on both ends of the branch (:code:`marginal`) (slower, experimental).

         **kwargs :
            Keyword arguments

        Keyword Args
        ------------

         verbose : int
            Output level

         store_old : bool
            If True, the old lengths will be saved in :code:`node._old_dist` attribute.
            Useful for testing, and special post-processing.


        """

        self.logger("TreeAnc.optimize_branch_length: running branch length optimization in mode %s..."%mode,1)
        if (self.tree is None) or (self.aln is None):
            self.logger("TreeAnc.optimize_branch_length: ERROR, alignment or tree are missing", 0)
            return ttconf.ERROR

        store_old_dist = False

        if 'store_old' in kwargs:
            store_old_dist = kwargs['store_old']

        if mode=='marginal':
            # a marginal ancestral reconstruction is required for
            # marginal branch length inference
            if not hasattr(self.tree.root, "marginal_profile"):
                self.infer_ancestral_sequences(marginal=True)

        max_bl = 0
        for node in self.tree.find_clades(order='postorder'):
            if node.up is None: continue # this is the root
            if store_old_dist:
                node._old_length = node.branch_length

            if mode=='marginal':
                new_len = self.optimal_marginal_branch_length(node)
            elif mode=='joint':
                new_len = self.optimal_branch_length(node)
            else:
                self.logger("treeanc.optimize_branch_length: unsupported optimization mode",4, warn=True)
                new_len = node.branch_length

            if new_len < 0:
                continue

            self.logger("Optimization results: old_len=%.4e, new_len=%.4e, naive=%.4e"
                   " Updating branch length..."%(node.branch_length, new_len, len(node.mutations)*self.one_mutation), 5)

            node.branch_length = new_len
            node.mutation_length=new_len
            max_bl = max(max_bl, new_len)

        # as branch lengths changed, the params must be fixed
        self.tree.root.up = None
        self.tree.root.dist2root = 0.0
        if max_bl>0.15 and mode=='joint':
            self.logger("TreeAnc.optimize_branch_length: THIS TREE HAS LONG BRANCHES."
                        " \n\t ****TreeTime IS NOT DESIGNED TO OPTIMIZE LONG BRANCHES."
                        " \n\t ****PLEASE OPTIMIZE BRANCHES WITH ANOTHER TOOL AND RERUN WITH"
                        " \n\t ****branch_length_mode='input'", 0, warn=True)
        self._prepare_nodes()
        return ttconf.SUCCESS
def optimize_branch_length(self, mode='joint', **kwargs):
        """
        Perform optimization for the branch lengths of the entire tree.
        This method only does a single path and needs to be iterated.

        **Note** this method assumes that each node stores information
        about its sequence as numpy.array object (node.sequence attribute).
        Therefore, before calling this method, sequence reconstruction with
        either of the available models must be performed.

        Parameters
        ----------

         mode : str
            Optimize branch length assuming the joint ML sequence assignment
            of both ends of the branch (:code:`joint`), or trace over all possible sequence
            assignments on both ends of the branch (:code:`marginal`) (slower, experimental).

         **kwargs :
            Keyword arguments

        Keyword Args
        ------------

         verbose : int
            Output level

         store_old : bool
            If True, the old lengths will be saved in :code:`node._old_dist` attribute.
            Useful for testing, and special post-processing.


        """

        self.logger("TreeAnc.optimize_branch_length: running branch length optimization in mode %s..."%mode,1)
        if (self.tree is None) or (self.aln is None):
            self.logger("TreeAnc.optimize_branch_length: ERROR, alignment or tree are missing", 0)
            return ttconf.ERROR

        store_old_dist = False

        if 'store_old' in kwargs:
            store_old_dist = kwargs['store_old']

        if mode=='marginal':
            # a marginal ancestral reconstruction is required for
            # marginal branch length inference
            if not hasattr(self.tree.root, "marginal_profile"):
                self.infer_ancestral_sequences(marginal=True)

        max_bl = 0
        for node in self.tree.find_clades(order='postorder'):
            if node.up is None: continue # this is the root
            if store_old_dist:
                node._old_length = node.branch_length

            if mode=='marginal':
                new_len = self.optimal_marginal_branch_length(node)
            elif mode=='joint':
                new_len = self.optimal_branch_length(node)
            else:
                self.logger("treeanc.optimize_branch_length: unsupported optimization mode",4, warn=True)
                new_len = node.branch_length

            if new_len < 0:
                continue

            self.logger("Optimization results: old_len=%.4e, new_len=%.4e, naive=%.4e"
                   " Updating branch length..."%(node.branch_length, new_len, len(node.mutations)*self.one_mutation), 5)

            node.branch_length = new_len
            node.mutation_length=new_len
            max_bl = max(max_bl, new_len)

        # as branch lengths changed, the params must be fixed
        self.tree.root.up = None
        self.tree.root.dist2root = 0.0
        if max_bl>0.15 and mode=='joint':
            self.logger("TreeAnc.optimize_branch_length: THIS TREE HAS LONG BRANCHES."
                        " \n\t ****TreeTime IS NOT DESIGNED TO OPTIMIZE LONG BRANCHES."
                        " \n\t ****PLEASE OPTIMIZE BRANCHES WITH ANOTHER TOOL AND RERUN WITH"
                        " \n\t ****branch_length_mode='input'", 0, warn=True)
        self._prepare_nodes()
        return ttconf.SUCCESS
Perform optimization for the branch lengths of the entire tree.
        This method only does a single path and needs to be iterated.

        **Note** this method assumes that each node stores information
        about its sequence as numpy.array object (node.sequence attribute).
        Therefore, before calling this method, sequence reconstruction with
        either of the available models must be performed.

        Parameters
        ----------

         mode : str
            Optimize branch length assuming the joint ML sequence assignment
            of both ends of the branch (:code:`joint`), or trace over all possible sequence
            assignments on both ends of the branch (:code:`marginal`) (slower, experimental).

         **kwargs :
            Keyword arguments

        Keyword Args
        ------------

         verbose : int
            Output level

         store_old : bool
            If True, the old lengths will be saved in :code:`node._old_dist` attribute.
            Useful for testing, and special post-processing.