pints-team / pints Goto Github PK

need a function/class that takes a model, a set of parameters, prior information and a data set and returns distributions on the parameters using adaptive covariance MCMC

Following Gary's suggestion on slack

"Gradient Profiling (Giles Hooker) CollocInfer in R - a mixture of Nonlinear least squares and Gradient matching I think… need to read up"

@jonc says: do this!

Maybe look into ABC?

See:

1. Brockwell (2006) Parallel Markov chain Monte Carlo simulation by pre-fetching
2. Strid (2009) Efficient parallelisation of Metropolis-Hastings algorithms using a pre-fetching approach

Due to the large filesize, the data should either not be stored at all, or store it in a directory with a suitable .gitignore file.

This can then be used to test some optimisation methods on echem data

https://github.com/mirams/PyHillFit/blob/master/python/PyHillFit.py

Share interface with #2

Build 1 and then add the rest later

Need a function/class that takes a model, some data, a set of parameters and some bounds and gives the best-fit parameters for that data

@mirams says: R has lots of libraries. We should set up a first method that uses on of those, test how well it works

Test different filtering methods from basic thinning to more sophisticated filters.

need a test script that compares STAN versus our own adaptive covariance MCMC, using

1. Electrochem FeIII dataset
2. Kylie's dataset (or other, as I don't have this yet...)

For example, explore each parameter individually (param x, score y) or plot any two parameters against each other (param 1 x, param 2 z, score y)

Use evaluator interface to parallelise

• Uniform
• Latin hypercube
• Sobol

Martin wrote:

"[C]an you create a score function that depends on the parameters (i.e. incorporate a prior into the score function)?"

Libraries shouldn't affect other libraries, it would be pretty rubbish if people imported pints and then found matplotlib's behavior had changed because of the cma and seaborn modules.

Ideally, this would automatically update priors too, somehow...?

Two options:

1. Implement a periodic transform on the parameter space (strategy used by CMA-ES, for example)
2. Implement a transform on the score function

Need to implement each model in the STAN language. Perhaps have a member function for each model that returns the STAN language representation as a string?

Following Gary's suggestion on slack

We need

• gamma
• beta
• exponential
• student-t

And tests for each

create a cvsin_type_1 reader for electrochemistry data

Generate data and store in CSV, add different levels of gaussian noise

It would be nice to have a class that implements a Model concept, which takes a cellml file (or string?) defining what the model is.

@MichaelClerx: you have some cellml conversion routines don't you. Would this be useful here?

The current cmaes method has some unused ipop code:

Once someone figures out a good way to get random samples in the parameter space we can either add an ipop setting to the CMAES class or rename the class IPOP_CMAES and create a wrapper called CMAES that disables it

Methods like IPOP_CMAES use multiple restarts from random positions in the search space to improve chances of finding optima and reducing chances of getting stuck.

We could add some code that does this automatically, maybe using the Boundaries class to generate new starting points or perhaps a Prior class.

Following @mirams suggestion

Not the actual results, just the examples of MCMC

Gary wrote: "Another 'whilst I remember' type thing! It would be good to get boring-old-Fisher Information / Hessian at the MLE and the covariance matrix that that implies, so we could compare max likelihood with Bayesian for some of these problems. Some of our peaks are so unimodal I suspect it may be an excellent approximation for a lot of our problems, and a zillion times faster."

• CMA-ES uses a NaN returned from the score function to trigger a resampling
• PSO uses Inf, because it wants a ranking of scores, and inf is the worst possible value (note that x < NaN --> false, x > NaN --> false, etc)
• We could use exceptions instead (but I quite like returning NaN, as it might happen undetected anyway!)

@martinjrobins Jonathan Cooper suggested we set up this repo to have automated testing with Travis (travis-ci.org).
I had a look but it tells me I don't have the authority. Would you like to give this a go?

Sanmitra says it's better :-)

1. Are all the algorithms happy with this?
2. Should this only happen under the hood? I imagine users would prefer to specify a prior rather than a log prior... We could even think about giving the Prior class a log() method?

Burn-in should always happen first
But after that, stop either when max iterations is reached or when converged (if it works well, maybe set default max-iter to something very big?)

Maybe update the names of the fixed/inferred ones to make the distinction clearer

Martin wrote:

[C]an you create an efficient score function that depends on the distance between experiment and model in the frequency domain, rather than time domain? I guess the score class can just take an FFT of the values when its created and re-use this?

One simple option is repeated restarts: #29

Another option might be to use PSO to find N good starting points (by running a search with N+M, M>=0 particles and returning the best N results) and then starting a search from each of these

• Aslanidi 2009
• Clancy 2001
• Courtemanche 1998
• Di Veroli 2013 a
• Di Veroli 2013 b
• Fink 2008
• Fox 2002
• Grandi 2010
• Hund 2004
• Inada 2009
• Kurata 2002
• Lindblad 1996
• Lin 1996
• Lu 2001
• Matsuoka 2003
• Mazhari 2001
• Noble 1998
• Nygren 1998
• Oehmen 2002
• O'Hara 2011
• Priebe 1998
• Ramirez 2000
• Seemann 2003
• Severi 2012
• Shannon 2004
• Ten Tusscher 2004
• Wang 1997
• Winslow 1999
• Zeng 1995
• Zhang 2000

Skipping Kiehn 1999 because it doesn't have equations for the rate constants but a look-up table instead

Possibly based on GPy or GPflow

Need a function/class that takes a model and a data set, passes this into STAN to be solved (using HMC) and returns distributions on the parameters

I think @MichaelClerx changes to CMA-ES to get it working in the new infrastructure removed the parallel aspects of CMA-ES? This is still in there as a comment, so should just need to integrate it with the new code

Chris Gill wrote:

I had an idea last summer about how one might go about doing model selection and parameter fitting in one go using mcmc but didn’t have time to get the details working enough to share it with anyone. It turns out someone has already developed the idea in quite a general framework, and it is useful - it’s called reversible jump mcmc. Essentially you can jump between different parameter spaces provided you have a suitable map between them. The wikipedia page has a fairly short intro to it. I wondered if that might be an interesting direction to try out with the electrochemistry stuff, e.g. determining mechanism of action and the parameters in one (admittedly computationally expensive) go? Just a thought, and I’m sure it will depend on how the different reaction models are specified, but I’ve been meaning to email you about it for some time now.