diff --git a/HISTORY.md b/HISTORY.md index aba2af7841..78d89beb1e 100644 --- a/HISTORY.md +++ b/HISTORY.md @@ -41,7 +41,18 @@ (k1, k2), A <--> B end ``` - +- Catalyst's network visualization capability has shifted from using Graphviz to [GraphMakie.jl](https://graph.makie.org/stable/). To use this functionality, load the GraphMakie extension by installing `Catalyst` and `GraphMakie`, along with a Makie backend like `GLMakie`. There are two new methods for visualizing graphs: `plot_network` and `plot_complexes`, which respectively display the species-reaction graph and complex graph. + ```julia + using Catalyst, GraphMakie, GLMakie + brusselator = @reaction_network begin + A, ∅ --> X + 1, 2X + Y --> 3X + B, X --> Y + 1, X --> ∅ + end + plot_network(brusselator) + ``` + ## Catalyst 14.4.1 - Support for user-defined functions on the RHS when providing coupled equations for CRNs using the @equations macro. For example, the following now works: diff --git a/Project.toml b/Project.toml index 4e6035d919..b39476030b 100644 --- a/Project.toml +++ b/Project.toml @@ -75,7 +75,6 @@ julia = "1.10" [extras] DiffEqCallbacks = "459566f4-90b8-5000-8ac3-15dfb0a30def" DomainSets = "5b8099bc-c8ec-5219-889f-1d9e522a28bf" -Graphviz_jll = "3c863552-8265-54e4-a6dc-903eb78fde85" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Logging = "56ddb016-857b-54e1-b83d-db4d58db5568" NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec" @@ -99,4 +98,4 @@ Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" [targets] -test = ["DiffEqCallbacks", "DomainSets", "Graphviz_jll", "Logging", "NonlinearSolve", "OrdinaryDiffEqBDF", "OrdinaryDiffEqDefault", "OrdinaryDiffEqRosenbrock", "OrdinaryDiffEqTsit5", "OrdinaryDiffEqVerner", "Pkg", "Plots", "Random", "SafeTestsets", "SciMLBase", "SciMLNLSolve", "StableRNGs", "StaticArrays", "Statistics", "SteadyStateDiffEq", "StochasticDiffEq", "Test", "Unitful"] +test = ["DiffEqCallbacks", "DomainSets", "Logging", "NonlinearSolve", "OrdinaryDiffEqBDF", "OrdinaryDiffEqDefault", "OrdinaryDiffEqRosenbrock", "OrdinaryDiffEqTsit5", "OrdinaryDiffEqVerner", "Pkg", "Plots", "Random", "SafeTestsets", "SciMLBase", "SciMLNLSolve", "StableRNGs", "StaticArrays", "Statistics", "SteadyStateDiffEq", "StochasticDiffEq", "Test", "Unitful"] diff --git a/docs/Project.toml b/docs/Project.toml index 73163c65ef..10e6eeda3f 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -18,6 +18,7 @@ Latexify = "23fbe1c1-3f47-55db-b15f-69d7ec21a316" LinearSolve = "7ed4a6bd-45f5-4d41-b270-4a48e9bafcae" Logging = "56ddb016-857b-54e1-b83d-db4d58db5568" ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78" +NetworkLayout = "46757867-2c16-5918-afeb-47bfcb05e46a" NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec" Optim = "429524aa-4258-5aef-a3af-852621145aeb" Optimization = "7f7a1694-90dd-40f0-9382-eb1efda571ba" diff --git a/docs/make.jl b/docs/make.jl index a877fb5d5e..f732463793 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -1,5 +1,7 @@ using Documenter using Catalyst, ModelingToolkit +# Add packages for plotting +using GraphMakie, CairoMakie docpath = Base.source_dir() assetpath = joinpath(docpath, "src", "assets") @@ -37,7 +39,9 @@ makedocs(sitename = "Catalyst.jl", collapselevel = 1, assets = ["assets/favicon.ico"], canonical = "https://docs.sciml.ai/Catalyst/stable/"), - modules = [Catalyst, ModelingToolkit], + modules = [Catalyst, ModelingToolkit, + isdefined(Base, :get_extension) ? Base.get_extension(Catalyst, :CatalystGraphMakieExtension) : + Catalyst.CatalystGraphMakieExtension], doctest = false, clean = true, pages = pages, diff --git a/docs/pages.jl b/docs/pages.jl index bc77aaf459..e17fbbe7b4 100644 --- a/docs/pages.jl +++ b/docs/pages.jl @@ -40,7 +40,7 @@ pages = Any[ "Steady state analysis" => Any[ "steady_state_functionality/homotopy_continuation.md", "steady_state_functionality/nonlinear_solve.md", - "steady_state_functionality/steady_state_stability_computation.md", + #"steady_state_functionality/steady_state_stability_computation.md", "steady_state_functionality/bifurcation_diagrams.md", "steady_state_functionality/dynamical_systems.md" ], diff --git a/docs/src/api.md b/docs/src/api.md index 573076851e..f78a5f70e4 100644 --- a/docs/src/api.md +++ b/docs/src/api.md @@ -271,7 +271,7 @@ isequivalent ==(rn1::ReactionSystem, rn2::ReactionSystem) ``` -## Network visualization +## [Network visualization](@id network_visualization) [Latexify](https://korsbo.github.io/Latexify.jl/stable/) can be used to convert networks to LaTeX equations by ```julia @@ -292,13 +292,10 @@ displayed as the ODE form) Finally, another optional argument (`expand_functions=true`) automatically expands functions defined by Catalyst (such as `mm`). To disable this, set `expand_functions=false`. -If [Graphviz](https://graphviz.org/) is installed and commandline accessible, it -can be used to create and save network diagrams using [`Graph`](@ref) and -[`savegraph`](@ref). +Reaction networks can be plotted using the `GraphMakie` extension, which is loaded whenever all of `Catalyst`, `GraphMakie`, and `NetworkLayout` are loaded (note that a Makie backend, like `CairoMakie`, must be loaded as well). The two functions for plotting networks are `plot_network` and `plot_complexes`, which are two distinct representations. ```@docs -Graph -complexgraph -savegraph +plot_network(::ReactionSystem) +plot_complexes(::ReactionSystem) ``` ## [Rate laws](@id api_rate_laws) diff --git a/docs/src/assets/network_graphs/brusselator_complex_graph_makie.png b/docs/src/assets/network_graphs/brusselator_complex_graph_makie.png new file mode 100644 index 0000000000..b22f527866 Binary files /dev/null and b/docs/src/assets/network_graphs/brusselator_complex_graph_makie.png differ diff --git a/docs/src/assets/network_graphs/brusselator_graph_makie.png b/docs/src/assets/network_graphs/brusselator_graph_makie.png new file mode 100644 index 0000000000..0c1ec1101c Binary files /dev/null and b/docs/src/assets/network_graphs/brusselator_graph_makie.png differ diff --git a/docs/src/assets/network_graphs/repressilator_graph_makie.png b/docs/src/assets/network_graphs/repressilator_graph_makie.png new file mode 100644 index 0000000000..bf151a0de1 Binary files /dev/null and b/docs/src/assets/network_graphs/repressilator_graph_makie.png differ diff --git a/docs/src/index.md b/docs/src/index.md index 03977c127f..c30cd42ca3 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -29,7 +29,7 @@ etc). - Models can be [coupled with events](@ref constraint_equations_events) that affect the system and its state during simulations. - By leveraging ModelingToolkit, users have a variety of options for generating optimized system representations to use in solvers. These include construction of [dense or sparse Jacobians](@ref ode_simulation_performance_sparse_jacobian), [multithreading or parallelization of generated derivative functions](@ref ode_simulation_performance_parallelisation), [automatic classification of reactions into optimized jump types for Gillespie type simulations](https://docs.sciml.ai/JumpProcesses/stable/jump_types/#jump_types), [automatic construction of dependency graphs for jump systems](https://docs.sciml.ai/JumpProcesses/stable/jump_types/#Jump-Aggregators-Requiring-Dependency-Graphs), and more. - [Symbolics.jl](https://github.com/JuliaSymbolics/Symbolics.jl) symbolic expressions and Julia `Expr`s can be obtained for all rate laws and functions determining the deterministic and stochastic terms within resulting ODE, SDE, or jump models. -- [Steady states](@ref homotopy_continuation) (and their [stabilities](@ref steady_state_stability)) can be computed for model ODE representations. +- [Steady states](@ref homotopy_continuation) can be computed for model ODE representations. #### [Features of Catalyst composing with other packages](@id doc_index_features_composed) - [OrdinaryDiffEq.jl](https://github.com/SciML/OrdinaryDiffEq.jl) Can be used to numerically solve generated reaction rate equation ODE models. @@ -37,7 +37,7 @@ etc). - [JumpProcesses.jl](https://github.com/SciML/JumpProcesses.jl) can be used to numerically sample generated Stochastic Chemical Kinetics Jump Process models. - Support for [parallelization of all simulations](@ref ode_simulation_performance_parallelisation), including parallelization of [ODE simulations on GPUs](@ref ode_simulation_performance_parallelisation_GPU) using [DiffEqGPU.jl](https://github.com/SciML/DiffEqGPU.jl). - [Latexify](https://korsbo.github.io/Latexify.jl/stable/) can be used to [generate LaTeX expressions](@ref visualisation_latex) corresponding to generated mathematical models or the underlying set of reactions. -- [Graphviz](https://graphviz.org/) can be used to generate and [visualize reaction network graphs](@ref visualisation_graphs) (reusing the Graphviz interface created in [Catlab.jl](https://algebraicjulia.github.io/Catlab.jl/stable/)). +- [GraphMakie](https://docs.makie.org/stable/) recipes are provided that can be used to generate and [visualize reaction network graphs](@ref visualisation_graphs) - Model steady states can be [computed through homotopy continuation](@ref homotopy_continuation) using [HomotopyContinuation.jl](https://github.com/JuliaHomotopyContinuation/HomotopyContinuation.jl) (which can find *all* steady states of systems with multiple ones), by [forward ODE simulations](@ref steady_state_solving_simulation) using [SteadyStateDiffEq.jl)](https://github.com/SciML/SteadyStateDiffEq.jl), or by [numerically solving steady-state nonlinear equations](@ref steady_state_solving_nonlinear) using [NonlinearSolve.jl](https://github.com/SciML/NonlinearSolve.jl). [BifurcationKit.jl](https://github.com/bifurcationkit/BifurcationKit.jl) can be used to compute bifurcation diagrams of model steady states (including finding periodic orbits). - [DynamicalSystems.jl](https://github.com/JuliaDynamics/DynamicalSystems.jl) can be used to compute model [basins of attraction](@ref dynamical_systems_basins_of_attraction), [Lyapunov spectrums](@ref dynamical_systems_lyapunov_exponents), and other dynamical system properties. diff --git a/docs/src/introduction_to_catalyst/introduction_to_catalyst.md b/docs/src/introduction_to_catalyst/introduction_to_catalyst.md index 8f1dc1d5c5..a766468add 100644 --- a/docs/src/introduction_to_catalyst/introduction_to_catalyst.md +++ b/docs/src/introduction_to_catalyst/introduction_to_catalyst.md @@ -81,17 +81,18 @@ latexify(rn) ```@example tut1 rn #hide ``` -Assuming [Graphviz](https://graphviz.org/) is installed and command line -accessible, within a Jupyter notebook we can also graph the reaction network by -```julia -g = Graph(rn) +Catalyst also has functionality for visualizing networks using the [Makie](https://docs.makie.org/stable/) +plotting ecosystem. The relevant packages to load are Catalyst, GraphMakie, NetworkLayout, and a Makie backend +such as CairoMakie. Doing so and then using the `plot_network` function allows us to +visualize the network: +```@example tut1 +using Catalyst +import CairoMakie, GraphMakie, NetworkLayout +g = plot_network(rn) ``` -giving - -![Repressilator solution](../assets/repressilator.svg) The network graph shows a variety of information, representing each species as a -blue node, and each reaction as an orange dot. Black arrows from species to +blue node, and each reaction as an green node. Black arrows from species to reactions indicate reactants, and are labelled with their input stoichiometry. Similarly, black arrows from reactions to species indicate products, and are labelled with their output stoichiometry. In contrast, red arrows from a species @@ -105,8 +106,11 @@ hillr(P₂,α,K,n), ∅ --> m₃ have rates that depend on the proteins, and hence lead to red arrows from each `Pᵢ`. -Note, from the REPL or scripts one can always use [`savegraph`](@ref) to save -the graph (assuming `Graphviz` is installed). +Note, from the REPL or scripts one can always use Makie's `save` function to save +the graph. +```julia +save("repressilator_graph.png", g) +``` ## Mass action ODE models Let's now use our `ReactionSystem` to generate and solve a corresponding mass @@ -176,7 +180,7 @@ At this point we are all set to solve the ODEs. We can now use any ODE solver from within the [OrdinaryDiffEq.jl](https://docs.sciml.ai/DiffEqDocs/stable/solvers/ode_solve/) package. We'll use the recommended default explicit solver, `Tsit5()`, and then -plot the solutions: +plot the solutions: ```@example tut1 sol = solve(oprob, Tsit5(), saveat=10.0) diff --git a/docs/src/model_creation/model_visualisation.md b/docs/src/model_creation/model_visualisation.md index 96768ecae0..d378deb3a6 100644 --- a/docs/src/model_creation/model_visualisation.md +++ b/docs/src/model_creation/model_visualisation.md @@ -36,21 +36,26 @@ If you wish to copy the output to your [clipboard](https://en.wikipedia.org/wiki For a model to be nicely displayed you have to use an IDE that actually supports this (such as a [notebook](https://jupyter.org/)). Other environments (such as [the Julia REPL](https://docs.julialang.org/en/v1/stdlib/REPL/)) will simply return the full LaTeX code which would generate the desired expression. ## [Displaying model networks](@id visualisation_graphs) -A network graph showing a Catalyst model's species and reactions can be displayed using the `Graph` function. This first requires [Graphviz](https://graphviz.org/) to be installed and command line accessible. Here, we first declare a [Brusselator model](@ref basic_CRN_library_brusselator) and then displays its network topology: +Catalyst uses `GraphMakie` to display representations of chemical reaction networks, including the complex graph and the species-reaction graph (which is similar to the [Petri net](https://en.wikipedia.org/wiki/Petri_net) representation). To get started, import Catalyst, GraphMakie, and NetworkLayout to load the `CatalystGraphMakieExtension` extension, and then load a Makie backend (`CairoMakie` is a good lightweight choice). + +```@example visualisation_graphs +using Catalyst, GraphMakie, NetworkLayout +using CairoMakie +nothing # hide +``` + +Let's declare a [Brusselator model](@ref basic_CRN_library_brusselator) to see this plotting functionality. The functions `plot_network` and `plot_complexes` are used to create the species-reaction and complex graphs, respectively. For a more thorough description of these two representations, please see the [network visualization](@ref network_visualization) section of the API, but the gist is that the species-reaction graph has species and reactions as nodes, and the complex graph has reaction complexes as nodes. Below we will plot the species-reaction graph using `plot_network`. ```@example visualisation_graphs -using Catalyst brusselator = @reaction_network begin A, ∅ --> X 1, 2X + Y --> 3X B, X --> Y 1, X --> ∅ end -Graph(brusselator) -nothing # hide +plot_network(brusselator) ``` -!["Brusselator Graph"](../assets/network_graphs/brusselator_graph.png) -The network graph represents species as blue nodes and reactions as orange dots. Black arrows from species to reactions indicate substrates, and are labelled with their respective stoichiometries. Similarly, black arrows from reactions to species indicate products (also labelled with their respective stoichiometries). If there are any reactions where a species affect the rate, but does not participate as a reactant, this is displayed with a dashed red arrow. This can be seen in the following [Repressilator model](@ref basic_CRN_library_repressilator): +The species-reaction graph (or network graph) represents species as blue nodes and reactions as green dots. Black arrows from species to reactions indicate substrates, and are labelled with their respective stoichiometries. Similarly, black arrows from reactions to species indicate products (also labelled with their respective stoichiometries). If there are any reactions where a species affect the rate, but does not participate as a reactant, this is displayed with a dashed red arrow. This can be seen in the following [Repressilator model](@ref basic_CRN_library_repressilator): ```@example visualisation_graphs repressilator = @reaction_network begin hillr(Z,v,K,n), ∅ --> X @@ -58,21 +63,78 @@ repressilator = @reaction_network begin hillr(Y,v,K,n), ∅ --> Z d, (X, Y, Z) --> ∅ end -Graph(repressilator) -nothing # hide +plot_network(repressilator) ``` -!["Repressilator Graph"](../assets/network_graphs/repressilator_graph.png) -A generated graph can be saved using the `savegraph` function: +A generated graph can be saved using Makie's `save` function. ```julia -repressilator_graph = Graph(repressilator) -savegraph(repressilator_graph, "repressilator_graph.png") +repressilator_graph = plot_network(repressilator) +save("repressilator_graph.png", repressilator_graph) ``` -Finally, a [network's reaction complexes](@ref network_analysis_reaction_complexes) (and the reactions in between these) can be displayed using the `complexgraph(brusselator)` function: +Finally, a [network's reaction complexes](@ref network_analysis_reaction_complexes) (and the reactions in between these) can be displayed using the `plot_complexes(brusselator)` function: ```@example visualisation_graphs -complexgraph(brusselator) -nothing # hide +plot_complexes(brusselator) +``` +Here, reaction complexes are displayed as blue nodes, and reactions between complexes are displayed as black arrows. Red arrows indicate that the rate constantof a reaction has a species-dependence. Edges can be optionally labeled with their rate expressions by calling with the option `show_rate_labels`. +```@example visualisation_graphs +plot_complexes(brusselator, show_rate_labels = true) +``` + +## Customizing Plots +In this section we demonstrate some of the ways that plot objects can be manipulated to give nicer images. Let's start with our brusselator plot once again. Note that the `plot` function returns three objects: the `Figure`, the `Axis`, and the `Plot`, which can each be customized independently. See the general [Makie documentation](https://docs.makie.org/stable/) for more information. + +```@example visualisation_graphs +f, ax, p = plot_complexes(brusselator, show_rate_labels = true) +``` + +It seems like a bit of the top node is cut off. Let's hide the tick marks and grid and increase the top and bottom margins by increasing `yautolimitmargin`. +```@example visualisation_graphs +hidedecorations!(ax) +ax.yautolimitmargin = (0.1, 0.1) # defaults to (0.05, 0.05) +ax.aspect = DataAspect() +``` + +There are many keyword arguments that can be passed to `plot_network` or `plot_complexes` to change the look of the graph (which get passed to the `graphplot` Makie recipe). Let's change the color of the nodes and make the inner labels a bit smaller. As before, we hide the tick marks and grid. Let's also give the plot a title. +```@example visualisation_graphs +f, ax, p = plot_complexes(brusselator, show_rate_labels = true, node_color = :yellow, ilabels_fontsize = 10) +hidedecorations!(ax) +ax.yautolimitmargin = (0.1, 0.1) # defaults to (0.05, 0.05) +ax.aspect = DataAspect() +ax.title = "Brusselator" +``` + +Most of the kwargs that modify the nodes or edges will also accept a vector with the same length as the number of nodes or edges, respectively. See [here](https://graph.makie.org/stable/#The-graphplot-Recipe) for a full list of keyword arguments to `graph_plot`. Note that `plot_complexes` and `plot_network` default to `layout = Stress()` rather than `layout = Spring()`, since `Stress()` is better at generating plots with fewer edge crossings. More layout options and customizations (such as pinning nodes to certain positions) can be found in the [`NetworkLayout` documentation](https://juliagraphs.org/NetworkLayout.jl/stable/). + +Once a graph is already created we can also change the keyword arguments by modifying the fields of the `Plot` object `p`. +```@example visualisation_graphs +p.node_color = :orange +``` + +Custom node positions can also be given, if the automatic layout is unsatisfactory. +```@example visualisation_graphs +fixedlayout = [(0,0), (1,0), (0,1), (1,1), (2,0)] +p.layout = fixedlayout +autolimits!(ax) +``` + +Makie graph plots can be made to be interactive, allowing one to drag nodes and edges. To do this, we retrieve the axis from the GraphMakie plot, and then register the interactions. **Note that this can only be done if `GLMakie` is the installed Makie backend. See the [GraphMakie docs](https://graph.makie.org/stable/#Predefined-Interactions) for more information about the types of interactions one can register.** Below is a non-interactive code example that shows how to do this: + +```julia +using GLMakie +f, ax, p = plot_network(brusselator) +deregister_interaction!(ax, :rectanglezoom) +register_interaction!(ax, :ndrag, NodeDrag(p)) +register_interaction!(ax, :edrag, EdgeDrag(p)) +``` + +The equivalent of `show` for Makie plots is the `display` function. +```julia +f = plot_network(brusselator) +display(f) +``` + +Once you are happy with the graph plot, you can save it using the `save` function. +```julia +save("fig.png", f) ``` -!["Repressilator Complex Graph"](../assets/network_graphs/repressilator_complex_graph.png) -Here, reaction complexes are displayed as blue nodes, and reactions in between these as black arrows. diff --git a/docs/src/model_creation/network_analysis.md b/docs/src/model_creation/network_analysis.md index 3f77a99f55..8085de1b60 100644 --- a/docs/src/model_creation/network_analysis.md +++ b/docs/src/model_creation/network_analysis.md @@ -31,10 +31,10 @@ In the [Introduction to Catalyst](@ref introduction_to_catalyst) tutorial we showed how the above network could be visualized as a species-reaction graph. There, species are represented by the nodes of the graph and edges show the reactions in which a given species is a substrate or product. -```julia -g = Graph(repressilator) +```@example s1 +using Catalyst, GraphMakie, NetworkLayout, CairoMakie +g = plot_network(repressilator) ``` -![Repressilator solution](../assets/repressilator.svg) We also showed in the [Introduction to Catalyst](@ref introduction_to_catalyst) tutorial that the reaction rate equation ODE model for the repressilator is @@ -187,22 +187,18 @@ N == Z*B ``` Reaction complexes give an alternative way to visualize a reaction network -graph. Catalyst's [`complexgraph`](@ref) command will calculate the complexes of +graph. Catalyst's [`plot_complexes`](@ref) command will calculate the complexes of a network and then show how they are related. For example, -```julia -complexgraph(rn) +```@example s1 +using Catalyst, GraphMakie, NetworkLayout, CairoMakie +plot_complexes(rn) ``` -gives - -![Simple example complex graph](../assets/simple_complexgraph.svg) while for the repressilator we find -```julia -complexgraph(repressilator) +```@example s1 +plot_complexes(repressilator) ``` -![Repressilator complex](../assets/repressilator_complexgraph.svg) - Here ∅ represents the empty complex, black arrows show reactions converting substrate complexes into product complexes where the rate is just a number or parameter, and red arrows indicate the conversion of substrate complexes into @@ -228,11 +224,10 @@ rn = @reaction_network begin end ``` with graph -```julia -complexgraph(rn) +```@example s1 +plot_complexes(rn) ``` -![network_1](../assets/complex_rn.svg) #### [Linkage classes and sub-networks of the reaction network](@id network_analysis_structural_aspects_linkage) The preceding reaction complex graph shows that `rn` is composed of two @@ -260,19 +255,15 @@ subnets = subnetworks(rn) reactions.(subnets) ``` The graphs of the reaction complexes in the two sub-networks are then -```julia - complexgraph(subnets[1]) +```@example s1 + plot_complexes(subnets[1]) ``` -![subnetwork_1](../assets/complex_subnets1.svg) - and, -```julia - complexgraph(subnets[2]) +```@example s1 + plot_complexes(subnets[2]) ``` -![subnetwork_2](../assets/complex_subnets2.svg) - #### [Deficiency of the network](@id network_analysis_structural_aspects_deficiency) A famous theorem in Chemical Reaction Network Theory, the Deficiency Zero Theorem [^1], allows us to use knowledge of the net stoichiometry matrix and the @@ -357,12 +348,10 @@ end reactioncomplexes(rn) isreversible(rn) ``` -```julia -complexgraph(rn) +```@example s1 +plot_complexes(rn) ``` -![reversibility](../assets/complex_reversibility.svg) - It is evident from the preceding graph that the network is not reversible. However, it satisfies a weaker property in that there is a path from each reaction complex back to itself within its associated subgraph. This is known as diff --git a/ext/CatalystGraphMakieExtension/rn_graph_plot.jl b/ext/CatalystGraphMakieExtension/rn_graph_plot.jl index 1ed92aa440..26271e2b38 100644 --- a/ext/CatalystGraphMakieExtension/rn_graph_plot.jl +++ b/ext/CatalystGraphMakieExtension/rn_graph_plot.jl @@ -191,9 +191,9 @@ function Catalyst.plot_network(rn::ReactionSystem; kwargs...) end layout = if !haskey(kwargs, :layout) - is_connected(srg) ? Stress() : Spring() + Stress() end - graphplot(srg; + f = graphplot(srg; layout, edge_color = edgecolors, elabels = edgelabels, @@ -206,23 +206,30 @@ function Catalyst.plot_network(rn::ReactionSystem; kwargs...) curve_distance = gen_distances(srg), kwargs... ) + + f.axis.xautolimitmargin = (0.15, 0.15) + f.axis.yautolimitmargin = (0.15, 0.15) + + f end """ - plot_complexes(rn::ReactionSystem; kwargs...) + plot_complexes(rn::ReactionSystem; show_rate_labels = false, kwargs...) - Creates a GraphMakie plot of the [`ReactionComplex`](@ref)s in `rn`. Reactions - correspond to arrows and reaction complexes to blue circles. +Creates a GraphMakie plot of the [`Catalyst.ReactionComplex`](@ref)s in `rn`. Reactions +correspond to arrows and reaction complexes to blue circles. - Notes: - - Black arrows from complexes to complexes indicate reactions whose rate is a - parameter or a `Number`. i.e. `k, A --> B`. - - Red arrows from complexes to complexes indicate reactions whose rate - depends on species. i.e. `k*C, A --> B` for `C` a species. +Notes: +- Black arrows from complexes to complexes indicate reactions whose rate is a + parameter or a `Number`. i.e. `k, A --> B`. +- Red arrows from complexes to complexes indicate reactions whose rate constants +depends on species. i.e. `k*C, A --> B` for `C` a species. +- The `show_rate_labels` keyword, if set to `true`, will annotate each edge +with the rate constant for the reaction. For a list of accepted keyword arguments to the graph plot, please see the [GraphMakie documentation](https://graph.makie.org/stable/#The-graphplot-Recipe). """ -function Catalyst.plot_complexes(rn::ReactionSystem; kwargs...) +function Catalyst.plot_complexes(rn::ReactionSystem; show_rate_labels = false, kwargs...) rxs = reactions(rn) specs = species(rn) edgecolors = [:black for i in 1:length(rxs)] @@ -238,13 +245,13 @@ function Catalyst.plot_complexes(rn::ReactionSystem; kwargs...) # Get complex graph and reaction order for edgecolors and edgelabels. rxorder gives the order of reactions(rn) that would match the edge order in edges(cg). cg, rxorder = ComplexGraphWrap(rn) - layout = if !haskey(kwargs, :layout) - is_connected(cg) ? Stress() : Spring() + layout = if !haskey(kwargs, :layout) + Stress() end - graphplot(cg; + f = graphplot(cg; layout, edge_color = edgecolors[rxorder], - elabels = edgelabels[rxorder], + elabels = show_rate_labels ? edgelabels[rxorder] : [], ilabels = complexlabels(rn), node_color = :skyblue3, elabels_rotation = 0, @@ -253,6 +260,10 @@ function Catalyst.plot_complexes(rn::ReactionSystem; kwargs...) curve_distance = gen_distances(cg), kwargs... ) + f.axis.xautolimitmargin = (0.15, 0.15) + f.axis.yautolimitmargin = (0.15, 0.15) + + f end function complexelem_tostr(e::Catalyst.ReactionComplexElement, specstrs) diff --git a/src/Catalyst.jl b/src/Catalyst.jl index 1d809cca60..52daa0cec5 100644 --- a/src/Catalyst.jl +++ b/src/Catalyst.jl @@ -60,19 +60,6 @@ const DEFAULT_IV = default_t() const DEFAULT_IV_SYM = Symbol(DEFAULT_IV) export default_t, default_time_deriv -# as used in Catlab -const USE_GV_JLL = Ref(false) -function __init__() - @require Graphviz_jll="3c863552-8265-54e4-a6dc-903eb78fde85" begin - USE_GV_JLL[] = true - let cfg = joinpath(Graphviz_jll.artifact_dir, "lib", "graphviz", "config6") - if !isfile(cfg) - Graphviz_jll.dot(path -> run(`$path -c`)) - end - end - end -end - ### Package Constants ### # Union type of types that can occur in expressions. @@ -145,8 +132,6 @@ include("latexify_recipes.jl") # for making and saving graphs/plots include("plotting.jl") -include("graphs.jl") -export Graph, savegraph, complexgraph # for creating compounds include("chemistry_functionality.jl") @@ -168,7 +153,7 @@ export hc_steady_states function make_si_ode end export make_si_ode -# GraphMakie +# GraphMakie: functionality for plotting species-reaction graphs and complexes function plot_network end function plot_complexes end export plot_network, plot_complexes diff --git a/src/graphs.jl b/src/graphs.jl deleted file mode 100644 index 2f1fb490c4..0000000000 --- a/src/graphs.jl +++ /dev/null @@ -1,428 +0,0 @@ -####################################################################### -# Contains code from Catlab.jl: -# https://raw.githubusercontent.com/AlgebraicJulia/Catlab.jl/master/src/graphics/Graphviz.jl -# -# That license for that code is: -# -# The MIT License -# -# Copyright (c) 2017-2020: Evan Patterson. -# -# Permission is hereby granted, free of charge, to any person obtaining a copy -# of this software and associated documentation files (the "Software"), to deal -# in the Software without restriction, including without limitation the rights -# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell -# copies of the Software, and to permit persons to whom the Software is -# furnished to do so, subject to the following conditions: -# -# The above copyright notice and this permission notice shall be included in -# all copies or substantial portions of the Software. -# -# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR -# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, -# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE -# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER -# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, -# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN -# THE SOFTWARE. -####################################################################### -""" AST and pretty printer for Graphviz's DOT language. - -References: - -- DOT grammar: http://www.graphviz.org/doc/info/lang.html -- DOT language guide: http://www.graphviz.org/pdf/dotguide.pdf -""" - -### AST ### - -abstract type Expression end -abstract type Statement <: Expression end - -""" AST type for Graphviz's "HTML-like" node labels. - -For now, the HTML content is just a string. -""" -struct Html - content::String -end -Base.print(io::IO, html::Html) = print(io, html.content) - -const AttributeValue = Union{String, Html} -const Attributes = OrderedDict{Symbol, AttributeValue} - -as_attributes(attrs::Attributes) = attrs -as_attributes(d::OrderedDict) = Attributes(Symbol(k) => d[k] for k in keys(d)) -function as_attributes(d::AbstractDict) - Attributes(Symbol(k) => d[k] for k in sort!(collect(keys(d)))) -end - -@with_kw_noshow struct Graph <: Expression - name::String - directed::Bool - prog::String = "dot" - stmts::Vector{Statement} = Statement[] - graph_attrs::Attributes = Attributes() - node_attrs::Attributes = Attributes() - edge_attrs::Attributes = Attributes() -end - -function Graph(name::String, stmts::Vector{Statement}; kw...) - Graph(; name = name, directed = false, stmts = stmts, kw...) -end -function Graph(name::String, stmts::Vararg{Statement}; kw...) - Graph(; name = name, directed = false, stmts = collect(stmts), kw...) -end -function Digraph(name::String, stmts::Vector{Statement}; kw...) - Graph(; name = name, directed = true, stmts = stmts, kw...) -end -function Digraph(name::String, stmts::Vararg{Statement}; kw...) - Graph(; name = name, directed = true, stmts = collect(stmts), kw...) -end - -@with_kw_noshow struct Subgraph <: Statement - name::String = "" # Subgraphs can be anonymous - stmts::Vector{Statement} = Statement[] - graph_attrs::Attributes = Attributes() - node_attrs::Attributes = Attributes() - edge_attrs::Attributes = Attributes() -end - -Subgraph(stmts::Vector{Statement}; kw...) = Subgraph(; stmts = stmts, kw...) -Subgraph(stmts::Vararg{Statement}; kw...) = Subgraph(; stmts = collect(stmts), kw...) -function Subgraph(name::String, stmts::Vector{Statement}; kw...) - Subgraph(; name = name, stmts = stmts, kw...) -end -function Subgraph(name::String, stmts::Vararg{Statement}; kw...) - Subgraph(; name = name, stmts = collect(stmts), kw...) -end - -struct Node <: Statement - name::String - attrs::Attributes -end -Node(name::String, attrs::AbstractDict) = Node(name, as_attributes(attrs)) -Node(name::String; attrs...) = Node(name, attrs) - -struct NodeID <: Expression - name::String - port::String - anchor::String - NodeID(name::String, port::String = "", anchor::String = "") = new(name, port, anchor) -end - -struct Edge <: Statement - path::Vector{NodeID} - attrs::Attributes -end -Edge(path::Vector{NodeID}, attrs::AbstractDict) = Edge(path, as_attributes(attrs)) -Edge(path::Vector{NodeID}; attrs...) = Edge(path, attrs) -Edge(path::Vararg{NodeID}; attrs...) = Edge(collect(path), attrs) -Edge(path::Vector{String}, attrs::AbstractDict) = Edge(map(NodeID, path), attrs) -Edge(path::Vector{String}; attrs...) = Edge(map(NodeID, path), attrs) -Edge(path::Vararg{String}; attrs...) = Edge(map(NodeID, collect(path)), attrs) - -### Bindings ### - -""" Run a Graphviz program. - -Invokes Graphviz through its command-line interface. If the `Graphviz_jll` -package is installed and loaded, it is used; otherwise, Graphviz must be -installed on the local system. - -For bindings to the Graphviz C API, see the the package -[GraphViz.jl](https://github.com/Keno/GraphViz.jl). At the time of this writing, -GraphViz.jl is unmaintained. -""" -function run_graphviz(io::IO, graph::Graph; prog::Union{String, Nothing} = nothing, - format::String = "json0") - if isnothing(prog) - prog = graph.prog - end - @assert prog in ("dot", "neato", "fdp", "sfdp", "twopi", "circo") - if USE_GV_JLL[] - fun = getfield(Graphviz_jll, Symbol(prog)) - fun() do exe - open(`$exe -T$format`, io, write = true) do gv - pprint(gv, graph) - end - end - else - open(`$prog -T$format`, io, write = true) do gv - pprint(gv, graph) - end - end -end -function run_graphviz(graph::Graph; kw...) - io = IOBuffer() - run_graphviz(io, graph; kw...) - seekstart(io) -end - -function Base.show(io::IO, ::MIME"image/svg+xml", graph::Graph) - run_graphviz(io, graph, format = "svg") -end - -### Pretty-print ### - -""" Pretty-print the Graphviz expression. -""" -pprint(expr::Expression) = pprint(stdout, expr) -pprint(io::IO, expr::Expression) = pprint(io, expr, 0) - -function pprint(io::IO, graph::Graph, n::Int) - indent(io, n) - print(io, graph.directed ? "digraph " : "graph ") - print(io, graph.name) - println(io, " {") - pprint_attrs(io, graph.graph_attrs, n + 2; pre = "graph", post = ";\n") - pprint_attrs(io, graph.node_attrs, n + 2; pre = "node", post = ";\n") - pprint_attrs(io, graph.edge_attrs, n + 2; pre = "edge", post = ";\n") - for stmt in graph.stmts - pprint(io, stmt, n + 2, directed = graph.directed) - println(io) - end - indent(io, n) - println(io, "}") -end - -function pprint(io::IO, subgraph::Subgraph, n::Int; directed::Bool = false) - indent(io, n) - if isempty(subgraph.name) - println(io, "{") - else - print(io, "subgraph ") - print(io, subgraph.name) - println(io, " {") - end - pprint_attrs(io, subgraph.graph_attrs, n + 2; pre = "graph", post = ";\n") - pprint_attrs(io, subgraph.node_attrs, n + 2; pre = "node", post = ";\n") - pprint_attrs(io, subgraph.edge_attrs, n + 2; pre = "edge", post = ";\n") - for stmt in subgraph.stmts - pprint(io, stmt, n + 2, directed = directed) - println(io) - end - indent(io, n) - print(io, "}") -end - -function pprint(io::IO, node::Node, n::Int; directed::Bool = false) - indent(io, n) - print(io, node.name) - pprint_attrs(io, node.attrs) - print(io, ";") -end - -function pprint(io::IO, node::NodeID, n::Int) - print(io, node.name) - if !isempty(node.port) - print(io, ":") - print(io, node.port) - end - if !isempty(node.anchor) - print(io, ":") - print(io, node.anchor) - end -end - -function pprint(io::IO, edge::Edge, n::Int; directed::Bool = false) - indent(io, n) - for (i, node) in enumerate(edge.path) - if i > 1 - print(io, directed ? " -> " : " -- ") - end - pprint(io, node, n) - end - pprint_attrs(io, edge.attrs) - print(io, ";") -end - -function pprint_attrs(io::IO, attrs::Attributes, n::Int = 0; - pre::String = "", post::String = "") - if !isempty(attrs) - indent(io, n) - print(io, pre) - print(io, " [") - for (i, (key, value)) in enumerate(attrs) - if (i > 1) - print(io, ",") - end - print(io, key) - print(io, "=") - print(io, value isa Html ? "<" : "\"") - print(io, value) - print(io, value isa Html ? ">" : "\"") - end - print(io, "]") - print(io, post) - end -end - -indent(io::IO, n::Int) = print(io, " "^n) - -####################################################################### -# following is adapted from Petri.jl -# https://github.com/mehalter/Petri.jl -####################################################################### - -const graph_attrs = Attributes(:rankdir => "LR") -const node_attrs = Attributes(:shape => "plain", :style => "filled", :color => "white") -const edge_attrs = Attributes(:splines => "splines") - -function edgify(δ, i, reverse::Bool) - attr = Attributes() - return map(δ) do p - val = String(getname(p[1])) - weight = "$(p[2])" - attr = Attributes(:label => weight, :labelfontsize => "6") - return Edge(reverse ? ["rx_$i", "$val"] : - ["$val", "rx_$i"], attr) - end -end - -# make distinguished edge based on rate constant -function edgifyrates(rxs, specs) - es = Edge[] - for (i, rx) in enumerate(rxs) - deps = rx.rate isa Number ? Any[] : get_variables(rx.rate, specs) - for dep in deps - val = String(getname(dep)) - attr = Attributes(:color => "#d91111", :style => "dashed") - e = Edge(["$val", "rx_$i"], attr) - push!(es, e) - end - end - es -end - -# create edges from one reaction complex to another reaction complex -function edgifycomplex(δ, attr) - return map(δ) do p - return Edge([p[1], p[2]], attr) - end -end - -# modify vector of string of complexes into graphviz compatible strings -function modifystrcomp(strcomp::Vector{String}) - for i in 1:length(strcomp) - # cannot allow (t) as per graphviz - if occursin("(t)", strcomp[i]) - strcomp[i] = replace(strcomp[i], "(t)" => "") - end - # prettify NUll with ∅ - if occursin("0", strcomp[i]) - strcomp[i] = replace(strcomp[i], "0" => "∅") - end - end - strcomp = "<" .* strcomp .* ">" -end - -### Public-facing API ### - -""" - complexgraph(rn::ReactionSystem; complexdata=reactioncomplexes(rn)) - -Creates a Graphviz graph of the [`ReactionComplex`](@ref)s in `rn`. Reactions -correspond to arrows and reaction complexes to blue circles. - -Notes: -- Black arrows from complexes to complexes indicate reactions whose rate is a - parameter or a `Number`. i.e. `k, A --> B`. -- Red dashed arrows from complexes to complexes indicate reactions whose rate - depends on species. i.e. `k*C, A --> B` for `C` a species. -- Requires the Graphviz jll to be installed, or Graphviz to be installed and - commandline accessible. -""" -function complexgraph(rn::ReactionSystem; complexdata = reactioncomplexes(rn)) - rxs = reactions(rn) - specs = species(rn) - complexes, B = complexdata - fun = rcel -> specs[rcel.speciesid] * rcel.speciesstoich - compfun(rc) = rc == Catalyst.ReactionComplex{Int64}[] ? 0 : sum(fun, rc) - - strcomp = [string(compfun(rc)) for rc in complexes] - - newstrcomp = modifystrcomp(strcomp) - compnodes = [Node(str, Attributes(:shape => "circle", :color => "#6C9AC3")) - for str in newstrcomp] - - edges = [] - for (i, r) in enumerate(rxs) - subcomp = newstrcomp[argmin(@view B[:, i])] - prodcomp = newstrcomp[argmax(@view B[:, i])] - deps = get_variables(r.rate, specs) - if deps != Any[] - attr = Attributes(:color => "#d91111", :style => "dashed") - push!(edges, edgifycomplex(zip([subcomp], [prodcomp]), attr)) - else - attr = Attributes() - push!(edges, edgifycomplex(zip([subcomp], [prodcomp]), attr)) - end - end - stmts2 = Vector{Statement}() - append!(stmts2, compnodes) - append!(stmts2, collect(Iterators.flatten(edges))) - g = Digraph("G", stmts2; graph_attrs = graph_attrs, node_attrs = node_attrs, - edge_attrs = edge_attrs) - return g -end - -""" - Graph(rn::ReactionSystem) - -Converts a [`ReactionSystem`](@ref) into a Graphviz graph. -Reactions correspond to small green circles, and species to blue circles. - -Notes: -- Black arrows from species to reactions indicate reactants, and are labelled - with their input stoichiometry. -- Black arrows from reactions to species indicate products, and are labelled - with their output stoichiometry. -- Red arrows from species to reactions indicate that species is used within the - rate expression. For example, in the reaction `k*A, B --> C`, there would be a - red arrow from `A` to the reaction node. In `k*A, A+B --> C`, there would be - red and black arrows from `A` to the reaction node. -- Requires the Graphviz jll to be installed, or Graphviz to be installed and - commandline accessible. -""" -function Graph(rn::ReactionSystem) - rxs = reactions(rn) - specs = species(rn) - statenodes = [Node(string(getname(s)), - Attributes(:shape => "circle", :color => "#6C9AC3")) for s in specs] - transnodes = [Node(string("rx_$i"), - Attributes(:shape => "point", :color => "#E28F41", :width => ".1")) - for (i, r) in enumerate(rxs)] - - stmts = vcat(statenodes, transnodes) - edges = map(enumerate(rxs)) do (i, r) - vcat(edgify(zip(r.substrates, r.substoich), i, false), - edgify(zip(r.products, r.prodstoich), i, true)) - end - es = edgifyrates(rxs, specs) - (!isempty(es)) && push!(edges, es) - - stmts2 = Vector{Statement}() - append!(stmts2, stmts) - append!(stmts2, collect(Iterators.flatten(edges))) - g = Digraph("G", stmts2; graph_attrs = graph_attrs, node_attrs = node_attrs, - edge_attrs = edge_attrs) - return g -end - -""" - savegraph(g::Graph, fname, fmt="png") - -Given a `Graph` generated by [`Graph`](@ref), save the graph to the file with -name `fname` and extension `fmt`. - -Notes: -- `fmt="png"` is the default output format. -- Requires the Graphviz jll to be installed, or Graphviz to be installed and commandline accessible. -""" -function savegraph(g::Graph, fname, fmt = "png") - open(fname, "w") do io - run_graphviz(io, g, format = fmt) - end - nothing -end diff --git a/test/runtests.jl b/test/runtests.jl index 9b3016cd0c..29c2d0ca74 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -60,12 +60,6 @@ end @time @safetestset "ReactionSystem Serialisation" begin include("miscellaneous_tests/reactionsystem_serialisation.jl") end # BROKEN # @time @safetestset "Latexify" begin include("visualisation/latexify.jl") end - - # Tests network visualisation. - # Disable on Macs as can't install GraphViz via jll - if !Sys.isapple() - @time @safetestset "Graphs Visualisations" begin include("visualisation/graphs.jl") end - end end if GROUP == "All" || GROUP == "Spatial" diff --git a/test/visualisation/graphs.jl b/test/visualisation/graphs.jl deleted file mode 100644 index 2c0199a4c2..0000000000 --- a/test/visualisation/graphs.jl +++ /dev/null @@ -1,29 +0,0 @@ -# Fetch packages. -using Catalyst, Graphviz_jll - -### Basic Tests ### - -# Check basic functionality. -let - rn = @reaction_network begin - α, S + I --> 2I - β, I --> R - S^2, R --> 0 - end - - # Check can make a graph. - gr = Graph(rn) - - # Check can save a graph. - fname = Base.Filesystem.tempname() - savegraph(gr, "$fname.svg", "svg") - - rcgr = complexgraph(rn) - fname = Base.Filesystem.tempname() - savegraph(rcgr, "$fname.svg", "svg") -end - -# these are broken in the jll, see -# https://github.com/JuliaPackaging/Yggdrasil/issues/1428 -# savegraph(gr, "$fname.pdf", "pdf") -# savegraph(gr, "$fname.png", "png")