Skip to content

Traits

The trait level evolves phenotypes: a body size, a habitat, the presence or absence of a structure. A trait evolves along the species tree like everything else here: the complete tree, extinct lineages included. Hand it the species run itself or its .complete_tree; every level takes either and runs on the complete tree. There are two kinds, continuous and discrete.

The two kinds, on the same tree. A **continuous** trait starts at one value and drifts along every branch, so every node ends at a number of its own, with each branch painted with the value at its far end, and close relatives end up close. A **discrete** trait sits in one state and switches to another now and then; the dots mark the two switches, and every branch below a switch is in the new state until it switches again. The chips at the tips are the state each one ended in. The two are different functions, simulate_continuous and simulate_discrete, because they answer different questions.

Continuous traits

A continuous trait diffuses. Give it a starting value and a rate, and it wanders down every branch, its variance growing in proportion to elapsed time. On its own that is Brownian motion, the null model of continuous trait evolution:

from zombi2 import species, traits
from zombi2.params import PerLineage

tree = species.simulate_species_tree(birth=1.0, death=0.3, n_extant=20, seed=1)

# BM: a body size diffusing from 0 at variance-rate σ² = 1.0
size = traits.simulate_continuous(tree, start=0.0, rate=1.0, seed=1)

size.values                # each extant tip's value, keyed by tip name: {'n19': -0.39, ...}
size.node_values           # every node of the complete tree, ancestors and extinct included
size.write("out/traits/")  # trait_values.tsv · trait_tree.nwk · trait_summary.json (Appendix B)

A discrete run gives the same shape back, states instead of numbers, and writes trait_events.tsv besides: the switch log, which is the file a conditioned run takes as its driver (Chapter 8).

Here rate is the Brownian variance-rate σ², the trait level's reading of "how fast". It is counted per lineage, the only scope here, so PerLineage(1.0) is the bare 1.0 in wrapper form, and like every rate in ZOMBI2 it takes the verbs. Two variations, the first by two more arguments, the second by a verb on the rate:

# OU: the same diffusion, pulled toward an optimum value
traits.simulate_continuous(tree, start=0.0, rate=1.0,
                           reverts_to=2.0, pull=0.5, seed=1)

# early burst: the diffusion rate itself decays through time
traits.simulate_continuous(tree, start=0.0,
                           rate=PerLineage(1.0).changing_at({0: 1.0, 5: 0.2}), seed=1)

The Ornstein–Uhlenbeck process is Brownian motion with a rubber band: reverts_to is the optimum it is pulled back toward, and pull is how hard. Early burst (or ACDC) is a diffusion rate that decays as the tree ages, so most of the divergence happens near the root (Tr3); it is written with the same changing_at that gives the species tree its skyline, a schedule of steps, so a finer schedule is a finer approximation of the decay.

Two more arguments sit alongside rate. regimes= paints a multi-optimum OU, where clades pull toward different optima (Tr6): a discrete run supplies the painting, and reverts_to becomes a dict keyed by its states, as in regimes=habitat, reverts_to={"cave": -1.0, "surface": 3.0}. And at_speciation= adds a jump at each split rather than along the branches, so change concentrates at branching (Tr7). The value is the jump's variance: at_speciation=0.5 draws each jump from a normal with variance 0.5, standard deviation about 0.7. None of these is a separate model with its own function and its own parameters, which is why they combine: a trait that bursts early and reverts to an optimum is one rate with one verb and two arguments.

The optimum can also follow another continuous trait, grown first on the same tree. reverts_to=set_by(group_size, lambda v: 1.0 + 0.8 * v, step=0.01) makes a larger group keep pulling the trait upward. The function turns the other trait's value into the optimum, and the optimum can be any number, negative included. ZOMBI2 knows the other trait only at the nodes, so along a branch it draws that trait's path between them and reads it in stretches of at most step. Halving step halves the error and also changes the numbers for a given seed, so convergence is checked across seeds rather than by rerunning one; Chapter 8 gives the recipe. Write path=False to read the straight line between the two node values instead, which is the path's average rather than a path. The effect runs one way: the other trait steers this one, and this one does not steer it back.

regimes= is the one argument that asks you to give things up, and it refuses loudly rather than ignoring what you passed: with it the σ² is a plain number (not a modified rate), the jump variance is one number shared across regimes, and the run is one trait (no correlation=).

Discrete traits

A discrete trait takes a finite set of states and switches between them along the branches: a continuous-time Markov chain, the field's Mk model (Lewis 2001).

# Mk: habitat flips between two states at rate 0.1
traits.simulate_discrete(tree, states=["marine", "terrestrial"],
                         switch=0.1, start="marine", seed=1)

Leave start out and the root state is drawn uniformly from states. When the flips are not symmetric, replace the single rate with directed rates, one per ordered pair: the dict below, or a full k × k matrix whose entry at row i, column j is the rate i → j (Tr9):

# asymmetric: gains are commoner than losses
traits.simulate_discrete(tree, states=["absent", "present"],
                         switch={"absent->present": 0.2, "present->absent": 0.05},
                         seed=1)

Two more discrete models live on this same function. The threshold model reads the states off a continuous liability: simulate_discrete(tree, states=["absent", "present"], liability=1.0, threshold=0.0) diffuses a Brownian liability from start (a number here, 0.0 by default) and cuts it at the thresholds, k − 1 increasing cuts for k states (Tr10; --liability and --threshold on the command line). Give liability as a dict and add a correlation={("wings", "eyes"): 0.6} overlay and two discrete traits evolve together, their liabilities diffusing jointly and each cut by the shared thresholds, which is the correlated-binary model of the Literature table (Tr13).

Correlated traits

Two traits that evolve independently are two separate calls, in either order. Two traits that drift together cannot be simulated one before the other, because each is entangled with the other as it unfolds (Tr11). Correlation is specified as per-trait rates plus a correlation overlay:

traits.simulate_continuous(tree,
    start={"size": 0.0, "limb": 0.0},
    rate={"size": 1.0, "limb": 0.8},        # one variance-rate per trait
    correlation={("size", "limb"): 0.6},    # the overlay, ∈ [−1, 1]
    seed=1)

The overlay is a dimensionless number in [−1, 1], not a covariance matrix. Under correlation= the per-trait rates are plain numbers.

The overlay carries the reversion and the speciation jumps as well. Add reverts_to and pull and each trait reverts to its own optimum at its own strength; add at_speciation and the jump at each split is drawn under the same correlation the diffusion uses. Each argument takes one value shared across the traits, or a dict giving one per trait:

traits.simulate_continuous(tree,
    start={"size": 0.0, "limb": 0.0},
    rate={"size": 1.0, "limb": 0.8},
    correlation={("size", "limb"): 0.6},
    reverts_to={"size": 3.0, "limb": -1.0},   # each trait its own optimum
    pull={"size": 1.5, "limb": 0.4},          # and its own strength
    at_speciation=0.5, seed=1)

The traits can also steer each other. Give pull as a full matrix, keyed by pairs of traits, and one trait's distance from its optimum moves the other:

traits.simulate_continuous(tree,
    start={"size": 0.0, "limb": 0.0},
    rate={"size": 1.0, "limb": 0.8},
    reverts_to={"size": 3.0, "limb": -1.0},
    pull={("size", "size"): 1.5, ("limb", "limb"): 0.4,
          ("limb", "size"): -0.6},            # size above its optimum pushes limb up
    seed=1)

The entry at ("limb", "size") is the effect of size on limb. A negative entry pushes limb up while size sits above its optimum, a positive one pushes it down, and an entry left out is zero. Each trait's own entry is its pull toward its own optimum, and it must be positive. The effect can run both ways: add ("size", "limb") and limb steers size too. The run is exact, with no step to tune. The matrix is minus the selection matrix of coevolve (Ringen et al. 2026), which writes the drift with the opposite sign.

Feedback with no exact solution

The drift matrix above is exact because each optimum is a straight line in the other trait. Bend that line, or make one of the two traits discrete, and the exact solution goes. simulate_traits runs those two models instead. It takes a list of trait specs, each with a name, and joint=True, because neither trait can be simulated before the other.

A curved optimum is the first. Each trait's optimum is a curve of the other's value, written with set_by on reverts_to (Jo5):

import math
from zombi2.params.connection import set_by

traits.simulate_traits(tree, [
    traits.continuous(name="brain", start=0.0, rate=0.1, pull=1.0,
        reverts_to=set_by("traits:group", lambda v: 2.0 * math.tanh(v), step=0.01)),
    traits.continuous(name="group", start=0.0, rate=0.1, pull=0.5,
        reverts_to=set_by("traits:brain", lambda v: 1.5 * math.tanh(0.7 * v), step=0.01)),
    ], joint=True, seed=1)

A bigger group pulls the brain upward, a bigger brain pulls the group upward, and the tanh flattens both effects at the extremes, so the pair settles instead of running away. A trait reads the other by "traits:<name>", the same spelling a joint run uses everywhere else.

A continuous trait with a discrete one is the second. The discrete trait paints the continuous trait's optimum, and the continuous trait's value sets how fast the discrete trait switches (Jo6):

from zombi2.params import Scalar

traits.simulate_traits(tree, [
    traits.discrete(name="habitat", states=["surface", "cave"], start="surface",
        switch=PerLineage(0.6).scaled_by("traits:size", Scalar(0.4), step=0.02)),
    traits.continuous(name="size", start=0.0, rate=0.4, pull=1.0,
        regimes="traits:habitat", reverts_to={"surface": 0.0, "cave": 2.0}),
    ], joint=True, seed=1)

A painted optimum is regimes=, exactly as it is when the discrete trait was grown first (Chapter 8). The only new thing is that the driver is named, because it does not exist yet. reverts_to then gives one optimum per state, and it needs one for every state the painting trait has.

Each run returns one TraitsResult per name, the same object a single-trait run returns. Both models hold every driver still over slices of at most step, so both are approximate, and step is written on the connection rather than on the call. Chapter 8 says what that costs and how to check it.

Two refusals are worth knowing. Traits that set each other's rate — one trait changing how fast the other wanders — are not built. And a pair of straight-line optima is refused rather than run here, because it has an exact answer: the message names the simulate_continuous call to write, drift matrix and optima included.

Literature

Trait models arrive under a thicket of names, and a reader who wants "an OU model" or "a threshold model" should be able to find it. The names live here, in one table, each beside the example that shows it; the example carries the run that made it, so the table does not spell the call out a second time. It organises nothing else in the chapter.

What it does From the literature Gallery
a value diffusing Brownian motion (BM) (Felsenstein 1985) Tr1
diffusion pulled to an optimum Ornstein–Uhlenbeck (OU) (Hansen 1997; Butler & King 2004) Tr2
diffusion rate decays through time Early burst (EB / ACDC) (Harmon et al. 2010) Tr3
diffusion rate drifts between lineages Variable-rates BM (Maliet et al. 2019) Tr4
diffusion rate slows as the clade fills Diversity-dependent / ecological limits (Etienne et al. 2012) Tr5
the optimum differs between painted clades Multi-optimum OU (OUM) (Beaulieu et al. 2012) Tr6
the value jumps at each split Cladogenetic / punctuational change Tr7
a discrete state switching Mk (k-state Markov) Tr8
gains and losses at different rates All-rates-different Mk (ARD) Tr9
discrete driven by continuous liability Threshold / liability (Wright–Felsenstein) (Felsenstein 2012) Tr10
traits evolving together Multivariate BM Tr11
traits reverting together, each to its own optimum Multivariate OU, diagonal drift (Clavel et al. 2015) Tr12
discrete traits evolving together Correlated binary / Pagel (Pagel 1994) Tr13

On the command line

--kind is required, continuous or discrete, and it decides which of the other flags apply, since --rate and the OU flags belong to a continuous trait while --states and --switch belong to a discrete one. Any default would silently pick a model for you, so there is none:

# a continuous (Brownian) trait along a species tree
zombi2 traits out/ --kind continuous \
    --start 0.0 --rate 1.0 --seed 1

# a discrete two-state trait
zombi2 traits out/ --kind discrete \
    --states marine,terrestrial --switch 0.1 --seed 1

# the states the next command drives from — a discrete run writes trait_events.tsv by default,
# the driver file a conditioned genome, sequence or trait run takes (Chapter 8)
zombi2 traits out/ --kind discrete \
    --states cave,surface --switch 0.1 --seed 1

# a second trait, kept apart in out/traits/diet/ by --name,
# whose switch rate depends on the first (Chapter 8)
zombi2 traits out/ --kind discrete --name diet --states plant,fish --seed 2 \
    --switch "PerLineage(0.2).scaled_by('out/traits/trait_events.tsv', {'cave': 5.0, 'surface': 1.0})"

An unnamed trait run writes to out/traits/ and replaces whatever an earlier unnamed run left there, so the file the last command reads belongs to the cave/surface run just above it. --name gives a run its own directory, out/traits/diet/, which is what lets two traits sit side by side.

Every rate flag takes a rate in its written form, --switch as much as --rate, so the expression above is the same text the Python API takes. --switch also takes the keyword's other two shapes: the {'a->b': rate} dict and the k × k matrix.

Two keywords have no flag. correlation= grows several traits in one call, and the command line grows one trait per run. regimes= takes a discrete result object, the painting, used as input to the call. Both stay in the Python API.

The trait evolves on the complete tree, extinct lineages included, so species_complete.nwk is the file to use. An external tree goes in with --from (a Newick file, or another run's directory); if it is not ultrametric you must declare each tip's fate with --tip-fates, a TSV of tip<TAB>extant|extinct|unsampled, because ZOMBI2 will not guess which early-ending tips are extinct.