Smoothed Collocation

Smoothed collocation, also referred to as the two-stage method, allows for fitting differential equations to time series data without relying on a numerical differential equation solver by building a smoothed collocating polynomial and using this to estimate the true (u',u) pairs, at which point u'-f(u,p,t) can be directly estimated as a loss to determine the correct parameters p. This method can be extremely fast and robust to noise, though, because it does not accumulate through time, is not as exact as other methods.

Note

This is one of many methods for calculating the collocation coefficients for the training process. For a more comprehensive set of collocation methods, see JuliaSimModelOptimizer.

DiffEqFlux.CollocationKernel — Type
CollocationKernel

Abstract interface for kernels used by collocate_data when estimating smoothed state and derivative values from sampled time-series data.

Interface

Concrete kernels are immutable marker types with zero fields. A kernel implementation must provide DiffEqFlux.calckernel(kernel, t) for an offset t, returning the kernel weight at that offset with the same numeric type as t when possible.

Compact-support kernels may implement DiffEqFlux.calckernel(kernel, t, abs(t)); the generic two-argument method handles the support check and calls the three-argument method only when abs(t) <= 1.

Implementations

The public kernel choices are EpanechnikovKernel, UniformKernel, TriangularKernel, QuarticKernel, TriweightKernel, TricubeKernel, GaussianKernel, CosineKernel, LogisticKernel, SigmoidKernel, and SilvermanKernel.

Examples

using DiffEqFlux

kernel = TriangularKernel()
du, u = collocate_data(rand(2, 10), range(0, 1; length = 10), kernel)
source
DiffEqFlux.collocate_data — Function
u′, u = collocate_data(data, tpoints, kernel = TriangularKernel(), bandwidth=nothing)
u′, u = collocate_data(data, tpoints, tpoints_sample, interp, args...)

Computes a non-parametrically smoothed estimate of u' and u given the data, where each column is a snapshot of the timeseries at tpoints[i].

Arguments

  • data: Array of observed state values. For matrix data, each column is one observation at the corresponding entry of tpoints.
  • tpoints: Sample times for data.
  • kernel: A CollocationKernel used for local regression smoothing. Defaults to TriangularKernel.
  • bandwidth: Smoothing bandwidth. If nothing, a rule-of-thumb bandwidth is used.
  • tpoints_sample: Time points where interpolation-based collocation should be sampled.
  • interp: DataInterpolations.jl interpolation constructor used by the extension method.
  • args...: Additional positional arguments forwarded to the interpolation constructor.

Returns

  • u′: Estimated time derivatives at the requested points.
  • u: Smoothed state estimates at the requested points.

Kernel Choices

The following kernel constructors are provided:

  • EpanechnikovKernel
  • UniformKernel
  • TriangularKernel
  • QuarticKernel
  • TriweightKernel
  • TricubeKernel
  • GaussianKernel
  • CosineKernel
  • LogisticKernel
  • SigmoidKernel
  • SilvermanKernel

https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2631937/

Additionally, we can use interpolation methods from DataInterpolations.jl to generate data from intermediate timesteps. In this case, pass any of the methods like QuadraticInterpolation as interp, and the timestamps to sample from as tpoints_sample.

Examples

using DiffEqFlux

tpoints = range(0, 1; length = 20)
data = reduce(hcat, ([sin(t), cos(t)] for t in tpoints))
du, u = collocate_data(data, tpoints, EpanechnikovKernel())
source

Developer Interface

Custom kernels extend the following internal developer interface. It is intended for packages implementing additional CollocationKernel types, not for ordinary collocation calls.

DiffEqFlux.calckernel — Function
DiffEqFlux.calckernel(kernel, t)

Evaluate a collocation kernel at the normalized offset t. This is the developer extension point for custom CollocationKernel implementations.

Rules

  • Implement either the two-argument method for a kernel with non-compact support, or the three-argument method calckernel(kernel, t, abs_t) for a kernel supported on [-1, 1].
  • Return a scalar with a numeric type compatible with t.
  • The generic two-argument method applies the compact-support check before calling the three-argument method.
source

Kernel Choice

Note that the kernel choices of DataInterpolations.jl, such as CubicSpline(), are exact, i.e. go through the data points, while the smoothed kernels are regression splines. Thus CubicSpline() is preferred if the data is not too noisy or is relatively sparse. If data is sparse and very noisy, a BSpline() can be the best regression spline, otherwise one of the other kernels such as as EpanechnikovKernel.

Non-Allocating Forward-Mode L2 Collocation Loss

The following is an example of a loss function over the collocation that is non-allocating and compatible with forward-mode automatic differentiation:

using PreallocationTools
du = PreallocationTools.dualcache(similar(prob.u0))
preview_est_sol = [@view estimated_solution[:, i] for i in 1:size(estimated_solution, 2)]
preview_est_deriv = [@view estimated_derivative[:, i]
                     for i in 1:size(estimated_solution, 2)]

function construct_iip_cost_function(f, du, preview_est_sol, preview_est_deriv, tpoints)
    function (p)
        _du = PreallocationTools.get_tmp(du, p)
        vecdu = vec(_du)
        cost = zero(first(p))
        for i in 1:length(preview_est_sol)
            est_sol = preview_est_sol[i]
            f(_du, est_sol, p, tpoints[i])
            vecdu .= vec(preview_est_deriv[i]) .- vec(_du)
            cost += sum(abs2, vecdu)
        end
        sqrt(cost)
    end
end
cost_function = construct_iip_cost_function(
    f, du, preview_est_sol, preview_est_deriv, tpoints)