From d3d947e497cd592de406afbc1ac22a74a3c2bc9c Mon Sep 17 00:00:00 2001 From: Alexander Rodin Date: Sun, 21 May 2017 16:52:15 +0300 Subject: [PATCH 1/7] Make spacing more compliant to PEP8 --- pymc3/distributions/dist_math.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/pymc3/distributions/dist_math.py b/pymc3/distributions/dist_math.py index 074b459335..def866eeb1 100644 --- a/pymc3/distributions/dist_math.py +++ b/pymc3/distributions/dist_math.py @@ -365,6 +365,7 @@ def conjugate_solve_triangular(outer, inner): grad = tt.triu(s + s.T) - tt.diag(tt.diagonal(s)) return [tt.switch(ok, grad, floatX(np.nan))] + class SplineWrapper (theano.Op): """ Creates a theano operation from scipy.interpolate.UnivariateSpline @@ -381,6 +382,7 @@ def perform(self, node, inputs, output_storage): x, = inputs output_storage[0][0] = np.asarray(self.spline(x)) + class DifferentiableSplineWrapper (SplineWrapper): """ Creates a theano operation with defined gradient from From 4d510f4569129d5cc07f3fd81b7b1d084d93459c Mon Sep 17 00:00:00 2001 From: Alexander Rodin Date: Sun, 21 May 2017 16:56:18 +0300 Subject: [PATCH 2/7] Raise NotImplementedError from SplineWrapper grad operation Fixes #2209 --- pymc3/distributions/dist_math.py | 3 +++ 1 file changed, 3 insertions(+) diff --git a/pymc3/distributions/dist_math.py b/pymc3/distributions/dist_math.py index def866eeb1..4648979496 100644 --- a/pymc3/distributions/dist_math.py +++ b/pymc3/distributions/dist_math.py @@ -382,6 +382,9 @@ def perform(self, node, inputs, output_storage): x, = inputs output_storage[0][0] = np.asarray(self.spline(x)) + def grad(self, inputs, grads): + raise NotImplementedError + class DifferentiableSplineWrapper (SplineWrapper): """ From a85bd3e0b7bf7b565154f38f1853142a5463aeaf Mon Sep 17 00:00:00 2001 From: Alexander Rodin Date: Mon, 22 May 2017 16:22:07 +0300 Subject: [PATCH 3/7] Instantiate SplineWrapper recursively to simplify the architecture --- pymc3/distributions/continuous.py | 4 ++-- pymc3/distributions/dist_math.py | 28 ++++++++++------------------ 2 files changed, 12 insertions(+), 20 deletions(-) diff --git a/pymc3/distributions/continuous.py b/pymc3/distributions/continuous.py index b4d4ddec92..bbd1ea104a 100644 --- a/pymc3/distributions/continuous.py +++ b/pymc3/distributions/continuous.py @@ -16,7 +16,7 @@ from pymc3.theanof import floatX from . import transforms -from .dist_math import bound, logpow, gammaln, betaln, std_cdf, i0, i1, alltrue_elemwise, DifferentiableSplineWrapper +from .dist_math import bound, logpow, gammaln, betaln, std_cdf, i0, i1, alltrue_elemwise, SplineWrapper from .distribution import Continuous, draw_values, generate_samples, Bound __all__ = ['Uniform', 'Flat', 'Normal', 'Beta', 'Exponential', 'Laplace', @@ -1430,7 +1430,7 @@ def __init__(self, x_points, pdf_points, transform='interval', Z = interp.integral(x_points[0], x_points[-1]) self.Z = tt.as_tensor_variable(Z) - self.interp_op = DifferentiableSplineWrapper(interp) + self.interp_op = SplineWrapper(interp, n_derivatives=1) self.x_points = x_points self.pdf_points = pdf_points / Z self.cdf_points = interp.antiderivative()(x_points) / Z diff --git a/pymc3/distributions/dist_math.py b/pymc3/distributions/dist_math.py index 4648979496..ae4d516514 100644 --- a/pymc3/distributions/dist_math.py +++ b/pymc3/distributions/dist_math.py @@ -371,33 +371,25 @@ class SplineWrapper (theano.Op): Creates a theano operation from scipy.interpolate.UnivariateSpline """ - __props__ = ('spline',) + __props__ = ('spline', 'grad_op') itypes = [tt.dscalar] otypes = [tt.dscalar] - def __init__(self, spline): + def __init__(self, spline, n_derivatives): self.spline = spline + if n_derivatives > 0: + self.grad_op = SplineWrapper(spline.derivative(), n_derivatives - 1) + else: + self.grad_op = None + def perform(self, node, inputs, output_storage): x, = inputs output_storage[0][0] = np.asarray(self.spline(x)) def grad(self, inputs, grads): - raise NotImplementedError - - -class DifferentiableSplineWrapper (SplineWrapper): - """ - Creates a theano operation with defined gradient from - scipy.interpolate.UnivariateSpline - """ - - def __init__(self, spline): - super(DifferentiableSplineWrapper, self).__init__(spline) - self.spline_grad = SplineWrapper(spline.derivative()) - self.__props__ += ('spline_grad',) - - def grad(self, inputs, grads): + if self.grad_op is None: + raise NotImplementedError x, = inputs x_grad, = grads - return [x_grad * self.spline_grad(x)] + return [x_grad * self.grad_op(x)] From 532235aa7b47d00a06e84bb06bf2a838cf38213c Mon Sep 17 00:00:00 2001 From: Alexander Rodin Date: Mon, 22 May 2017 16:22:07 +0300 Subject: [PATCH 4/7] Create spline derivatives lazily --- pymc3/distributions/continuous.py | 2 +- pymc3/distributions/dist_math.py | 23 +++++++++++++---------- 2 files changed, 14 insertions(+), 11 deletions(-) diff --git a/pymc3/distributions/continuous.py b/pymc3/distributions/continuous.py index bbd1ea104a..49a2246f59 100644 --- a/pymc3/distributions/continuous.py +++ b/pymc3/distributions/continuous.py @@ -1430,7 +1430,7 @@ def __init__(self, x_points, pdf_points, transform='interval', Z = interp.integral(x_points[0], x_points[-1]) self.Z = tt.as_tensor_variable(Z) - self.interp_op = SplineWrapper(interp, n_derivatives=1) + self.interp_op = SplineWrapper(interp) self.x_points = x_points self.pdf_points = pdf_points / Z self.cdf_points = interp.antiderivative()(x_points) / Z diff --git a/pymc3/distributions/dist_math.py b/pymc3/distributions/dist_math.py index ae4d516514..d0fdf48a66 100644 --- a/pymc3/distributions/dist_math.py +++ b/pymc3/distributions/dist_math.py @@ -371,25 +371,28 @@ class SplineWrapper (theano.Op): Creates a theano operation from scipy.interpolate.UnivariateSpline """ - __props__ = ('spline', 'grad_op') + __props__ = ('spline',) itypes = [tt.dscalar] otypes = [tt.dscalar] - def __init__(self, spline, n_derivatives): + def __init__(self, spline): self.spline = spline - if n_derivatives > 0: - self.grad_op = SplineWrapper(spline.derivative(), n_derivatives - 1) - else: - self.grad_op = None - def perform(self, node, inputs, output_storage): x, = inputs output_storage[0][0] = np.asarray(self.spline(x)) def grad(self, inputs, grads): - if self.grad_op is None: - raise NotImplementedError x, = inputs x_grad, = grads - return [x_grad * self.grad_op(x)] + + if not hasattr(self, 'grad_op'): + try: + self.grad_op = SplineWrapper(self.spline.derivative()) + except ValueError: + self.grad_op = None + + if self.grad_op is None: + raise NotImplementedError('Spline of order 0 is not differentiable') + else: + return [x_grad * self.grad_op(x)] From 8fab09536ced1026224c9de3c24aae9307e3d0ea Mon Sep 17 00:00:00 2001 From: Alexander Rodin Date: Sat, 27 May 2017 09:34:20 +0300 Subject: [PATCH 5/7] Move grad_op to a separate property --- pymc3/distributions/dist_math.py | 23 +++++++++++++---------- 1 file changed, 13 insertions(+), 10 deletions(-) diff --git a/pymc3/distributions/dist_math.py b/pymc3/distributions/dist_math.py index d0fdf48a66..a56005dc9f 100644 --- a/pymc3/distributions/dist_math.py +++ b/pymc3/distributions/dist_math.py @@ -378,6 +378,18 @@ class SplineWrapper (theano.Op): def __init__(self, spline): self.spline = spline + @property + def grad_op(self): + if not hasattr(self, '_grad_op'): + try: + self._grad_op = SplineWrapper(self.spline.derivative()) + except ValueError: + self._grad_op = None + + if self._grad_op is None: + raise NotImplementedError('Spline of order 0 is not differentiable') + return self._grad_op + def perform(self, node, inputs, output_storage): x, = inputs output_storage[0][0] = np.asarray(self.spline(x)) @@ -386,13 +398,4 @@ def grad(self, inputs, grads): x, = inputs x_grad, = grads - if not hasattr(self, 'grad_op'): - try: - self.grad_op = SplineWrapper(self.spline.derivative()) - except ValueError: - self.grad_op = None - - if self.grad_op is None: - raise NotImplementedError('Spline of order 0 is not differentiable') - else: - return [x_grad * self.grad_op(x)] + return [x_grad * self.grad_op(x)] From b31d6f3b6dfc22b07de3d1bb63451161933fd75c Mon Sep 17 00:00:00 2001 From: Alexander Rodin Date: Sat, 27 May 2017 16:46:46 +0300 Subject: [PATCH 6/7] Add tests for SplineWrapper --- pymc3/tests/test_dist_math.py | 22 ++++++++++++++++++++-- 1 file changed, 20 insertions(+), 2 deletions(-) diff --git a/pymc3/tests/test_dist_math.py b/pymc3/tests/test_dist_math.py index 277e70556d..eeb4731317 100644 --- a/pymc3/tests/test_dist_math.py +++ b/pymc3/tests/test_dist_math.py @@ -4,13 +4,13 @@ import theano import theano.tests.unittest_tools as utt import pymc3 as pm -from scipy import stats +from scipy import stats, interpolate import pytest from ..theanof import floatX from ..distributions import Discrete from ..distributions.dist_math import ( - bound, factln, alltrue_scalar, MvNormalLogp) + bound, factln, alltrue_scalar, MvNormalLogp, SplineWrapper) def test_bound(): @@ -176,3 +176,21 @@ def test_hessian(self): logp = MvNormalLogp()(cov, delta) g_cov, g_delta = tt.grad(logp, [cov, delta]) tt.grad(g_delta.sum() + g_cov.sum(), [delta, cov]) + + +class TestSplineWrapper: + + def test_grad(self): + x = np.linspace(0, 1, 100) + y = x * x + spline = SplineWrapper(interpolate.InterpolatedUnivariateSpline(x, y, k=1)) + utt.verify_grad(spline, [0.5]) + + def test_hessian(self): + x = np.linspace(0, 1, 100) + y = x * x + spline = SplineWrapper(interpolate.InterpolatedUnivariateSpline(x, y, k=1)) + x_var = tt.dscalar('x') + g_x, = tt.grad(spline(x_var), [x_var]) + with pytest.raises(NotImplementedError): + tt.grad(g_x, [x_var]) From 5da1ce2424df5cc94c33fca4eaa238cbb36a979a Mon Sep 17 00:00:00 2001 From: Alexander Rodin Date: Sat, 27 May 2017 18:07:25 +0300 Subject: [PATCH 7/7] Fix style issues --- pymc3/distributions/continuous.py | 5 ++++- pymc3/distributions/dist_math.py | 2 +- pymc3/tests/test_dist_math.py | 3 +-- 3 files changed, 6 insertions(+), 4 deletions(-) diff --git a/pymc3/distributions/continuous.py b/pymc3/distributions/continuous.py index 49a2246f59..b9d2addde1 100644 --- a/pymc3/distributions/continuous.py +++ b/pymc3/distributions/continuous.py @@ -16,7 +16,10 @@ from pymc3.theanof import floatX from . import transforms -from .dist_math import bound, logpow, gammaln, betaln, std_cdf, i0, i1, alltrue_elemwise, SplineWrapper +from .dist_math import ( + bound, logpow, gammaln, betaln, std_cdf, i0, + i1, alltrue_elemwise, SplineWrapper +) from .distribution import Continuous, draw_values, generate_samples, Bound __all__ = ['Uniform', 'Flat', 'Normal', 'Beta', 'Exponential', 'Laplace', diff --git a/pymc3/distributions/dist_math.py b/pymc3/distributions/dist_math.py index a56005dc9f..9b3a95b064 100644 --- a/pymc3/distributions/dist_math.py +++ b/pymc3/distributions/dist_math.py @@ -366,7 +366,7 @@ def conjugate_solve_triangular(outer, inner): return [tt.switch(ok, grad, floatX(np.nan))] -class SplineWrapper (theano.Op): +class SplineWrapper(theano.Op): """ Creates a theano operation from scipy.interpolate.UnivariateSpline """ diff --git a/pymc3/tests/test_dist_math.py b/pymc3/tests/test_dist_math.py index eeb4731317..ac84a5424f 100644 --- a/pymc3/tests/test_dist_math.py +++ b/pymc3/tests/test_dist_math.py @@ -178,8 +178,7 @@ def test_hessian(self): tt.grad(g_delta.sum() + g_cov.sum(), [delta, cov]) -class TestSplineWrapper: - +class TestSplineWrapper(object): def test_grad(self): x = np.linspace(0, 1, 100) y = x * x