AdityaSavara / PEUQSE

Parameter estimation for complex physical problems often suffers from finding ‘solutions’ that are not physically realistic. The PEUQSE software provides tools for finding physically realistic parameter estimates, graphs of the parameter estimate positions within parameter space, and plots of the final simulation results.
13 stars 5 forks source link

change TPR case that has the integral to Better Method #33

Open AdityaSavara opened 4 years ago

AdityaSavara commented 4 years ago

As described in this paper: https://www.sciencedirect.com/science/article/pii/S0039602816300632 A similar issue occurs for the posterior because the likelihood has a shape similar to SSR.

So we should make a function that takes the integral for TPR, and can use Euler method or Runge Kutta with uncertainty propagation (perhaps with uncertainties package).

Interestingly, the same idea was had by John Kitchin, independently: http://kitchingroup.cheme.cmu.edu/blog/2018/10/11/A-differentiable-ODE-integrator-for-sensitivity-analysis/ https://kitchingroup.cheme.cmu.edu/blog/category/uncertainty/ https://rosettacode.org/wiki/Runge-Kutta_method

To use this with cantera we'll have to make a custom ODE integration function rather than using the build in time stepper: https://cantera.org/examples/python/reactors/custom.py.html

AdityaSavara commented 4 years ago

Example 3 in this pull request does that using Euler's method: https://github.com/AdityaSavara/ODE-KIN-BAYES-SG-EW/commit/0b7e62b4f9a316a13b921d749bcf2c55bb472bd9 So now the next thing is to make an example 4 that does that for the Cantera case.

AdityaSavara commented 4 years ago

Eventually we should get a cantera pyrolysis example https://groups.google.com/forum/?utm_medium=email&utm_source=footer#!msg/cantera-users/ZNirpzpdsEc/ItT-QAtNCAAJ

AdityaSavara commented 4 years ago

Perhaps should copy NLREG and use Romberg Method: http://www.nlreg.com/integral.htm https://en.wikipedia.org/wiki/Romberg%27s_method There is a Romberg method in scipy, but it's not for the same purpose. https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.romberg.html

There is also numpy trapz, which could be compatible with the uncertainties module. https://github.com/lebigot/uncertainties/ https://docs.scipy.org/doc/numpy/reference/generated/numpy.trapz.html

Another practical choice is to use Simpson's rule, since it should not be too hard to propagate the uncertainty that way, and also it could still be compatible with the uncertainties module noted above. https://github.com/jkitchin/f18-06623/blob/master/lectures/02-integration-1.ipynb

AdityaSavara commented 4 years ago

In simple testing, I found that the "optimum" TPR peak fit is shifting to the right with Euler method. This means that either we need to do something more accurate (e.g., runge-kutta) and/or that we need to do an untransformed optimization again once the map is found.

AdityaSavara commented 4 years ago

Evidently, uncertainties gets installed with pycse. It would be not so bad if we make pycse an optional dependency. https://pypi.org/project/pycse/