Parsing a Newick string

Table of Contents

• Introduction
  • Background
  • Grammar
  • Objectives
• To a Quetzal k-ary tree
  • Default properties
  • Custom properties
• Interfacing legacy code

Introduction

---

Background

Newick tree format is a way of representing graph-theoretical trees with edge lengths using parentheses and commas.

Yet not very efficient, Newick is a simple and popular format for representing trees, and provides a rather useful abstraction in coalescence theory for representing the shared history of populations (species trees, population trees) or gene genealogies that coalesce into these trees.

Quetzal provides some utilities to parse, generate and manipulate such format.

---

Grammar

Whitespace here refer to any of the following: spaces, tabs, carriage returns, and linefeeds.

• Whitespace within number is prohibited.
• Whitespace outside of node names is ignored.
• Grammar characters (semicolon, parentheses, comma, and colon) are prohibited
• Comments are enclosed in square brackets.

---

Objectives

In this tutorial section you will learn how to:

• Parse a Newick string to a Quetzal k-ary tree class,
• Customize the default Quetzal tree properties to describe your data and problem correctly,
• Parse a Newick string to your own legacy tree class

---

To a Quetzal k-ary tree

---

Default properties

In this example we will use the default (and simplistic) properties of the graph synthetized by the parsing. By properties we mean the type the data structures used to store arbitrary information along the vertices and edges of a graph.

For general graphs it could be anything, but since we are parsing Newick formats, we can reasonably expect ad minima:

• a std::string to describe a vertex (the name of the node)
• a double to describe an edge (the distance to a parent node).

Quetzal defaults to these minima when it tries to build a tree from a Newick string.

Remarks

In the next sections we will see that you can actually parse into more complex data structures - as long as vertices (resp. edges) remain constructible from a simple std::string (resp. double)!

Input

#include "quetzal/quetzal.hpp"
#include <cassert>
 
// Define a namespace alias to shorten notation
namespace newick = quetzal::format::newick;
 
int main()
{
    std::string s1 = "(A:0.1,B:0.2,(C:0.3,D:0.4)E:0.5)F;";
    std::string s2 = "(,,(,));";
 
    // For default graph properties, leave <> empty
    auto [tree1, root1] = newick::to_k_ary_tree<>(s1);
    auto [tree2, root2] = newick::to_k_ary_tree<>(s2);
 
    std::cout << "These graphs are " << (tree1.is_isomorphic(tree2) ? "toootally" : "not") << " isomorphic!"
              << std::endl;
 
    return 0;
}

Output

3
2
3
2
These graphs are toootally isomorphic!

---

Custom properties

In the previous example we used default types to represent the properties of a vertex and an edge. These defaults are useful, but not always sufficient. For example, you may want to visit your tree graph and trigger specific events when specific nodes or edges are visited.

In this case, you would need to define your own small structures to represent the graph properties you want:

• the vertex should be constructible from a std::string
• the edge should be constructible from a double
• the event would be to print to std::cout with a pretty display

Remarks

You can think of pretty much any type of events here:

• convert the branch length to coalescence units, generations, etc
• update a mutational state data member,
• construct a slightly modified copy of the tree while visiting its clone
• etc ...

Input

#include "quetzal/quetzal.hpp"
#include <iostream>
 
// Your custom little vertex class
struct vertex_info
{
    std::string name;
};
 
// Your custom little edge class
struct edge_info
{
    double length;
};
 
// Make your custom vertex streamable for the visitor
static inline auto &operator<<(std::ostream &os, vertex_info const &v)
{
    return os << quoted(v.name);
}
 
// we want to trigger a custom event when a vertex is discovered by the DFS
template <class Graph> struct MyVisitor
{
    Graph &_tree;
    MyVisitor(Graph &tree) : _tree(tree)
    {
    }
    void operator()(auto stage, auto vertex)
    {
        switch (stage)
        {
        case boost::visit::pre:
            // vertex_info g[u] is streamable!
            std::cout << "Discover vertex " << vertex << " (" << _tree[vertex] << ")\n";
            break;
        case boost::visit::in:
            break;
        case boost::visit::post:
            break;
        default:
            throw std::invalid_argument("received invalid DFS order");
        }
    }
};
 
// Define a namespace alias to shorten notation
namespace newick = quetzal::format::newick;
 
int main()
{
    std::string s = "(A:0.1,B:0.2,(C:0.3,D:0.4)E:0.5)F;";
 
    // Inject your custom littles classes here: the parser will generate these
    // structures and embed them into the graph.
    auto [tree, root] = newick::to_k_ary_tree<vertex_info, edge_info>(s);
 
    // depth-first traversal of the vertices in the directed graph
    MyVisitor visitor(tree);
    tree.depth_first_search(root, visitor);
 
    return 0;
}

Output

3
2
Discover vertex 0 ("F")
Discover vertex 1 ("A")
Discover vertex 2 ("B")
Discover vertex 3 ("E")
Discover vertex 4 ("C")
Discover vertex 5 ("D")

---

Interfacing legacy code

If you have been coding for a while, you may already depend heavily on your own class to describe a tree graph, and not feel very excited about refactoring your whole project to switch to Quetzal trees. Fair enough.

This section shows how to use Quetzal parser to populate your own legacy class and then forget about Quetzal.

1. You will need to parse Newick strings into an AST (Abstract Syntax Tree): this is a temporary data structure that represents the syntactic structure of the tree.
2. This AST can then be iterated over to convert it into the data structure of your choice.

Note

1. Quetzal does not know anything about your own classes, so you will need a bit of extra effort to write a recursion on the AST data
2. But at least you don't have to worry about the parsing logic anymore!
3. If the task becomes tedious, don't hesitate to ask for help !
4. On most compilers but Apple Clang, you could simplify the example code by writing: subtree.left_child = std::make_unique<MyNode>(ast.name, ast.distance_to_parent); in the recursion function.

Input

#include "quetzal/quetzal.hpp"
#include <iostream> // std::cout
#include <memory>   // std::unique_ptr
 
// Assume this binary tree class pervades your code
struct MyNode
{
    std::string id;                     // the vertex name
    double length;                      // the edge length
    std::unique_ptr<MyNode> left_child; // never use bare pointers in interfaces
    std::unique_ptr<MyNode> right_child;
};
 
// Define a namespace alias to shorten notation
namespace newick = quetzal::format::newick;
 
int main()
{
    std::string s = "(A:0.1,(C:0.3,D:0.4)B:0.5)E;";
 
    // We define a shorter type alias
    using ast_type = newick::ast::node;
 
    // We initialize the root of Abstract Syntax Tree (AST)
    ast_type ast;
 
    // We parse the input string into the AST
    auto ok = quetzal::parse(begin(s), end(s), newick::parser::tree, ast);
 
    // We check by printing the AST if the string was successfully parsed
    std::cout << quoted(s) << "\t " << (ok ? "Parsed" : "Failed");
    if (ok)
        std::cout << "\n\nAbstract Tree Syntax:\n\n" << ast << std::endl;
 
    // Let's use a recursive lambda until C++23 makes this much simpler
    static const auto populate = [](auto &tree_root, const auto &ast_root) {
        static auto recurse = [](auto &subtree, const auto &ast, auto &self) mutable -> void {
            // the idea is simple: traverse the AST ...
            if (ast.children.size() == 2)
            {
                // ... by instantiating the children and propagating to the next subtree
                subtree.left_child = std::unique_ptr<MyNode>(new MyNode{ast.name, ast.distance_to_parent});
                self(*subtree.left_child, ast.children.front(), self);
 
                subtree.right_child = std::unique_ptr<MyNode>(new MyNode{ast.name, ast.distance_to_parent});
                self(*subtree.right_child, ast.children.back(), self);
            }
        };
        recurse(tree_root, ast_root, recurse);
    };
 
    // Initialize your root using the AST data members
    MyNode my_root{.id = ast.name, .length = ast.distance_to_parent};
 
    // And then call the recursion to populate the children
    populate(my_root, ast);
 
    std::cout << (my_root.left_child == nullptr ? "Failed" : "Populated") << std::endl;
    std::cout << (my_root.right_child == nullptr ? "Failed" : "Populated") << std::endl;
 
    return 0;
}

Output

"(A:0.1,(C:0.3,D:0.4)B:0.5)E;"   Parsed
 
Abstract Tree Syntax:
 
"E"
├──0.1──"A"
└──0.5──"B"
    ├──0.3──"C"
    └──0.4──"D"
 
Populated
Populated