C++ TopologyNode::setAge方法代码示例

/**
 * Perform the proposal.
 *
 * A Uniform-simplex proposal randomly changes some values of a simplex, although the other values
 * change too because of the renormalization.
 * First, some random indices are drawn. Then, the proposal draws a new somplex
 *   u ~ Uniform(val[index] * alpha)
 * where alpha is the tuning parameter.The new value is set to u.
 * The simplex is then renormalized.
 *
 * \return The hastings ratio.
 */
double RootTimeSlideUniformProposal::doProposal( void )
{
    
    // Get random number generator
    RandomNumberGenerator* rng     = GLOBAL_RNG;
    
    Tree& tau = variable->getValue();
    
    // pick a random node which is not the root and neithor the direct descendant of the root
    TopologyNode* node = &tau.getRoot();
    
    // we need to work with the times
    double my_age      = node->getAge();
    double child_Age   = node->getChild( 0 ).getAge();
    if ( child_Age < node->getChild( 1 ).getAge())
    {
        child_Age = node->getChild( 1 ).getAge();
    }
    
    // now we store all necessary values
    storedAge = my_age;
    
    // draw new ages and compute the hastings ratio at the same time
    double my_new_age = (origin->getValue() - child_Age) * rng->uniform01() + child_Age;
    
    // set the age
    node->setAge( my_new_age );
    
    return 0.0;
    
}

void AbstractRootedTreeDistribution::simulateTree( void )
{

    // the time tree object (topology & times)
    Tree *psi = new Tree();

    // internally we treat unrooted topologies the same as rooted
    psi->setRooted( true );

    // create the tip nodes
    std::vector<TopologyNode*> nodes;
    for (size_t i=0; i<num_taxa; ++i)
    {

        // create the i-th taxon
        TopologyNode* node = new TopologyNode( taxa[i], i );

        // set the age of this tip node
        node->setAge( taxa[i].getAge() );

        if (node->getAge() > 0)
        {
            node->setFossil(true);
        }

        // add the new node to the list
        nodes.push_back( node );

    }


    double ra = root_age->getValue();
    double present = ra;

    // we need a sorted vector of constraints, starting with the smallest
    std::vector<Clade> sorted_clades;

    std::vector<Clade> constraints;

    for (size_t i = 0; i < constraints.size(); ++i)
    {
        if (constraints[i].getAge() > ra)
        {
            throw RbException("Cannot simulate tree: clade constraints are older than the root age.");
        }

        // set the ages of each of the taxa in the constraint
        for (size_t j = 0; j < constraints[i].size(); ++j)
        {
            for (size_t k = 0; k < num_taxa; k++)
            {
                if ( taxa[k].getName() == constraints[i].getTaxonName(j) )
                {
                    constraints[i].setTaxonAge(j, taxa[k].getAge());
                    break;
                }
            }
        }

        if ( constraints[i].size() > 1 && constraints[i].size() < num_taxa )
        {
            sorted_clades.push_back( constraints[i] );
        }
    }

    // create a clade that contains all species
    Clade all_species = Clade(taxa);
    all_species.setAge( ra );
    sorted_clades.push_back(all_species);

    // next sort the clades
    std::sort(sorted_clades.begin(),sorted_clades.end());

    // remove duplicates
    std::vector<Clade> tmp;
    tmp.push_back( sorted_clades[0] );
    for (size_t i = 1; i < sorted_clades.size(); ++i)
    {
        Clade &a = tmp[tmp.size()-1];
        Clade &b = sorted_clades[i];

        if ( a.size() != b.size() )
        {
            tmp.push_back( sorted_clades[i] );
        }
        else
        {
            bool equal = true;
            for (size_t j = 0; j < a.size(); ++j)
            {
                if ( a.getTaxon(j) != b.getTaxon(j) )
                {
                    equal = false;
                    break;
                }
            }
            if ( equal == false )
            {
                tmp.push_back( sorted_clades[i] );
            }
//.........这里部分代码省略.........

void AbstractRootedTreeDistribution::simulateClade(std::vector<TopologyNode *> &n, double age, double present)
{

    // Get the rng
    RandomNumberGenerator* rng = GLOBAL_RNG;

    // get the minimum age
    double current_age = n[0]->getAge();
    for (size_t i = 1; i < n.size(); ++i)
    {

        if ( current_age > n[i]->getAge() )
        {
            current_age = n[i]->getAge();
        }

    }
    

    while ( n.size() > 2 && current_age < age )
    {

        // get all the nodes before the current age
        std::vector<TopologyNode*> active_nodes;
        for (size_t i = 0; i < n.size(); ++i)
        {

            if ( current_age >= n[i]->getAge() )
            {
                active_nodes.push_back( n[i] );
            }

        }

        // we need to get next age of a node larger than the current age
        double next_node_age = age;
        for (size_t i = 0; i < n.size(); ++i)
        {

            if ( current_age < n[i]->getAge() && n[i]->getAge() < next_node_age )
            {
                next_node_age = n[i]->getAge();
            }

        }

        // only simulate if there are at least to valid/active nodes
        if ( active_nodes.size() < 2 )
        {
            current_age = next_node_age;
        }
        else
        {

            // now we simulate new ages
            double next_sim_age = simulateNextAge(active_nodes.size()-1, present-age, present-current_age, present);

            if ( next_sim_age < next_node_age )
            {

                // randomly pick two nodes
                size_t index_left = static_cast<size_t>( floor(rng->uniform01()*active_nodes.size()) );
                TopologyNode* left_child = active_nodes[index_left];
                active_nodes.erase(active_nodes.begin()+long(index_left));
                size_t index_right = static_cast<size_t>( floor(rng->uniform01()*active_nodes.size()) );
                TopologyNode* right_right = active_nodes[index_right];
                active_nodes.erase(active_nodes.begin()+long(index_right));

                // erase the nodes also from the origin nodes vector
                n.erase(std::remove(n.begin(), n.end(), left_child), n.end());
                n.erase(std::remove(n.begin(), n.end(), right_right), n.end());


                // create a parent for the two
                TopologyNode *parent = new TopologyNode();
                parent->addChild( left_child );
                parent->addChild( right_right );
                left_child->setParent( parent );
                right_right->setParent( parent );
                parent->setAge( next_sim_age );

                // insert the parent to our list
                n.push_back( parent );
                
                current_age = next_sim_age;
            }
            else
            {
                current_age = next_node_age;
            }

        }
        
    }


    if ( n.size() == 2 )
    {

        // pick two nodes
//.........这里部分代码省略.........