Simulations

In addition to inference, the Gaussian process modelling framework allows simulations and predictions of the underlying process. This is done by conditioning the Gaussian process on some observations and then sampling from the conditioned distribution.

Sampling from the Gaussian process

Assuming a SingleBendingPowerLaw we can draw realisations from the Gaussian process with a mean. First, we define the power spectral density function and plot it.

using Random
using Plots
using Pioran

rng = MersenneTwister(1234)
# Define power spectral density function
𝓟 = SingleBendingPowerLaw(0.3,1e-2,2.9)
f_min, f_max = 1e-3, 1e3
f0,fM = f_min/20.,f_max*20.
f = 10 .^ range(log10(f0), log10(fM), length=1000)
plot(f, 𝓟(f), xlabel="Frequency (day^-1)",ylabel="Power Spectral Density",legend=false,framestyle = :box,xscale=:log10,yscale=:log10,lw=2)
Example block output

Then we approximate it to form a covariance function. We then build the Gaussian process and draw realisations from it using rand.

variance = 15.2
# Approximation of the PSD to form a covariance function
𝓡 = approx(𝓟, f_min, f_max, 20, variance, basis_function="SHO")
μ = 1.3
# Build the GP
f = ScalableGP(μ, 𝓡) # Define the GP

T = 450 # days
t = range(0, stop=T, length=1000)
σ² = ones(length(t))*0.25

fx = f(t, σ²) # Gaussian process
realisations = [rand(rng,fx) for _ in 1:5]
plot(t, realisations, label="Realisations", xlabel="Time [days]", framestyle = :box,ylabel="Value")
Example block output
Note

Sampling from a Gaussian process built with semi-separable covariance functions is very efficient. The time complexity is O(N) where N is the number of data points.

Conditioning the Gaussian process

We can compute the conditioned or posterior distribution of the Gaussian process given some observations. Let's use a subset of the realisations to condition the Gaussian process and then sample from the conditioned distribution.

using StatsBase
idx = sort(sample(1:length(t), 50, replace = false));
t_obs = t[idx]
y_obs = realisations[1][idx]
yerr = 0.25*ones(length(t_obs))
fx = f(t_obs, yerr.^2) # Gaussian process
AbstractGPs.FiniteGP{ScalableGP{AbstractGPs.GP{AbstractGPs.ConstMean{Float64}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}, Vector{Float64}, LinearAlgebra.Diagonal{Float64, Vector{Float64}}}(
f: ScalableGP{AbstractGPs.GP{AbstractGPs.ConstMean{Float64}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}(AbstractGPs.GP{AbstractGPs.ConstMean{Float64}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}(AbstractGPs.ConstMean{Float64}(1.3), Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}(Celerite[Celerite(0.21179581122061306, 0.21179581122061306, 0.0002221441469079183, 0.0002221441469079183), Celerite(0.38542244831627853, 0.38542244831627853, 0.0006300497947253256, 0.0006300497947253256), Celerite(0.9038882037020084, 0.9038882037020084, 0.0017869601758986303, 0.0017869601758986303), Celerite(1.6067166095810141, 1.6067166095810141, 0.0050682131745472155, 0.0050682131745472155), Celerite(4.964666291074675, 4.964666291074675, 0.014374570362059993, 0.014374570362059993), Celerite(22.755812872815056, 22.755812872815056, 0.04076945187142277, 0.04076945187142277), Celerite(5.378025805869529, 5.378025805869529, 0.11563115724719714, 0.11563115724719714), Celerite(0.5417886639897328, 0.5417886639897328, 0.32795546451037977, 0.32795546451037977), Celerite(0.08070461417658467, 0.08070461417658467, 0.9301540282286335, 0.9301540282286335), Celerite(0.010848075938726597, 0.010848075938726597, 2.638121970376775, 2.638121970376775), Celerite(0.0015090806302273192, 0.0015090806302273192, 7.482295748198312, 7.482295748198312), Celerite(0.0002076692891940062, 0.0002076692891940062, 21.221440969050757, 21.221440969050757), Celerite(2.8676529647079397e-5, 2.8676529647079397e-5, 60.18868699641365, 60.18868699641365), Celerite(3.955543296634202e-6, 3.955543296634202e-6, 170.70839099171252, 170.70839099171252), Celerite(5.458037196298439e-7, 5.458037196298439e-7, 484.16664674402534, 484.16664674402534), Celerite(7.53042020169437e-8, 7.53042020169437e-8, 1373.203393562152, 1373.203393562152), Celerite(1.0389899331678754e-8, 1.0389899331678754e-8, 3894.70768540476, 3894.70768540476), Celerite(1.4338973791019396e-9, 1.4338973791019396e-9, 11046.249977144691, 11046.249977144691), Celerite(1.9677308815711412e-10, 1.9677308815711412e-10, 31329.601195702533, 31329.601195702533), Celerite(3.012647957419647e-11, 3.012647957419647e-11, 88857.65876316732, 88857.65876316732)], [0.21179581122061306, 0.38542244831627853, 0.9038882037020084, 1.6067166095810141, 4.964666291074675, 22.755812872815056, 5.378025805869529, 0.5417886639897328, 0.08070461417658467, 0.010848075938726597, 0.0015090806302273192, 0.0002076692891940062, 2.8676529647079397e-5, 3.955543296634202e-6, 5.458037196298439e-7, 7.53042020169437e-8, 1.0389899331678754e-8, 1.4338973791019396e-9, 1.9677308815711412e-10, 3.012647957419647e-11], [0.21179581122061306, 0.38542244831627853, 0.9038882037020084, 1.6067166095810141, 4.964666291074675, 22.755812872815056, 5.378025805869529, 0.5417886639897328, 0.08070461417658467, 0.010848075938726597, 0.0015090806302273192, 0.0002076692891940062, 2.8676529647079397e-5, 3.955543296634202e-6, 5.458037196298439e-7, 7.53042020169437e-8, 1.0389899331678754e-8, 1.4338973791019396e-9, 1.9677308815711412e-10, 3.012647957419647e-11], [0.0002221441469079183, 0.0006300497947253256, 0.0017869601758986303, 0.0050682131745472155, 0.014374570362059993, 0.04076945187142277, 0.11563115724719714, 0.32795546451037977, 0.9301540282286335, 2.638121970376775, 7.482295748198312, 21.221440969050757, 60.18868699641365, 170.70839099171252, 484.16664674402534, 1373.203393562152, 3894.70768540476, 11046.249977144691, 31329.601195702533, 88857.65876316732], [0.0002221441469079183, 0.0006300497947253256, 0.0017869601758986303, 0.0050682131745472155, 0.014374570362059993, 0.04076945187142277, 0.11563115724719714, 0.32795546451037977, 0.9301540282286335, 2.638121970376775, 7.482295748198312, 21.221440969050757, 60.18868699641365, 170.70839099171252, 484.16664674402534, 1373.203393562152, 3894.70768540476, 11046.249977144691, 31329.601195702533, 88857.65876316732])), Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}(Celerite[Celerite(0.21179581122061306, 0.21179581122061306, 0.0002221441469079183, 0.0002221441469079183), Celerite(0.38542244831627853, 0.38542244831627853, 0.0006300497947253256, 0.0006300497947253256), Celerite(0.9038882037020084, 0.9038882037020084, 0.0017869601758986303, 0.0017869601758986303), Celerite(1.6067166095810141, 1.6067166095810141, 0.0050682131745472155, 0.0050682131745472155), Celerite(4.964666291074675, 4.964666291074675, 0.014374570362059993, 0.014374570362059993), Celerite(22.755812872815056, 22.755812872815056, 0.04076945187142277, 0.04076945187142277), Celerite(5.378025805869529, 5.378025805869529, 0.11563115724719714, 0.11563115724719714), Celerite(0.5417886639897328, 0.5417886639897328, 0.32795546451037977, 0.32795546451037977), Celerite(0.08070461417658467, 0.08070461417658467, 0.9301540282286335, 0.9301540282286335), Celerite(0.010848075938726597, 0.010848075938726597, 2.638121970376775, 2.638121970376775), Celerite(0.0015090806302273192, 0.0015090806302273192, 7.482295748198312, 7.482295748198312), Celerite(0.0002076692891940062, 0.0002076692891940062, 21.221440969050757, 21.221440969050757), Celerite(2.8676529647079397e-5, 2.8676529647079397e-5, 60.18868699641365, 60.18868699641365), Celerite(3.955543296634202e-6, 3.955543296634202e-6, 170.70839099171252, 170.70839099171252), Celerite(5.458037196298439e-7, 5.458037196298439e-7, 484.16664674402534, 484.16664674402534), Celerite(7.53042020169437e-8, 7.53042020169437e-8, 1373.203393562152, 1373.203393562152), Celerite(1.0389899331678754e-8, 1.0389899331678754e-8, 3894.70768540476, 3894.70768540476), Celerite(1.4338973791019396e-9, 1.4338973791019396e-9, 11046.249977144691, 11046.249977144691), Celerite(1.9677308815711412e-10, 1.9677308815711412e-10, 31329.601195702533, 31329.601195702533), Celerite(3.012647957419647e-11, 3.012647957419647e-11, 88857.65876316732, 88857.65876316732)], [0.21179581122061306, 0.38542244831627853, 0.9038882037020084, 1.6067166095810141, 4.964666291074675, 22.755812872815056, 5.378025805869529, 0.5417886639897328, 0.08070461417658467, 0.010848075938726597, 0.0015090806302273192, 0.0002076692891940062, 2.8676529647079397e-5, 3.955543296634202e-6, 5.458037196298439e-7, 7.53042020169437e-8, 1.0389899331678754e-8, 1.4338973791019396e-9, 1.9677308815711412e-10, 3.012647957419647e-11], [0.21179581122061306, 0.38542244831627853, 0.9038882037020084, 1.6067166095810141, 4.964666291074675, 22.755812872815056, 5.378025805869529, 0.5417886639897328, 0.08070461417658467, 0.010848075938726597, 0.0015090806302273192, 0.0002076692891940062, 2.8676529647079397e-5, 3.955543296634202e-6, 5.458037196298439e-7, 7.53042020169437e-8, 1.0389899331678754e-8, 1.4338973791019396e-9, 1.9677308815711412e-10, 3.012647957419647e-11], [0.0002221441469079183, 0.0006300497947253256, 0.0017869601758986303, 0.0050682131745472155, 0.014374570362059993, 0.04076945187142277, 0.11563115724719714, 0.32795546451037977, 0.9301540282286335, 2.638121970376775, 7.482295748198312, 21.221440969050757, 60.18868699641365, 170.70839099171252, 484.16664674402534, 1373.203393562152, 3894.70768540476, 11046.249977144691, 31329.601195702533, 88857.65876316732], [0.0002221441469079183, 0.0006300497947253256, 0.0017869601758986303, 0.0050682131745472155, 0.014374570362059993, 0.04076945187142277, 0.11563115724719714, 0.32795546451037977, 0.9301540282286335, 2.638121970376775, 7.482295748198312, 21.221440969050757, 60.18868699641365, 170.70839099171252, 484.16664674402534, 1373.203393562152, 3894.70768540476, 11046.249977144691, 31329.601195702533, 88857.65876316732]))
x: [1.3513513513513513, 1.8018018018018018, 15.765765765765765, 44.14414414414414, 54.95495495495496, 58.55855855855856, 62.16216216216216, 85.58558558558559, 98.64864864864865, 100.45045045045045  …  396.39639639639637, 400.9009009009009, 401.8018018018018, 403.15315315315314, 415.31531531531533, 418.4684684684685, 422.52252252252254, 422.97297297297297, 425.22522522522524, 444.14414414414415]
Σy: [0.0625 0.0 … 0.0 0.0; 0.0 0.0625 … 0.0 0.0; … ; 0.0 0.0 … 0.0625 0.0; 0.0 0.0 … 0.0 0.0625]
)

We can compute the posterior distribution of the Gaussian process given the observations.

fp = posterior(fx, y_obs) # Posterior distribution
Pioran.PosteriorGP{AbstractGPs.FiniteGP{ScalableGP{AbstractGPs.GP{AbstractGPs.ConstMean{Float64}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}, Vector{Float64}, LinearAlgebra.Diagonal{Float64, Vector{Float64}}}, Vector{Float64}}(AbstractGPs.FiniteGP{ScalableGP{AbstractGPs.GP{AbstractGPs.ConstMean{Float64}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}, Vector{Float64}, LinearAlgebra.Diagonal{Float64, Vector{Float64}}}(
f: ScalableGP{AbstractGPs.GP{AbstractGPs.ConstMean{Float64}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}(AbstractGPs.GP{AbstractGPs.ConstMean{Float64}, Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}}(AbstractGPs.ConstMean{Float64}(1.3), Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}(Celerite[Celerite(0.21179581122061306, 0.21179581122061306, 0.0002221441469079183, 0.0002221441469079183), Celerite(0.38542244831627853, 0.38542244831627853, 0.0006300497947253256, 0.0006300497947253256), Celerite(0.9038882037020084, 0.9038882037020084, 0.0017869601758986303, 0.0017869601758986303), Celerite(1.6067166095810141, 1.6067166095810141, 0.0050682131745472155, 0.0050682131745472155), Celerite(4.964666291074675, 4.964666291074675, 0.014374570362059993, 0.014374570362059993), Celerite(22.755812872815056, 22.755812872815056, 0.04076945187142277, 0.04076945187142277), Celerite(5.378025805869529, 5.378025805869529, 0.11563115724719714, 0.11563115724719714), Celerite(0.5417886639897328, 0.5417886639897328, 0.32795546451037977, 0.32795546451037977), Celerite(0.08070461417658467, 0.08070461417658467, 0.9301540282286335, 0.9301540282286335), Celerite(0.010848075938726597, 0.010848075938726597, 2.638121970376775, 2.638121970376775), Celerite(0.0015090806302273192, 0.0015090806302273192, 7.482295748198312, 7.482295748198312), Celerite(0.0002076692891940062, 0.0002076692891940062, 21.221440969050757, 21.221440969050757), Celerite(2.8676529647079397e-5, 2.8676529647079397e-5, 60.18868699641365, 60.18868699641365), Celerite(3.955543296634202e-6, 3.955543296634202e-6, 170.70839099171252, 170.70839099171252), Celerite(5.458037196298439e-7, 5.458037196298439e-7, 484.16664674402534, 484.16664674402534), Celerite(7.53042020169437e-8, 7.53042020169437e-8, 1373.203393562152, 1373.203393562152), Celerite(1.0389899331678754e-8, 1.0389899331678754e-8, 3894.70768540476, 3894.70768540476), Celerite(1.4338973791019396e-9, 1.4338973791019396e-9, 11046.249977144691, 11046.249977144691), Celerite(1.9677308815711412e-10, 1.9677308815711412e-10, 31329.601195702533, 31329.601195702533), Celerite(3.012647957419647e-11, 3.012647957419647e-11, 88857.65876316732, 88857.65876316732)], [0.21179581122061306, 0.38542244831627853, 0.9038882037020084, 1.6067166095810141, 4.964666291074675, 22.755812872815056, 5.378025805869529, 0.5417886639897328, 0.08070461417658467, 0.010848075938726597, 0.0015090806302273192, 0.0002076692891940062, 2.8676529647079397e-5, 3.955543296634202e-6, 5.458037196298439e-7, 7.53042020169437e-8, 1.0389899331678754e-8, 1.4338973791019396e-9, 1.9677308815711412e-10, 3.012647957419647e-11], [0.21179581122061306, 0.38542244831627853, 0.9038882037020084, 1.6067166095810141, 4.964666291074675, 22.755812872815056, 5.378025805869529, 0.5417886639897328, 0.08070461417658467, 0.010848075938726597, 0.0015090806302273192, 0.0002076692891940062, 2.8676529647079397e-5, 3.955543296634202e-6, 5.458037196298439e-7, 7.53042020169437e-8, 1.0389899331678754e-8, 1.4338973791019396e-9, 1.9677308815711412e-10, 3.012647957419647e-11], [0.0002221441469079183, 0.0006300497947253256, 0.0017869601758986303, 0.0050682131745472155, 0.014374570362059993, 0.04076945187142277, 0.11563115724719714, 0.32795546451037977, 0.9301540282286335, 2.638121970376775, 7.482295748198312, 21.221440969050757, 60.18868699641365, 170.70839099171252, 484.16664674402534, 1373.203393562152, 3894.70768540476, 11046.249977144691, 31329.601195702533, 88857.65876316732], [0.0002221441469079183, 0.0006300497947253256, 0.0017869601758986303, 0.0050682131745472155, 0.014374570362059993, 0.04076945187142277, 0.11563115724719714, 0.32795546451037977, 0.9301540282286335, 2.638121970376775, 7.482295748198312, 21.221440969050757, 60.18868699641365, 170.70839099171252, 484.16664674402534, 1373.203393562152, 3894.70768540476, 11046.249977144691, 31329.601195702533, 88857.65876316732])), Pioran.SumOfCelerite{StructArrays.StructArray{<:Pioran.SemiSeparable}}(Celerite[Celerite(0.21179581122061306, 0.21179581122061306, 0.0002221441469079183, 0.0002221441469079183), Celerite(0.38542244831627853, 0.38542244831627853, 0.0006300497947253256, 0.0006300497947253256), Celerite(0.9038882037020084, 0.9038882037020084, 0.0017869601758986303, 0.0017869601758986303), Celerite(1.6067166095810141, 1.6067166095810141, 0.0050682131745472155, 0.0050682131745472155), Celerite(4.964666291074675, 4.964666291074675, 0.014374570362059993, 0.014374570362059993), Celerite(22.755812872815056, 22.755812872815056, 0.04076945187142277, 0.04076945187142277), Celerite(5.378025805869529, 5.378025805869529, 0.11563115724719714, 0.11563115724719714), Celerite(0.5417886639897328, 0.5417886639897328, 0.32795546451037977, 0.32795546451037977), Celerite(0.08070461417658467, 0.08070461417658467, 0.9301540282286335, 0.9301540282286335), Celerite(0.010848075938726597, 0.010848075938726597, 2.638121970376775, 2.638121970376775), Celerite(0.0015090806302273192, 0.0015090806302273192, 7.482295748198312, 7.482295748198312), Celerite(0.0002076692891940062, 0.0002076692891940062, 21.221440969050757, 21.221440969050757), Celerite(2.8676529647079397e-5, 2.8676529647079397e-5, 60.18868699641365, 60.18868699641365), Celerite(3.955543296634202e-6, 3.955543296634202e-6, 170.70839099171252, 170.70839099171252), Celerite(5.458037196298439e-7, 5.458037196298439e-7, 484.16664674402534, 484.16664674402534), Celerite(7.53042020169437e-8, 7.53042020169437e-8, 1373.203393562152, 1373.203393562152), Celerite(1.0389899331678754e-8, 1.0389899331678754e-8, 3894.70768540476, 3894.70768540476), Celerite(1.4338973791019396e-9, 1.4338973791019396e-9, 11046.249977144691, 11046.249977144691), Celerite(1.9677308815711412e-10, 1.9677308815711412e-10, 31329.601195702533, 31329.601195702533), Celerite(3.012647957419647e-11, 3.012647957419647e-11, 88857.65876316732, 88857.65876316732)], [0.21179581122061306, 0.38542244831627853, 0.9038882037020084, 1.6067166095810141, 4.964666291074675, 22.755812872815056, 5.378025805869529, 0.5417886639897328, 0.08070461417658467, 0.010848075938726597, 0.0015090806302273192, 0.0002076692891940062, 2.8676529647079397e-5, 3.955543296634202e-6, 5.458037196298439e-7, 7.53042020169437e-8, 1.0389899331678754e-8, 1.4338973791019396e-9, 1.9677308815711412e-10, 3.012647957419647e-11], [0.21179581122061306, 0.38542244831627853, 0.9038882037020084, 1.6067166095810141, 4.964666291074675, 22.755812872815056, 5.378025805869529, 0.5417886639897328, 0.08070461417658467, 0.010848075938726597, 0.0015090806302273192, 0.0002076692891940062, 2.8676529647079397e-5, 3.955543296634202e-6, 5.458037196298439e-7, 7.53042020169437e-8, 1.0389899331678754e-8, 1.4338973791019396e-9, 1.9677308815711412e-10, 3.012647957419647e-11], [0.0002221441469079183, 0.0006300497947253256, 0.0017869601758986303, 0.0050682131745472155, 0.014374570362059993, 0.04076945187142277, 0.11563115724719714, 0.32795546451037977, 0.9301540282286335, 2.638121970376775, 7.482295748198312, 21.221440969050757, 60.18868699641365, 170.70839099171252, 484.16664674402534, 1373.203393562152, 3894.70768540476, 11046.249977144691, 31329.601195702533, 88857.65876316732], [0.0002221441469079183, 0.0006300497947253256, 0.0017869601758986303, 0.0050682131745472155, 0.014374570362059993, 0.04076945187142277, 0.11563115724719714, 0.32795546451037977, 0.9301540282286335, 2.638121970376775, 7.482295748198312, 21.221440969050757, 60.18868699641365, 170.70839099171252, 484.16664674402534, 1373.203393562152, 3894.70768540476, 11046.249977144691, 31329.601195702533, 88857.65876316732]))
x: [1.3513513513513513, 1.8018018018018018, 15.765765765765765, 44.14414414414414, 54.95495495495496, 58.55855855855856, 62.16216216216216, 85.58558558558559, 98.64864864864865, 100.45045045045045  …  396.39639639639637, 400.9009009009009, 401.8018018018018, 403.15315315315314, 415.31531531531533, 418.4684684684685, 422.52252252252254, 422.97297297297297, 425.22522522522524, 444.14414414414415]
Σy: [0.0625 0.0 … 0.0 0.0; 0.0 0.0625 … 0.0 0.0; … ; 0.0 0.0 … 0.0625 0.0; 0.0 0.0 … 0.0 0.0625]
)
, [5.049625421463632, 5.788046102707981, 3.253813167202233, -6.434451728276888, -2.153775909938405, 1.0264990265324576, 4.751459874105995, 4.31344119799946, 9.309942440653156, 9.295111340218552  …  2.5120861273426445, 0.07640695112513551, 0.7411102972888532, -0.19778064914068172, -8.3981535447054, -8.348219512676872, -8.300576243904722, -8.754491553645984, -8.714588428544848, 0.7835131871531228])

The mean and standard deviation of this distribution, can be computed using the mean and std functions. The posterior covariance matrix can be computed using the cov function.

m = mean(fp,t);
s = std(fp,t);
1000-element Vector{Float64}:
 0.8128885705176537
 0.6108530914946403
 0.4006365343352816
 0.2120747150299907
 0.21080303025566322
 0.3912678115025849
 0.5858127957039118
 0.7648927567581641
 0.9277798746859259
 1.0747445973130927
 ⋮
 1.2244987023111926
 1.405812711700324
 1.5784983878280654
 1.7434754428873536
 1.9013910043560451
 2.052690740509442
 2.1977093081230383
 2.336733400913537
 2.4700354318778466

We can plot the realisations, the observations and the posterior distribution.

plot(t, realisations[1], label="Realisation", xlabel="Time [days]", framestyle = :box,ylabel="Value")
plot!(t_obs, y_obs,yerr=yerr, label="Observations", seriestype=:scatter)
plot!(t,m, ribbon=s,label="Posterior distribution", lw=2)
plot!(t,m,ribbon=2*s, label=nothing)
Example block output
Note

The computation of the mean of the distribution is very efficient. The time complexity is O(N) where N is the number of data points. However, the computation of the covariance is very inefficient as the posterior covariance is not semi-separable.

Sampling from the conditioned distribution

We can draw realisations from the conditioned distribution using rand.

samples_cond = rand(rng,fp,t,5);
1000×5 Matrix{Float64}:
 5.53396    5.10929     6.40639   6.23271  4.27956
 5.71575    4.95406     6.0142    6.00398  4.52487
 5.63491    4.84978     5.81612   5.73105  4.71497
 5.63823    4.9363      5.81306   5.62295  5.20891
 5.39262    5.44513     5.93852   5.44746  5.3652
 5.05934    5.54439     5.62429   5.56458  4.9885
 4.78084    5.68507     5.5569    5.20994  4.46708
 4.51815    5.60471     5.6951    4.94701  4.20013
 4.01389    5.99095     6.0906    4.56868  3.97646
 3.40165    6.32814     5.91974   4.38783  3.47988
 ⋮                                         
 0.189211   1.31306    -0.315965  2.20256  2.43656
 0.498496   1.01516    -0.329997  2.28257  2.53727
 0.82524    0.245869   -0.430343  2.72871  2.77738
 0.994561   0.0734135  -0.297912  3.10073  2.47235
 1.38801   -0.199538   -0.313239  2.93133  1.75012
 1.75062   -0.584455   -0.355962  3.10093  1.47234
 2.12952   -0.374646   -0.299121  3.11443  1.27404
 2.65848   -0.501399   -0.178592  3.26339  1.48661
 3.0586    -0.472663   -0.842167  3.62231  1.34132

We can plot the realisations, the observations and the posterior distribution.

plot(t_obs, y_obs,yerr=yerr, label="Observations", seriestype=:scatter, xlabel="Time [days]", framestyle = :box,ylabel="Value",color=:black,lw=2)
plot!(t, samples_cond, label=nothing)
Example block output
Note

Sampling from the conditioned Gaussian process is very inefficient as the posterior distribution is not semi-separable.