Fitting

The fit function can be used to fit a model to a recorded series with either the Whittle or debiased Whittle likelihoods using the interior point Newton method as implemented in Optim.jl.

First load the package:

using WhittleLikelihoodInference

As a basic example, consider a univariate Gaussian process with JONSWAP spectral density function. We can simulate such a process with the following:

n = 2048
Δ = 1.0
nreps = 1
σ = 0.8
θ = 0.6
ts = simulate_gp(OU(σ,θ),n,Δ,nreps)[1]

WhittleLikelihoodInference.TimeSeries{Float64, 1}([0.8957923382764735, -0.2869095098623978, -0.5843731274042048, -1.1002520862595118, 0.14245954716497766, 1.9878330267833317, 1.5509562964594539, 0.616577021066969, 0.7720808384681657, -0.02261673450707613  …  0.7219220841823171, 0.6494898047231441, 0.7373143779725513, 0.883422919316776, -0.2419594447544285, 0.07526589891407498, 0.009206327103899757, 0.6409224395179699, 0.0933146308439764, 0.38942079856662287], 1.0)

In this case, we have simulated a Gaussian process with OU spectral density function with parameters σ = 0.8, θ = 0.6 of length 2048 sampled every 1.0 seconds. The function simulate_gp will simulate a vector of nreps series, which is why we recover the first of these to get one series.

We can now fit a model with the following:

x₀ = [σ,θ] .+ 0.1
res = fit(ts,model=OU,x₀=x₀)

 * Status: success

 * Candidate solution
    Final objective value:     -2.293723e+03

 * Found with
    Algorithm:     Interior Point Newton

 * Convergence measures
    |x - x'|               = 0.00e+00 ≤ 0.0e+00
    |x - x'|/|x'|          = 0.00e+00 ≤ 0.0e+00
    |f(x) - f(x')|         = 0.00e+00 ≤ 0.0e+00
    |f(x) - f(x')|/|f(x')| = 0.00e+00 ≤ 0.0e+00
    |g(x)|                 = 3.19e-05 ≰ 1.0e-08

 * Work counters
    Seconds run:   1  (vs limit Inf)
    Iterations:    25
    f(x) calls:    119
    ∇f(x) calls:   119

Here x₀ is the vector of initial parameter guesses. The estimated parameter vector can be recovered by doing:

x̂ = res.minimizer

2-element Vector{Float64}:
 0.8077832308235766
 0.6100036495790229

The full example is:

import Random
Random.seed!(1234)
using WhittleLikelihoodInference
n = 2048
Δ = 1.0
nreps = 1
σ = 0.8
θ = 0.6
ts = simulate_gp(OU(σ,θ),n,Δ,nreps)[1]
x₀ = [σ,θ] .+ 0.1
res = fit(ts,model=OU,x₀=x₀)
x̂ = res.minimizer

2-element Vector{Float64}:
 0.8077832308235766
 0.6100036495790229

For more options, see below:

fit(ts,Δ;model::Type{<:TimeSeriesModel},x₀,lowerΩcutoff,upperΩcutoff,x_lowerbounds,x_upperbounds,method,taper)
fit(timeseries::TimeSeries;model::Type{<:TimeSeriesModel},x₀,lowerΩcutoff,upperΩcutoff,x_lowerbounds,x_upperbounds,method,taper)