Contributed Discussion to Bayesian Solution Uncertainty Quantification for Differential Equations by Oksana A. Chkrebtii, David A. Campbell, Ben Calderhead, and Mark A. Girolami
We begin by introducing the main ideas of the paper, and we give a brief description of the method proposed. Next, we discuss an alternative approach based on B-spline expansion, and lastly we make some comments on the method’s convergence rate.
2016 \volumeTBA \issueTBA \doi0000 \firstpage1 \lastpage2 \startlocaldefs\endlocaldefs
Differential equation \kwddiscretization uncertainty \kwdB-splines \kwdtensor product B-splines \kwdconvergence rate
I would like to congratulate the authors for such an interesting research. The Bayesian method with the probabilistic solver introduced is highly innovative and practical. The various examples presented in the paper show the wide applicability of the proposed method. However, I do find the title a bit of a misnomer, as I initially thought that the authors are constructing credible sets for the fixed but unknown solution of the differential equation.
The inverse problem that the authors are trying to solve, in its most basic form is this: Suppose you have observations , where is some normal errors and follows . Here, is a known transformation from the state space to the observation space , is the first order derivative with respect to its argument , is the known form of the differential equation, and ’s are the equation’s parameters. The method proposed consists of two steps, with one nested within the other. First, solve for probabilistically to obtain a discretized solution at some grid points. Then we embed these discretized version of in a Bayesian hierarchical framework to estimate . To model discretization uncertainty associated with using only evaluated at grid points, the authors endow priors based on Gaussian process jointly on and , where the covariance function is constructed by convolving kernels.
There is an alternative and perhaps a conceptually easier way to achieve the same result. We can first represent by a B-spline series, i.e., with denoting the th B-spline of order , and we endow the coefficients ’s with normal priors. Here, the number of basis plays the role of , where is the length-scale parameter defined in the paper. It turns out that the first derivative of this is another B-spline series where is some weighted first order finite difference of the ’s ((4.23) of Schumaker (2007)). Therefore, and are jointly normal and their associated covariance matrices are banded due the support separation property of B-splines. To enforce the given initial condition, we can condition the joint prior on .
Moreover, this approach can be generalized to the partial differential equation case, where we take tensor product of B-splines to model both the spatial and temporal components, i.e., . As in the univariate case, partial derivatives of tensor product B-splines will be another tensor-product B-splines ((3.2) of Yoo and Ghosal (2016)). Hence we will obtain the same Gaussian process prior for and all its mixed partial derivatives if we endow normal priors on the coefficients. As before, we enjoy some simplification in computing the covariance matrices because they are banded.
In addition, I would like to comment on the effect of grid point distribution on the convergence rate of the proposed algorithm. Intuitively, one would expect that the grid points should be chosen roughly uniformly across the domain . Suppose we choose grid points and we further assume that they are quasi-uniform, i.e., for some constant with . In other words, the max grid increment is of the same order as the min grid increment. Then it follows that is of the order of and by Theorem 1 of Chkrebtii et al. (2016), the rate of convergence is . Therefore for quasi-uniform grids (which includes uniform discrete grids), increasing the number of grid points will result in more accurate solution.
The paper under discussion Chkrebtii et al. (2016) makes significant contribution to the new field of probabilistic numerics. I have learnt a great deal by reading this paper, and I hope that there will be more papers in uncertainty quantification for differential equation models in the future.
- Chkrebtii et al. (2016) Chkrebtii, O. A., Campbell, D. A., Calderhead, B., and Girolami, M. A. (2016). “Bayesian solution uncertainty quantification for differential equations.” Bayesian Anal., 29 pages. Advance publication.
- Schumaker (2007) Schumaker, L. (2007). Spline Functions: Basic Theory. Cambridge University Press, New York, third edition.
- Yoo and Ghosal (2016) Yoo, W. W. and Ghosal, S. (2016). “Supremum norm posterior contraction and credible sets for nonparametric multivariate regression.” Ann. Statist., 44(3): 1069–1102.