Bayesian inversion analysis of nonlinear dynamics in surface heterogeneous reactions

It is essential to extract nonlinear dynamics from time-series data as an inverse problem in natural sciences. We propose a Bayesian statistical framework for extracting nonlinear dynamics of surface heterogeneous reactions from sparse and noisy observable data. Surface heterogeneous reactions are chemical reactions with conjugation of multiple phases, and they have the intrinsic nonlinearity of their dynamics caused by the effect of surface-area between different phases. We adapt a belief propagation method and an expectation-maximization (EM) algorithm to partial observation problem, in order to simultaneously estimate the time course of hidden variables and the kinetic parameters underlying dynamics. The proposed belief propagation method is performed by using sequential Monte Carlo algorithm in order to estimate nonlinear dynamical system. Using our proposed method, we show that the rate constants of dissolution and precipitation reactions, which are typical examples of surface heterogeneous reactions, as well as the temporal changes of solid reactants and products, were successfully estimated only from the observable temporal changes in the concentration of the dissolved intermediate product.

Figure 1

Target experimental system for surface heterogeneous reactions. A solid product $n^{\left(\mathrm{p}\right)}$ is generated from a solid reactant $n^{\left(\mathrm{r}\right)}$ via an intermediate product $C^{\left(\mathrm{i}\right)}$. The reactions $n^{\left(\mathrm{r}\right)}\to C^{\left(\mathrm{i}\right)}$ and $C^{\left(\mathrm{i}\right)}\to n^{\left(\mathrm{p}\right)}$ depend on the surface areas of the solid reactant $S_{\mathrm{r}}\left(n^{\left(\mathrm{r}\right)}\right)$ and the solid product $S_{\mathrm{p}}\left(n^{\left(\mathrm{p}\right)}\right)$, respectively. In the proposed method, the rate constants $k_{\mathrm{r}}$ and $k_{\mathrm{p}}$ underlying surface heterogeneous reactions between reactant $n^{\left(\mathrm{r}\right)}$ and product $n^{\left(\mathrm{p}\right)}$ via intermediate product $C^{\left(\mathrm{i}\right)}$ are estimated by using only the intermediate product concentration $C^{\left(\mathrm{i}\right)}$ observed every time step interval ${\mathrm{\Delta }}_{\mathrm{obs}}$.

Figure 2

Typical time course of surface heterogeneous reaction and noisy observable data. The solid product $n^{\left(\mathrm{p}\right)}$ (blue dotted line) is generated from the solid reactant $n^{\left(\mathrm{r}\right)}$ (red solid line) through the intermediate product dissolved in the liquid phase, $C^{\left(\mathrm{i}\right)}$ (green dashed line). Sparse and noisy data of the intermediate product $y$ (circles) are observed with low time-resolution.

Figure 3

Graphical model of the proposed nonlinear state-space model for surface heterogeneous reactions. Time evolution of three hidden variables ${\mathit{x}}_t=\{n_t^{\left(\mathrm{r}\right)},C_t^{\left(\mathrm{i}\right)},n_t^{\left(\mathrm{p}\right)}\}$ at time step $t$ is described by the system model $p\left({\mathit{x}}_{t+1}\right|{\mathit{x}}_t)$, whereas the observation process from hidden variables ${\mathit{x}}_t$ to observed data $y_t$ can be described by the observation model $p\left(y_t\right|{\mathit{x}}_t)$. Within three hidden variables ${\mathit{x}}_t=\{n_t^{\left(\mathrm{r}\right)},C_t^{\left(\mathrm{i}\right)},n_t^{\left(\mathrm{p}\right)}\}$, only noisy data of intermediate product, $y_t$, is assumed to be observed.

Figure 4

(a) Time courses of true hidden variables: solid reactant $n_t^{\left(\mathrm{r}\right)}$ (solid line), intermediate product $C_t^{\left(\mathrm{i}\right)}$ (dashed line), and solid product $n_t^{\left(\mathrm{p}\right)}$ (dotted line). (b)–(d) Observed data $y_t$ generated by adding Gaussian noise to intermediate product $C_t^{\left(\mathrm{i}\right)}$. The standard deviation of the observation noise is set to be $\sigma ={10}^{-2.5}$ for (b), $\sigma ={10}^{-2}$ for (c), and $\sigma ={10}^{-1.5}$ for (d), whereas the number of observation points is set to be $N_{\mathrm{obs}}=31$ for (b)–(d).

Figure 5

Estimated hidden variables (solid lines): solid reactant $n_t^{\left(\mathrm{r}\right)}$ (red line), intermediate product $C_t^{\left(\mathrm{i}\right)}$ (green line), and solid product $n_t^{\left(\mathrm{p}\right)}$ (blue line). Dashed lines show the true time courses of the hidden variables. Parameters are set to be $N_{\mathrm{obs}}=31$ and ${\sigma }_y=1.0\times {10}^{-2}$.

Figure 6

Estimation of rate constants $k_{\mathrm{r}}$ and $k_{\mathrm{p}}$ from partially observable data. Estimated values for different levels of noise $[{\sigma }_y={10}^{-2.5}$ (black), ${10}^{-2}$ (red), and ${10}^{-1.5}$ (blue)] are shown. (a), (b) The estimated rate constants (solid lines) converge to true rate constants as the iteration steps proceed, even when we set the initial values are far from the true values (dashed lines). (c) Trajectory of estimated rate constants in two-dimensional space $(k_{\mathrm{r}},k_{\mathrm{p}})$. The square shows initial values and the circle shows true values.

Figure 7

Dependence of discrepancy between true and estimated kinetic parameters $k_{\mathrm{r}}$ and $k_{\mathrm{p}}$ on the number of observation $N_{\mathrm{obs}}$ and the observation noise ${\sigma }_y$. The noisy data of intermediate products are set to be observable at every time step interval ${\mathrm{\Delta }}_{\mathrm{obs}}\sim T/N_{\mathrm{obs}}$ within total time steps $T=500$.