JuliaStats / TimeModels.jl

Modeling time series in Julia
Other
57 stars 28 forks source link

Switch Kalman smoothing algorithm to Durbin-Koopman #50

Closed GordStephen closed 8 years ago

GordStephen commented 8 years ago

As mentioned in #49, this provides performance advantages, particularily when working with large state vectors. Some unscientific benchmarks from my initial tests:

master (RTS smoothing):

julia> @time kalman_smooth(y_noisy, lm)
  1.497366 seconds (17.18 M allocations: 875.542 MB, 10.39% gc time)
KalmanSmoothed{Float64}
100001 observations, 2-D process x 1-D observations
Negative log-likelihood: 225675.43959948

DK Smoothing

julia> @time smooth(y_noisy, lm)
  1.268853 seconds (17.18 M allocations: 752.029 MB, 8.73% gc time)
KalmanSmoothed{Float64}
100001 observations, 2-D process x 1-D observations
Negative log-likelihood: 126623.43103532399

For this small system the DK smoother is consistently a little faster, and also returns slightly better estimates (maybe from avoided numerical imprecisions from the matrix inversions?) which results in the lower negative-log-likelihood.

The performance gains become more apparent with larger state vectors:

master:

julia> @time kalman_smooth(y_noisy, blm)
 19.162492 seconds (17.18 M allocations: 20.914 GB, 6.62% gc time)
KalmanSmoothed{Float64}
100001 observations, 30-D process x 1-D observations
Negative log-likelihood: 228215.11686602756

DK Smoothing

julia> @time smooth(y_noisy, blm)
  9.845013 seconds (17.18 M allocations: 16.056 GB, 8.84% gc time)
KalmanSmoothed{Float64}
100001 observations, 30-D process x 1-D observations
Negative log-likelihood: 130779.03074016975
rob-luke commented 8 years ago

Nice improvements! You can't argue with faster and more accurate.

As in #49 I would suggest keeping stand alone filtering too. And maybe a selection of smoothing methods, but the default being what is fastest and most accurate.

GordStephen commented 8 years ago

Done in #52