Need help with handling ODE model. Error: unhashable type: 'slice'

Questions version agnostic
,
1

Hello,

I am following the example ODE Lotka-Volterra With Bayesian Inference in Multiple Ways and trying to apply it with my ODE model.

The problem I have is that my rhs (right hand side function defining the system of ODEs) requires an input of 17 species concentration, C and it outputs 17 ODEs (or my dC/dt). I only have observed data for 6 of the 17 species therefore I only need to take 6 concentration profiles or 6 ode solutions.

I tried doing the following:

@as_op(itypes=[pt.dvector], otypes=[pt.dmatrix])
def pytensor_forward_model_matrix(theta):
    answer = solve_ivp(fun=rhs, t_span = [0, 1080], y0=y0_data[0], t_eval=time, args=(theta,), method='Radau')
    return answer[:,[5, 9, 11, 6, 12, 13]]

where [5,9,11,6,12,13] are the column indices of the 6 species.

I got the following error:

Traceback (most recent call last):
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pytensor/compile/function/types.py", line 970, in __call__
    self.vm()
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pytensor/graph/op.py", line 543, in rval
    r = p(n, [x[0] for x in i], o)
        ^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pytensor/compile/ops.py", line 259, in perform
    outs = self.__fn(*inputs)
           ^^^^^^^^^^^^^^^^^^
  File "/home/admin/pymc-projects/sunovion_sampler.py", line 291, in pytensor_forward_model_matrix
    return answer[:,[5, 9, 11, 6, 12, 13]]
           ~~~~~~^^^^^^^^^^^^^^^^^^^^^^^^^
TypeError: unhashable type: 'slice'

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/home/admin/pymc-projects/sunovion_sampler.py", line 329, in <module>
    trace_DEMZ = pm.sample(step=[pm.DEMetropolisZ(vars_list)], tune=tune, draws=draws, cores=1)
                 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pymc/sampling/mcmc.py", line 608, in sample
    model.check_start_vals(ip)
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pymc/model.py", line 1776, in check_start_vals
    initial_eval = self.point_logps(point=elem)
                   ^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pymc/model.py", line 1810, in point_logps
    self.compile_fn(factor_logps_fn)(point),
    ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pymc/pytensorf.py", line 764, in __call__
    return self.f(**state)
           ^^^^^^^^^^^^^^^
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pytensor/compile/function/types.py", line 983, in __call__
    raise_with_op(
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pytensor/link/utils.py", line 535, in raise_with_op
    raise exc_value.with_traceback(exc_trace)
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pytensor/compile/function/types.py", line 970, in __call__
    self.vm()
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pytensor/graph/op.py", line 543, in rval
    r = p(n, [x[0] for x in i], o)
        ^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/admin/anaconda3/envs/pymc_env/lib/python3.11/site-packages/pytensor/compile/ops.py", line 259, in perform
    outs = self.__fn(*inputs)
           ^^^^^^^^^^^^^^^^^^
  File "/home/admin/pymc-projects/sunovion_sampler.py", line 291, in pytensor_forward_model_matrix
    return answer[:,[5, 9, 11, 6, 12, 13]]
           ~~~~~~^^^^^^^^^^^^^^^^^^^^^^^^^
TypeError: unhashable type: 'slice'
Apply node that caused the error: FromFunctionOp{pytensor_forward_model_matrix}(MakeVector{dtype='float64'}.0)
Toposort index: 22
Inputs types: [TensorType(float64, (17,))]
Inputs shapes: [(17,)]
Inputs strides: [(8,)]
Inputs values: ['not shown']
Outputs clients: [[Elemwise{Composite}(Y_obs{[[4.200321..9168e-01]]}, FromFunctionOp{pytensor_forward_model_matrix}.0, InplaceDimShuffle{x,x}.0, TensorConstant{(1, 1) of -0.5}, TensorConstant{(1, 1) of ..5332046727}, Elemwise{log,no_inplace}.0, Elemwise{gt,no_inplace}.0, TensorConstant{(1, 1) of -inf})]]

Backtrace when the node is created (use PyTensor flag traceback__limit=N to make it longer):
  File "/home/admin/pymc-projects/sunovion_sampler.py", line 316, in <module>
    ode_solution = pytensor_forward_model_matrix(

HINT: Use the PyTensor flag `exception_verbosity=high` for a debug print-out and storage map footprint of this Apply node.

Any help would be greatly appreciated. Thanks.

2

Could you slice before handing the ode_solution to the likelihood, rather than inside the pytensor_forward_model_matrix function? For example:

    # ... priors ...
    # Ode solution function
    ode_solution = pytensor_forward_model_matrix(theta)

    # Likelihood
    pm.Normal("Y_obs", mu=ode_solution[:, [5, 9, 11, 6, 12, 13]], sigma=sigma, observed=data.values)
3

Hi. Thank you for your reply.

Your suggestion was actually the first method I did. I receive a “Segmentation fault” error when I do that.

4

A segfault is very bad, and definitely shouldn’t be happening when you slice a matrix. Some thoughts:

  1. Can you post a whole model (with some data, even artificial) that reproduces the error?

  2. Also, did you install PyMC using the suggested method in the wiki?

  3. Can you execute this code? If so, is the shape/type of obs as expected?

    obs = pm.Normal("Y_obs", mu=ode_solution, sigma=sigma)
    print(obs.eval())
  1. Can you use an indicator matrix to “slice” ode_solution, i.e.:
Z = np.zeros((6, 17))
Z[np.arange(6), [5, 9, 11, 6, 12, 13]] = 1
mu = Z @ ode_solution
obs = pm.Normal("Y_obs", mu=mu, sigma=sigma, observed=data.values)
  1. Can you solve your model using the least_squares method described in the notebook you linked without issue?
5

Hi,

I appreciate the help.

  1. This is the whole model with artificial data:
    toshare.py (7.3 KB)
  2. Yes. Here’s my package list:
    package_list.txt (13.4 KB)
  3. I can’t seem to execute the code.
  4. This also returns errors.
  5. Yes, the model works using least_squares and I use the least squares result as an initial value for the sampling just like in the notebook.
6

It seems the problem is that solve_ivp, unlike odeint, returns a results object instead of a numpy array. You need your pytensor_forward_model_matrix function to return an array, so you can make the following change:

@as_op(itypes=[pt.dvector], otypes=[pt.dmatrix])
def pytensor_forward_model_matrix(theta):
    answer = solve_ivp(fun=rhs, t_span = [0, 1080], y0=y0_data, t_eval=time, args=(theta,), method='Radau')
    return answer.y

In general when you hit an error, make sure you check the type and shape of everything in your model. The .eval() method is very useful here. I uncovered the problem by doing ode_solution.eval() and seeing that it was not an array.

1 Like