Skip to content

Commit

Permalink
Minor style nitpicks.
Browse files Browse the repository at this point in the history
- Vertical whitespace.
- Multiline strings.
- `; arg = arg` elisions.
- ..
  • Loading branch information
iago-lito committed Apr 18, 2023
1 parent c4fc87f commit d5b316c
Show file tree
Hide file tree
Showing 6 changed files with 135 additions and 170 deletions.
34 changes: 7 additions & 27 deletions src/measures/functioning.jl
Original file line number Diff line number Diff line change
@@ -1,6 +1,5 @@
#=
Quantifying functions
Adapted from BioenergeticFoodWeb.jl
Quantifying functions.
=#

"""
Expand Down Expand Up @@ -167,14 +166,12 @@ https://en.wikipedia.org/wiki/Diversity_index#Shannon_index
"""
function shannon_diversity(solution; threshold = 0, kwargs...)
measure_on = extract_last_timesteps(solution; kwargs...)

shan = shannon_diversity.(eachcol(measure_on); threshold)
mean(shan)
end

function shannon_diversity(n::AbstractVector; threshold = 0)
x = filter(>(threshold), n)

if length(x) >= 1
p = x ./ sum(x)
p_ln_p = p .* log.(p)
Expand Down Expand Up @@ -203,14 +200,12 @@ https://en.wikipedia.org/wiki/Diversity_index#Simpson_index
"""
function simpson(solution; threshold = 0, kwargs...)
measure_on = extract_last_timesteps(solution; kwargs...)

simp = simpson.(eachcol(measure_on); threshold)
mean(simp)
end

function simpson(n::AbstractVector; threshold = 0)
x = filter(>(threshold), n)

if length(x) >= 1
p = x ./ sum(x)
p2 = 2 .^ p
Expand All @@ -237,7 +232,6 @@ https://en.wikipedia.org/wiki/Species_evenness
"""
function evenness(solution; threshold = 0, kwargs...)
measure_on = extract_last_timesteps(solution; kwargs...)

piel = evenness.(eachcol(measure_on); threshold)
mean(piel)
end
Expand Down Expand Up @@ -371,10 +365,10 @@ julia> foodweb = FoodWeb([0 0 0; 0 1 0; 1 1 0]; quiet = true);
"""
function trophic_structure(solution; threshold = 0, idxs = nothing, kwargs...)

isnothing(idxs) || throw(ArgumentError("`trophic_structure()` operates at the whole \
network level, so it makes no sense to ask for \
particular species with anything other than \
`idxs = nothing`."))
isnothing(idxs) ||
throw(ArgumentError("`trophic_structure()` operates at the whole network level, \
so it makes no sense to ask for particular species \
with anything other than `idxs = nothing`."))


# Measure trophic structure over last timesteps
Expand All @@ -388,7 +382,6 @@ function trophic_structure(solution; threshold = 0, idxs = nothing, kwargs...)
alive = alive_trophic_network(measure_on[:, i], net; threshold)
tlvl = alive.trophic_level
bm = alive.species_biomass

push!(maxl, max_trophic_level(tlvl))
push!(avgl, mean_trophic_level(tlvl))
push!(wavg, weighted_average_trophic_level(bm, tlvl))
Expand Down Expand Up @@ -568,12 +561,7 @@ end


"""
living_species(
solution::Solution;
threshold = 0,
idxs = nothing,
kwargs...,
)
living_species(solution::Solution; threshold = 0, idxs = nothing, kwargs...)
Returns the vectors of alive species and their indices in the original network.
Living species are the ones having, in average, a biomass above `threshold` over
Expand All @@ -600,17 +588,10 @@ julia> B0 = [0, 0.5, 0.5];
(species = ["s2", "s3"], idxs = [2, 3])
```
"""
function living_species(
solution::Solution;
threshold = 0,
idxs = nothing,
kwargs...,
)
function living_species(solution::Solution; threshold = 0, idxs = nothing, kwargs...)

measure_on = extract_last_timesteps(solution; idxs, kwargs...)

alive_sp = living_species(measure_on; threshold)

sp = get_parameters(solution).network.species

tmp_idxs = process_idxs(solution; idxs)
Expand All @@ -622,7 +603,6 @@ end

living_species(mat::AbstractMatrix; threshold = 0) =
findall(>(threshold), biomass(mat).species)

living_species(n::AbstractVector; threshold = 0) = findall(>(threshold), n)

"""
Expand Down
19 changes: 7 additions & 12 deletions src/measures/stability.jl
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
#=
Various measures of stability
Various measures of stability.
=#

"""
Expand All @@ -11,6 +11,7 @@ See [`coefficient_of_variation`](@ref) for the details.
"""
function species_cv(solution::Solution; threshold = 0, last = "10%", kwargs...)
measure_on = extract_last_timesteps(solution; last, kwargs...)

# Fetch species that are alive, whose avg biomass is > threshold
living_sp = living_species(measure_on; threshold)

Expand All @@ -21,7 +22,6 @@ function species_cv(solution::Solution; threshold = 0, last = "10%", kwargs...)
end

function species_cv(mat::AbstractMatrix; corrected = true)

if any(size(mat) .== 0)
average = NaN
species = NaN
Expand Down Expand Up @@ -65,19 +65,18 @@ function synchrony(
corrected = true,
kwargs...,
)

measure_on = extract_last_timesteps(solution; last, kwargs...)

# Fetch species that are alive, whose avg biomass is > threshold
living_sp = living_species(measure_on; threshold)

# Transpose to get the time x species matrix
mat = transpose(measure_on[living_sp, :])

synchrony(mat; corrected = corrected)

end
function synchrony(mat::AbstractMatrix; corrected = true)

function synchrony(mat::AbstractMatrix; corrected = true)
if any(size(mat) .== 0)
phi = NaN
else
Expand All @@ -88,9 +87,7 @@ function synchrony(mat::AbstractMatrix; corrected = true)

phi = com_var / std_sp^2
end

phi

end

"""
Expand All @@ -113,19 +110,18 @@ function community_cv(
corrected = true,
kwargs...,
)
measure_on = extract_last_timesteps(solution; last, kwargs...)

measure_on = extract_last_timesteps(solution; last = last, kwargs...)
# Fetch species that are alive, whose avg biomass is > threshold
living_sp = living_species(measure_on; threshold = threshold, kwargs...)
living_sp = living_species(measure_on; threshold, kwargs...)

# Transpose to get the time x species matrix
mat = transpose(measure_on[living_sp, :])

community_cv(mat; corrected)

end
function community_cv(mat::AbstractMatrix; corrected = true)

function community_cv(mat::AbstractMatrix; corrected = true)
if any(size(mat) .== 0)
cv_com = NaN
else
Expand Down Expand Up @@ -187,7 +183,6 @@ function coefficient_of_variation(
corrected = true,
kwargs...,
)

measure_on = extract_last_timesteps(solution; last, kwargs...)
# Fetch species that are alive, whose avg biomass is > threshold
alive_sp = living_species(measure_on; threshold)
Expand Down
104 changes: 50 additions & 54 deletions src/measures/utils.jl
Original file line number Diff line number Diff line change
Expand Up @@ -29,10 +29,9 @@ julia> fw = FoodWeb([0 0; 1 0]);
[0.219659439; 0.188980349;;]
true
julia> sim = extract_last_timesteps(m; last = 1, idxs = [1, 2]);
sim ==
extract_last_timesteps(m; last = 1, idxs = ["s1", "s2"]) ≈
[0.188980349; 0.219659439;;]
julia> sim1 = extract_last_timesteps(m; last = 1, idxs = [1, 2]);
sim2 = extract_last_timesteps(m; last = 1, idxs = ["s1", "s2"]);
sim1 ≈ sim2 ≈ [0.188980349; 0.219659439;;]
true
julia> sim = extract_last_timesteps(m; last = 1, idxs = [2]);
Expand All @@ -41,17 +40,13 @@ true
```
"""
function extract_last_timesteps(solution; idxs = nothing, quiet = false, kwargs...)

last = process_last_timesteps(solution; quiet = quiet, kwargs...)


last = process_last_timesteps(solution; quiet, kwargs...)
out = solution[:, end-(last-1):end]

idxs = process_idxs(solution; idxs)
# Extract species indices:
idxs = process_idxs(solution; idxs)
out = out[idxs, :]

quiet || check_last_extinction(solution; idxs = idxs, last = last)
quiet || check_last_extinction(solution; idxs, last)

deepcopy(out)
end
Expand All @@ -78,24 +73,27 @@ function process_idxs(solution; idxs = nothing)
idxs_in = indexin(idxs, sp)
# Handle missing species
absent_sp = isnothing.(idxs_in)
any(absent_sp) && throw(ArgumentError("Species $(idxs[absent_sp]) are not \
found in the network. Any mispelling?"))

any(absent_sp) && throw(
ArgumentError("Species $(idxs[absent_sp]) are not found in the network. \
Any mispelling?"),
)
# Get the index of the species names in solution matrix
idxs = idxs_in[.!absent_sp]
# remove Union{nothing, ...}
idxs = something.(idxs)

elseif eltype(idxs) <: Integer
check_bound_idxs = 1 .<= idxs .<= length(sp)
absent_sp = idxs[findall(.!check_bound_idxs)]
all(check_bound_idxs) || throw(ArgumentError("Cannot extract idxs $(absent_sp) \
when there are $(length(sp)) species."))

all(check_bound_idxs) || throw(
ArgumentError(
"Cannot extract idxs $(absent_sp) when there are $(length(sp)) species.",
),
)
else
throw(ArgumentError("`idxs` should be a vector of Integer (species indices) or \
String (species names)"))
throw(ArgumentError("`idxs` should be a vector of integers (species indices) \
or strings (species names)"))
end

idxs
end

Expand All @@ -104,66 +102,67 @@ function process_last_timesteps(solution; last = 1, quiet = false)
n_timesteps = length(solution.t)

if last isa AbstractString
endswith(last, "%") || throw(ArgumentError("The `last` argument, when given as a \
string, should end with character '%'"))

endswith(last, "%") ||
throw(ArgumentError("The `last` argument, when given as a string, \
should end with character '%'"))
perc = parse(Float64, last[1:(end-1)])
is_valid_perc = 0.0 < perc <= 100.0
is_valid_perc || throw(ArgumentError("Cannot extract $(perc)% of the solution's \
timesteps: 0% < `last` <= 100% must hold."))
is_valid_perc ||
throw(ArgumentError("Cannot extract $(perc)% of the solution's timesteps: \
0% < `last` <= 100% must hold."))
last = round(Int, n_timesteps * perc / 100)
last > 0 || quiet || @warn "$perc% of $n_timesteps \
timesteps correspond to $last output lines: an empty table has been extracted."
(last > 0 || quiet) ||
@warn("$perc% of $n_timesteps timesteps correspond to $last output lines: \
an empty table has been extracted.")
elseif last isa Integer
last > 0 || throw(ArgumentError("Cannot extract $last timesteps. `last` should be \
a positive integer."))
last > 0 || throw(ArgumentError("Cannot extract $last timesteps. \
`last` should be a positive integer."))
elseif last isa Float64
throw(ArgumentError("Cannot extract `last` from a floating point number. \
Did you mean \"$last%\"?"))
Did you mean \"$last%\"?"))
else
throw(ArgumentError("Cannot extract timesteps with `last=$last` \
of type $(typeof(last)). \
`last` should be a positive integer \
or a string representing a percentage."))
of type $(typeof(last)). \
`last` should be a positive integer \
or a string representing a percentage."))
end

last > n_timesteps && throw(ArgumentError("Cannot extract $last timesteps from a \
trajectory solution with only \
$(n_timesteps) timesteps. \
Consider decreasing the `last` argument value \
and/or specifying it as a percentage instead \
(e.g. `\"10%\"`)."))
last > n_timesteps && throw(
ArgumentError("Cannot extract $last timesteps from a trajectory solution \
with only $(n_timesteps) timesteps. \
Consider decreasing the `last` argument value \
and/or specifying it as a percentage instead (e.g. `\"10%\"`)."),
)

last
end

function check_last_extinction(solution; idxs = nothing, last = 1)
ext = get_extinction_timesteps(solution; idxs)
ext_t = ext.extinction_timestep

n_timesteps = length(solution.t)

check_last_extinction(n_timesteps; t = ext_t, species = ext.species, last = last)
check_last_extinction(n_timesteps; t = ext_t, species = ext.species, last)

end
function check_last_extinction(n_timesteps::Integer; t, species, last)
extinct = t .!== nothing
if any(extinct)
check_last = findall(>(n_timesteps - (last - 1)), t[extinct])

isempty(check_last) || @warn "With `last` = $last, a table has been extracted with \
the species $(species[extinct][check_last]), \
that went extinct at timesteps = $(t[extinct][check_last]). Set `last` <= \
$(n_timesteps - (maximum(t[extinct]) - 1)) to get rid of them."
sp = species[extinct][check_last]
ts = t[extinct][check_last]
max = n_timesteps - (maximum(t[extinct]) - 1)
isempty(check_last) ||
@warn("With `last` = $last, a table has been extracted with the species $sp, \
that went extinct at timesteps = $ts. \
Set `last` <= $max to get rid of them.")
end
end

function get_extinction_timesteps(solution; idxs = nothing)
idxs = process_idxs(solution; idxs)
sp = get_parameters(solution).network.species[idxs]

ext_t = findfirst.(isequal(0), eachrow(solution[idxs, :]))
extinct = ext_t .!== nothing

(
species = sp[extinct],
idxs = idxs[extinct],
Expand Down Expand Up @@ -205,12 +204,9 @@ julia> sim = simulate(params, [0, 0]; tmax = 20);
function get_alive_species(solution; idxs = nothing, threshold = 0)
idxs = process_idxs(solution; idxs)
sp = get_parameters(solution).network.species[idxs]

alive_idxs = get_alive_species(solution[idxs, end]; threshold = threshold)

(species = sp[alive_idxs], idxs = idxs[alive_idxs])
alive = get_alive_species(solution[idxs, end]; threshold = threshold)
(species = sp[alive], idxs = idxs[alive])
end

get_extinct_species(m::AbstractVector; threshold = 0) = findall(<=(threshold), m)

get_alive_species(m::AbstractVector; threshold = 0) = findall(>(threshold), m)
Loading

0 comments on commit d5b316c

Please sign in to comment.