In this notebook we consider the latent force model paradigm (Alvarez et al, 2013, Alvarez et al, 2011, Alvarez et al, 2009) to drive spatio-temporal differential equation with a Gaussian process. Latent force models are differential equations whose initial conditions or driving forces are given by a stochastic process. Linear latent force models are linear (partial or ordinary) differential equations. If a linear latent force model is driven by a Gaussian process latent force, this describes a joint Gaussian process distribution across all the variables of interest. The Gaussian process covariance function encodes the relationships between the variables, as proscribed by the differential equation, through the covariance function.
This direction of research is part of our ideas as to how to merge mechanistic and statistical models of data. It also maps onto our ideas about 'Gaussian processes over everything'. The covariance function interelates the protein and mRNA concentrations in the same model.
A model of post-transcriptional processing is formulated to describe the spatio-temporal Drosophila protein expression data (Becker et al, 2013, Alvarez et al, 2011). Protein production is considered to be linearly dependent on the concentration of mRNA at an earlier time point. The model also allows for diffusion of protein between nuclei and linear protein decay. These processes are dependent on the diffusion parameter and the degradation rate of protein respectively.
The coefficients $a$, $b$ and $c$ are unknown. In this study, we use Gaussian process with an exponentiated quadratic kernel as a prior over $y_{x,t}$ (protein). The kernel of $f_{x,t}$ (mRNA) is derived by applying the partial differential operator on the spatial-temporal kernel of protein. The multi-output Gaussian process are developed by combining the covariance matrix of mRNA and protein and their cross covariance.
%matplotlib inline
%config InlineBackend.figure_format = 'svg'
import numpy as np
import pylab as pb
import GPy
import pandas
from pandas import read_csv
The spatio-temporal multi-output partial differential equation covariance functions have been developed with the kernel name GPy.kern.ODE_st
in GPy. The inputs are one dimension spatial data, one dimension temporal data and one dimension index which is used to indicate $f$ and $y$.
data = GPy.util.datasets.drosophila_knirps()
The next thing to do is to compose the data ready for presentation to GPy. Here we need to use the time field as the input and a concatanation of the expression levels as the output. We need to provide a corresponding index for to each input to describe what output is being represented.
Now we set up the covariance function with the relevant parameters.
kern = GPy.kern.ODE_st(input_dim=3,
a=1., b=1., c=1.,
variance_Yx=1., variance_Yt=1.,
lengthscale_Yx=15.,
lengthscale_Yt=15.)
With the data correctly presented and the covariance function defined, we are ready to proceed with the Gaussian process regression.
data['X']
data['Y']
array([[ 5.90000000e-01], [ 5.00000000e-01], [ 9.00000000e-02], [ 2.50000000e-01], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 8.00000000e-02], [ 2.30000000e-01], [ 0.00000000e+00], [ 1.90000000e-01], [ 0.00000000e+00], [ 2.70000000e-01], [ 1.40000000e-01], [ 1.00000000e-01], [ 2.70000000e-01], [ 8.20000000e-01], [ 2.04000000e+00], [ 2.01000000e+00], [ 3.19000000e+00], [ 4.74000000e+00], [ 8.47000000e+00], [ 9.74000000e+00], [ 1.44700000e+01], [ 2.02400000e+01], [ 2.74100000e+01], [ 3.54200000e+01], [ 4.33400000e+01], [ 5.29000000e+01], [ 6.03400000e+01], [ 6.57900000e+01], [ 7.28200000e+01], [ 7.45700000e+01], [ 7.55900000e+01], [ 7.54700000e+01], [ 7.48900000e+01], [ 7.24400000e+01], [ 6.73000000e+01], [ 6.15700000e+01], [ 5.56400000e+01], [ 4.92700000e+01], [ 4.11100000e+01], [ 3.35900000e+01], [ 2.76600000e+01], [ 2.22200000e+01], [ 1.69400000e+01], [ 1.26200000e+01], [ 1.03800000e+01], [ 7.10000000e+00], [ 5.81000000e+00], [ 4.27000000e+00], [ 3.05000000e+00], [ 2.55000000e+00], [ 1.76000000e+00], [ 1.35000000e+00], [ 1.44000000e+00], [ 1.22000000e+00], [ 4.00000000e-02], [ 1.00000000e-02], [ 0.00000000e+00], [ 7.00000000e-02], [ 0.00000000e+00], [ 8.00000000e-02], [ 0.00000000e+00], [ 1.60000000e-01], [ 8.00000000e-01], [ 4.90000000e-01], [ 9.70000000e-01], [ 1.10000000e+00], [ 1.28000000e+00], [ 1.50000000e+00], [ 1.15000000e+00], [ 1.74000000e+00], [ 2.31000000e+00], [ 2.45000000e+00], [ 3.40000000e+00], [ 4.53000000e+00], [ 5.91000000e+00], [ 9.42000000e+00], [ 1.37000000e+01], [ 2.05500000e+01], [ 2.95600000e+01], [ 4.14100000e+01], [ 5.40700000e+01], [ 6.88200000e+01], [ 8.37000000e+01], [ 9.77600000e+01], [ 1.09540000e+02], [ 1.16130000e+02], [ 1.20800000e+02], [ 1.21370000e+02], [ 1.19050000e+02], [ 1.14270000e+02], [ 1.06430000e+02], [ 9.58400000e+01], [ 8.20200000e+01], [ 7.07700000e+01], [ 5.85700000e+01], [ 4.64300000e+01], [ 3.57600000e+01], [ 2.91800000e+01], [ 2.30300000e+01], [ 1.71900000e+01], [ 1.42700000e+01], [ 1.05000000e+01], [ 8.81000000e+00], [ 7.22000000e+00], [ 6.12000000e+00], [ 5.20000000e+00], [ 4.17000000e+00], [ 3.32000000e+00], [ 3.02000000e+00], [ 2.76000000e+00], [ 1.97000000e+00], [ 4.00000000e-01], [ 2.60000000e-01], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 4.10000000e-01], [ 1.57000000e+00], [ 3.16000000e+00], [ 5.27000000e+00], [ 8.68000000e+00], [ 1.48000000e+01], [ 2.42800000e+01], [ 3.78200000e+01], [ 5.36300000e+01], [ 7.41300000e+01], [ 9.58900000e+01], [ 1.17020000e+02], [ 1.36130000e+02], [ 1.49680000e+02], [ 1.57060000e+02], [ 1.59200000e+02], [ 1.56380000e+02], [ 1.49960000e+02], [ 1.36460000e+02], [ 1.18340000e+02], [ 9.93600000e+01], [ 8.16100000e+01], [ 6.28500000e+01], [ 4.85500000e+01], [ 3.77600000e+01], [ 2.90100000e+01], [ 2.17400000e+01], [ 1.62400000e+01], [ 1.18200000e+01], [ 9.01000000e+00], [ 6.69000000e+00], [ 5.40000000e+00], [ 4.44000000e+00], [ 3.80000000e+00], [ 3.08000000e+00], [ 2.46000000e+00], [ 2.06000000e+00], [ 1.63000000e+00], [ 1.60000000e+00], [ 1.44000000e+00], [ 3.80000000e-01], [ 1.60000000e-01], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 2.60000000e-01], [ 4.40000000e-01], [ 1.55000000e+00], [ 1.38000000e+00], [ 1.16000000e+00], [ 1.54000000e+00], [ 1.06000000e+00], [ 9.30000000e-01], [ 1.47000000e+00], [ 2.04000000e+00], [ 3.61000000e+00], [ 4.61000000e+00], [ 6.02000000e+00], [ 8.15000000e+00], [ 1.25200000e+01], [ 1.92400000e+01], [ 3.05000000e+01], [ 4.59700000e+01], [ 6.48700000e+01], [ 9.24300000e+01], [ 1.23410000e+02], [ 1.49660000e+02], [ 1.71480000e+02], [ 1.85620000e+02], [ 1.92350000e+02], [ 1.93720000e+02], [ 1.86770000e+02], [ 1.75700000e+02], [ 1.58920000e+02], [ 1.37310000e+02], [ 1.12660000e+02], [ 8.71900000e+01], [ 6.75500000e+01], [ 5.13100000e+01], [ 3.89600000e+01], [ 2.91500000e+01], [ 2.24100000e+01], [ 1.70300000e+01], [ 1.31100000e+01], [ 9.82000000e+00], [ 8.22000000e+00], [ 7.00000000e+00], [ 6.37000000e+00], [ 5.53000000e+00], [ 5.27000000e+00], [ 4.75000000e+00], [ 4.11000000e+00], [ 3.31000000e+00], [ 2.78000000e+00], [ 2.11000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 1.40000000e-01], [ 3.00000000e-01], [ 3.20000000e-01], [ 7.00000000e-02], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 4.90000000e-01], [ 1.63000000e+00], [ 2.71000000e+00], [ 3.17000000e+00], [ 4.86000000e+00], [ 8.34000000e+00], [ 1.60300000e+01], [ 2.77100000e+01], [ 4.26200000e+01], [ 6.24900000e+01], [ 8.68800000e+01], [ 1.22690000e+02], [ 1.57310000e+02], [ 1.81910000e+02], [ 1.98100000e+02], [ 2.04160000e+02], [ 2.07460000e+02], [ 2.04180000e+02], [ 1.96550000e+02], [ 1.80970000e+02], [ 1.56840000e+02], [ 1.24870000e+02], [ 9.50400000e+01], [ 6.91700000e+01], [ 5.02900000e+01], [ 3.68100000e+01], [ 2.77600000e+01], [ 2.00500000e+01], [ 1.53700000e+01], [ 1.21600000e+01], [ 9.60000000e+00], [ 7.65000000e+00], [ 6.37000000e+00], [ 5.76000000e+00], [ 5.11000000e+00], [ 4.38000000e+00], [ 3.70000000e+00], [ 2.94000000e+00], [ 2.46000000e+00], [ 1.85000000e+00], [ 1.33000000e+00], [ 6.50000000e-01], [ 1.70000000e-01], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 1.60000000e-01], [ 2.70000000e-01], [ 2.10000000e-01], [ 1.00000000e-01], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 3.00000000e-01], [ 1.56000000e+00], [ 2.69000000e+00], [ 3.37000000e+00], [ 4.87000000e+00], [ 9.06000000e+00], [ 1.78500000e+01], [ 2.90600000e+01], [ 4.27800000e+01], [ 6.01100000e+01], [ 8.42900000e+01], [ 1.16780000e+02], [ 1.48300000e+02], [ 1.71600000e+02], [ 1.88430000e+02], [ 1.97980000e+02], [ 1.99570000e+02], [ 1.95680000e+02], [ 1.86220000e+02], [ 1.68300000e+02], [ 1.42590000e+02], [ 1.16470000e+02], [ 8.96500000e+01], [ 6.17400000e+01], [ 4.44800000e+01], [ 3.31100000e+01], [ 2.39100000e+01], [ 1.75700000e+01], [ 1.30500000e+01], [ 9.73000000e+00], [ 7.60000000e+00], [ 6.15000000e+00], [ 5.53000000e+00], [ 5.03000000e+00], [ 4.36000000e+00], [ 3.90000000e+00], [ 3.06000000e+00], [ 2.52000000e+00], [ 1.74000000e+00], [ 1.58000000e+00], [ 1.05000000e+00], [ 7.50000000e-01], [ 1.26000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 5.30000000e-01], [ 6.20000000e-01], [ 2.40000000e-01], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 1.39000000e+00], [ 2.45000000e+00], [ 3.01000000e+00], [ 3.79000000e+00], [ 6.01000000e+00], [ 1.32400000e+01], [ 2.62000000e+01], [ 4.06700000e+01], [ 5.35400000e+01], [ 7.02200000e+01], [ 9.75800000e+01], [ 1.32250000e+02], [ 1.61920000e+02], [ 1.84690000e+02], [ 1.96490000e+02], [ 1.93580000e+02], [ 1.87730000e+02], [ 1.80960000e+02], [ 1.65470000e+02], [ 1.42400000e+02], [ 1.14660000e+02], [ 8.52300000e+01], [ 5.92000000e+01], [ 4.41300000e+01], [ 3.23600000e+01], [ 2.28700000e+01], [ 1.60200000e+01], [ 1.15200000e+01], [ 8.47000000e+00], [ 6.32000000e+00], [ 5.36000000e+00], [ 4.78000000e+00], [ 4.19000000e+00], [ 3.67000000e+00], [ 3.09000000e+00], [ 1.98000000e+00], [ 1.10000000e+00], [ 5.40000000e-01], [ 7.00000000e-02], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 4.80000000e+00], [ 9.30000000e-01], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 1.30000000e-01], [ 2.11000000e+00], [ 2.49000000e+00], [ 1.71000000e+00], [ 3.30000000e-01], [ 2.40000000e-01], [ 1.80000000e-01], [ 1.00000000e-01], [ 4.40000000e-01], [ 3.31000000e+00], [ 4.53000000e+00], [ 4.35000000e+00], [ 4.78000000e+00], [ 8.90000000e+00], [ 2.19300000e+01], [ 3.68400000e+01], [ 4.94300000e+01], [ 6.58800000e+01], [ 8.36400000e+01], [ 1.08710000e+02], [ 1.35760000e+02], [ 1.61310000e+02], [ 1.78790000e+02], [ 1.75260000e+02], [ 1.65930000e+02], [ 1.64850000e+02], [ 1.53460000e+02], [ 1.34780000e+02], [ 1.10040000e+02], [ 8.71300000e+01], [ 6.10500000e+01], [ 4.57600000e+01], [ 3.38100000e+01], [ 2.53000000e+01], [ 1.76600000e+01], [ 1.24000000e+01], [ 8.71000000e+00], [ 6.55000000e+00], [ 5.65000000e+00], [ 5.49000000e+00], [ 4.93000000e+00], [ 3.95000000e+00], [ 3.25000000e+00], [ 2.31000000e+00], [ 1.63000000e+00], [ 1.10000000e+00], [ 1.70000000e-01], [ 5.60000000e-01], [ 1.50000000e-01], [ 2.40000000e-01], [ 0.00000000e+00], [ 7.49000000e+00], [ 6.96000000e+00], [ 4.73000000e+00], [ 4.20000000e+00], [ 6.12000000e+00], [ 7.15000000e+00], [ 6.70000000e+00], [ 6.91000000e+00], [ 5.57000000e+00], [ 6.08000000e+00], [ 7.77000000e+00], [ 5.54000000e+00], [ 4.36000000e+00], [ 4.27000000e+00], [ 5.92000000e+00], [ 2.99000000e+00], [ 3.46000000e+00], [ 4.01000000e+00], [ 3.55000000e+00], [ 5.34000000e+00], [ 4.55000000e+00], [ 8.77000000e+00], [ 1.36900000e+01], [ 1.92600000e+01], [ 2.82600000e+01], [ 4.04600000e+01], [ 5.21300000e+01], [ 6.87400000e+01], [ 8.51200000e+01], [ 1.01550000e+02], [ 1.09970000e+02], [ 1.21430000e+02], [ 1.25950000e+02], [ 1.36690000e+02], [ 1.21220000e+02], [ 1.17650000e+02], [ 1.06490000e+02], [ 9.43300000e+01], [ 7.57700000e+01], [ 6.22000000e+01], [ 5.63100000e+01], [ 3.78600000e+01], [ 2.97000000e+01], [ 2.27400000e+01], [ 2.01100000e+01], [ 1.45500000e+01], [ 1.08300000e+01], [ 9.30000000e+00], [ 7.59000000e+00], [ 5.80000000e+00], [ 2.81000000e+00], [ 3.60000000e+00], [ 3.06000000e+00], [ 1.35000000e+00], [ 2.64000000e+00], [ 2.58000000e+00], [ 4.29000000e+00], [ 2.27000000e+00], [ 4.11000000e+00], [ 8.40000000e-01], [ 2.92000000e+00], [ 2.92000000e+00], [ 1.76000000e+00], [ 1.23000000e+00], [ 4.90000000e-01], [ 1.72000000e+00], [ 1.46000000e+00], [ 1.88000000e+00], [ 5.10000000e-01], [ 0.00000000e+00], [ 2.30000000e-01], [ 1.67000000e+00], [ 0.00000000e+00], [ 9.50000000e-01], [ 7.00000000e-01], [ 2.71000000e+00], [ 2.41000000e+00], [ 5.45000000e+00], [ 1.19900000e+01], [ 1.65900000e+01], [ 3.23400000e+01], [ 4.86000000e+01], [ 7.45900000e+01], [ 9.83700000e+01], [ 1.31450000e+02], [ 1.57250000e+02], [ 1.86570000e+02], [ 1.85670000e+02], [ 1.88920000e+02], [ 1.93460000e+02], [ 1.89980000e+02], [ 1.66270000e+02], [ 1.42150000e+02], [ 1.15260000e+02], [ 8.94400000e+01], [ 6.21500000e+01], [ 4.91100000e+01], [ 3.18100000e+01], [ 2.20200000e+01], [ 1.30800000e+01], [ 9.26000000e+00], [ 5.29000000e+00], [ 3.32000000e+00], [ 2.41000000e+00], [ 1.30000000e+00], [ 1.55000000e+00], [ 1.69000000e+00], [ 0.00000000e+00], [ 1.90000000e+00], [ 3.90000000e-01], [ 1.37000000e+00], [ 1.65000000e+00], [ 4.27000000e+00], [ 4.36000000e+00], [ 7.17000000e+00], [ 7.45000000e+00], [ 5.27000000e+00], [ 7.56000000e+00], [ 4.34000000e+00], [ 6.87000000e+00], [ 9.60000000e+00], [ 8.12000000e+00], [ 8.26000000e+00], [ 6.87000000e+00], [ 6.70000000e+00], [ 7.15000000e+00], [ 7.15000000e+00], [ 6.01000000e+00], [ 5.38000000e+00], [ 6.57000000e+00], [ 4.76000000e+00], [ 5.34000000e+00], [ 5.75000000e+00], [ 7.10000000e+00], [ 4.57000000e+00], [ 1.01600000e+01], [ 1.69600000e+01], [ 3.11300000e+01], [ 4.74200000e+01], [ 6.79300000e+01], [ 9.62600000e+01], [ 1.28880000e+02], [ 1.64790000e+02], [ 1.93420000e+02], [ 1.91820000e+02], [ 2.02190000e+02], [ 2.07870000e+02], [ 1.91080000e+02], [ 1.66440000e+02], [ 1.34880000e+02], [ 1.07760000e+02], [ 8.03400000e+01], [ 5.57500000e+01], [ 3.55900000e+01], [ 2.61900000e+01], [ 1.94900000e+01], [ 1.12300000e+01], [ 9.05000000e+00], [ 6.43000000e+00], [ 7.01000000e+00], [ 6.64000000e+00], [ 4.57000000e+00], [ 4.20000000e+00], [ 4.92000000e+00], [ 4.94000000e+00], [ 6.52000000e+00], [ 4.62000000e+00], [ 4.59000000e+00], [ 5.24000000e+00], [ 6.80000000e+00], [ 5.31000000e+00], [ 6.29000000e+00], [ 4.41000000e+00], [ 6.12000000e+00], [ 7.15000000e+00], [ 6.84000000e+00], [ 1.02800000e+01], [ 8.70000000e+00], [ 1.10400000e+01], [ 1.11400000e+01], [ 1.04600000e+01], [ 6.45000000e+00], [ 5.24000000e+00], [ 7.40000000e+00], [ 6.87000000e+00], [ 9.35000000e+00], [ 6.17000000e+00], [ 6.24000000e+00], [ 5.59000000e+00], [ 3.71000000e+00], [ 4.52000000e+00], [ 4.27000000e+00], [ 1.00000000e+01], [ 2.00000000e+01], [ 3.27400000e+01], [ 4.02800000e+01], [ 5.48900000e+01], [ 7.31300000e+01], [ 1.01500000e+02], [ 1.30410000e+02], [ 1.66530000e+02], [ 1.85670000e+02], [ 1.83470000e+02], [ 1.89470000e+02], [ 1.75320000e+02], [ 1.53000000e+02], [ 1.24000000e+02], [ 9.26600000e+01], [ 6.00200000e+01], [ 3.71700000e+01], [ 2.88400000e+01], [ 2.07400000e+01], [ 1.61500000e+01], [ 1.23700000e+01], [ 8.96000000e+00], [ 9.91000000e+00], [ 9.19000000e+00], [ 8.54000000e+00], [ 1.03500000e+01], [ 9.65000000e+00], [ 8.31000000e+00], [ 8.86000000e+00], [ 5.73000000e+00], [ 7.38000000e+00], [ 4.41000000e+00], [ 5.92000000e+00], [ 5.36000000e+00], [ 3.80000000e+00], [ 4.92000000e+00], [ 2.51000000e+00], [ 3.69000000e+00], [ 3.76000000e+00], [ 2.60000000e+00], [ 3.94000000e+00], [ 4.80000000e+00], [ 7.40000000e+00], [ 7.10000000e+00], [ 8.33000000e+00], [ 3.57000000e+00], [ 5.78000000e+00], [ 6.17000000e+00], [ 6.66000000e+00], [ 6.40000000e+00], [ 4.59000000e+00], [ 5.57000000e+00], [ 5.89000000e+00], [ 4.96000000e+00], [ 7.08000000e+00], [ 9.58000000e+00], [ 1.60300000e+01], [ 2.89800000e+01], [ 3.74700000e+01], [ 4.70500000e+01], [ 6.58000000e+01], [ 7.91400000e+01], [ 9.78600000e+01], [ 1.32150000e+02], [ 1.43070000e+02], [ 1.44420000e+02], [ 1.33890000e+02], [ 1.23240000e+02], [ 1.12130000e+02], [ 9.55400000e+01], [ 6.61700000e+01], [ 4.82600000e+01], [ 3.18300000e+01], [ 2.68900000e+01], [ 1.82600000e+01], [ 1.42000000e+01], [ 1.04900000e+01], [ 1.10900000e+01], [ 9.28000000e+00], [ 7.17000000e+00], [ 9.47000000e+00], [ 8.84000000e+00], [ 9.12000000e+00], [ 9.07000000e+00], [ 7.28000000e+00], [ 1.03000000e+01], [ 9.74000000e+00], [ 8.51000000e+00], [ 9.02000000e+00], [ 5.71000000e+00], [ 8.65000000e+00], [ 9.05000000e+00], [ 4.08000000e+00], [ 7.10000000e+00], [ 4.50000000e+00], [ 4.45000000e+00], [ 3.62000000e+00], [ 5.17000000e+00], [ 6.17000000e+00], [ 4.20000000e+00], [ 5.57000000e+00], [ 9.50000000e-01], [ 5.15000000e+00], [ 3.02000000e+00], [ 2.58000000e+00], [ 2.76000000e+00], [ 2.81000000e+00], [ 4.59000000e+00], [ 3.22000000e+00], [ 4.57000000e+00], [ 4.06000000e+00], [ 1.67000000e+01], [ 2.84200000e+01], [ 4.03200000e+01], [ 5.13400000e+01], [ 6.14600000e+01], [ 7.26400000e+01], [ 8.06700000e+01], [ 1.00970000e+02], [ 1.16020000e+02], [ 1.10130000e+02], [ 1.04790000e+02], [ 9.14100000e+01], [ 7.64400000e+01], [ 6.24100000e+01], [ 5.10900000e+01], [ 3.87000000e+01], [ 3.18500000e+01], [ 2.00200000e+01], [ 1.72100000e+01], [ 1.35700000e+01], [ 1.02500000e+01], [ 8.51000000e+00], [ 1.02500000e+01], [ 7.93000000e+00], [ 7.54000000e+00], [ 8.75000000e+00], [ 9.19000000e+00], [ 7.70000000e+00], [ 6.68000000e+00], [ 7.66000000e+00], [ 8.44000000e+00], [ 5.38000000e+00], [ 8.00000000e+00], [ 1.53000000e+00], [ 1.44000000e+00], [ 2.48000000e+00], [ 2.32000000e+00], [ 1.90000000e-01], [ 4.59000000e+00], [ 6.01000000e+00], [ 6.22000000e+00], [ 1.06000000e+01], [ 9.60000000e+00], [ 9.21000000e+00], [ 8.89000000e+00], [ 1.14400000e+01], [ 8.05000000e+00], [ 1.32900000e+01], [ 8.54000000e+00], [ 8.75000000e+00], [ 4.11000000e+00], [ 7.40000000e+00], [ 4.50000000e+00], [ 8.98000000e+00], [ 4.64000000e+00], [ 5.80000000e+00], [ 4.22000000e+00], [ 8.82000000e+00], [ 1.36900000e+01], [ 2.51300000e+01], [ 2.87200000e+01], [ 3.26200000e+01], [ 3.50600000e+01], [ 4.06000000e+01], [ 5.13400000e+01], [ 5.24600000e+01], [ 7.27600000e+01], [ 5.97600000e+01], [ 4.51500000e+01], [ 3.50600000e+01], [ 3.06500000e+01], [ 3.00700000e+01], [ 2.65600000e+01], [ 2.31800000e+01], [ 1.65200000e+01], [ 1.51700000e+01], [ 1.00000000e+01], [ 1.05300000e+01], [ 1.04900000e+01], [ 7.54000000e+00], [ 7.35000000e+00], [ 8.84000000e+00], [ 7.52000000e+00], [ 9.65000000e+00], [ 9.47000000e+00], [ 9.67000000e+00], [ 7.47000000e+00], [ 1.08100000e+01], [ 6.80000000e+00], [ 5.71000000e+00], [ 5.87000000e+00], [ 4.92000000e+00], [ 9.84000000e+00], [ 6.59000000e+00], [ 7.12000000e+00], [ 1.46000000e+00], [ 2.80000000e-01], [ 4.62000000e+00], [ 4.27000000e+00], [ 2.13000000e+00], [ 3.06000000e+00], [ 1.74000000e+00], [ 1.02000000e+00], [ 0.00000000e+00], [ 3.04000000e+00], [ 1.48000000e+00], [ 2.13000000e+00], [ 1.21000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 0.00000000e+00], [ 9.30000000e-01], [ 1.00000000e+00], [ 1.62000000e+00], [ 1.16000000e+00], [ 2.00000000e+00], [ 3.18000000e+00], [ 9.51000000e+00], [ 1.60500000e+01], [ 1.50800000e+01], [ 1.84700000e+01], [ 1.93700000e+01], [ 1.98400000e+01], [ 2.76500000e+01], [ 1.97200000e+01], [ 1.86100000e+01], [ 1.26000000e+01], [ 1.07400000e+01], [ 1.37600000e+01], [ 1.01600000e+01], [ 2.18000000e+00], [ 4.66000000e+00], [ 3.55000000e+00], [ 5.13000000e+00], [ 1.88000000e+00], [ 4.41000000e+00], [ 4.08000000e+00], [ 3.00000000e-01], [ 5.36000000e+00], [ 5.45000000e+00], [ 0.00000000e+00], [ 4.48000000e+00], [ 2.51000000e+00], [ 4.71000000e+00], [ 3.36000000e+00], [ 4.13000000e+00], [ 3.99000000e+00], [ 4.50000000e+00], [ 4.41000000e+00], [ 7.68000000e+00], [ 6.73000000e+00], [ 7.08000000e+00]])
m = GPy.models.GPRegression(data['X'],data['Y'],kern)
The initial value of $a$, $b$ and $c$ are 1. For these choices of covariance function parameters, we can plot the random field of $f$ and $y$ separately.
leng = data['X'].shape[0]
m.plot(fixed_inputs=[(2,0)], which_data_rows = slice(0,leng*2/2))
m.plot(fixed_inputs=[(2,1)], which_data_rows = slice(leng*2/2,leng*2))
{'contour': <matplotlib.contour.QuadContourSet instance at 0x10f9a32d8>, 'dataplot': <matplotlib.collections.PathCollection at 0x10f9b97d0>}
Now we optimize the model, this will take a few minutes.
m.optimize(messages=True)
I F Scale |g| 0003 7.614997e+04 2.500000e-01 2.504914e+09 0008 1.133084e+04 7.812500e-03 8.476310e+07 0016 3.600478e+03 3.051758e-05 6.694209e+04 0022 2.830839e+03 4.768372e-07 5.181429e+02 0031 2.768527e+03 9.313226e-10 7.879161e+01 0053 2.747934e+03 1.000000e-15 7.722679e+00 0069 2.746164e+03 1.000000e-15 5.889059e-02 0082 2.746153e+03 1.000000e-15 7.165007e-03 0089 2.746153e+03 1.000000e-15 1.479437e-04 0089 2.746153e+03 1.000000e-15 1.479437e-04 converged - relative reduction in objective
After optimization, the estimated value of $a$, $b$ and $c$ can be printed.
print m
GP_regression. | Value | Constraint | Prior | Tied to ode_st.a | 0.438768542606 | +ve | | ode_st.b | 10.0747613738 | +ve | | ode_st.c | 0.620385966622 | +ve | | ode_st.variance_Yt | 32.8716774332 | +ve | | ode_st.variance_Yx | 32.8716774332 | +ve | | ode_st.lengthscale_Yt | 18.3997301624 | +ve | | ode_st.lengthscale_Yx | 14.3560437898 | +ve | | Gaussian_noise.variance | 3.26998881271 | +ve | |
The plot of the random fields are plotted below.
m.plot(fixed_inputs=[(2,0)], which_data_rows = slice(0,leng*2/2))
pb.savefig("gene.pdf")
m.plot(fixed_inputs=[(2,1)], which_data_rows = slice(leng*2/2,leng*2))
pb.savefig("protein.pdf")
In his diploma thesis Kolja Becker estimated the parameters of the partial differential equation where $a = 0.159$, $b =12.77$ and $c =0.983$. The results from the GP approach is $a =0.439$, $b= 10.06$ and $c=0.62$. The GP results are different but they are still within the range of the confidence interval defined in Becker's paper. Two further issues may lead to the difference. One reason could be that the original partial differential equation had a delay parameter $\tau$ for the mRNA ($f$). In our GP model, we did not include this parameter. However, since the Protein-mRNA partial differential equation is linear and Becker's estimate of $\tau$ is small comparing with the time step of the data. The delay impact should not be too big. The estimation algorithm used in Becker's thesis is based on least squares optimization. The GP approach considers protein and mRNA as a nonlinear function of space and time. And treat the PDE as a linear function to link them. The different modelling methods could also lead to different estimation of the parameters.
Finally, we didn't perform a sensitivity analysis in the above notebook. One thing we can do next is to run a Hamiltonian Monte-Carlo sampler on the model to form error bars for our own analysis. It may be that they are not very well determined given the data.