diff --git a/src/measures/functioning.jl b/src/measures/functioning.jl index 52ad5a297..03ef84998 100644 --- a/src/measures/functioning.jl +++ b/src/measures/functioning.jl @@ -439,18 +439,20 @@ function alive_trophic_network(n::AbstractVector, A::AbstractMatrix; kwargs...) (; species, species_biomass, trophic_level, A) end -""" -max_trophic_level(solution; threshold, kwargs...) +docstring = """ + max_trophic_level(solution; threshold, kwargs...) + mean_trophic_level(solution; threshold, kwargs...) + weighted_average_trophic_level(solution; threshold, kwargs...) -Returns the maximum trophic level averaged over the `last` timesteps. +Return the aggregated trophic level over the `last` timesteps, +either with `max` or `mean` aggregation, +or by the mean trophic level weighted by species biomasses. kwargs... arguments are forwarded to [`extract_last_timesteps`](@ref). See [`extract_last_timesteps`](@ref) for the argument details. -`max_trophic_level()` also handles biomass vectors associated with a network, as well as +These functions also handle biomass vectors associated with a network, as well as a vector of trophic levels (See examples). -`mean_trophic_level()` and `weighted_average_trophic_level()`, respectively computing the -mean trophic level and the mean trophic level weighted by species biomasses work the same. # Examples @@ -497,58 +499,58 @@ julia> w = weighted_average_trophic_level(sol); 1.54 ``` """ -function max_trophic_level(solution::SciMLBase.ODESolution; threshold = 0, kwargs...) - measure_on = extract_last_timesteps(solution; kwargs...) - net = get_parameters(solution).network.A - out = [] - for i in 1:size(measure_on, 2) - tmp = max_trophic_level(measure_on[:, i], net; threshold) - push!(out, tmp) +# Generate the three functions in a row with macros. +function on_matrices(op_name, aggregate_function) + op_trophic_level = Symbol(op_name, :_trophic_level) + return quote + function $op_trophic_level(n::AbstractVector, A::AbstractMatrix; kwargs...) + out = alive_trophic_network(n, A; kwargs...) + tlvl = out.trophic_level + $op_trophic_level(tlvl) + end + $op_trophic_level(A::AbstractMatrix;) = $op_trophic_level(trophic_levels(A)) + $op_trophic_level(tlvl::AbstractVector) = + isempty(tlvl) ? NaN : $aggregate_function(tlvl) end - mean(out) -end -function max_trophic_level(n::AbstractVector, A::AbstractMatrix; kwargs...) - out = alive_trophic_network(n, A; kwargs...) - tlvl = out.trophic_level - max_trophic_level(tlvl) end -max_trophic_level(A::AbstractMatrix;) = max_trophic_level(trophic_levels(A)) -max_trophic_level(tlvl::AbstractVector) = isempty(tlvl) ? NaN : maximum(tlvl) -function mean_trophic_level(solution::SciMLBase.ODESolution; threshold = 0, kwargs...) - measure_on = extract_last_timesteps(solution; kwargs...) - net = get_parameters(solution).network.A +function aggregate_trophic_level(op_name, aggregate_function) + op_trophic_level = Symbol(op_name, :_trophic_level) + if isnothing(aggregate_function) + on_matrices_code = :() + else + on_matrices_code = on_matrices(op_name, aggregate_function) + end + return quote + + @doc $docstring function $op_trophic_level( + solution::SciMLBase.ODESolution; + threshold = 0, + kwargs..., + ) + measure_on = extract_last_timesteps(solution; kwargs...) + net = get_parameters(solution).network.A + + out = [] + for i in 1:size(measure_on, 2) + tmp = $op_trophic_level(measure_on[:, i], net; threshold) + push!(out, tmp) + end + mean(out) + end + + $on_matrices_code - out = [] - for i in 1:size(measure_on, 2) - tmp = mean_trophic_level(measure_on[:, i], net; threshold) - push!(out, tmp) end - mean(out) end -function mean_trophic_level(n::AbstractVector, A::AbstractMatrix; kwargs...) - out = alive_trophic_network(n, A; kwargs...) - mean_trophic_level(out.trophic_level) -end -mean_trophic_level(A::AbstractMatrix;) = mean_trophic_level(trophic_levels(A)) -mean_trophic_level(tlvl::AbstractVector) = isempty(tlvl) ? NaN : mean(tlvl) -function weighted_average_trophic_level( - solution::SciMLBase.ODESolution; - threshold = 0, - kwargs..., -) - measure_on = extract_last_timesteps(solution; kwargs...) - net = get_parameters(solution).network.A +# Generation happens here. +eval(aggregate_trophic_level(:max, :maximum)) +eval(aggregate_trophic_level(:mean, :mean)) +eval(aggregate_trophic_level(:weighted_average, nothing)) - out = [] - for i in 1:size(measure_on, 2) - tmp = weighted_average_trophic_level(measure_on[:, i], net; threshold) - push!(out, tmp) - end - mean(out) -end +# The code differs slightly in the weighted_average case. function weighted_average_trophic_level(n::AbstractVector, A::AbstractMatrix; kwargs...) out = alive_trophic_network(n, A; kwargs...) bm = out.species_biomass