Merge pull request #612 from molpopgen/discrete_crossover_positions

Add support for discrete crossover positions

doc/examples/precapitation.rst

@@ -3,6 +3,10 @@
 Starting a simulation with a tree sequence from msprime
 ================================================================================
 
+.. todo::
+
+   Document what to do with discrete crossover positions
+
 The following command line uses :mod:`msprime` to simulate under the discrete
 time Wright-Fisher model using the methods described in [Nelson2020]_.  Then,
 `fwdpy11` simulates for another 1,000 generations.  From the end of each simulation,

doc/examples/recapitation.rst

@@ -3,6 +3,10 @@
 Finishing a simulation with a tree sequence from msprime
 ================================================================================
 
+.. todo::
+
+   Document what to do with discrete crossover positions
+
 This example is complementary to :ref:`precapitation`.  Rather than starting
 with a tree sequence from :mod:`msprime`, we instead finish a simulation by "coalescing
 back" the first generation of the simulation using :mod:`msprime`.  [Haller2019]_

doc/pages/tutorial.rst

@@ -397,6 +397,28 @@ one-half, meaning that each region is assorting independently (50 `cM` between e
 
 As an aside, this example is not creating objects in order by their positions.  Such ordering is not required.
 
+Beginning in version ``0.12.0``, it is possible to restrict crossover positions to integer values.
+For the examples given above, crossover positions are floating-point values sampled uniformly from :math:`[beg, end)`.
+To restrict positions to integer values, we pass ``discrete=True`` when creating object instances:
+
+.. jupyter-execute::
+
+    recregions = [
+        fwdpy11.PoissonInterval(beg=0, end=5, mean=2e-3 / 3, discrete=True),
+        fwdpy11.PoissonInterval(beg=5, end=10, mean=1e-3 / 3, discrete=True),
+    ]
+
+Now, breakpoints from the first region will only take on values of ``0``, ``1``, ``2``, ``3``, or ``4``.
+
+Setting ``discrete=True`` requires the following:
+
+* Values for ``beg`` and ``end`` must be :class:`int`.  Thus, ``1`` is valid but ``1.0`` will raise a ``TypeError``.
+* ``end - beg`` must be ``> 1``.  This requirement prevents you from using ``beg=0`` and ``end=1``, for example,
+  which would result in the only possible crossover position being ``0``.
+* You must be more careful when using :mod:`msprime` to start/finish a simulation.
+  See :ref:`here <precapitation>` and :ref:`here <recapitation>` for details.
+
+
 The following classes are available:
 
 * :class:`fwdpy11.PoissonInterval`

fwdpy11/CMakeLists.txt

@@ -11,13 +11,7 @@ set(FWDPP_TYPES_SOURCES src/fwdpp_types/init.cc ;
     src/fwdpp_types/EdgeTable.cc
     src/fwdpp_types/MutationTable.cc
     src/fwdpp_types/SiteTable.cc
-    src/fwdpp_types/TableCollection.cc
-    src/fwdpp_types/GeneticMapUnit.cc 
-    src/fwdpp_types/PoissonInterval.cc 
-    src/fwdpp_types/BinomialPoint.cc
-    src/fwdpp_types/PoissonPoint.cc
-    src/fwdpp_types/FixedCrossovers.cc
-    src/fwdpp_types/BinomialInterval.cc)
+    src/fwdpp_types/TableCollection.cc)
 
 
 set(FWDPY11_TYPES_SOURCES src/fwdpy11_types/init.cc
@@ -32,7 +26,8 @@ set(FWDPY11_TYPES_SOURCES src/fwdpy11_types/init.cc
     src/fwdpy11_types/DiploidPopulation.cc
     src/fwdpy11_types/ts_from_tskit.cc
     src/fwdpy11_types/tsrecorders.cc
-    src/fwdpy11_types/RecordNothing.cc)
+    src/fwdpy11_types/RecordNothing.cc
+    src/fwdpy11_types/GeneticMapUnit.cc)
 
 set(REGION_SOURCES src/regions/init.cc src/regions/Region.cc src/regions/Sregion.cc src/regions/GammaS.cc src/regions/ConstantS.cc
     src/regions/ExpS.cc src/regions/UniformS.cc src/regions/GaussianS.cc src/regions/MutationRegions.cc

fwdpy11/_types/model_params.py

@@ -238,6 +238,13 @@ def validate_recregions(self, attribute, value):
                     attr.validators.instance_of(fwdpy11.GeneticMapUnit)(
                         self, attribute, i
                     )
+                if not all([i.discrete for i in value]) and not all(
+                    [not i.discrete for i in value]
+                ):
+                    warnings.warn(
+                        "genetic map has a mix of discrete=True and discrete=False"
+                    )
+
             except TypeError:
                 raise
 

fwdpy11/genetic_map_unit.py

@@ -17,6 +17,8 @@
 # along with fwdpy11.  If not, see <http://www.gnu.org/licenses/>.
 #
 
+import typing
+
 import attr
 import numpy as np
 
@@ -28,6 +30,11 @@
 _common_attr_attribs = {"frozen": True, "auto_attribs": True, "repr_ns": "fwdpy11"}
 
 
+def _is_integer_if_discrete(self, attribute, value):
+    if self.discrete is True:
+        attr.validators.instance_of(int)(self, attribute, value)
+
+
 @attr_add_asblack
 @attr_class_pickle_with_super
 @attr_class_to_from_dict
@@ -41,11 +48,14 @@ class PoissonInterval(fwdpy11._fwdpy11._ll_PoissonInterval):
     also determine the order of positional arguments:
 
     :param beg: The beginning of the region
-    :type beg: float
+    :type beg: int or float
     :param end: The end of the region
-    :type end: float
+    :type end: int or float
     :param mean: The mean number of breakpoints per meiosis
     :type mean: float
+    :param discrete: If ``False``, positions are continuous and uniform from ``[beg, end)``.
+                     If ``True``, positions take integer values uniformly from ``[beg, end)``.
+    :type discrete: bool
 
     .. versionadded:: 0.3.0
 
@@ -61,14 +71,21 @@ class PoissonInterval(fwdpy11._fwdpy11._ll_PoissonInterval):
 
         Refactored to use attrs and inherit from
         low-level C++ class
+
+    .. versionchanged:: 0.12.0
+
+        Added ``discrete`` option to initializer.
     """
 
-    beg: float
-    end: float
+    beg: typing.Union[int, float] = attr.ib(validator=_is_integer_if_discrete)
+    end: typing.Union[int, float] = attr.ib(validator=_is_integer_if_discrete)
     mean: float
+    discrete: bool = attr.ib(kw_only=True, default=False)
 
     def __attrs_post_init__(self):
-        super(PoissonInterval, self).__init__(self.beg, self.end, self.mean)
+        super(PoissonInterval, self).__init__(
+            beg=self.beg, end=self.end, mean=self.mean, discrete=self.discrete
+        )
 
 
 @attr_add_asblack
@@ -85,9 +102,12 @@ class PoissonPoint(fwdpy11._fwdpy11._ll_PoissonPoint):
     also determine the order of positional arguments:
 
     :param position: The position of the crossover
-    :type position: float
+    :type position: int or float
     :param mean: The mean number of breakpoints per meiosis
     :type mean: float
+    :param discrete: If ``False``, positions are continuous and uniform from ``[beg, end)``.
+                     If ``True``, positions take integer values uniformly from ``[beg, end)``.
+    :type discrete: bool
 
     .. versionadded:: 0.3.0
 
@@ -103,13 +123,20 @@ class PoissonPoint(fwdpy11._fwdpy11._ll_PoissonPoint):
 
         Refactored to use attrs and inherit from
         low-level C++ class
+
+    .. versionchanged:: 0.12.0
+
+        Added ``discrete`` option to initializer.
     """
 
-    position: float
+    position: typing.Union[int, float] = attr.ib(validator=_is_integer_if_discrete)
     mean: float
+    discrete: bool = attr.ib(kw_only=True, default=False)
 
     def __attrs_post_init__(self):
-        super(PoissonPoint, self).__init__(self.position, self.mean)
+        super(PoissonPoint, self).__init__(
+            position=self.position, mean=self.mean, discrete=self.discrete
+        )
 
     def __getstate__(self):
         return self.asdict()
@@ -131,11 +158,14 @@ class BinomialInterval(fwdpy11._fwdpy11._ll_BinomialInterval):
     also determine the order of positional arguments:
 
     :param beg: The beginning of the region
-    :type beg: float
+    :type beg: int or float
     :param end: The end of the region
-    :type end: float
+    :type end: int or float
     :param probability: The probability of a recombination (per meiosis).
     :type probability: float
+    :param discrete: If ``False``, positions are continuous and uniform from ``[beg, end)``.
+                     If ``True``, positions take integer values uniformly from ``[beg, end)``.
+    :type discrete: bool
 
     .. versionadded:: 0.5.2
 
@@ -147,14 +177,24 @@ class BinomialInterval(fwdpy11._fwdpy11._ll_BinomialInterval):
 
         Refactored to use attrs and inherit from
         low-level C++ class
+
+    .. versionchanged:: 0.12.0
+
+        Added ``discrete`` option to initializer.
     """
 
-    beg: float
-    end: float
+    beg: typing.Union[int, float] = attr.ib(validator=_is_integer_if_discrete)
+    end: typing.Union[int, float] = attr.ib(validator=_is_integer_if_discrete)
     probability: float
+    discrete: bool = attr.ib(kw_only=True, default=False)
 
     def __attrs_post_init__(self):
-        super(BinomialInterval, self).__init__(self.beg, self.end, self.probability)
+        super(BinomialInterval, self).__init__(
+            beg=self.beg,
+            end=self.end,
+            probability=self.probability,
+            discrete=self.discrete,
+        )
 
 
 @attr_add_asblack
@@ -172,9 +212,12 @@ class BinomialPoint(fwdpy11._fwdpy11._ll_BinomialPoint):
     also determine the order of positional arguments:
 
     :param position: The beginning of the region
-    :type position: float
+    :type position: int or float
     :param probability: The probability of a recombination (per meiosis).
     :type probability: float
+    :param discrete: If ``False``, positions are continuous and uniform from ``[beg, end)``.
+                     If ``True``, positions take integer values uniformly from ``[beg, end)``.
+    :type discrete: bool
 
     .. versionadded:: 0.3.0
 
@@ -190,13 +233,20 @@ class BinomialPoint(fwdpy11._fwdpy11._ll_BinomialPoint):
 
         Refactored to use attrs and inherit from
         low-level C++ class
+
+    .. versionchanged:: 0.12.0
+
+        Added ``discrete`` option to initializer.
     """
 
-    position: float
-    probability: float
+    position: typing.Union[int, float] = attr.ib(validator=_is_integer_if_discrete)
+    probability: typing.Union[int, float]
+    discrete: bool = attr.ib(kw_only=True, default=False)
 
     def __attrs_post_init__(self):
-        super(BinomialPoint, self).__init__(self.position, self.probability)
+        super(BinomialPoint, self).__init__(
+            position=self.position, probability=self.probability, discrete=self.discrete
+        )
 
 
 @attr_add_asblack
@@ -212,11 +262,14 @@ class FixedCrossovers(fwdpy11._fwdpy11._ll_FixedCrossovers):
     also determine the order of positional arguments:
 
     :param beg: The beginning of the region
-    :type beg: float
+    :type beg: int or float
     :param end: The end of the region
-    :type end: float
+    :type end: int or float
     :param num_xovers: The number of breakpoints per meiosis
     :type num_xovers: float
+    :param discrete: If ``False``, positions are continuous and uniform from ``[beg, end)``.
+                     If ``True``, positions take integer values uniformly from ``[beg, end)``.
+    :type discrete: bool
 
     .. versionadded:: 0.3.0
 
@@ -232,11 +285,21 @@ class FixedCrossovers(fwdpy11._fwdpy11._ll_FixedCrossovers):
 
         Refactored to use attrs and inherit from
         low-level C++ class
+
+    .. versionchanged:: 0.12.0
+
+        Added ``discrete`` option to initializer.
     """
 
-    beg: float
-    end: float
+    beg: typing.Union[int, float] = attr.ib(validator=_is_integer_if_discrete)
+    end: typing.Union[int, float] = attr.ib(validator=_is_integer_if_discrete)
     num_xovers: int
+    discrete: bool = attr.ib(kw_only=True, default=False)
 
     def __attrs_post_init__(self):
-        super(FixedCrossovers, self).__init__(self.beg, self.end, self.num_xovers)
+        super(FixedCrossovers, self).__init__(
+            beg=self.beg,
+            end=self.end,
+            num_xovers=self.num_xovers,
+            discrete=self.discrete,
+        )