Skip to content

Commit

Permalink
Merge pull request #40 from lrennels/override
Browse files Browse the repository at this point in the history
Add keyword args; update tests
  • Loading branch information
lrennels authored Aug 7, 2020
2 parents e4340fe + 3bc8436 commit cee1d1e
Show file tree
Hide file tree
Showing 5 changed files with 69 additions and 33 deletions.
2 changes: 1 addition & 1 deletion .appveyor.yml
Original file line number Diff line number Diff line change
@@ -1,8 +1,8 @@
environment:
matrix:
- julia_version: 1.2
- julia_version: 1.3
- julia_version: 1.4
- julia_version: 1.5
- julia_version: nightly

platform:
Expand Down
2 changes: 1 addition & 1 deletion .travis.yml
Original file line number Diff line number Diff line change
Expand Up @@ -3,9 +3,9 @@ os:
- linux
- osx
julia:
- 1.2
- 1.3
- 1.4
- 1.5
- nightly
notifications:
email: false
Expand Down
8 changes: 6 additions & 2 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -46,12 +46,16 @@ specific problem using Sobol Analysis:
After sampling with `sample`, use the resulting of matrix of parameter combinations to run your model, producing a vector of results. The next and final step is to analyze the results with your `model_output` using the `analyze` function with the signature below. This function takes the same `SobolData` as `sample`, as well as the `model_output` vector and produces a dictionary of results. This dictionary will include the `:firstorder`, `:totalorder`, and (optionally) `:secondorder` indices for each parameter.

```julia
analyze(data::SobolData, model_output::AbstractArray{<:Number, S}; num_resamples::Int = 10_000, conf_level::Number = 0.95)
function analyze(data::SobolData, model_output::AbstractArray{<:Number, S}; num_resamples::Union{Nothing, Int} = 10_000, conf_level::Union{Nothing, Number} = 0.95, progress_meter::Bool = true, N_override::Union{Nothing, Integer}=nothing)

Performs a Sobol Analysis on the `model_output` produced with the problem
defined by the information in `data` and returns the a dictionary of results
with the sensitivity indices and respective confidence intervals for each of the
parameters.
parameters defined using the `num_resamples` and `conf_level` keyword args. If these
are Nothing than no confidence intervals will be calculated. The `progress_meter`
keyword argument indicates whether a progress meter will be displayed and defaults
to true. The `N_override` keyword argument allows users to override the `N` used in
a specific `analyze` call to analyze just a subset (useful for convergence graphs).
```

An example of the basic flow can be found in `src/main.jl` using the Ishigami test function in `src/test_functions/ishigami.jl`, and is copied and commented below for convenience.
Expand Down
24 changes: 20 additions & 4 deletions src/analyze_sobol.jl
Original file line number Diff line number Diff line change
Expand Up @@ -20,17 +20,18 @@ References
=#

"""
analyze(data::SobolData, model_output::AbstractArray{<:Number, S}; num_resamples::Union{Nothing, Int} = 10_000, conf_level::Union{Nothing, Number} = 0.95, progress_meter::Bool = true) where S
analyze(data::SobolData, model_output::AbstractArray{<:Number, S}; num_resamples::Union{Nothing, Int} = 10_000, conf_level::Union{Nothing, Number} = 0.95, progress_meter::Bool = true, N_override::Union{Nothing, Integer}=nothing) where S
Performs a Sobol Analysis on the `model_output` produced with the problem
defined by the information in `data` and returns the a dictionary of results
with the sensitivity indices and respective confidence intervals for each of the
parameters defined using the `num_resamples` and `conf_level` keyword args. If these
are Nothing than no confidence intervals will be calculated. The `progress_meter`
keyword argument indicates whether a progress meter will be displayed and defaults
to true.
to true. The `N_override` keyword argument allows users to override the `N` used in
a specific `analyze` call to analyze just a subset (useful for convergence graphs).
"""
function analyze(data::SobolData, model_output::AbstractArray{<:Number, S}; num_resamples::Union{Nothing, Int} = 10_000, conf_level::Union{Nothing, Number} = 0.95, progress_meter::Bool = true) where S
function analyze(data::SobolData, model_output::AbstractArray{<:Number, S}; num_resamples::Union{Nothing, Int} = 10_000, conf_level::Union{Nothing, Number} = 0.95, progress_meter::Bool = true, N_override::Union{Nothing, Integer}=nothing) where S

# handle confidence interval flag
num_nothings = (num_resamples === nothing) + (conf_level === nothing)
Expand All @@ -45,7 +46,22 @@ function analyze(data::SobolData, model_output::AbstractArray{<:Number, S}; num_
# define constants
calc_second_order = data.calc_second_order
D = length(data.params) # number of uncertain parameters in problem
N = data.N # number of samples

# deal with overriding N
if N_override === nothing
N = data.N # number of samples
else
N_override > data.N ? error("N_override ($N_override) cannot be greater than original N used in sampling ($(data.N))") : nothing
N = N_override # number of samples

# reduce the output to just what should be considered for this N
if data.calc_second_order
lastrow = N * ((2 * D) + 2)
else
lastrow = N * (D + 2)
end
model_output = model_output[1:lastrow]
end

# values for CI calculations
if conf_flag
Expand Down
66 changes: 41 additions & 25 deletions test/unit_tests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ data5 = SobolData(N = 100)
# scale_sobol_seq!
seq = rand(100, 8)
original_seq = copy(seq)
dists = [Normal(1, 0.2), Uniform(0.75, 1.25), LogNormal(20,4), TriangularDist(0, 4, 1)]
dists = [Normal(1, 0.2), Uniform(0.75, 1.25), LogNormal(0, 0.5), TriangularDist(0, 4, 1)]
scale_sobol_seq!(seq, dists)
@test size(seq) == size(original_seq)
@test seq != original_seq
Expand Down Expand Up @@ -69,54 +69,70 @@ end
@test_throws ErrorException sample(data5) # params are nothing

##
## 4. Analyze Sobol
## 4a. Analyze Sobol
##

Y1 = ishigami(samples)
results1 = analyze(data1, Y1)
for Si in results1[:firstorder]
results = analyze(data1, Y1)
for Si in results[:firstorder]
@test Si <= 1
end
for CI in results1[:firstorder_conf]
for CI in results[:firstorder_conf]
@test CI > 0
end
for St in results1[:totalorder]
for St in results[:totalorder]
@test St <= 1
end
for CI in results1[:totalorder_conf]
for CI in results[:totalorder_conf]
@test CI > 0
end
@test sum(results1[:totalorder]) > sum(results1[:firstorder])
@test sum(results[:totalorder]) > sum(results[:firstorder])

Y3 = ishigami(samples3)
results3 = analyze(data3, Y3)
for Si in results3[:firstorder]
results = analyze(data3, Y3)
for Si in results[:firstorder]
@test Si <= 1
end
for CI in results3[:firstorder_conf]
for CI in results[:firstorder_conf]
@test ismissing(CI) || CI > 0
end
for S2i in results3[:secondorder]
for S2i in results[:secondorder]
@test ismissing(S2i) || S2i <= 1
end
for CI in results3[:secondorder_conf]
for CI in results[:secondorder_conf]
@test ismissing(CI) || CI > 0
end
for St in results3[:totalorder]
for St in results[:totalorder]
@test St <= 1
end
for CI in results3[:totalorder_conf]
for CI in results[:totalorder_conf]
@test ismissing(CI) || CI > 0
end
@test sum(results3[:totalorder]) > sum(results3[:firstorder])
@test sum(results[:totalorder]) > sum(results[:firstorder])

if(VERSION < v"1.3.0")
errortype = ArgumentError
else
errortype = ErrorException
end
@test_throws errortype results4 = analyze(data3, Y3; num_resamples = nothing)
@test_throws errortype results4 = analyze(data3, Y3; conf_level = nothing)
##
## 4b. Analyze Sobol Optional Keyword Args
##

results4 = analyze(data3, Y3; num_resamples = nothing, conf_level = nothing)
@test length(results4) == 3
data = SobolData(
params = OrderedDict(:x1 => Normal(1, 0.2),
:x2 => Uniform(0.75, 1.25),
:x3 => LogNormal(0, 0.5)),
N = 1000
)
samples = sample(data)
Y = ishigami(samples)
results = analyze(data, Y)

@test_throws ErrorException analyze(data, Y; num_resamples = nothing)
@test_throws ErrorException analyze(data, Y; conf_level = nothing)

@test length(analyze(data, Y; num_resamples = nothing, conf_level = nothing)) == 3 # no confidence intervals
results = analyze(data, Y; progress_meter = false) # no progress bar should show

@test length(analyze(data, Y; N_override = 10)) == 6
results_override = analyze(data, Y, N_override = data.N)
results_original = analyze(data, Y)
@test results_override[:firstorder] == results_original[:firstorder]
@test results_override[:totalorder] == results_original[:totalorder]
@test_throws ErrorException analyze(data1, Y1; N_override = data.N + 1) # N_override > N

0 comments on commit cee1d1e

Please sign in to comment.