diff --git a/src/decision/dMCDA.jl b/src/decision/dMCDA.jl index 1105ed721f..bf138cedcf 100644 --- a/src/decision/dMCDA.jl +++ b/src/decision/dMCDA.jl @@ -74,15 +74,15 @@ end """ DMCDA_vars(domain::Domain, criteria::NamedDimsArray, - site_ids::AbstractArray, leftover_space::AbstractArray, area_to_seed::Float64, - waves::AbstractArray, dhws::AbstractArray)::DMCDA_vars + site_ids::AbstractArray, leftover_space::AbstractArray, area_to_seed::Float64, + waves::AbstractArray, dhws::AbstractArray)::DMCDA_vars DMCDA_vars(domain::Domain, criteria::NamedDimsArray, site_ids::AbstractArray, - leftover_space::AbstractArray, area_to_seed::Float64)::DMCDA_vars + leftover_space::AbstractArray, area_to_seed::Float64)::DMCDA_vars + DMCDA_vars(domain::Domain, criteria::DataFrameRow, site_ids::AbstractArray, + leftover_space::AbstractArray, area_to_seed::Float64)::DMCDA_vars DMCDA_vars(domain::Domain, criteria::DataFrameRow, site_ids::AbstractArray, - leftover_space::AbstractArray, area_to_seed::Float64)::DMCDA_vars - DMCDA_vars(domain::Domain, criteria::DataFrameRow, site_ids::AbstractArray, - leftover_space::AbstractArray, area_to_seed::Float64, - waves::AbstractArray, dhw::AbstractArray)::DMCDA_vars + leftover_space::AbstractArray, area_to_seed::Float64, + waves::AbstractArray, dhw::AbstractArray)::DMCDA_vars Constuctors for DMCDA variables. """ @@ -161,7 +161,7 @@ function DMCDA_vars( waves::AbstractArray, dhw::AbstractArray, )::DMCDA_vars - criteria_vec::NamedDimsArray = NamedDimsArray(collect(criteria); rows = names(criteria)) + criteria_vec::NamedDimsArray = NamedDimsArray(collect(criteria); rows=names(criteria)) return DMCDA_vars( domain, criteria_vec, site_ids, leftover_space, area_to_seed, waves, dhw ) @@ -173,7 +173,7 @@ function DMCDA_vars( leftover_space::AbstractArray, area_to_seed::Float64, )::DMCDA_vars - criteria_vec::NamedDimsArray = NamedDimsArray(collect(criteria); rows = names(criteria)) + criteria_vec::NamedDimsArray = NamedDimsArray(collect(criteria); rows=names(criteria)) return DMCDA_vars(domain, criteria_vec, site_ids, leftover_space, area_to_seed) end @@ -192,7 +192,7 @@ end Normalize a Matrix (SE/SH) for MCDA. """ function mcda_normalize(x::Matrix)::Matrix - return x ./ sqrt.(sum(x .^ 2; dims = 1)) + return x ./ sqrt.(sum(x .^ 2; dims=1)) end """ @@ -201,7 +201,7 @@ end Normalize weights for a set of scenarios (wse/wsh) for MCDA. """ function mcda_normalize(x::DataFrame)::DataFrame - return x ./ sum(Matrix(x); dims = 2) + return x ./ sum(Matrix(x); dims=2) end """ @@ -240,7 +240,7 @@ function rank_sites!( mcda_func::Union{Function, Type{<:MCDMMethod}}, rank_col)::Tuple{Vector{Int64}, Matrix{Union{Float64, Int64}}} # Filter out all non-preferred sites - selector = vec(.!all(S[:, 2:end] .== 0; dims = 1)) + selector = vec(.!all(S[:, 2:end] .== 0; dims=1)) # weights in order of: in_conn, out_conn, wave, heat, predecessors, low cover weights = weights[selector] @@ -302,7 +302,7 @@ function retrieve_ranks( scores::Vector, maximize::Bool, )::Matrix{Union{Float64, Int64}} - s_order::Vector{Int64} = sortperm(scores; rev = maximize) + s_order::Vector{Int64} = sortperm(scores; rev=maximize) return Union{Float64, Int64}[Int64.(site_ids[s_order]) scores[s_order]] end @@ -315,15 +315,15 @@ Sites are then filtered based on heat and wave stress risk. Where no sites are filtered, size of ``A := n_sites × 9 criteria``. Columns indicate: -1. Site ID -2. Incoming node connectivity centrality -3. Outgoing node connectivity centrality -4. Wave stress -5. Heat stress -6. Priority predecessors -7. GBRMPA zoning criteria -8. Available area (relative to max cover) -9. Location depth + 1. Site ID + 2. Incoming node connectivity centrality + 3. Outgoing node connectivity centrality + 4. Wave stress + 5. Heat stress + 6. Priority predecessors + 7. GBRMPA zoning criteria + 8. Available area (relative to max cover) + 9. Location depth # Arguments - `site_ids` : Unique site ids. @@ -382,7 +382,7 @@ function create_decision_matrix( rule = (A[:, 4] .<= risk_tol) .& (A[:, 5] .> risk_tol) A[rule, 5] .= NaN - filtered = vec(.!any(isnan.(A); dims = 2)) + filtered = vec(.!any(isnan.(A); dims=2)) # Remove rows with NaNs A = A[filtered, :] @@ -460,7 +460,7 @@ function create_seed_matrix( SE[SE[:, 8] .<= 0.0, 8] .= NaN # Filter out sites with no space # Filter out identified locations - SE = SE[vec(.!any(isnan.(SE); dims = 2)), :] + SE = SE[vec(.!any(isnan.(SE); dims=2)), :] return SE, wse end @@ -552,10 +552,10 @@ end # Returns Tuple : - - `pref_seed_locs` : Vector, Indices of preferred seeding locations - - `pref_fog_locs` : Vector, Indices of preferred shading locations - - `rankings` : Matrix[n_sites ⋅ 3] where columns are site_id, seeding_rank, shading_rank - Values of 0 indicate sites that were not considered +- `pref_seed_locs` : Vector, Indices of preferred seeding locations +- `pref_fog_locs` : Vector, Indices of preferred shading locations +- `rankings` : Matrix[n_sites ⋅ 3] where columns are site_id, seeding_rank, shading_rank +Values of 0 indicate sites that were not considered """ function guided_site_selection( d_vars::DMCDA_vars, @@ -568,8 +568,10 @@ function guided_site_selection( in_conn::Vector{Float64}, out_conn::Vector{Float64}, strong_pred::Vector{Int64}; - methods_mcda = mcda_methods() -)::Tuple{Vector{T}, Vector{T}, Matrix{T}} where { + methods_mcda=mcda_methods() +)::Tuple{ + Vector{T}, Vector{T}, Matrix{T} +} where { T <: Int64, IA <: AbstractArray{<:Int64}, IB <: AbstractArray{<:Int64}, B <: Bool } site_ids = copy(d_vars.site_ids) @@ -585,7 +587,9 @@ function guided_site_selection( end n_iv_locs::Int64 = d_vars.n_site_int - priority_sites::Array{Int64} = d_vars.priority_sites[in.(d_vars.priority_sites, [site_ids])] + priority_sites::Array{Int64} = d_vars.priority_sites[in.( + d_vars.priority_sites, [site_ids] + )] priority_zones::Array{String} = d_vars.priority_zones zones = d_vars.zones[site_ids] @@ -612,7 +616,9 @@ function guided_site_selection( zone_preds_temp::Vector{Int64} = strong_pred[zones .== z_name] for s::Int64 in unique(zone_preds_temp) # for each predecessor site, add zone_weights * (no. of zone sites the site is a strongest predecessor for) - zone_preds[site_ids .== s] .= zone_preds[site_ids .== s] .+ (zone_weights[k] .* sum(zone_preds_temp .== s)) + zone_preds[site_ids .== s] .= + zone_preds[site_ids .== s] .+ + (zone_weights[k] .* sum(zone_preds_temp .== s)) end # add zone_weights for sites in the zone (whether a strongest predecessor of a zone or not) zone_sites[zones .== z_name] .= zone_weights[k] @@ -732,7 +738,7 @@ end Tuple : - `pref_locs` : Vector, Indices of preferred intervention locations - `rankings` : Matrix[n_sites ⋅ 3] where columns are site_id, seeding_rank, shading_rank - Values of 0 indicate sites that were not considered +Values of 0 indicate sites that were not considered """ function constrain_reef_cluster( reefs::Union{Vector{String}, Vector{Float64}}, @@ -753,7 +759,9 @@ function constrain_reef_cluster( local num_locs::Int64 for _ in 1:max_iters # If enough space for seeding corals, keep n_site_int, else expand as needed - num_locs = max(findfirst(>=(area_to_seed), cumsum(available_space[loc_ordered_ids])), n_iv_locs) + num_locs = max( + findfirst(>=(area_to_seed), cumsum(available_space[loc_ordered_ids])), n_iv_locs + ) pref_locs = loc_ordered_ids[1:num_locs] @@ -766,7 +774,7 @@ function constrain_reef_cluster( pref_reefs = reefs[pref_locs] # Reefs that selected locations sit within # Number of times a reef appears within each location - reef_occurances = vec(sum(pref_reefs .== unique_reefs; dims = 1)) + reef_occurances = vec(sum(pref_reefs .== unique_reefs; dims=1)) # If more than n_reefs locations in a reef, swap out the worst locations reefs_swap = unique_reefs[(reef_occurances .> max_members)] @@ -779,9 +787,11 @@ function constrain_reef_cluster( # Find locations in reefs which need replacement, and find the ids of lowest # ranked locations in this set - locs_to_replace = vcat([ - pref_locs[pref_reefs .== reef][replace_start:end] for reef in reefs_swap - ]...) + locs_to_replace = vcat( + [ + pref_locs[pref_reefs .== reef][replace_start:end] for reef in reefs_swap + ]..., + ) # Acceptable reefs to switch out for reef_switch_ids = unique_reefs[(reef_occurances .+ 1) .<= max_members] @@ -797,10 +807,13 @@ function constrain_reef_cluster( # Indices of the subset of locations which can be added which also sit within an # allowed reef - add_locs_ind = findall(dropdims(any( - reshape(reefs[alternate_loc_ids], 1, length(reefs[alternate_loc_ids])) - .== - reef_switch_ids; dims = 1); dims = 1)) + add_locs_ind = findall( + dropdims( + any( + reshape(reefs[alternate_loc_ids], 1, length(reefs[alternate_loc_ids])) + .== + reef_switch_ids; dims=1); dims=1), + ) # New preferred location set locs_to_add_inds = add_locs_ind[1:length(locs_to_replace)] @@ -864,13 +877,15 @@ function unguided_site_selection( if seed_years pref_seed_locs = zeros(Int64, n_site_int) - pref_seed_locs[1:s_n_site_int] .= StatsBase.sample(candidate_sites, s_n_site_int; replace = false) + pref_seed_locs[1:s_n_site_int] .= StatsBase.sample( + candidate_sites, s_n_site_int; replace=false + ) end if fog_years pref_fog_locs = zeros(Int64, n_site_int) pref_fog_locs[1:s_n_site_int] .= StatsBase.sample( - candidate_sites, s_n_site_int; replace = false + candidate_sites, s_n_site_int; replace=false ) end @@ -890,14 +905,16 @@ Calculates mean over specified dimensions plus half the standard deviation. # Returns Weighted combination of mean and standard deviation of the projected environmental conditions (e.g., DHWs, wave stress, etc): - (μ * w) + (σ * (1 - w)) +(μ * w) + (σ * (1 - w)) """ function summary_stat_env( env_layer::AbstractArray, dims::Union{Int64, Symbol, Tuple{Symbol, Symbol}}; - w = 0.5, + w=0.5, )::Vector{Float64} - return vec((mean(env_layer; dims = dims) .* w) .+ (std(env_layer; dims = dims) .* (1.0 - w))) + return vec( + (mean(env_layer; dims=dims) .* w) .+ (std(env_layer; dims=dims) .* (1.0 - w)) + ) end """ diff --git a/test/site_selection.jl b/test/site_selection.jl index 321242ce0c..e50e650ca4 100644 --- a/test/site_selection.jl +++ b/test/site_selection.jl @@ -9,90 +9,91 @@ if !@isdefined(ADRIA_DIR) end function test_site_ranks( - weights_set::Vector{Vector{Float64}}, - A::Matrix{Float64}, - rankings::Matrix{Int64}, - n_site_int::Int64, - mcda_func::Function, - inv::Float64, + weights_set::Vector{Vector{Float64}}, + A::Matrix{Float64}, + rankings::Matrix{Int64}, + n_site_int::Int64, + mcda_func::Function, ) - S = ADRIA.decision.mcda_normalize(A) - S[:, 1] .= A[:, 1] - criteria_names = [ - "heat stress", - "wave stress", - "median depth", - "coral cover space", - "in connectivity", - "out connectivity", - ] - for weights in weights_set - crit_inds = findall(weights .> 0.0) - - prefsites, s_order = ADRIA.decision.rank_sites!( - S, weights, rankings, n_site_int, mcda_func, inv - ) - - names_temp = criteria_names[crit_inds] - names_string = string(["$(name), " for name in names_temp[1:(end - 1)]]...) - - # Check that 2 best sites are selected (5 and 6) - @test any(prefsites .== 5) & any(prefsites .== 6) || string( - "For the ", - names_string, - "and $(names_temp[end]) criteria, the best sites (5 and 6) were not selected.", - ) - # Check that 2 worst sites aren't selected (9 and 10) - @test any(prefsites .!= 9) & any(prefsites .!= 10) || string( - "For the ", - names_string, - "and $(names_temp[end]) criteria, the worst sites (9 and 10) were selected.", - ) - end + S = ADRIA.decision.mcda_normalize(A) + S[:, 1] .= A[:, 1] + criteria_names = [ + "heat stress", + "wave stress", + "median depth", + "coral cover space", + "in connectivity", + "out connectivity", + ] + for weights in weights_set + # Find criteria being used for decision + crit_inds = findall(weights .> 0.0) + + # Get site preference order + prefsites, s_order = ADRIA.decision.rank_sites!( + S, weights, rankings, n_site_int, mcda_func, 2 + ) + + # Get names of criteria being used for error message + names_temp = criteria_names[crit_inds] + names_string = string(["$(name), " for name in names_temp[1:(end - 1)]]...) + + # Check that 2 best sites are selected (5 and 6) + @test any(prefsites .== 5) & any(prefsites .== 6) || string( + "For the ", + names_string, + "and $(names_temp[end]) criteria, the best sites (5 and 6) were not selected.", + ) + # Check that 2 worst sites aren't selected (9 and 10) + @test any(prefsites .!= 9) & any(prefsites .!= 10) || string( + "For the ", + names_string, + "and $(names_temp[end]) criteria, the worst sites (9 and 10) were selected.", + ) + end end function test_mcda_funcs(rankings::Matrix{Int64}, S::Matrix{Float64}, - weights::Vector{Float64}, - mcda_funcs, n_site_int::Int64) - for mf in mcda_funcs - prefsites, s_order = ADRIA.decision.rank_sites!( - S, weights, rankings, n_site_int, mf, 2 - ) - # Check that 2 best sites are selected (5 and 6) - @test in(5, prefsites) .& in(6, prefsites) || - "The best overall sites (5 and 6) were not chosen by method $mf." - - # Check that 2 worst sites aren't selected (9 and 10) - @test any(prefsites .!= 9) & any(prefsites .!= 10) || - "The worst overall sites (9 and 10) were chosen by method $mf." - end + weights::Vector{Float64}, mcda_funcs, n_site_int::Int64) + for mf in mcda_funcs + prefsites, s_order = ADRIA.decision.rank_sites!( + S, weights, rankings, n_site_int, mf, 2 + ) + # Check that 2 best sites are selected (5 and 6) + @test in(5, prefsites) .& in(6, prefsites) || + "The best overall sites (5 and 6) were not chosen by method $mf." + + # Check that 2 worst sites aren't selected (9 and 10) + @test any(prefsites .!= 9) & any(prefsites .!= 10) || + "The worst overall sites (9 and 10) were chosen by method $mf." + end end function get_test_decision_matrix(dom::Domain)::Matrix{Float64} - cover = sum(dom.init_coral_cover; dims = :species)[species = 1] - leftover_space = 1 .- cover - k_area = dom.site_data.area .* dom.site_data.k - dhw_av = ADRIA.decision.summary_stat_env(dom.dhw_scens, (:timesteps, :scenarios)) - wave_av = ADRIA.decision.summary_stat_env(dom.wave_scens, (:timesteps, :scenarios)) - depth_med = dom.site_data.depth_med - - TP_data = ADRIA.connectivity_strength( - dom.TP_data .* ADRIA.site_k_area(dom), collect(cover), dom.TP_data - ) - - site_ids = dom.site_data.site_id - - heat_stress = - 1 .- vec((dhw_av .- minimum(dhw_av)) ./ (maximum(dhw_av) - minimum(dhw_av))) - wave_stress = - 1 .- vec((wave_av .- minimum(wave_av)) ./ (maximum(wave_av) - minimum(wave_av))) - space_area = leftover_space .* k_area - in_conn = TP_data.in_conn - out_conn = TP_data.out_conn - - A = hcat(site_ids, heat_stress, wave_stress, depth_med, space_area, in_conn, out_conn) - - return A + cover = sum(dom.init_coral_cover; dims=:species)[species=1] + leftover_space = 1 .- cover + k_area = dom.site_data.area .* dom.site_data.k + dhw_av = ADRIA.decision.summary_stat_env(dom.dhw_scens, (:timesteps, :scenarios)) + wave_av = ADRIA.decision.summary_stat_env(dom.wave_scens, (:timesteps, :scenarios)) + depth_med = dom.site_data.depth_med + + TP_data = ADRIA.connectivity_strength( + dom.TP_data .* ADRIA.site_k_area(dom), collect(cover), dom.TP_data + ) + + heat_stress = + 1 .- vec((dhw_av .- minimum(dhw_av)) ./ (maximum(dhw_av) - minimum(dhw_av))) + wave_stress = + 1 .- vec((wave_av .- minimum(wave_av)) ./ (maximum(wave_av) - minimum(wave_av))) + space_area = leftover_space .* k_area + in_conn = TP_data.in_conn + out_conn = TP_data.out_conn + + site_ids = dom.site_data.site_id + # Simplified decision matrix for testing + A = hcat(site_ids, heat_stress, wave_stress, depth_med, space_area, in_conn, out_conn) + + return A end @testset "site selection" begin @@ -127,9 +128,11 @@ end dhw_scens, ) n_sites = length(mcda_vars.site_ids) - @test (size(mcda_vars.conn, 1) == n_sites) && (size(mcda_vars.conn, 2) == n_sites) || "Connectivity input is incorrect size." + @test (size(mcda_vars.conn, 1) == n_sites) && (size(mcda_vars.conn, 2) == n_sites) || + "Connectivity input is incorrect size." @test length(mcda_vars.dam_prob) == n_sites || "Wave damage input is incorrect size." - @test length(mcda_vars.heat_stress_prob) == n_sites || "Heat stress input is incorrect size." + @test length(mcda_vars.heat_stress_prob) == n_sites || + "Heat stress input is incorrect size." @test length(mcda_vars.leftover_space) == n_sites || "Initial cover input is incorrect size." end @@ -161,9 +164,10 @@ end end @testset "Guided site selection without ADRIA ecological model" begin - dom = ADRIA.load_domain(EXAMPLE_DOMAIN_PATH, 45) + dom = ADRIA.load_domain(TEST_DOMAIN_PATH, 45) site_ids = collect(1:length(dom.site_data.site_id)) - A = get_test_decision_matrix(dom) + mcda_funcs = ADRIA.decision.mcda_methods() + A = get_test_decision_matrix(dom) n_sites = length(site_ids) n_site_int = 5 @@ -171,27 +175,27 @@ end # Test 0.0 or 1.0 weights in all combinations weights = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0] - weights_choice = collect(0.0:0.1:1.0) for num_crit in 1:6 weights[num_crit] = 1.0 weights_set = collect(distinct(permutations(weights, 6))) test_site_ranks( - weights_set, A, rankings, criteria_names, n_site_int, mcda_funcs[1], 1 + weights_set, A, rankings, n_site_int, mcda_funcs[1] ) end + + # Test randomised weights in all combinations weights = rand(100, 6) weights = weights ./ sum(weights; dims=2) for ww in eachrow(weights) weights_set = collect(distinct(permutations(ww, 6))) test_site_ranks( - weights_set, A, rankings, criteria_names, n_site_int, mcda_funcs[1], 1 + weights_set, A, rankings, n_site_int, mcda_funcs[1] ) end - end -@testset "Test ranks line up with ordering" begin +@testset "Ranks line up with ordering" begin supported_methods = ADRIA.decision.mcda_methods() mcda_func = supported_methods[rand(1:length(supported_methods))] n_sites = 20 @@ -209,15 +213,17 @@ end S, weights, rankings, n_site_int, mcda_func, 2 ) - @test all([(rankings[rankings[:, 1].==s_order[rank, 1], 2].==rank)[1] for rank in 1:size(s_order, 1)]) || "Ranking does not match mcda score ordering" - + @test all([ + (rankings[rankings[:, 1] .== s_order[rank, 1], 2] .== rank)[1] for + rank in 1:size(s_order, 1) + ]) || "Ranking does not match mcda score ordering" end -@testset "Test each mcda method is working" begin +@testset "MCDA methods" begin mcda_funcs = ADRIA.decision.mcda_methods() - dom = ADRIA.load_domain(EXAMPLE_DOMAIN_PATH, 45) - A = get_test_decision_matrix(dom) + dom = ADRIA.load_domain(TEST_DOMAIN_PATH, 45) + A = get_test_decision_matrix(dom) site_ids = dom.site_data.site_id n_sites = length(site_ids) @@ -231,4 +237,4 @@ end weights = [1.0, 1.0, 1.0, 1.0, 0.0, 0.0] test_mcda_funcs(rankings, S, weights, mcda_funcs, n_site_int) -end \ No newline at end of file +end