Skip to content

Broadcasting with metropolis sampler leads to issues where NUTS works perfectly. #4491

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom
Closed
ExpectationMax opened this issue Feb 27, 2021 · 1 comment
Closed

Broadcasting with metropolis sampler leads to issues where NUTS works perfectly. #4491

ExpectationMax opened this issue Feb 27, 2021 · 1 comment

Comments

@ExpectationMax
Copy link
Contributor

Using broadcasting together with the Metropolis sampler often leads to problems as indicated in these previous issues:

These are still present even though the linked issues are already closed.
Interesting these problems only occur when using metropolis and not when using NUTS (see examples below).

Minimal examples

Successful with NUTS sampler:

import theano
import pymc3 as pm

theano.config.mode = "FAST_COMPILE"
theano.config.exception_verbosity = "high"


def main():
    n_sample = 5
    n_seq = 13

    with pm.Model() as simple_model:
        N = pm.Uniform('N', lower=0, upper=100000, shape=(1, n_seq))
        p = pm.Beta('p', 1, 1, shape=(n_sample, 1))
        pm.Normal('dummy', mu=N, sigma=p, shape=(n_sample, n_seq))
        idata = pm.sample(200, tune=100, return_inferencedata=True)


if __name__ == '__main__':
    main()

Failure with metropolis:

import theano
import pymc3 as pm

theano.config.mode = "FAST_COMPILE"
theano.config.exception_verbosity = "high"

def main():
    n_sample = 5
    n_seq = 13

    with pm.Model() as simple_model:
        N = pm.Uniform('N', lower=0, upper=100000, shape=(1, n_seq))
        p = pm.Beta('p', 1, 1, shape=(n_sample, 1))
        pm.Normal('dummy', mu=N, sigma=p, shape=(n_sample, n_seq))
        step = pm.Metropolis()
        idata = pm.sample(200, step=step, tune=100, return_inferencedata=True)

if __name__ == '__main__':
    main()

Traceback.

Only 200 samples in chain.
Multiprocess sampling (2 chains in 2 jobs)
CompoundStep
>Metropolis: [dummy]
>Metropolis: [p]
>Metropolis: [N]
pymc3.parallel_sampling.RemoteTraceback: -----------------------------------------------------------| 0.00% [0/600 00:00<00:00 Sampling 2 chains, 0 divergences]
"""
Traceback (most recent call last):
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/theano/link/vm.py", line 313, in __call__
    thunk()
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/theano/graph/op.py", line 476, in rval
    r = p(n, [x[0] for x in i], o)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/theano/tensor/elemwise.py", line 776, in perform
    raise ValueError(base_exc_str)
ValueError: Dimension mismatch; shapes are (5, 13), (1, 13)

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/parallel_sampling.py", line 137, in run
    self._start_loop()
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/parallel_sampling.py", line 191, in _start_loop
    point, stats = self._compute_point()
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/parallel_sampling.py", line 216, in _compute_point
    point, stats = self._step_method.step(self._point)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/step_methods/compound.py", line 42, in step
    point, state = method.step(point)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/step_methods/arraystep.py", line 196, in step
    apoint, stats = self.astep(self.bij.map(point))
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/step_methods/metropolis.py", line 217, in astep
    accept = self.delta_logp(q, q0)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/theano/compile/function/types.py", line 975, in __call__
    if output_subset is None
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/theano/link/vm.py", line 317, in __call__
    raise_with_op(self.fgraph, node, thunk)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/theano/link/utils.py", line 508, in raise_with_op
    raise exc_value.with_traceback(exc_trace)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/theano/link/vm.py", line 313, in __call__
    thunk()
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/theano/graph/op.py", line 476, in rval
    r = p(n, [x[0] for x in i], o)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/theano/tensor/elemwise.py", line 776, in perform
    raise ValueError(base_exc_str)
ValueError: Dimension mismatch; shapes are (5, 13), (1, 13)
Apply node that caused the error: Elemwise{sub,no_inplace}(Reshape{2}.0, Elemwise{add,no_inplace}.0)
Toposort index: 17
Inputs types: [TensorType(float64, matrix), TensorType(float64, matrix)]
Inputs shapes: [(5, 13), (1, 13)]
Inputs strides: [(104, 8), (104, 8)]
Inputs values: ['not shown', 'not shown']
Inputs type_num: [12, 12]
Outputs clients: [[Elemwise{pow,no_inplace}(Elemwise{sub,no_inplace}.0, TensorConstant{(1, 1) of 2})]]

Backtrace when the node is created(use Theano flag traceback__limit=N to make it longer):
  File "test.py", line 22, in <module>
    main()
  File "test.py", line 17, in main
    step = pm.Metropolis()
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/step_methods/arraystep.py", line 86, in __new__
    step.__init__([var], *args, **kwargs)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/step_methods/metropolis.py", line 187, in __init__
    self.delta_logp = delta_logp(model.logpt, vars, shared)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/step_methods/metropolis.py", line 883, in delta_logp
    [logp0], inarray0 = pm.join_nonshared_inputs([logp], vars, shared)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/theanof.py", line 278, in join_nonshared_inputs
    xs_special = [theano.clone(x, replace, strict=False) for x in xs]
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/theanof.py", line 278, in <listcomp>
    xs_special = [theano.clone(x, replace, strict=False) for x in xs]

Debugprint of the apply node:
Elemwise{sub,no_inplace} [id A] <TensorType(float64, matrix)> ''
 |Reshape{2} [id B] <TensorType(float64, matrix)> ''
 | |Subtensor{int64:int64:} [id C] <TensorType(float64, vector)> ''
 | | |inarray [id D] <TensorType(float64, vector)>
 | | |Constant{0} [id E] <int64>
 | | |Constant{65} [id F] <int64>
 | |TensorConstant{[ 5 13]} [id G] <TensorType(int64, vector)>
 |Elemwise{add,no_inplace} [id H] <TensorType(float64, matrix)> ''
   |Elemwise{mul,no_inplace} [id I] <TensorType(float64, matrix)> ''
   | |sigmoid [id J] <TensorType(float64, matrix)> ''
   | | |N_interval___shared [id K] <TensorType(float64, matrix)>
   | |TensorConstant{(1, 1) of 100000.0} [id L] <TensorType(float64, (True, True))>
   |Elemwise{mul,no_inplace} [id M] <TensorType(float64, matrix)> ''
     |Elemwise{sub,no_inplace} [id N] <TensorType(float64, matrix)> ''
     | |TensorConstant{(1, 1) of 1} [id O] <TensorType(int8, (True, True))>
     | |sigmoid [id J] <TensorType(float64, matrix)> ''
     |TensorConstant{(1, 1) of 0.0} [id P] <TensorType(float64, (True, True))>

HINT: Use the Theano flag 'exception_verbosity=high' for a debugprint and storage map footprint of this apply node.
"""

The above exception was the direct cause of the following exception:

ValueError: Dimension mismatch; shapes are (5, 13), (1, 13)
Apply node that caused the error: Elemwise{sub,no_inplace}(Reshape{2}.0, Elemwise{add,no_inplace}.0)
Toposort index: 17
Inputs types: [TensorType(float64, matrix), TensorType(float64, matrix)]
Inputs shapes: [(5, 13), (1, 13)]
Inputs strides: [(104, 8), (104, 8)]
Inputs values: ['not shown', 'not shown']
Inputs type_num: [12, 12]
Outputs clients: [[Elemwise{pow,no_inplace}(Elemwise{sub,no_inplace}.0, TensorConstant{(1, 1) of 2})]]

Backtrace when the node is created(use Theano flag traceback__limit=N to make it longer):
  File "test.py", line 22, in <module>
    main()
  File "test.py", line 17, in main
    step = pm.Metropolis()
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/step_methods/arraystep.py", line 86, in __new__
    step.__init__([var], *args, **kwargs)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/step_methods/metropolis.py", line 187, in __init__
    self.delta_logp = delta_logp(model.logpt, vars, shared)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/step_methods/metropolis.py", line 883, in delta_logp
    [logp0], inarray0 = pm.join_nonshared_inputs([logp], vars, shared)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/theanof.py", line 278, in join_nonshared_inputs
    xs_special = [theano.clone(x, replace, strict=False) for x in xs]
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/theanof.py", line 278, in <listcomp>
    xs_special = [theano.clone(x, replace, strict=False) for x in xs]

Debugprint of the apply node:
Elemwise{sub,no_inplace} [id A] <TensorType(float64, matrix)> ''
 |Reshape{2} [id B] <TensorType(float64, matrix)> ''
 | |Subtensor{int64:int64:} [id C] <TensorType(float64, vector)> ''
 | | |inarray [id D] <TensorType(float64, vector)>
 | | |Constant{0} [id E] <int64>
 | | |Constant{65} [id F] <int64>
 | |TensorConstant{[ 5 13]} [id G] <TensorType(int64, vector)>
 |Elemwise{add,no_inplace} [id H] <TensorType(float64, matrix)> ''
   |Elemwise{mul,no_inplace} [id I] <TensorType(float64, matrix)> ''
   | |sigmoid [id J] <TensorType(float64, matrix)> ''
   | | |N_interval___shared [id K] <TensorType(float64, matrix)>
   | |TensorConstant{(1, 1) of 100000.0} [id L] <TensorType(float64, (True, True))>
   |Elemwise{mul,no_inplace} [id M] <TensorType(float64, matrix)> ''
     |Elemwise{sub,no_inplace} [id N] <TensorType(float64, matrix)> ''
     | |TensorConstant{(1, 1) of 1} [id O] <TensorType(int8, (True, True))>
     | |sigmoid [id J] <TensorType(float64, matrix)> ''
     |TensorConstant{(1, 1) of 0.0} [id P] <TensorType(float64, (True, True))>

HINT: Use the Theano flag 'exception_verbosity=high' for a debugprint and storage map footprint of this apply node.

The above exception was the direct cause of the following exception:

Traceback (most recent call last):
  File "test.py", line 22, in <module>
    main()
  File "test.py", line 18, in main
    idata = pm.sample(200, step=step, tune=100, return_inferencedata=True)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/sampling.py", line 558, in sample
    trace = _mp_sample(**sample_args, **parallel_args)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/sampling.py", line 1476, in _mp_sample
    for draw in sampler:
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/parallel_sampling.py", line 479, in __iter__
    draw = ProcessAdapter.recv_draw(self._active)
  File "/Users/hornm/Library/Caches/pypoetry/virtualenvs/grass-collaboration-i2H8sdro-py3.7/lib/python3.7/site-packages/pymc3/parallel_sampling.py", line 359, in recv_draw
    raise error from old_error
RuntimeError: Chain 0 failed.

Potential fix

As was indicated in #3337, this might well be due to loss of broadcasting information when calling make_shared_replacements. Theanos error message and traceback then indicate that the source of the problem is in pm.join_nonshared_inputs, which is slightly misleading. The error occurs here as this is the location where the original variable ist replaced with the shared version of itself where broadcasting information is missing.

I will open a pull request to fix this, I do not anticipate any negative consequences for other parts of pymc3.

Versions and main components

  • PyMC3 Version: 3.11.1
  • Theano Version: 1.1.2
  • Python Version: 3.7.9
  • Operating system: MacOS Catalina 10.15.7
  • How did you install PyMC3: pip (to be exact via poetry)
@ExpectationMax ExpectationMax mentioned this issue Feb 27, 2021
Merged
3 tasks
@michaelosthege
Copy link
Member

closed by #4492

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Assignees
No one assigned
Labels
None yet
Projects
None yet
Milestone
No milestone
Development

No branches or pull requests
2 participants
@ExpectationMax @michaelosthege