Physics Informed Deep Kernel Learning
Authors: Zheng Wang
Wei Xing
Robert Kirby
Shandian Zhe
Venue:  (2021)
Keywords: Deep Kernel
Bayesian Learning

Abstract: Deep kernel learning is a promising combination of deep neural networks and nonparametric function estimation. However, as a data-driven approach, the performance of deep kernel learning can still be restricted by scarce or insufficient data, especially in extrapolation tasks.  To address these limitations, we propose Physics Informed Deep Kernel Learning (PI-DKL) that exploits physics  knowledge represented by differential equations with latent sources. Specifically, we use the posterior function sample of the  Gaussian process as the surrogate for the solution of the differential equation, and construct  a generative component to integrate the equation in a principled Bayesian hybrid framework. For efficient and effective inference, we marginalize out the latent variables in the joint probability and derive a simple model evidence lower bound (ELBO), based on which we develop a stochastic collapsed inference algorithm. Our ELBO can be viewed as a nice, interpretable posterior regularization objective. On synthetic datasets and real-world applications, we show the advantage of our approach in both prediction accuracy and uncertainty quantification. 

Under review as a conference paper at ICLR 2021
PHYSICS INFORMED DEEP KERNEL LEARNING
Anonymous authors
Paper under double-blind review
ABSTRACT
Deep kernel learning is a promising combination of deep neural networks and nonparametric function estimation. However, as data driven approach, the per- formance of deep kernel learning can still be restricted by scarce or insufficient data, especially in extrapolation tasks. To address these limitations, we propose Physics Informed Deep Kernel Learning (PI-DKL) that exploits physics knowl- edge represented by differential equations with latent sources. Specifically, we use the posterior function sample of the Gaussian process as the surrogate for the solution of the differential equation, and construct generative component integrate the equation in a principled Bayesian hybrid framework. For efficient and effective inference, we marginalize out the latent variables in the joint probability and derive simple model evidence lower bound (ELBO), based which we develop a stochastic collapsed inference algorithm. Our ELBO can be viewed as nice, interpretable posterior regularization objective. On synthetic datasets and real-world applications, we show the advantage of our approach in both prediction accuracy and uncertainty quantification.
Introduction
Deep kernel learning (Wilson et al., 2016a) uses deep neural networks to construct kernels for nonparametric function estimation (e.g., Gaussian processes (Williams and Rasmussen (006)) and unifies both the expressive power of neural networks and self-adaptation of nonparametric function learning. Many applications have shown that deep kernel learning substantially outperforms the conventional shallow kernel learning (e.g., RBF). Compared to standard neural networks, deep kernel learning enjoys closed-form posterior distributions and hence is more convenient for uncertainty quantification and reasoning, which is important for decision making.
Nonetheless, as data driven approach, the performance of deep kernel learning can still be restricted by scarce data, especially when the training samples are insufficient to reflect the complexity of the system (that produced the data) or the test points are far away from the training set, i.e., extrapolation. On the other hand, physies knowledge, expressed as differential equations, are used to build physical models for various science and engineering applications (Lapidus and Pinder 2011). These models meant to characterize the underlying mechanism (i.e., physical processes) that drives the system (e.g,, how the heat diffuses across the spatial and temporal domains) and are much less restricted by data availability: they can make accurate predictions even without training data, e.g., the landing of Curiosity on Mars and flight of Voyager 1.
Therefore, we consider integrating physics knowledge into deep kernel learning to further improve its performance in prediction and uncertainty quantification, especially for scarce data and extrapolation tasks. Our work enlightened by the recent Physics Informed Neural Networks (PINNs) (Raissi et al. 2019). However, there are two substantial differences. First, PINNs require the form of the differential equations be fully specified. We allow the equations include unknown latent sources (functions), which is of often the case in practice. Second, we integrate the differential equations with a principled Bayesian manner to pursue better calibrated posterior estimations.
Specifically, we use the posterior sample of the Gaussian process (GP), which is a random function, as the surrogate of the solution of the differential equation. We then apply the differential operators in the equation to obtain the sample of the latent source (function), for which we assign another GP prior. To ensure the sampling procedure is valid, we use the symmetric property of the Gaussian distribution to sample set of virtual observations {0). which is computationally equivalent to placing the GP prior with zero mean function over the latent source. The sampling procedure constitutes a generative component and ties to the original deep kernel model in the Bayesian hybrid framework (Lasserre

---

Under review as a conference paper at ICLR 2021
al. 2006). For efficient and high-quality inference, we marginalize out all the latent variables in the joint distribution to avoid approximating their complex posteriors. Then we use Jensen's 's inequality to derive simple model evidence lower bound (ELBO), based on which we develop stochastic collapsed inference algorithm. The ELBO can be further explained as a soft posterior regularization objective (Ganchev et al.] 20T0), regularized by physics.
For evaluation, we examined our physics informed deep kernel learning (PI-DKL) in both simulation and real-world applications. On synthetic datasets based on two commonly used differential equations, PI-DKL outperforms the standard deep kernel learning, shallow kernel learning, and latent force models (LFM) that combine the physics via kernel convolution, in both ground-trut function recovery and prediction uncertainty, especially in the case of extrapolation. We then examined PI-DKL in four real-world applications. PI-DKI. consistently improves upon the competing approaches prediction error and test log-likelihood. We applied PI-DKL. for nonlinear differential equation where LFM is infeasible. PI-DKL significantly outperforms standard deep/shallow kernel learning methods.
2 Background
Gaussian Process and Kernel Learning, The Gaussian process (GP) is the most commonly used nonparametric function prior for kernel learning. Suppose we aim to learn a function f Rd - I from training set D = (X,y). where X [x,.... xN [T y [y/ yN] each x,, is a d dimensional input vector and yn the observed output. To avoid muder-fitting and over-fitting, we do not want to assume any parametric form of f. Instead, we want the complexity of f-)) to automatically adapt to the data. To this end, we introduce a kernel function (=;:,) that measures the similarity of the function values in terms of their inputs. The similarity only brings in smoothness assumption about
the target function. For example, the commonly used RBF kernel, (x, xj) = exp( -xx,,x;1""
implies the function is infinitely differentiable. We then use the kernel to construct GP prior, f GP (mn(-), (m.,),((.,))) ()--)) where (.))) the mean function that is usually set to constant O. According to GP definition, the finite projection of (-) the training inputs X, namely f [f(x)..
(xxx)
follow a multivariate Gaussian distribution, Xf = (f O, K) where K is the kernel matrix on X and each [K],, x,). Given the function values f, the observed outputs y are sampled from noisy model. For example, when y are continuous, we can use the isotropic Gaussian noise model, p(y f)) ((ylffr-1l where is the inverse variance. We can then integrate out f to obtain the marginal likelihood,
ply|X) N(y|0,
To learn the model, we can maximize the likelihood to estimate the kernel parameters and the inverse variance т. According to the GP prior, given new input x", the posterior (or predictive) distribution of the output f(x") is conditional Gaussian,
(((xx)))((x)))))))) v(x))) where u(x*) k. (K + - Iy- v(x*) x(xxxx")) KK + II' and k. [x(x",, x1). ... k(x' xN)]
p(f(x") J|x", X,y)
(2)
Deep Kernel Learning. While GP priors with shallow kernels (e.g., RBF and Matern) have achieved great success in many applications, these shallow structures can limit the expressiveness of the kernels estimating highly complicated functions, e.g., sharp discontinuities and high curvatures. To address this problem, Wilson et al. (2016a) proposed to construct deep kernels with neural networks. Specifically, they first choose a shallow kernel as the base kernel. Each input is first fed into a neural network (NN), and the NN outputs are then fed into the base kernel to compute the final kernel function value. Take RBF as an example of the base kernel, we can construct a deep kernel by
kDEEP X,, x,) KREF NN(x)))
Note that the NN weights now become the kernel parameters, We can then use the deep kernel to construct GP prior for nonparametric function estimation. The likelihood and predictive distribution have the same forms as in and
3 Model
By using deep neural networks construct highly expressive kernels, deep kernel learning greatly en- hances the capability of estimating complicated functions, and meanwhile inherits the self-adaptation of the nonparametric function learning and convenient posterior inference. However, purely
2

---

Under review as a conference paper at ICLR 2021
data-drive approach, deep kernel learning can still suffer from data searcity, especially when the training examples are inadequate to reflect the complexity of the underlying mapping and when the test points are distant from all the training samples, i.e., extrapolation. To overcome this limitation, propose PI-DKL., a physics informed deep kernel learning model that exploit physics prior knowledge to improve the function learning and uncertainty reasoning. Our model is presented as follows.
3.1 Physics Informed Deep Kernel Learning
We assume that in general, the physics is described by differential equation of the following form,
+[f(x]] = gxx)
(4)
where t is a functional that combines a set of differential operators, [x) is the target (or solution) function we want to estimate from the training dataset whose form is unknown, Note that the functional x5[-] may include unknown parameters as well. One example is ([[())] df(a)
yX, and g(x) is latent source
de + a f - 3. where the input x scalar, and and 3 are unknown parameters. This functional represents a linear operator. Another commonly seen example is from the 0f(x)
viscous version of Burger's equation (Olsen Kettle 2011). ([](x)]
8r: fxx 0f(x) 0 f(x)
dag
0ug
where x [2] Ii the spatial variable, 22 the time variable, and the unknown viscosity parameter. This functional includes a nonlinear operator, f((x) 0f(x)
To incorporate the physics knowledge in (4), we propose a hybrid of conditional and generative models based on the general framework proposed by Lasserre et al. (2006). The conditional component is the standard deep-kemnll GP that given the training inputs X, samples the (noisy) output observations y. and the probability p(y X) is given in ① The generative component fulfills another GP prior over the latent source (()), but avoids the double prior problem to ensure a valid joint Bayesian model for posterior inference. Coupled with the differential operators, the generative component regularizes and guides the deep kernel learning of the f(-). The graphical illustration of PI-DKL is shown in Fig. 1 of the supplementary document.
Specifically, to consider a GP prior over g(.), we first sample a finite set of input locations Z [z...... Zm] (we will discuss the choice of p(Z) later). Then the projection of g(-) on Z follows a multivariate Gaussian distribution,
p(g Z) = N( Σ)
(5)
where g =
[ = x( 2,) and x(.,)) -) is another kernel.
Next, we link the GP model of the target ((-) to the latent source g(-) via the differential equation Our key idea is that from the GP posterior (2). we can construct sample of the target function "
f(.) = (a-) + + where e ~ [c(o,, 1), pa-)) and vv(.) are the posterior mean and standard deviation functions. While this is a random function (due to E), it has a closed form and we apply the functional y to obtain the sample of g()),
g(-) h(,cc = (h(p(-) + cv)())
(6)
Therefore, to sample g - the values of g(-) on Z. we can first sample standard Gaussian white noise c, and then sample from
p(g c, x, = II ( h(2,,)))
where 3, = 922)) and 5(-) the Dirac delta prior. Note that we can also directly view g as transformation of the Gaussian noise and derive the marginal distribution p'g/X, y) (see the discussion in the supplementary material), which, however, much more difficult to compute.
Now, we want to tie the GP prior for .(-) in (5) to the samples g generated from the GP model of the target function f()). i.e., through (7). In this way, the learning of f(-) can be guided or regularized
by the differential equation (4). However, directly multiplying (5) and (7) is problematic, because g will have double priors and the sampling procedure is invalid if g is sampled from it cannot
In computational physics, this is viewed a surrogate for the solution function of the differential equation.
3

---

We can see that placing a (finite) GP prior over g(-) is equivalent to sampling a set of virtual observations O, due to the symmetry of the Gaussian distribution. Therefore, can turn the GP prior of the latent source to a generative component that samples a set of virtual observations From the computational perspective, they totally equivalent. However, the sampling procedure Under review as a conference paper at ICLR 2021
be sampled again from 77, and vice versa. To ensure model is valid probabilistic model for posterior inference, we utilize the symmetric property of the Gaussian distribution,
p(g = Ng 0, 2)) = No >) - p(0)g, 2).
(8)
now becomes valid - we first sample g from 7 and then sample 0 from (8). Note that the virtual observations 0 come from the zero-mean function of the GP prior of g(-). We can use different virtual observations by choosing a nonzero mean function.
Finally, we combine the conditional model and the generative model (see T and 8 obtain joint probability distribution,
ply, 0, Z, e, x) = (()()))(()((pppl x, y)p(0 Z) = -, +T ())(() (Tyyo,K++*TT(2))(((()) 1) II [0p - h(zz, ()) Σ)
The choice of p(Z) is flexible. If we have knowledge about the input distribution, we can use a uniform distribution for the bounded domain, and for unbounded domains we can use a wide Gaussian distribution with zero mean or uniform distribution on a region large enough to cover our interested predictions.
Algorithm
4.1 Stochastic Collapsed Inference
We now present the model inference algorithm. The exact posterior of the latent random variables Z, and in (9) are infeasible to calculate because they are coupled in kernels and differential operators. While we can use variational approximations, they will introduce extra variational parameters, complicate the optimization and affect the integration of the physics knowledge. Therefore, we marginalize out all the latent variables to conduct collapsed inference to avoid approximating their complex posteriors. Specifically, we derive that p(y, 0 p(y ( X), where
/ p( y)p(0 Z Z)dZdedg
p(O X)
- Eptell Eptell /
- h)W E)dg Ep( EAr(Ejo,)) W (h(Z.,))0, E)].
(10)
where h(Z,e) = " T)
Note that ((-,.) is defined in
To allow to adjust the importance of the generative component and the influence of the physics during training, we weight the likelihood of the generative component by a free hyper-parameter
y 0. The weighted marginal likelihood (Warm 1989 Hu and Zidek 2002) is
Py 0|X)= = p(y X))
(11)
Our inference is to maximize the log weighted marginal likelihood to optimize the kernel parameters
in (DEEP(-.)) and x(...), the inverse noise variance T and unknown parameters in the differential equation, log P- (y, 0|X) log (N(y)o, K +T)) + 1 log (Epcz) However, the log likelihood is infeasible compute due the intractable expectation inside the logarithm. To address this problem, we use Jensen's is inequality on the log function to obtain a model evidence lower bound (ELBO), L < log p, (y, 0 where
[W (h(Z,e) e) Σ)] ).
== log (N K +++--11 -T) +
Σ)) Epc E [log (N (h(zz ,)
(12)
While L is still intractable, it is straightforward to maximize & with stochastic optimization. Each time, we generate a sample of the input locations from pZZ) and the noise from [N 1), denoted
by Z and ट. We then obtain Σ log (N(y) Wyy|0, - "D) + log (N(hhZ, 2) Σ) unbiased stochastic estimation of C. We calculate Vă as an unbiased stochastic gradient of C, with which any stochastic optimization to estimate the model parameters. While h(:,)) couples the
4

---

Under review as a conference paper at ICLR 2021
deep kernels and complex operators in v. it is differentiable and we can use automatic differentiation libraries to calculate the stochastic gradient conveniently.
The ELBO C in 12 is the GP log marginal likelihood plus an extra term, Ep(z)EN(ejo.1 [log (N ( =) Each element of h is obtained by applying the functional ei on the posterior sample of f(-) (see (6)) Jointly maximizing this term in L encourages that all the possible latent source values (at locations) obtained from the GP posterior function f()) (through the equation) should be considered as the samples of another GP. This can be viewed as soft constraint over the posterior function of the GP. Therefore, our ELBO is also posterior regularization objective (Ganchev et al.) 2010), and our inference estimates the standard deep-kernel GP model with soft regularization on its posterior distribution.
4.2 Algorithm Complexity
The time complexity for the inference of our model is O(N ((N++m"), +mn), because it involves the calculation for two GPs: one is the standard model, and the other is in the generative component. The time complexity for prediction is still ((N3), The space complexity is O(N2 m²), including the storage of the kernel matrices of the two GPs.
5 Related Work
An influential work, physics informed neural networks (PINNs) (Raissi et al.) 2019). were recently proposed to train neural networks that respect physical laws. The key idea is to use neural networks surrogate for the solution of the (partial) differential equation, and minimize the NN loos plus the residual error on a set of randomly collected collocation points in the input domain. Research along this line are quickly growing: (Mao et al. 2020 Jagtap et al.) 2020 Zhang et al., 2020 Chen et al. 2020 Pang et al. 2019), to name a few. While our work is enlightened by PINNs, there are several substantial differences. First, PINNs demand the form of the PDE is fully specified, i.e., (x)) = 0. while we assume there can be some unknown source (function), [x)] = g(x). Thus, our work to exploit incomplete physics knowledge. Second, we use the posterior of the deep-kernel GP to construct random surrogate for the PDE solution, and cast the integration of the physics into principled Bayesian framework to enable posterior inference and uncertainty quantification, while PINNs only conduct point estimations. Our experiments show that the incomplete physics knowledge can also improve uncertainty quantification. Note that in the mean time, Zhhnn et al. (2019) combined
polynomial chaos (Xiu and Karniadakis] 2002) and dropout (Gal and Ghahramani the total uncertainty for PINNs with stochastic PDEs
2016 to estimate
Lawrence al. 2007: Gao et al. 2008 Alvarez et al. 2009; 2013 Raissi et al.) 2017). For example, Graepel (2003) uses GPs solve the linear equation given observed noisy sources. He first defines the kernel for the solution function with which to derive the kernel for the source function. The kernel parameters are then estimated from the noisy source data, given which the solution can be predicted. Raissi iet al. (217)) assume both the noisy forces and solutions are observed, and they jointly model these examples in one single GP with a heterogeneous block covariance matrix. Latent force models (LFM) (Alvarez et al. 2009) make the same assumption about the differential equations as our work. LFMs convolve the Green's function of the equation with the kernel for the latent source to obtain the kernel for the target, and then learn the kernel parameters from data. While LFMs enable a hard encoding of physics, they rely on an analytical Green's 's function, which is not available for many equations. In addition, LFMs construct shallow kernels, which can be less expressive as deep kernels. Other works include (Calderhead et al.) 2009 [Barber and Wang] 2014 Macdonald et al.) 2015) Heinonen et al. 2018 Wenk et al. (020)) etc. They mainly focus estimating parameters/oppeattrss in ODEs without latent funetions/sources.
Many works have used GPs to model or learn physical systems (Graepel 2003
6 Experiments
6.1 Simulation
We first examined if PI-DKL can improve extrapolation with right physics knowledge. We generated two synthetic datasets. The first dataset, ISODE, was simulated from a first-order Ordinary Differ-
f(t) - D = g(t) where B B=1 D = 1, g(t) sin(2xt)exp(-1) and the initial condition f()) = 0.1. We set the time domain t E [o, 1]. We ran the finite dif- ference algorithm (Mitchell and Griffiths| 1980) to obtain the accurate solution. We chose 1. 001 equally spaced time points (fn = 0, 1000 = 1) and their solution values as the dataset. The second dataset, IdDiffusion, was simulated from diffusion equation with one dimensional spatial domain,
ential Equation (ODE), af()) B 2r
5

---

Under review as a conference paper at ICLR 2021
0f(2.)) Ør a 2x = g(x,1) 1) where a = 10, g00,(1 () = and the initial condition [ar 0) is a square wave. We set the domain (s, E [0.1] x [0. We ran a numerical solver to obtain the accurate solution. Then we discretized the entire spatial and time domain into 48 101 grid with equal spacing in each dimension. We retrieved the grid points and their solution values our dataset.
0f(2.1)
Competing methods. We compared with (1) shallow kemnel learning (SKL) with SE-ARD kenrel, (2) deep kernel learning (DKL), and (3) LFM, which uses SE-ARD for the latent source, and then convolves it with Green's 's function to obtain the kernel for the target function. To construct a deep kemel, we followed (Wilson et al. 2016a) to feed the input variables to (deep) neural network (NN) and calculated the RBF kernel over the neural network outputs (see ( Across our experiments, we used 5-layer NN, with 20 nodes in each hidden and output layer. We used tanh(.) as the activation function. For our method PI-DKL, we used the same deep kernel for the target function. As in LFM, we used SE-ARD kernel for the latent source. We set the number of virtual observations m 10 for the generative component, and uniformly sampled the input locations from the entire domain (see (I2). We chose the weight of the generative component y from f0.01, 0.05, 0.1, 0.5, 1, 5, 10). For both LFM and PI-DKL, the parameters of differential equations are unknown. All the methods were implemented with TensorFlow (Abadi et al. 2016). For our method, we used ADAM (Kingma and Ba 2014) for stochastic inference. We ran 10K epochs to ensure convergence. For the other methods, we used L-BFGS for optimization and set the maximum number of iterations to 5K.
For IstODE, we used the first 101 samples (t, E 0.1)) for training, and the remaining 900 samples (t, E (0. for test. We show the posterior distribution of the functions learned by all the methods We can that the predictions of SKL and DKL are largely biased when the test points are far from the training region [0. 0.1]. On average. DKL obtains better accuracy than SKL. The not-mean-square errors (RMSEs) are (DKL:0.21, SKL:0.25). As a comparison, the posterior means of LFM and PI-DKL. are much closer to the ground-truth in the test region, and the RMSEs are (LFM: 0.09, PI-DKL.: 0.04), showing the benefit of the physies. However, LFM is quite unstable in extrapolation: the farther away the test area, the more fluctuating the prediction. By contrast, PI-DKL obtains much smoother curves that are even closer to the ground-truth and smaller posterior variances in the test region. Hence, it shows that the LFM kernel obtained from shallow kernel convolution is less expressiveppowerful than the regularized deep kernel in PI-DKL. Note that unlike SKL/DKL, both LFM and PI-DKL estimated nontrivial posterior variances (i.e., not extremely close to 0) in the training region, implying that the physics also helps prevent overfitting.
and the ground-trute truth in Fig. [
Since for diffuusion equations, LFM cannot derive the kernel for time variable i,, for a fair comparison on IdDiffission, we fixed t = 0.5 and used the 48 spatial points as the training inputs. We then evaluated the posterior distribution of the function values at all the grid points (48 x 101) in the entire domain. We report the absolute difference between the posterior mean and ground-truth in Fig. - We that the prediction errors of SKL/DKL are close to 0 (dark colors) in regions close to the training data (t- = 0.5). However, when the test points are getting far away, say, close to the boundary 0 or 1), the error grows significantly (see the bright colors). Overall, DKL still achieves smaller extrapolation error than SKL, implying an advantage of using more flexible deep kernels. From Fig. 2c. we can see that while LFM misses the time information, it still exhibits better extrapolation results, as compared with SKL/DKL, showing the benefit of the physics PI-DKL achieves even smaller prediction error (i.e. darker) when t away from the training time point and exhibits even best extrapolation performance. The RMSES of all the methods are (SKL: 0.18, DKL: 0.11, LFM: 0.09, P1.DKL:0.07). We also report the predictive standard deviation (PSD) of each method in Fig. - We can see that the PSDs of SKL/DKL are both close to 0 in the training region, and quickly increase when the inputs move away (on average DKL shows smaller PSDs and smoother changes). By contrast, LFM and PI-DKL obtain PSDs quite uniformly across the entire domain and less than SKL/DKL. means that the physics knowledge help inhibit overfitting and reduce the uncertainty in extrapolation. Compared with LFM, PI-DKL. obtains even smaller PSDs (darker color) across the domain, showing even smaller uncertainty in extrapolation.
6.2 Real-World Applications
Metal Pollution in Swiss Jura. Next, we evaluated PI-DKL in real-world applications. We examined the predictive performance in terms of normalized RMSE (nRMSE) and test log-likelihood (LL). Due to the space limit, the test LL results are provided in the supplementary material. We first considered predicting the metal concentration in Swiss Jura. The data were collected from 300 locations in a 14.5 km2 region tpt / /ddrr. 10/ /ran/ gtttt////////t/ /tr////aanaaaaaa///a/aaaaaaaaaaa jura. html). The diffusion of the metal concentration is naturally modelled by diffusion equation with the two-dimensional

---

Under review as a conference paper at ICLR 2021
0.5
0.5
0.5
0.5
0.1
0.5
0.1
0.5
0.1
0.5
0.1
0.5
(a) SKL
(b) DKL
(c) LFM
(d) PI-DKL
Figure I: The posterior distribution of the learned solution functions on IsODE. The red lines in the middle are the posterior means and the red dashed lines on the boundary of the shaded region the posterior mean plus/minus posterior standard deviation. The black line is the ground-truth solution. The training inputs stay in [0, 0.1] (left to the green line).
0.5
0.25
x 0.5
x 0.5
0.5
0.5
0.5
t
(a) SKL
(b) DKL
(c) LFM
(d) PI-DKL
0.15
0.1
0.5
0.5
x0.5
0.05
0
0
0,5
0.5
0
0.5
0
0.5
(e) SKL (PSD)
(f) DKL (PSD)
(g) LFM (PSD)
(h) PI-DKL (PSD)
Figure 2: The absolute value of the difference between the postcrior mean the ground-truth (Ist row) and posterior standard deviation (2nd row) on IdDiffusion. The training examples stay on = 0.5 (the green linc).
x10
0.9
nRMSE
ORMSE
ORMSE
0.6
0.7
0.85
7M5
BN
19151
THS
DA
(a) Copper
(b) Cadmium
(c) Joint Motion
Figure 3: Metal concentration prediction in Swiss Jura (a. b) and joint angle prediction motion capture (c).
results are averaged over runs. The normalized root-mean-squaa error (nRMSE) in each run computed by normalizing the RMSE by the mean of the test outputs.
spatial domain, 0f(m).m2.t) a( (85 (m
(f(11.22.1)) where f(s - is the concentration of 82
the metal at particular location and time point. However, the dataset do not include the time t, when
these concentrations were measured. LFM considers the initial condition f 222 as the latent function and obtains a kernel of the locations where l, can be viewed as a kernel parameter learned from data. In our approach,
estimated the solution function f.. h( x22) [( 22, f.)) ((=1,=22) where the latent function
Hence, the equation can be viewed as ah"(z) (h2(2)
duf
dag
g( .2) 05 3f(m) (a(() lt=t,- We were interested in predicting the concentration of cadium and 35 copper. The input variables include the coordinates of the location (21, the concentrations of (nickel, zinc] for cadmium, and [lead, nickel, zinc} for copper. For PI-DKL, we selected m from {10, 50, 100, 200, 500} for the generative component and 7 from {0.01, 0.05, 0.1, 0.5, 2, 5, 10). We normalized the training inputs and then sampled latent inputs from A ( I) in model estimation. For LFM, we varied the number of latent forces from 1133. 5). We randomly selected 50 samples for training, and used the remaining 250 samples for test. We repeated the experiments for times, and report the average aRMSE and its standard deviation of each method in Fig. 3 and b. PI-DKI.
7

---

Under review as a conference paper at ICLR 2021
outperforms all the competing approaches for both prediction tasks. PI-DKL always significantly improves upon SKL and DKL (p < 0.05). In addition, PI-DKL significantly outperforms LFM in predicting Cadium concentration (Fig. Ep). Note that LFM does improve upon SKL in predicting Copper concentration (Fig. 3 but not as significant as PI-DKL..
Motion Capture. We then looked into predicting trajectories of joints in the motion capture appli- cation. To this end, we used CMU motion capture database (fatt::/mmoc from which we used the samples collected from subject 35 in the walk and jog motion lasting for 2,644 seconds. We trained all the models to predict the angles of Joint 60 along with time. We used the first order ODE in simulation to represent the physical model, based on which we ran LFM and PI-DKL. Note this physical system might be oversimplified (Alvarez et al.) (2009). For LFM, we varied the number of latent forces from (1,3, 5). Again, we randomly selected 500 samples for training and 2, 000 samples for test. We repeated the experiments for 5 times and report the average nRMSE and its standard deviation in Fig. 33. As we can see, PI-DKL improves upon all the competing methods by large margin. Note that LFM is even far worse than SKL. This might because LFM over-exploits the simp-simllified physics, which harms the prediction. By contrast, PI-DKL allows us tune the number of virtual observations and the likelihood weight (y in 12 and hence can consistently improve upon DKL.
cs. cmu edu/),
PM2.5 in Salt Lake City. Second, consid- ered predicting the Particulate Matter (PM2.5) levels across Salt Lake City. The dataset were collected from sensors' reads at differ- ent time and locations pectps://açann org/). We chose the time range from 07/04/2018 to 07/06/2018. Following (Wang et al.] 2018). we used the diffusion equation plus source term (i.e,, the latent function) to
07
nANGE
MAMGE
06
06
Tap
Ded
(a) PM2.5
(b) Traffic flow
Figure 4: PM2.5 and traffic flow prediction. (0) g( 22, l). where f is the con-
represent the physical model, 2f(2).22.1) a Σz=
8x)
centration level and g the source term. The input variables include both the location coordinates and detailed time points. Since LFM cannot construct full kernel of the input variables from the physics, we did not test it to avoid unfair comparisons. We trained SKL and DKL with both the spatial and time imputs. We randomly selected 500 samples for training and 000 samples for test. We repeated the experiments for 5 times and report the average nRMSE and its standard deviation in Fig. 4a. As we can see, with a more expressive kernel, DKL. improves upon SKI. significantly, and with the incorporation of the physics, PI-DKL in turn outperforms DKL significantly (p < 0.05).
Wayh-Way Traffic Flow Prediction. Finally, we applied PI-DKL to predict the traffic flow in the interstate highway 215 across Utah state. The Utah Department of Transportation (UDOT) has installed sensors every a few miles along the high way. Each sensor counts the number of vehicles passed every minute, and sends the data back to central database. The real time data and road conditions are available at at/htpps: /udott teri s pemsp com/ We used the data collected by 20 sensors continuously installed in a segment of 30 miles, and the time was chosen from 08/05/2019 to 08/1 1/2019. The input variables include the location coordinates of each sensor and the time of each read. Following (Nagatani 2000), we used the Burger's is equation plus source term to describe the system,
Øf σi + J -ド Jz. - j== as² g(0, T2, !), where f is the traffic flow, v the 2zy
unknown viscous coefficient, and y the source term, i.e., the latent function. Note that the equation is nonlinear and we do not have an analytical form of Green's 's function. Hence we cannot use LFM to incorporate the physics to enhance GP training. Hence we compared with SKL and DKL only. We randomly selected 500 and 2, 000 samples for training and test, respectively. and repeated for 5 times. The average nRMSES and the standard deviations are reported in Fig. 4p. As we can see, DKL significantly outperforms SKL, which demonstrates the advantage of the more expressive, deep kernel. More important, PRGP further improves upon DKL, showing that the physics incorporated by our approach indeed promotes the predictive performance.
7 Conclusion
We have presented PI-DKL, a physics informed deep kernel learning model that can flexibly incor- porate physics knowledge from incomplete differential equations to improve function learning and uncertainty quantification. In the future work, we will extend our model with sparse approxima- tions (Hensman et al. 2013] Wilson et al. 2016b) exploit physics in large-scale applications.
8

---

Under review as a conference paper at ICLR 2021
References
Abadi, M., Barham, P., Chen, J., Chen, Z., Davis, A., Dean, J., Devin, M., Ghemawat, S., Irving, G., Isard, M., et al. (2016). Tensorflow; A system for large-scale machine learning. In 12th USENIX Symposium Operating Systems Design and Implementation (OSDI 16), pages 283.
Alvarez, M., Luengo, D., and Lawrence, N. D. (2009). Latent force models. In Artificial Intelligence and Statistics, pages 9-16.
Alvarez, M. A., Luengo, D., and Lawrence, N. D. (2013). Linear latent force models using Gaussian processes, IEEE transactions on pattern analysis and machine intelligence, 35691)266332275
Barber, D. and Wang, Y. (2014). Gaussian processes for bayesian estimation ordinary differential equations. In International Conference Machine Learning, pages 1485-1493.
Calderhead, B., Girolami, M., and Lawrence, N. D. (2009). Accelerating bayesian inference over nonlinear differential equations with gaussian processes. In Advances in neural information processing systems, pages 217-224.
Chen, Y., Lu, L., Karniadakis, G. E., and Dal Negro, L. (2020). Physics-informed neural networks for inverse problems in nano-optics and metamaterials. Optics Express, 28(8):1
Gal, Y. and Ghahramani, Z. (2016). Dropout as bayesian approximation: Representing model uncertainty in deep learning. In international conference on machine learning, pages 1050-1059.
Ganchev, K., Gillenwater, J., Taskar, B., et al. (2010). Posterior regularization for structured latent variable models. Journal Machine Learning Research, 11(Jul)2201-2200
Gao, P., Honkela, A., Rattray, M., and Lawrence, N. D. (2008). Gaussian process modelling of latent chemical species: applications to inferring transcription factor activities. Bioinformatics, 24(16):
Graepel, T. (2003). Solving noisy linear operator equations by gaussian processes: Application to ordinary and partial differential equations. ICML, pages 234-241.
Heinonen, M., Yildiz, C., Mannerström, H., Intosalmi, J., and Lähdesmäki, H. (2018). Learning unknown ode models with gaussian processes. In International Conference on Machine Learning, pages 1959-1968.
Hensman, J., Fusi, N., and Lawrence, N. D. (2013). Gaussian processes for big data. Proceedings of the Conference on Uncertainty in Artificial Intelligence (UAI).
Hu, F. and Zidek, V. (2002). The weighted likelihood. Canadian Journal of Statistics, 30(3):347- 371.
Jagtap, A. D., Kawaguchi, K., and Karniadakis, G. E. (2020). Adaptive activation functions accelerate convergence in deep and inysics-informed neural networks. Journal of Computational Physics, 404:109136.
Kingma, D. and Ba, J. (2014). Adam: method for stochastic optimization. arXiv preprint 41X.v:1412.6980.
Lapidus, L. and Pinder, G. F. (2011). Numerical solution of partial differential equations in science and engineering. John Wiley & Sons.
Lasserre, J. A., Bishop, C. M., and Minka, P. (2006). Principled hybrids of generative and discriminative models. In 2006 IEEE Computer Society Conference Computer Vision and Pattern Recognition (CVPR 06), volume I, pages 87-94. IEEE.
Lawrence, N. D., Sanguinetti, G., and Rattray, M. (2007). Modelling transcriptional regulation using gaussian processes. In Advances in Neural Information Processing Systems, pages 785-792.
Macdonald, B., Higham, C., and Husmeier, D. (2015). Controversy in mechanistic modelling with gaussian processes. Proceedings of Machine Learning Research, 37:1539-155 1547.

---

Under review as a conference paper at ICLR 2021
Mao, Z., Jagtap, A. D., and Karniadakis, G. E. (2020). informed neural networks for high-speed flows. Computer Methods in Applied Mechanics and Engineering, 360:112789.
Mitchell, A. R. and Griffiths, D. F. (1980). The finite difference method in partial differential equations. Number BOOK. John Wiley.
Nagatani, T. (2000). Density waves in traffic flow. Physical Review E, (44:3366.
Olsen-Kettle, (2011). Numerical solution of partial differential equations. Lecture notes at University of Queensland, Australia.
Pang. G., Lu, L., and Kamiaaakii, G. E. (2019). fpinns: Fractional physics-informed neural networks. SIAM Journal on Scientific Computing. 41(4):A2033-22266
Raissi, M., Perdikaris, P., and Karniadakis, G. E. (2017). Machine learning of linear differential
equations using gaussian processes. Journal of Computational Physics, 348:683-699 693.
Raissi, M., Perdikaris, P., and Karniadakis, G. (2019). Physics-informet neural networks:
deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations. Journal Computational Physics, 378:686-707.
Wang. Y., Wang, H., Chang, S., and Avram, A. (2018). Prediction of daily pm2. 5 concentration in china using partial differential equations. PloS one, 13(6):e0197666
Warm, A. (1989). Weighted likelihood estimation of ability item response theory, Psychometrika, 54(3):427-40..
Wenk, P.. Abbati, G., Osborne, M. A., Schölkopf, B., Krause, A., and Bauer, S. (2020). Odin: Ode-informed regression for parameter and state inference time-continuous dynamical systems. In AAAL pages 6364-6371.
Williams, C. K. and Rasmussen, C. E. (2006). Gaussian processes for machine learning, volume 2.
MIT press Cambridge, MA.
Wilson, A. G., Hu, Z., Salakhutdinov, R., and Xing, E. P. (2016a). Deep kernel learning. In Artificial Intelligence and Statistics, pages 370-378.
Wilson, A. G., Hu, Z., Salakhutdinov, R., and Xing, E. (2016b). Stochastic variational deep kernel learning. In Advances in Neural Information Processing Systems, pages 2586-2594,
Xiu, D. and Karniadakis, G. E. (2002). The wiener-asker polynomial chaos for stochastic differential equations. SIAM journal on scientific computing, 24(2):619 644644
Zhang, D., Guo, L., and Karniadakis, G. E. (2020). Learning in modal space: Solving time-dependent stochastic pdes using thysics-informed neural networks. SIAM Journal on Scientific Computing. A2(2):A639-A66)
Zhang, D., Lu, L., Guo, L., and Karniadakis, G. E. (2019). Quantifying total uncertainty in physics-informed neural networks for solving forward and inverse stochastic problems. Journal of Computational Physics, 397:108850.
10


Reviews:
[{"review_id": "vuQKw3jioQv", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "This work proposes a deep Gaussian Process (GP) framework for data modeling informed by dynamical systems. \nThe analogous to the one of gradient matching, where a GP regression problem is penalised by constraints taking the form of differential equations acting on the GP itself. This study merges this framework to the scalable deep-GP framework of Wilson et al, and proposes an approximated inference setting in which the problem can be optimised through stochastic variational inference. \nThe proposed method is tested in several experimental scenarios and compared to standard kernel-based methods, and to the Latent Force Model of Alvarez and colleagues. The results show promising performances in terms of prediction and stability of extrapolation.\n\nThe proposed methodology contributes to the domain of gradient-matching and extends the classical approaches to allow for scalability in deep models. In this sense, the comparison with respect to the state of the art appears weak. The comparisons proposed in the study are with respect to either shallow models, or deep-models not allowing the integration of physics-informed contraints. Recent works (e.g. [1,2,3,4]) could provide a fairer benchmark for the proposed approach. In this sense, the feeling is that proposed methodology largely overlaps with these more recent studies, for example for what concerns the use of deep models and variational inference schemes. In particular, the idea of soft-regularisation of deep-GP models has been already explored in [2], and the authors may want to compare the method with respect to this approach. \n\nIt is not clear whether the parameters of the dynamical systems in the term h(Z,\\epsilon) of (12) are inferred. If this is the case, the ability of the framework in estimating the systems’ should be assessed and compared to other moment matching approaches. This aspect is overlooked in the paper, although it is quite relevant for interpretability purposes. ", "extra_scores": {}, "extra_text": {}}, {"review_id": "gTo_Ff5nbw", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "Post-discussion update: The authors have clarified their work considerably, and I believe the work is probably correct. However, the paper still suffers from poor presentation and poorly-motivated or justified modelling choices. The current version of the paper has not been updated, and all below issues thus stand. The paper overall presents a good idea, but I believe the authors made poor modelling choices, which led the kludgy math. \n\n----\n\nThe paper proposes to combine deep kernel learning (DKL) with PDE/ODE prior knowledge for learning spatiotemporal systems. The idea is sound and sensible, although it represents an incremental combination of two existing techniques. The differential assumption clearly improves in spatiotemporal experiments, which also motivates the idea. \n\nThe paper has unfortunately three major mathematical errors. \n\nFirst, the very first equation of the method (eq 4) is already incorrect: g(x) is not a function of input “x” but instead it’s a function of solution “f”. The paper itself acknowledges this by giving series of examples which all are functions of “f”. The following GP prior for g(x) is then placed on incorrect inputs and seems misguided or at least insufficient to model the differential. Intuitively this is also easy to see: the differential of the solution is obviously not just a function of space and time (but of solution as well). \n\nSecond, the eq 6 states that a non-linear differential of a GP is some other Gaussian process. This is incorrect: non-linear transforms of a Gaussian process does not retain Gaussianity. Since this is a central definition of the paper, the whole method is likely incorrect.\n\nThird, the paper introduces a bizarre concept of constant zero vectors as random variables, ie. p(0|g). This is clearly wrong, and the probabilistic model is then wrong as well. While this could be easily fixed, the model is an example of GP-matching (which has known pathologies), and this has been extensively studied by many authors (eg. Wenk’19, Wenk’18, Gorbach’17, MacDonalds’15).\n\nGiven these obvious errors, the paper needs a major revision.\n\n\nTechnical comments:\no Eq 4 should be g(f,x), since later the \\varphi is a function of f and its x-derivatives. Currently none of the examples (eg. burger) follow eq 4. Also eq6 is also wrong with same argument\no Please define (mathematically) the tackled problem and the problem domain. The paper starts by discussing classification, but at sec 3 the context suddenly changes into PDEs without the reader being informed about this.\no The paragraph “to incorporate..” is difficult to follow since its technical but does not open up the math yet. It would help to write this in more conceptual way\no It would greatly help the reader understand the method to include the sup-fig1 in the main paper\no It seems that both f and g are assumed to be separate GPs. This is unclear from text, please clarify\no eq7 is strange, since it defines a dirac between a random variable g and a sample of h. This is nonsensical.\no The 0-sampling is obviously wrong and unnecessary. Double priors are perfectly fine in Bayesian modelling, and eq9 could be defined without the 0-stuff. The double priors are called product of expert priors, see eg. MacDonalds'15 or Wenk’19. They do have their own pathologies, but these papers discuss them extensively.  \n ", "extra_scores": {}, "extra_text": {}}, {"review_id": "eK3tye_1mf", "reviewer": "AnonReviewer4", "summary": "", "questions": "", "limitations": "", "rating": 8, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "The paper is clearly written, technically sound and innovative.\n\n### Impact:\nThe paper is bringing an important contribution in the domain of Physics-Informed Machine Learning. It proposes an interesting and technically sounds (and not obvious) way to combine multiple successful components of recent literature in the field.\n\n### Clarity and technical soundness:\nThe paper is clear, technically sound and manipulates tools in an advanced way.\nIt would be useful to provide more intuition and high-level interpretation on the most technical aspects of Deep Kernel inference and on the usage of ELBO methods. That would go at the benefit of understanding for the non-specialists of Kernel-based learning who are interested in Physics Informed ML.\n\n### Results:\nThe results are extensive, are exploiting both simulation-type data and real-world data.\nThe benchmark used for comparison are relevant.\n\n### Applicability:\nIt would be interesting to have a more high-level interpretation and analysis of the applicability of the method to 3D simulation data. \nIn 4.2, it was not clear to me what N and M stand for. Could you clarify it and provide a more high-level analysis with it too?", "extra_scores": {}, "extra_text": {}}, {"review_id": "0W4eVJrGqdM", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 7, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "### Summary of my understanding\n\nThe authors address the problem of function estimation given noisy observations of the function's values. The estimation is done in a non-parametric manner as posterior inference of Gaussian processes. The authors' proposal is to incorporate physics knowledge into this process. The knowledge here refers to a differential equation on the function to be estimated ($f$), and it may include unknown parameters and an unknown external force term $g$. They derive a joint distribution of 1) observation of $f$, i.e., $y$, 2) \"virtual\" observations $0$ of $g$, 3) latent variable $Z$ for latent force $g$, 4) randomness between $f$ and $g$, and 5) values of $g$. Afterward, they propose to perform a collapsed inference based on a lower bound of the log marginal likelihood (of observations $y$ and virtual observations $0$).\n\n### Evaluation\n\nThe problem setting is clear. The proposed method is technically sound, and the experiments are enough convincing to show its superiority to baselines. Regardless of the common points with the latent force models [Alvarez+ 09] in terms of the problem setting, this paper proposes a somewhat new technique to perform inference under the presence of differential equation-based constraints on a target function. I think this paper is a sound work.\n\n### Note\n\nI think a concern lies in the paper's clear violation of the formatting instruction of ICLR 2021. Especially, the too narrow margins before and after headings seem to compress the main text, which might seem quite unfair. I delegate judgment on this regard to the chairs and for now, I decided the rating ignoring this matter.\n", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: The paper presents a framework for incorporating physics knowledge (through, potentially incomplete, differential equations) into the deep kernel learning approach of Wilson et al. The reviewers found the paper addresses an important problem and presents good results.  However, one of the main issues raised by R1 is that, although the proposed method can be applied to broader settings such as that of incomplete differential equations, there are still regimes where the comparison is not only possible but perhaps insightful. An example baseline is the work of Lorenzi and Filippone, “Constraining the Dynamics of Deep Probabilistic Models” (ICML, 2018). Another critical issue, raised by R4, is the insufficient clarity in the presentation. Many of the concerns raised by this reviewer were clarified in the discussion and I thank the authors for their engagement. However, the AC believes some of the points raised by R4 in this regard were left unaddressed in the paper and the manuscript does indeed require at least one more iteration.

The format violation concerns raised during the reviewing process did not affect the decision on this paper, as the PCs confirmed that they did not meet the bar for desk rejection and recommended to assess the paper on its technical merits.
Decision: Reject