Skip to content

Commit

Permalink
Parameterized sample rate for risk index caculation; Added UT
Browse files Browse the repository at this point in the history
  • Loading branch information
yihuicai committed Aug 23, 2023
1 parent d83fce6 commit 7fe8168
Show file tree
Hide file tree
Showing 3 changed files with 59 additions and 10 deletions.
22 changes: 14 additions & 8 deletions simglucose/analysis/report.py
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,8 @@
# from pandas.plotting import lag_plot
import logging

from simglucose.sensor.cgm import CGMSensor

logger = logging.getLogger(__name__)


Expand Down Expand Up @@ -92,18 +94,18 @@ def percent_stats(BG, ax=None):
return p_stats, fig, ax


def risk_index_trace(df_BG, visualize=False):
chunk_BG = [df_BG.iloc[i:i + 20, :] for i in range(0, len(df_BG), 20)]
def risk_index_trace(df_BG, sample_rate, visualize=False):
chunk_BG = [df_BG.iloc[i:i + sample_rate, :] for i in range(0, len(df_BG), sample_rate)]

if len(chunk_BG[-1]) != 20: # Remove the last chunk which is not full
if len(chunk_BG[-1]) != sample_rate: # Remove the last chunk which is not full
chunk_BG.pop()

fBG = [
1.509 * (np.log(BG[BG > 0]) ** 1.084 - 5.381) for BG in chunk_BG
]

rl = [(10 * (fbg * (fbg < 0))**2).mean() for fbg in fBG]
rh = [(10 * (fbg * (fbg > 0))**2).mean() for fbg in fBG]
rl = [(10 * (fbg * (fbg < 0)) ** 2).mean() for fbg in fBG]
rh = [(10 * (fbg * (fbg > 0)) ** 2).mean() for fbg in fBG]

LBGI = pd.concat(rl, axis=1).transpose()
HBGI = pd.concat(rh, axis=1).transpose()
Expand Down Expand Up @@ -238,7 +240,7 @@ def CVGA(BG_list, label=None):
edgecolors='k',
zorder=4,
label='%s (A: %d%%, B: %d%%, C: %d%%, D: %d%%, E: %d%%)' %
(l, 100 * A, 100 * B, 100 * C, 100 * D, 100 * E))
(l, 100 * A, 100 * B, 100 * C, 100 * D, 100 * E))
zone_stats.append((A, B, C, D, E))

zone_stats = pd.DataFrame(zone_stats, columns=['A', 'B', 'C', 'D', 'E'])
Expand All @@ -247,12 +249,16 @@ def CVGA(BG_list, label=None):
return zone_stats, fig, ax


def report(df, save_path=None):
def report(df, cgm_sensor=None, save_path=None):
BG = df.unstack(level=0).BG

sample_rate = 20 # Use 20 as default sample rate in one hour
if cgm_sensor is not None:
sample_rate = int(60 / cgm_sensor.sample_time)

fig_ensemble, ax1, ax2, ax3 = ensemblePlot(df)
pstats, fig_percent, ax4 = percent_stats(BG)
ri_per_hour, ri_mean, fig_ri, ax5 = risk_index_trace(BG, visualize=False)
ri_per_hour, ri_mean, fig_ri, ax5 = risk_index_trace(BG, sample_rate, visualize=False)
zone_stats, fig_cvga, ax6 = CVGA(BG, label='')
axes = [ax1, ax2, ax3, ax4, ax5, ax6]
figs = [fig_ensemble, fig_percent, fig_ri, fig_cvga]
Expand Down
5 changes: 3 additions & 2 deletions simglucose/simulation/user_interface.py
Original file line number Diff line number Diff line change
Expand Up @@ -361,9 +361,10 @@ def simulate(sim_time=None,
if controller is None:
controller = pick_controller()

cgm_sensor = CGMSensor.withName(cgm_name, seed=cgm_seed)

def local_build_env(pname):
patient = T1DPatient.withName(pname)
cgm_sensor = CGMSensor.withName(cgm_name, seed=cgm_seed)
insulin_pump = InsulinPump.withName(insulin_pump_name)
scen = copy.deepcopy(scenario)
env = T1DSimEnv(patient, cgm_sensor, insulin_pump, scen)
Expand All @@ -380,7 +381,7 @@ def local_build_env(pname):
results = batch_sim(sim_instances, parallel=parallel)

df = pd.concat(results, keys=[s.env.patient.name for s in sim_instances])
results, ri_per_hour, zone_stats, figs, axes = report(df, save_path)
results, ri_per_hour, zone_stats, figs, axes = report(df, cgm_sensor, save_path)

return results

Expand Down
42 changes: 42 additions & 0 deletions tests/test_report.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,42 @@
from simglucose.analysis.report import risk_index_trace
from simglucose.sensor.cgm import CGMSensor
from simglucose.simulation.rendering import Viewer
from datetime import datetime
import pandas as pd
import unittest
import logging
import os

logger = logging.getLogger(__name__)
TESTDATA_FILENAME = os.path.join(os.path.dirname(__file__), 'sim_results.csv')


class TestReport(unittest.TestCase):
def setUp(self):
self.df = pd.concat([pd.read_csv(TESTDATA_FILENAME, index_col=0)], keys=['test'])

def test_risk_index_trace(self):
BG = self.df.unstack(level=0).BG
sample_time = CGMSensor.withName("Dexcom").sample_time
sample_rate = int(60 / sample_time)
ri_per_hour, ri_mean, fig, axes = risk_index_trace(BG, sample_rate)

LBGI = ri_per_hour.transpose().LBGI
HBGI = ri_per_hour.transpose().HBGI
RI = ri_per_hour.transpose()["Risk Index"]

self.assertEqual(LBGI.size, 48)
self.assertEqual(LBGI.iloc[-1].test, 0.8429957158900777)
self.assertEqual(LBGI.iloc[0].test, 0.0)

self.assertEqual(HBGI.size, 48)
self.assertEqual(HBGI.iloc[-1].test, 0.0)
self.assertEqual(HBGI.iloc[0].test, 2.755277346918188)

self.assertEqual(RI.size, 48)
self.assertEqual(RI.iloc[-1].test, 0.8429957158900777)
self.assertEqual(RI.iloc[0].test, 2.755277346918188)


if __name__ == '__main__':
unittest.main()

0 comments on commit 7fe8168

Please sign in to comment.