From 0e1f7a8210aa19d6a0b2310e152b68e05f5fe806 Mon Sep 17 00:00:00 2001 From: Magnus Urquhart Date: Mon, 30 Mar 2020 13:00:21 +0200 Subject: [PATCH] Fixed bug when creating periodic subLHC --- Project.toml | 2 +- src/AudzeEglaisObjective.jl | 24 +++++++++++++----------- src/LatinHypercubeSampling.jl | 22 ++++++++++++++++++---- 3 files changed, 32 insertions(+), 16 deletions(-) diff --git a/Project.toml b/Project.toml index b5a9238..e52b886 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "LatinHypercubeSampling" uuid = "a5e1c1ea-c99a-51d3-a14d-a9a37257b02d" authors = ["Magnus Urquhart "] -version = "1.6.3" +version = "1.6.4" [deps] Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" diff --git a/src/AudzeEglaisObjective.jl b/src/AudzeEglaisObjective.jl index cf2135f..16d9359 100644 --- a/src/AudzeEglaisObjective.jl +++ b/src/AudzeEglaisObjective.jl @@ -1,4 +1,4 @@ -@inline function _AudzeEglaisDist(LHC,periodic_ae,ae_power) +@inline function _AudzeEglaisDist(LHC,periodic_ae,ae_power,periodic_n) n,d = size(LHC) dist = 0.0 @@ -9,7 +9,7 @@ for k = 1:d if periodic_ae @inbounds dist_comp = abs(LHC[i,k]-LHC[j,k]) - dist_comp = min(dist_comp,n-dist_comp) + dist_comp = min(dist_comp,periodic_n-dist_comp) else @inbounds dist_comp = LHC[i,k]-LHC[j,k] end @@ -28,13 +28,13 @@ return output end -function _AudzeEglaisObjective(dim::Continuous,LHC,periodic_ae,ae_power) - output = _AudzeEglaisDist(LHC,periodic_ae,ae_power) +function _AudzeEglaisObjective(dim::Continuous,LHC,periodic_ae,ae_power,periodic_n) + output = _AudzeEglaisDist(LHC,periodic_ae,ae_power,periodic_n) return output end -function _AudzeEglaisObjective(dim::Categorical,LHC,periodic_ae,ae_power) - output = _AudzeEglaisDist(LHC,periodic_ae,ae_power) +function _AudzeEglaisObjective(dim::Categorical,LHC,periodic_ae,ae_power,periodic_n) + output = _AudzeEglaisDist(LHC,periodic_ae,ae_power,periodic_n) output == Inf ? 0 : output end @@ -45,20 +45,21 @@ distance as the objective function. Note this is the inverse of the typical Audze-Eglais distance which normally is minimized. """ function AudzeEglaisObjective(LHC::T; dims::Array{V,1} =[Continuous() for i in 1:size(LHC,2)], - interSampleWeight::Float64=1.0,periodic_ae::Bool=false,ae_power::Union{Int,Float64}=2 + interSampleWeight::Float64=1.0,periodic_ae::Bool=false,ae_power::Union{Int,Float64}=2, + periodic_n::Int = size(LHC,1) ) where T <: AbstractArray where V <: LHCDimension out = 0.0 #Compute the objective function among all points - out += _AudzeEglaisObjective(Continuous(),LHC,periodic_ae,ae_power)*interSampleWeight + out += _AudzeEglaisObjective(Continuous(),LHC,periodic_ae,ae_power,periodic_n)*interSampleWeight #Compute the objective function within each categorical dimension categoricalDimInds = findall(x->typeof(x)==Categorical,dims) for i in categoricalDimInds for j = 1:dims[i].levels subLHC = @view LHC[LHC[:,i] .== j,:] - out += _AudzeEglaisObjective(dims[i],subLHC,periodic_ae,ae_power)*dims[i].weight + out += _AudzeEglaisObjective(dims[i],subLHC,periodic_ae,ae_power,periodic_n)*dims[i].weight end end @@ -67,8 +68,9 @@ end # Remove depwarning in release 2.x.x function AudzeEglaisObjective!(dist,LHC::T; dims::Array{V,1} =[Continuous() for i in 1:size(LHC,2)], - interSampleWeight::Float64=1.0,periodic_ae::Bool=false,ae_power::Union{Int,Float64}=2 + interSampleWeight::Float64=1.0,periodic_ae::Bool=false,ae_power::Union{Int,Float64}=2, + periodic_n::Int = size(LHC,1) ) where T <: AbstractArray where V <: LHCDimension @warn "AudzeEglaisObjective!(dist,LHC) is deprecated and does not differ from AudzeEglaisObjective(LHC)" - AudzeEglaisObjective(LHC; dims = dims, interSampleWeight = interSampleWeight, periodic_ae=periodic_ae, ae_power=ae_power) + AudzeEglaisObjective(LHC; dims = dims, interSampleWeight = interSampleWeight, periodic_ae=periodic_ae, ae_power=ae_power, periodic_n=periodic_n) end \ No newline at end of file diff --git a/src/LatinHypercubeSampling.jl b/src/LatinHypercubeSampling.jl index 53b8a40..16186cf 100644 --- a/src/LatinHypercubeSampling.jl +++ b/src/LatinHypercubeSampling.jl @@ -266,12 +266,20 @@ function LHCoptim!(X::Array{Int,2},gens; popsize::Int=100, end """ - function subLHCoptim(X,n::Int,gens;popsize::Int=100,ntour::Int=2,ptour::Float64=0.8,periodic_ae::Bool=false,ae_power::Union{Int,Float64}=2) + function subLHCoptim(X,n::Int,gens; popsize::Int=100, + ntour::Int=2, + ptour::Float64=0.8, + periodic_ae::Bool=false, + ae_power::Union{Int,Float64}=2) Produce an optimized Latin Hyper Cube with `n` sample points from a subset of points in `X`. Optimization is run for `gens` generations. Returns a tuple of the sample plan and the optimization fitness history. """ -function subLHCoptim(X,n::Int,gens;popsize::Int=100,ntour::Int=2,ptour::Float64=0.8,periodic_ae::Bool=false,ae_power::Union{Int,Float64}=2) +function subLHCoptim(X,n::Int,gens; popsize::Int=100, + ntour::Int=2, + ptour::Float64=0.8, + periodic_ae::Bool=false, + ae_power::Union{Int,Float64}=2) #preallocate memory nLarge, d = size(X) @@ -298,7 +306,10 @@ function subLHCoptim(X,n::Int,gens;popsize::Int=100,ntour::Int=2,ptour::Float64= for i = 1:popsize+1 subInds = sample(1:nLarge, n, replace = false) pop[i] = X[subInds,:] - fitness[i] = AudzeEglaisObjective(pop[i];periodic_ae=periodic_ae,ae_power=ae_power) + fitness[i] = AudzeEglaisObjective(pop[i]; + periodic_ae=periodic_ae, + ae_power=ae_power, + periodic_n=nLarge) end @@ -334,7 +345,10 @@ function subLHCoptim(X,n::Int,gens;popsize::Int=100,ntour::Int=2,ptour::Float64= #evaluate fitness for i = 1:popsize+1 - fitness[i] = AudzeEglaisObjective(nextpop[i];periodic_ae=periodic_ae,ae_power=ae_power) + fitness[i] = AudzeEglaisObjective(nextpop[i]; + periodic_ae=periodic_ae, + ae_power=ae_power, + periodic_n=nLarge) end #set the first individual to the best and save the fitness