Skip to content

Commit 786f6ea

Browse files
committed
Reduced graph sizes in default config
1 parent 2403287 commit 786f6ea

4 files changed

Lines changed: 90 additions & 54 deletions

File tree

paper.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -50,7 +50,7 @@ models and contact networks for particular initial conditions) is provided in th
5050
Our software answers an open question from [@kiss:2017] by facilitating algorithmic generation and solving of equations
5151
describing compartmental models with underlying contact networks. Others have written their own software to generate
5252
equations for systems up to a certain length of terms (usually terms on one and two vertices), which approximates the
53-
system dynamics - see, for example, [@sharkey:2011],. Instead, our software generates and solves equations up to the
53+
system dynamics - see, for example, [@sharkey:2011]. Instead, our software generates and solves equations up to the
5454
full system size, yielding full, complete and deterministic representations of model dynamics.
5555

5656
The more usual approach to obtaining modelling results for compartmental models on contact networks is to use a

src/equation/generation.py

Lines changed: 12 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -9,7 +9,7 @@
99
from equation.Term import Vertex, Term
1010
from equation.closing import can_be_closed, replace_with_closures
1111
from model_params.cmodel import CModel
12-
from model_params.helpers import dynamically_relevant, Coupling, coupling_types
12+
from model_params.helpers import dynamically_relevant, Coupling, coupling_types, parse_coupling_descriptor
1313

1414
t = sym.symbols('t')
1515

@@ -136,6 +136,8 @@ def add_terms(v: Vertex, term: Term, transition: tuple, model: CModel, neighbour
136136
neighbours_of_v_not_in_tuple = list(set(neighbours_of_v) - set(copy.deepcopy(term).node_list()))
137137
neighbours_of_v_in_tuple = list(set(copy.deepcopy(term).vertices) - {v})
138138
terms = []
139+
_, _, _, target_state = parse_coupling_descriptor(transition[1])
140+
target_state = target_state or v.state
139141
# Examples relate to vanilla SIR model
140142
if transition[0] == Coupling.NEIGHBOUR_ENTER:
141143
# E.G. I state requires something in S coming into contact with something in I
@@ -144,7 +146,7 @@ def add_terms(v: Vertex, term: Term, transition: tuple, model: CModel, neighbour
144146
elif transition[0] == Coupling.NEIGHBOUR_EXIT:
145147
# Converse of neighbour enter - S contacts I, transitions to I
146148
neighbour_exit(model, neighbours_of_v_in_tuple, neighbours_of_v_not_in_tuple, transition[2],
147-
sym.symbols(f'{transition[3]}'), copy.deepcopy(term), terms, transition, v, transition[1][-1])
149+
sym.symbols(f'{transition[3]}'), copy.deepcopy(term), terms, transition, v, target_state)
148150
elif transition[0] == Coupling.ISOLATED_ENTER:
149151
# E.G. I -> R without input of neighbours
150152
isolated_enter(transition[2], sym.symbols(f'{transition[3]}'), copy.deepcopy(term), terms, transition, v)
@@ -161,9 +163,10 @@ def find_neighbour_entries(model, neighbours_of_v_in_tuple, neighbours_of_v_not_
161163
# e.g. v is in state I, so change v to S and each neighbour in turn to I
162164
rate_of_transition = transition[2]
163165
symbol_for_rate_of_transition = sym.symbols(f'{transition[3]}')
164-
v_transitions_to = transition[1].split(':')[1][-1]
165-
v_starts_as = transition[1].split(':')[1][0] # the state v would be in before transitioning to current state
166-
other_state_for_neighbours = transition[1].split('*')[1][0]
166+
_, neighbour_state, source_state, target_state = parse_coupling_descriptor(transition[1])
167+
other_state_for_neighbours = neighbour_state or v.state
168+
v_starts_as = source_state or v.state # the state v would be in before transitioning to current state
169+
v_transitions_to = target_state or v.state
167170
# First, look at all external vertices that could lead to entering this state
168171
for n in neighbours_of_v_not_in_tuple:
169172
# Make sure new term contains all same terms as original
@@ -182,7 +185,8 @@ def find_neighbour_entries(model, neighbours_of_v_in_tuple, neighbours_of_v_not_
182185
def neighbour_exit(model, neighbours_of_v_in_tuple, neighbours_of_v_not_in_tuple, rate_of_transition,
183186
symbol_for_rate_of_transition, term_clone, terms, transition, v, v_transitions_to):
184187
# e.g. v in state S so can exit if any neighbour in I
185-
other_state_for_neighbours = transition[1].split('*')[1][0]
188+
_, neighbour_state, _, _ = parse_coupling_descriptor(transition[1])
189+
other_state_for_neighbours = neighbour_state or v.state
186190
for n in neighbours_of_v_not_in_tuple:
187191
vertices = set(term_clone.vertices).union({v, Vertex(other_state_for_neighbours, n)})
188192
term = Term(list(vertices)).function()(sym.symbols('t'))
@@ -203,7 +207,8 @@ def isolated_exit(rate_of_transition, symbol_for_rate_of_transition, term_clone,
203207

204208
def isolated_enter(rate_of_transition, symbol_for_rate_of_transition, term_clone, terms, transition, v):
205209
# e.g. v in state R, gets there through recovery after being in I
206-
other_state_for_v = transition[1].split(':')[1][0]
210+
_, _, source_state, _ = parse_coupling_descriptor(transition[1])
211+
other_state_for_v = source_state or v.state
207212
vertices = set(term_clone.vertices).union({Vertex(other_state_for_v, v.node)})
208213
term = Term(list(vertices)).function()(sym.symbols('t'))
209214
return append_term_to_terms(rate_of_transition, symbol_for_rate_of_transition, term, terms)

src/experiments/default_config.json

Lines changed: 21 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -2,7 +2,7 @@
22
"experiments": [
33
{
44
"graph_type": "erdos_renyi",
5-
"num_vertices": 100,
5+
"num_vertices": 50,
66
"graph_params": {"p": 0.05},
77
"method": "both",
88
"iterations": 12,
@@ -15,7 +15,7 @@
1515
},
1616
{
1717
"graph_type": "erdos_renyi",
18-
"num_vertices": 100,
18+
"num_vertices": 50,
1919
"graph_params": {"p": 0.1},
2020
"method": "both",
2121
"iterations": 12,
@@ -28,7 +28,7 @@
2828
},
2929
{
3030
"graph_type": "erdos_renyi",
31-
"num_vertices": 100,
31+
"num_vertices": 50,
3232
"graph_params": {"p": 0.15},
3333
"method": "both",
3434
"iterations": 12,
@@ -41,7 +41,7 @@
4141
},
4242
{
4343
"graph_type": "erdos_renyi",
44-
"num_vertices": 100,
44+
"num_vertices": 50,
4545
"graph_params": {"p": 0.2},
4646
"method": "both",
4747
"iterations": 12,
@@ -54,7 +54,7 @@
5454
},
5555
{
5656
"graph_type": "erdos_renyi",
57-
"num_vertices": 100,
57+
"num_vertices": 50,
5858
"graph_params": {"p": 0.25},
5959
"method": "both",
6060
"iterations": 12,
@@ -67,7 +67,7 @@
6767
},
6868
{
6969
"graph_type": "barabasi_albert",
70-
"num_vertices": 100,
70+
"num_vertices": 50,
7171
"graph_params": {"m": 2},
7272
"method": "both",
7373
"iterations": 12,
@@ -80,7 +80,7 @@
8080
},
8181
{
8282
"graph_type": "barabasi_albert",
83-
"num_vertices": 100,
83+
"num_vertices": 50,
8484
"graph_params": {"m": 3},
8585
"method": "both",
8686
"iterations": 12,
@@ -93,7 +93,7 @@
9393
},
9494
{
9595
"graph_type": "barabasi_albert",
96-
"num_vertices": 100,
96+
"num_vertices": 50,
9797
"graph_params": {"m": 4},
9898
"method": "both",
9999
"iterations": 12,
@@ -106,7 +106,7 @@
106106
},
107107
{
108108
"graph_type": "watts_strogatz",
109-
"num_vertices": 100,
109+
"num_vertices": 50,
110110
"graph_params": {"k": 4, "p": 0.1},
111111
"method": "both",
112112
"iterations": 12,
@@ -119,7 +119,7 @@
119119
},
120120
{
121121
"graph_type": "watts_strogatz",
122-
"num_vertices": 100,
122+
"num_vertices": 50,
123123
"graph_params": {"k": 8, "p": 0.1},
124124
"method": "both",
125125
"iterations": 12,
@@ -132,7 +132,7 @@
132132
},
133133
{
134134
"graph_type": "watts_strogatz",
135-
"num_vertices": 100,
135+
"num_vertices": 50,
136136
"graph_params": {"k": 12, "p": 0.1},
137137
"method": "both",
138138
"iterations": 12,
@@ -145,7 +145,7 @@
145145
},
146146
{
147147
"graph_type": "watts_strogatz",
148-
"num_vertices": 100,
148+
"num_vertices": 50,
149149
"graph_params": {"k": 16, "p": 0.1},
150150
"method": "both",
151151
"iterations": 12,
@@ -158,7 +158,7 @@
158158
},
159159
{
160160
"graph_type": "watts_strogatz",
161-
"num_vertices": 100,
161+
"num_vertices": 50,
162162
"graph_params": {"k": 20, "p": 0.1},
163163
"method": "both",
164164
"iterations": 12,
@@ -171,7 +171,7 @@
171171
},
172172
{
173173
"graph_type": "random_regular",
174-
"num_vertices": 100,
174+
"num_vertices": 50,
175175
"graph_params": {"d": 2},
176176
"method": "both",
177177
"iterations": 12,
@@ -184,7 +184,7 @@
184184
},
185185
{
186186
"graph_type": "random_regular",
187-
"num_vertices": 100,
187+
"num_vertices": 50,
188188
"graph_params": {"d": 3},
189189
"method": "both",
190190
"iterations": 12,
@@ -197,7 +197,7 @@
197197
},
198198
{
199199
"graph_type": "random_regular",
200-
"num_vertices": 100,
200+
"num_vertices": 50,
201201
"graph_params": {"d": 4},
202202
"method": "both",
203203
"iterations": 12,
@@ -210,7 +210,7 @@
210210
},
211211
{
212212
"graph_type": "random_regular",
213-
"num_vertices": 100,
213+
"num_vertices": 50,
214214
"graph_params": {"d": 5},
215215
"method": "both",
216216
"iterations": 12,
@@ -223,7 +223,7 @@
223223
},
224224
{
225225
"graph_type": "random_geometric",
226-
"num_vertices": 100,
226+
"num_vertices": 50,
227227
"graph_params": {"r": 0.2},
228228
"method": "both",
229229
"iterations": 12,
@@ -236,7 +236,7 @@
236236
},
237237
{
238238
"graph_type": "random_geometric",
239-
"num_vertices": 100,
239+
"num_vertices": 50,
240240
"graph_params": {"r": 0.4},
241241
"method": "both",
242242
"iterations": 12,
@@ -249,7 +249,7 @@
249249
},
250250
{
251251
"graph_type": "random_geometric",
252-
"num_vertices": 100,
252+
"num_vertices": 50,
253253
"graph_params": {"r": 0.6},
254254
"method": "both",
255255
"iterations": 12,
@@ -262,7 +262,7 @@
262262
},
263263
{
264264
"graph_type": "random_geometric",
265-
"num_vertices": 100,
265+
"num_vertices": 50,
266266
"graph_params": {"r": 0.8},
267267
"method": "both",
268268
"iterations": 12,

src/model_params/helpers.py

Lines changed: 56 additions & 25 deletions
Original file line numberDiff line numberDiff line change
@@ -8,32 +8,63 @@ class Coupling(Enum):
88
ISOLATED_EXIT = 3
99

1010

11+
def parse_coupling_descriptor(descriptor):
12+
"""Parse a coupling descriptor into its component states.
13+
14+
Returns a tuple ``(s1, s2, s3, s4)`` representing the states involved in
15+
the driving interaction (``s1`` and optional neighbour ``s2``) together with
16+
the source (``s3``) and destination (``s4``) states of the transition. When
17+
the descriptor omits any of these components they are inferred using the
18+
same defaults as :meth:`CModel.set_coupling_rate`.
19+
"""
20+
21+
transition_part, _, remainder = descriptor.partition(':')
22+
transition_part = transition_part.strip()
23+
remainder = remainder.strip()
24+
25+
lhs_states = [token for token in transition_part.split('*') if token]
26+
s1 = lhs_states[0] if lhs_states else None
27+
s2 = lhs_states[1] if len(lhs_states) > 1 else None
28+
29+
s3 = s4 = None
30+
if remainder:
31+
if '=>' in remainder:
32+
src_part, dst_part = remainder.split('=>', 1)
33+
src_part = src_part.strip()
34+
dst_part = dst_part.strip()
35+
s3 = src_part or None
36+
s4 = dst_part or None
37+
else:
38+
s4 = remainder or None
39+
40+
if s3 is None and s4 is None:
41+
# No explicit target provided; mirror the driving states.
42+
if s1 is not None and s2 is not None:
43+
s3, s4 = s1, s2
44+
elif s1 is not None:
45+
s3 = s4 = s1
46+
elif s3 is None and s4 is not None:
47+
# Target provided but no explicit source; default to the first driver.
48+
s3 = s1
49+
50+
return s1, s2, s3, s4
51+
52+
1153
def coupling_types(model):
12-
_coupling_map = {}
13-
# Initialise the dictionary keys (one for each model state)
14-
for state in model.states:
15-
_coupling_map[state] = []
16-
for state in model.states:
17-
for couple in model.couplings:
18-
# Get the situation under which a transition occurs
19-
transition = model.couplings[couple][0].split(':')[0]
20-
# Does the transition contain a state we are currently interested in?
21-
if model.couplings[couple][0].count(state) > 0:
22-
if transition[0] == state: # EXIT TRANSITION
23-
if transition.count('*') > 0: # NEEDS A NEIGHBOUR
24-
_coupling_map[state].append((Coupling.NEIGHBOUR_EXIT, model.couplings[couple][0],
25-
model.couplings[couple][1], couple))
26-
else:
27-
_coupling_map[state].append((Coupling.ISOLATED_EXIT, model.couplings[couple][0],
28-
model.couplings[couple][1], couple))
29-
else: # ENTRY TRANSITION
30-
if transition.count('*') > 0: # NEEDS A NEIGHBOUR
31-
if state in model.couplings[couple][0].split(':')[1]:
32-
_coupling_map[state].append((Coupling.NEIGHBOUR_ENTER, model.couplings[couple][0],
33-
model.couplings[couple][1], couple))
34-
else:
35-
_coupling_map[state].append((Coupling.ISOLATED_ENTER, model.couplings[couple][0],
36-
model.couplings[couple][1], couple))
54+
_coupling_map = {state: [] for state in model.states}
55+
56+
for couple, (descriptor, rate) in model.couplings.items():
57+
s1, s2, s3, s4 = parse_coupling_descriptor(descriptor)
58+
uses_neighbour = s2 is not None
59+
60+
if s3 in _coupling_map:
61+
exit_type = Coupling.NEIGHBOUR_EXIT if uses_neighbour else Coupling.ISOLATED_EXIT
62+
_coupling_map[s3].append((exit_type, descriptor, rate, couple))
63+
64+
if s4 in _coupling_map:
65+
entry_type = Coupling.NEIGHBOUR_ENTER if uses_neighbour else Coupling.ISOLATED_ENTER
66+
_coupling_map[s4].append((entry_type, descriptor, rate, couple))
67+
3768
return _coupling_map
3869

3970

0 commit comments

Comments
 (0)