From 32e1e6de92d9906ab4a5b6194de8900caf1f40f6 Mon Sep 17 00:00:00 2001 From: Iago Bonnici Date: Mon, 22 Apr 2024 18:16:29 +0200 Subject: [PATCH] =?UTF-8?q?=F0=9F=9A=A7=20=E2=8C=9BStart=20designing=20`To?= =?UTF-8?q?pology`=20to=20represent=20whole=20Model.?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/EcologicalNetworksDynamics.jl | 3 + src/Internals/Internals.jl | 1 + src/Internals/model/model_parameters.jl | 1 + src/components/foodweb_graph.jl | 2 + .../nontrophic_components_utils.jl | 5 +- src/methods/graphs.jl | 38 +++ src/topology.jl | 286 ++++++++++++++++++ 7 files changed, 335 insertions(+), 1 deletion(-) create mode 100644 src/topology.jl diff --git a/src/EcologicalNetworksDynamics.jl b/src/EcologicalNetworksDynamics.jl index 60c59c9d8..7c0b97cbe 100644 --- a/src/EcologicalNetworksDynamics.jl +++ b/src/EcologicalNetworksDynamics.jl @@ -40,6 +40,9 @@ using .AliasingDicts include("./multiplex_api.jl") using .MultiplexApi +# Types to represent the model under a pure topological perspective. +include("./topology.jl") # Will be part of the internals after their refactoring. + #------------------------------------------------------------------------------------------- # "Inner" parts: legacy internals. diff --git a/src/Internals/Internals.jl b/src/Internals/Internals.jl index c0a9bb2ca..f174d5d7a 100644 --- a/src/Internals/Internals.jl +++ b/src/Internals/Internals.jl @@ -58,6 +58,7 @@ const Option{T} = Union{Nothing,T} # Since parts of the API is being extracted out of this module to survive, # authorize using it here. using ..EcologicalNetworksDynamics +const Topology = EcologicalNetworksDynamics.Topology # Part of future refactoring here. const equal_fields = EcologicalNetworksDynamics.equal_fields include("./macros.jl") diff --git a/src/Internals/model/model_parameters.jl b/src/Internals/model/model_parameters.jl index 443b42b0a..f5d17c3e6 100644 --- a/src/Internals/model/model_parameters.jl +++ b/src/Internals/model/model_parameters.jl @@ -13,6 +13,7 @@ mutable struct ModelParameters biorates::BioRates # (this one does but all values are initially 'nothing' inside) functional_response::Option{FunctionalResponse} producer_growth::Option{ProducerGrowth} + topology::Topology # This will actually be part of the future refactoring. # Since 'foodweb' is still a mandatory input to construct interaction layers, # keep this artificial reference to it to make the system/components API work # until we refactor all the internals. diff --git a/src/components/foodweb_graph.jl b/src/components/foodweb_graph.jl index 679ecd587..7e53197e2 100644 --- a/src/components/foodweb_graph.jl +++ b/src/components/foodweb_graph.jl @@ -1,5 +1,7 @@ # Wrap a SimpleDiGraph into a custom type indexed by stable labels # to represent either a foodweb or its biomass-foodweb restriction. +# HERE: rather a whole model graph with all nodes/edges layers, +# and make this only a view into it. """ Representation of a restriction of the graph model diff --git a/src/components/nontrophic_layers/nontrophic_components_utils.jl b/src/components/nontrophic_layers/nontrophic_components_utils.jl index ac3f54f15..cf4dd4cae 100644 --- a/src/components/nontrophic_layers/nontrophic_components_utils.jl +++ b/src/components/nontrophic_layers/nontrophic_components_utils.jl @@ -139,10 +139,13 @@ end # ========================================================================================== # Check/expand full layer components. +has_nontrophic_layers(model) = model.network isa MultiplexNetwork +export has_nontrophic_layers + # The application procedure differs # whether the NTI layer is the first to be set or not. function set_layer!(model, interaction, layer) - if model.network isa Internals.FoodWeb + if !has_nontrophic_layers(model) # First NTI component to be added. # Switch from plain foodweb to a multiplex network. S = model.richness diff --git a/src/methods/graphs.jl b/src/methods/graphs.jl index 18ee84df1..df5556165 100644 --- a/src/methods/graphs.jl +++ b/src/methods/graphs.jl @@ -1,3 +1,10 @@ +# Retrieve underlying model topology. +@expose_data graph begin + property(topology, model_graph) + ref(m -> m.topology) + get(m -> deepcopy(m.topology)) +end + """ restrict_to_live!(g::TrophicGraph, biomasses; threshold = 0) @@ -118,3 +125,34 @@ function starving_consumers(graph::TrophicGraph) consumers end export starving_consumers + +# ========================================================================================== +# Expose as checked model methods. + +function biomass_foodweb( + m::InnerParms, + biomasses; + threshold = 0, + warn_nontrophic_links = true, +) + g = m.trophic_graph + if warn_nontrophic_links && has_nontrophic_layers(m) + @warn "Restricting ecological model to only a trophic graph \ + with extinct species nodes removed \ + may yield inconsistent expectations regarding \ + the meaning of remaining non-trophic interactions." + end + restrict_to_live!(g, biomasses; threshold) +end +@method biomass_foodweb depends(Foodweb) + +function disconnected_components(m::InnerParms, biomasses; threshold = 0) + g = biomass(m, biomasses; threshold) + # HERE: this raises too many complications + # just because non-trophic interactions + # are not represented by the trophic graph. + # Make the trophic graph a simple view into the whole model topology, + # reified as a true `ModelGraph` type, + # with proper layers for nodes and edges. + # Drop Graphs.jl in this respect: only locally useful for very particular algorithms. +end diff --git a/src/topology.jl b/src/topology.jl new file mode 100644 index 000000000..baf707fd9 --- /dev/null +++ b/src/topology.jl @@ -0,0 +1,286 @@ +module Topologies + +argerr(mess) = throw(ArgumentError(mess)) + +# Mark removed nodes. +struct Tombstone end + +""" +Values of this type are constructed from a model value, +to represent its pure topology: + + - Nodes identity and types: species, nutrients, patches.. + - Edges types: trophic interaction, non-trophic interactions, migrations.. + +They are supposed not to be mutated, +as they carry faithful topological information reflecting the model +at the moment it has been extracted from it. +Nodes and edge information can be queried from it using labels or indices. +They *may* be removed to represent *e.g.* species extinction, +and study topological consequences or removing them. +But the indices and labels remain stable, +always consistent with their indices from the model value when extracted. +As a consequence, every node (at least) in the topology +can be queried for having been 'removed' or not: tombstones remain. +No tombstone remain for edges: once removed, there is no trace left of them. +""" +struct Topology + # List/index possible types for nodes and edges. + node_types_labels::Vector{Symbol} + node_types_index::Dict{Symbol,Int} + edge_types_labels::Vector{Symbol} + edge_types_index::Dict{Symbol,Int} + # List nodes and their associated types. + # Sorted are *sorted by type*: + # so all nodes of the same type are stored contiguously in this array. + nodes_labels::Vector{Symbol} + nodes_index::Dict{Symbol,Int} + nodes_types::Vector{UnitRange{Int}} # {type -> (start, end) of nodes with this type} + # Topological information: paired, layered adjacency lists. + outgoing::Vector{Union{Tombstone,Vector{Vector{Int}}}} + incoming::Vector{Union{Tombstone,Vector{Vector{Int}}}} + # [node: [edgetype: [nodeid]]] + # ^--------------------------^- : Adjacency list: one entry per node 'N'. + # ^-------------------^- : Tombstone (removed node) or one entry per edge type. + # ^--------^- : One entry per neighbour of 'N': its index. + # Cached redundant information. + n_edges::Vector{Int} # Per edge type. + n_nodes::Vector{Int} # Per node type, not counting tombstones. + Topology() = new([], Dict(), [], Dict(), [], Dict(), [], [], [], []) +end +export Topology + +#------------------------------------------------------------------------------------------- +# Construction primitives. + +# Only push whole slices of nodes of a new type at once. +function push_nodes!(top::Topology, labels, type::Symbol) + haskey(top.node_types_index, type) && + argerr("Node type :$type has already been pushed.") + haskey(top.edge_types_index, type) && + argerr("Node type :$type could be confused with edge type :$type.") + push!(top.node_types_labels, type) + top.node_types_index[type] = length(top.node_types_labels) + n_before = length(top.nodes_labels) + nl = top.nodes_labels + for lab::Symbol in labels + haskey(nl, lab) && argerr("Label :$lab was already given \ + to a node of type :$(node_type(top, lab)).") + push!(nl, lab) + top.nodes_index[lab] = length(nl) + for adj in (top.outgoing, top.incoming) + push!(adj, Vector{Vector{Int}}()) + end + end + n_after = length(nl) + push!(top.nodes_types, n_before+1:n_after) + top +end +export push_nodes! + +function new_edge_type!(top::Topology, type::Symbol) + haskey(top.edge_types_index, type) && argerr("Edge type :$type has already been added.") + haskey(top.node_types_index, type) && + argerr("Edge type :$type would be confused with node type :$type.") + push!(top.edge_types_labels, type) + top.edge_types_index[type] = length(top.edge_types_labels) + for adj in (top.outgoing, top.incoming) + for node in adj + node isa Tombstone && continue + push!(node, []) + end + end + push!(top.n_edges, 0) + top +end +export new_edge_type! + +#------------------------------------------------------------------------------------------- +# Basic unchecked queries. +# (ref-leaking methods are still _protected) + +const imap = Iterators.map +const ifilter = Iterators.filter +idmap(x) = imap(identity, x) # Useful to not leak refs to private collections. + +# Information about types. +n_node_types(top::Topology) = length(top.node_types_labels) +n_edge_types(top::Topology) = length(top.edge_types_labels) +node_type_label(top::Topology, i::Int) = top.node_types_labels[i] +node_type_index(top::Topology, lab::Symbol) = top.node_types_index[lab] +node_type_label(::Topology, lab::Symbol) = lab +node_type_index(::Topology, i::Int) = i +edge_type_label(top::Topology, i::Int) = top.edge_types_labels[i] +edge_type_index(top::Topology, lab::Symbol) = top.edge_types_index[lab] +edge_type_label(::Topology, lab::Symbol) = lab +edge_type_index(::Topology, i::Int) = i +node_types(top::Topology) = idmap(identity, top.node_types_labels) +edge_types(top::Topology) = idmap(identity, top.edge_types_labels) + +# General information about nodes. +n_nodes(top::Topology, type) = top.n_nodes[node_type_index(type)] +n_nodes_and_removed(top::Topology, type) = length(top.nodes_types[node_type_index(type)]) +nodes_indices(top::Topology, type) = top.nodes_types[node_type_index(type)] # Okay to leak: immutable. +_nodes_labels(top::Topology, type) = top.nodes_labels[nodes_indices(top, type)] +node_labels(top::Topology, type) = idmap(identity, _nodes_labels(top, type)) + +# Particular information about nodes. +node_label(top::Topology, i::Int) = top.nodes_labels[i] +node_index(top::Topology, label::Symbol) = top.nodes_index[label] +node_label(::Topology, lab::Symbol) = lab +node_index(::Topology, i::Int) = i +node_type_index(top::Topology, id) = + findfirst(range -> node_index(top, id) in range, top.nodes_types) +node_type(top::Topology, id) = node_type_label(node_type_index(top, id)) + +# Information about edges. +n_edges(top::Topology, type) = top.n_edges[edge_type_index(type)] + +# Direct neighbourhood when querying particular edge type. +# (assuming target is not a tombstone) +function _outgoing_indices(top::Topology, node, edge_type) + i_type = edge_type_index(top, edge_type) + _outgoing_indices(top, node)[i_type] +end +function _incoming_indices(top::Topology, node, edge_type) + i_type = edge_type_index(top, edge_type) + _incoming_indices(top, node)[i_type] +end +outgoing_indices(top::Topology, node, type) = idmap(_outgoing_indices(top, node, type)) +incoming_indices(top::Topology, node, type) = idmap(_incoming_indices(top, node, type)) +outgoing_labels(top::Topology, node, type) = + imap(i -> top.nodes_labels[i], _outgoing_indices(top, node, type)) +incoming_labels(top::Topology, node, type) = + imap(i -> top.nodes_labels[i], _incoming_indices(top, node, type)) + +# Direct neighbourhood: return twolevel iterator: +# first iterate on edge types, then neighbours. +# (assuming target is not a tombstone) +function _outgoing_indices(top::Topology, node) + i_node = node_index(top, node) + top.outgoing[i_node] +end +function _incoming_indices(top::Topology, node) + i_node = node_index(top, node) + top.incoming[i_node] +end +outgoing_indices(top::Topology, node) = imap( + (i_edge_type, _neighbours) -> (i_edge_type, idmap(_neighbours)), + enumerate(_outgoing_indices(top, node)), +) +incoming_indices(top::Topology, node) = imap( + (i_edge_type, _neighbours) -> (i_edge_type, idmap(_neighbours)), + enumerate(_incoming_indices(top, node)), +) +outgoing_labels(top::Topology, node) = imap( + (i_edge, _neighbours) -> ( + top.edge_types_labels[i_edge], + imap(i_node -> top.nodes_labels[i_node], _neighbours), + ), + enumerate(_outgoing_indices(top, node)), +) +incoming_labels(top::Topology, node) = imap( + (i_edge, _neighbours) -> ( + top.edge_types_labels[i_edge], + imap(i_node -> top.nodes_labels[i_node], _neighbours), + ), + enumerate(_incoming_indices(top, node)), +) + +# Filter adjacency iterators given one particular edge type. +# Also return twolevel iterators: focal node, then its neighbours. +function _outgoing_edges_indices(top::Topology, edge_type) + i_type = edge_type_index(top, edge_type) + imap(ifilter(enumerate(top.outgoing)) do (_, node) + !(node isa Tombstone) + end) do (i, _neighbours) + (i, _neighbours[i_type]) + end +end +function _incoming_edges_indices(top::Topology, edge_type) + i_type = edge_type_index(top, edge_type) + imap(ifilter(enumerate(top.incoming)) do (_, node) + !(node isa Tombstone) + end) do (i, _neighbours) + (i, _neighbours[i_type]) + end +end +outgoing_edges_indices(top::Topology, edge_type) = imap( + (i_node, _neighbours) -> (i_node, idmap(_neighbours)), + _outgoing_edges_indices(top, edge_type), +) +incoming_edges_indices(top::Topology, edge_type) = imap( + (i_node, _neighbours) -> (i_node, idmap(_neighbours)), + _incoming_edges_indices(top, edge_type), +) +outgoing_edges_labels(top::Topology, edge_type) = imap( + (i_node, _neighbours) -> + (node_label(top, i_node), idmap(i -> node_label(top, i), _neighbours)), + _outgoing_edges_indices(top, edge_type), +) +incoming_edges_labels(top::Topology, edge_type) = imap( + (i_node, _neighbours) -> + (node_label(top, i_node), idmap(i -> node_label(top, i), _neighbours)), + _incoming_edges_indices(top, edge_type), +) + + +#------------------------------------------------------------------------------------------- +# Display. + +function Base.show(io::IO, top::Topology) + nnt = n_node_types(top) + net = new_edge_type!(top) + s = n -> n > 1 ? "s" : "" + print( + io, + "Topology for $nnt node type$(s(nnt)) \ + and $net edge type$s(net):\n", + ) + tombs = [] # Collect removed nodes to display at the end. + for (i_type, type) in enumerate(node_types(top)) + tomb = [] + push!(tombs, tomb) + print(io, " Nodes '$type':") + empty = true + for node in _nodes_labels(i_type) + if node isa Tombstone + push!(tomb, node) + continue + end + print(io, " :$node") + empty = false + end + if empty + print(io, " ") + end + println(io) + end + for (i_type, type) in enumerate(edge_types(top)) + print(io, " Edges of type '$type':") + empty = true + for (i_source, _neighbours) in _outgoing_edges_indices(top, i_type) + isempty(_neighbours) && continue + source = node_label(top, i_source) + print(io, " $source:") + for target in _neighbours + print(io, " :$target") + end + empty = false + end + if empty + print(io, " ") + end + println(io) + end +end + +end + +# ========================================================================================== +# DEBUG +using .Topologies +top = Topology() +push_nodes!(top, Symbol.(collect("abc")), :species) +new_edge_type!(top, :trophic) +top # HERE: debug all the above at least until displaying works.