Skip to content

Commit 53dba72

Browse files
gcappongcappon
committed
v 1.1.0
1 parent 15ac716 commit 53dba72

16 files changed

+167
-54
lines changed

‎dist/py_replay_bg-1.0.5-py3-none-any.whl‎

102 KB
Binary file not shown.

‎dist/py_replay_bg-1.0.5.tar.gz‎

75.1 KB
Binary file not shown.

‎dist/py_replay_bg-1.1.0-py3-none-any.whl‎

104 KB
Binary file not shown.

‎dist/py_replay_bg-1.1.0.tar.gz‎

75.2 KB
Binary file not shown.

‎py_replay_bg.egg-info/PKG-INFO‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
Metadata-Version: 2.2
22
Name: py_replay_bg
3-
Version: 1.0.4
3+
Version: 1.1.0
44
Summary: ReplayBG is a digital twin-based methodology to assess new strategies for type 1 diabetes management.
55
Author-email: Giacomo Cappon <cappongiacomo@gmail.com>
66
Project-URL: Homepage, https://github.com/gcappon/py_replay_bg

‎py_replay_bg.egg-info/SOURCES.txt‎

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -26,7 +26,9 @@ py_replay_bg/example/code/replay_mcmc_dss.py
2626
py_replay_bg/example/code/twin_intervals_map.py
2727
py_replay_bg/example/code/twin_intervals_mcmc.py
2828
py_replay_bg/example/code/twin_map.py
29+
py_replay_bg/example/code/twin_map_extended.py
2930
py_replay_bg/example/code/twin_mcmc.py
31+
py_replay_bg/example/code/twin_mcmc_extended.py
3032
py_replay_bg/example/code/utils.py
3133
py_replay_bg/input_validation/__init__.py
3234
py_replay_bg/input_validation/input_validator_init.py
@@ -47,12 +49,16 @@ py_replay_bg/tests/test_replay_intervals_map.py
4749
py_replay_bg/tests/test_replay_intervals_mcmc.py
4850
py_replay_bg/tests/test_replay_map.py
4951
py_replay_bg/tests/test_replay_map_dss.py
52+
py_replay_bg/tests/test_replay_map_extended.py
5053
py_replay_bg/tests/test_replay_mcmc.py
5154
py_replay_bg/tests/test_replay_mcmc_dss.py
55+
py_replay_bg/tests/test_replay_mcmc_extended.py
5256
py_replay_bg/tests/test_twin_intervals_map.py
5357
py_replay_bg/tests/test_twin_intervals_mcmc.py
5458
py_replay_bg/tests/test_twin_map.py
59+
py_replay_bg/tests/test_twin_map_extended.py
5560
py_replay_bg/tests/test_twin_mcmc.py
61+
py_replay_bg/tests/test_twin_mcmc_extended.py
5662
py_replay_bg/twinning/__init__.py
5763
py_replay_bg/twinning/map.py
5864
py_replay_bg/twinning/mcmc.py

‎py_replay_bg/model/logpriors_t1d.py‎

Lines changed: 18 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -39,10 +39,10 @@ def log_prior_single_meal(
3939
"""
4040

4141
# unpack the model parameters
42-
Gb, SG, ka2, kd, kempt, SI, kabs, beta = theta
43-
42+
Gb, SG, p2, ka2, kd, kempt, SI, kabs, beta = theta
4443
# compute each log prior
4544
logprior_SI = log_gamma(SI * VG, 3.3, 1 / 5e-4)
45+
logprior_p2 = log_norm(np.sqrt(p2), mu=0.11, sigma=0.004) if 0 < p2 < 1 else -np.inf
4646
logprior_Gb = log_norm(Gb, mu=119.13, sigma=7.11) if 70 <= Gb <= 180 else -np.inf
4747
logprior_SG = log_lognorm(SG, mu=-3.8, sigma=0.5) if 0 < SG < 1 else -np.inf
4848
logprior_ka2 = log_lognorm(ka2, mu=-4.2875, sigma=0.4274) if 0 < ka2 < kd and ka2 < 1 else -np.inf
@@ -57,6 +57,7 @@ def log_prior_single_meal(
5757
return (logprior_SI +
5858
logprior_Gb +
5959
logprior_SG +
60+
logprior_p2 +
6061
logprior_ka2 +
6162
logprior_kd +
6263
logprior_kempt +
@@ -97,11 +98,10 @@ def log_prior_single_meal_exercise(
9798
"""
9899

99100
# unpack the model parameters
100-
Gb, SG, ka2, kd, kempt, SI, kabs, beta, e1, e2 = theta
101-
101+
Gb, SG, p2, ka2, kd, kempt, SI, kabs, beta, e1, e2 = theta
102102
# compute each log prior
103103
logprior_SI = log_gamma(SI * VG, 3.3, 1 / 5e-4)
104-
104+
logprior_p2 = log_norm(np.sqrt(p2), mu=0.11, sigma=0.004) if 0 < p2 < 1 else -np.inf
105105
logprior_Gb = log_norm(Gb, mu=119.13, sigma=7.11) if 70 <= Gb <= 180 else -np.inf
106106
logprior_SG = log_lognorm(SG, mu=-3.8, sigma=0.5) if 0 < SG < 1 else -np.inf
107107
logprior_ka2 = log_lognorm(ka2, mu=-4.2875, sigma=0.4274) if 0 < ka2 < kd and ka2 < 1 else -np.inf
@@ -120,6 +120,7 @@ def log_prior_single_meal_exercise(
120120
return (logprior_SI +
121121
logprior_Gb +
122122
logprior_SG +
123+
logprior_p2 +
123124
logprior_ka2 +
124125
logprior_kd +
125126
logprior_kempt +
@@ -235,9 +236,9 @@ def log_prior_multi_meal(
235236
"""
236237

237238
# unpack the model parameters
238-
# SI, Gb, SG, p2, ka2, kd, kempt, kabs, beta = theta
239+
# Gb, SG, p2, ka2, kd, kempt, kabs, beta = theta
239240

240-
Gb, SG, ka2, kd, kempt = theta[0:5]
241+
Gb, SG, p2, ka2, kd, kempt = theta[0:6]
241242

242243
SI_B = theta[pos_SI_B] if pos_SI_B else SI_B
243244
SI_L = theta[pos_SI_L] if pos_SI_L else SI_L
@@ -262,7 +263,8 @@ def log_prior_multi_meal(
262263

263264
logprior_Gb = log_norm(Gb, mu=119.13, sigma=7.11) if 70 <= Gb <= 180 else -np.inf
264265
logprior_SG = log_lognorm(SG, mu=-3.8, sigma=0.5) if 0 < SG < 1 else -np.inf
265-
# logprior_p2 = np.log(stats.norm.pdf(np.sqrt(p2), 0.11, 0.004)) if 0 < p2 < 1 else -np.inf
266+
267+
logprior_p2 = log_norm(np.sqrt(p2), mu=0.11, sigma=0.004) if 0 < p2 < 1 else -np.inf
266268
logprior_ka2 = log_lognorm(ka2, mu=-4.2875, sigma=0.4274) if 0 < ka2 < kd and ka2 < 1 else -np.inf
267269
logprior_kd = log_lognorm(kd, mu=-3.5090, sigma=0.6187) if 0 < ka2 < kd and kd < 1 else -np.inf
268270
logprior_kempt = log_lognorm(kempt, mu=-1.9646, sigma=0.7069) if 0 < kempt < 1 else -np.inf
@@ -289,6 +291,7 @@ def log_prior_multi_meal(
289291
logprior_SI_D +
290292
logprior_Gb +
291293
logprior_SG +
294+
logprior_p2 +
292295
logprior_ka2 +
293296
logprior_kd +
294297
logprior_kempt +
@@ -423,7 +426,7 @@ def log_prior_multi_meal_exercise(
423426
# unpack the model parameters
424427
# SI, Gb, SG, p2, ka2, kd, kempt, kabs, beta = theta
425428

426-
Gb, SG, ka2, kd, kempt = theta[0:5]
429+
Gb, SG, p2, ka2, kd, kempt = theta[0:6]
427430

428431
SI_B = theta[pos_SI_B] if pos_SI_B else SI_B
429432
SI_L = theta[pos_SI_L] if pos_SI_L else SI_L
@@ -451,7 +454,7 @@ def log_prior_multi_meal_exercise(
451454

452455
logprior_Gb = log_norm(Gb, mu=119.13, sigma=7.11) if 70 <= Gb <= 180 else -np.inf
453456
logprior_SG = log_lognorm(SG, mu=-3.8, sigma=0.5) if 0 < SG < 1 else -np.inf
454-
# logprior_p2 = np.log(stats.norm.pdf(np.sqrt(p2), 0.11, 0.004)) if 0 < p2 < 1 else -np.inf
457+
logprior_p2 = log_norm(np.sqrt(p2), mu=0.11, sigma=0.004) if 0 < p2 < 1 else -np.inf
455458
logprior_ka2 = log_lognorm(ka2, mu=-4.2875, sigma=0.4274) if 0 < ka2 < kd and ka2 < 1 else -np.inf
456459
logprior_kd = log_lognorm(kd, mu=-3.5090, sigma=0.6187) if 0 < ka2 < kd and kd < 1 else -np.inf
457460
logprior_kempt = log_lognorm(kempt, mu=-1.9646, sigma=0.7069) if 0 < kempt < 1 else -np.inf
@@ -481,6 +484,7 @@ def log_prior_multi_meal_exercise(
481484
logprior_SI_D +
482485
logprior_Gb +
483486
logprior_SG +
487+
logprior_p2 +
484488
logprior_ka2 +
485489
logprior_kd +
486490
logprior_kempt +
@@ -645,9 +649,9 @@ def log_prior_multi_meal_extended(
645649
"""
646650

647651
# unpack the model parameters
648-
# SI, Gb, SG, p2, ka2, kd, kempt, kabs, beta = theta
652+
# Gb, SG, p2, ka2, kd, kempt, kabs, beta = theta
649653

650-
Gb, SG, ka2, kd, kempt = theta[0:5]
654+
Gb, SG, p2, ka2, kd, kempt = theta[0:6]
651655

652656
SI_B = theta[pos_SI_B] if pos_SI_B else SI_B
653657
SI_L = theta[pos_SI_L] if pos_SI_L else SI_L
@@ -683,7 +687,7 @@ def log_prior_multi_meal_extended(
683687

684688
logprior_Gb = log_norm(Gb, mu=119.13, sigma=7.11) if 70 <= Gb <= 180 else -np.inf
685689
logprior_SG = log_lognorm(SG, mu=-3.8, sigma=0.5) if 0 < SG < 1 else -np.inf
686-
# logprior_p2 = np.log(stats.norm.pdf(np.sqrt(p2), 0.11, 0.004)) if 0 < p2 < 1 else -np.inf
690+
logprior_p2 = log_norm(np.sqrt(p2), mu=0.11, sigma=0.004) if 0 < p2 < 1 else -np.inf
687691
logprior_ka2 = log_lognorm(ka2, mu=-4.2875, sigma=0.4274) if 0 < ka2 < kd and ka2 < 1 else -np.inf
688692
logprior_kd = log_lognorm(kd, mu=-3.5090, sigma=0.6187) if 0 < ka2 < kd and kd < 1 else -np.inf
689693
logprior_kempt = log_lognorm(kempt, mu=-1.9646, sigma=0.7069) if 0 < kempt < 1 else -np.inf
@@ -721,6 +725,7 @@ def log_prior_multi_meal_extended(
721725
logprior_SI_D +
722726
logprior_Gb +
723727
logprior_SG +
728+
logprior_p2 +
724729
logprior_ka2 +
725730
logprior_kd +
726731
logprior_kempt +

‎py_replay_bg/model/model_step_equations_t1d.py‎

Lines changed: 15 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -142,10 +142,10 @@ def model_step_equations_single_meal(A, I, cho, hour_of_the_day, xkm1, B,
142142
risk = 1
143143

144144
# Compute the risk
145-
if 119.13 > xkm1[0] >= 60:
146-
risk = risk + 10 * r1 * (np.log(xkm1[0]) ** r2 - np.log(119.13) ** r2) ** 2
145+
if Gb > xkm1[0] >= 60:
146+
risk = risk + 10 * r1 * (np.log(xkm1[0]) ** r2 - np.log(Gb) ** r2) ** 2
147147
elif xkm1[0] < 60:
148-
risk = risk + 10 * r1 * (np.log(60) ** r2 - np.log(119.13) ** r2) ** 2
148+
risk = risk + 10 * r1 * (np.log(60) ** r2 - np.log(Gb) ** r2) ** 2
149149

150150
# Compute the model state at time k using backward Euler method
151151
B[:] = [cho / (1 + kgri), 0, 0,
@@ -176,10 +176,10 @@ def model_step_equations_single_meal_exercise(A, I, cho, vo2, hour_of_the_day, x
176176
risk = 1
177177

178178
# Compute the risk
179-
if 119.13 > xkm1[0] >= 60:
180-
risk = risk + 10 * r1 * (np.log(xkm1[0]) ** r2 - np.log(119.13) ** r2) ** 2
179+
if Gb > xkm1[0] >= 60:
180+
risk = risk + 10 * r1 * (np.log(xkm1[0]) ** r2 - np.log(Gb) ** r2) ** 2
181181
elif xkm1[0] < 60:
182-
risk = risk + 10 * r1 * (np.log(60) ** r2 - np.log(119.13) ** r2) ** 2
182+
risk = risk + 10 * r1 * (np.log(60) ** r2 - np.log(Gb) ** r2) ** 2
183183

184184
if vo2 == 0:
185185
xk[8] = 0
@@ -227,10 +227,10 @@ def model_step_equations_multi_meal(A, I, cho_b, cho_l, cho_d, cho_s, cho_h, hou
227227
risk = 1
228228

229229
# Compute the risk
230-
if 119.13 > xkm1[0] >= 60:
231-
risk = risk + 10 * r1 * (np.log(xkm1[0]) ** r2 - np.log(119.13) ** r2) ** 2
230+
if Gb > xkm1[0] >= 60:
231+
risk = risk + 10 * r1 * (np.log(xkm1[0]) ** r2 - np.log(Gb) ** r2) ** 2
232232
elif xkm1[0] < 60:
233-
risk = risk + 10 * r1 * (np.log(60) ** r2 - np.log(119.13) ** r2) ** 2
233+
risk = risk + 10 * r1 * (np.log(60) ** r2 - np.log(Gb) ** r2) ** 2
234234

235235
# Compute the model state at time k using backward Euler method
236236

@@ -276,10 +276,10 @@ def model_step_equations_multi_meal_extended(A, I, cho_b, cho_l, cho_d, cho_s, c
276276
risk = 1
277277

278278
# Compute the risk
279-
if 119.13 > xkm1[0] >= 60:
280-
risk = risk + 10 * r1 * (np.log(xkm1[0]) ** r2 - np.log(119.13) ** r2) ** 2
279+
if Gb > xkm1[0] >= 60:
280+
risk = risk + 10 * r1 * (np.log(xkm1[0]) ** r2 - np.log(Gb) ** r2) ** 2
281281
elif xkm1[0] < 60:
282-
risk = risk + 10 * r1 * (np.log(60) ** r2 - np.log(119.13) ** r2) ** 2
282+
risk = risk + 10 * r1 * (np.log(60) ** r2 - np.log(Gb) ** r2) ** 2
283283

284284
# Compute the model state at time k using backward Euler method
285285

@@ -330,10 +330,10 @@ def model_step_equations_multi_meal_exercise(A, I, cho_b, cho_l, cho_d, cho_s, c
330330
risk = 1
331331

332332
# Compute the risk
333-
if 119.13 > xkm1[0] >= 60:
334-
risk = risk + 10 * r1 * (np.log(xkm1[0]) ** r2 - np.log(119.13) ** r2) ** 2
333+
if Gb > xkm1[0] >= 60:
334+
risk = risk + 10 * r1 * (np.log(xkm1[0]) ** r2 - np.log(Gb) ** r2) ** 2
335335
elif xkm1[0] < 60:
336-
risk = risk + 10 * r1 * (np.log(60) ** r2 - np.log(119.13) ** r2) ** 2
336+
risk = risk + 10 * r1 * (np.log(60) ** r2 - np.log(Gb) ** r2) ** 2
337337

338338
if vo2 == 0:
339339
xk[20] = 0

‎py_replay_bg/model/t1d_model_multi_meal.py‎

Lines changed: 9 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -138,15 +138,15 @@ def __init__(self,
138138
self.model_parameters = ModelParametersT1DMultiMeal(data, bw, u2ss, extended)
139139

140140
# Unknown parameters
141-
self.unknown_parameters = ['Gb', 'SG', 'ka2', 'kd', 'kempt']
141+
self.unknown_parameters = ['Gb', 'SG', 'p2', 'ka2', 'kd', 'kempt']
142142

143143
# initial guess for unknown parameter
144144
self.start_guess = np.array(
145-
[self.model_parameters.Gb, self.model_parameters.SG,
145+
[self.model_parameters.Gb, self.model_parameters.SG, self.model_parameters.p2,
146146
self.model_parameters.ka2, self.model_parameters.kd, self.model_parameters.kempt])
147147

148148
# initial guess for the SD of each parameter
149-
self.start_guess_sigma = np.array([1, 5e-4, 1e-3, 1e-3, 1e-3])
149+
self.start_guess_sigma = np.array([1, 5e-4, 1e-3, 1e-3, 1e-3, 1e-3])
150150

151151
# Get the hour of the day for each data point
152152
t = np.array(data.t.dt.hour.values).astype(int)
@@ -938,9 +938,10 @@ def __log_likelihood(self, theta: np.ndarray, rbg_data: ReplayBGData):
938938
# Set model parameters to current guess
939939
(self.model_parameters.Gb,
940940
self.model_parameters.SG,
941+
self.model_parameters.p2,
941942
self.model_parameters.ka2,
942943
self.model_parameters.kd,
943-
self.model_parameters.kempt) = theta[0:5]
944+
self.model_parameters.kempt) = theta[0:6]
944945

945946
self.model_parameters.SI_B = theta[self.pos_SI_B] if self.pos_SI_B else self.model_parameters.SI_B
946947
self.model_parameters.SI_L = theta[self.pos_SI_L] if self.pos_SI_L else self.model_parameters.SI_L
@@ -1002,9 +1003,10 @@ def __log_likelihood_exercise(self, theta: np.ndarray, rbg_data: ReplayBGData):
10021003
# Set model parameters to current guess
10031004
(self.model_parameters.Gb,
10041005
self.model_parameters.SG,
1006+
self.model_parameters.p2,
10051007
self.model_parameters.ka2,
10061008
self.model_parameters.kd,
1007-
self.model_parameters.kempt) = theta[0:5]
1009+
self.model_parameters.kempt) = theta[0:6]
10081010

10091011
self.model_parameters.SI_B = theta[self.pos_SI_B] if self.pos_SI_B else self.model_parameters.SI_B
10101012
self.model_parameters.SI_L = theta[self.pos_SI_L] if self.pos_SI_L else self.model_parameters.SI_L
@@ -1069,9 +1071,10 @@ def __log_likelihood_extended(self, theta: np.ndarray, rbg_data: ReplayBGData):
10691071
# Set model parameters to current guess
10701072
(self.model_parameters.Gb,
10711073
self.model_parameters.SG,
1074+
self.model_parameters.p2,
10721075
self.model_parameters.ka2,
10731076
self.model_parameters.kd,
1074-
self.model_parameters.kempt) = theta[0:5]
1077+
self.model_parameters.kempt) = theta[0:6]
10751078

10761079
self.model_parameters.SI_B = theta[self.pos_SI_B] if self.pos_SI_B else self.model_parameters.SI_B
10771080
self.model_parameters.SI_L = theta[self.pos_SI_L] if self.pos_SI_L else self.model_parameters.SI_L

‎py_replay_bg/model/t1d_model_single_meal.py‎

Lines changed: 4 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -134,15 +134,15 @@ def __init__(self,
134134
self.model_parameters = ModelParametersT1DSingleMeal(data, bw, u2ss)
135135

136136
# Unknown parameters
137-
self.unknown_parameters = ['Gb', 'SG', 'ka2', 'kd', 'kempt']
137+
self.unknown_parameters = ['Gb', 'SG', 'p2', 'ka2', 'kd', 'kempt']
138138

139139
# initial guess for unknown parameter
140140
self.start_guess = np.array(
141-
[self.model_parameters.Gb, self.model_parameters.SG,
141+
[self.model_parameters.Gb, self.model_parameters.SG, self.model_parameters.p2,
142142
self.model_parameters.ka2, self.model_parameters.kd, self.model_parameters.kempt])
143143

144144
# initial guess for the SD of each parameter
145-
self.start_guess_sigma = np.array([1, 5e-4, 1e-3, 1e-3, 1e-3])
145+
self.start_guess_sigma = np.array([1, 5e-4, 1e-3, 1e-3, 1e-3, 1e-3])
146146

147147
if is_twin:
148148
# Attach SI
@@ -598,6 +598,7 @@ def __log_likelihood(
598598
# Set model parameters to current guess
599599
(self.model_parameters.Gb,
600600
self.model_parameters.SG,
601+
self.model_parameters.p2,
601602
self.model_parameters.ka2,
602603
self.model_parameters.kd,
603604
self.model_parameters.kempt,

0 commit comments

Comments
 (0)