Skip to content

Commit

Permalink
Merge pull request #21 from isaacsas/updates
Browse files Browse the repository at this point in the history
updates
  • Loading branch information
isaacsas authored Aug 16, 2020
2 parents d99c430 + e2e712b commit d27ba3a
Show file tree
Hide file tree
Showing 4 changed files with 65 additions and 20 deletions.
26 changes: 26 additions & 0 deletions .github/workflows/CompatHelper.yml
Original file line number Diff line number Diff line change
@@ -0,0 +1,26 @@
name: CompatHelper

on:
schedule:
- cron: '00 * * * *'
issues:
types: [opened, reopened]

jobs:
build:
runs-on: ${{ matrix.os }}
strategy:
matrix:
julia-version: [1.3.0]
julia-arch: [x86]
os: [ubuntu-latest]
steps:
- uses: julia-actions/setup-julia@latest
with:
version: ${{ matrix.julia-version }}
- name: Pkg.add("CompatHelper")
run: julia -e 'using Pkg; Pkg.add("CompatHelper")'
- name: CompatHelper.main()
env:
GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }}
run: julia -e 'using CompatHelper; CompatHelper.main()'
41 changes: 30 additions & 11 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -3,20 +3,29 @@
[![Build Status](https://travis-ci.org/isaacsas/ReactionNetworkImporters.jl.svg?branch=master)](https://travis-ci.org/isaacsas/ReactionNetworkImporters.jl)
[![Build status](https://ci.appveyor.com/api/projects/status/wqq5flk2w8asad78/branch/master?svg=true)](https://ci.appveyor.com/project/isaacsas/reactionnetworkimporters-jl/branch/master)

*Note, v0.2.0 is a breaking release in that ReactionNetworkImporters now
generates a Catalyst.jl `ReactionSystem`. Do not upgrade if you rely on the
previous DiffEqBiological.jl reaction network representation.*

This package provides importers to load reaction networks into
[Catalyst.jl](https://github.com/SciML/Catalyst.jl)
[`ReactionSystem`s](https://catalyst.sciml.ai/dev/api/catalyst_api/#ModelingToolkit.ReactionSystem)
from several file formats. Currently it supports loading networks in the
following formats:
1. A *subset* of the BioNetGen .net file format.
2. Networks represented by dense or sparse substrate and product stoichiometric matrices.
2. Networks represented by dense or sparse substrate and product stoichiometric
matrices.
<!-- 3. The basic format used by the [RSSA](https://www.cosbi.eu/research/prototypes/rssa) group at COSBI in their [model collection](https://www.cosbi.eu/prototypes/jLiexDeBIgFV4zxwnKiW97oc4BjTtIoRGajqdUz4.zip). -->

----
## Examples

### Loading a BioNetGen .net file
A simple network from the builtin BioNetGen bngl examples is the [repressilator](data/repressilator/Repressilator.bngl). The `generate_network` command in the bngl file outputs a reduced network description, i.e. a [.net](data/repressilator/Repressilator.net) file, which can be loaded into a Catalyst `ReactionSystem` as:
A simple network from the builtin BioNetGen bngl examples is the
[repressilator](data/repressilator/Repressilator.bngl). The `generate_network`
command in the bngl file outputs a reduced network description, i.e. a
[.net](data/repressilator/Repressilator.net) file, which can be loaded into a
Catalyst `ReactionSystem` as:
```julia
using ReactionNetworkImporters
fname = "PATH/TO/Repressilator.net"
Expand All @@ -27,11 +36,19 @@ Here `BNGNetwork` is a type specifying the file format that is being loaded.
- `rn`, a Catalyst `reaction_network`
- `u₀`, the initial condition (as a `Vector{Float64}`)
- `p`, the parameter vector (as a `Vector{Float64}`)
- `paramexprs`, the parameter vector as a mix of `Numbers`, `Symbols` and `Exprs`. `p` is generated by evaluation of these expressions and symbols.
- `varstonames`, a `Dict` mapping from the internal `Symbol` of a species used in the generated `ReactionSystem` to a `Symbol` generated from the name in the .net file. This is necessary as BioNetGen can generate exceptionally long species names, involving characters that lead to malformed species names when used with `Catalyst`.
- `groupstoids`, a `Dict` mapping the `Symbol`s (i.e. names) for any species groups defined in the .net file to a vector of indices into `u₀` where the corresponding species are stored.
- `paramexprs`, the parameter vector as a mix of `Numbers`, `Symbols` and
`Exprs`. `p` is generated by evaluation of these expressions and symbols.
- `varstonames`, a `Dict` mapping from the internal `Symbol` of a species used
in the generated `ReactionSystem` to a `Symbol` generated from the name in the
.net file. This is necessary as BioNetGen can generate exceptionally long
species names, involving characters that lead to malformed species names when
used with `Catalyst`.
- `groupstoids`, a `Dict` mapping the `Symbol`s (i.e. names) for any species
groups defined in the .net file to a vector of indices into `u₀` where the
corresponding species are stored.

Given `prnbng`, we can construct and solve the corresponding ODE model for the reaction system by
Given `prnbng`, we can construct and solve the corresponding ODE model for the
reaction system by
```julia
using OrdinaryDiffEq, Catalyst
rn = prnbng.rn
Expand Down Expand Up @@ -85,8 +102,8 @@ prn = loadrxnetwork(MatrixNetwork(),
rateexprs::AbstractVector,
substoich::AbstractMatrix,
prodstoich::AbstractMatrix;
species=Symbol[],
params=Symbol[])
species::AbstractVector=Operation[],
params::AbstractVector=Operation[])
```
Here `MatrixNetwork()` is the dispatch type, which selects that we are
constructing a matrix-based stoichiometric representation as input. The other
Expand All @@ -102,8 +119,10 @@ parameters are:
- `prodstoich` - A number of species by number of reactions matrix with entry
`(i,j)` giving the stoichiometric coefficient of species `i` as a product in
reaction `j`.
- `species` - Optional `ModelingToolkit.Operation`s representing each species in the network.
- `parameters` - Optional `ModelingToolkit.Operation`s representing each parameter in the network.
- `species` - Optional `ModelingToolkit.Operation`s representing each species in
the network.
- `parameters` - Optional `ModelingToolkit.Operation`s representing each
parameter in the network.

`prn` is again a `ParsedReactionNetwork`, with only the `reaction_network`
field, `prn.rn`, defined.
Expand All @@ -114,7 +133,7 @@ A dispatch is added if `substoich` and `prodstoich` both have the type

If the keyword argument `species` is not set, the resulting reaction network
will simply name the species `S1`, `S2`,..., `SN` for a system with `N` total
species. `params` defaults to an empty vector of `Symbol`s, so that it does not
species. `params` defaults to an empty vector of `Operation`s, so that it does not
need to be set for systems with no parameters.

<!-- ### Loading a RSSA format network file
Expand Down
14 changes: 7 additions & 7 deletions examples/test_bcr_odes.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@ using ReactionNetworkImporters
using TimerOutputs

# parameters
doplot = true
doplot = false
networkname = "testbcrbng"
tf = 10000.
build_jac = false
Expand All @@ -31,7 +31,7 @@ show(to)
rnbng = prnbng.rn; u0 = prnbng.u₀; p = prnbng.p;
@timeit to "bODESystem" bosys = convert(ODESystem, rnbng)
show(to)
@timeit to "bODEProb" boprob = ODEProblem(bosys, Pair.(species(rnbng,u0), (0.,tf), Pair.(params(rnbng),p))
@timeit to "bODEProb" boprob = ODEProblem(bosys, Pair.(species(rnbng),u0), (0.,tf), Pair.(params(rnbng),p))
show(to)
u = copy(u0);
du = similar(u);
Expand All @@ -47,8 +47,8 @@ show(to)
println()

# BNG simulation results for Activated Syk
asykgroups = prnbng.groupstoids[!,:Activated_Syk]
asyksyms = findall(x -> x asykgroups, rnbng.syms_to_ints)
asykgroups = prnbng.groupstoids[:Activated_Syk]
#asyksyms = findall(x -> x ∈ asykgroups, rnbng.syms_to_ints)
# asynbng = zeros(length(gdatdf[:time]))
# for sym in asyksyms
# global asynbng
Expand All @@ -72,14 +72,14 @@ basyk = sum(bsol[asykgroups,:], dims=1)

if doplot
plotlyjs()
plot(gdatdf[!,:time][2:end], gdatdf[!,:Activated_Syk][2:end], xscale=:log10, label=:AsykGroup, linestyle=:dot)
plot(gdatdf[!,:time][2:end], gdatdf[!,:Activated_Syk][2:end], xscale=:log10, label="AsykGroup", linestyle=:dot)
# # plot!(cdatdf[:time][2:end], asynbng[2:end], xscale=:log10, label=:AsykSum)
plot!(bsol.t[2:end], basyk'[2:end], label=:AsykDEBio, xscale=:log10)
plot!(bsol.t[2:end], basyk'[2:end], label="AsykDEBio", xscale=:log10)
end

# test the error, note may be large in abs value though small relatively
# #norm(gdatdf[:Activated_Syk] - asynbng, Inf)
# #norm(asynbng - basyk', Inf)
norm(gdatdf[!,:Activated_Syk] - basyk', Inf)

@assert all(abs.(gdatdf[:Activated_Syk] - asynbng) .< 1e-6 * abs.(gdatdf[:Activated_Syk]))
#@assert all(abs.(gdatdf[!,:Activated_Syk] - basyk') .< 1e-6 * abs.(gdatdf[:Activated_Syk]))
4 changes: 2 additions & 2 deletions src/parsing_routines_matrixnetworks.jl
Original file line number Diff line number Diff line change
Expand Up @@ -70,8 +70,8 @@ function loadrxnetwork(ft::MatrixNetwork,
rateexprs::AbstractVector,
substoich::SparseMatrixCSC,
prodstoich::SparseMatrixCSC;
species=Symbol[],
params=Symbol[])
species::AbstractVector=Operation[],
params::AbstractVector=Operation[])

sz = size(substoich)
@assert sz == size(prodstoich)
Expand Down

0 comments on commit d27ba3a

Please sign in to comment.