jonathf / chaospy

Chaospy - Toolbox for performing uncertainty quantification.
https://chaospy.readthedocs.io/
MIT License
437 stars 87 forks source link

Can't calculate moments with Trunc and TruncNormal distributions #371

Open spinjet opened 2 years ago

spinjet commented 2 years ago

I'm using chaospy version 4.3.4 installed with pip, on Python 3.8.12 I can't calculate the Expected value with the function cp.E(), as it throws the error UnsupportedFeature.

Here's the input for minimal testing.

import chaospy as cp

dist = cp.TruncNormal(upper=2, mu=1, sigma=0.5)
cp.E(dist)

With the following output:

Traceback (most recent call last):
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\distributions\baseclass\distribution.py", line 618, in mom
    out = [self._get_mom(kdata) for kdata in K.T]
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\distributions\baseclass\distribution.py", line 618, in <listcomp>
    out = [self._get_mom(kdata) for kdata in K.T]
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\distributions\baseclass\distribution.py", line 641, in _get_mom
    ret_val = float(self._mom(kdata, **parameters))
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\distributions\baseclass\shift_scale.py", line 108, in _mom
    out = sum([dist._get_mom(key)*coeff
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\distributions\baseclass\shift_scale.py", line 108, in <listcomp>
    out = sum([dist._get_mom(key)*coeff
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\distributions\baseclass\distribution.py", line 641, in _get_mom
    ret_val = float(self._mom(kdata, **parameters))
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\distributions\baseclass\simple.py", line 109, in _mom
    raise chaospy.UnsupportedFeature(
chaospy.UnsupportedFeature: trunc_normal(lower=-inf, upper=2, mu=1, sigma=0.5): does not support analytical raw moments.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\descriptives\expected.py", line 42, in E
    moments = dist.mom(poly.exponents.T, **kws)
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\distributions\baseclass\distribution.py", line 625, in mom
    out = [chaospy.approximate_moment(self, kdata) for kdata in K.T]
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\distributions\baseclass\distribution.py", line 625, in <listcomp>
    out = [chaospy.approximate_moment(self, kdata) for kdata in K.T]
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\distributions\approximation.py", line 206, in approximate_moment
    MOMENTS_QUADS[distribution, order] = chaospy.generate_quadrature(
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\quadrature\frontend.py", line 152, in generate_quadrature
    abscissas, weights = _generate_quadrature(
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\quadrature\frontend.py", line 248, in _generate_quadrature
    abscissas, weights = quad_function(order, dist, **parameters)
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\quadrature\clenshaw_curtis.py", line 67, in clenshaw_curtis
    return hypercube_quadrature(
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\quadrature\hypercube.py", line 85, in hypercube_quadrature
    return quad_func(**kwargs)
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\quadrature\hypercube.py", line 286, in distribution_to_domain
    abscissas, weights = quad_func(lower=lower, upper=upper, **kwargs)
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\quadrature\hypercube.py", line 237, in univariate_to_multivariate
    abscissas, weights = quad_func(**sizable)
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\quadrature\hypercube.py", line 405, in scale_quadrature
    abscissas, weights = quad_func(order=order, **kwargs)
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\quadrature\hypercube.py", line 342, in split_into_segments
    return quad_func(order=order, **kwargs)
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\quadrature\hypercube.py", line 171, in ensure_output
    abscissas, weights = quad_func(**kwargs)
  File "C:\Users\<username>\Anaconda3\envs\<envname>\lib\site-packages\chaospy\quadrature\clenshaw_curtis.py", line 107, in clenshaw_curtis_simple
    assert numpy.isclose(numpy.sum(weights), 1)
AssertionError

The same error is raised if I pass a distribution in the cp.Trunc() distribution to limit its range.

Many thanks for your help in advance.

jonathf commented 2 years ago

I'm not able to replicate the bug on my system. Mind posting your current version of numpoly and numpy by doing pip freeze?

Also if you have only tried to install through conda, do you mind trying out with pip?

spinjet commented 2 years ago
numpoly==1.2.3
numpy==1.21.3

from the pip freeze output

And yes, I had already installed chaospy using pip, conda usually has outdated packages.

spinjet commented 2 years ago

May I add I also tested the same minimal code by installing a fresh Python 3.10.2 distribution (from the official website, so not from Anaconda), and installed chaospy with pip.

Python 3.10.2 (tags/v3.10.2:a58ebcc, Jan 17 2022, 14:12:15) [MSC v.1929 64 bit (AMD64)] on win32
Type "help", "copyright", "credits" or "license" for more information.
>>> import chaospy as cp
>>> dist = cp.TruncNormal(upper=2, mu=1, sigma=0.5)
>>> cp.E(dist)

Traceback (most recent call last):
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 618, in mom
    out = [self._get_mom(kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 618, in <listcomp>
    out = [self._get_mom(kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 641, in _get_mom
    ret_val = float(self._mom(kdata, **parameters))
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\shift_scale.py", line 108, in _mom
    out = sum([dist._get_mom(key)*coeff
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\shift_scale.py", line 108, in <listcomp>
    out = sum([dist._get_mom(key)*coeff
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 641, in _get_mom
    ret_val = float(self._mom(kdata, **parameters))
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\simple.py", line 109, in _mom
    raise chaospy.UnsupportedFeature(
chaospy.UnsupportedFeature: trunc_normal(lower=-inf, upper=2, mu=1, sigma=0.5): does not support analytical raw moments.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\descriptives\expected.py", line 42, in E
    moments = dist.mom(poly.exponents.T, **kws)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 625, in mom
    out = [chaospy.approximate_moment(self, kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 625, in <listcomp>
    out = [chaospy.approximate_moment(self, kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\approximation.py", line 206, in approximate_moment
    MOMENTS_QUADS[distribution, order] = chaospy.generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 152, in generate_quadrature
    abscissas, weights = _generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 248, in _generate_quadrature
    abscissas, weights = quad_function(order, dist, **parameters)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\clenshaw_curtis.py", line 67, in clenshaw_curtis
    return hypercube_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 85, in hypercube_quadrature
    return quad_func(**kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 286, in distribution_to_domain
    abscissas, weights = quad_func(lower=lower, upper=upper, **kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 237, in univariate_to_multivariate
    abscissas, weights = quad_func(**sizable)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 405, in scale_quadrature
    abscissas, weights = quad_func(order=order, **kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 342, in split_into_segments
    return quad_func(order=order, **kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 171, in ensure_output
    abscissas, weights = quad_func(**kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\clenshaw_curtis.py", line 107, in clenshaw_curtis_simple
    assert numpy.isclose(numpy.sum(weights), 1)
AssertionError

pip freeze returns these installed packages:

chaospy==4.3.4
numpoly==1.2.3
numpy==1.22.2
scipy==1.7.3
jonathf commented 2 years ago

My best guess is that on your windows machine the machine precision is lower than on mine.

I've made the assert more lenient now in version 4.3.5. Place check it out to see if it solves you problem.

spinjet commented 2 years ago

Hi, did you push the update through pip? It still installs the 4.3.4 version and not the updated one.

jonathf commented 2 years ago

Hm... the building of the docs depend on the existence of indexes for external docs to numpy, scipy etc. Looks like scipy reference guide is down and its index with it, making the pipeline fail.

Hopefully scipy docs will be up soon. In the mean time I've uploaded 4.3.5 manually now.

spinjet commented 2 years ago

Hi jonathf

So I successfully installed the updated version and quickly re-tested the minimal code using the Python REPL. Sadly it seems the error is still present. I tested TruncNormal() with the truncation and without, same error.

Here's the output:

Python 3.10.2 (tags/v3.10.2:a58ebcc, Jan 17 2022, 14:12:15) [MSC v.1929 64 bit (AMD64)] on win32
Type "help", "copyright", "credits" or "license" for more information.
>>> import chaospy as cp
>>> cp.__version__
'4.3.5'
>>> dist = cp.TruncNormal(lower=-0.5, upper=0.5, mu=0, sigma=1)
>>> cp.E(dist)
Traceback (most recent call last):
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 618, in mom
    out = [self._get_mom(kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 618, in <listcomp>
    out = [self._get_mom(kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 641, in _get_mom
    ret_val = float(self._mom(kdata, **parameters))
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\shift_scale.py", line 108, in _mom
    out = sum([dist._get_mom(key)*coeff
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\shift_scale.py", line 108, in <listcomp>
    out = sum([dist._get_mom(key)*coeff
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 641, in _get_mom
    ret_val = float(self._mom(kdata, **parameters))
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\simple.py", line 109, in _mom
    raise chaospy.UnsupportedFeature(
chaospy.UnsupportedFeature: trunc_normal(lower=-0.5, upper=0.5, mu=0, sigma=1): does not support analytical raw moments.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\descriptives\expected.py", line 42, in E
    moments = dist.mom(poly.exponents.T, **kws)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 625, in mom
    out = [chaospy.approximate_moment(self, kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 625, in <listcomp>
    out = [chaospy.approximate_moment(self, kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\approximation.py", line 206, in approximate_moment
    MOMENTS_QUADS[distribution, order] = chaospy.generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 152, in generate_quadrature
    abscissas, weights = _generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 248, in _generate_quadrature
    abscissas, weights = quad_function(order, dist, **parameters)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\clenshaw_curtis.py", line 67, in clenshaw_curtis
    return hypercube_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 85, in hypercube_quadrature
    return quad_func(**kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 286, in distribution_to_domain
    abscissas, weights = quad_func(lower=lower, upper=upper, **kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 237, in univariate_to_multivariate
    abscissas, weights = quad_func(**sizable)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 407, in scale_quadrature
    assert numpy.sum(weights) <= 1+1e-10
AssertionError
>>> dist2 = cp.TruncNormal(mu=0, sigma=1)
>>> cp.E(dist2)
Traceback (most recent call last):
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 618, in mom
    out = [self._get_mom(kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 618, in <listcomp>
    out = [self._get_mom(kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 641, in _get_mom
    ret_val = float(self._mom(kdata, **parameters))
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\shift_scale.py", line 108, in _mom
    out = sum([dist._get_mom(key)*coeff
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\shift_scale.py", line 108, in <listcomp>
    out = sum([dist._get_mom(key)*coeff
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 641, in _get_mom
    ret_val = float(self._mom(kdata, **parameters))
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\simple.py", line 109, in _mom
    raise chaospy.UnsupportedFeature(
chaospy.UnsupportedFeature: trunc_normal(lower=-inf, upper=inf, mu=0, sigma=1): does not support analytical raw moments.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\descriptives\expected.py", line 42, in E
    moments = dist.mom(poly.exponents.T, **kws)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 625, in mom
    out = [chaospy.approximate_moment(self, kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 625, in <listcomp>
    out = [chaospy.approximate_moment(self, kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\approximation.py", line 206, in approximate_moment
    MOMENTS_QUADS[distribution, order] = chaospy.generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 152, in generate_quadrature
    abscissas, weights = _generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 248, in _generate_quadrature
    abscissas, weights = quad_function(order, dist, **parameters)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\clenshaw_curtis.py", line 67, in clenshaw_curtis
    return hypercube_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 85, in hypercube_quadrature
    return quad_func(**kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 286, in distribution_to_domain
    abscissas, weights = quad_func(lower=lower, upper=upper, **kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 237, in univariate_to_multivariate
    abscissas, weights = quad_func(**sizable)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 407, in scale_quadrature
    assert numpy.sum(weights) <= 1+1e-10
AssertionError
>>>

I took a look into the traceback and added a print statement for numpy.sum(weights). Here's what I found:

>>> import chaospy as cp
>>> dist = cp.TruncNormal(mu=0, sigma=1, lower=-1.5, upper=1.5)
>>> cp.E(dist)
1.000019999120172
Traceback (most recent call last):
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 618, in mom
    out = [self._get_mom(kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 618, in <listcomp>
    out = [self._get_mom(kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 641, in _get_mom
    ret_val = float(self._mom(kdata, **parameters))
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\shift_scale.py", line 108, in _mom
    out = sum([dist._get_mom(key)*coeff
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\shift_scale.py", line 108, in <listcomp>
    out = sum([dist._get_mom(key)*coeff
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 641, in _get_mom
    ret_val = float(self._mom(kdata, **parameters))
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\simple.py", line 109, in _mom
    raise chaospy.UnsupportedFeature(
chaospy.UnsupportedFeature: trunc_normal(lower=-1.5, upper=1.5, mu=0, sigma=1): does not support analytical raw moments.

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\descriptives\expected.py", line 42, in E
    moments = dist.mom(poly.exponents.T, **kws)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 625, in mom
    out = [chaospy.approximate_moment(self, kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\baseclass\distribution.py", line 625, in <listcomp>
    out = [chaospy.approximate_moment(self, kdata) for kdata in K.T]
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\approximation.py", line 206, in approximate_moment
    MOMENTS_QUADS[distribution, order] = chaospy.generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 152, in generate_quadrature
    abscissas, weights = _generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 248, in _generate_quadrature
    abscissas, weights = quad_function(order, dist, **parameters)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\clenshaw_curtis.py", line 67, in clenshaw_curtis
    return hypercube_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 85, in hypercube_quadrature
    return quad_func(**kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 286, in distribution_to_domain
    abscissas, weights = quad_func(lower=lower, upper=upper, **kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 237, in univariate_to_multivariate
    abscissas, weights = quad_func(**sizable)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 408, in scale_quadrature
    assert numpy.sum(weights) <= 1+1e-10
AssertionError

Since it's using the cp.approximate_moment() function, I tested it with different quadrature orders. It seems the assertion does not trigger with smaller orders.

>>> cp.approximate_moment(dist, (1,), order=1000)
0.9999999999999999
5.551115123125783e-17
>>> cp.approximate_moment(dist, (1,), order=10000)
1.0
5.551115123125783e-17
>>> cp.approximate_moment(dist, (1,), order=100000)
1.000019999120172
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\approximation.py", line 206, in approximate_moment
    MOMENTS_QUADS[distribution, order] = chaospy.generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 152, in generate_quadrature
    abscissas, weights = _generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 248, in _generate_quadrature
    abscissas, weights = quad_function(order, dist, **parameters)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\clenshaw_curtis.py", line 67, in clenshaw_curtis
    return hypercube_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 85, in hypercube_quadrature
    return quad_func(**kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 286, in distribution_to_domain
    abscissas, weights = quad_func(lower=lower, upper=upper, **kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 237, in univariate_to_multivariate
    abscissas, weights = quad_func(**sizable)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 408, in scale_quadrature
    assert numpy.sum(weights) <= 1+1e-10
AssertionError

You might want to test this behavior more thoroughly. Perhaps also add in the code some rounding if the order is incredibly small (in this case it was e-17, practically zero as it should be).

Finally, I also tested passing in Trunc() a random distribution (in this case a Weibull) with the respective approximated mean from sampling:

>>> w = cp.Weibull(shape=2)
>>> dist = cp.Trunc(w, lower=0.25, upper=1)
>>> cp.approximate_moment(dist, (1,),order=1000)
0.9999999999999999
0.6454739602268855
>>> cp.approximate_moment(dist, (1,),order=10000)
1.0
0.6454739602268855
>>> cp.approximate_moment(dist, (1,),order=100000)
1.000019999120172
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\approximation.py", line 206, in approximate_moment
    MOMENTS_QUADS[distribution, order] = chaospy.generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 152, in generate_quadrature
    abscissas, weights = _generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 248, in _generate_quadrature
    abscissas, weights = quad_function(order, dist, **parameters)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\clenshaw_curtis.py", line 67, in clenshaw_curtis
    return hypercube_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 85, in hypercube_quadrature
    return quad_func(**kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 286, in distribution_to_domain
    abscissas, weights = quad_func(lower=lower, upper=upper, **kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 237, in univariate_to_multivariate
    abscissas, weights = quad_func(**sizable)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 408, in scale_quadrature
    assert numpy.sum(weights) <= 1+1e-10
AssertionError
>>> cp.approximate_moment(dist, (1,),order=1000000)
1.000029999370875
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\distributions\approximation.py", line 206, in approximate_moment
    MOMENTS_QUADS[distribution, order] = chaospy.generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 152, in generate_quadrature
    abscissas, weights = _generate_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\frontend.py", line 248, in _generate_quadrature
    abscissas, weights = quad_function(order, dist, **parameters)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\clenshaw_curtis.py", line 67, in clenshaw_curtis
    return hypercube_quadrature(
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 85, in hypercube_quadrature
    return quad_func(**kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 286, in distribution_to_domain
    abscissas, weights = quad_func(lower=lower, upper=upper, **kwargs)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 237, in univariate_to_multivariate
    abscissas, weights = quad_func(**sizable)
  File "C:\Users\<username>\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.10_qbz5n2kfra8p0\LocalCache\local-packages\Python310\site-packages\chaospy\quadrature\hypercube.py", line 408, in scale_quadrature
    assert numpy.sum(weights) <= 1+1e-10
AssertionError
>>> samples = dist.sample(1e5)
>>> samples.mean()
0.6586637468320664
>>> dist.sample(1e6, rule='latin_hypercube').mean()
0.6587808374697979
>>> dist.sample(1e6, rule='sobol').mean()
0.6587800828647302
jonathf commented 2 years ago

Still challanging for me to diagnose.

Could you perhaps try using python -O to disable the safety checks? Does it produce resultss similar to the Monte Carlo results?

spinjet commented 2 years ago

Still challanging for me to diagnose.

Could you perhaps try using python -O to disable the safety checks? Does it produce resultss similar to the Monte Carlo results?

Here's the full output. I also used scipy.stats.truncnorm to compare the moments that were calculated. The -O flag prevented any errors apparently.

PS C:\Users\username> python -O
Python 3.10.2 (tags/v3.10.2:a58ebcc, Jan 17 2022, 14:12:15) [MSC v.1929 64 bit (AMD64)] on win32
Type "help", "copyright", "credits" or "license" for more information.
>>> import chaospy as cp
>>> from scipy.stats import truncnorm
>>> cp.__version__
'4.3.5'
>>> myclip_a, myclip_b, my_mean, my_std = -0.5, 0.5, 0, 1
>>> tn_cp = cp.TruncNormal(lower=myclip_a, upper=myclip_b, mu=my_mean, sigma=my_std)
>>> a, b = (myclip_a - my_mean) / my_std, (myclip_b - my_mean) / my_std
>>> tn_sp = truncnorm(a,b)
>>> cp.E(tn_cp)
array(2.32561366e-17)
>>> cp.Std(tn_cp)
0.28389061510439345
>>> cp.Var(tn_cp)
array(0.08059388)
>>> tn_sp.mean()
0.0
>>> tn_sp.std()
0.2838822900443274
>>> tn_sp.var()
0.08058915460081162
>>> tn_cp.sample(1e5).mean()
-0.0008270271748199187
>>> tn_cp.sample(1e5, rule='sobol').mean()
-6.042816383210119e-06
>>> tn_cp.sample(1e5, rule='latin_hypercube').mean()
-9.949751367432924e-10
>>> tn_cp.sample(1e5).std()
0.2848249303005903
>>> tn_cp.sample(1e5, rule='latin_hypercube').std()
0.28388227838058516
>>>

The mean is very close to the true mean. You may want to check how scipy calculates the moments, as they seem an analytical result and not approximate.

Here's some calculation of the error, both from the scipy values and a few Monte-Carlo runs.

>>> abs(tn_sp.mean() - cp.E(tn_cp))
2.325613659981407e-17
>>> abs(tn_sp.var() - cp.Var(tn_cp))
4.72674353924063e-06
>>> abs(tn_sp.std() - cp.Std(tn_cp))
8.325060066038947e-06
>>> abs(tn_cp.sample(1e6, rule='latin_hypercube').mean() - cp.E(tn_cp))
7.780774649919844e-11
>>> abs(tn_cp.sample(1e6, rule='latin_hypercube').mean() - cp.E(tn_cp))
1.758722960732602e-10
>>> abs(tn_cp.sample(1e6, rule='latin_hypercube').mean() - cp.E(tn_cp))
2.354790927511612e-10
>>> abs(tn_cp.sample(1e6, rule='latin_hypercube').mean() - cp.E(tn_cp))
2.839611536308677e-10
>>> abs(tn_cp.sample(1e6, rule='latin_hypercube').mean() - cp.E(tn_cp))
5.504369003564609e-10
>>> abs(tn_cp.sample(1e6, rule='latin_hypercube').mean() - cp.E(tn_cp))
1.0836413074784626e-10
>>> abs(tn_cp.sample(1e6, rule='latin_hypercube').std() - cp.Std(tn_cp))
8.324944964777092e-06
>>> abs(tn_cp.sample(1e6, rule='latin_hypercube').std() - cp.Std(tn_cp))
8.32488275115395e-06
>>> abs(tn_cp.sample(1e6, rule='latin_hypercube').std() - cp.Std(tn_cp))
8.325408364151521e-06
>>> abs(tn_cp.sample(1e6, rule='latin_hypercube').std() - cp.Std(tn_cp))
8.324838066897655e-06
>>> abs(tn_cp.sample(1e6, rule='latin_hypercube').std() - cp.Std(tn_cp))
8.325140821885402e-06
>>>
jonathf commented 2 years ago

That is an excellent idea. Once upon a time I did that for many of the distributions, but scipy did not have that particular distribution long ago and I have not kept updated on the progression of scipy.

I've just released 4.3.6 with a fix. Let me know if it solves the problem.