1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
|
#! /usr/bin/env python
##############################################################################
## DendroPy Phylogenetic Computing Library.
##
## Copyright 2010-2015 Jeet Sukumaran and Mark T. Holder.
## All rights reserved.
##
## See "LICENSE.rst" for terms and conditions of usage.
##
## If you use this work or any portion thereof in published work,
## please cite it as:
##
## Sukumaran, J. and M. T. Holder. 2010. DendroPy: a Python library
## for phylogenetic computing. Bioinformatics 26: 1569-1571.
##
##############################################################################
"""
Tests of birth-death model fitting.
"""
import unittest
import dendropy
from dendropy.test.support.mockrandom import MockRandom
from dendropy.test.support import pathmap
from dendropy.model import birthdeath
class PureBirthModelEstimationTestCase(unittest.TestCase):
def testBasicEst(self):
# list of tuples, (birth-rate, log-likelihood)
expected_results = (
# birth rate # log-likelihood
(0.02879745490817826186758, -59.41355682054444287132355),
(0.03074708092192806122012, -57.38280732060526645454956),
(0.02539588437187430269848, -63.31025321526630023072357),
(0.02261951969802362960582, -66.89924384677527768872096),
(0.02804607815688910446572, -60.23314120509648716961237),
(0.02748663302756114423797, -60.85775993426526042640035),
(0.02816256618562208019485, -60.10465085978295007862471),
(0.03592126646048716259729, -52.56123967307649991198559),
(0.02905144990609926855529, -59.14133401672411594063306),
(0.02703739196351075124714, -61.36860953277779628933786),
(0.01981322730236481297061, -71.00561162515919022553135),
)
trees = dendropy.TreeList.get_from_path(
pathmap.tree_source_path("pythonidae.reference-trees.newick"), "newick")
self.assertEqual(len(trees), len(expected_results))
for tree, expected_result in zip(trees, expected_results):
obs_result1 = birthdeath.fit_pure_birth_model(tree=tree, ultrametricity_precision=1e-5)
obs_result2 = birthdeath.fit_pure_birth_model(internal_node_ages=tree.internal_node_ages(ultrametricity_precision=1e-5))
for obs_result in (obs_result1, obs_result2):
self.assertAlmostEqual(obs_result["birth_rate"], expected_result[0], 5)
self.assertAlmostEqual(obs_result["log_likelihood"], expected_result[1], 5)
class BirthDeathTreeTest(unittest.TestCase):
def testGSABD(self):
"""test that the birth-death process produces the correct number of tips with GSA."""
_RNG = MockRandom()
for num_leaves in range(2, 15):
t = birthdeath.birth_death_tree(birth_rate=1.0, death_rate=0.2, ntax=num_leaves, gsa_ntax=3*num_leaves, rng=_RNG)
self.assertTrue(t._debug_tree_is_valid())
self.assertEqual(num_leaves, len(t.leaf_nodes()))
def testYule(self):
"""test that the pure-birth process produces the correct number of tips."""
_RNG = MockRandom()
for num_leaves in range(2, 20):
t = birthdeath.birth_death_tree(birth_rate=1.0, death_rate=0.0, ntax=num_leaves, rng=_RNG)
self.assertTrue(t._debug_tree_is_valid())
self.assertEqual(num_leaves, len(t.leaf_nodes()))
def testGSA(self):
"""test that the pure-birth process produces the correct number of tips with GSA."""
_RNG = MockRandom()
for num_leaves in range(2, 20):
t = birthdeath.birth_death_tree(birth_rate=1.0, death_rate=0.0, ntax=num_leaves, gsa_ntax=4*num_leaves, rng=_RNG)
self.assertTrue(t._debug_tree_is_valid())
self.assertEqual(num_leaves, len(t.leaf_nodes()))
def testBDTree(self):
"""PureCoalescentTreeTest -- tree generation without checking [TODO: checks]"""
_RNG = MockRandom()
for num_leaves in range(2, 20):
t = birthdeath.birth_death_tree(birth_rate=1.0, death_rate=0.2, ntax=num_leaves, rng=_RNG)
self.assertTrue(t._debug_tree_is_valid())
self.assertEqual(num_leaves, len(t.leaf_nodes()))
if __name__ == "__main__":
unittest.main()
|
