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.38458076; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0743635 … 0.0 0.0]), states_modifiers = (), readout = (weight = [155.28373476559958 -135.78488253111578 … 155.81065440534564 -17.073510821132764; 269.8674304337835 -308.7893343814745 … 308.2896055814805 -132.12494440168194; … ; 278.42285642481767 -303.7029763829606 … 244.19644626190083 -197.56890846440453; -3.5290205093834364 37.19219131890375 … -1.255283593537792 31.084410740022005],)), (reservoir = (cell = (rng = Random.MersenneTwister(42, (0, 25050, 24048, 870)),), carry = (Float32[0.14156413; -0.036176555; … ; -0.38727778; 0.09736929;;],)), states_modifiers = (), readout = NamedTuple(), states = Float32[0.6755669 0.16517994 … 0.16209146 0.14156413; 0.00037146357 -0.30405423 … 0.0013690992 -0.036176555; … ; -0.2935045 0.2056673 … -0.34522423 -0.38727778; -0.7566877 -0.59265745 … 0.04711806 0.09736929]))

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.38458076; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0743635 … 0.0 0.0]), states_modifiers = (), readout = (weight = Float32[-0.4062664 0.85807514 … -8.48817 4.1298366; -1.7177653 10.826405 … -3.3635206 -7.9875226; … ; 0.47673988 14.949505 … 4.5631976 -1.9362493; 5.3215685 6.846592 … -10.189642 -2.0689604],)), (reservoir = (cell = (rng = Random.MersenneTwister(42, (0, 25050, 24048, 870)),), carry = (Float32[0.14502546; -0.030527089; … ; -0.3896055; 0.09629535;;],)), states_modifiers = (), readout = NamedTuple(), states = Float32[0.24012244 0.09647706 … 0.15855107 0.14502546; 0.070600435 -0.07172653 … -0.0055307266 -0.030527089; … ; -0.3597555 -0.3799423 … -0.3430129 -0.3896055; 0.038221136 0.031300407 … 0.048109192 0.09629535]))

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.38458076; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0743635 … 0.0 0.0]), states_modifiers = (), readout = (weight = Float32[24.169455 77.279396 … 113.73274 305.27502; 184.7662 -159.7874 … -480.27 -178.68349; … ; -400.65005 469.0911 … 1360.4556 1019.7482; -90.6014 85.59617 … 211.81766 401.6941],)), (reservoir = (cell = (rng = Random.MersenneTwister(42, (0, 25050, 24048, 870)),), carry = (Float32[0.14514823; -0.030314906; … ; -0.38968092; 0.0962635;;],)), states_modifiers = (), readout = NamedTuple(), states = Float32[0.2366644 0.09987559 … 0.1584202 0.14514823; 0.064417884 -0.06638202 … -0.005780449 -0.030314906; … ; -0.35722053 -0.3824757 … -0.3429394 -0.38968092; 0.03926979 0.030110914 … 0.048138287 0.0962635]))

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.38458076; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0743635 … 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.14515251; -0.03030742; … ; -0.38968354; 0.09626235;;],)), states_modifiers = (), readout = NamedTuple(), states = Float32[0.23654222 0.09999258 … 0.15841565 0.14515251; 0.06419541 -0.0661864 … -0.0057892012 -0.03030742; … ; -0.35713533 -0.38255987 … -0.34293678 -0.38968354; 0.03930358 0.03007245 … 0.048139296 0.09626235]))

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.05689015 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 -0.0641678; 0.0 0.0 … 0.00967149 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.05689015 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 -0.0641678; 0.0 0.0 … 0.00967149 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.07993371; -0.056284487; … ; -0.29779446; -0.016760059;;],)), states_modifiers = (), readout = NamedTuple(), states = Float32[0.77943015 -0.6124105 … -0.038226303 -0.07993371; -0.8753577 0.0026252898 … 0.0798458 -0.056284487; … ; -0.8539995 -0.89032215 … -0.25453463 -0.29779446; 0.27237403 -0.120416425 … -0.07133219 -0.016760059]))