SMC ABC Inference doesn't support 1D parameter space

[BryanRumsey]

If this is an intended feature we should raise an error early.

Error

Statistics
Compute Environment: Local
2023-04-07 19:12:35,614 - GillesPy2 - ERROR - Job errors: tuple index out of range
Traceback (most recent call last):
File "/stochss/stochss/handlers/util/scripts/start_job.py", line 120, in
job.run(verbose=args.verbose)
File "/stochss/stochss/handlers/util/model_inference.py", line 943, in run
results = smc_abc_inference.infer(*infer_args, **infer_kwargs)
File "/opt/conda/lib/python3.8/site-packages/sciope/inference/smc_abc.py", line 210, in infer
abc_results = abc_instance.infer(num_samples=t,
File "/opt/conda/lib/python3.8/site-packages/sciope/inference/abc_inference.py", line 354, in infer
return self.rejection_sampling(num_samples, batch_size, chunk_size, ensemble_size, normalize)
File "/opt/conda/lib/python3.8/site-packages/sciope/inference/abc_inference.py", line 210, in rejection_sampling
for f, dist in as_completed(futures_dist, with_results=True):
File "/opt/conda/lib/python3.8/site-packages/distributed/client.py", line 5199, in __next__
return self._get_and_raise()
File "/opt/conda/lib/python3.8/site-packages/distributed/client.py", line 5188, in _get_and_raise
raise exc.with_traceback(tb)
File "/opt/conda/lib/python3.8/site-packages/sciope/inference/smc_abc.py", line 72, in _weighted_draw_perturb
if self.ref_prior.pdf(sz) > 0:
File "/opt/conda/lib/python3.8/site-packages/sciope/utilities/priors/uniform_prior.py", line 82, in pdf
for i in range(z.shape[0]):
IndexError: tuple index out of range

Model

def create_michaelis_menten(parameter_values=None):
    # Initialize Model
    model = gillespy2.Model(name="Michaelis_Menten")

    # Define Variables (GillesPy2.Species)
    A = gillespy2.Species(name='Substrate', initial_value=301)
    B = gillespy2.Species(name='Enzyme', initial_value=120)
    C = gillespy2.Species(name='Enzyme_Substrate_Complex', initial_value=0)
    D = gillespy2.Species(name='Product', initial_value=0)
    
    # Add Variables to Model
    model.add_species([A, B, C, D])

    # Define Parameters
    rate1 = gillespy2.Parameter(name='rate1', expression=0.0017)
    rate2 = gillespy2.Parameter(name='rate2', expression=0.5)
    rate3 = gillespy2.Parameter(name='rate3', expression=0.1)
    
    # Add Parameters to Model
    model.add_parameter([rate1, rate2, rate3])

    # Define Reactions
    r1 = gillespy2.Reaction(
        name="r1", reactants={'Substrate': 1, 'Enzyme': 1}, products={'Enzyme_Substrate_Complex': 1}, rate='rate1'
    )
    r2 = gillespy2.Reaction(
        name="r2", reactants={'Enzyme_Substrate_Complex': 1}, products={'Substrate': 1, 'Enzyme': 1}, rate='rate2'
    )
    r3 = gillespy2.Reaction(
        name="r3", reactants={'Enzyme_Substrate_Complex': 1}, products={'Enzyme': 1, 'Product': 1}, rate='rate3'
    )
    
    # Add Reactions to Model
    model.add_reaction([r1, r2, r3])
    
    # Define Timespan
    tspan = gillespy2.TimeSpan.linspace(t=100, num_points=101)
    
    # Set Model Timespan
    model.timespan(tspan)
    return model

Inference

def configure_simulation(trajectories=1):
    solver = gillespy2.SSACSolver(model=model, delete_directory=False)
    kwargs = {
        'number_of_trajectories': trajectories,
        'seed': None,
        'solver': solver,
    }
    return kwargs

kwargs = configure_simulation()
results = model.run(**kwargs)

unshaped_obs_data = results.to_array()
shaped_obs_data = unshaped_obs_data.swapaxes(1, 2)
obs_data = shaped_obs_data[:,1:, :]
print(obs_data.shape)

def process(raw_results):
    return raw_results.to_array().swapaxes(1, 2)[:,1:, :]

def simulator(parameter_point):
    model = create_michaelis_menten()

    labels = [
        'rate1'
    ]
    for ndx, parameter in enumerate(parameter_point):
        model.listOfParameters[labels[ndx]].expression = str(parameter)

    kwargs = configure_simulation()
    raw_results = model.run(**kwargs)

    return process(raw_results)

values = numpy.array([
    0.0017
])
parameter_names = [
    'rate1'
]

lower_bounds = numpy.array([
    0.00085
])
upper_bounds = numpy.array([
    0.00255
])

uni_prior = uniform_prior.UniformPrior(lower_bounds, upper_bounds)

summ_func = auto_tsfresh.SummariesTSFRESH()

eps_selector = RelativeEpsilonSelector(20, max_rounds=2)

c = Client()

smc_abc_inference = smc_abc.SMCABC(
    obs_data, sim=simulator, prior_function=uni_prior, summaries_function=summ_func.compute
)

smc_abc_results = smc_abc_inference.infer(
    num_samples=100, batch_size=10, eps_selector=eps_selector, chunk_size=10
)