Python msprime.simulate函数代码示例

 def test_two_populations_no_migration_one_locus(self):
     seed = 1234
     ts1 = msprime.simulate(
         population_configurations=[
             msprime.PopulationConfiguration(10),
             msprime.PopulationConfiguration(10)],
         migration_matrix=np.zeros((2, 2)),
         __tmp_max_time=0.1,
         random_seed=seed)
     ts2 = msprime.simulate(
         population_configurations=[
             msprime.PopulationConfiguration(),
             msprime.PopulationConfiguration()],
         migration_matrix=np.zeros((2, 2)),
         from_ts=ts1,
         demographic_events=[
             msprime.MassMigration(100, 0, 1, 1.0)],
         random_seed=seed)
     tree = ts2.first()
     # We should have two children at the root, and every node below
     # should be in one population.
     root_children = tree.children(tree.root)
     self.assertEqual(len(root_children), 2)
     populations = {ts2.node(u).population: u for u in root_children}
     self.assertEqual(len(populations), 2)
     for pop in [0, 1]:
         for node in tree.nodes(populations[pop]):
             self.assertEqual(ts2.node(node).population, pop)

    def test_gaps(self):
        # SMC simulations should never have adjacent edgesets with
        # a non-zero distance between them and the same parent.
        # First we do a simulation with the standard model to make sure
        # we have plausible parameter values.
        sample_size = 10
        recombination_rate = 20
        random_seed = 1

        ts = msprime.simulate(
            sample_size=sample_size, recombination_rate=recombination_rate,
            random_seed=random_seed)
        edgesets = sorted(ts.edgesets(), key=lambda e: (e.parent, e.left))
        num_found = 0
        for j in range(1, len(edgesets)):
            r = edgesets[j - 1]
            s = edgesets[j]
            if r.right != s.left and r.parent == s.parent:
                num_found += 1
        self.assertGreater(num_found, 10)  # Make a reasonable threshold

        # Now do the same for SMC and SMC'.
        for model in ["smc", "smc_prime"]:
            ts = msprime.simulate(
                sample_size=sample_size, recombination_rate=recombination_rate,
                random_seed=random_seed, model=model)
            edgesets = sorted(ts.edgesets(), key=lambda e: (e.parent, e.left))
            num_found = 0
            for j in range(1, len(edgesets)):
                r = edgesets[j - 1]
                s = edgesets[j]
                if r.right != s.left and r.parent == s.parent:
                    num_found += 1
            self.assertEqual(num_found, 0)

 def test_simple_cases(self):
     for n in range(2, 10):
         st = next(msprime.simulate(n).trees())
         self.verify_sparse_tree(st)
     for n in [11, 13, 19, 101]:
         st = next(msprime.simulate(n).trees())
         self.verify_sparse_tree(st)

 def test_two_populations_migration(self):
     n = 10
     seed = 1234
     ts1 = msprime.simulate(
         population_configurations=[
             msprime.PopulationConfiguration(n),
             msprime.PopulationConfiguration(0)],
         migration_matrix=[[0, 1], [1, 0]],
         random_seed=seed)
     tables = msprime.TableCollection(1)
     tables.populations.add_row()
     tables.populations.add_row()
     for _ in range(n):
         tables.nodes.add_row(
             flags=msprime.NODE_IS_SAMPLE, time=0, population=0)
     ts2 = msprime.simulate(
         from_ts=tables.tree_sequence(), start_time=0,
         population_configurations=[
             msprime.PopulationConfiguration(),
             msprime.PopulationConfiguration()],
         migration_matrix=[[0, 1], [1, 0]],
         random_seed=seed)
     tables1 = ts1.dump_tables()
     tables2 = ts2.dump_tables()
     tables1.provenances.clear()
     tables2.provenances.clear()
     self.assertEqual(tables1, tables2)

 def test_sequence_length(self):
     from_ts = msprime.simulate(
         5, __tmp_max_time=0.1, random_seed=5, length=5)
     start_time = from_ts.tables.nodes.time.max()
     final_ts = msprime.simulate(
         from_ts=from_ts, start_time=start_time, random_seed=2, length=5)
     self.verify_from_tables(from_ts, final_ts, start_time)
     self.verify_simulation_completed(final_ts)

 def test_single_locus_max_time(self):
     from_ts = msprime.simulate(20, __tmp_max_time=1, random_seed=5)
     self.assertGreater(max(tree.num_roots for tree in from_ts.trees()), 1)
     start_time = from_ts.tables.nodes.time.max()
     final_ts = msprime.simulate(
         from_ts=from_ts, start_time=start_time, random_seed=2)
     self.verify_from_tables(from_ts, final_ts, start_time)
     self.verify_simulation_completed(final_ts)

 def test_from_single_locus_decapitated(self):
     ts = msprime.simulate(10, random_seed=5)
     from_ts = tsutil.decapitate(ts, ts.num_edges // 2)
     start_time = from_ts.tables.nodes.time.max()
     final_ts = msprime.simulate(
         from_ts=from_ts, start_time=start_time, random_seed=2)
     self.verify_from_tables(from_ts, final_ts, start_time)
     self.verify_simulation_completed(final_ts)

 def test_provenance(self):
     ts = msprime.simulate(10)
     self.assertEqual(ts.num_provenances, 1)
     self.verify_provenance(ts.provenance(0))
     # TODO check the form of the dictionary
     for ts in msprime.simulate(10, num_replicates=10):
         self.assertEqual(ts.num_provenances, 1)
         self.verify_provenance(ts.provenance(0))

 def test_sequence_length_mismatch(self):
     base_ts = self.get_example_base(length=5)
     for bad_length in [1, 4.99, 5.01, 100]:
         with self.assertRaises(ValueError):
             msprime.simulate(from_ts=base_ts, start_time=100, length=bad_length)
         recomb_map = msprime.RecombinationMap.uniform_map(bad_length, 1)
         with self.assertRaises(ValueError):
             msprime.simulate(
                 from_ts=base_ts, start_time=100, recombination_map=recomb_map)

 def test_decapitated_mutations(self):
     ts = msprime.simulate(10, random_seed=5, mutation_rate=10)
     from_ts = tsutil.decapitate(ts, ts.num_edges // 2)
     self.assertGreater(from_ts.num_mutations, 0)
     start_time = from_ts.tables.nodes.time.max()
     final_ts = msprime.simulate(
         from_ts=from_ts, start_time=start_time, random_seed=2)
     self.verify_from_tables(from_ts, final_ts, start_time)
     self.verify_simulation_completed(final_ts)

 def test_individuals(self):
     from_ts = msprime.simulate(25, random_seed=5, __tmp_max_time=0.5)
     self.assertTrue(any(tree.num_roots > 1 for tree in from_ts.trees()))
     from_ts = tsutil.insert_random_ploidy_individuals(from_ts, seed=2)
     start_time = from_ts.tables.nodes.time.max()
     final_ts = msprime.simulate(
         from_ts=from_ts, start_time=start_time, random_seed=2)
     self.verify_from_tables(from_ts, final_ts, start_time)
     self.verify_simulation_completed(final_ts)

 def test_from_multi_locus_old_recombination(self):
     ts = msprime.simulate(10, recombination_rate=2, random_seed=5)
     self.assertGreater(ts.num_trees, 1)
     from_ts = tsutil.decapitate(ts, ts.num_edges // 2)
     start_time = from_ts.tables.nodes.time.max()
     final_ts = msprime.simulate(
         from_ts=from_ts, start_time=start_time, random_seed=2,
         recombination_rate=2)
     self.verify_from_tables(from_ts, final_ts, start_time)
     self.verify_simulation_completed(final_ts)

 def test_models(self):
     # Exponential growth of 0 and constant model should be identical.
     for n in [2, 10, 100]:
         m1 = msprime.PopulationConfiguration(n, growth_rate=0)
         m2 = msprime.PopulationConfiguration(n, initial_size=1.0)
         st1 = next(msprime.simulate(
             random_seed=1, population_configurations=[m1]).trees())
         st2 = next(msprime.simulate(
             random_seed=1, population_configurations=[m2]).trees())
         self.assertEqual(st1.parent_dict, st2.parent_dict)

 def test_population_mismatch(self):
     for N in range(1, 5):
         base_ts = self.get_example_base(num_populations=N)
         start_time = max(node.time for node in base_ts.nodes())
         for k in range(1, N):
             if k != N:
                 with self.assertRaises(ValueError):
                     msprime.simulate(
                         from_ts=base_ts, start_time=start_time,
                         population_configurations=[
                             msprime.PopulationConfiguration() for _ in range(k)])

 def test_zero_recombination_rate(self):
     from_ts = msprime.simulate(
         sample_size=4, __tmp_max_time=1, random_seed=5, recombination_rate=1)
     self.assertGreater(max(tree.num_roots for tree in from_ts.trees()), 1)
     start_time = from_ts.tables.nodes.time.max()
     with self.assertRaises(_msprime.InputError):
         # Raises:
         # The specified recombination map is too coarse to translate...
         msprime.simulate(
             from_ts=from_ts, start_time=start_time, random_seed=2,
             recombination_rate=0)