Operads and the Tree of Life

This week Lisa and I are visiting her 90-year-old mother in Montréal. Friday I’m giving a talk at the Université du Québec à Montréal. The main person I know there is André Joyal, an expert on category theory and algebraic topology. So, I decided to give a talk explaining how some ideas from these supposedly ‘pure’ branches of math show up in biology.

My talk is called ‘Operads and the Tree of Life’.

Trees

In biology, trees are very important:

So are trees of a more abstract sort: phylogenetic trees describe the history of evolution. The biggest phylogenetic tree is the ‘Tree of Life’. It includes all the organisms on our planet, alive now or anytime in the past. Here’s a rough sketch of this enormous tree:

Its structure is far from fully understood. So, biologists typically study smaller phylogenetic trees, like this tree of dog-like species made by Elaine Ostrander:

Abstracting still further, we can also think of a tree as a kind of purely mathematical structure, like this:

Trees are important in combinatorics, but also in algebraic topology. The reason is that in algebraic topology we get pushed into studying spaces equipped with enormous numbers of operations. We’d get hopelessly lost without a good way of drawing these operations. We can draw an operation with n inputs and one output as a little tree like this:

We can also draw the various ways of composing these operations. Composing them is just like building a big tree out of little trees!

An operation with n inputs and one output is called an n-ary operation. In the late 1960s, various mathematicians including Boardmann and Vogt realized that spaces with tons of n-ary operations were crucial to algebraic topology. To handle all these operations, Peter May invented the concept of an operad. This formalizes the way operations can be drawn as trees. By now operads are a standard tool, not just in topology, but also in algebraic geometry, string theory and many other subjects.

But how do operads show up in biology?

When attending a talk by Susan Holmes on phylogenetic trees, I noticed that her work on phylogenetic trees was closely related to a certain operad. And when I discussed her work here, James Griffin pointed out that this operad can be built using a slight variant of a famous construction due to Boardman and Vogt: their so-called ‘W construction’!

I liked the idea that trees and operads in topology might be related to phylogenetic trees. And thinking further, I found that the relation was real, and far from a coincidence. In fact, phylogenetic trees can be seen as operations in a certain operad… and this operad is closely related to the way computational biologists model DNA evolution as a branching sort of random walk.

That’s what I’d like to explain now.

I’ll be a bit sketchy, because I’d rather get across the basic ideas than the technicalities. I could even be wrong about some fine points, and I’d be glad to talk about those in the comments. But the overall picture is solid.

Phylogenetic trees

First, let’s ponder the mathematical structure of a phylogenetic tree. First, it’s a tree: a connected graph with no circuits. Second, it’s a rooted tree, meaning it has one vertex which is designated the root. And third, the leaves are labelled.

I should explain the third part! For any rooted tree, the vertices with just one edge coming out of them are called leaves. If the root is drawn at the bottom of the tree, the leaves are usually drawn at the top. In biology, the leaves are labelled by names of species: these labels matter. In mathematics, we can label the leaves by numbers where is the number of leaves.

Summarizing all this, we can say a phylogenetic tree should at least be a leaf-labelled rooted tree.

That’s not all there is to it. But first, a comment. When you see a phylogenetic tree drawn by a biologist, it’ll pretty much always a binary tree, meaning that as we move up any edge, away from the root, it either branches into two new edges or ends in a leaf. The reason is that while species often split into two as they evolve, it is less likely for a species to split into three or more new species all at once.

So, the phylogenetic trees we see in biology are usually leaf-labeled rooted binary trees. However, we often want to guess such a tree from some data. In this game, trees that aren’t binary become important too!

Why? Well, here another fact comes into play. In a phylogenetic tree, typically each edge can be labeled with a number saying how much evolution occurred along that edge. But as this number goes to zero, we get a tree that’s not binary anymore. So, we think of non-binary trees as conceptually useful ‘borderline cases’ between binary trees.

So, it’s good to think about phylogenetic trees that aren’t necessarily binary… and have edges labelled by numbers. Let’s make this into a formal definition:

Definition A phylogenetic tree is a leaf-labeled rooted tree where each edge not touching a leaf is labeled by a positive real number called its length.

By the way, I’m not claiming that biologists actually use this definition. I’ll write for the set of phylogenetic trees with leaves. This becomes a topological space in a fairly obvious way, where we can trace out a continuous path by continuously varying the edge lengths of a tree. But when some edge lengths approach zero, our graph converges to one where the vertices at ends of these edges ‘fuse into one’, leaving us with a graph with fewer vertices.

Here’s an example for you to check your understanding of what I just said. With the topology I’m talking about, there’s a continuous path in that looks like this:

These trees are upside-down, but don’t worry about that. You can imagine this path as a process where biologists slowly change their minds about a phylogenetic tree as new data dribbles in. As they change their minds, the tree changes shape in a continuous way.

For more on the space of phylogenetic trees, see:

• Louis Billera, Susan Holmes and Karen Vogtmann, Geometry of the space of phylogenetic trees, Advances in Applied Mathematics 27 (2001), 733-767.

Operads

How are phylogenetic trees related to operads? I have three things to say about this. First, they are the operations of an operad:

Theorem 1. There is an operad called the phylogenetic operad, or whose space of n-ary operations is

If you don’t know what an operad is, I’d better tell you now. They come in different flavors, and technically I’ll be using ‘symmetric topological operads’. But instead of giving the full definition, which you can find on the nLab, I think it’s better if I sketch some of the key points.

For starters, an operad consists of a topological space for each . The point in are called the n-ary operations of You can visualize an n-ary operation as a black box with input wires and one output wire:

Of course, this also looks like a tree.

We can permute the inputs of an n-ary operation and get a new n-ary operation, so we have an action of the permutation group on . You visualize this as permuting input wires:

More importantly, we can compose operations! If we have an n-ary operation , and n more operations, say , we can compose with all the rest and get an operation called

Here’s how you should imagine it:

Composition and permutation must obey some laws, all of which are completely plausible if you draw them as pictures. For example, the associative law makes a composite of composites like this well-defined:

Now, these pictures look a lot like trees. So it shouldn’t come as a shock that phylogenetic trees are the operations of some operad But let’s sketch why it’s true.

First, we can permute the ‘inputs’—meaning the labels on the leaves—of any phylogenetic tree and get a new phylogenetic tree. This is obvious.

Second, and more importantly, we can ‘compose’ phylogenetic trees. How do we do this? Simple: we glue the roots of a bunch of phylogenetic trees to the leaves of another and get a new one!

More precisely, suppose we have a phylogenetic tree with leaves, say . And suppose we have more, say Then we can glue the roots of to the leaves of to get a new phylogenetic tree called

Third and finally, all the operad laws hold. Since these laws all look obvious when you draw them using pictures, this is really easy to show.

If you’ve been paying careful attention, you should be worrying about something now. In operad theory, we think of an operation as having n inputs and one output. For example, this guy has 3 inputs and one output:

But in biology, we think of a phylogenetic tree as having one input and n outputs. We start with one species (or other grouping of organisms) at the bottom of the tree, let it evolve and branch, and wind up with n of them!

In other words, operad theorists read a tree from top to bottom, while biologists read it from bottom to top.

Luckily, this isn’t a serious problem. Mathematicians often use a formal trick where they take an operation with n inputs and one output and think of it as having one input and n outputs. They use the prefix ‘co-‘ to indicate this formal trick.

So, we could say that phylogenetic trees stand for ‘co-operations’ rather than operations. Soon this trick will come in handy. But not just yet!

The W construction

Boardman and Vogt had an important construction for getting new operads for old, called the ‘W construction’. Roughly speaking, if you start with an operad , this gives a new operad whose operations are leaf-labelled rooted trees where:

1) all vertices except leaves are labelled by operations of and a vertex with n input edges must be labelled by an n-ary operation of

and

2) all edges except those touching the leaves are labelled by numbers in .

If you think about it, the operations of are strikingly similar to phylogenetic trees, except that:

1) in phylogenetic trees the vertices don’t seem to be labelled by operations of a operad,

and

2) we use arbitrary nonnegative numbers to label edges, instead of numbers in .

The second point is a real difference, but it doesn’t matter much: if Boardman and Vogt had used nonnegative numbers instead of numbers in to label edges in the W construction, it would have worked just as well. Technically, they’d get a ‘weakly equivalent’ operad.

The first point is not a real difference. You see, there’s an operad called which has exactly one operation of each arity. So, labelling vertices by operations of is a completely trivial process.

As a result, we conclude:

Theorem 2. The phylogenetic operad is weakly equivalent to .

If you’re not an expert on operads (such a person is called an ‘operadchik’), you may be wondering what stands for. The point is that operads have ‘algebras’, where the abstract operations of the operad are realized as actual operations on some topological space. And the algebras of are precisely commutative topological monoids: that is, topological spaces equipped with a commutative associative product!

Branching Markov processes and evolution

By now, if you haven’t fallen asleep, you should be brimming with questions, such as:

1) What does it mean that phylogenetic trees are the operations of some operad ? Why should we care?

2) What does it mean to apply the W construction to the operad ? What’s the significance of doing this?

3) What does it mean that is weakly equivalent to ? You can see the definition of weak equivalence here, but it’s pretty technical, so it needs some explanation.

The answers to questions 2) and 3) take us quickly into fairly deep waters of category theory and algebraic topology—deep, that is, if you’ve never tried to navigate them. However, these waters are well-trawled by numerous experts, and I have little to say about questions 2) and 3) that they don’t already know. So given how long this talk already is, I’ll instead try to answer question 1). This is where some ideas from biology come into play.

I’ll summarize my answer in a theorem, and then explain what the theorem means:

Theorem 3. Given any continuous-time Markov process on a finite set , the vector space whose basis is naturally becomes a coalgebra of the phylogenetic operad.

Impressive, eh? But this theorem is really just saying that biologists are already secretly using the phylogenetic operad.

Biologists who try to infer phylogenetic trees from present-day genetic data often use simple models where the genotype of each species follows a ‘random walk’. Also, species branch in two at various times. These models are called Markov models.

The simplest Markov model for DNA evolution is the Jukes–Cantor model. Consider a genome of fixed length: that is, one or more pieces of DNA having a total of base pairs. For example, this tiny genome has base pairs, just enough to illustrate the 4 possible choices, which are called A, T, C and G:

Since there are 4 possible choices for each base pair, there are possible genotypes with base pairs. In the human genome, is about . So, there are about

genotypes of this length. That’s a lot!

As time passes, the Jukes–Cantor model says that the human genome randomly walks through this enormous set of possibilities, with each base pair having the same rate of randomly flipping to any other base pair.

Biologists have studied many ways to make this model more realistic in many ways, but in a Markov model of DNA evolution we’ll typically have some finite set of possible genotypes, together with some random walk on this set. But the term ‘random walk’ is a bit imprecise: what I really mean is a ‘continuous-time Markov process’. So let me define that.

Fix a finite set . For each time and pair of points i, j in , a continuous-time Markov process gives a number saying the probability that starting at the point i at time zero, the random walk will go to the point j at time . We can think of these numbers as forming an square matrix at each time . We demand that four properties hold:

1) depends continuously on .

2) For all s, t we have .

3) is the identity matrix.

4) For all j and t we have:

.

All these properties make a lot of sense if you think a bit, though condition 2) says that the random walk does not change character with the passage of time, which would be false given external events like, say, ice ages. As far as math jargon goes, conditions 1)-3) say that is a continuous one-parameter semigroup, while condition 4) together with the fact that says that at each time, is a stochastic matrix.

Let be the vector space whose basis is . To avoid getting confused, let’s write for the basis vector corresponding to . Any probability distribution on gives a vector in . Why? Because it gives a probability for each , and we can think of these as the components of a vector .

Similarly, for any time , we can think of the matrix as a linear operator

So, if we start with some probability distribution of genotypes, and let them evolve for a time according to our continuous-time Markov process, by the end the probability distribution will be .

But species do more than evolve this way: they also branch! A phylogenetic tree describes a way for species to evolve and branch.

So, you might hope that any phylogenetic tree gives a ‘co-operation’ that takes one probability distribution as input and returns n probability distributions as output.

That’s true. But these n probability distributions will be correlated, so it’s better to think of them as a single probability distribution on the set This can be seen as a vector in the vector space , the tensor product of n copies of

So, any phylogenetic tree gives a linear operator from to . We’ll call it

because we’ll build it starting from the Markov process

Here’s a sketch of how we build it—I’ll give a more precise account in the next and final section. A phylogenetic tree is made of a bunch of vertices and edges. So, I just need to give you an operator for each vertex and each edge, and you can compose them and tensor them to get the operator :

1) For each vertex with one edge coming in and n coming out:

we need an operator

that describes what happens when one species branches into n species. This operator takes the probability distribution we put in and makes n identical and perfectly correlated copies. To define this operator, we use the fact that the vector space has a basis labelled by the genotypes Here’s how the operator is defined:

2) For each edge of length , we need an operator that describes a random walk of length This operator is provided by our continuous-time Markov process: it’s

And that’s it! By combining these two kinds of operators, one for ‘branching’ and one for ‘random walking’, we get a systematic way to take any phylogenetic tree and get an operator

In fact, these operators obey just the right axioms to make into what’s called a ‘coalgebra’ of the phylogenetic operad. But to see this—that is, to prove Theorem 3—it helps to use a bit more operad technology.

The proof

I haven’t even defined coalgebras of operads yet. And I don’t think I’ll bother. Why not? Well, while the proof of Theorem 3 is fundamentally trivial, it’s sufficiently sophisticated that only operadchiks would enjoy it without a lengthy warmup. And you’re probably getting tired by now.

So, to most of you reading this: bye! It was nice seeing you! And I hope you sensed the real point of this talk:

Some of the beautiful structures used in algebraic topology are also lurking in biology. These structures may or may not be useful in biology… but we’ll never know if we don’t notice them and say what they are! So, it makes sense for mathematicians to spend some time looking for them.

Now, let me sketch a proof of Theorem 3. It follows from a more general theorem:

Theorem 4. Suppose is an object in some symmetric monoidal topological category . Suppose that is equipped with an action of the additive monoid . Suppose also that is a cocommutative coalgebra. Then naturally becomes a coalgebra of the phylogenetic operad.

How does this imply Theorem 3? In Theorem 3, is the category of finite-dimensional real vector space. The action of on is the continuous-time Markov process. And becomes a cocommutative coalgebra because it’s a vector space with a distinguished basis, namely the finite set . This makes into a cocommutative coalgebra in the usual way, where the comultiplication:

‘duplicates’ basis vectors:

while the counit:

‘deletes’ them:

These correspond to species splitting in two and species going extinct, respectively. (Biologists trying to infer phylogenetic trees often ignore extinction, but it’s mathematically and biologically natural to include it.) So, all the requirements are met to apply Theorem 4 and make into coalgebra of the phylogenetic operad.

But how do we prove Theorem 4? It follows immediately from Theorem 5:

Theorem 5. The phylogenetic operad is the coproduct of the operad and the additive monoid , viewed as an operad with only 1-ary operations.

Given how coproducts works, this means that an algebra of both and is automatically an algebra of . In other words, any commutative algebra with an action of is an algebra of . Dualizing, it follows that any cocommutative coalgebra with an action of is an coalgebra of And that’s Theorem 4!

But why is Theorem 5 true? First of all, I should emphasize that the idea of using it was suggested by Tom Leinster in our last blog conversation on the phylogenetic operad. And in fact, Tom proved a result very similar to Theorem 5 here:

• Tom Leinster, Coproducts of operads, and the W-construction, 14 September 2000.

He gives an explicit description of the coproduct of an operad and a monoid, viewed as an operad with only unary operations. He works with non-symmetric, non-topological operads, but his ideas also work for symmetric, topological ones. Applying his ideas to the coproduct of and , we see that we get the phylogenetic operad!

And so, phylogenetic trees turn out to be related to coproducts of operads. Who’d have thought it? But we really don’t have as many fundamentally different ideas as you might think: it’s hard to have new ideas. So if you see biologists and algebraic topologists both drawing pictures of trees, you should expect that they’re related.

Earth System Research for Global Sustainability

Some good news!

The International Mathematical Union or IMU is mainly famous for running the International Congress of Mathematicians every four years. But they do other things, too. The new vice-president of the IMU is Christiane Rousseau. Rousseau was already spearheading the Mathematics of Planet Earth 2013 project. Now she’s trying to get the IMU involved in a ten-year research initiative on sustainability.

As you can see from this editorial, she treats climate change and sustainability with the seriousness they deserve. Let’s hope more mathematicians join in!

I would like to get involved somehow, but I’m not exactly sure how.

Editorial

I had the privilege of being elected Vice-president of the IMU at the last General Assembly, and it is now five months that I am following the activities of the IMU. The subjects discussed at the Executive Committee are quite diverse, from the establishment of the permanent office to the ranking and pricing of journals, to mathematics in developing countries and the future ICM, and the members of the
Executive Committee tend to specialize on one or two dossiers. Although I am a pure mathematician myself, I am becoming more and more interested in the science of sustainability, so let me talk to you of this.

IMU is one of the international unions inside the International Council of Science (ICSU). At the Executive we regularly receive messages from ICSU asking for input from its members. While it is not new that scientists are involved in the study of climate change and sustainability issues, a new feeling of emergency has developed. The warning signs are becoming more numerous that urgent action is needed if we want to save the planet from a disastrous future, since we may not be far from a point of no return: climate change with more extreme weather events, rising of the sea level with the melting of glaciers, shortage of food and water in the near future because of the increase of the world population and the climate change, loss of biodiversity, new epidemics or invasive species, etc. This explains why ICSU is starting a new 10-year research initiative: EARTH SYSTEM RESEARCH FOR GLOBAL SUSTAINABILITY, and a Steering Committee for this initiative is presently nominated. The goals of the Initiative are to:

1. Deliver at global and regional scales the knowledge that societies need to effectively respond to global change while meeting economic and social goals;

2. Coordinate and focus international scientific research to address the Grand Challenges and Belmont Challenge;

3. Engage a new generation of researchers in the social, economic, natural, health, and engineering sciences in global sustainability research.

In the same spirit, ICSU is preparing a strong scientific presence at the next United Nations Conference on Sustainable Development (Rio+20) that will take place on June 4-6, 2012 in Rio de Janeiro. For this, ICSU is organizing a number of preparatory regional and global meetings. It is clear that mathematical sciences have an essential role in the interdisciplinary research that needs to take place in order to achieve significant impact. The other scientific disciplines concerned are numerous from physics, to biology, to economics, etc.

Let me quote Graciela Chichilnisky, the author of the carbon market of the UN Kyoto Protocol: “It is the physicists that study the climate change, but it is the economists who advise the politicians that take the decisions.” Considering the importance of the contribution of mathematical sciences in sustainability issues, IMU has asked to participate actively in these preparatory meetings and be represented at Rio+20. This should be an occasion to build partnerships with the other scientific unions inside ICSU. More and more mathematicians and research institutes around the world become interested in sustainable development as is acknowledged by the large participation in Mathematics of Planet Earth 2013 which was recently endorsed by IMU. But the world needs more than a one year initiative. The science of sustainability is full of challenging problems which are very interesting mathematically. Many of these problems require new mathematical techniques. We could hope that these initiatives will allow training a new generation of researchers in mathematical sciences who will be able to work in interdisciplinary teams to address these issues.

Christiane Rousseau
Vice-President of Executive Committee of IMU

Conferences on Math and Climate Change

Here are some conferences on climate change and related issues, specially designed to get mathematicians interacting with scientists who work on these things! If you know of any more coming up, please let me know. These ones are sponsored by the Mathematics and Climate Research Network, a US-based organization, but there are probably others.

• Society for Industrial and Applied Mathematics (SIAM) Conference on Applications of Dynamical Systems, Snowbird, Utah, US, 22-26 May, 2011. Organized by Jonathan Dawes and Vivien Kirk.

Climate modeling and data assimilation are among the themes of this conference, which is aimed at starting communication between mathematicians who develop dynamical systems techniques and the applied scientists who use them.

• Mathematical Biosciences Institute (MBI) Workshop on Ocean Ecologies and their Physical Habitats in a Changing Climate, Columbus, Ohio, US, June 2011. Organized by Ken Golden, Chris Jones, Hans Kaper, and Mary Lou Zeeman.

The goal of this workshop is to bring together biologists studying ocean and polar ecologies; oceanographers, biogeochemists, and climate scientists studying the changing physical habitats; and mathematicians with ecological and physical expertise. The interactions between ocean ecological systems and their physical environments may dramatically impact both marine biodiversity and the planetary response to the changing atmosphere. The types of mathematics used to model ecological and physical processes are typically quite different. The team organizing this workshop anticipates interesting new mathematical challenges arising from combining these different approaches. The workshop will focus on two main themes:

1) polar and sea ice ecologies;

2) phytoplankton and the carbon cycle.

• Minisymposium on the Dynamics of the Earth’s Climate, as part of the International Congress on Industrial and Applied Mathematics (ICIAM), Vancouver, British Columbia, July 2011.
Organized by Hans G. Kaper, Mary C. Silber and Mary Lou Zeeman.

The speakers in this mini-symposium will highlight some interesting mathematical problems that have come from climate science and can be addressed with techniques developed in the dynamical systems community.

• Institute of Mathematics and its Applications (IMA) Conference on Mathematics of the Climate System, University of Reading, United Kingdom, 12-15 September, 2011. Organized by Paul Williams, Colin Cotter, Mike Cullen, Mike Davey, Christopher Ferro, John Huthnance and David Stainforth.

This conference is about the construction and use of mathematical models of the climate system. The conference will focus on three related topics:

1) the extraction of mathematical models from climate data and climate-model output (homogenisation, stochastic model reduction, bistability and metastable states, low frequency variability, data-driven coarse-graining, set-oriented methods, trend identification, time-series analysis);

2) reduced models and their dynamics (linear response theory, bifurcations, extreme events, uncertainty);

3) testing hypotheses about the climate system using statistical frameworks (emulators, Bayesian methods, nonparametric methods, equitability).

Networks and Population Biology (Part 3)

I think today I’ll focus on one aspect of the talks Susan Holmes gave today: the space of phylogenetic trees. Her talks were full of interesting applications to genetics, but I’m afraid my summary will drift off into a mathematical daydream inspired by what she said! Luckily you can see her actual talk uncontaminated by my reveries here:

• Susan Holmes, Treespace: distances between trees.

It’s based on this paper:

• Louis Billera, Susan Holmes and Karen Vogtmann, Geometry of the space of phylogenetic trees, Advances in Applied Mathematics 27 (2001), 733-767.

As I mentioned last time, a phylogenetic tree is something like this:

More mathematically, what we see here is a tree (a connected graph with no circuits), with a distinguished vertex called the root, and vertices of degree 1, called leaves, that are labeled with elements from some -element set. We shall call such a thing a leaf-labelled rooted tree.

Now, the tree above is actually a binary tree, meaning that as we move up an edge, away from the root, it either branches into two new edges or ends in a leaf. (More precisely: each vertex that doesn’t have degree 1 has degree 3.) This makes sense in biology because while species often split into two as they evolve, it is less likely for a species to split into three all at once.

So, the phylogenetic trees we see in biology are usually leaf-labeled rooted binary trees. However, we often want to guess such a tree from some data. In this game, trees that aren’t binary become important too!

Why? Well, each edge of the tree can be labeled with a number saying how much evolution occurred along that edge: for example, how many DNA base pairs changed. But as this number goes to zero, we get a tree that’s not binary anymore. So, we think of non-binary trees as conceptually useful ‘intermediate cases’ between binary trees.

This idea immediately leads us to consider a topological space consisting of phylogenetic trees which are not necessarily binary. And at this point in the lecture I drifted off into a daydream about ‘operads’, which are a nice piece of mathematics that’s closely connected to this idea.

So, I will deviate slightly from Holmes and define a phylogenetic tree to be a leaf-labeled rooted tree where each edge is labeled by a number called its length. This length must be positive for every edge except the edge incident to the root; for that edge any nonnegative length is allowed.

Let’s write for the set of phylogenetic trees with leaves. This becomes a topological space in a fairly obvious way. For example, there’s a continuous path in that looks like this:

Moreover we have this fact:

Theorem. There is a topological operad called the phylogenetic operad, or whose space of n-ary operations is for and the empty set for .

If you don’t know what an operad is, don’t be scared. This mainly just means that you can glue a bunch of phylogenetic trees to the top of another one and get a new phylogenetic tree! More precisely, suppose you have a phylogenetic tree with leaves, say . And suppose you have more, say . Then you can glue the roots of to the leaves of to get a new phylogenetic tree called . Furthermore, this gluing operation obeys some rules which look incredibly intimidating when you write them out using symbols, but pathetically obvious when you draw them using pictures of trees. And these rules are the definition of an operad.

I would like to know if mathematicians have studied the operad It’s closely related to Stasheff’s associahedron operad, but importantly different. Operads have ‘algebras’, and the algebras of the associahedron operad are topological spaces with a product that’s ‘associative up to coherent homotopy’. I believe algebras of the phylogenetic operad are topological spaces with a commutative product that’s associative up to coherent homotopy. Has someone studied these?

In their paper Holmes and her coauthors discuss the associahedron in relation to their own work, but they don’t mention operads. I’ve found another paper that mentions ‘the space of phylogenetic trees’:

• David Speyer and Bernd Sturmfels, The tropical Grassmannian, Adv. Geom. 4 (2004), 389–411.

but they don’t seem to study the operad aspect either.

Perhaps one reason is that Holmes and her coauthors deliberately decide to ignore the labellings on the edges incident to the leaves. So, they get a space of phylogenetic trees with leaves whose product with is the space I’m calling . As they mention, this simplifies the geometry a bit. However, it’s not so nice if you want an operad that accurately describes how you can build a big phylogenetic tree from smaller ones.

They don’t care about operads; they do some wonderful things with the geometry of their space of phylogenetic trees. They construct a natural metric on it, and show it’s a CAT(0) space in the sense of Gromov. This means that the triangles in this space are more skinny than those in Euclidean space—more like triangles in hyperbolic space:

They study geodesics in this space—even though it’s not a manifold, but something more singular. And so on!

There’s a lot of great geometry here. But for Holmes, all this is just preparation for doing some genomics— for example, designing statistical tests to measure how reliable the phylogenetic trees guessed from data actually are. And for that aspect, try this:

• Susan Holmes, Statistical approach to tests involving phylogenies, in O. Gascuel, editor, Mathematics of Evolution and Phylogeny, Oxford U. Press, Oxford, 2007.

Networks and Population Biology (Part 2)

Yesterday I was too sick to attend this conference:

• Tutorials on discrete mathematics and probability in networks and population biology, Institute of Mathematical Sciences, National University of Singapore, 2-6 May 2011. Organized by Andrew Barbour, Malwina Luczak, Gesine Reinert and Rongfeng Sun.

But today I bounced back, and I want to tell you about the first lectures by Persi Diaconis and his wife Susan Holmes, both statisticians at Stanford University.

Since Persi Diaconis is one of those people who make it really obvious that being a mathematician is more fun than any other profession, I should say a bit about him. He left home at the age of 14 to become a professional stage magician, and only returned to school so that he could learn all the math needed to understand Feller’s famous book An Introduction to Probability Theory and Its Applications. Since then he has worked on the mathematics of shuffling cards and the physics of flipping coins. He also works on random matrix theory, the statistics of Riemann zeta zeroes, applications of Hopf algebras to probability theory, and all sorts of other things that don’t sound fun unless you know enough to realize that yes, they are.

In my last blog post on this conference Tim van Beek asked if Persi could really flip a coin and have it land with whatever side up he wanted. So, I asked him about that. He said yes, he could. He didn’t demonstrate it at the time, since we were just about to go to a talk. But he said you can see videos of other people doing it on the web.

Unfortunately, the videos I’ve found so far mainly convince me, not that there are people who have mastered the physical skill of controlling a coin flip, but that there are lots of tricksters out there! Here are two:

• Derrin Brown, The system – coin toss.

• EricSurf6, How to win a coin toss.

Can you figure out how they do it? If you don’t get the second one, you can watch this.

So, I have to wonder: can Persi really flip a coin and have it land the way he wants, or is just fooling people into thinking he can? Mind you, I wouldn’t consider the latter a bad thing, at least insofar as he’s still a practicing stage magician: a basic part of the craft is trickery of this sort. I should ask if he’s sworn the Magician’s Oath:

“As a magician I promise never to reveal the secret of any illusion to a non-magician, unless that one swears to uphold the Magician’s Oath in turn. I promise never to perform any illusion for any non-magician without first practicing the effect until I can perform it well enough to maintain the illusion of magic.”

Anyway, he is giving a series of talks on “Graph limits and exponential models”. There are lots of very large graphs in biology, so he plans to talk about ways to study large graphs by taking limits where the graphs become infinite. He began by describing a result so beautiful that I would like to state it here.

‘Graph’ means many things, but in his talk a graph is simple (no self-loops or multiple edges between vertices), undirected and labeled. In other words, something like this:

A graph homomorphism

is a map sending vertices of to vertices of , such that whenever vertices of have an edge between them, the vertices have an edge between them in .

Given this, we can define to be the fraction of the functions from vertices of to vertices of that are actually graph homomorphisms.

The famous graph theorist Lovasz said that a sequence of graphs converges if converges for all graphs .

Using this definition, we have:

Theorem (Aldous-Hoover, Lovasz et al). There is a compact metric space such that each graph gives a point in , and a sequence of graphs converges iff

for some point .

The space is called the space of graph limits, and the important thing is that one can construct it rather explicitly! Persi explained how. For a written explanation, see:

• Miklos Racz, Dense graph limits, 31 October 2010.

After a break, Susan Holmes began her series of talks on “Phylogeny, trees and the human microbiome”. I can’t easily summarize all she said, but here are the slides for the first talk:

• Susan Holmes, Molecular evolution: phylogenetic tree building.

Her basic plan today was to teach us a little genomics, a little statistics, and a little bit about Markov chains, so that we could start to understand how people attempt to reconstruct phylogenetic trees from DNA data.

A phylogenetic tree is something like this:

or something less grand where, for example, you only consider species of frogs, or even HIV viruses in a single patient. These days there are programs to generate such trees given copies of the ‘same gene’ in the different organisms you’re considering. These programs have names like PhyML, FastML, RaxML and Mr. Bayes. In case you’re wondering, ‘ML’ stands for ‘maximum likelihood’ a standard technique in statistics, while Bayes was the originator of some very important ideas in statistical reasoning.

But even though there are programs to do these things, there are still lots of fascinating problems to solve in this area! And indeed, even understanding the programs is a bit of a challenge.

Since we were talking about the genetic code here recently, I’ll wrap up by mentioning one thing learned about this from Susan’s talk.

It’s common to describe the accumulation of changes in DNA using substitution models. These are continuous time Markov chains where each base pair has a given probability of mutating to another one. Often people assume this probability for each base pair is independent of what its neighbors do. The last assumption is known to be false, but that’s another story for another day! What I wanted to say is that there are two famous models. The simplest is the Jukes-Cantor model, where each of the four bases—A, T, C, and G—has an equal probability of mutating into any other. But a somewhat better model is the Kimura model, where the transitions

A ↔ T
C ↔ G

have a different rate of occuring than all the remaining possibilities. If you look at the pictures of the A, T, C, and G molecules here you’ll instantly see why:

Since I’m busy learning about continuous-time Markov chains, it’s nice to see more good examples!

Networks and Population Biology (Part 1)

There are some tutorials starting today here:

• Tutorials on discrete mathematics and probability in networks and population biology, Institute of Mathematical Sciences, National University of Singapore, 2-6 May 2011. Organized by Andrew Barbour, Malwina Luczak, Gesine Reinert and Rongfeng Sun.

Rick Durrett is speaking on “Cancer modelling”. For his slides, see here, here and here. But here’s a quick taste:

Back in 1954, Armitage and Doll noticed that log-log plots of cancer incidence as a function of age are close to linear, except for breast cancer, which slows down in older women. They suggested an explanation: a chain of independent random events have to occur before cancer can start. A simple model based on a Markov process gives a simple formula for how many events it must take—see the first batch of slides for details. This work was the first of a series of ever more sophisticated multi-stage models of carcinogenesis.

One of the first models Durrett explained was the Moran process: a stochastic model of a finite population of constant size in which things of two types, say and are competing for dominance. I believe this model can be described by a stochastic Petri net with two states, and , and two transitions:

and

Since I like stochastic Petri nets, I’d love to add this to my collection.

Chris Cannings is talking about “Evolutionary conflict theory” and the concept of ‘evolutionary stable strategies’ for two-party games. Here’s the basic idea, in a nutshell.

Suppose a population of animals roams around randomly and whenever two meet, they engage in some sort of conflict… or more generally, any sort of ‘game’. Suppose each can choose from some set of strategies. Suppose that if one chooses strategy and the other chooses strategy , the expected ‘payoff’ to the one is , while for the other it’s .

More generally, the animals might choose their strategies probabilistically. If the first chooses the ith strategy with probability and the second chooses it with probability then the expected payoff to the first player is

where the angle brackets are the usual inner product in I’m saying this in an overly fancy way, and making it look like quantum mechanics, in the hope that some bright kid out there will get some new ideas. But it’s not rocket science; the angle bracket is just a notation for this sum:

Let me tell you what it means for a probabilistic strategy to be ‘evolutionarily stable’. Suppose we have a ‘resident’ population of animals with strategy and we add a few ‘invaders’ with some other strategy, say . Say the fraction of animals who are invaders is some small number , while the fraction of residents is .

If a resident plays the game against a randomly chosen animal, its expected payoff will be

Indeed, it’s just as if the resident was playing the game against an animal with probabilistic strategy ! On the other hand, if an invader plays the game against a randomly chosen animal, its expected payoff will be

The strategy is evolutionarily stable if the residents do better:

for all probability distributions and sufficiently small .

Canning showed us how to manipulate this condition in various ways and prove lots of nice theorems. His slides will appear online later, and then I’ll include a link to them. Naturally, I’m hoping we’ll see that a dynamical model, where animals with greater payoff get to reproduce more, has the evolutionary stable strategies as stable equilibria. And I’m hoping that some model of this sort can be described using a stochastic Petri net—though I’m not sure I see how.

On another note, I was happy to see Persi Diaconis and meet his wife Susan Holmes. Both will be speaking later in the week. Holmes is a statistician who specializes in “large, messy datasets” from biology. Lately she’s been studying ant networks! Using sophisticated image analysis to track individual ants over long periods of time, she and her coauthors have built up networks showing who meets who in ant ant colony. They’ve found, for example, that some harvester ants interact with many more of their fellows than the average ant. However, this seems to be due to their location rather than any innate proclivity. They’re the ants who hang out near the entrance of the nest!

That’s my impression from a short conversation, anyway. I should read her brand-new paper:

• Noa Pinter-Wollman, Roy Wollman, Adam Guetz, Susan Holmes and Deborah M. Gordon, The effect of individual variation on the structure and function of interaction networks in harvester ants, Journal of the Royal Society Interface, 13 April 2011.

She said this is a good book to read:

• Deborah M. Gordon, Ant Encounters: Interaction Networks and Colony Behavior, Princeton U. Press, Princeton New Jersey, 2010.

There are also lots of papers available at Gordon’s website.

Network Theory (Part 7)

guest post by Jacob Biamonte

This post is part of a series on what John and I like to call Petri net field theory. Stochastic Petri nets can be used to model everything from vending machines to chemical reactions. Chemists have proven some powerful theorems about when these systems have equilibrium states. We’re trying to bind these old ideas into our fancy framework, in hopes that quantum field theory techniques could also be useful in this deep subject. We’ll describe the general theory later; today we’ll do an example from population biology.

Those of you following this series should know that I’m the calculation bunny for this project, with John playing the role of the wolf. If I don’t work quickly, drawing diagrams and trying to keep up with John’s space-bending quasar of information, I’ll be eaten alive! It’s no joke, so please try to respond and pretend to enjoy anything you read here. This will keep me alive for longer. If I did not take notes during our meetings, lots of this stuff would have never made it here, so hope you enjoy.

Amoeba reproduction and competition

Here’s a stochastic Petri net:

It shows a world with one state, amoeba, and two transitions:

• reproduction, where one amoeba turns into two. Let’s call the rate constant for this transition .

• competition, where two amoebas battle for resources and only one survives. Let’s call the rate constant for this transition .

We are going to analyse this example in several ways. First we’ll study the deterministic dynamics it describes: we’ll look at its rate equation, which turns out to be the logistic equation, familiar in population biology. Then we’ll study the stochastic dynamics, meaning its master equation. That’s where the ideas from quantum field theory come in.

The rate equation

If is the population of amoebas at time , we can follow the rules explained in Part 3 and crank out this rate equation:

We can rewrite this as

where

What’s the meaning of and ?

• is the carrying capacity, that is, the maximum sustainable population the environment can support.

• is the growth rate describing the approximately exponential growth of population when is small.

It’s a rare treat to find such an important differential equation that can be solved by analytical methods. Let’s enjoy solving it.

We start by separating variables and integrating both sides:

We need to use partial fractions on the left side above, resulting in

and so we pick up a constant , let

and rearrange things as

so the population as a function of time becomes

At we can determine uniquely. We write and becomes

The model now becomes very intuitive. Let’s set and make a plot for various values of :

We arrive at three distinct cases:

• equilibrium (). The horizontal blue line corresponds to the case where the initial population exactly equals the carrying capacity. In this case the population is constant.

• dieoff (). The three decaying curves above the horizontal blue line correspond to cases where initial population is higher than the carrying capacity. The population dies off over time and approaches the carrying capacity.

• growth (). The four increasing curves below the horizontal blue line represent cases where the initial population is lower than the carrying capacity. Now the population grows over time and approaches the carrying capacity.

The master equation

Next, let us follow the rules explained in Part 6 to write down the master equation for our example. Remember, now we write:

where is the probability of having amoebas at time , and is a formal variable. The master equation says

where is an operator on formal power series called the Hamiltonian. To get the Hamiltonian we take each transition in our Petri net and build an operator built from creation and annihilation operators, as follows. Reproduction works like this:

while competition works like this:

Here is the annihilation operator, is the creation operator and is the number operator. Last time John explained precisely how the ‘s arise. So the theory is already in place, and we arrive at this Hamiltonian:

Remember, is the rate constant for reproduction, while is the rate constant for competition.

The master equation can be solved: it’s equivalent to

so that is constant, and so

and that’s it! We can calculate the time evolution starting from any initial probability distribution of populations. Maybe everyone is already used to this, but I find it rather remarkable.

Here’s how it works. We pick a population, say amoebas at This would mean . We then evolve this state using . We expand this operator as

This operator contains the full information for the evolution of the system. It contains the histories of all possible amoeba populations—an amoeba mosaic if you will. From this, we can construct amoeba Feynman diagrams.

To do this, we work out each of the terms in the expansion above. The first-order terms correspond to the Hamiltonian acting once. These are proportional to either or . The second-order terms correspond to the Hamiltonian acting twice. These are proportional to either , or . And so on.

This is where things start to get interesting! To illustrate how it works, we will consider two possibilities for the second-order terms:

1) We start with a lone amoeba, so . It reproduces and splits into two. In the battle of the century, the resulting amoebas compete and one dies. At the end we have:

We can draw this as a Feynman diagram:

You might find this tale grim, and you may not like the odds either. It’s true, the odds could be better, but people are worse off than amoebas! The great Japanese swordsman Miyamoto Musashi quoted the survival odds of fair sword duels as 1/3, seeing that 1/3 of the time both participants die. A remedy is to cheat, but these amoeba are competing honestly.

2) We start with two amoebas, so the initial state is . One of these amoebas splits into two. One of these then gets into an argument with the original amoeba over the Azimuth blog. The amoeba who solved all John’s puzzles survives. At the end we have

with corresponding Feynman diagram:

This should give an idea of how this all works. The exponential of the Hamiltonian gives all possible histories, and each of these can be translated into a Feynman diagram. In a future blog entry, we might explain this theory in detail.

An equilibrium state

We’ve seen the equilibrium solution for the rate equation; now let’s look for equilibrium solutions of the master equation. This paper:

• D. F. Anderson, G. Craciun and T.G. Kurtz, Product-form stationary distributions for deficiency zero chemical reaction networks, arXiv:0803.3042.

proves that for a large class of stochastic Petri nets, there exists an equilibrium solution of the master equation where the number of things in each state is distributed according to a Poisson distribution. Even more remarkably, these probability distributions are independent, so knowing how many things are in one state tells you nothing about how many are in another!

Here’s a nice quote from this paper:

The surprising aspect of the deficiency zero theorem is that the assumptions of the theorem are completely related to the network of the system whereas the conclusions of the theorem are related to the dynamical properties of the system.

The ‘deficiency zero theorem’ is a result of Feinberg, which says that for a large class of stochastic Petri nets, the rate equation has a unique equilibrium solution. Anderson showed how to use this fact to get equilibrium solutions of the master equation!

We will consider this in future posts. For now, we need to talk a bit about ‘coherent states’.

These are all over the place in quantum theory. Legend (or at least Wikipedia) has it that Erwin Schrödinger himself discovered coherent states when he was looking for states of a quantum system that look ‘as classical as possible’. Suppose you have a quantum harmonic oscillator. Then the uncertainty principle says that

where is the uncertainty in the momentum and is the uncertainty in position. Suppose we want to make as small as possible, and suppose we also want . Then we need our particle to be in a ‘coherent state’. That’s the definition. For the quantum harmonic oscillator, there’s a way to write quantum states as formal power series

where is the amplitude for having quanta of energy. A coherent state then looks like this:

where can be any complex number. Here we have omitted a constant factor necessary to normalize the state.

We can also use coherent states in classical stochastic systems like collections of amoebas! Now the coefficient of tells us the probability of having amoebas, so had better be real. And probabilities should sum to 1, so we really should normalize as follows:

Now, the probability distribution

is called a Poisson distribution. So, for starters you can think of a ‘coherent state’ as an over-educated way of talking about a Poisson distribution.

Let’s work out the expected number of amoebas in this Poisson distribution. In the answers to the puzzles in Part 6, we started using this abbreviation:

We also saw that the expected number of amoebas in the probability distribution is

What does this equal? Remember that . The annihilation operator is just , so

and we get

But we saw in Part 5 that is stochastic, meaning

for any Furthermore, our here has

since it’s a probability distribution. So:

The expected number of amoebas is just .

Puzzle 1. This calculation must be wrong if is negative: there can’t be a negative number of amoebas. What goes wrong then?

Puzzle 2. Use the same tricks to calculate the standard deviation of the number of amoebas in the Poisson distribution .

Now let’s return to our problem and consider the initial amoeba state

Here aren’t bothering to normalize it, because we’re going to look for equilibrium solutions to the master equation, meaning solutions where doesn’t change with time. So, we want to solve

Since this equation is linear, the normalization of doesn’t really matter.

Remember,

Let’s work this out. First consider the two terms:

and

Likewise for the terms we find

and

Here I’m using something John showed in Part 6: the product equals the ‘falling power’ .

The sum of all four terms must vanish. This happens whenever

which is satisfied for

Yipee! We’ve found an equilibrium solution, since we found a value for that makes . Even better, we’ve seen that the expected number of amoebas in this equilibrium state is

This is just the same as the equilibrium population we saw for the rate equation—that is, the carrying capacity for the logistic equation! That’s pretty cool, but it’s no coincidence: in fact, Anderson proved it works like this for lots of stochastic Petri nets.

I’m not sure what’s up next or what’s in store, since I’m blogging at gun point from inside a rabbit cage:

I’d imagine we’re going to work out the theory behind this example and prove the existence of equilibrium solutions for master equations more generally. One idea John had was to have me start a night shift—that way you’ll get Azimuth posts 24 hours a day.