
Research
/Security News
Toptal’s GitHub Organization Hijacked: 10 Malicious Packages Published
Threat actors hijacked Toptal’s GitHub org, publishing npm packages with malicious payloads that steal tokens and attempt to wipe victim systems.
Simulation of rooted phylogenetic trees under a given Multi-Type Birth–Death model (with or without Contact Tracing) (with or without a Skyline).
Simulation of rooted phylogenetic trees under a given Multi-Type Birth–Death (MTBD) model, with or without contact tracing (CT), and with or without Skyline.
treesimulator is described in the S1 Appendix of the following article:
Anna Zhukova, Olivier Gascuel. Accounting for contact tracing in epidemiological birth-death models. medRxiv 2024.09.09.24313296; doi:10.1101/2024.09.09.24313296
The MTBD models were introduced by Stadler & Bonhoeffer [Philos. Trans. R. Soc. B 2013].
An MTBD model with m states has
m(m-1) state transition rate parameters:
m2 transmission rate parameters:
m removal (becoming non-infectious) rate parameters:
m sampling probability upon removal parameters:
Using these probabilities, one can calculate the equilibrium frequencies $πi$ of the model's states, where 0 ≤ πi ≤ 1 and π1 + ... + πm = 1.
The MTBD model has the following epidemiological parameters:
Contact tracing extension was introduced by Zhukova & Gascuel [medRxiv 2024]. It adds two parameters to the initial MTBD model:
and a meta-parameter κ, which defines how many most recent contacts can be notified by each index case. Each is notified independently, with the probability υ.
CT extension adds m notified contact states (1C, ..., mC) to its original MTBD model, where iC is a notified version of the state i. Transition rates for a notified state iC are analogous to those of i (μiCjC=μij), while transitions from non-notified states to notified ones and vice versa are not allowed (μiCj=μjiC=0). Transmission rates for a notified state iC are the same as those of i (λiCj=λij), where the recipients are always in a non-notified state (λiCjC=λijC=0). The removal rate for a notified state iC is φ.
For CT models exit times from notified contact states iC is calculated as:
Skyline was introduced by Stadler et al. [PNAS 2013] and extended to MTBD by Kühnert et al. [MBE 2016]. It enables piece-wise constant parameter value changes. To use a skyline with k models, one needs to specify k-1 model change times t1, ..., k-1, and k sets of model parameters (see above). At time 0 the simulation starts with model 1, it switches to models 2 at time t1, etc. All the models in the Skyline must have the same states. CT-related parameters can also change at skyline changing times, in that case if some skyline intervals do not have CT, υ=0 and any value for φ must be specified for them. The same κ value is shared among all the skyline -CT models.
We pay particular interest to the classical BD model, the BD Exposed-Infectious (BDEI) model, and BD with super-spreading (BDSS), as they are described in [Voznica et al. 2021], and to their -CT(κ) versions.
1 state:
3 parameters:
Epidemiological parameters:
2 states:
5 parameters:
2 states:
4 parameters:
BDEI-specific epidemiological parameter:
4 states:
6 parameters:
2 states:
5(+1) parameters:
λnn = λII -- transmission rate from I to I
λns = λIS -- transmission rate from I to S
λsn = λSI -- transmission rate from S to I
λss = λSS -- transmission rate from S to S
(with a constraint that λss/λns=λsn/λnn)
ψ = ψI = ψS -- removal rate
p = pI = pS -- sampling probability upon removal
BDSS-specific epidemiological parameters:
4 states:
7(+1) parameters:
λnn = λII = λICI -- transmission rate from a normal spreader (notified or not) to I
λns = λIS = λICS -- transmission rate from a normal spreader (notified or not) to S
λsn = λSI = λSCI -- transmission rate from a superspreader (notified or not) to I
λss = λSS = λSCS -- transmission rate from a superspreader (notified or not) to S
(with a constraint that λss/λns=λsn/λnn)
ψ = ψI = ψS -- removal rate of a non-notified state
p = pI = pS -- sampling probability upon removal of a non-notified state
υ -- probability to notify contacts upon sampling
φ -- notified contact removal and sampling rate: φ >> ψ
3 states:
6 parameters:
BDEISS-specific epidemiological parameters:
6 states:
8 parameters:
There are 4 alternative ways to run treesimulator on your computer: with docker, apptainer, in Python3, or via command line (requires installation with Python3).
You could either install python (version 3.6 or higher) system-wide and then install treesimulator via pip:
sudo apt install -y python3 python3-pip python3-setuptools python3-distutils
pip3 install treesimulator
or alternatively, you could install python (version 3.6 or higher) and treesimulator via conda (make sure that conda is installed first).
(Optional) to install treesimulator in a new conda environment (e.g., called phyloenv below), first create and activate the environment:
conda create --name phyloenv python=3.6
conda activate phyloenv
Install treesimulator with conda
conda install treesimulator
If you installed treesimulator in a conda environment (here named phyloenv), do not forget to first activate it, e.g.
conda activate phyloenv
The following command simulates a tree with 200-500 tips under the BD model, with λ=0.5, ψ=0.25, p=0.5, and saves it to the file tree.nwk, while saving the parameters to the comma-separated file params.csv:
generate_bd --min_tips 200 --max_tips 500 \
--la 0.5 --psi 0.25 --p 0.5 \
--nwk tree.nwk --log params.csv
The following command simulates a tree with 200-500 tips under the BD-CT(1) model, with λ=0.5, ψ=0.25, p=0.5, φ=2.5, υ=0.2, and allowing to notify only the most recent contact of each sampled index case. The simulated tree is saved to the file tree.nwk, while the model parameters are saved to the comma-separated file params.csv:
generate_bd --min_tips 200 --max_tips 500 \
--la 0.5 --psi 0.25 --p 0.5 \
--phi 2.5 --upsilon 0.2 --max_notified_contacts 1 \
--nwk tree.nwk --log params.csv
The following command simulates a tree with 200-500 tips under the BD-CT(1)-Skyline model with two time intervals, with λ=0.5, ψ=0.25, p=0.5, φ=2.5, υ=0 between t=0 and t=3, and λ=1, ψ=0.25, p=0.75, φ=2.5, υ=0.2 starting at t=3, and allowing to notify only the most recent contact of each sampled index case. The simulated tree is saved to the file tree.nwk, while the model parameters are saved to the comma-separated file params.csv:
generate_bd --min_tips 200 --max_tips 500 \
--la 0.5 1 --psi 0.25 0.25 --p 0.5 0.75 \
--phi 2.5 2.5 --upsilon 0 0.2 --max_notified_contacts 1 \
--skyline_times 3 \
--nwk tree.nwk --log params.csv
To see detailed options, run:
generate_bd --help
The following command simulates a tree with 200-500 tips under the BDEI model, with μ=1, λ=0.5, ψ=0.25, p=0.5, and saves it to the file tree.nwk, while saving the parameters to the comma-separated file params.csv:
generate_bdei --min_tips 200 --max_tips 500 \
--mu 1 --la 0.5 --psi 0.25 --p 0.5 \
--nwk tree.nwk --log params.csv
The following command simulates a tree with 200-500 tips under the BDEI-CT(2) model, with μ=1, λ=0.5, ψ=0.25, p=0.5, φ=2.5, υ=0.2, and allowing to notify last two contacts of each sampled index case. The simulated tree is saved to the file tree.nwk, while the model parameters are saved to the comma-separated file params.csv:
generate_bdei --min_tips 200 --max_tips 500 \
--mu 1 --la 0.5 --psi 0.25 --p 0.5 \
--phi 2.5 --upsilon 0.2 --max_notified_contacts 2 \
--nwk tree.nwk --log params.csv
The following command simulates a tree with 200-500 tips under the BDEI-CT(2)-Skyline model with three time intervals, with μ=1, λ=0.5, ψ=0.25, p=0.2, φ=2.5, υ=0.2, between t=0 and t=2, with μ=1, λ=0.5, ψ=0.3, p=0.3, φ=2.5, υ=0.3, between t=2 and t=3, and μ=1, λ=0.5, ψ=0.5, p=0.5, φ=5, υ=0.3 starting at t=3, and allowing to notify last two contacts of each sampled index case. The simulated tree is saved to the file tree.nwk, while the model parameters are saved to the comma-separated file params.csv:
generate_bdei --min_tips 200 --max_tips 500 \
--mu 1 1 1 --la 0.5 0.5 0.5 --psi 0.25 0.3 0.5 --p 0.2 0.3 0.5 \
--phi 2.5 2.5 5 --upsilon 0.2 0.3 0.3 --max_notified_contacts 2 \
--skyline_times 2 3 \
--nwk tree.nwk --log params.csv
To see detailed options, run:
generate_bdei --help
The following command simulates a tree with 200-500 tips under the BDSS model, with λnn=0.1, λns=0.3, λsn=0.5, λss=1.5, ψ=0.25, p=0.5, and saves it to the file tree.nwk, while saving the parameters to the comma-separated file params.csv:
generate_bdss --min_tips 200 --max_tips 500 \
--la_nn 0.1 --la_ns 0.3 --la_sn 0.5 --la_ss 1.5 --psi 0.25 --p 0.5 \
--nwk tree.nwk --log params.csv
The following command simulates a tree with 200-500 tips under the BDSS-CT(3) model, with λnn=0.1, λns=0.3, λsn=0.5, λss=1.5, ψ=0.25, p=0.5, φ=2.5, υ=0.2, and allowing to notify last three contacts of each sampled index case. The simulated tree is saved to the file tree.nwk, while the model parameters are saved to the comma-separated file params.csv:
generate_bdss --min_tips 200 --max_tips 500 \
--la_nn 0.1 --la_ns 0.3 --la_sn 0.5 --la_ss 1.5 --psi 0.25 --p 0.5 \
--phi 2.5 --upsilon 0.2 --max_notified_contacts 3 \
--nwk tree.nwk --log params.csv
The following command simulates a tree with 200-500 tips under the BDSS-CT(3)-Skyline model with two time intervals, with λnn=0.1, λns=0.3, λsn=0.5, λss=1.5, ψ=0.25, p=0.5, φ=2.5, υ=0.2 between t=0 and t=2, and λnn=0.1, λns=0.3, λsn=1, λss=3, ψ=0.25, p=0.5, φ=5, υ=0.5 starting at t=2, and allowing to notify last three contacts of each sampled index case. The simulated tree is saved to the file tree.nwk, while the model parameters are saved to the comma-separated file params.csv:
generate_bdss --min_tips 200 --max_tips 500 \
--la_nn 0.1 0.1 --la_ns 0.3 0.3 --la_sn 0.5 1 --la_ss 1.5 3 --psi 0.25 0.25 --p 0.5 0.5 \
--phi 2.5 5 --upsilon 0.2 0.5 --max_notified_contacts 3 \
--skyline_times 2 \
--nwk tree.nwk --log params.csv
To see detailed options, run:
generate_bdss --help
The following command simulates a tree with 200-500 tips under the BDEISS model, with μn=0.1, μs=0.3, λn=0.5, λs=1.5, ψ=0.25, p=0.5, and saves it to the file tree.nwk, while saving the parameters to the comma-separated file params.csv:
generate_bdeiss --min_tips 200 --max_tips 500 \
--mu_n 0.1 --mu_s 0.3 --la_n 0.5 --la_s 1.5 --psi 0.25 --p 0.5 \
--nwk tree.nwk --log params.csv
The following command simulates a tree with 200-500 tips under the BDEISS-CT(1) model, with μn=0.1, μs=0.3, λn=0.5, λs=1.5, ψ=0.25, p=0.5, φ=2.5, υ=0.2, and allowing to notify the last contact of each sampled index case. The simulated tree is saved to the file tree.nwk, while the model parameters are saved to the comma-separated file params.csv:
generate_bdeiss --min_tips 200 --max_tips 500 \
--mu_n 0.1 --mu_s 0.3 --la_n 0.5 --la_s 1.5 --psi 0.25 --p 0.5 \
--phi 2.5 --upsilon 0.2 --max_notified_contacts 1 \
--nwk tree.nwk --log params.csv
The following command simulates a tree with 200-500 tips under the BDSS-CT(1)-Skyline model with two time intervals, with μn=0.1, μs=0.3, λn=0.5, λs=1.5, ψ=0.25, p=0.5, φ=2.5, υ=0.2 between t=0 and t=2, and μn=0.1, μs=0.3, λn=0.5, λs=1.5, ψ=0.25, p=0.2, φ=2.5, υ=0.5 starting at t=2, and allowing to notify the last contact of each sampled index case. The simulated tree is saved to the file tree.nwk, while the model parameters are saved to the comma-separated file params.csv:
generate_bdeiss --min_tips 200 --max_tips 500 \
--mu_n 0.1 0.1 --mu_s 0.3 0.3 --la_n 0.5 0.5 --la_s 1.5 1.5 --psi 0.25 0.25 --p 0.5 0.2 \
--phi 2.5 2.5 --upsilon 0.2 0.5 --max_notified_contacts 1 \
--skyline_times 2 \
--nwk tree.nwk --log params.csv
To see detailed options, run:
generate_bdss --help
The following command simulates a tree with 200-500 tips under a generic MTBD model, with two states A and B, with μab=0.6, μba=0.7 (note that μaa=μbb=0 as only transitions between different states are possible!), λaa=0.1, λab=0.2, λba=0.3, λbb=0.4, ψa=0.05, ψb=0.08, p=a0.15, p=b0.65, and saves it to the file tree.nwk, while saving the parameters to the comma-separated file params.csv:
generate_mtbd --min_tips 200 --max_tips 500 \
--states A B \
--transition_rates 0 0.6 0.7 0 \
--transmission_rates 0.1 0.2 0.3 0.4 \
--removal_rates 0.05 0.08 \
--sampling_probabilities 0.15 0.65 \
--nwk tree.nwk --log params.csv
The following command simulates a tree with 200-500 tips under a generic MTBD-CT(1) model, with two states A and B, with μab=0.6, μba=0.7, λaa=0.1, λab=0.2, λba=0.3, λbb=0.4, ψa=0.05, ψb=0.08, p=a0.15, p=b0.65, φ=2.5, υ=0.2, and allowing to notify only the most recent contact of each sampled index case. The simulated tree is saved to the file tree.nwk, while the model parameters are saved to the comma-separated file params.csv:
generate_mtbd --min_tips 200 --max_tips 500 \
--states A B \
--transition_rates 0 0.6 0.7 0 \
--transmission_rates 0.1 0.2 0.3 0.4 \
--removal_rates 0.05 0.08 \
--sampling_probabilities 0.15 0.65 \
--phi 2.5 --upsilon 0.2 --max_notified_contacts 1 \
--nwk tree.nwk --log params.csv
The following command simulates a tree with 200-500 tips under a generic MTBD-CT(1)-Skyline model, with two states A and B, with μab=0.6, μba=0.7, λaa=0.1, λab=0.2, λba=0.3, λbb=0.4, ψa=0.05, ψb=0.08, p=a0.15, p=b0.65, φ=2.5, υ=0.2 between t=0 and t=8, and μab=1.6, μba=1.7, λaa=1.1, λab=1.2, λba=1.3, λbb=1.4, ψa=1.05, ψb=1.08, p=a0.1, p=b0.6, φ=3.5, υ=0.4 starting at t=8, and allowing to notify only the most recent contact of each sampled index case. The simulated tree is saved to the file tree.nwk, while the model parameters are saved to the comma-separated file params.csv:
generate_mtbd --min_tips 200 --max_tips 500 \
--states A B \
--transition_rates 0 0.6 0.7 0 0 1.6 1.7 0 \
--transmission_rates 0.1 0.2 0.3 0.4 1.1 1.2 1.3 1.4 \
--removal_rates 0.05 0.08 1.05 1.08 \
--sampling_probabilities 0.15 0.65 0.1 0.6 \
--phi 2.5 3.5 --upsilon 0.2 0.4 --max_notified_contacts 1 \
--skyline_times 8 \
--nwk tree.nwk --log params.csv
To see detailed options, run:
generate_mtbd --help
To simulate trees with 200-500 tips under the above models and settings:
from treesimulator.generator import generate
from treesimulator import save_forest
from treesimulator.mtbd_models import Model, BirthDeathModel, BirthDeathExposedInfectiousModel,
BirthDeathWithSuperSpreadingModel, BirthDeathExposedInfectiousWithSuperSpreadingModel, CTModel
# 1. BD, BD-CT(1) and BD-CT(1)-Skyline
## BD model
bd_model = BirthDeathModel(p=0.5, la=0.5, psi=0.25)
print(bd_model.get_epidemiological_parameters())
[bd_tree], _, _ = generate([bd_model], min_tips=200, max_tips=500)
save_forest([bd_tree], 'BD_tree.nwk')
## Adding -CT to the model above
bdct_model = CTModel(model=bd_model, upsilon=0.2, phi=2.5)
[bdct_tree], _, _ = generate([bdct_model], min_tips=200, max_tips=500, max_notified_contacts=1)
save_forest([bdct_tree], 'BDCT_tree.nwk')
## BD-CT(1)-Skyline models
bdct_model_1 = CTModel(BirthDeathModel(p=0.5, la=0.5, psi=0.25),
upsilon=0, phi=2.5)
bdct_model_2 = CTModel(BirthDeathModel(p=0.75, la=1, psi=0.25),
upsilon=0.2, phi=2.5)
[bdct_skyline_tree], _, _ = generate([bdct_model_1, bdct_model_2], skyline_times=[3],
min_tips=200, max_tips=500, max_notified_contacts=1)
save_forest([bdct_skyline_tree], 'BDCTSkyline_tree.nwk')
# BDEI, BDEI-CT(2) and BDEI-CT(2)-Skyline
## BDEI model
bdei_model = BirthDeathExposedInfectiousModel(p=0.5, mu=1, la=0.5, psi=0.25)
print(bdei_model.get_epidemiological_parameters())
[bdei_tree], _, _ = generate([bdei_model], min_tips=200, max_tips=500)
save_forest([bdei_tree], 'BDEI_tree.nwk')
## Adding -CT to the model above
bdeict_model = CTModel(model=bdei_model, upsilon=0.2, phi=2.5)
[bdeict_tree], _, _ = generate([bdeict_model], min_tips=200, max_tips=500, max_notified_contacts=2)
save_forest([bdeict_tree], 'BDEICT_tree.nwk')
## BDEI-CT(2)-Skyline with three time intervals
bdeict_model_1 = CTModel(model=BirthDeathExposedInfectiousModel(p=0.2, mu=1, la=0.5, psi=0.25), upsilon=0.2,
phi=2.5)
bdeict_model_2 = CTModel(model=BirthDeathExposedInfectiousModel(p=0.3, mu=1, la=0.5, psi=0.3), upsilon=0.3, phi=2.5)
bdeict_model_3 = CTModel(model=BirthDeathExposedInfectiousModel(p=0.5, mu=1, la=0.5, psi=0.5), upsilon=0.3, phi=5)
[bdeict_skyline_tree], _, _ = generate([bdeict_model_1, bdeict_model_2, bdeict_model_3], skyline_times=[2, 3],
min_tips=200, max_tips=500, max_notified_contacts=2)
save_forest([bdeict_skyline_tree], 'BDEICTSkyline_tree.nwk')
# BDSS, BDSS-CT(3) and BDSS-CT(3)-Skyline
## BDSS model
bdss_model = BirthDeathWithSuperSpreadingModel(p=0.5, la_nn=0.1, la_ns=0.3, la_sn=0.5, la_ss=1.5, psi=0.25)
print(bdss_model.get_epidemiological_parameters())
[bdss_tree], _, _ = generate([bdss_model], min_tips=200, max_tips=500)
save_forest([bdss_tree], 'BDSS_tree.nwk')
## Adding -CT to the model above
bdssct_model = CTModel(model=bdss_model, upsilon=0.2, phi=2.5)
[bdssct_tree], _, _ = generate([bdssct_model], min_tips=200, max_tips=500, max_notified_contacts=3)
save_forest([bdssct_tree], 'BDSSCT_tree.nwk')
## BDSS-CT(3)-Skyline with two time intervals, using the model above for the first interval
bdssct_model_2 = CTModel(
model=BirthDeathWithSuperSpreadingModel(p=0.5, la_nn=0.1, la_ns=0.3, la_sn=1, la_ss=3, psi=0.25),
upsilon=0.5, phi=5)
[bdssct_skyline_tree], _, _ = generate([bdssct_model, bdssct_model_2], skyline_times=[2], min_tips=200, max_tips=500,
max_notified_contacts=3)
save_forest([bdssct_skyline_tree], 'BDSSCTSkyline_tree.nwk')
# BDEISS, BDEISS-CT(1) and BDEISS-CT(1)-Skyline
## BDEISS model
bdeiss_model = BirthDeathExposedInfectiousWithSuperSpreadingModel(p=0.5, mu_n=0.1, mu_s=0.3, la_n=0.5, la_s=1.5,
psi=0.25)
print(bdeiss_model.get_epidemiological_parameters())
[bdeiss_tree], _, _ = generate([bdeiss_model], min_tips=200, max_tips=500)
save_forest([bdeiss_tree], 'BDEISS_tree.nwk')
## Adding -CT to the model above
bdeissct_model = CTModel(model=bdeiss_model, upsilon=0.2, phi=2.5)
[bdeissct_tree], _, _ = generate([bdeissct_model], min_tips=200, max_tips=500, max_notified_contacts=1)
save_forest([bdeissct_tree], 'BDEISSCT_tree.nwk')
## BDEISS-CT(1)-Skyline with two time intervals, using the model above for the first interval
bdeissct_model_2 = CTModel(
model=BirthDeathExposedInfectiousWithSuperSpreadingModel(p=0.2, mu_n=0.1, mu_s=0.3, la_n=0.5, la_s=1.5, psi=0.25),
upsilon=0.5, phi=5)
[bdeissct_skyline_tree], _, _ = generate([bdeissct_model, bdeissct_model_2], skyline_times=[2], min_tips=200,
max_tips=500,
max_notified_contacts=1)
save_forest([bdeissct_skyline_tree], 'BDEISSCTSkyline_tree.nwk')
# MTBD, MTBD-CT(1) and MTBD-CT(1)-Skyline
## MTBD model with two states: A and B
mtbd_model = Model(states=['A', 'B'], transition_rates=[[0, 0.6], [0.7, 0]],
transmission_rates=[[0.1, 0.2], [0.3, 0.4]],
removal_rates=[0.05, 0.08], ps=[0.15, 0.65])
[mtbd_tree], _, _ = generate([mtbd_model], min_tips=200, max_tips=500)
save_forest([mtbd_tree], 'MTBD_tree.nwk')
## Adding -CT to the model above
mtbdct_model = CTModel(model=mtbd_model, upsilon=0.2, phi=2.5)
[mtbdct_tree], _, _ = generate([mtbdct_model], min_tips=200, max_tips=500, max_notified_contacts=1)
save_forest([mtbdct_tree], 'MTBDCT_tree.nwk')
## MTBD-CT(1)-Skyline with two time intervals, using the model above for the first interval
mtbdct_model_2 = CTModel(model=Model(states=['A', 'B'], transition_rates=[[0, 1.6], [1.7, 0]],
transmission_rates=[[1.1, 1.2], [1.3, 1.4]],
removal_rates=[1.05, 1.08], ps=[0.1, 0.6]),
upsilon=0.4, phi=3.5)
[mtbdct_skyline_tree], _, _ = generate([mtbdct_model, mtbdct_model_2], skyline_times=[8],
min_tips=200, max_tips=500, max_notified_contacts=1)
save_forest([mtbdct_skyline_tree], 'MTBDCTSkyline_tree.nwk')
Once apptainer is installed, run the following command:
apptainer run docker://evolbioinfo/treesimulator
This will launch a terminal session within the container, in which you can run treesimulator following the instructions for the command line ("Basic usage in a command line") above.
FAQs
Simulation of rooted phylogenetic trees under a given Multi-Type Birth–Death model (with or without Contact Tracing) (with or without a Skyline).
We found that treesimulator demonstrated a healthy version release cadence and project activity because the last version was released less than a year ago. It has 1 open source maintainer collaborating on the project.
Did you know?
Socket for GitHub automatically highlights issues in each pull request and monitors the health of all your open source dependencies. Discover the contents of your packages and block harmful activity before you install or update your dependencies.
Research
/Security News
Threat actors hijacked Toptal’s GitHub org, publishing npm packages with malicious payloads that steal tokens and attempt to wipe victim systems.
Research
/Security News
Socket researchers investigate 4 malicious npm and PyPI packages with 56,000+ downloads that install surveillance malware.
Security News
The ongoing npm phishing campaign escalates as attackers hijack the popular 'is' package, embedding malware in multiple versions.