#!/usr/bin/env python # coding: utf-8 # [![img](..\ep8_2018-57x57.png)](http://endlesspint.com/) # # # Sleep No More - Caffeine # # **post @** [endlesspint.com](http://endlesspint.com/2019-05-17-sleep-no-more-caffeine/) # In[1]: import numpy as np import pandas as pd import matplotlib.pyplot as plt import datetime plt.style.use('ggplot') get_ipython().run_line_magic('matplotlib', 'inline') # In[2]: half_life = 5.0 half_hour_decay = np.exp(-np.log(2) * (0.5 / half_life)) half_hour_decay # 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) # In[4]: intake = pd.read_excel("11day_intake.xlsx") intake.head(10) # In[5]: intake_ts_day = pd.pivot_table(intake, values="Caffeine (mg)", index="Time", columns="Day").fillna(0) intake_ts_day # In[6]: intake_ts_day[2][datetime.time(9,0)] # 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() # 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) # 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 # 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) # ## different attempts at plotting with time on x-axis # In[11]: df_caff_down2 = df_caff_down.reset_index().set_index('day_time_split') df_caff_down2.head() # In[12]: df_caff_down2[df_caff_down2.columns[1:-2]].plot(figsize=(18,10)) # In[13]: df_caff_down.reset_index()[df_caff_down.columns[0:11]].plot(figsize=(18,10)) # 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) # 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() # 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") # ## the caffeine comedown when extending half-life to 7 hours # In[22]: half_life = 7.0 half_hour_decay = np.exp(-np.log(2) * (0.5 / half_life)) half_hour_decay # 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) # 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() # ## self v others # 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[ ]: