Van der Pol Oscillator

The Van der Pol oscillator is a nonlinear self-sustaining oscillator originally developed to model vacuum tube circuits (van der Pol, 1926). Its limit-cycle dynamics bear structural similarities to cardiac rhythm generation.

Theory

The Van der Pol system is governed by a second-order ODE, written as two first-order equations:

\[\frac{dV}{dt} = W\]

\[\frac{dW}{dt} = \mu(1 - V^2)W - V\]

The nonlinear damping term $\mu(1 - V^2)$ acts as negative damping near $V = 0$ (amplifying small oscillations) and positive damping far from $V = 0$ (suppressing large excursions). This produces a limit cycle: a stable closed orbit that trajectories converge to regardless of initial conditions.

Effect of μ

The parameter $\mu$ controls the character of the oscillation:

μ range	Oscillation type	HRV analogy
0.1 – 1.0	Near-sinusoidal, weakly nonlinear	Regular, low-variability rhythm
1.5 – 2.5	Strong relaxation oscillations	Normal cardiac range
> 3.0	Sharp peaks with slow recovery	Extreme nonlinearity

IBI Generation

The simulation scales the oscillation period to the desired heart_rate (in BPM), with $\mu$ controlling the amplitude and shape of beat-to-beat modulation. The mean IBI is $60000 / \text{heart\_rate}$ ms.

Parameters

Parameter	Range	Default	Description
`μ`	0.1 – 3.0	1.0	Nonlinearity; controls oscillation strength
`heart_rate`	40 – 120 BPM	70	Base heart rate in beats per minute
`σ_noise`	1 – 50	(fitted)	Observation-noise scale (ms); prior `Exponential(10)`, inferred during fitting

Fitting Methods

Method	Algorithm	Loss function
`:gradient`	LBFGS with box constraints (Optim.jl)	Normalized feature-space distance over `mean`, `sdnn`, `rmssd`
`:bayesian`	NUTS MCMC (Turing.jl)	Likelihood: $\text{IBI} \sim \mathcal{N}(\hat{\text{IBI}}, \sigma)$

Examples

Create and Simulate

julia> using HeartRateLab

julia> vdp = VanDerPol()
VanDerPol()

julia> params = (μ = 1.5, heart_rate = 70.0);

julia> ibis = simulate(vdp, params, 10);

julia> length(ibis)
10

Fit with Gradient Method

using HeartRateLab

vdp = VanDerPol()
ibis = read_txt("data.txt")

result = fit(vdp, ibis; method=:gradient)

println("Fitted μ:          ", round(result.params.μ; digits=3))
println("Fitted heart rate: ", round(result.params.heart_rate; digits=1), " bpm")
println("Converged:         ", result.diagnostics["converged"])

# Generate synthetic
synthetic = simulate(result.model, result.params, length(ibis))

Fit with Bayesian Inference

using HeartRateLab

vdp = VanDerPol()
ibis = read_txt("data.txt")

result = fit(vdp, ibis; method=:bayesian, chains=4, samples=1000)

println("Posterior mean μ:          ", round(result.params.μ; digits=3))
println("Posterior mean heart rate: ", round(result.params.heart_rate; digits=1), " bpm")
println("R-hat μ:                   ", round(result.diagnostics["rhat_mu"]; digits=3))

# Posterior samples
μ_samples = parameter_series(result, :μ)

API Reference

HeartRateLab.Models.VanDerPol — Type

VanDerPol <: AbstractHRVModel

Van der Pol oscillator model for HRV simulation. Simple mechanistic model for cardiac oscillations.

source

HeartRateLab.Models.parameter_space — Method

parameter_space(model::VanDerPol) -> NamedTuple

Return the parameter space for Van der Pol model with priors for Bayesian fitting.

Returns

NamedTuple with keys: μ, heart_rate, σ_noise Each parameter has: lower, upper, prior (TruncatedNormal or Exponential)

source

HeartRateLab.Models.simulate — Method

simulate(model::VanDerPol, params::NamedTuple, n_beats::Int) -> Vector{Float64}

Simulate IBI time series using Van der Pol oscillator.

Parameters

μ: Non-linearity parameter (amplitude of oscillation)
heart_rate: Base heart rate in BPM

Returns

Vector of inter-beat intervals in milliseconds

source

HeartRateLab.Models.fit — Method

fit(model::VanDerPol, data::Vector{Float64}; method=:bayesian, kwargs...) -> ModelFitResult

Fit Van der Pol model to IBI data using Bayesian or gradient-based optimization.

Arguments

model::VanDerPol: Van der Pol model instance
data::Vector{Float64}: IBI time series in milliseconds
method::Symbol: :bayesian (NUTS MCMC) or :gradient (LBFGS)
chains::Int: Number of MCMC chains (default: 4, for Bayesian only)
samples::Int: Samples per chain (default: 1000, for Bayesian only)

Returns

ModelFitResult with posterior samples, MAP estimates, and diagnostics

source