In [1]:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import datetime
plt.style.use('ggplot')
%matplotlib inline
In [2]:
half_life = 5.0
half_hour_decay = np.exp(-np.log(2) * (0.5 / half_life))
half_hour_decay
Out[2]:
In [3]:
# toy sample
t = np.zeros(10)
c = np.array([0, 0, 0, 30, 0, 0, 10, 0, 0, 0])
for i in range(2, len(t)):
t[i] += t[i-1] * half_hour_decay + c[i-1]
plt.plot(t)
Out[3]:
In [4]:
intake = pd.read_excel("11day_intake.xlsx")
intake.head(10)
Out[4]:
In [5]:
intake_ts_day = pd.pivot_table(intake, values="Caffeine (mg)", index="Time", columns="Day").fillna(0)
intake_ts_day
Out[5]:
In [6]:
intake_ts_day[2][datetime.time(9,0)]
Out[6]:
In [7]:
# dummy Y/M/D
six_am = datetime.datetime(2019, 1, 1, 6, 0)
time_idx = []
for i in range(48):
step = 30 * i
time_idx.append((six_am + datetime.timedelta(minutes = step)).time())
print(time_idx[:10])
df_time_idx = pd.DataFrame(np.zeros(len(time_idx)), index=time_idx)
df_time_idx.head()
Out[7]:
In [8]:
df = pd.concat([intake_ts_day, df_time_idx], axis=1, sort=True).fillna(0)
df_intake_ts_day = df.loc[time_idx, :11]
df_intake_ts_day.columns = ['day_' + str(i) for i in df_intake_ts_day.columns]
df_intake_ts_day.head(20)
Out[8]:
In [9]:
caff_down = {}
for col in df_intake_ts_day.columns:
dayy = np.zeros(len(time_idx))
caff = df_intake_ts_day[col].values
for i in range(2, len(dayy)):
dayy[i] += dayy[i-1] * half_hour_decay + caff[i-1]
caff_down[col] = dayy
df_caff_down = pd.DataFrame(caff_down, index=time_idx)
df_caff_down
Out[9]:
In [10]:
df_caff_down['day_split'] = 'A - '
df_caff_down.loc[time_idx[-12:], 'day_split'] = 'B - '
df_caff_down['day_time'] = [str(i)[:5] for i in df_caff_down.index]
df_caff_down['day_time_split'] =df_caff_down.day_split + df_caff_down.day_time
df_caff_down.tail(15)
Out[10]:
In [11]:
df_caff_down2 = df_caff_down.reset_index().set_index('day_time_split')
df_caff_down2.head()
Out[11]:
In [12]:
df_caff_down2[df_caff_down2.columns[1:-2]].plot(figsize=(18,10))
Out[12]:
In [13]:
df_caff_down.reset_index()[df_caff_down.columns[0:11]].plot(figsize=(18,10))
Out[13]:
In [14]:
plt.figure(figsize=(18, 10))
plt.xticks(rotation='vertical')
days = ['day_1', 'day_2', 'day_3', 'day_4', 'day_5', 'day_6', 'day_7', 'day_8',
'day_9', 'day_10', 'day_11']
for col in days:
plt.plot( 'day_time_split', col, data=df_caff_down)
# plt.axvline(x='A - 22:00')
# plt.axhline(y=50)
plt.legend()
plt.title("Caffeine Half-Life", size=20)
plt.savefig("caff_hl.jpg")
In [15]:
df_caff_down['day_q25'] = df_caff_down.quantile(q=0.25, axis=1)
df_caff_down['day_q50'] = df_caff_down.quantile(q=0.50, axis=1)
df_caff_down['day_q75'] = df_caff_down.quantile(q=0.75, axis=1)
df_caff_down.head(15)
Out[15]:
In [16]:
plt.figure(figsize=(18, 10))
plt.xticks(rotation='vertical')
for col in zip(['day_q25', 'day_q50', 'day_q75'],['skyblue', 'darkblue', 'skyblue']):
plt.plot( 'day_time_split', col[0], data=df_caff_down, color=col[1], linewidth=4)
# plt.axvline(x='A - 22:00')
# plt.axhline(y=50)
plt.title("Caffeine Half-Life, IQR", size=20)
plt.savefig("caff_hl_iqr.jpg")
# plt.legend()
In [17]:
intake_ts_day.sum().describe()
Out[17]:
In [18]:
intake_ts_day.sum().hist()
plt.title("Hist")
plt.show()
In [19]:
import matplotlib.pyplot as plt
import numpy as np
import scipy.stats as stats
import math
mu = intake_ts_day.sum().describe()['mean']
sigma = intake_ts_day.sum().describe()['std']
x = np.linspace(mu - 2*sigma, mu + 2*sigma, 100)
plt.plot(x, stats.norm.pdf(x, mu, sigma))
plt.title("Norm")
plt.show()
In [20]:
from scipy.stats import poisson
mu = intake_ts_day.sum().describe()['mean']
x_p = np.arange(poisson.ppf(0.01, mu),
poisson.ppf(0.99, mu))
plt.plot(x_p, poisson.pmf(x_p, mu), 'bo', ms=8, label='poisson pmf')
plt.title("Poisson")
plt.show()
In [21]:
# summary charts, of a sort
f, (ax1, ax2, ax3) = plt.subplots(1, 3, figsize=(18,6))
ax1.hist(intake_ts_day.sum())
ax1.set_title("Hist")
ax2.plot(x, stats.norm.pdf(x, mu, sigma))
ax2.set_title("Norm")
ax3.plot(x_p, poisson.pmf(x_p, mu), 'bo', ms=8, label='poisson pmf')
ax3.set_title("Poisson")
plt.tight_layout()
plt.savefig("caff_intake_summaries.jpg")
In [22]:
half_life = 7.0
half_hour_decay = np.exp(-np.log(2) * (0.5 / half_life))
half_hour_decay
Out[22]:
In [23]:
caff_down = {}
for col in df_intake_ts_day.columns:
dayy = np.zeros(len(time_idx))
caff = df_intake_ts_day[col].values
for i in range(2, len(dayy)):
dayy[i] += dayy[i-1] * half_hour_decay + caff[i-1]
caff_down[col] = dayy
df_caff_down7 = pd.DataFrame(caff_down, index=time_idx)
df_caff_down7['day_split'] = 'A - '
df_caff_down7.loc[time_idx[-12:], 'day_split'] = 'B - '
df_caff_down7['day_time'] = [str(i)[:5] for i in df_caff_down7.index]
df_caff_down7['day_time_split'] =df_caff_down7.day_split + df_caff_down7.day_time
df_caff_down7['day_q25'] = df_caff_down7.quantile(q=0.25, axis=1)
df_caff_down7['day_q50'] = df_caff_down7.quantile(q=0.50, axis=1)
df_caff_down7['day_q75'] = df_caff_down7.quantile(q=0.75, axis=1)
df_caff_down7.tail(15)
Out[23]:
In [41]:
plt.figure(figsize=(18, 10))
plt.xticks(rotation='vertical')
plt.plot( 'day_time_split', 'day_q50', data=df_caff_down, color='skyblue', linewidth=4, label='5hr HL')
plt.plot( 'day_time_split', 'day_q50', data=df_caff_down7, color='deepskyblue', linewidth=4, label='7hr HL')
plt.title("Median Caffeine Intake, 5 v 7 hrs Half-Life", size=20)
plt.legend()
plt.savefig("caff_diff_hl.jpg")
plt.show()
In [42]:
mu = intake_ts_day.sum().describe()['mean']
sigma = intake_ts_day.sum().describe()['std']
x = np.linspace(mu - 2*sigma, mu + 5*sigma, 100)
plt.plot(x, stats.norm.pdf(x, mu, sigma), label="ep8")
mu_jcsm = 319.0
sigma_jcsm = 181.0
plt.plot(x, stats.norm.pdf(x, mu_jcsm, sigma_jcsm), label="JCSM")
plt.title("Comp")
plt.legend()
plt.savefig("caff_comp.jpg")
plt.show()
In [ ]: