Bugfix: rounding of discrete distribution (#305)

* Bugfix: rounding of discrete distribution

Rounding at 0.5 is a touchy subject as the rules differs between
implementation. The wrong assumption makes the edge cases for discrete
distributions a bit problematic as NN.5 numbers are used to construct
piecewise linear distributions. To fix this at the limits a nudge is
added to keep the rounding correct.

* Bugfix: mixing of discrete/continuous vars

CHANGELOG.rst

@@ -1,6 +1,12 @@
 Master Branch
 =============
 
+Version 4.2.1 (2020-11-24)
+==========================
+
+CHANGED:
+  * Bugfix in rounding for discrete distributions.
+
 Version 4.2.0 (2020-11-23)
 ==========================
 

README.rst

@@ -61,6 +61,7 @@ you find fit. The only requirement is that the output is compatible with
 .. code-block:: python
 
     >>> coordinates = numpy.linspace(0, 10, 100)
+
     >>> def forward_solver(coordinates, parameters):
     ...     """Function to do uncertainty quantification on."""
     ...     param_init, param_rate = parameters
@@ -72,8 +73,8 @@ probability distribution. For example:
 
 .. code-block:: python
 
-    >>> distribution = chaospy.J(
-    ...     chaospy.Uniform(1, 2), chaospy.Normal(0, 2))
+    >>> distribution = chaospy.J(chaospy.Uniform(1, 2), chaospy.Normal(0, 2))
+
     >>> print(distribution)
     J(Uniform(lower=1, upper=2), Normal(mu=0, sigma=2))
 
@@ -83,6 +84,7 @@ polynomials. These can be automatically constructed:
 .. code-block:: python
 
     >>> expansion = chaospy.generate_expansion(8, distribution)
+
     >>> print(expansion[:5].round(8))
     [1.0 q1 q0-1.5 q0*q1-1.5*q1 q0**2-3.0*q0+2.16666667]
 
@@ -95,6 +97,7 @@ low-discrepancy sequences. For example to create Sobol sequence samples:
 .. code-block:: python
 
     >>> samples = distribution.sample(1000, rule="sobol")
+
     >>> print(samples[:, :4].round(8))
     [[ 1.5         1.75        1.25        1.375     ]
      [ 0.         -1.3489795   1.3489795  -0.63727873]]
@@ -103,8 +106,9 @@ We can evaluating the forward solver using these samples:
 
 .. code-block:: python
 
-    >>> evaluations = numpy.array([
-    ...     forward_solver(coordinates, sample) for sample in samples.T])
+    >>> evaluations = numpy.array([forward_solver(coordinates, sample)
+    ...                            for sample in samples.T])
+
     >>> print(evaluations[:3, :5].round(8))
     [[1.5        1.5        1.5        1.5        1.5       ]
      [1.75       2.00546578 2.29822457 2.63372042 3.0181921 ]
@@ -116,8 +120,8 @@ of ``forward_solver``:
 
 .. code-block:: python
 
-    >>> approx_solver = chaospy.fit_regression(
-    ...     expansion, samples, evaluations)
+    >>> approx_solver = chaospy.fit_regression(expansion, samples, evaluations)
+
     >>> print(approx_solver[:2].round(4))
     [q0 -0.0002*q0*q1**3+0.0051*q0*q1**2-0.101*q0*q1+q0]
 
@@ -127,9 +131,10 @@ directly. For example:
 .. code-block:: python
 
     >>> expected = chaospy.E(approx_solver, distribution)
+    >>> deviation = chaospy.Std(approx_solver, distribution)
+
     >>> print(expected[:5].round(8))
     [1.5        1.53092356 1.62757217 1.80240142 2.07915608]
-    >>> deviation = chaospy.Std(approx_solver, distribution)
     >>> print(deviation[:5].round(8))
     [0.28867513 0.43364958 0.76501802 1.27106355 2.07110879]
 

chaospy/distributions/baseclass/distribution.py

@@ -441,29 +441,21 @@ def sample(self, size=(), rule="random", antithetic=None, include_axis_dim=False
 
         Changing the sampling scheme, use the following ``rule`` flag:
 
-        +----------------------+-----------------------------------------------+
-        | key                  | Description                                   |
-        +======================+===============================================+
-        | ``chebyshev``        | Roots of first order Chebyshev polynomials.   |
-        +----------------------+-----------------------------------------------+
-        | ``nested_chebyshev`` | Chebyshev nodes adjusted to ensure nested.    |
-        +----------------------+-----------------------------------------------+
-        | ``korobov``          | Korobov lattice.                              |
-        +----------------------+-----------------------------------------------+
-        | ``random``           | Classical (Pseudo-)Random samples.            |
-        +----------------------+-----------------------------------------------+
-        | ``grid``             | Regular spaced grid.                          |
-        +----------------------+-----------------------------------------------+
-        | ``nested_grid``      | Nested regular spaced grid.                   |
-        +----------------------+-----------------------------------------------+
-        | ``latin_hypercube``  | Latin hypercube samples.                      |
-        +----------------------+-----------------------------------------------+
-        | ``sobol``            | Sobol low-discrepancy sequence.               |
-        +----------------------+-----------------------------------------------+
-        | ``halton``           | Halton low-discrepancy sequence.              |
-        +----------------------+-----------------------------------------------+
-        | ``hammersley``       | Hammersley low-discrepancy sequence.          |
-        +----------------------+-----------------------------------------------+
+        ----------------------  -------------------------------------------
+        key                     description
+        ----------------------  -------------------------------------------
+        ``additive_recursion``  Modulus of golden ratio samples.
+        ``chebyshev``           Roots of first order Chebyshev polynomials.
+        ``grid``                Regular spaced grid.
+        ``halton``              Halton low-discrepancy sequence.
+        ``hammersley``          Hammersley low-discrepancy sequence.
+        ``korobov``             Korobov lattice.
+        ``latin_hypercube``     Latin hypercube samples.
+        ``nested_chebyshev``    Chebyshev nodes adjusted to ensure nested.
+        ``nested_grid``         Nested regular spaced grid.
+        ``random``              Classical (Pseudo-)Random samples.
+        ``sobol``               Sobol low-discrepancy sequence.
+        ----------------------  -------------------------------------------
 
         All samples are created on the ``[0, 1]``-hypercube, which then is
         mapped into the domain of the distribution using the inverse Rosenblatt
@@ -499,6 +491,9 @@ def sample(self, size=(), rule="random", antithetic=None, include_axis_dim=False
         out = sampler.generator.generate_samples(
             order=size_, domain=self, rule=rule, antithetic=antithetic)
 
+        for idx, dist in enumerate(self):
+            if dist.interpret_as_integer:
+                out[idx] = numpy.round(out[idx])
         if self.interpret_as_integer:
             out = numpy.round(out).astype(int)
         out = out.reshape(shape)

chaospy/distributions/collection/binomial.py

@@ -41,10 +41,10 @@ def _pdf(self, x_data, size, prob):
         return special.comb(size, x_data)*prob**x_data*(1-prob)**(size-x_data)
 
     def _lower(self, size, prob):
-        return -0.5
+        return -0.5+1e-10
 
     def _upper(self, size, prob):
-        return numpy.round(size)+0.5
+        return numpy.round(size)+0.5-1e-10
 
     def _mom(self, k_data, size, prob):
         x_data = numpy.arange(int(size)+1, dtype=int)
@@ -76,7 +76,7 @@ class Binomial(J):
         >>> distribution.inv(uloc).round(2)
         array([-0.5 ,  0.55,  0.93,  1.31,  1.69,  2.07,  2.45,  3.5 ])
         >>> distribution.sample(10)
-        array([2, 1, 0, 2, 2, 2, 2, 3, 3, 0])
+        array([1, 3, 2, 1, 0, 1, 1, 0, 2, 1])
         >>> distribution.mom([1, 2, 3]).round(4)
         array([1.5 , 3.  , 6.75])
 

chaospy/distributions/collection/discrete_uniform.py

@@ -53,11 +53,11 @@ def _cdf(self, x_data, lower, upper):
 
     def _lower(self, lower, upper):
         """Lower bounds."""
-        return numpy.round(lower)-0.5
+        return numpy.round(lower)-0.5+1e-10
 
     def _upper(self, lower, upper):
         """Upper bounds."""
-        return numpy.round(upper)+0.5
+        return numpy.round(upper)+0.5-1e-10
 
     def _pdf(self, x_data, lower, upper):
         """Probability density function."""

pyproject.toml

@@ -4,7 +4,7 @@ build-backend = "poetry.masonry.api"
 
 [tool.poetry]
 name = "chaospy"
-version = "4.2.0"
+version = "4.2.1"
 description = "Numerical tool for perfroming uncertainty quantification"
 license = "MIT"
 authors = ["Jonathan Feinberg"]