From e6b997b44b632f2ba9d29fd9d151daa507aa5be0 Mon Sep 17 00:00:00 2001 From: Art Wild Date: Wed, 8 Dec 2021 17:50:03 -0500 Subject: [PATCH] allow better control over used RNG (part 1) - added `rng` parameter to `Options` - wired calls that use RNG to use the option `rng` parameter in GA & GP --- src/api/types.jl | 1 + src/ga.jl | 14 +++++++------- src/gp.jl | 20 +++++++++++++------- test/gp.jl | 26 +++++++++++++------------- test/rosenbrock.jl | 10 +++++----- 5 files changed, 39 insertions(+), 32 deletions(-) diff --git a/src/api/types.jl b/src/api/types.jl index 7c88f1e..18ad065 100644 --- a/src/api/types.jl +++ b/src/api/types.jl @@ -36,6 +36,7 @@ There are following options available: callback::TCallback = nothing time_limit::Float64 = NaN parallelization::Symbol = :serial + rng::AbstractRNG = Random.GLOBAL_RNG end function show(io::IO, o::Options) for k in fieldnames(typeof(o)) diff --git a/src/ga.jl b/src/ga.jl index e471173..4fa723e 100644 --- a/src/ga.jl +++ b/src/ga.jl @@ -61,19 +61,19 @@ end function update_state!(objfun, constraints, state, population::AbstractVector{IT}, method::GA, options, itr) where {IT} @unpack populationSize,crossoverRate,mutationRate,ɛ,selection,crossover,mutation = method evaltype = options.parallelization - + rng = options.rng offspring = similar(population) # Select offspring - selected = selection(state.fitpop, populationSize) + selected = selection(state.fitpop, populationSize, rng=rng) # Perform mating - offidx = randperm(populationSize) + offidx = randperm(rng, populationSize) offspringSize = populationSize - state.eliteSize for i in 1:2:offspringSize j = (i == offspringSize) ? i-1 : i+1 - if rand() < crossoverRate - offspring[i], offspring[j] = crossover(population[selected[offidx[i]]], population[selected[offidx[j]]]) + if rand(rng) < crossoverRate + offspring[i], offspring[j] = crossover(population[selected[offidx[i]]], population[selected[offidx[j]]], rng=rng) else offspring[i], offspring[j] = population[selected[i]], population[selected[j]] end @@ -88,8 +88,8 @@ function update_state!(objfun, constraints, state, population::AbstractVector{IT # Perform mutation for i in 1:offspringSize - if rand() < mutationRate - mutation(offspring[i]) + if rand(rng) < mutationRate + mutation(offspring[i], rng=rng) end end diff --git a/src/gp.jl b/src/gp.jl index c498dab..1313d74 100644 --- a/src/gp.jl +++ b/src/gp.jl @@ -55,7 +55,7 @@ function randterm(rng::AbstractRNG, t::TreeGP) if isa(term, Symbol) || isa(term, Real) term elseif isa(term, Function) - term() + term(rng) # terminal functions must accept RNG as an argument else # Code shouldn't reach branch but left as a catchall dim = t.terminals[term] @@ -97,10 +97,11 @@ rand(t::TreeGP, maxdepth::Int=2; kwargs...) = Initialize a random population of expressions derived from `expr`. """ -function initial_population(m::TreeGP, expr::Union{Expr,Nothing}=nothing) +function initial_population(m::TreeGP, expr::Union{Expr,Nothing}=nothing; + rng::AbstractRNG=Random.GLOBAL_RNG) n = population_size(m) return [ - expr === nothing ? rand(m, m.maxdepth, mindepth=m.mindepth) : deepcopy(expr) + expr === nothing ? rand(rng, m, m.maxdepth, mindepth=m.mindepth) : deepcopy(expr) for i in 1:n ] end @@ -129,8 +130,13 @@ function update_state!(objfun, constraints, state::GPState, population::Abstract end # Custom optimization call -function optimize(f, method::TreeGP, options::Options = Options(;default_options(method)...)) - method.optimizer.mutation = subtree(method) - method.optimizer.crossover = crosstree - optimize(f, NoConstraints(), nothing, method, options) +function optimize(f, mthd::TreeGP, options::Options = Options(;default_options(mthd)...)) + optimize(f, mthd, initial_population(mthd, rng=options.rng), options) end + +function optimize(f, mthd::TreeGP, population, options::Options = Options(;default_options(mthd)...)) + mthd.optimizer.mutation = subtree(mthd) + mthd.optimizer.crossover = crosstree + optimize(f, NoConstraints(), mthd, population, options) +end + diff --git a/test/gp.jl b/test/gp.jl index 7f6f84b..38b8d85 100644 --- a/test/gp.jl +++ b/test/gp.jl @@ -1,5 +1,6 @@ @testset "Genetic Programming" begin - Random.seed!(9874984737486); + rng = StableRNG(42) + pop = 10 terms = Terminal[:x, :y, rand] funcs = Function[+,-,*,/] @@ -16,18 +17,18 @@ @test_skip summary(t) == "TreeGP[P=10,Parameter[x,y],Function[*, +, /, -]]" # population initialization - popexp = Evolutionary.initial_population(t); + popexp = Evolutionary.initial_population(t, rng=rng); @test length(popexp) == pop popexp = Evolutionary.initial_population(t, :(x + 1)); @test popexp[1] == :(x + 1) # recursive helper functions - Random.seed!(9874984737482) - gt = rand(TreeGP(pop, terms, funcs, maxdepth=2, initialization=:grow), 3) + Random.seed!(rng, 1) + gt = rand(rng, TreeGP(pop, terms, funcs, maxdepth=2, initialization=:grow), 3) @test Evolutionary.nodes(gt) < 15 @test Evolutionary.height(gt) <= 3 @test length(gt) < 15 - ft = rand(TreeGP(pop, terms, funcs, maxdepth=2, initialization=:full), 3) + ft = rand(rng, TreeGP(pop, terms, funcs, maxdepth=2, initialization=:full), 3) @test Evolutionary.nodes(ft) == 15 @test Evolutionary.height(ft) == 3 @test length(ft) == 15 @@ -57,7 +58,7 @@ # evaluation ex = Expr(:call, +, 1, :x) |> Evolutionary.Expression - xs = rand(10) + xs = rand(rng, 10) @test ex(xs[1]) == xs[1]+1 @test ex.(xs) == xs.+1 io = IOBuffer() @@ -67,24 +68,23 @@ depth = 5 fitfun(x) = x*x + x + 1.0 function fitobj(expr) - rg = -5.0:0.5:5.0 + rg = -5.0:0.1:5.0 ex = Evolutionary.Expression(expr) - sum(v->isnan(v) ? 1.0 : v, abs2.(fitfun.(rg) - ex.(rg)) )/length(rg) |> sqrt + sum(v->isnan(v) ? 1.0 : v, abs2.(fitfun.(rg) - ex.(rg)) )/length(rg) end - Random.seed!(9874984737426); res = Evolutionary.optimize(fitobj, TreeGP(50, Terminal[:x, randn], Function[+,-,*,Evolutionary.pdiv], mindepth=1, maxdepth=depth, optimizer = GA( - selection = tournament(5), + selection = tournament(4), ɛ = 0.1, - mutationRate = 0.85, + mutationRate = 0.8, crossoverRate = 0.2, ), - ) + ), Evolutionary.Options(rng=StableRNG(1)) ) - @test minimum(res) < 1.1 + @test minimum(res) < 1 end diff --git a/test/rosenbrock.jl b/test/rosenbrock.jl index 8259642..8113da1 100644 --- a/test/rosenbrock.jl +++ b/test/rosenbrock.jl @@ -52,12 +52,12 @@ c = PenaltyConstraints(100.0, fill(0.0, 5N), Float64[], [1.0], [1.0], con_c!) result = Evolutionary.optimize(rosenbrock, c, (() -> rand(5N)), CMAES(mu = 40, lambda = 100)) println("(5/5,100)-CMA-ES [penalty] => F: $(minimum(result)), C: $(Evolutionary.iterations(result))") - @test Evolutionary.minimizer(result) |> sum ≈ 1.0 atol=1e-1 + @test Evolutionary.minimizer(result) |> sum ≈ 1.0 atol=0.1 c = PenaltyConstraints(100.0, fill(0.0, 2N), fill(0.5, 2N), [1.0], [1.0], con_c!) result = Evolutionary.optimize(rosenbrock, c, (() -> rand(2N)), CMAES(mu = 8, lambda = 100)) println("(5/5,100)-CMA-ES [penalty] => F: $(minimum(result)), C: $(Evolutionary.iterations(result))") - @test Evolutionary.minimizer(result) |> sum ≈ 1.0 atol=1e-1 + @test Evolutionary.minimizer(result) |> sum ≈ 1.0 atol=0.1 @test all(0.0 <= x+0.01 && x-0.01 <= 0.5 for x in abs.(Evolutionary.minimizer(result))) # Testing: GA @@ -65,15 +65,15 @@ populationSize = 100, ɛ = 0.1, selection = rouletteinv, - crossover = intermediate(0.25), + crossover = IC(0.2), mutation = BGA(fill(0.5,N)) ) result = Evolutionary.optimize(rosenbrock, (() -> rand(N)), m) println("GA(p=100,x=0.8,μ=0.1,ɛ=0.1) => F: $(minimum(result)), C: $(Evolutionary.iterations(result))") - test_result(result, N, 1e-1) + test_result(result, N, 0.1) result = Evolutionary.optimize(rosenbrock, BoxConstraints(0.0, 0.5, N), (() -> rand(N)), m) println("GA(p=100,x=0.8,μ=0.1,ɛ=0.1)[box] => F: $(minimum(result)), C: $(Evolutionary.iterations(result))") - @test Evolutionary.minimizer(result) ≈ [0.5, 0.25] atol=1e-1 + @test Evolutionary.minimizer(result) ≈ [0.5, 0.25] atol=0.1 # Testing: DE result = Evolutionary.optimize(rosenbrock, (() -> rand(N)), DE(populationSize = 100))