Lotka-Volterra Parameter Estimation Benchmarks

Vaibhav Dixit, Chris Rackauckas

Parameter estimation of Lotka Volterra model using optimisation methods

using ParameterizedFunctions, OrdinaryDiffEq, DiffEqParamEstim
using BlackBoxOptim, NLopt, Plots, RecursiveArrayTools, QuadDIRECT
gr(fmt=:png)
Plots.GRBackend()
loc_bounds = Tuple{Float64, Float64}[(0, 5), (0, 5), (0, 5), (0, 5)]
glo_bounds = Tuple{Float64, Float64}[(0, 10), (0, 10), (0, 10), (0, 10)]
loc_init = [1,0.5,3.5,1.5]
glo_init = [5,5,5,5]
f = @ode_def LotkaVolterraTest begin
    dx = a*x - b*x*y
    dy = -c*y + d*x*y
end a b c d
u0 = [1.0,1.0]                          #initial values
tspan = (0.0,10.0)
p = [1.5,1.0,3.0,1,0]                   #parameters used, these need to be estimated from the data
tspan = (0.0, 30.0)                     # sample of 3000 observations over the (0,30) timespan
prob = ODEProblem(f, u0, tspan,p)
tspan2 = (0.0, 3.0)                     # sample of 3000 observations over the (0,30) timespan
prob_short = ODEProblem(f, u0, tspan2,p)
dt = 30.0/3000
tf = 30.0
tinterval = 0:dt:tf
t  = collect(tinterval)
h = 0.01
M = 300
tstart = 0.0
tstop = tstart + M * h
tinterval_short = 0:h:tstop
t_short = collect(tinterval_short)
#Generate Data
data_sol_short = solve(prob_short,Tsit5(),saveat=t_short,reltol=1e-9,abstol=1e-9)
data_short = convert(Array, data_sol_short)
data_sol = solve(prob,Tsit5(),saveat=t,reltol=1e-9,abstol=1e-9)
data = convert(Array, data_sol)

Plot of the solution

Short Solution
p1 = plot(data_sol_short)
Longer Solution
p2 = plot(data_sol)

Local Solution from the short data set

obj_short = build_loss_objective(prob_short,Tsit5(),L2Loss(t_short,data_short),tstops=t_short)
res1 = bboptimize(obj_short;SearchRange = glo_bounds, MaxSteps = 7e3)
Starting optimization with optimizer BlackBoxOptim.DiffEvoOpt{BlackBoxOptim
.FitPopulation{Float64},BlackBoxOptim.RadiusLimitedSelector,BlackBoxOptim.A
daptiveDiffEvoRandBin{3},BlackBoxOptim.RandomBound{BlackBoxOptim.RangePerDi
mSearchSpace}}
0.00 secs, 0 evals, 0 steps
0.50 secs, 1807 evals, 1696 steps, improv/step: 0.203 (last = 0.2034), fitn
ess=245.789636091
1.00 secs, 3281 evals, 3170 steps, improv/step: 0.189 (last = 0.1716), fitn
ess=13.445992864
1.50 secs, 4849 evals, 4738 steps, improv/step: 0.184 (last = 0.1754), fitn
ess=0.227123561
2.00 secs, 6305 evals, 6194 steps, improv/step: 0.186 (last = 0.1902), fitn
ess=0.002694708

Optimization stopped after 7001 steps and 2.35699200630188 seconds
Termination reason: Max number of steps (7000) reached
Steps per second = 2970.311304103474
Function evals per second = 3017.4052270795464
Improvements/step = 0.187
Total function evaluations = 7112


Best candidate found: [1.50062, 1.00013, 2.99721, 0.999021]

Fitness: 0.000346021
# Lower tolerance could lead to smaller fitness (more accuracy)
obj_short = build_loss_objective(prob_short,Tsit5(),L2Loss(t_short,data_short),tstops=t_short,reltol=1e-9)
res1 = bboptimize(obj_short;SearchRange = glo_bounds, MaxSteps = 7e3)
Starting optimization with optimizer BlackBoxOptim.DiffEvoOpt{BlackBoxOptim
.FitPopulation{Float64},BlackBoxOptim.RadiusLimitedSelector,BlackBoxOptim.A
daptiveDiffEvoRandBin{3},BlackBoxOptim.RandomBound{BlackBoxOptim.RangePerDi
mSearchSpace}}
0.00 secs, 0 evals, 0 steps
0.50 secs, 1584 evals, 1460 steps, improv/step: 0.214 (last = 0.2144), fitn
ess=423.685215058
1.00 secs, 3296 evals, 3172 steps, improv/step: 0.178 (last = 0.1478), fitn
ess=25.337585437
1.50 secs, 5005 evals, 4882 steps, improv/step: 0.178 (last = 0.1760), fitn
ess=0.248299253
2.00 secs, 6539 evals, 6416 steps, improv/step: 0.176 (last = 0.1721), fitn
ess=0.005103412

Optimization stopped after 7001 steps and 2.2087929248809814 seconds
Termination reason: Max number of steps (7000) reached
Steps per second = 3169.60450259376
Function evals per second = 3225.291026493065
Improvements/step = 0.17514285714285716
Total function evaluations = 7124


Best candidate found: [1.49895, 0.998787, 3.0043, 1.00127]

Fitness: 0.002280114
# Change in tolerance makes it worse
obj_short = build_loss_objective(prob_short,Vern9(),L2Loss(t_short,data_short),tstops=t_short,reltol=1e-9,abstol=1e-9)
res1 = bboptimize(obj_short;SearchRange = glo_bounds, MaxSteps = 7e3)
Starting optimization with optimizer BlackBoxOptim.DiffEvoOpt{BlackBoxOptim
.FitPopulation{Float64},BlackBoxOptim.RadiusLimitedSelector,BlackBoxOptim.A
daptiveDiffEvoRandBin{3},BlackBoxOptim.RandomBound{BlackBoxOptim.RangePerDi
mSearchSpace}}
0.00 secs, 0 evals, 0 steps
0.50 secs, 1538 evals, 1396 steps, improv/step: 0.184 (last = 0.1841), fitn
ess=696.750948315
1.00 secs, 3078 evals, 2936 steps, improv/step: 0.146 (last = 0.1123), fitn
ess=130.810573523
1.50 secs, 4677 evals, 4535 steps, improv/step: 0.150 (last = 0.1557), fitn
ess=2.455411701
2.00 secs, 6253 evals, 6111 steps, improv/step: 0.159 (last = 0.1859), fitn
ess=0.025589192

Optimization stopped after 7001 steps and 2.2848219871520996 seconds
Termination reason: Max number of steps (7000) reached
Steps per second = 3064.1336784080704
Function evals per second = 3126.28294027551
Improvements/step = 0.16157142857142856
Total function evaluations = 7143


Best candidate found: [1.49631, 0.999809, 3.02062, 1.00567]

Fitness: 0.019885287
# using the moe accurate Vern9() reduces the fitness marginally and leads to some increase in time taken

Using NLopt

Global Optimisation first

obj_short = build_loss_objective(prob_short,Vern9(),L2Loss(t_short,data_short),tstops=t_short,reltol=1e-9,abstol=1e-9)
opt = Opt(:GN_ORIG_DIRECT_L, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[10.0,10.0,10.0,10.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,glo_init)
1.216149 seconds (5.24 M allocations: 479.339 MiB, 11.61% gc time)
(368.38768828452805, [1.7284, 2.22222, 3.58025, 1.11721], :XTOL_REACHED)
opt = Opt(:GN_CRS2_LM, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[10.0,10.0,10.0,10.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,glo_init)
2.754464 seconds (11.91 M allocations: 1.065 GiB, 11.82% gc time)
(1.6661812282399235e-16, [1.5, 1.0, 3.0, 1.0], :XTOL_REACHED)
opt = Opt(:GN_ISRES, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[10.0,10.0,10.0,10.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,glo_init)
3.083688 seconds (13.41 M allocations: 1.199 GiB, 11.85% gc time)
(129.27677895286664, [1.21892, 0.987869, 5.15834, 1.6624], :MAXEVAL_REACHED
)
opt = Opt(:GN_ESCH, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[10.0,10.0,10.0,10.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,glo_init)
3.240287 seconds (13.41 M allocations: 1.199 GiB, 11.91% gc time)
(250.62230628672717, [1.09494, 0.956526, 6.8101, 2.36687], :MAXEVAL_REACHED
)

Now local optimization algorithms are used to check the global ones, these use the local constraints, different intial values and time step

opt = Opt(:LN_BOBYQA, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
0.145208 seconds (427.79 k allocations: 39.158 MiB, 17.74% gc time)
(1.6661750310234297e-16, [1.5, 1.0, 3.0, 1.0], :SUCCESS)
opt = Opt(:LN_NELDERMEAD, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
0.200497 seconds (671.86 k allocations: 61.498 MiB, 14.77% gc time)
(1.6661765163255584e-16, [1.5, 1.0, 3.0, 1.0], :XTOL_REACHED)
opt = Opt(:LD_SLSQP, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
0.174350 seconds (463.94 k allocations: 38.670 MiB, 15.49% gc time)
(4.1924574575045983e-16, [1.5, 1.0, 3.0, 1.0], :XTOL_REACHED)
opt = Opt(:LN_COBYLA, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
3.550159 seconds (13.41 M allocations: 1.199 GiB, 11.82% gc time)
(4.3894137711062486e-10, [1.5, 1.0, 3.0, 1.0], :MAXEVAL_REACHED)
opt = Opt(:LN_NEWUOA_BOUND, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
0.221720 seconds (382.20 k allocations: 34.984 MiB, 3.12% gc time)
(1.6746576191653103e-9, [1.5, 1.0, 3.0, 0.999998], :SUCCESS)
opt = Opt(:LN_PRAXIS, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
0.074558 seconds (285.65 k allocations: 26.146 MiB, 18.42% gc time)
(1.6729183094709276e-16, [1.5, 1.0, 3.0, 1.0], :SUCCESS)
opt = Opt(:LN_SBPLX, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
3.085757 seconds (13.41 M allocations: 1.199 GiB, 11.70% gc time)
(3.857624965591935e-12, [1.5, 1.0, 3.0, 1.0], :MAXEVAL_REACHED)
opt = Opt(:LD_MMA, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
5.927612 seconds (26.04 M allocations: 2.332 GiB, 11.69% gc time)
(4.632433645015572e-15, [1.5, 1.0, 3.0, 1.0], :XTOL_REACHED)
opt = Opt(:LD_TNEWTON_PRECOND_RESTART, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj_short.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
0.113176 seconds (469.54 k allocations: 43.051 MiB, 11.80% gc time)
(4.192581262191412e-16, [1.5, 1.0, 3.0, 1.0], :SUCCESS)

Now the longer problem is solved for a global solution

Vern9 solver with reltol=1e-9 and abstol=1e-9 is used and the dataset is increased to 3000 observations per variable with the same integration time step of 0.01.

obj = build_loss_objective(prob,Vern9(),L2Loss(t,data),tstops=t,reltol=1e-9,abstol=1e-9)
res1 = bboptimize(obj;SearchRange = glo_bounds, MaxSteps = 4e3)
Starting optimization with optimizer BlackBoxOptim.DiffEvoOpt{BlackBoxOptim
.FitPopulation{Float64},BlackBoxOptim.RadiusLimitedSelector,BlackBoxOptim.A
daptiveDiffEvoRandBin{3},BlackBoxOptim.RandomBound{BlackBoxOptim.RangePerDi
mSearchSpace}}
0.00 secs, 0 evals, 0 steps
0.51 secs, 160 evals, 98 steps, improv/step: 0.490 (last = 0.4898), fitness
=24393.830576126
1.01 secs, 325 evals, 236 steps, improv/step: 0.403 (last = 0.3406), fitnes
s=22582.541396785
1.51 secs, 488 evals, 392 steps, improv/step: 0.324 (last = 0.2051), fitnes
s=13765.676102501
2.01 secs, 650 evals, 548 steps, improv/step: 0.285 (last = 0.1859), fitnes
s=12722.438059839
2.51 secs, 816 evals, 712 steps, improv/step: 0.272 (last = 0.2317), fitnes
s=12722.438059839
3.01 secs, 977 evals, 869 steps, improv/step: 0.250 (last = 0.1465), fitnes
s=12722.438059839
3.52 secs, 1142 evals, 1034 steps, improv/step: 0.220 (last = 0.0606), fitn
ess=12722.438059839
4.02 secs, 1310 evals, 1202 steps, improv/step: 0.205 (last = 0.1190), fitn
ess=12722.438059839
4.52 secs, 1475 evals, 1367 steps, improv/step: 0.189 (last = 0.0727), fitn
ess=12722.438059839
5.02 secs, 1639 evals, 1531 steps, improv/step: 0.178 (last = 0.0793), fitn
ess=12722.438059839
5.53 secs, 1804 evals, 1696 steps, improv/step: 0.167 (last = 0.0667), fitn
ess=12119.173247369
6.03 secs, 1962 evals, 1854 steps, improv/step: 0.156 (last = 0.0380), fitn
ess=12119.173247369
6.53 secs, 2125 evals, 2017 steps, improv/step: 0.147 (last = 0.0429), fitn
ess=6290.592388272
7.03 secs, 2288 evals, 2180 steps, improv/step: 0.141 (last = 0.0675), fitn
ess=6290.592388272
7.53 secs, 2451 evals, 2343 steps, improv/step: 0.133 (last = 0.0307), fitn
ess=6290.592388272
8.03 secs, 2615 evals, 2507 steps, improv/step: 0.129 (last = 0.0671), fitn
ess=6290.592388272
8.53 secs, 2784 evals, 2676 steps, improv/step: 0.123 (last = 0.0414), fitn
ess=6290.592388272
9.03 secs, 2949 evals, 2841 steps, improv/step: 0.119 (last = 0.0424), fitn
ess=6290.592388272
9.53 secs, 3110 evals, 3002 steps, improv/step: 0.115 (last = 0.0435), fitn
ess=6290.592388272
10.04 secs, 3275 evals, 3167 steps, improv/step: 0.111 (last = 0.0424), fit
ness=6290.592388272
10.54 secs, 3436 evals, 3328 steps, improv/step: 0.109 (last = 0.0683), fit
ness=5302.520129161
11.04 secs, 3600 evals, 3492 steps, improv/step: 0.106 (last = 0.0427), fit
ness=5302.520129161
11.54 secs, 3766 evals, 3658 steps, improv/step: 0.104 (last = 0.0663), fit
ness=5302.520129161
12.04 secs, 3929 evals, 3821 steps, improv/step: 0.104 (last = 0.1104), fit
ness=3644.281165843
12.54 secs, 4086 evals, 3978 steps, improv/step: 0.102 (last = 0.0573), fit
ness=2022.747734424

Optimization stopped after 4001 steps and 12.616769075393677 seconds
Termination reason: Max number of steps (4000) reached
Steps per second = 317.11763733578186
Function evals per second = 325.67767353479826
Improvements/step = 0.10225
Total function evaluations = 4109


Best candidate found: [1.58137, 1.3244, 2.68362, 1.00624]

Fitness: 2022.747734424
opt = Opt(:GN_ORIG_DIRECT_L, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[10.0,10.0,10.0,10.0])
min_objective!(opt, obj.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,glo_init)
6.386443 seconds (25.43 M allocations: 2.094 GiB, 10.11% gc time)
(23525.885834893048, [8.2716, 7.42112, 7.40283, 3.7037], :XTOL_REACHED)
opt = Opt(:GN_CRS2_LM, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[10.0,10.0,10.0,10.0])
min_objective!(opt, obj.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 20000)
@time (minf,minx,ret) = NLopt.optimize(opt,glo_init)
31.580591 seconds (125.90 M allocations: 10.368 GiB, 9.99% gc time)
(1.2705965301453002e-14, [1.5, 1.0, 3.0, 1.0], :XTOL_REACHED)
opt = Opt(:GN_ISRES, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[10.0,10.0,10.0,10.0])
min_objective!(opt, obj.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 50000)
@time (minf,minx,ret) = NLopt.optimize(opt,glo_init)
145.946593 seconds (607.55 M allocations: 50.034 GiB, 10.22% gc time)
(20604.40797301213, [0.847035, 1.00954, 8.33068, 2.79177], :MAXEVAL_REACHED
)
opt = Opt(:GN_ESCH, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[10.0,10.0,10.0,10.0])
min_objective!(opt, obj.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 20000)
@time (minf,minx,ret) = NLopt.optimize(opt,glo_init)
58.578064 seconds (243.02 M allocations: 20.014 GiB, 9.87% gc time)
(4680.147063624253, [0.911551, 0.734183, 5.80582, 2.33091], :MAXEVAL_REACHE
D)
opt = Opt(:LN_BOBYQA, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
1.490528 seconds (6.49 M allocations: 547.189 MiB, 8.15% gc time)
(1.2705773526988974e-14, [1.5, 1.0, 3.0, 1.0], :SUCCESS)
opt = Opt(:LN_NELDERMEAD, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj.cost_function2)
xtol_rel!(opt,1e-9)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
1.349663 seconds (5.78 M allocations: 487.756 MiB, 8.33% gc time)
(1.7646805842304758e-14, [1.5, 1.0, 3.0, 1.0], :XTOL_REACHED)
opt = Opt(:LD_SLSQP, 4)
lower_bounds!(opt,[0.0,0.0,0.0,0.0])
upper_bounds!(opt,[5.0,5.0,5.0,5.0])
min_objective!(opt, obj.cost_function2)
xtol_rel!(opt,1e-12)
maxeval!(opt, 10000)
@time (minf,minx,ret) = NLopt.optimize(opt,loc_init)
2.877662 seconds (8.37 M allocations: 703.304 MiB, 5.48% gc time)
(21569.713608934053, [3.25971, 2.72384, 0.85817, 0.404497], :XTOL_REACHED)

Using QuadDIRECT

obj_short = build_loss_objective(prob_short,Tsit5(),L2Loss(t_short,data_short),tstops=t_short)
lower = [0.0,0.0,0.0,0.0]
upper = [5.0,5.0,5.0,5.0]
splits = ([0.0,1.0,3.0],[0.0,1.0,3.0],[0.0,1.0,3.0],[0.0,1.0,3.0])
root, x0 = analyze(obj_short,splits,lower,upper)
minimum(root)
Box1.1284001360402185@[1.49017, 0.994299, 2.99232, 1.0]
obj = build_loss_objective(prob,Vern9(),L2Loss(t,data),tstops=t,reltol=1e-9,abstol=1e-9)
lower = [0.0,0.0,0.0,0.0]
upper = [10.0,10.0,10.0,10.0]
splits = ([0.0,3.0,6.0],[0.0,3.0,6.0],[0.0,3.0,6.0],[0.0,3.0,6.0])
root, x0 = analyze(obj,splits,lower,upper)
minimum(root)
Box23222.23744589307@[2.84639, 3.0, 5.78411, 3.0]

Parameter estimation on the longer sample proves to be extremely challenging for some of the global optimizers. A few give the accurate values, BlacBoxOptim also performs quite well while others seem to struggle with accuracy a lot.

Conclusion

In general we observe that lower tolerance lead to higher accuracy but too low tolerance could affect the convergance time drastically. Also fitting a shorter timespan seems to be easier in comparision (quite intutively). NLOpt methods seem to give great accuracy in the shorter problem with a lot of the algorithms giving 0 fitness, BBO performs very well on it with marginal change with tol values. In case of global optimization of the longer problem there is some difference in the perfomance amongst the algorithms with LD_SLSQP GN_ESCH GN_ISRES GN_ORIG_DIRECT_L performing among the worse, BBO also gives a bit high fitness in comparison. QuadDIRECT gives accurate results in the case of the shorter problem but doesn't perform very well in the longer problem case.