starttree.h

//
//  Copyright (C) 2020, James Barbetti.
//
//  LICENSE:
//* This program is free software; you can redistribute it and/or modify
//* it under the terms of the GNU General Public License as published by
//* the Free Software Foundation; either version 2 of the License, or
//* (at your option) any later version.
//*
//* This program is distributed in the hope that it will be useful,
//* but WITHOUT ANY WARRANTY; without even the implied warranty of
//* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//* GNU General Public License for more details.
//*
//* You should have received a copy of the GNU General Public License
//* along with this program; if not, write to the
//* Free Software Foundation, Inc.,
//* 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
//

#ifndef starttree_h
#define starttree_h

#include <string>
#include <map>
#include <iostream>
#include <sstream>                //for std::stringstream
#include <vector>
#include <utils/timeutil.h>       //for getRealTime()
#include <utils/vectortypes.h>    //for StrVector template class

namespace StartTree
{
    class BuilderInterface
    {
    public:
        BuilderInterface              () = default;
        virtual ~BuilderInterface     () = default;
        /**
         * @brief indicate if the tree building algorithm does benchmarking
         * @return true  - if this is a benchmarking algorithm
         *                 (which shouldn't report progress, because 
         *                 it is doing the same thing over and over
         *                 with different algorithms).
         * @return false - otherwise
         */
        virtual bool isBenchmark      () const = 0;

        /**
         * @brief  specify whether output (to newick format files)
         *         should be compressed.
         * @param  zipIt  true to turn on compression, false to turn it off
         * @return true  if the setting will be honoured
         * @return false false if won't
         */
        virtual bool setZippedOutput  ( bool zipIt) = 0;
        /**
         * @brief Indicate whether the tree to be generated should
         *        have a root node (normally tree construction 
         *        algorithms generated unrooted trees, with the 
         *        first node in the output file having three 
         *        neighbours): setting this to true asks that the
         *        tree construction algorithm output arooted tree
         *        whose first node has only two neighbours
         * @param rootIt true to request a rooted tree
         * @return true  if the request will be honoured.
         * @return false if not.
         */
        virtual bool setIsRooted      ( bool rootIt) = 0;

        /**
         * @brief Indicate whether the tree to be generated next
         *        is a subtree of another tree (in which case 
         *        it won't need enclosing parentheses and a 
         *        terminating semi-colon).
         *        This is used when tree files are being 
         *        generated by appending the results of 
         *        phylogenetic inferences, of multiple subtrees, 
         *        to a single newick format file.
         * 
         * @param wantSubtree 
         * @return true   if the tree generation algorithm will
         *                do as it is asked.          
         * @return false  if not.
         */
        virtual bool setSubtreeOnly   ( bool wantSubtree) = 0;

        /**
         * @brief Construct a phylogenetic tree, by loading sequence
         *        names, and distances (in) a distance matrix, from a
         *        Phylip distance matrix file, generating a 
         *        phylogenetic tree, based on that input, and writing
         *        it to the file indicated by newickTreeFilePath.
         * 
         * @param distanceMatrixFilePath path of input distance matrix file
         * @param newickTreeFilePath     path of output newick format file
         * @return true  on success
         * @return false on failure
         */
        virtual bool constructTree
            ( const std::string & distanceMatrixFilePath
            , const std::string & newickTreeFilePath) = 0;
        
        /**
         * @brief construct a tree, from sequences and a distance matrix
         *        and write the tree to the file indicated by 
         *        newickTreeFilePath.
         * 
         * @param sequenceNames    (vector of string) the names of the n sequences
         * @param distanceMatrix   pointer to an n*n matrix of distances 
         *                         between sequences.
         * @param newickTreeFilePath path of output newick format file
         * @return true  success
         * @return false an error occurred
         */
        virtual bool constructTreeInMemory
            ( const StrVector&   sequenceNames
            , const double*      distanceMatrix
            , const std::string& newickTreeFilePath) = 0;


        /**
         * @brief construct a tree, from sequences and a distance matrix
         *        and write the tree to the file indicated by 
         *        newickTreeFilePath.
         * 
         * @param sequenceNames    (vector of string) the names of the n sequences
         * @param distanceMatrix   pointer to an n*n matrix of distances 
         *                         between sequences.
         * @param output_string    (output) string, in which to return the 
         *                         newick representation of the phylogenetic tree
         * @return true  success
         * @return false an error occurred
         */
        virtual bool constructTreeStringInMemory
            ( const StrVector& sequenceNames
            , const double*    distanceMatrix
            , std::string&     output_string) = 0;

        /**
         * @brief construct a tree, from sequences and a distance matrix
         *        and write the tree to the file indicated by 
         *        newickTreeFilePath.
         * 
         * @param sequenceNames  (vector of string) the names of the n sequences
         * @param distanceMatrix pointer to an n*n matrix of distances 
         *                       between sequences.
         * @param treeStream     (output) a stream, to which to write the
         *                       newick format representation of the
         *                       phylogenetic tree.
         * @return true  success
         * @return false an error occurred
         */
        virtual bool constructTreeAndAppendToStream
            ( const StrVector &sequenceNames
            , const double *distanceMatrix
            , std::iostream& treeStream) = 0;

        /**
         * @brief Return the name of the algorithm
         * @return const std::string& (output) the name
         */
        virtual const std::string& getName() const = 0;

        /**
         * @brief Return a description of the algorithm
         * @return const std::string& (output) the description
         */
        virtual const std::string& getDescription() = 0;

        /**
         * @brief indicate whether the algorithm is to append its output file
         *        (rather than truncating it before writing to it)
         * @param appendIt - true if the algorithm is to append rather than 
         *                   truncate and replace the existing content of the
         *                   file (if any).
         * @return true  if the request will be honoured
         * @return false if not
         */
        virtual bool setAppendFile(bool appendIt) = 0; //returns true if supported

        /**
         * @brief tell the distance matrix algorithm to avoid writing anything
         *        to standard output (it may still write to standard error if 
         *        something goes wrong).
         */
        virtual void beSilent() = 0;

        /**
         * @brief set the precision (number of digits after the decimal point)
         *        that the algorithm is to use when writing out distances 
         *        between nodes in the phylogenetic tree, that it calculates.
         * @param precision how many digits of precision
         * @return true   if the algorithm will use that precision
         * @return false  otherwise
         */
        virtual bool setPrecision(int precision) = 0;
    };

    class BenchmarkingTreeBuilder;

    class Registry
    {
        friend class BenchmarkingTreeBuilder;
        
    private:
        std::map<std::string, BuilderInterface*> mapOfTreeBuilders;
        //Note: Owned by the Registry, and will be deleted in ~Registry.
        std::string nameOfDefaultTreeBuilder;
        
    protected:
        Registry();
        ~Registry();
        /**
         * @brief  return the number of tree-construction algorithms that
         *         are avaiable in the registry.
         * @return size_t how many
         */
        size_t getBuilderCount();
        /**
         * @brief Add a tree construction algorithm implementation to
         *        the registry.
         * @param name     the name under which to register the 
         *                 distance matrix phylogenetic tree 
         *                 inference implementation
         * @param builder  a pointer to the phylogenetic tree
         *                 inference implementation.
         */
        void addBuilder(const std::string& name, BuilderInterface* builder);

        /**
         * @brief Look up a distance-matrix phylogenetic inference 
         *        implementation by name.
         * @param name the name
         * @return BuilderInterface* 
         */
        BuilderInterface* getBuilder(const std::string& name) const;
        BuilderInterface* getBuilder(const char* name) const;

        /**
         * @brief  Return the default distance matrix phylogenetic 
         *         inference implementation.
         * @return BuilderInterface* a pointer to the requested 
         *         implementation (if it is registered, with that 
         *         name, in the registry).
         */
        BuilderInterface* getDefaultTreeBuilder() const;

    public:
        /**
         * @brief  return the singleton Registry instance.
         * @return Registry& a reference to the singleton instance.
         */
        static Registry& getInstance();

        /**
         * @brief Add a distance matrix phylogenetic inference implementation
         *        to the registry
         * @param builder The implementation
         * @note  The registry takes ownership of the implementation.
         *        If another implementation, p, is already registered, with the 
         *        same name, builder replaces it.  If p!=builder, p is deleted.
         *        (as it is owned by the registry).
         */
        void   advertiseTreeBuilder(BuilderInterface* builder);

        /**
         * @brief set the name for the default tree builder
         * @param name 
         */
        void   setNameOfDefaultTreeBuilder(const char* name);

        /**
         * @brief  returns the name of the default distance matrix phylogenetic
         *         inference implementation, in the registry (or blank if that
         *         has not been set).
         * @return const std::string& 
         */
        static const std::string& getNameOfDefaultTreeBuilder();

        /**
         * @brief  look up a distance matrix phylogenetic inference
         *         implementation by name.
         * @param  name 
         * @return BuilderInterface* - the implementation registered
         *         with that name, if there is one.  nullptr if there
         *         isn't one.
         */
        static BuilderInterface* getTreeBuilderByName(const std::string& name);
        static BuilderInterface* getTreeBuilderByName(const char* name);

        /**
         * @brief  return a list of the names and 
         *         descriptions of the distance matrix phylogenetic 
         *         inference implementations that are registered 
         *         (in lexicographic order).
         * @return std::string
         * @note   format is [name]: [description]<LineFeed>
         *         for each implementation.
         *         
         */
        std::string getListOfTreeBuilders() const;

        /**
         * @brief  Return a vector, containing the names of the registered 
         *         distance matrix phylogenetic inference implementations,
         *         in lexicographic order.
         * @param  withDescriptions if true, indicates that each string in
         *         the vector is to be in the format "[name]: [description]"
         * @return StrVector - the vector
         */
        StrVector getVectorOfTreeBuilderNames(bool withDescriptions) const;
    };

    /**
     * @brief  Wrapper for distance matrix phylogenetic inference
     *         implementations (which remembers name, description,
     *         and other state information, for the implementation
     *         it wraps.
     * @tparam B - the implementation that is to be wrapped so. 
     * 
     * @note   The wrapper provides diagnostic information (like, 
     *         how long did phylogenetic inference take, how much CPU
     *         did it use) (unless it is told to be quiet!)
     */
    template <class B> class Builder: public BuilderInterface
    {
        //Note: B must have:
        //      1. a constructor that takes the name of an ".mldist"
        //         distance matrix file as a parameter;
        //      2. a constructTree()   member function; 
        //      3. a writeTreeFile()   member function;
        //      4. a setZippedOutput() member function;
        //      5. a setAppendFile()   member function;
        //      6. a setIsRooted()     member function;
        //      7. a setSubtreeOnly()  member function; and
        //      8. a calculateRMSOfTMinusD() member function.
        //
    protected:
        const std::string name;        //name under which B will be registered
        const std::string description; //description
        bool  silent;                  //whether diagnostic and timing output is to 
                                       //be suppressed (true if so)
        bool  isOutputToBeAppended;    //whether output files are to
        bool  isOutputToBeZipped;      //whether output files are to be compressed
        int   precision;               //precision to pass down
        bool  isRooted;                //indicates whether a rooted tree (where the 
                                       //"root" node has two neighbours), or an
                                       //unrooted tree (where the "root" has three)
                                       //is wanted.
        bool  subtreeOnly;             //indicates whethere enclosing ( ) ; 
                                       //characters are wanted in the Newick 
                                       //tree string (or file content) that will be
                                       //calculated.

        /**
         * @brief  wraps up phylogenetic tree construction (adding the
         *         logging of timing information)
         * @param  builder 
         * @return true (always)
         */
        bool  constructTreeWith(B& builder) {
            double buildStart    = getRealTime();
            double buildStartCPU = getCPUTime();
            if (silent) {
                builder.beSilent();
            }
            builder.setIsRooted(isRooted);
            builder.constructTree();
            double buildElapsed = getRealTime() - buildStart;
            double buildCPU     = getCPUTime() - buildStartCPU;
            if (!silent) {
                std::cout.precision(6);
                std::cout << "Computing "
                          << name << " tree took " 
                          << buildElapsed << " sec (of wall-clock time) " 
                          << buildCPU << " sec (of CPU time)" ;
                if (0<buildElapsed) {
                    double percentCPU = (buildCPU/buildElapsed)*100.0;
                    std::cout << "(" << floor(percentCPU) << "%)";
                }
                std::cout << std::endl;
                std::cout.precision(3);
            }
            return true;
        }
    public:
        Builder(const char* nameToUse, const char *descriptionToGive)
        : name(nameToUse), description(descriptionToGive), silent(false)
        , isOutputToBeAppended(false), isOutputToBeZipped(false)
        , precision(6), subtreeOnly(false) {
        }

        /**
         * @brief Turn on silent mode (where timing and other
         *        diagnostic information will not be logged).
         */
        virtual void beSilent() override {
            silent = true;
        }

        /**
         * @brief  return the registered name of the algorithm
         * @return const std::string& - takes the name
         */
        virtual const std::string& getName() const override{
            return name;
        }

        /**
         * @brief  get the algorithn's description
         * @return const std::string& - takes the description
         */
       virtual const std::string& getDescription() override {
            return description;
        }

        /**
         * @brief  indicate whether this is a multiple-algorithm 
         *         benchmarking distance matrix phylogenetic inference
         *         implementation (it never is!).
         * @return false (this class isn't ever used to
         *         benchmark multiple algorithms)
         */
        virtual bool isBenchmark() const override {
            return false;
        }

        /**
         * @brief  set whether compression is to be used when writing
         *         Newick format output files.
         * @param  zipIt boolean (true if compression is to be used,
         *         false if not))
         * @return true 
         * @note   compression is passed to the B instance by calling
         *         its setZippedOutput() method after a tree is
         *         determined.
         */
        virtual bool setZippedOutput(bool zipIt) override {
            isOutputToBeZipped = zipIt;
            return true;
        }

        /**
         * @brief  set whether files are to be appended 
         *         (rather than truncated)
         *         when B is asked to write to them 
         * @param appendIt true if they are, false if not
         * @return always returns true
         * @note   this is passed on to the B instance by calling
         *         its setZippedOutput() method after a tree is
         *         determined.
         */
        virtual bool setAppendFile(bool appendIt) override {
            isOutputToBeAppended = appendIt;
            return true;
        }

        /**
         * @brief  set whether rooted trees (with the root having
         *         two neighboours, rather than three) are to be
         *         constructed, by B
         * @param  rootIt true if trees are to be rooted, 
         *         false if they are to be unrooted.
         * @return always returns true
         * @note   this is passed on to the B instance by calling
         *         its setIsRooted() method just before a 
         *         tree is determined.
         */
        virtual bool setIsRooted(bool rootIt) override {
            isRooted = rootIt;
            return true;
        }

        /**
         * @brief  set whether B is to write trees (which are enclosed
         *         in brackets, and have a trailing semi-colon), or
         *         subtrees (which are not so enclosed, and lack a 
         *         semi-colon) to its output files, streams, or strings.
         * @param  wantSubtree true for subtrees, false for trees
         * @return always returns true
         * @note   this is passed on to the B instance by calling
         *         its setSubtreeOnly() method just before a 
         *         tree output is written to a file, stream, or string.
         */
        virtual bool setSubtreeOnly(bool wantSubtree) override {
            subtreeOnly = wantSubtree;
            return true;
        }

        /**
         * @brief Generate (or append!) an output Newick Tree file
         *        with a representation of the phylogenetic tree,
         *        constructed by an instance of B, from the sequences
         *        and distance matrix, read from the specified 
         *        phylip format distanc matrix (input) file.
         * 
         * @param distanceMatrixFilePath file path of phylip format 
         *                               distance matrix file (input)
         * @param newickTreeFilePath     file path of newick format
         *                               phylogenetic tree file (output)
         * @return true if the distance matrix could be loaded,
         *         the phylogenetic tree inferred, and the tree 
         *         written to the output file.  false otherwise.
         * @note   if the newickTreeFilePath is empty, no tree 
         *         is written to a file, and true is returned 
         *         (if the matrix loaded and a tree was constructed).
         */
        virtual bool constructTree
            ( const std::string &distanceMatrixFilePath
            , const std::string & newickTreeFilePath) override {
            B builder;
            if (silent) {
                builder.beSilent();
            }
            if (!builder.loadMatrixFromFile(distanceMatrixFilePath)) {
                return false;
            }
            builder.setIsRooted(isRooted);
            if (!constructTreeWith(builder)) {
                return false;
            }
            if (newickTreeFilePath.empty()) {
                return true;
            }
            builder.setZippedOutput(isOutputToBeZipped);
            builder.setAppendFile(isOutputToBeAppended);
            builder.setSubtreeOnly(subtreeOnly);
            return builder.writeTreeFile(precision, newickTreeFilePath);
        }

        /**
         * @brief  from a string vector of n sequence names, and a
         *         distance matrix (a pointer to an array of n*n 
         *         doubles, in row-major order), have B determine
         *         a phylogenetic tree.
         * @param  sequenceNames   a vector of sequence names
         *                         (these don't actually have to be
         *                         unique! but: from ambiguous input,
         *                         you will get an ambiguous output!)
         * @param  distanceMatrix  a pointer to the first element in
         *                         an n*n-element const double array
         *                         of pairwise distances between taxa
         *                         (in row major order) (n is
         *                         sequenceNames.size()).
         * @param  builder         a reference to an instance of B.
         * @return true if the matrix could be loaded and the tree
         *         constructed.
         * @note   this is the only constructTree() - or similar - 
         *         member function that can display the root mean square
         *         error of the difference, between the input distance matrix
         *         and the distance matrix implicit in the output tree
         *         (calculated by adding distances along the path from each
         *         between each pair of leaf taxa).  The others do not have 
         *         the input distance matrix, and so cannot "ask" for it.
         */
        bool constructTree(const StrVector &sequenceNames
                          , const double *distanceMatrix, B& builder) {
            if (silent) {
                builder.beSilent();
            }
            if (!builder.loadMatrix(sequenceNames, distanceMatrix)) {
                return false;
            }
            if (!constructTreeWith(builder)) {
                return false;
            }
            if (!silent) {
                double rms = 0;
                if (builder.calculateRMSOfTMinusD(distanceMatrix, 
                                                  sequenceNames.size(), rms)) {
                    std::cout << "Root Mean Square Error was " << rms << std::endl;
                }
            }
            builder.setZippedOutput(isOutputToBeZipped);
            builder.setAppendFile(isOutputToBeAppended);
            builder.setSubtreeOnly(subtreeOnly);
            return true;        
        }

        /**
         * @brief  from a string vector of n sequence names, and a
         *         distance matrix (a pointer to an array of n*n 
         *         doubles, in row-major order), have a (local 
         *         instance of) B determine a phylogenetic tree,
         *         and write it to the specified newick format 
         *         output file.
         * 
         * @param  sequenceNames   a vector of sequence names
         *                         (these don't actually have to be
         *                         unique! but: from ambiguous input,
         *                         you will get an ambiguous output!)
         * @param  distanceMatrix  a pointer to the first element in
         *                         an n*n-element const double array
         *                         of pairwise distances between taxa
         *                         (in row major order) (n is
         *                         sequenceNames.size()).
         * @param  newickTreeFilePath  the path of the output file to write
         * @return true if the instance of B could construct the tree,
         *         and write it to the specified output file. false otherwise.
         */
        virtual bool constructTreeInMemory
            ( const StrVector &sequenceNames
            , const double *distanceMatrix
            , const std::string & newickTreeFilePath) override {
            B builder;
            if (!constructTree(sequenceNames, distanceMatrix, builder)) {
                return false;
            }
            if (newickTreeFilePath.empty()) {
                return true;
            }
            builder.setZippedOutput(isOutputToBeZipped);
            builder.setAppendFile(isOutputToBeAppended);
            builder.setSubtreeOnly(subtreeOnly);
            return builder.writeTreeFile(precision, newickTreeFilePath);
        }

        /**
         * @brief  from a string vector of n sequence names, and a
         *         distance matrix (a pointer to an array of n*n 
         *         doubles, in row-major order), have a (local 
         *         instance of) B determine a phylogenetic tree,
         *         and write a newick-format description of it to a
         *         string.
         * 
         * @param  sequenceNames   a vector of sequence names
         *                         (these don't actually have to be
         *                         unique! but: from ambiguous input,
         *                         you will get an ambiguous output!)
         * @param  distanceMatrix  a pointer to the first element in
         *                         an n*n-element const double array
         *                         of pairwise distances between taxa
         *                         (in row major order) (n is
         *                         sequenceNames.size()).
         * @param  output_string   the string to hold the Newick format 
         *                         description of the phylogenetic tree.
         * @return true if the instance of B could construct the tree,
         *         and write it to output_string.  false otherwise.
         * @note   this was added so that pydecenttree could return
         *         strings to python clients.
         */
         virtual bool constructTreeStringInMemory
            ( const StrVector& sequenceNames
            , const double*    distanceMatrix
            , std::string&     output_string) override {
            B builder;
            constructTree(sequenceNames, distanceMatrix, builder);
            std::stringstream stream;
            stream.precision(precision);
            builder.setSubtreeOnly(subtreeOnly);
            bool rc =  builder.writeTreeToOpenFile(stream);
            output_string = stream.str();
            return rc;
        }

        /**
         * @brief  from a string vector of n sequence names, and a
         *         distance matrix (a pointer to an array of n*n 
         *         doubles, in row-major order), have a (local 
         *         instance of) B determine a phylogenetic tree,
         *         and append a newick-format description of it to
         *         an open output stream.
         * 
         * @param  sequenceNames   a vector of sequence names
         *                         (these don't actually have to be
         *                         unique! but: from ambiguous input,
         *                         you will get an ambiguous output!)
         * @param  distanceMatrix  a pointer to the first element in
         *                         an n*n-element const double array
         *                         of pairwise distances between taxa
         *                         (in row major order) (n is
         *                         sequenceNames.size()).
         * @param  stream          the stream, to which the Newick format 
         *                         description of the phylogenetic tree
         *                         is to be appended.
         * @return true if the instance of B could construct the tree,
         *         and append it to stream.  false otherwise.
         */
        virtual bool constructTreeAndAppendToStream
            ( const StrVector &sequenceNames
            , const double* distanceMatrix
            , std::iostream& stream) override {
            B builder;
            constructTree(sequenceNames, distanceMatrix, builder);
            stream.precision(precision);
            builder.setSubtreeOnly(subtreeOnly);
            return builder.writeTreeToOpenFile(stream);
        }

        /**
         * @brief  Set the precision (number of decimal places to have
         *         in distances, in the Newick format tree).
         * @param  precision_to_use 
         * @return true - always works
         * @note   precision is either set by changing precision on the
         *         output stream (when the output is being written), or
         *         (less often) by passing it to B's writeTreeFile()
         *         member function, which is expected to take precision
         *         as its first parameter.
         */
        virtual bool setPrecision(int precision_to_use) override {
            precision = precision_to_use;
            return true;
        }
    };

    class BenchmarkingTreeBuilder: public BuilderInterface
    {
    protected:
        const std::string name;
        const std::string description;
        std::vector<BuilderInterface*> builders;
        bool isOutputToBeZipped;
        bool silent;
        int  precision;
        bool subtreeOnly;
        bool isRooted;
    public:
        BenchmarkingTreeBuilder(Registry& f, const char* nameToUse, 
                                const char *descriptionToGive);
        virtual const std::string& getName() const override;
        virtual const std::string& getDescription() override;
        virtual bool isBenchmark() const override;
        virtual bool constructTree
            ( const std::string& distanceMatrixFilePath
            , const std::string& newickTreeFilePath) override;
        virtual bool constructTreeInMemory
            ( const StrVector &sequenceNames
            , const double *distanceMatrix
            , const std::string& newickTreeFilePath) override;
        virtual bool constructTreeAndAppendToStream
            ( const StrVector& sequenceNames
            , const double*    distanceMatrix
            , std::iostream&   stream) override;
        virtual bool constructTreeStringInMemory
            ( const StrVector& sequenceNames
            , const double*    distanceMatrix
            , std::string&     output_string) override;
        virtual bool setZippedOutput (bool zipIt) override;
        virtual void beSilent        () override;
        virtual bool setPrecision    (int precisionToUse) override;
        virtual bool setAppendFile   (bool appendIt)      override;
        virtual bool setIsRooted     (bool rootIt)        override;
        virtual bool setSubtreeOnly  (bool wantSubtree)   override;
    };
}

#define START_TREE_RECOGNIZED(name) \
    ( StartTree::Registry::getTreeBuilderByName(name) != nullptr )

#endif /* starttree_h */