Merge pull request #25 from fbaumdicker/master

added ecoli model from Lapierre et al

ecoli_example.py

@@ -0,0 +1,23 @@
+"""
+Example of using the stdpopsim library with msprime.
+"""
+import msprime
+
+
+import stdpopsim
+from stdpopsim import e_coli
+
+chrom = e_coli.genome.chromosomes["chr"]
+#recomb_map = chrom.recombination_map()
+
+model = e_coli.LapierreConstant()
+model.debug()
+
+samples=[msprime.Sample(population=0, time=0)]*100
+
+ts = msprime.simulate(
+    samples=samples,
+    mutation_rate=chrom.mean_mutation_rate,
+    **model.asdict())
+print("simulated:", ts.num_trees, ts.num_sites)
+print("mean mut rate", chrom.mean_mutation_rate)

stdpopsim/e_coli.py

@@ -0,0 +1,54 @@
+"""
+Genome and demographic model definitions for Escherichia coli.
+"""
+import msprime
+
+import stdpopsim.models as models
+import stdpopsim.genomes as genomes
+
+
+###########################################################
+#
+# Genome definition
+#
+###########################################################
+
+_chromosomes = []
+_chromosomes.append(genomes.Chromosome(
+        name="chr",
+        length=4641652,
+        mean_mutation_rate=1e-5+2e-4,
+        mean_recombination_rate=0.0))
+# mean_conversion_rate=8.9e-11 # not implemented yet!
+# mean_conversion_length=542 # not implemented yet!
+
+#: :class:`stdpopsim.Genome` definition for E. Coli.
+# Chromosome length data is based on strain K-12.
+genome = genomes.Genome(
+    species="e_coli",
+    chromosomes=_chromosomes,
+    default_genetic_map=None)
+
+
+###########################################################
+#
+# Demographic models
+#
+###########################################################
+
+
+class LapierreConstant(models.Model):
+    """
+    The constant population size model from Lapierre et al. 2016
+    https://doi.org/10.1093/molbev/msw048
+    """
+
+    def __init__(self):
+
+        N_e = 1.8e8
+        # Single population
+        self.population_configurations = [
+            msprime.PopulationConfiguration(initial_size=N_e),
+        ]
+        self.migration_matrix = [[0]]
+        self.demographic_events = []