Browse Source

Merge branch 'master' into jpf/bilayernetworks

pull/50/head
James 5 months ago
committed by GitHub
parent
commit
9ff654b1fa
No known key found for this signature in database GPG Key ID: 4AEE18F83AFDEB23
  1. 2
      docs/src/api.md
  2. 200
      src/AlgebraicPetri.jl
  3. 38
      src/Epidemiology.jl

2
docs/src/api.md

@ -1,5 +1,5 @@
# Library Reference
```@autodocs
Modules = [AlgebraicPetri]
Modules = [AlgebraicPetri, AlgebraicPetri.Epidemiology]
```

200
src/AlgebraicPetri.jl

@ -29,6 +29,10 @@ state_dict(n) = Dict(s=>i for (i, s) in enumerate(n))
# Petri Nets
############
""" ACSet definition for a Petri net.
See Catlab.jl documentation for description of the @present syntax.
"""
@present TheoryPetriNet(FreeSchema) begin
T::Ob
S::Ob
@ -45,12 +49,39 @@ const AbstractPetriNet = AbstractACSetType(TheoryPetriNet)
const PetriNet = CSetType(TheoryPetriNet,index=[:it,:is,:ot,:os])
const OpenPetriNetOb, OpenPetriNet = OpenCSetTypes(PetriNet,:S)
""" Open(p::AbstractPetriNet)
Converts a PetriNet to an OpenPetriNet where each state is exposed as a leg of
the cospan. The OpenPetriNet can be composed over an undirected wiring diagram
(see this
[blog post](https://www.algebraicjulia.org/blog/post/2020/11/structured-cospans-2/)
for a description of this compositional tooling)
"""
Open(p::AbstractPetriNet) = OpenPetriNet(p, map(x->FinFunction([x], ns(p)), 1:ns(p))...)
""" Open(p::AbstractPetriNet, legs...)
Generates on OpenPetriNet with legs bundled as described by `legs`
"""
Open(p::AbstractPetriNet, legs...) = OpenPetriNet(p, map(l->FinFunction(l, ns(p)), legs)...)
""" Open(n, p::AbstractPetriNet, m)
Generates on OpenPetriNet with two legs, `n` and `m`
"""
Open(n, p::AbstractPetriNet, m) = Open(p, n, m)
# PetriNet([:S, :I, :R], :infection=>((1, 2), 3))
""" PetriNet(n::Int, ts::Vararg{Union{Pair,Tuple}})
Constructs a PetriNet object with `n` states and transitions described by `ts`.
Transitions are given as `(input_states)=>(output_states)`.
A PetriNet modelling the SIR model with 3 states and 2 transitions can be
constructed as follows:
```@example
PetriNet(3, (1,2)=>(2,2), 2=>3)
```
"""
PetriNet(n::Int, ts::Vararg{Union{Pair,Tuple}}) = begin
p = PetriNet()
add_species!(p, n)
@ -64,34 +95,122 @@ PetriNet(n::Int, ts::Vararg{Union{Pair,Tuple}}) = begin
p
end
""" Number of states in a Petri net
"""
ns(p::AbstractPetriNet) = nparts(p,:S)
""" Number of transitions in a Petri net
"""
nt(p::AbstractPetriNet) = nparts(p,:T)
""" Number of input relationships in a Petri net
"""
ni(p::AbstractPetriNet) = nparts(p,:I)
""" Number of output relationships in a Petri net
"""
no(p::AbstractPetriNet) = nparts(p,:O)
""" Add a species to the Petri net. Label and concentration can be provided
depending on the kind of Petri net.
Returns the ID of the species
"""
add_species!(p::AbstractPetriNet;kw...) = add_part!(p,:S;kw...)
""" Add `n` species to the Petri net. Label and concentration can be provided
depending on the kind of Petri net.
Returns the ID of the species
"""
add_species!(p::AbstractPetriNet,n;kw...) = add_parts!(p,:S,n;kw...)
""" Add a transition to the Petri net. Label and rate can be provided
depending on the kind of Petri net.
Returns the ID of the transition
"""
add_transition!(p::AbstractPetriNet;kw...) = add_part!(p,:T;kw...)
""" Add `n` transitions to the Petri net. Label and rate can be provided
depending on the kind of Petri net.
Returns the ID of the transition
"""
add_transitions!(p::AbstractPetriNet,n;kw...) = add_parts!(p,:T,n;kw...)
""" add_input!(p::AbstractPetriNet,t,s;kw...)
Add an input relationship to the Petri net between the transition `t` and species `s`.
Returns the ID of the input relationship
"""
add_input!(p::AbstractPetriNet,t,s;kw...) = add_part!(p,:I;it=t,is=s,kw...)
""" add_inputs!(p::AbstractPetriNet,n,t,s;kw...)
Add input relationships to the Petri net between the transitions `t` and species `s`.
Returns the ID of the input relationship
"""
add_inputs!(p::AbstractPetriNet,n,t,s;kw...) = add_parts!(p,:I,n;it=t,is=s,kw...)
""" add_output!(p::AbstractPetriNet,t,s;kw...)
Add an output relationship to the Petri net between the transition `t` and species `s`.
Returns the ID of the input relationship
"""
add_output!(p::AbstractPetriNet,t,s;kw...) = add_part!(p,:O;ot=t,os=s,kw...)
""" add_outputs!(p::AbstractPetriNet,n,t,s;kw...)
Add output relationships to the Petri net between the transitions `t` and species `s`.
Returns the ID of the input relationship
"""
add_outputs!(p::AbstractPetriNet,n,t,s;kw...) = add_parts!(p,:O,n;ot=t,os=s,kw...)
""" Name of species
Note that this returns an index if labels are not present in the PetriNet
"""
sname(p::AbstractPetriNet,s) = (1:ns(p))[s]
""" Name of transition
Note that this returns an index if labels are not present in the PetriNet
"""
tname(p::AbstractPetriNet,t) = (1:nt(p))[t]
""" Names of species in a Petri net
Note that this returns indices if labels are not present in the PetriNet
"""
snames(p::AbstractPetriNet) = map(s->sname(p, s), 1:ns(p))
""" Names of transitions in a Petri net
Note that this returns indices if labels are not present in the PetriNet
"""
tnames(p::AbstractPetriNet) = map(t->tname(p, t), 1:nt(p))
# Note: although indexing makes this pretty fast, it is often faster to bulk-convert
# the PetriNet net into a transition matrix, if you are working with all of the transitions
""" Input relationships for a transition
"""
inputs(p::AbstractPetriNet,t) = subpart(p,incident(p,t,:it),:is)
""" Output relationships for a transition
"""
outputs(p::AbstractPetriNet,t) = subpart(p,incident(p,t,:ot),:os)
""" TransitionMatrices
This data structure stores the transition matrix of an AbstractPetriNet object.
This is primarily used for constructing the vectorfield representation of the
Petri net.
"""
struct TransitionMatrices
input::Matrix{Int}
output::Matrix{Int}
@ -107,6 +226,10 @@ struct TransitionMatrices
end
end
""" PetriNet(tm::TransitionMatrices)
Constructs a PetriNet from its TransitionMatrices representation.
"""
PetriNet(tm::TransitionMatrices) = begin
(m,n) = size(tm.input)
p = PetriNet()
@ -124,6 +247,14 @@ end
valueat(x::Number, u, t) = x
valueat(f::Function, u, t) = try f(u,t) catch e f(t) end
""" vectorfield(pn::AbstractPetriNet)
Generates a Julia function which calculates the vectorfield of the Petri net
being simulated under the law of mass action.
The resulting function has a signature of the form `f!(du, u, p, t)` and can be
passed to the DifferentialEquations.jl solver package.
"""
vectorfield(pn::AbstractPetriNet) = begin
tm = TransitionMatrices(pn)
dt = tm.output - tm.input
@ -142,6 +273,10 @@ vectorfield(pn::AbstractPetriNet) = begin
return f
end
""" ACSet definition for a Petri net with labels on transitions and states.
See Catlab.jl documentation for description of the @present syntax.
"""
@present TheoryLabelledPetriNet <: TheoryPetriNet begin
Name::Data
@ -155,6 +290,7 @@ const LabelledPetriNet = LabelledPetriNetUntyped{Symbol}
const OpenLabelledPetriNetObUntyped, OpenLabelledPetriNetUntyped = OpenACSetTypes(LabelledPetriNetUntyped,:S)
const OpenLabelledPetriNetOb, OpenLabelledPetriNet = OpenLabelledPetriNetObUntyped{Symbol}, OpenLabelledPetriNetUntyped{Symbol}
Open(p::AbstractLabelledPetriNet) = OpenLabelledPetriNet(p, map(x->FinFunction([x], ns(p)), 1:ns(p))...)
Open(p::AbstractLabelledPetriNet, legs...) = begin
s_idx = Dict(sname(p, s)=>s for s in 1:ns(p))
@ -162,6 +298,18 @@ Open(p::AbstractLabelledPetriNet, legs...) = begin
end
Open(n, p::AbstractLabelledPetriNet, m) = Open(p, n, m)
""" LabelledPetriNet(n, ts::Vararg{Union{Pair,Tuple}})
Constructs a LabelledPetriNet object with state names as elements of `n` and
labelled transitions described by `ts`.
Transitions are given as `transition_name=>((input_states)=>(output_states))`.
A LabelledPetriNet modelling the SIR model with 3 states and 2 transitions can be
constructed as follows:
```@example
LabelledPetriNet([:S, :I, :R], :inf=>((:S,:I)=>(:I,:I)), :rec=>(:I=>:R))
```
"""
LabelledPetriNet(n, ts::Vararg{Union{Pair,Tuple}}) = begin
p = LabelledPetriNet()
n = vectorify(n)
@ -180,6 +328,11 @@ end
# Reaction Nets
###############
""" ACSet definition for a Petri net with rates on transitions and
concentrations on states.
See Catlab.jl documentation for description of the @present syntax.
"""
@present TheoryReactionNet <: TheoryPetriNet begin
Rate::Data
Concentration::Data
@ -196,6 +349,22 @@ Open(p::AbstractReactionNet{R,C}) where {R,C} = OpenReactionNet{R,C}(p, map(x->F
Open(p::AbstractReactionNet{R,C}, legs...) where {R,C} = OpenReactionNet{R,C}(p, map(l->FinFunction(l, ns(p)), legs)...)
Open(n, p::AbstractReactionNet, m) = Open(p, n, m)
""" ReactionNet{R,C}(n, ts::Vararg{Union{Pair,Tuple}}) where {R,C}
Constructs a ReactionNet object with state concentrations as elements of `n`
and transitions described by `ts`. `R` is the data type used to store rates and
`C` is the data type used to store concentrations.
Transitions are given as `transition_rate=>((input_states)=>(output_states))`.
A ReactionNet modelling the SIR model with 3 states and 2 transitions, an
initial population of 10 susceptible, 1 infected, 0 recovered and an infection
rate of 0.5 and recovery rate of 0.1 can be
constructed as follows:
```@example
ReactionNet{Float64, Float64}([10,1,0], 0.5=>((1,2)=>(2,2)), 0.1=>(2=>3))
```
"""
ReactionNet{R,C}(n, ts::Vararg{Union{Pair,Tuple}}) where {R,C} = begin
p = ReactionNet{R,C}()
add_species!(p, length(n), concentration=n)
@ -215,6 +384,10 @@ rate(p::AbstractReactionNet,t) = subpart(p,t,:rate)
concentrations(p::AbstractReactionNet) = map(s->concentration(p, s), 1:ns(p))
rates(p::AbstractReactionNet) = map(t->rate(p, t), 1:nt(p))
""" ACSet definition for a ReactionNet with labels on transitions and states.
See Catlab.jl documentation for description of the @present syntax.
"""
@present TheoryLabelledReactionNet <: TheoryReactionNet begin
Name::Data
@ -237,6 +410,22 @@ end
Open(n, p::AbstractLabelledReactionNet, m) = Open(p, n, m)
# Ex. LabelledReactionNet{Number, Int}((:S=>990, :I=>10, :R=>0), (:inf, .3/1000)=>((:S, :I)=>(:I,:I)), (:rec, .2)=>(:I=>:R))
""" LabelledReactionNet{R,C}(n, ts::Vararg{Union{Pair,Tuple}}) where {R,C}
Constructs a LabelledReactionNet object with labelled state concentrations as elements of `n`
and labelled transitions described by `ts`. `R` is the data type used to store rates and
`C` is the data type used to store concentrations.
Transitions are given as `(t_name=>t_rate)=>((input_states)=>(output_states))`.
A LabelledReactionNet modelling the SIR model with 3 states and 2 transitions, an
initial population of 10 susceptible, 1 infected, 0 recovered and an infection
rate of 0.5 and recovery rate of 0.1 can be
constructed as follows:
```@example
ReactionNet{Float64, Float64}([:S=>10,:I=>1,:R=>0], (:inf=>0.5)=>((1,2)=>(2,2)), (:rec=>0.1)=>(2=>3))
```
"""
LabelledReactionNet{R,C}(n, ts::Vararg{Union{Pair,Tuple}}) where {R,C} = begin
p = LabelledReactionNet{R,C}()
n = vectorify(n)
@ -323,14 +512,23 @@ LabelledReactionNet{R,C}(pn::Union{AbstractPetriNet}, states, transitions) where
pn′
end
""" Concentration of a ReactionNet
"""
concentration(p::AbstractLabelledReactionNet,s) = subpart(p,s,:concentration)
""" Rate of a RectionNet
"""
rate(p::AbstractLabelledReactionNet,t) = subpart(p,t,:rate)
""" All concentrations of a ReactionNet
"""
concentrations(p::AbstractLabelledReactionNet) = begin
snames = [sname(p, s) for s in 1:ns(p)]
LVector(;[(snames[s]=>concentration(p, s)) for s in 1:ns(p)]...)
end
""" All rates of a ReactionNet
"""
rates(p::AbstractLabelledReactionNet) = begin
tnames = [tname(p, s) for s in 1:nt(p)]
LVector(;[(tnames[t]=>rate(p, t)) for t in 1:nt(p)]...)

38
src/Epidemiology.jl

@ -1,3 +1,6 @@
""" Specific generators and useful tools for constructing epidemiological
systems
"""
module Epidemiology
using AlgebraicPetri
@ -8,19 +11,52 @@ using Catlab.CategoricalAlgebra.FinSets
export oapply_epi, infection, exposure, illness, recovery, death
#spontaneous_petri(x::Symbol, y::Symbol, z::Symbol)
#
#Generates an OpenLabelledPetriNet which represents a transition named `z` which
#changes a token from state `x` to state `y`.
spontaneous_petri(x::Symbol, y::Symbol, z::Symbol) = Open(LabelledPetriNet(unique([x,y]), z=>(x, y)))
#exposure_petri(x::Symbol, y::Symbol, z::Symbol, transition::Symbol)
#
#Generates an OpenLabelledPetriNet which represents an exposure transition where
#a token from state `y` "exposes" a token from state `x`, converting the token
#from state `x` to state `z`.
exposure_petri(x::Symbol, y::Symbol, z::Symbol, transition::Symbol) =
Open(LabelledPetriNet(unique([x,y,z]), transition=>((x,y)=>(z,y))))
""" LabelledPetriNet which describes the infection process of tokens in state S
by tokens in state I
"""
infection = exposure_petri(:S, :I, :I, :inf)
""" LabelledPetriNet which describes the exposure process where tokens in I
"expose" tokens in S, changing them from S to E
"""
exposure = exposure_petri(:S, :I, :E, :exp)
""" LabelledPetriNet which describes the illness process which moves tokens from E to I.
"""
illness = spontaneous_petri(:E,:I,:ill)
""" LabelledPetriNet which describes the recovery process which moves tokens from I to R
"""
recovery = spontaneous_petri(:I,:R,:rec)
""" LabelledPetriNet which describes the death process which moves tokens from I to D
"""
death = spontaneous_petri(:I,:D,:death)
epi_dict = Dict(:infection=>infection, :exposure=>exposure, :illness=>illness, :recovery=>recovery, :death=>death)
""" oapply_epi
""" oapply_epi(ex, args...)
Generates a LabelledPetriNet under a composition pattern described by the
undirected wiring diagram `ex`. This requires that the nodes in `ex` are only
labelled with labels from the following set:
```
[:infection, :exposure, :illness, :recovery, :death]
```
"""
oapply_epi(ex, args...) = oapply(ex, epi_dict, args...)

Loading…
Cancel
Save