Symbolic Regression of a Nonlinear System via Lifting

To infer more complex examples, EQSearch also can be called with a Basis to use predefined features. Let's look at the well-known pendulum model

using DataDrivenDiffEq
using LinearAlgebra
using OrdinaryDiffEq
using SciMLBase: NoSpecialize
using DataDrivenSR
using Plots

function pendulum!(du, u, p, t)
    du[1] = u[2]
    return du[2] = -9.81 * sin(u[1])
end

u0 = [0.1, π / 2]
tspan = (0.0, 10.0)
sys = ODEProblem{true, NoSpecialize}(pendulum!, u0, tspan)
sol = solve(sys, Tsit5());

We will use the data provided by our problem.

prob = DataDrivenProblem(sol)

Continuous DataDrivenProblem{Float64} ##DDProblem#469 in 2 dimensions and 43 samples

And plot the problems data.

plot(prob)

To solve our problem, we will use EQSearch, which provides a wrapper for the symbolic regression interface. We will stick to simple operations, use a L1DistLoss, and limit the verbosity of the algorithm. Note that we do not include sin, but rather lift the search space of variables.

@variables u[1:2]
u = collect(u)

basis = Basis([polynomial_basis(u, 2); sin.(u)], u)

eqsearch_options = SymbolicRegression.Options(
    binary_operators = [+, *],
    loss = L1DistLoss(),
    verbosity = -1, progress = false, npop = 30,
    timeout_in_seconds = 60.0
)

alg = EQSearch(eq_options = eqsearch_options)

EQSearch
  weights: Nothing nothing
  numprocs: Nothing nothing
  procs: Nothing nothing
  addprocs_function: Nothing nothing
  parallelism: Symbol serial
  runtests: Bool true
  eq_options: Options{SymbolicRegression.CoreModule.OptionsStructModule.ComplexityMapping{Int64, Int64}, OperatorEnum, Node, Expression, @NamedTuple{}, MutationWeights, false, false, nothing, Nothing, 5}

Again, we solve the problem to obtain a DataDrivenSolution with similar options as the previous example but provide a Basis along the arguments.

res = solve(prob, basis, alg, options = DataDrivenCommonOptions(maxiters = 100))
println(res)

"DataDrivenSolution{Float64}" with 2 equations and 0 parameters.
Returncode: Success
Residual sum of squares: 0.0

Currently, the parameters of the result found by EQSearch are not turned into symbolic parameters. This affects some functions like dof, aicc, bic.

system = get_basis(res)

Model ##Basis#475 with 2 equations
States : u[1] u[2]
Independent variable: t
Equations
Differential(t, 1)(u[1]) = u[2]
Differential(t, 1)(u[2]) = -9.81sin(u[1])

Copy-Pasteable Code

using DataDrivenDiffEq
using LinearAlgebra
using OrdinaryDiffEq
using SciMLBase: NoSpecialize
using DataDrivenSR

function pendulum!(du, u, p, t)
    du[1] = u[2]
    return du[2] = -9.81 * sin(u[1])
end

u0 = [0.1, π / 2]
tspan = (0.0, 10.0)
sys = ODEProblem{true, NoSpecialize}(pendulum!, u0, tspan)
sol = solve(sys, Tsit5());

prob = DataDrivenProblem(sol)

@variables u[1:2]
u = collect(u)

basis = Basis([polynomial_basis(u, 2); sin.(u)], u)

eqsearch_options = SymbolicRegression.Options(
    binary_operators = [+, *],
    loss = L1DistLoss(),
    verbosity = -1, progress = false, npop = 30,
    timeout_in_seconds = 60.0
)

alg = EQSearch(eq_options = eqsearch_options)

res = solve(prob, basis, alg, options = DataDrivenCommonOptions(maxiters = 100))

system = get_basis(res)

# This file was generated using Literate.jl, https://github.com/fredrikekre/Literate.jl

This page was generated using Literate.jl.