Training Reservoir Computing Models

Training reservoir computing (RC) models usually means solving a linear regression problem. ReservoirComputing.jl offers multiple stratedies to provide a readout; in this page we will show the basics, while also pointing out the possible extensions.

Training in ReservoirComputing.jl: Ridge Regression

The most simple training of RC models is through ridge regression. Given the widepread adoption of this training mechanism, ridge regression is the default training algorithm for RC models in the library.

using ReservoirComputing
using Random
Random.seed!(42)
rng = MersenneTwister(42)

input_data = rand(Float32, 3, 100)
target_data = rand(Float32, 5, 100)

model = ESN(3, 100, 5)
ps, st = setup(rng, model)
ps, st = train!(model, input_data, target_data, ps, st,
    StandardRidge(); # default
    solver = QRSolver()) # default

((reservoir = (input_matrix = Float32[0.086408615 0.0998888 -0.013271046; -0.06170206 -0.06359978 0.0921186; … ; -0.027093245 -0.04336114 -0.04771917; 0.07098396 0.015740562 0.019439578], reservoir_matrix = Float32[0.0 0.0 … 0.0 0.38458103; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.07436356 … 0.0 0.0]), states_modifiers = (), readout = (weight = [155.35108756323686 -135.8712657819953 … 155.88074896713223 -17.109907041437385; 269.94220128782416 -308.8849300590367 … 308.3685883855184 -132.17113965976264; … ; 278.4173811801401 -303.6919995502757 … 244.19595651464243 -197.58585602683507; -3.5393982571778024 37.21017225452454 … -1.2648540501081846 31.08003819647461],)), (reservoir = (cell = (rng = Random.MersenneTwister(42, (0, 25050, 24048, 870)),), carry = (Float32[0.14156468; -0.036177248; … ; -0.38727927; 0.09736941;;],)), states_modifiers = (), readout = NamedTuple(), states = Float32[0.67556715 0.16518015 … 0.1620921 0.14156468; 0.0003713667 -0.30405453 … 0.0013688962 -0.036177248; … ; -0.2935047 0.20566761 … -0.345226 -0.38727927; -0.7566881 -0.592658 … 0.047118302 0.09736941]))

In this call you can see that there are two possible knobs to be modified: the loss function, in this case ridge, and the solver, in this case the build in QR factorization. In the remaining part of this tutorial we will see how it is possible to change either.

Changing Ridge Regression Solver

Building on SciML's LinearSolve.jl, it is possible to leverage multiple solvers for the ridge problem. For instance, building on the previous example:

using LinearSolve

ps, st = train!(model, input_data, target_data, ps, st,
    StandardRidge(); # default
    solver = SVDFactorization()) # from LinearSolve

((reservoir = (input_matrix = Float32[0.086408615 0.0998888 -0.013271046; -0.06170206 -0.06359978 0.0921186; … ; -0.027093245 -0.04336114 -0.04771917; 0.07098396 0.015740562 0.019439578], reservoir_matrix = Float32[0.0 0.0 … 0.0 0.38458103; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.07436356 … 0.0 0.0]), states_modifiers = (), readout = (weight = Float32[-0.42554864 2.830525 … -10.93533 9.562259; -1.7472765 11.085349 … -2.2991953 -6.7392015; … ; 0.51117396 16.642477 … 3.350666 1.195735; 4.609715 -1.8165002 … -16.057856 6.02251],)), (reservoir = (cell = (rng = Random.MersenneTwister(42, (0, 25050, 24048, 870)),), carry = (Float32[0.14502607; -0.030527564; … ; -0.38960707; 0.09629541;;],)), states_modifiers = (), readout = NamedTuple(), states = Float32[0.24012297 0.09647721 … 0.15855159 0.14502607; 0.070599996 -0.07172731 … -0.0055311597 -0.030527564; … ; -0.3597568 -0.3799437 … -0.34301454 -0.38960707; 0.038221236 0.031300593 … 0.048109535 0.09629541]))

or

ps, st = train!(model, input_data, target_data, ps, st,
    StandardRidge(); # default
    solver = QRFactorization()) # from LinearSolve

((reservoir = (input_matrix = Float32[0.086408615 0.0998888 -0.013271046; -0.06170206 -0.06359978 0.0921186; … ; -0.027093245 -0.04336114 -0.04771917; 0.07098396 0.015740562 0.019439578], reservoir_matrix = Float32[0.0 0.0 … 0.0 0.38458103; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.07436356 … 0.0 0.0]), states_modifiers = (), readout = (weight = Float32[122.6651 11.0168705 … -262.47482 122.44317; 170.44348 -230.26588 … -110.79341 272.9332; … ; 46.673862 70.66685 … -362.86142 471.12817; -93.76712 82.558914 … 130.86298 -78.306595],)), (reservoir = (cell = (rng = Random.MersenneTwister(42, (0, 25050, 24048, 870)),), carry = (Float32[0.14514889; -0.030315312; … ; -0.38968256; 0.09626353;;],)), states_modifiers = (), readout = NamedTuple(), states = Float32[0.2366648 0.099875934 … 0.15842067 0.14514889; 0.06441728 -0.0663825 … -0.0057809157 -0.030315312; … ; -0.3572218 -0.38247722 … -0.34294105 -0.38968256; 0.039269898 0.030111087 … 0.048138604 0.09626353]))

For a detailed explanation of the different solvers, as well as a complete list of them, we suggest visiting the appropriate page in LinearSolve's documentation

Changing Linear Regression Problem

Linear regression is a general problem, which can be espressed through multiple different loss functions. While ridge regression is the most common in RC, due to its closed form, there are multiple other available. ReservoirComputing.jl leverages MLJLinearModels.jl to access all the methods available from that library.

The train function can be called as before, only this time you can specify different models and different solvers for the linear regression problem:

using MLJLinearModels

ps, st = train!(model, input_data, target_data, ps, st,
    LassoRegression(fit_intercept=false); # from MLJLinearModels
    solver = ProxGrad()) # from MLJLinearModels

((reservoir = (input_matrix = Float32[0.086408615 0.0998888 -0.013271046; -0.06170206 -0.06359978 0.0921186; … ; -0.027093245 -0.04336114 -0.04771917; 0.07098396 0.015740562 0.019439578], reservoir_matrix = Float32[0.0 0.0 … 0.0 0.38458103; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.07436356 … 0.0 0.0]), states_modifiers = (), readout = (weight = Float32[0.0 -0.0 … -0.0 0.0; 0.0 -0.0 … -0.0 0.0; … ; 0.0 -0.0 … -0.0 0.0; 0.0 -0.0 … -0.0 0.0],)), (reservoir = (cell = (rng = Random.MersenneTwister(42, (0, 25050, 24048, 870)),), carry = (Float32[0.1451532; -0.030307855; … ; -0.3896852; 0.09626243;;],)), states_modifiers = (), readout = NamedTuple(), states = Float32[0.23654257 0.09999292 … 0.15841615 0.1451532; 0.06419469 -0.066186875 … -0.0057896706 -0.030307855; … ; -0.35713658 -0.3825614 … -0.34293845 -0.3896852; 0.039303705 0.030072607 … 0.048139594 0.09626243]))

Make sure to check the MLJLinearModels documentation pages for the available models and solvers. Please note that not all solvers can be used on all the models.

Support Vector Regression

ReservoirComputing.jl also allows users to train RC models with support vector regression through LIBSVM.jl. However, the majority of builtin models in the library uses a LinearReadout by default, which can only be trained with linear regression. In order to use support vector regression, one needs to build a model with SVMReadout

using LIBSVM

model = ReservoirComputer(
    StatefulLayer(ESNCell(3=>100)),
    SVMReadout(100=>5)
)

ps, st = setup(rng, model)

((reservoir = (input_matrix = Float32[0.036005475 -0.07744928 0.00549283; -0.027602648 -0.03803761 0.029640054; … ; 0.020895792 -0.04784565 -0.09396603; -0.05864432 0.05747981 0.08333655], reservoir_matrix = Float32[0.0 0.0 … -0.056890085 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 -0.06416772; 0.0 0.0 … 0.009671479 0.0]), states_modifiers = (), readout = NamedTuple()), (reservoir = (cell = (rng = Random.MersenneTwister(42, (0, 50100, 49098, 357)),), carry = nothing), states_modifiers = (), readout = NamedTuple()))

We can now train our new model similarly to before:

ps, st = train!(model, input_data, target_data, ps, st,
    EpsilonSVR() # from LIBSVM
    )

((reservoir = (input_matrix = Float32[0.036005475 -0.07744928 0.00549283; -0.027602648 -0.03803761 0.029640054; … ; 0.020895792 -0.04784565 -0.09396603; -0.05864432 0.05747981 0.08333655], reservoir_matrix = Float32[0.0 0.0 … -0.056890085 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 -0.06416772; 0.0 0.0 … 0.009671479 0.0]), states_modifiers = (), readout = (models = Any[LIBSVM.EpsilonSVR(LIBSVM.Kernel.RadialBasis, 0.01, 0.1, 1.0, 3, 0.0, 0.001, true, false, LIBSVM.SVM{Float64, LIBSVM.Kernel.KERNEL}(LIBSVM.EpsilonSVR, LIBSVM.Kernel.RadialBasis, nothing, 100, 100, 2, Float64[], Int32[], Float64[], Int32[], LIBSVM.SupportVectors{Vector{Float32}, Matrix{Float32}}(0, Int32[], Float32[], Matrix{Float32}(undef, 100, 0), Int32[], LIBSVM.SVMNode[]), 0.0, Matrix{Float64}(undef, 0, 1), Float64[], Float64[], [-0.0036369473622437157], 3, 0.01, 200.0, 0.001, 1.0, 0.5, 0.1, true, false)), LIBSVM.EpsilonSVR(LIBSVM.Kernel.RadialBasis, 0.01, 0.1, 1.0, 3, 0.0, 0.001, true, false, LIBSVM.SVM{Float64, LIBSVM.Kernel.KERNEL}(LIBSVM.EpsilonSVR, LIBSVM.Kernel.RadialBasis, nothing, 100, 100, 2, Float64[], Int32[], Float64[], Int32[], LIBSVM.SupportVectors{Vector{Float32}, Matrix{Float32}}(0, Int32[], Float32[], Matrix{Float32}(undef, 100, 0), Int32[], LIBSVM.SVMNode[]), 0.0, Matrix{Float64}(undef, 0, 1), Float64[], Float64[], [-0.0036369473622437157], 3, 0.01, 200.0, 0.001, 1.0, 0.5, 0.1, true, false)), LIBSVM.EpsilonSVR(LIBSVM.Kernel.RadialBasis, 0.01, 0.1, 1.0, 3, 0.0, 0.001, true, false, LIBSVM.SVM{Float64, LIBSVM.Kernel.KERNEL}(LIBSVM.EpsilonSVR, LIBSVM.Kernel.RadialBasis, nothing, 100, 100, 2, Float64[], Int32[], Float64[], Int32[], LIBSVM.SupportVectors{Vector{Float32}, Matrix{Float32}}(0, Int32[], Float32[], Matrix{Float32}(undef, 100, 0), Int32[], LIBSVM.SVMNode[]), 0.0, Matrix{Float64}(undef, 0, 1), Float64[], Float64[], [-0.0036369473622437157], 3, 0.01, 200.0, 0.001, 1.0, 0.5, 0.1, true, false)), LIBSVM.EpsilonSVR(LIBSVM.Kernel.RadialBasis, 0.01, 0.1, 1.0, 3, 0.0, 0.001, true, false, LIBSVM.SVM{Float64, LIBSVM.Kernel.KERNEL}(LIBSVM.EpsilonSVR, LIBSVM.Kernel.RadialBasis, nothing, 100, 100, 2, Float64[], Int32[], Float64[], Int32[], LIBSVM.SupportVectors{Vector{Float32}, Matrix{Float32}}(0, Int32[], Float32[], Matrix{Float32}(undef, 100, 0), Int32[], LIBSVM.SVMNode[]), 0.0, Matrix{Float64}(undef, 0, 1), Float64[], Float64[], [-0.0036369473622437157], 3, 0.01, 200.0, 0.001, 1.0, 0.5, 0.1, true, false)), LIBSVM.EpsilonSVR(LIBSVM.Kernel.RadialBasis, 0.01, 0.1, 1.0, 3, 0.0, 0.001, true, false, LIBSVM.SVM{Float64, LIBSVM.Kernel.KERNEL}(LIBSVM.EpsilonSVR, LIBSVM.Kernel.RadialBasis, nothing, 100, 100, 2, Float64[], Int32[], Float64[], Int32[], LIBSVM.SupportVectors{Vector{Float32}, Matrix{Float32}}(0, Int32[], Float32[], Matrix{Float32}(undef, 100, 0), Int32[], LIBSVM.SVMNode[]), 0.0, Matrix{Float64}(undef, 0, 1), Float64[], Float64[], [-0.0036369473622437157], 3, 0.01, 200.0, 0.001, 1.0, 0.5, 0.1, true, false))],)), (reservoir = (cell = (rng = Random.MersenneTwister(42, (0, 50100, 49098, 457)),), carry = (Float32[-0.07993154; -0.056282256; … ; -0.2977938; -0.016758788;;],)), states_modifiers = (), readout = NamedTuple(), states = Float32[0.77942973 -0.61240965 … -0.038224187 -0.07993154; -0.87535733 0.0026252489 … 0.07984722 -0.056282256; … ; -0.85399914 -0.89032173 … -0.25453427 -0.2977938; 0.27237374 -0.12041608 … -0.07133142 -0.016758788]))