57 bool const fully_resolve,
58 std::string
const& name_prefix
73 std::string
const& name_prefix
86 throw std::runtime_error(
87 "Tree provided for labelled_tree() is not compatible " 88 "with the Tree of the provided Sample." 93 std::vector< tree::TreeEdge* > edge_list;
94 edge_list.reserve( result.edge_count() );
95 for(
auto& edge : result.edges() ) {
96 edge_list.push_back( &edge );
115 auto placement_pairs_per_edge = [](
Sample const& sample ){
116 std::vector< std::vector< PlacementPair >> place_map;
118 for(
auto const& pqry : sample.
pqueries() ) {
121 auto max_it = std::max_element(
122 pqry.placements().begin(),
123 pqry.placements().end(),
131 if( max_it != pqry.placements().end() ) {
132 place_map[ max_it->edge().index() ].emplace_back( &pqry, &*max_it );
150 auto place_map = placement_pairs_per_edge( sample );
151 assert( place_map.size() == edge_list.size() );
157 auto sort_placements = fully_resolve
158 ? []( PlacementPair
const& lhs, PlacementPair
const& rhs ){
159 return lhs.placement->proximal_length < rhs.placement->proximal_length;
161 : []( PlacementPair
const& lhs, PlacementPair
const& rhs ){
162 return lhs.placement->pendant_length < rhs.placement->pendant_length;
165 for(
auto& edge_l : place_map ) {
166 std::sort( edge_l.begin(), edge_l.end(), sort_placements );
174 auto add_lonely_placement = [ &name_prefix ](
177 PlacementPair
const& placement_pair
181 auto& pendant_edge = pendant_node.link().edge();
182 auto& proximal_edge = pendant_edge.primary_link().next().edge();
183 auto& distal_edge = pendant_edge.primary_link().next().next().edge();
186 assert(
degree( pendant_edge.primary_node() ) == 3 );
187 assert(
is_leaf( pendant_edge.secondary_node() ));
195 assert( pend_data.branch_length == 0.0 );
199 double const distal_len = prox_data.
branch_length - placement_pair.placement->proximal_length;
200 pend_data.branch_length = placement_pair.placement->pendant_length;
201 prox_data.
branch_length = placement_pair.placement->proximal_length;
206 assert( placement_pair.pquery->name_size() <= 1 );
207 if( placement_pair.pquery->name_size() == 1 ) {
209 pendant_node_data.
name = name_prefix + placement_pair.pquery->name_at(0).name;
218 auto process_edge_fully_resolved = [ &name_prefix ] (
221 std::vector<PlacementPair>
const& placement_pairs
227 double used_length = 0.0;
234 auto insertion_edge = &edge;
237 for(
auto const& placement_pair : placement_pairs ) {
239 for(
auto const& pquery_name : placement_pair.pquery->names() ) {
243 auto& pendant_edge = pendant_node.link().edge();
244 auto& proximal_edge = pendant_edge.primary_link().next().edge();
245 auto& distal_edge = pendant_edge.primary_link().next().next().edge();
248 assert(
degree( pendant_edge.primary_node() ) == 3 );
249 assert(
is_leaf( pendant_edge.secondary_node() ));
258 assert( pend_data.branch_length == 0.0 );
262 assert( placement_pair.placement->proximal_length >= used_length );
265 pend_data.branch_length = placement_pair.placement->pendant_length;
266 prox_data.
branch_length = placement_pair.placement->proximal_length - used_length;
270 node_data.
name = name_prefix + pquery_name.name;
273 used_length = placement_pair.placement->proximal_length;
274 insertion_edge = &distal_edge;
288 assert( end_edge_data.branch_length == 0.0 );
289 if( used_length < orig_length ) {
290 end_edge_data.branch_length = orig_length - used_length;
300 auto process_edge_multifurcating = [ &name_prefix ] (
303 std::vector<PlacementPair>
const& placement_pairs
309 auto& base_edge = new_node.link().edge();
310 auto& pri_edge = base_edge.primary_link().next().edge();
311 auto& sec_edge = base_edge.primary_link().next().next().edge();
314 assert(
degree( base_edge.primary_node() ) == 3 );
315 assert(
is_leaf( base_edge.secondary_node() ));
323 assert( base_data.branch_length == 0.0 );
327 auto const avg_prox_len = std::accumulate(
328 placement_pairs.begin(),
329 placement_pairs.end(),
331 [](
double in, PlacementPair
const& cur ){
332 return in + cur.placement->proximal_length;
334 ) / static_cast<double>( placement_pairs.size() );
349 assert( placement_pairs.size() > 0 );
350 auto const min_pen_len = std::max( 0.0, std::min_element(
351 placement_pairs.begin(),
352 placement_pairs.end(),
353 []( PlacementPair
const& lhs, PlacementPair
const& rhs ){
354 return lhs.placement->pendant_length < rhs.placement->pendant_length;
356 )->placement->pendant_length );
360 assert( min_pen_len == placement_pairs[0].placement->pendant_length );
363 base_data.branch_length = min_pen_len;
366 for(
auto const& placement_pair : placement_pairs ) {
368 for(
auto const& pquery_name : placement_pair.pquery->names() ) {
372 auto& p_edge = p_node.link().edge();
377 assert( placement_pair.placement->pendant_length >= min_pen_len );
378 p_data.branch_length = placement_pair.placement->pendant_length - min_pen_len;
382 p_node_data.
name = name_prefix + pquery_name.name;
392 for(
size_t edge_i = 0; edge_i < edge_list.size(); ++edge_i ) {
394 assert( edge_list[ edge_i ] );
397 if( place_map[ edge_i ].size() == 0 ) {
403 if( place_map[ edge_i ].size() == 1 && place_map[ edge_i ][0].pquery->name_size() <= 1 ) {
404 assert( place_map[ edge_i ][0].pquery && place_map[ edge_i ][0].placement );
405 add_lonely_placement( result, *edge_list[ edge_i ], place_map[ edge_i ][0] );
410 if( fully_resolve ) {
411 process_edge_fully_resolved(
412 result, *edge_list[ edge_i ], place_map[ edge_i ]
415 process_edge_multifurcating(
416 result, *edge_list[ edge_i ], place_map[ edge_i ]
PlacementTree & tree()
Get the PlacementTree of this Sample.
A pquery holds a set of PqueryPlacements and a set of PqueryNames.
tree::Tree labelled_tree(Sample const &sample, bool const fully_resolve, std::string const &name_prefix)
Produce a Tree where the most probable PqueryPlacement of each Pquery in a Sample is turned into an E...
size_t edge_count() const
Return the number of TreeEdges of the Tree.
One placement position of a Pquery on a Tree.
Container namespace for all symbols of genesis in order to keep them separate when used as a library...
utils::Range< iterator_pqueries > pqueries()
Return a Range iterator to the Pqueries .
Class for representing phylogenetic trees.
CommonTree convert_to_common_tree(Tree const &source_tree)
Convert a Tree to a CommonTree with CommonNodeData and CommonEdgeData.
TreeNode & add_new_leaf_node(Tree &tree, TreeEdge &target_edge, std::function< void(TreeEdge &target_edge, TreeEdge &new_edge)> adjust_edges)
Add a new Node as a leaf to an existing Edge, by also adding a new Node in the middle of that Edge...
std::string name
Name of the node.
TreeNode & add_new_node(Tree &tree, TreeNode &target_node)
Add a new Node as a leaf to an existing Node.
double branch_length
Branch length of the edge.
bool identical_topology(Tree const &lhs, Tree const &rhs, bool identical_indices)
Return whether both trees have an identical topology.
Manage a set of Pqueries along with the PlacementTree where the PqueryPlacements are placed on...
Common class containing the commonly needed data for tree nodes.
Common class containing the commonly needed data for tree edges.
size_t degree(TreeLink const &link)
Return the degree of the node for a given TreeLink, i.e. how many neighbouring nodes it has...
double like_weight_ratio
Likelihood weight ratio of this placement.
bool is_leaf(TreeLink const &link)
Return true iff the node of the given link is a leaf node.