Duplicate emulators when repeating `emulator = Emulator(); optimize_hyperparameters!(emulator)`

[nefrathenrici]

Repeating emulator = Emulator(gauss_proc, test_data); optimize_hyperparameters!(emulator) seems to create multiple copies of the emulator within the emulator object, rather than overwriting the emulator each time.

I'm not sure if this is a bug or if it is intended and I am using the API improperly.

I also see LinearAlgebra.PosDefException(3) popping up in the output of optimize_hyperparameters, I am not sure if that is related.

Simple reproducer:

julia> import CalibrateEmulateSample as CES

julia> using CalibrateEmulateSample.Emulators

julia> import EnsembleKalmanProcesses as EKP

julia> input_data = hcat(collect(range(1.0,50.0))...)
1×50 Matrix{Float64}:
 1.0  2.0  3.0  4.0  5.0  …  47.0  48.0  49.0  50.0

julia> output_data = hcat(collect(range(51.0,100.0))...)

1×50 Matrix{Float64}:
 51.0  52.0  53.0  54.0  …  97.0  98.0  99.0  100.0

julia> test_data = EKP.DataContainers.PairedDataContainer(input_data, output_data)
EnsembleKalmanProcesses.DataContainers.PairedDataContainer{Float64}(EnsembleKalmanProcesses.DataContainers.DataContainer{Float64}([1.0 2.0 … 49.0 50.0]), EnsembleKalmanProcesses.DataContainers.DataContainer{Float64}([51.0 52.0 … 99.0 100.0]))

julia> gppackage = Emulators.GPJL()
GPJL()

julia> gauss_proc = Emulators.GaussianProcess(gppackage, noise_learn = false)
GaussianProcess{GPJL, Float64}(Union{Nothing, PyCall.PyObject, GaussianProcesses.GPE}[], nothing, false, 1.0, YType())

julia> emulator = Emulator(gauss_proc,test_data); optimize_hyperparameters!(emulator)
Using default squared exponential kernel, learning length scale and variance parameters
Using default squared exponential kernel: Type: GaussianProcesses.SEArd{Float64}, Params: [-0.0, 0.0]
kernel in GaussianProcess:
Type: GaussianProcesses.SEArd{Float64}, Params: [-0.0, 0.0]
created GP: 1
LinearAlgebra.PosDefException(3)
optimized hyperparameters of GP: 1
Type: GaussianProcesses.SEArd{Float64}, Params: [3.0834863232540557, 5.290309213025948]

julia> emulator = Emulator(gauss_proc,test_data); optimize_hyperparameters!(emulator)
Using default squared exponential kernel, learning length scale and variance parameters
Using default squared exponential kernel: Type: GaussianProcesses.SEArd{Float64}, Params: [-0.0, 0.0]
kernel in GaussianProcess:
Type: GaussianProcesses.SEArd{Float64}, Params: [-0.0, 0.0]
created GP: 1
optimized hyperparameters of GP: 1
Type: GaussianProcesses.SEArd{Float64}, Params: [3.0834863232540557, 5.290309213025948]
LinearAlgebra.PosDefException(3)
optimized hyperparameters of GP: 2
Type: GaussianProcesses.SEArd{Float64}, Params: [3.0834863232540557, 5.290309213025948]

julia> emulator = Emulator(gauss_proc,test_data); optimize_hyperparameters!(emulator)
Using default squared exponential kernel, learning length scale and variance parameters
Using default squared exponential kernel: Type: GaussianProcesses.SEArd{Float64}, Params: [-0.0, 0.0]
kernel in GaussianProcess:
Type: GaussianProcesses.SEArd{Float64}, Params: [-0.0, 0.0]
created GP: 1
optimized hyperparameters of GP: 1
Type: GaussianProcesses.SEArd{Float64}, Params: [3.0834863232540557, 5.290309213025948]
optimized hyperparameters of GP: 2
Type: GaussianProcesses.SEArd{Float64}, Params: [3.0834863232540557, 5.290309213025948]
LinearAlgebra.PosDefException(3)
optimized hyperparameters of GP: 3
Type: GaussianProcesses.SEArd{Float64}, Params: [3.0834863232540557, 5.290309213025948]

Note that recreating the GaussianProcess resets the emulator var.

gauss_proc = Emulators.GaussianProcess(gppackage, noise_learn = false)
emulator = Emulator(gauss_proc, input_output_pairs)#, normalize_inputs = true,  obs_noise_cov = Γ)
optimize_hyperparameters!(emulator)

Repeated calls to just optimize_hyperparameters! do not have this effect either:

julia> optimize_hyperparameters!(emulator)
optimized hyperparameters of GP: 1
Type: GaussianProcesses.SEArd{Float64}, Params: [4.185063788518014, 5.533119621637727]

julia> optimize_hyperparameters!(emulator)
optimized hyperparameters of GP: 1
Type: GaussianProcesses.SEArd{Float64}, Params: [4.185063788518014, 5.533119621637727]

[odunbar]

Thanks for raising this:

Response

1. the posdef errors frequently occur during Optim's attempt to train GP so I think that can be ignored. ML training is a messy problem, and this happens inside the packages we are depending on.

2. As you note, the duplication happens because you are passing in the same machine_learning_tool i.e. gauss_proc into emulator. The Emulator is effective just a wrapper while gauss_proc is the object which is actually trained. In the end, there is a loop with a push!(models, m) which new trained GP(s) to be stored for later prediction.

   Solution

   I think in terms of resolution, there are two options,

3. We either don't add to a non-empty model vector

4. We throw a warning and replace the old GP in the emulator

I'll have a think about what is safest!

Update I went option 1. plus added a warning. This will alert users that they are about to overwrite information and encourage them to use a new GP instead