diff --git a/.github/workflows/Lint.yml b/.github/workflows/Lint.yml index 8ab8f34..77959dc 100644 --- a/.github/workflows/Lint.yml +++ b/.github/workflows/Lint.yml @@ -27,7 +27,7 @@ jobs: - name: Use Julia cache uses: julia-actions/cache@v2 - name: Install JuliaFormatter.jl - run: julia -e 'using Pkg; pkg"add JuliaFormatter@1"' + run: julia -e 'using Pkg; pkg"add JuliaFormatter@2"' - name: Setup Python uses: actions/setup-python@v6 with: diff --git a/Project.toml b/Project.toml index 0ff143f..c160581 100644 --- a/Project.toml +++ b/Project.toml @@ -1,43 +1,31 @@ name = "GraphDynamicalSystems" uuid = "13529e2e-ed53-56b1-bd6f-420b01fca819" +version = "0.0.7" authors = ["Reuben Gardos Reid <5456207+ReubenJ@users.noreply.github.com>"] -version = "0.0.6" [deps] AbstractTrees = "1520ce14-60c1-5f80-bbc7-55ef81b5835c" -AutoHashEquals = "15f4f7f2-30c1-5605-9d31-71845cf9641f" DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae" DynamicalSystemsBase = "6e36e845-645a-534a-86f2-f5d4aa5a06b4" Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6" -HerbConstraints = "1fa96474-3206-4513-b4fa-23913f296dfc" -HerbCore = "2b23ba43-8213-43cb-b5ea-38c12b45bd45" -HerbGrammar = "4ef9e186-2fe5-4b24-8de7-9f7291f24af7" -HerbSearch = "3008d8e8-f9aa-438a-92ed-26e9c7b4829f" -JSON = "682c06a0-de6a-54ab-a142-c8b1cf79cde6" MLStyle = "d8e11817-5142-5d16-987a-aa16d5891078" MacroTools = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09" MetaGraphsNext = "fa8bd995-216d-47f1-8a91-f3b68fbeb377" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" -StructUtils = "ec057cc2-7a8d-4b58-b3b3-92acb9f63b42" +TermInterface = "8ea1fca8-c5ef-4a55-8b96-4e9afe9c9a3c" [compat] AbstractTrees = "0.4.5" -AutoHashEquals = "2.2.0" DocStringExtensions = "0.9.3" DynamicalSystemsBase = "3.10" Graphs = "1.12" -HerbConstraints = "0.4" -HerbCore = "0.3.4" -HerbGrammar = "0.6" -HerbSearch = "0.4.1" -JSON = "1" MLStyle = "0.4.17" MacroTools = "0.5.16" MetaGraphsNext = "0.7" Random = "1.10" SciMLBase = "2.74.1" StaticArrays = "1.9.12" -StructUtils = "2.5.1" +TermInterface = "2.0.0" julia = "1.10" diff --git a/src/GraphDynamicalSystems.jl b/src/GraphDynamicalSystems.jl index 5262512..6f5394e 100644 --- a/src/GraphDynamicalSystems.jl +++ b/src/GraphDynamicalSystems.jl @@ -3,32 +3,9 @@ module GraphDynamicalSystems using DocStringExtensions include("gds_interface.jl") -export GraphDynamicalSystem, - GDS, - ScheduleStyle, - Asynchronous, - Synchronous, - get_n_entities, - entities, - get_schedule, - get_state, - get_graph include("qualitative_networks.jl") -export QualitativeNetwork, - QN, - Entity, - build_qn_grammar, - update_functions_to_interaction_graph, - sample_qualitative_network, - get_domain, - target_functions, - interpret, - create_qn_system, - default_target_function, - set_state!, - current_parameters -include("io/bma.jl") +# include("io/bma.jl") end diff --git a/src/gds_interface.jl b/src/gds_interface.jl index 8a58656..0472aaa 100644 --- a/src/gds_interface.jl +++ b/src/gds_interface.jl @@ -1,11 +1,27 @@ -import DynamicalSystemsBase.get_state +import DynamicalSystemsBase as DSB using MetaGraphsNext: labels +public ScheduleStyle, + Asynchronous, + Synchronous, + GraphDynamicalSystem, + GDS, + get_n_entities, + get_schedule, + get_graph, + get_domain, + get_fn + abstract type ScheduleStyle end struct Asynchronous <: ScheduleStyle end struct Synchronous <: ScheduleStyle end -abstract type GraphDynamicalSystem{N,S} end +# A GDS is a DiscreteTimeDynamicalSystem where the model is a MetaGraph with +# functions at each of the vertices in the graph. The parameters of the system +# include the ranges of each of the entities. The state of the system is just +# the current value of each of the entities + +abstract type GraphDynamicalSystem{N,S} <: DSB.DiscreteTimeDynamicalSystem end const GDS = GraphDynamicalSystem """ @@ -44,11 +60,34 @@ function entities(gds::GDS) return collect(labels(get_graph(gds))) end +abstract type AbstractEntity end + +function get_fn end + """ $(TYPEDSIGNATURES) Get the state of the GDS. """ -function get_state(gds::GDS) - return gds.state +function DSB.current_state(gds::GDS) + g = get_graph(gds) + l = labels(g) + return DSB.current_state.(getindex.((g,), l)) +end + +function DSB.current_state(gds::GDS, entity) + g = get_graph(gds) + return DSB.current_state(g[entity]) +end + +function get_fn(gds::GDS) + g = get_graph(gds) + l = labels(g) + return get_fn.(getindex.((g,), l)) +end + +function get_domain(gds::GDS) + g = get_graph(gds) + l = labels(g) + return get_domain.(getindex.((g,), l)) end diff --git a/src/io/bma.jl b/src/io/bma.jl deleted file mode 100644 index cc50be0..0000000 --- a/src/io/bma.jl +++ /dev/null @@ -1,393 +0,0 @@ -module BMA - -import GraphDynamicalSystems -import GraphDynamicalSystems: - Asynchronous, - EntityIdName, - QN, - QualitativeNetwork, - ScheduleStyle, - Synchronous, - default_target_function, - domain, - entities, - get_domain, - get_graph, - id, - name, - range_from, - range_to, - target_function, - target_functions -import Graphs: SimpleDiGraph, add_edge!, add_vertex! -import JSON -import MLStyle: @match -import MetaGraphsNext: MetaGraph, edge_labels, inneighbor_labels, labels -import StructUtils - -using DocStringExtensions -using MacroTools: @capture, postwalk - -StructUtils.@tags struct Relationship - id::Int & (json = (name = "id",),) - from::Int & (json = (name = "fromvariable",),) - to::Int & (json = (name = "tovariable",),) - type::String & (json = (name = "type",),) -end - -GraphDynamicalSystems.id(r::Relationship) = r.id -from(r::Relationship) = r.from -to(r::Relationship) = r.to -type(r::Relationship) = r.type - -StructUtils.@defaults struct Entity - target_function::Any & (json = (name = "formula",),) - id::Int & (json = (name = "id",),) - range_from::Int & (json = (name = "rangefrom",),) - range_to::Int & (json = (name = "rangeto",),) - name::String = "" & (json = (name = "name",),) -end -GraphDynamicalSystems.id(e::Entity) = e.id -target_function(e::Entity) = e.target_function -range_from(e::Entity) = e.range_from -range_to(e::Entity) = e.range_to -GraphDynamicalSystems.name(e::Entity) = e.name - -function GraphDynamicalSystems.Entity(e::Entity) - GraphDynamicalSystems.Entity( - (id(e), name(e)), - target_function(e), - range_from(e):range_to(e), - ) -end - -StructUtils.@tags struct Model - entities::Vector{Entity} & (json = (name = "variables",),) - relationships::Vector{Relationship} -end -entities(m::Model) = m.entities -relationships(m::Model) = m.relationships - -struct BMAInputFormat - model::Model - # layout is not needed here - # and thus ignored -end -model(bma::BMAInputFormat) = bma.model - -function bma_dict_to_qn(bma_model::Model) - bma_entities = GraphDynamicalSystems.Entity.(entities(bma_model)) - bma_relationships = relationships(bma_model) - - names = name.(bma_entities) - id_to_name = Dict(id.(bma_entities) .=> names) - mg = MetaGraph(SimpleDiGraph(), Int, Union{Expr,Integer,Symbol}, String) - - foreach(bma_entities) do v - # adding an empty expression: :() - # because we need to construct the interaction graph - # first before parsing the functions correctly - added = add_vertex!(mg, id(v), :()) - if !added - error( - """ - Failed to add the entity (\"$(name(v))\", id: #$(id(v))) from the input file while \ - constructing the QN. Check that there is only one entity in the model with \ - the id #$(id(v)). - """, - ) - end - end - - foreach(bma_relationships) do r - if from(r) ∉ labels(mg) || to(r) ∉ labels(mg) - error("Either the source or destination of the edge is not in the graph.") - end - added = add_edge!(mg, from(r), to(r), type(r)) - if !added - @warn """ - Could not create an edge between entities (from: \ - #$(from(r)); to: #$(to(r))) while constructing \ - the QN. - """ - end - end - - entities_with_functions = [ - GraphDynamicalSystems.Entity( - (id(e), name(e)), - create_target_function( - e, - collect(inneighbor_labels(mg, id(e))), - id_to_name, - mg, - ), - domain(e), - ) for e in bma_entities - ] - - return QualitativeNetwork(entities_with_functions; schedule = Synchronous) -end - -""" - $(SIGNATURES) - -Classify all symbols in `ex` as activators or inhibitors. - -## Examples - - -""" -function classify_activators_inhibitors(ex, sign = 1, activators = [], inhibitors = []) - (activators, inhibitors) = @match ex begin - ::Symbol => if sign == 1 - (push!(activators, ex), inhibitors) - else - (activators, push!(inhibitors, ex)) - end - ::Int => (activators, inhibitors) - Expr(:call, :(-), child) => - classify_activators_inhibitors(child, -sign, activators, inhibitors) - Expr(:call, :(-), left_child, right_child) => begin - (activators, inhibitors) = classify_activators_inhibitors( - left_child, - sign, - activators, - inhibitors, - ) - (activators, inhibitors) = classify_activators_inhibitors( - right_child, - -sign, - activators, - inhibitors, - ) - (activators, inhibitors) - end - Expr(:call, f, children...) => begin - for child in children - (activators, inhibitors) = classify_activators_inhibitors( - child, - sign, - activators, - inhibitors, - ) - end - (activators, inhibitors) - end - Expr(expr_type, _...) => error("Can't classify expression of type $expr_type") - end - - return activators, inhibitors -end - -function classify_activators_inhibitors(d::AbstractDict) - return Dict(e => classify_activators_inhibitors(fn) for (e, fn) in d) -end - -function swap_entity_names_to_var_ids(ex) - @match ex begin - ::Symbol && if (ex ∉ [:+, :-, :/, :*, :min, :max, :ceil, :floor]) - end => :(var($(parse(Int, last(rsplit(string(ex), "_"; limit = 2)))))) - Expr(:call, op, children...) => - Expr(:call, op, swap_entity_names_to_var_ids.(children)...) - _ => ex - end -end - -""" - stringify_fn(ex, lower_bound, upper_bound) - -Take an `ex` and if it's of the form of a default function, return "". -""" -function stringify_fn(ex, lower_bound, upper_bound) - if is_default_function(ex, lower_bound, upper_bound) - return "" - else - string(ex) - end -end - -function is_default_function(ex, lower_bound, upper_bound) - @match ex begin - # no inputs - -1 => true - # single activator - :(var($id)) => true - - # multiple activators - Expr(:call, :/, Expr(:call, :+, vars...), denom) && ( - if length(vars) == denom - end - ) => true - - # only inhibitor(s) - :($bound - $inh) && ( - if bound == upper_bound - end - ) => is_default_function(inh, lower_bound, upper_bound) - - # both inhibitor(s) and activator(s) - :(max($bound, $act - $inh)) && ( - if bound == lower_bound - end - ) => - is_default_function(act, lower_bound, upper_bound) && - is_default_function(inh, lower_bound, upper_bound) - _ => false - end -end - -""" - $(SIGNATURES) - -Write QN to a dictionary to output as JSON. - -Use `JSON.json(qn)` directly to convert to JSON. -""" -function qn_to_bma_dict(qn::QN{N,S,M}) where {N,S,C,G,L<:EntityIdName,M<:MetaGraph{C,G,L}} - lower_upper = extrema.(get_domain(qn)) - ids = id.(entities(qn)) - entity_names = name.(entities(qn)) - functions = [target_functions(qn)[e] for e in entities(qn)] - activator_inhibitor_pairs = classify_activators_inhibitors(target_functions(qn)) - functions = swap_entity_names_to_var_ids.(functions) - functions = stringify_fn.(functions, first.(lower_upper), last.(lower_upper)) - - variables = [ - Dict( - "RangeFrom" => d[1], - "RangeTo" => d[2], - "Id" => i, - "Formula" => f, - "Name" => n, - ) for (d, i, n, f) in zip(lower_upper, ids, entity_names, functions) - ] - relationships = [ - Dict( - "Id" => i, - "FromVariable" => id(src), - "ToVariable" => id(dst), - "Type" => let (activators, inhibitors) = activator_inhibitor_pairs[dst] - activators_transformed = EntityIdName.(activators) - inhibitors_transformed = EntityIdName.(inhibitors) - if src in activators_transformed - "Activator" - elseif src in inhibitors_transformed - "Inhibitor" - else - error( - "Malformed edge. $src not found in activators ($activators_transformed) or inhibitors ($inhibitors_transformed).", - ) - end - end, - ) for (i, (src, dst)) in enumerate(edge_labels(get_graph(qn))) - ] - output_dict = Dict( - "Model" => Dict("Variables" => variables, "Relationships" => relationships), - "Layout" => Dict( - "Variables" => - [Dict("Id" => v["Id"], "Name" => v["Name"]) for v in variables], - ), - ) - - return output_dict -end - -function sanitize_formula(f) - # surround variable names with quotes - return replace(f, r"var\(([^\)]+)\)" => s"var(\"\1\")") -end - -function entity_name_from_in_neighbors(entity, in_neighbors) - # the formulas can reference their incoming edges - # with either the name of the neighbor entity or - # its id - e_id = tryparse(Int, entity) - - entity_name = [ - Symbol("$(name)_$id") for - (id, name, _) in in_neighbors if isnothing(e_id) ? name == entity : id == e_id - ] - - if length(entity_name) != 1 - error( - """ - Error while constructing name for entity: $entity, with in neighbors: \ - $in_neighbors. There are more than one incoming neighbor entities with the same \ - name. To fix this error, remove the erroneous relationships from the JSON file, \ - or reference the entity by id (like `var(3)`). - """, - ) - end - return only(entity_name) -end - -function create_target_function( - variable::GraphDynamicalSystems.Entity{EntityIdName{S},Int}, - in_neighbor_ids::Vector{Int}, - id_to_name::Dict, - mg::MetaGraph, -) where {S} - formula = Meta.parse(sanitize_formula(target_function(variable))) - in_neighbor_names = getindex.((id_to_name,), in_neighbor_ids) - in_neighbor_types = getindex.((mg.edge_data,), in_neighbor_ids, (id(variable),)) - in_neighbors = zip(in_neighbor_ids, in_neighbor_names, in_neighbor_types) - - if isnothing(formula) # default target function - if length(in_neighbor_ids) == 0 - @warn "$(name(variable)) has no inputs, defaulting formula to -1" - return -1 - else - activators = [ - Symbol("$(name)_$id") for - (id, name, ty) in in_neighbors if ty == "Activator" - ] - inhibitors = [ - Symbol("$(name)_$id") for - (id, name, ty) in in_neighbors if ty == "Inhibitor" - ] - return default_target_function( - range_from(variable), - range_to(variable), - activators, - inhibitors, - ) - end - else # custom target function - return postwalk( - x -> - @capture(x, var(v_String)) ? - :($(entity_name_from_in_neighbors(v, in_neighbors))) : x, - formula, - ) - end -end - -function nested_dicts_keys_to_lowercase(d) - if d isa AbstractDict - return Dict([lowercase(k) => nested_dicts_keys_to_lowercase(v) for (k, v) in d]) - elseif d isa AbstractVector - return [nested_dicts_keys_to_lowercase(v) for v in d] - else - return d - end -end - -function QualitativeNetwork(bma_file_path::AbstractString) - bma_def_raw = JSON.parsefile(bma_file_path) - bma_def_raw = nested_dicts_keys_to_lowercase(bma_def_raw) - bma_def = JSON.parse(JSON.json(bma_def_raw), BMAInputFormat) - bma_model = model(bma_def) - - return bma_dict_to_qn(bma_model) -end - -function JSON.json(qn::QualitativeNetwork) - return JSON.json(qn_to_bma_dict(qn); omit_empty = false) -end -JSON.json(io_or_filename, qn::T) where {T<:QualitativeNetwork} = - JSON.json(io_or_filename, qn_to_bma_dict(qn); omit_empty = false) -JSON.json(io::IO, qn::T) where {T<:QualitativeNetwork} = - JSON.json(io, qn_to_bma_dict(qn); omit_empty = false) - -end diff --git a/src/qualitative_networks.jl b/src/qualitative_networks.jl index c97d0cf..85aebc8 100644 --- a/src/qualitative_networks.jl +++ b/src/qualitative_networks.jl @@ -1,612 +1,170 @@ -import AutoHashEquals: @auto_hash_equals -import DynamicalSystemsBase: current_parameters, get_state, set_state! -import JSON +import MetaGraphsNext as MG +import Graphs +import AbstractTrees as AT +import TermInterface as TI +import DynamicalSystemsBase as DSB +import MLStyle import SciMLBase -using AbstractTrees: Leaves, PostOrderDFS -using DynamicalSystemsBase: ArbitrarySteppable, initial_state -using Graphs: AbstractGraph, SimpleDiGraph, add_edge!, add_vertex!, ne -using HerbConstraints: DomainRuleNode, Forbidden, Ordered, Unique, VarNode, addconstraint! -using HerbCore: AbstractGrammar, RuleNode, get_rule -using HerbGrammar: @csgrammar, add_rule!, rulenode2expr -using HerbSearch: rand -using MLStyle: @match -using MetaGraphsNext: MetaGraph, add_edge!, edge_labels, inneighbor_labels, labels, nv -using StaticArrays: MVector, SVector - -const base_qn_grammar = @csgrammar begin - Val = Val + Val - Val = Val - Val - Val = Val / Val - Val = Val * Val - Val = min(Val, Val) - Val = max(Val, Val) - Val = ceil(Val) - Val = floor(Val) -end - -const default_qn_constants = [0, 1, 2] - -""" - $(TYPEDSIGNATURES) - -Builds a grammar based on the base QN grammar adding `entity_names` and `constants` -to the grammar. - -The following constraints are currently included - -1. removing symmetry due to commutativity of `+`/`*`/`min`/`max` -2. forbidding same arguments of two argument functions -3. forbidding constant arguments to 2-argument functions -4. forbidding constant arguments to 1-argument functions -5. using each of the entities only once per function -6. forbidding adding or subtracting zero -7. forbidding multiplication and division by 1 or 0 -8. forcing the first operator inside `ceil` and `floor` to be `÷` -9. forbidding `max(□, X)` and `min(□, X)` where X is either the max or min -constant in the grammar. - -""" -function build_qn_grammar( - entity_names, - constants = default_qn_constants; - unique_constr = true, -) - g = deepcopy(base_qn_grammar) - - for e in entity_names - add_rule!(g, :(Val = $e)) - end - - for c in constants - add_rule!(g, :(Val = $c)) - end - - add_rule!(g, :(Start = Val)) - - # +, *, min, max, are all commutative - domain = BitVector(zeros(length(g.rules))) - @. domain[[1, 4:6...]] = true - template_tree = DomainRuleNode(domain, [VarNode(:a), VarNode(:b)]) - order = [:a, :b] - - addconstraint!(g, Ordered(deepcopy(template_tree), order)) - - # Forbid same arguments for 2-argument functions - domain = BitVector(zeros(length(g.rules))) - @. domain[length(g.childtypes)==2] = true - template_tree = DomainRuleNode(domain, [VarNode(:a), VarNode(:a)]) - - addconstraint!(g, Forbidden(deepcopy(template_tree))) - - # Forbid constant arguments for 2-argument functions - domain = falses(length(g.rules)) - @. domain[length(g.childtypes)==2] = true - consts_domain = falses(length(g.rules)) - consts_domain[findall(x -> x isa Int, g.rules)] .= true - consts_domain_rn = DomainRuleNode(consts_domain) - template_tree = DomainRuleNode(domain, [consts_domain_rn, consts_domain_rn]) - - addconstraint!(g, Forbidden(deepcopy(template_tree))) - - # Forbid constant arguments for 1-argument functions - domain = falses(length(g.rules)) - @. domain[[7, 8]] = true - consts_domain = falses(length(g.rules)) - consts_domain[findall(x -> x isa Int, g.rules)] .= true - consts_domain_rn = DomainRuleNode(consts_domain) - template_tree = DomainRuleNode(domain, [consts_domain_rn]) - - addconstraint!(g, Forbidden(deepcopy(template_tree))) - - n_original_rules = length(base_qn_grammar.rules) - - # Only use each of the entities once per function - n_consts = length(constants) - entities = (n_original_rules+1):(length(g.rules)-n_consts) - - if unique_constr - addconstraint!.((g,), Unique.(entities)) - end - - # Forbid □ + 0, □ - 0 - plus_or_minus = falses(length(g.rules)) - plus_or_minus[[1, 2]] .= true - zero_rule = findfirst(==(0), g.rules) - if !isnothing(zero_rule) - template_tree = DomainRuleNode(plus_or_minus, [VarNode(:a), RuleNode(zero_rule)]) - - addconstraint!(g, Forbidden(deepcopy(template_tree))) - - # Both orderings, but only for plus. Allow 0 - □ - plus_or_minus[2] = false - template_tree = DomainRuleNode(plus_or_minus, [RuleNode(zero_rule), VarNode(:a)]) - addconstraint!(g, Forbidden(deepcopy(template_tree))) - end - - # Forbid □ * 1, □ / 1, □ * 0, □ / 0 - mult_or_div = falses(length(g.rules)) - mult_or_div[[3, 4]] .= true - one_zero_domain = falses(length(g.rules)) - one_zero_domain[findfirst(==(1), g.rules)] = true - if !isnothing(findfirst(==(0), g.rules)) - one_zero_domain[findfirst(==(0), g.rules)] = true - end - - template_tree = - DomainRuleNode(mult_or_div, [VarNode(:a), DomainRuleNode(one_zero_domain)]) - - addconstraint!(g, Forbidden(deepcopy(template_tree))) - - # Forbid ceil(X) and floor(X) unless X = □ ÷ □ - ceil_or_floor = BitVector(zeros(length(g.rules))) - ceil_or_floor[[7, 8]] .= true - all_except_div = trues(length(g.rules)) - all_except_div[3] = false - template_tree = DomainRuleNode(ceil_or_floor, [DomainRuleNode(all_except_div)]) - - addconstraint!(g, Forbidden(deepcopy(template_tree))) - - # Forbid max(□, X) and min(□, X) where X is either the largest or smallest constant in the grammar - min_max_rules = falses(length(g.rules)) - min_max_rules[[5, 6]] .= true - (min_const, max_const) = extrema(filter(x -> isa(x, Int), g.rules)) - extrema_domain = falses(length(g.rules)) - extrema_domain[findall(x -> x == min_const || x == max_const, g.rules)] .= true - rule_extrema_consts = DomainRuleNode(extrema_domain) - template_tree = DomainRuleNode(min_max_rules, [VarNode(:a), rule_extrema_consts]) - - addconstraint!(g, Forbidden(deepcopy(template_tree))) - - return g -end - -""" - $TYPEDSIGNATURES - -Construct a default target function for an entity in a QN from a list of -`activators` and `inhibitors`. - -Follows the definition given in Eq. 3 of ["Qualitative networks: a symbolic -approach to analyze biological signaling -networks"](https://doi.org/10.1186/1752-0509-1-4). - -## Examples - -Say we have a component `X` and it has an lower bound on its state value of 0, -an upper bound of 4, activators `A`, `B`, `C`, and inhibitors `D`, `E`, `F`, -then the following example constructs an expression for its default target -function. - -```jldoctest -julia> default_target_function(0, 4, [:A, :B, :C], [:D, :E, :F]) -:(max(0, (A + B + C) / 3 - (D + E + F) / 3)) -``` - -""" -function default_target_function( - lower_bound::Integer, - upper_bound::Integer, - activators::AbstractVector = [], - inhibitors::AbstractVector = [], -) - sum_only_or_nothing = x -> if length(x) == 0 - nothing - elseif length(x) == 1 - :($(only(x))) - elseif length(x) > 1 - :($(Expr(:call, :+, x...)) / $(length(x))) - end - - expr_activators = sum_only_or_nothing(activators) - expr_inhibitors = sum_only_or_nothing(inhibitors) - - if isnothing(expr_activators) && isnothing(expr_inhibitors) - error("Constructing a default target function for a QN with no \ - activators or inhibitors.") - elseif isnothing(expr_activators) # no activators, special case mentioned in paper - return :($upper_bound - $expr_inhibitors) - elseif isnothing(expr_inhibitors) - return :($expr_activators) - else - return :(max($lower_bound, $expr_activators - $expr_inhibitors)) - end -end - -abstract type EntityLabel end - -struct EntityId <: EntityLabel - id::Int -end -id(e::EntityId) = e.id - -struct EntityName{S} <: EntityLabel - name::S -end -name(e::EntityName) = e.name - -@auto_hash_equals struct EntityIdName{S} <: EntityLabel - id::Int - name::S -end -function EntityIdName(s::Symbol) - en_str = string(s) - name_id_str_split = rsplit(en_str, "_"; limit = 2) - if length(name_id_str_split) != 2 - error("""Failed to convert the Symbol $s to an EntityIdName. \ - Expecting an EntityName with a name in the form of "Name_00".""") - end - (name_str, id_str) = name_id_str_split - - id_val = tryparse(Int, id_str) - if isnothing(id_val) - error("""Entity name ($s) contained an underscore but the \ - content after the underscore ($id_str) could not be parsed as \ - an integer to convert it to an ID.""") - end - - return EntityIdName(id_val, string(name_str)) -end -EntityIdName{String}(s::Symbol) = EntityIdName(s) -id(e::EntityIdName) = e.id -name(e::EntityIdName) = e.name -combined_name(e::EntityIdName) = Symbol("$(name(e))_$(id(e))") - -struct Entity{I<:EntityLabel,D} - label::I - target_function::Any - domain::UnitRange{D} -end -Entity(name::Symbol, args...) = Entity(EntityName(name), args...) -Entity(id::Int, args...) = Entity(EntityId(id), args...) -Entity((id, name), args...) = Entity(EntityIdName(id, name), args...) - -label(e::Entity) = e.label -id(e::Entity) = id(label(e)) -name(e::Entity) = name(label(e)) -target_function(e::Entity) = e.target_function -domain(e::Entity) = e.domain -range_from(e::Entity) = first(domain(e)) -range_to(e::Entity) = last(domain(e)) - -function get_used_entities(fn, entities_in_model::Vector{<:Entity{<:EntityIdName}}) - filter(in(combined_name.(label.(entities_in_model))), collect(Leaves(fn))) -end - -function get_used_entities(fn, entities_in_model) - filter(in(name.(entities_in_model)), collect(Leaves(fn))) -end - -""" - $(TYPEDSIGNATURES) -""" -function update_functions_to_interaction_graph( - entities_in_model::AbstractVector{<:E}, - schedule = Synchronous, -) where {EntityLabelType,E<:Entity{EntityLabelType}} - graph = MetaGraph( - SimpleDiGraph(); - label_type = EntityLabelType, - vertex_data_type = E, - graph_data = schedule, - ) - if !allunique(label.(entities_in_model)) - val_counts = Dict() - for e in entities_in_model - val_counts[name(e)] = append!(get(val_counts, name(e), []), [e]) - end - duplicates = [v for v in values(val_counts) if length(v) > 1] - error("""The QN implementation only supports models with unique \ - entity name/id combinations. - - Duplicates: - - $duplicates""") - end - - for entity in entities_in_model - graph[label(entity)] = entity - end - - for dst in entities_in_model - input_entities = get_used_entities(target_function(dst), entities_in_model) - for src in EntityLabelType.(input_entities) - dst_label = label(dst) - l = collect(labels(graph)) - if !(src ∈ l && dst_label ∈ l) - error( - """Could not add edge from $src to $(dst_label). The vertex labels in the graph are currently $(collect(labels(graph))).""", - ) - end - add_edge!(graph, src, dst_label) - end - end - - return graph -end - -""" - $(TYPEDSIGNATURES) -""" -function sample_qualitative_network( - entities::AbstractVector{Symbol}, - domains::AbstractVector{UnitRange{Int}}, - max_eq_depth::Int; - schedule = Synchronous, -) - g = build_qn_grammar(entities, default_qn_constants) - update_fns = Union{Expr,Integer,Symbol}[ - rulenode2expr(rand(RuleNode, g, :Val, max_eq_depth), g) for _ in entities - ] - - qn = QualitativeNetwork(Entity.(entities, update_fns, domains); schedule = schedule) - - return qn -end - -sample_qualitative_network(N::Int, args...; kwargs...) = - sample_qualitative_network(Symbol.(('A':'Z')[1:N]), args...; kwargs...) +public QualitativeNetwork, QN, interpret """ $(TYPEDEF) -A qualitative network model as described in ["Qualitative networks: a symbolic approach to -analyze biological signaling networks"](https://doi.org/10.1186/1752-0509-1-4). - -This implementation encompasses both the synchronous and asynchonous cases. In the paper, it -is assumed that the synchronous case is used. As such, the default constructor uses a -synchronous schedule. - -$(FIELDS) - -Systems that include the model semantics wrap around this struct with an -[`ArbitrarySteppable`](https://juliadynamics.github.io/DynamicalSystems.jl/stable/tutorial/#DynamicalSystemsBase.ArbitrarySteppable) -from [`DynamicalSystems`](https://juliadynamics.github.io/DynamicalSystems.jl/stable/). See -[`create_qn_system`](@ref) for an example. +A graph dynamical system with a finite domain. State values of each entity are +limited to change by at most 1 per time step. """ -struct QualitativeNetwork{ - N, - Schedule, - M<:MetaGraph{Int,<:SimpleDiGraph,<:EntityLabel,<:Entity}, -} <: GraphDynamicalSystem{N,Schedule} - "Graph containing the topology and target functions of the network" - graph::M - "State of the network" - state::MVector{N,Int} - - function QualitativeNetwork(graph, state; schedule = Synchronous) - N = nv(graph) - if N != length(state) - error("""The number of entities in the model ($N) must match the \ - length of the provided state vector ($(length(state))).""") - end - - return new{N,schedule(),typeof(graph)}(graph, state) - end +struct QualitativeNetwork{N_Entities,Schedule,Graph<:MG.MetaGraph} <: + GraphDynamicalSystem{N_Entities,Schedule} + graph::Graph end -function QualitativeNetwork( - entities::AbstractVector{<:Entity}; - state = nothing, - schedule = Synchronous, -) - graph = update_functions_to_interaction_graph(entities, schedule) - - if isnothing(state) - state = rand.(domain.(entities)) - end - - return QualitativeNetwork(graph, state; schedule) -end - -QualitativeNetwork(entities::AbstractVector{<:AbstractString}, args...; kwargs...) = - QualitativeNetwork(Symbol.(entities), args...; kwargs...) - """ - $(TYPEDSIGNATURES) + $(TYPEDEF) -Shorthand for [`QualitativeNetwork`](@ref). +Alias for a [`QualitativeNetwork`](@ref). """ const QN = QualitativeNetwork -""" - $(TYPEDSIGNATURES) +DSB.isinplace(::QualitativeNetwork) = true +DSB.dynamic_rule(qn::QualitativeNetwork) = get_fn(qn) +DSB.current_parameters(qn::QualitativeNetwork) = () +DSB.current_time(::QualitativeNetwork) = 0 +DSB.reinit!(qn::QN, state::AbstractVector) = DSB.set_state!(qn, Int.(state)) -Get the domain of the entity `entity_label` in `qn`. -""" -function get_domain(qn::QN, entity_label) - graph = get_graph(qn) - entity = graph[entity_label] - return domain(entity) -end - -""" - $(TYPEDSIGNATURES) - -Get all of the domains of the entities in `qn`. -""" -function get_domain(qn::QN) - return get_domain.((qn,), labels(get_graph(qn))) -end - -function _get_entity_index(qn::QN, entity) - # Ugly, we shouldn't pass entities around as Symbols anymore - # but this works for now - # actually, doesn't because there are sometimes variables with underscores... time to rewrite more - entity_proc = if entity isa EntityName - split_res = rsplit(String(name(entity)), "_"; limit = 2) - entity_proc = if length(split_res) == 2 - (entity_str, id_str) = split_res - id_val = tryparse(Int, id_str) - if !isnothing(id_val) - EntityIdName(id_val, entity_str) - else - entity - end - else - entity - end - else - entity - end - i = findfirst(isequal(entity_proc), entities(qn)) - if isnothing(i) - error("""Tried to get the state of $entity_proc but could not retrieve it. \ - The entities in the model are $(entities(qn))""") - end - return i +mutable struct QNEntity{F,S,D} <: AbstractEntity + fn::F + state::S + domain::D end -""" - $(TYPEDSIGNATURES) -""" -function target_functions(qn::QN) - return Dict([ - c => target_function(entity) for (c, (_, entity)) in get_graph(qn).vertex_properties - ]) -end - -""" - $(TYPEDSIGNATURES) -""" -function get_state(qn::QN, entity) - i = _get_entity_index(qn, entity) - return qn.state[i] -end -get_state(qn::QN, entity::Symbol) = get_state(qn, EntityName(entity)) - -function _set_state!(qn::QN, entity, value::Integer) - i = _get_entity_index(qn::QN, entity) - qn.state[i] = value -end - -""" - $(TYPEDSIGNATURES) -""" -function set_state!(qn::QN, entity, value::Integer) - max_for_entity = maximum(get_domain(qn, entity)) - if value > max_for_entity +get_fn(qne::QNEntity) = qne.fn +DSB.current_state(qne::QNEntity) = qne.state +function DSB.set_state!(qne::QNEntity, s) + domain = get_domain(qne) + if s < minimum(domain) || s > maximum(domain) error( - "Value ($value) cannot be larger than the maximum level for $entity ($(max_for_entity))", + "New state value for entity must be within its domain (domain: $domain, new state: $s)", ) end - - _set_state!(qn, entity, value) + qne.state = s end +get_domain(qne::QNEntity) = qne.domain -function set_state!(qn::QN, entity::Symbol, value::Integer) - set_state!(qn, EntityName(entity), value) +function QualitativeNetwork(graph::G) where {G<:Graphs.AbstractGraph} + QualitativeNetwork{Graphs.nv(graph),Synchronous(),G}(graph) end -function set_state!(qn::QN, values) - set_state!.((qn,), entities(qn), values) +function QualitativeNetwork(update_functions::AbstractDict, domains) + QualitativeNetwork{Graphs.SimpleDiGraph}(update_functions, domains) end """ $(TYPEDSIGNATURES) -Interpret target functions from a [`QualitativeNetwork`](@ref). -""" -function interpret(e::Union{Expr,Symbol,Int}, qn::QN) - @match e begin - ::Symbol => get_state(qn, e) - ::Int => e - :($v1 + $v2) => interpret(v1, qn) + interpret(v2, qn) - :($v1 - $v2) => interpret(v1, qn) - interpret(v2, qn) - :($v1 / $v2) => interpret(v1, qn) / interpret(v2, qn) - :($v1 * $v2) => interpret(v1, qn) * interpret(v2, qn) - :(min($v1, $v2)) => min(interpret(v1, qn), interpret(v2, qn)) - :(max($v1, $v2)) => max(interpret(v1, qn), interpret(v2, qn)) - :(ceil($v)) => ceil(interpret(v, qn)) - :(floor($v)) => floor(interpret(v, qn)) - _ => error("Unhandled Expr in `interpret`: $e") - end +Create a [`QualitativeNetwork`](@ref) from the dictionary `update_functions` which should map from entities (`E`) to their functions (`F`) + +The entity names (`E` in the signature) can be anything, while the functions (`F`) are required to either + +- be a numerical constant (like `1`), or a reference to a single entity (like `A`) +- implement the `TermInterface.jl` interface. Any terminal nodes in the functions must be numerical constants or reference an entity. +""" +function QualitativeNetwork{GraphType}( + update_functions::AbstractDict{E,F}, + domains, +)::QualitativeNetwork where {E,F,GraphType<:Graphs.AbstractGraph} + entity_keys = collect(keys(update_functions)) + entity_fns = getindex.((update_functions,), entity_keys) + entity_domains = getindex.((domains,), entity_keys) + get_arguments_or_empty = x -> TI.isexpr(x) ? (x, TI.arguments(x)) : (x, ()) + collect_arguments = + x -> + AT.treemap(get_arguments_or_empty, x) |> + AT.Leaves .|> + AT.nodevalue |> + filter(in(entity_keys)) + + referenced_entities = union.(collect_arguments.(entity_fns)) + referenced_indices = + map(ref_for_e -> findfirst.(.==(ref_for_e), (entity_keys,)), referenced_entities) + edges = + Iterators.flatten( + map(((j, idxs),) -> tuple.(idxs, (j,)), enumerate(referenced_indices)), + ) |> collect + graph = GraphType() + Graphs.add_vertices!(graph, length(entity_keys)) + Graphs.add_edge!.((graph,), Graphs.Edge.(edges)) + vertices_description = Pair{Symbol,QNEntity}[ + (e => QNEntity(fn, 0, d)) for + (e, fn, d) in zip(entity_keys, entity_fns, entity_domains) + ] + edges_description = Pair{Tuple{E,E},Nothing}[ + (entity_keys[s], entity_keys[d]) => nothing for (s, d) in edges + ] + return QualitativeNetwork( + MG.MetaGraph(graph, vertices_description, edges_description, nothing), + ) end -""" - $(TYPEDSIGNATURES) - -Returns the limited value of `next_value` which is at most 1 different than `prev_value`. +function SciMLBase.step!(qn::QualitativeNetwork{N,S}) where {N,S} + SciMLBase.step!(S, qn) +end +SciMLBase.step!(qn::QN, n::Int, _...) = foreach(_ -> SciMLBase.step!(qn), 1:n) -It is also never negative, or larger than `N`. -""" -function limit_change( - prev_value::Integer, - next_value::Integer, - min_level::Integer, - max_level::Integer, -) - if next_value > prev_value - limited_value = min(prev_value + 1, max_level) - elseif next_value < prev_value - limited_value = max(prev_value - 1, min_level) +function limit_change(next, prev, lower, upper) + if next > prev + min(upper, prev + 1) + elseif next < prev + max(lower, prev - 1) else - limited_value = next_value + next end - - return limited_value end -function _compute_next_state!(qn::QN, entity) - (min_level, max_level) = extrema(get_domain(qn, entity)) - t = target_functions(qn)[entity] - old_state = get_state(qn, entity) - new_state = interpret(t, qn) - new_state = isnan(new_state) ? min_level : new_state - new_state = isinf(new_state) ? max_level : new_state - limited_state = limit_change(old_state, floor(Int, new_state), min_level, max_level) +function DSB.set_state!( + qn::QN{N,S,M}, + new_state::Int, + entity::L, +) where {N,S,I,G,L,M<:MG.MetaGraph{I,G,L}} + g = get_graph(qn) + DSB.set_state!(g[entity], new_state) end -""" - $(TYPEDSIGNATURES) -""" -function async_qn_step!(qn::QN) - entity_labels = entities(qn) - entity = rand(entity_labels) - next_state = _compute_next_state!(qn, entity) - set_state!(qn, entity, next_state) -end +function DSB.set_state!(qn::QN, new_state::AbstractVector) + g = get_graph(qn) -""" - $(TYPEDSIGNATURES) -""" -function sync_qn_step!(qn::QN) - next_states = _compute_next_state!.((qn,), entities(qn)) - set_state!.((qn,), entities(qn), next_states) + DSB.set_state!.((qn,), new_state, MG.labels(g)) end -extract_state(model::QN) = model.state -extract_parameters(model::QN) = model.graph -current_parameters(model::QN) = model.graph -reset_model!(model::QN, u, _) = model.state .= u +function SciMLBase.step!(::Synchronous, qn::QualitativeNetwork) + l = MG.labels(get_graph(qn)) + current_state = DSB.current_state(qn) + current_state_dict = Dict(l .=> current_state) + domains = get_domain(qn) + lower_bounds = minimum.(domains) + upper_bounds = maximum.(domains) + fns = get_fn(qn) + next_state_uncapped = interpret.(fns, (current_state_dict,)) + next_state_capped = + limit_change.(next_state_uncapped, current_state, lower_bounds, upper_bounds) + DSB.set_state!(qn, next_state_capped) -function SciMLBase.reinit!( - ds::ArbitrarySteppable{<:AbstractVector{<:Real},<:QualitativeNetwork}, - u::AbstractVector{<:Real} = initial_state(ds); - p = current_parameters(ds), - t0 = 0, # t0 is not used but required for downstream. -) - ds.reinit(ds.model, u, p) - ds.t[] = 0 - return ds + return qn end -""" - $(TYPEDSIGNATURES) - -Construct an asynchronous [`QualitativeNetwork`](@ref) system using the -[`async_qn_step!`](@ref) as a step function. -""" -function create_qn_system(qn::QN) - step_fn = get_schedule(qn) == Asynchronous() ? async_qn_step! : sync_qn_step! - - return ArbitrarySteppable( - qn, - step_fn, - extract_state, - extract_parameters, - reset_model!, - isdeterministic = get_schedule(qn) == Synchronous(), - ) +function interpret(fn, state) + MLStyle.@match fn begin + ::Int => fn + ::Symbol => state[fn] + :($a + $b) => interpret(a, state) + interpret(b, state) + :($a - $b) => interpret(a, state) - interpret(b, state) + :($a / $b) => interpret(a, state) / interpret(b, state) + :(min($a, $b)) => min(interpret(a, state), interpret(b, state)) + :(max($a, $b)) => max(interpret(a, state), interpret(b, state)) + :(ceil($a)) => ceil(interpret(a, state)) + :(floor($a)) => floor(interpret(a, state)) + _ => error("Unhandled expression: $fn") + end end +interpret(fn, qn::QN) = interpret(fn, Dict(entities(qn) .=> DSB.current_state(qn))) diff --git a/test/Project.toml b/test/Project.toml index 38d34ab..95a6b3d 100644 --- a/test/Project.toml +++ b/test/Project.toml @@ -3,12 +3,10 @@ Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595" Attractors = "f3fd9213-ca85-4dba-9dfd-7fc91308fec7" DynamicalSystemsBase = "6e36e845-645a-534a-86f2-f5d4aa5a06b4" Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6" -HerbCore = "2b23ba43-8213-43cb-b5ea-38c12b45bd45" IterTools = "c8e1da08-722c-5040-9ed9-7db0dc04731e" JET = "c3a54625-cd67-489e-a8e7-0a5a0ff4e31b" -JSON = "682c06a0-de6a-54ab-a142-c8b1cf79cde6" MetaGraphsNext = "fa8bd995-216d-47f1-8a91-f3b68fbeb377" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" ReTestItems = "817f1d60-ba6b-4fd5-9520-3cf149f6a823" -SoleLogics = "b002da8f-3cb3-4d91-bbe3-2953433912b5" +SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/test/io/bma_test.jl b/test/io/bma_test.jl deleted file mode 100644 index 326b2c2..0000000 --- a/test/io/bma_test.jl +++ /dev/null @@ -1,147 +0,0 @@ -@testitem "Construct default target functions" begin - lower_bound = 0 - upper_bound = 4 - activators = [:A, :B, :C] - inhibitors = [:D, :E, :F] - - @test default_target_function(lower_bound, upper_bound, activators, inhibitors) == - :(max($lower_bound, (A + B + C) / 3 - (D + E + F) / 3)) - - activators = [:A, :B] - inhibitors = [:D] - - @test default_target_function(lower_bound, upper_bound, activators, inhibitors) == - :(max($lower_bound, (A + B) / 2 - D)) - - activators = [] - inhibitors = [:D] - - @test default_target_function(lower_bound, upper_bound, activators, inhibitors) == - :($upper_bound - D) - - activators = [:A] - inhibitors = [] - - @test default_target_function(lower_bound, upper_bound, activators, inhibitors) == :(A) - - @test_throws r"no activators or inhibitors" default_target_function(0, 4) -end - -@testitem "Load from BMA" begin - using DynamicalSystemsBase: step! - bma_models_path = joinpath(@__DIR__, "..", "resources", "bma_models") - good_models = joinpath(bma_models_path, "well_formed_examples") - - for model_path in readdir(good_models; join = true) - if occursin("Skin1D", model_path) - @test_broken QN(model_path) isa GraphDynamicalSystem - else - qn = QN(model_path) - @test qn isa GraphDynamicalSystem - step!(create_qn_system(qn), 100) - end - end - - bad_models = joinpath(bma_models_path, "error_examples") - - @test_throws "Failed to add" QN(joinpath(bad_models, "duplicate_entity_ids.json")) - @test_throws "Error while constructing name for entity" QN( - joinpath(bad_models, "multiple_incoming_edges_same_name.json"), - ) -end - -@testitem "Save to BMA" begin - import MetaGraphsNext: edge_labels, labels - import GraphDynamicalSystems.BMA: is_default_function - using JSON - - function test_json_roundtrip(model_path::AbstractString) - qn = QN(model_path) - @test length(labels(get_graph(qn))) > 0 - @test length(edge_labels(get_graph(qn))) > 0 - - output_str = JSON.json(qn) - output_dict = JSON.parse(output_str) - orig_dict = JSON.parse(read(model_path, String)) - - if !haskey(orig_dict, "Model") - @warn "Skipping test for file $model_path for now because of the nonstandard key names" - return - end - - @test haskey(output_dict, "Model") - - model_dict = output_dict["Model"] - - @test haskey(model_dict, "Variables") - variables = model_dict["Variables"] - for (orig_v, v) in zip(orig_dict["Model"]["Variables"], variables) - if v["Name"] == "" - @test !haskey(orig_v, "Name") - else - @test v["Name"] == orig_v["Name"] - end - orig_f = Meta.parse(orig_v["Formula"]) - f = Meta.parse(v["Formula"]) - - if is_default_function(orig_f, orig_v["RangeFrom"], orig_v["RangeTo"]) - @test isnothing(f) - else - @test orig_f == f - end - end - orig_variables_no_f = [ - Dict(k => v for (k, v) in var if k != "Formula" && k != "Name") for - var in orig_dict["Model"]["Variables"] - ] - output_variables_no_f = [ - Dict(k => v for (k, v) in var if k != "Formula" && k != "Name") for - var in variables - ] - @test orig_variables_no_f == output_variables_no_f - - @test haskey(model_dict, "Relationships") - relationships = model_dict["Relationships"] - orig_relationships_no_id = [ - Dict(k => v for (k, v) in rel if k != "Id") for - rel in orig_dict["Model"]["Relationships"] - ] - output_relationships_no_id = - [Dict(k => v for (k, v) in rel if k != "Id") for rel in relationships] - - # At least check that we are not creating new relationships out of nothing - @test length(setdiff(output_relationships_no_id, orig_relationships_no_id)) == 0 - - if occursin("VPC.json", model_path) - @test_broken false # VPC has some weird relationships that are not used in the target functions - else - @test Set(orig_relationships_no_id) == Set(output_relationships_no_id) - end - end - bma_models_path = joinpath(@__DIR__, "..", "resources", "bma_models") - good_models = joinpath(bma_models_path, "well_formed_examples") - - # just another reminder that the "Skin1D" example isn't working with this test - @test_broken false - - for model_path in filter(!contains(r"Skin1D"), readdir(good_models; join = true)) - test_json_roundtrip(model_path) - end -end - -@testitem "is default function" begin - import IterTools: subsets - import GraphDynamicalSystems.BMA: - is_default_function, default_target_function, swap_entity_names_to_var_ids - combinations = Iterators.filter( - x -> !all(isempty.(x)), - Iterators.product(subsets([:A_1, :B_2, :X_5, :Y_6]), subsets([:C_3, :D_4, :Z_7])), - ) - activators = first.(combinations) - inhibitors = last.(combinations) - - fns = default_target_function.(0, 4, activators, inhibitors) - for f in swap_entity_names_to_var_ids.(fns) - @test is_default_function(f, 0, 4) - end -end diff --git a/test/qn_test.jl b/test/qn_test.jl index 71c203a..ba0932c 100644 --- a/test/qn_test.jl +++ b/test/qn_test.jl @@ -4,122 +4,125 @@ seed!(42) end @testsetup module ExampleQN -export qn_size, max_eq_depth, domains, qn -using GraphDynamicalSystems +export qn_size, qn_fns, qn_domains, qn, qn_entities +import GraphDynamicalSystems as GDS qn_size = 3 -max_eq_depth = 3 -domains = [1:5, 1:5, 1:5] -qn = sample_qualitative_network(qn_size, domains, max_eq_depth) +qn_entities = [:X, :Y, :Z] +qn_fns = Dict(:X => 1, :Y => :X, :Z => :(Y - X)) +qn_domains = Dict(e => 0:5 for e in qn_entities) +qn = GDS.QN(qn_fns, qn_domains) end -@testitem "QN Grammar Creation" begin - entities = [:a, :b, :c] - constants = [i for i = 0:10] - g = build_qn_grammar(entities, constants) - - @test issubset(Set(entities), Set(g.rules)) - @test issubset(Set(constants), Set(g.rules)) +@testitem "From a dictionary of functions" begin + import Graphs + import SciMLBase + import GraphDynamicalSystems as GDS + import DynamicalSystemsBase as DSB + + fns = Dict(:A => 1, :B => :A, :C => :(B - A)) + d = Dict(k => 0:3 for k in keys(fns)) + qn = GDS.QualitativeNetwork{Graphs.SimpleGraph}(fns, d) + SciMLBase.step!(qn) + @test DSB.current_state(qn) == [1, 0, 0] + DSB.set_state!(qn, [1, 2, 3]) + @test DSB.current_state(qn) == [1, 2, 3] + SciMLBase.step!(qn) + @test DSB.current_state(qn) == [1, 1, 2] + SciMLBase.step!(qn) + @test DSB.current_state(qn) == [1, 1, 1] + SciMLBase.step!(qn) + @test DSB.current_state(qn) == [1, 1, 0] end @testitem "QN Graph Correctness" begin - import GraphDynamicalSystems: EntityName + import GraphDynamicalSystems: QN, get_graph import MetaGraphsNext: edge_labels - entity_labels = [:a, :b, :c] - target_fns = Union{Expr,Integer,Symbol}[:(-c), :a, :b] - domains = [0:2 for _ = 1:3] + target_fns = Dict(:a => :(-c), :b => :a, :c => :b) + domains = Dict(:a => 0:2, :b => 0:2, :c => 0:2) - qn = QN(Entity.(entity_labels, target_fns, domains)) + qn = QN(target_fns, domains) g = get_graph(qn) - @test haskey(g, EntityName(:c), EntityName(:a)) - @test haskey(g, EntityName(:a), EntityName(:b)) - @test haskey(g, EntityName(:b), EntityName(:c)) -end - -@testitem "QN Sampling" setup = [RandomSetup, ExampleQN] begin - using Graphs: ne, nv - graph = get_graph(qn) - - @test nv(graph) == qn_size - @test ne(graph) > 0 + @test haskey(g, :c, :a) + @test haskey(g, :a, :b) + @test haskey(g, :b, :c) end @testitem "QN properties, fields" setup = [RandomSetup, ExampleQN] begin - import GraphDynamicalSystems: EntityName - using DynamicalSystemsBase: get_state, set_state!, step! - - set_state!(qn, :A, 1) + import GraphDynamicalSystems: entities, get_fn, get_domain + using DynamicalSystemsBase: + current_state, + set_state!, + step!, + isinplace, + dynamic_rule, + current_parameters, + current_time + set_state!(qn, 1, :X) @test length(entities(qn)) == qn_size - @test length(target_functions(qn)) == qn_size + @test length(get_fn(qn)) == qn_size - @test all(get_state(qn) .<= maximum.(get_domain(qn))) + @test all(current_state(qn) .<= maximum.(get_domain(qn))) - @test get_state(qn, :A) == 1 + @test current_state(qn, :X) == 1 - @test_throws r"max" set_state!(qn, :A, 6) + @test_throws r"domain" set_state!(qn, 6, :X) + @test isinplace(qn) + @test dynamic_rule(qn) == get_fn(qn) + @test current_parameters(qn) == () + @test current_time(qn) == 0 end @testitem "Target Function" setup = [RandomSetup, ExampleQN] begin - using DynamicalSystemsBase: step!, get_state, set_state! - set_state!(qn, :A, 1) - set_state!(qn, :B, 1) - set_state!(qn, :C, 1) + import GraphDynamicalSystems: interpret + using DynamicalSystemsBase: set_state! + set_state!(qn, 1, :X) + set_state!(qn, 1, :Y) + set_state!(qn, 1, :Z) # All state values should be 1, so adding two of them == 2, etc. # The size of the test network is 3, so there should be A, B, C # as available entities to work with. - @test interpret(:(A + B), qn) == 2 - @test interpret(:(A - B), qn) == 0 - set_state!(qn, :B, 2) - @test interpret(:(A / B), qn) == 0.5 - @test interpret(:(A / 2), qn) == 0.5 - @test interpret(:(min(A, B)), qn) == 1 - @test interpret(:(max(A, B)), qn) == 2 - @test interpret(:(ceil(A / B)), qn) == 1 - @test interpret(:(floor(A / B)), qn) == 0 + @test interpret(:(X + Y), qn) == 2 + @test interpret(:(X - Y), qn) == 0 + set_state!(qn, 2, :Y) + @test interpret(:(X / Y), qn) == 0.5 + @test interpret(:(X / 2), qn) == 0.5 + @test interpret(:(min(X, Y)), qn) == 1 + @test interpret(:(max(X, Y)), qn) == 2 + @test interpret(:(ceil(X / Y)), qn) == 1 + @test interpret(:(floor(X / Y)), qn) == 0 @test_throws r"Unhandled" interpret(:(nonexistent_function(A)), qn) end -@testitem "Async QN" setup = [RandomSetup] begin - using DynamicalSystemsBase: step!, get_state, set_state! +@testitem "Async QN" setup = [RandomSetup, ExampleQN] begin + using DynamicalSystemsBase: step!, current_state, set_state! + import GraphDynamicalSystems: Asynchronous, QN qn_size = 3 max_eq_depth = 3 for N = 2:5 # a few different levels of N for _ = 1:100 # 100 different initializations - domains = [1:N for _ = 1:qn_size] - async_qn = sample_qualitative_network( - qn_size, - domains, - max_eq_depth; - schedule = Asynchronous, - ) - async_qn_system = create_qn_system(async_qn) - step!(async_qn_system, 100) - @test all(get_state(async_qn_system.model) .<= maximum.(domains)) + domains = Dict(e => 0:N for e in qn_entities) + async_qn = QN(qn_fns, domains) + step!(async_qn, 100) + @test all(current_state(async_qn) .<= maximum.(values(domains))) end end end @testitem "Get attractors" setup = [RandomSetup, ExampleQN] begin + import GraphDynamicalSystems: Asynchronous using Attractors: AttractorsViaRecurrences, basins_of_attraction - qn_size = 3 - max_eq_depth = 3 - N = 3 - domains = [1:N for _ = 1:qn_size] - - async_qn = - sample_qualitative_network(qn_size, domains, max_eq_depth; schedule = Asynchronous) - async_qn_system = create_qn_system(async_qn) grid = Tuple(range(0, 1) for _ = 1:qn_size) - mapper = AttractorsViaRecurrences(async_qn_system, grid) + mapper = AttractorsViaRecurrences(qn, grid) basins = basins_of_attraction(mapper, grid) end diff --git a/test/quick.jl b/test/quick.jl index b99bb1b..455cba4 100644 --- a/test/quick.jl +++ b/test/quick.jl @@ -8,11 +8,5 @@ using TestEnv using ReTestItems using GraphDynamicalSystems -TestEnv.activate() do - runtests( - GraphDynamicalSystems, - name = r"^(?!Code).+$", - failfast = true, - failures_first = true, - ) -end +TestEnv.activate() +rt() = runtests(GraphDynamicalSystems, failfast = true, failures_first = true)