03_Models_replication_and_evaluation.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "|                                                                    title                                                                    |          authors          |  N |     age group    |     data     |  window |           Features           |    methods   |           results           |\n",
    "|:-------------------------------------------------------------------------------------------------------------------------------------------:|:-------------------------:|:--:|:----------------:|:------------:|:-------:|:----------------------------:|:------------:|:---------------------------:|\n",
    "| Improving energy expenditure estimates from wearable devices: A machine learning approach                                                   | Driscoll et al., 2020     | 59 |      44(14)      | accel+, demo |  1 min  |                              |   RF regres  | RMSE=1-1.37 METS,   R2=0.85 |\n",
    "| Comparison of linear and non-linear models for predicting energy expenditure from raw accelerometer data                                    | Montoye et al., 2017      |    |                  |              |  30sec  | 6 to 36-time domain and demo | LM, LMM, ANN |                             |\n",
    "| A random forest classifier for the prediction of energy expenditure and type of physical activity from wrist and hip accelerometers         | Elleis et al., 2014       | 40 |     35.8(12)     |    accel+    |  1 min  |      45-time-freq on svm     |   RF regres  |        RMSE=1.09 METS       |\n",
    "| An   artificial neural network to estimate physical activity energy expenditure   and identify physical activity type from an accelerometer | Staudenmayer et al., 2009 | 48 |        35        |  accel, demo |  1 min  |               6              |      ANN     |     RMSE=1.22(0.08) METS    |\n",
    "| Using   Deep Learning for Energy Expenditure Estimation with Wearable Sensors                                                               | Zhu et al., 2015          | 30 | 27.8(6.9)[19-45] | accel+, demo | 5.12sec |                              |      CNN     |                             |"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "os.getcwd()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "import os\n",
    "os.chdir(\"...\")\n",
    "os.listdir()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import datetime\n",
    "import pandas as pd\n",
    "import pickle\n",
    "import random\n",
    "import glob\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "# Display figures inline in Jupyter notebook\n",
    "import seaborn as sns\n",
    "# Use seaborn style defaults and set the default figure size\n",
    "sns.set(rc={'figure.figsize':(11, 4)})\n",
    "from sklearn import linear_model\n",
    "from sklearn.metrics import r2_score\n",
    "from sklearn.metrics import mean_squared_error\n",
    "from math import sqrt\n",
    "import csv\n",
    "import statsmodels.api as sm\n",
    "import statsmodels.formula.api as smf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def read_data(data, dataDir, withActivities= False): \n",
    "\n",
    "    ankle = ['ankle_x', 'ankle_y', 'ankle_z']\n",
    "    wrist = ['wrist_x', 'wrist_y', 'wrist_z'] \n",
    "\n",
    "    print('read data:', data)\n",
    "    df = pd.read_csv(dataDir+data+\".csv\", header = 0, index_col = None, low_memory=False)\n",
    "\n",
    "    df_cosmed = df[['EEm', 'time']]\n",
    "    df_cosmed.rename(columns={'time':'time_cosmed'})\n",
    "\n",
    "    # df = df.drop('labels', axis=1).dropna(axis=0, how='any')\n",
    "\n",
    "    time = df.time\n",
    "\n",
    "    # accel data\n",
    "    df_ankle = df[ankle]\n",
    "    # df_ankle = SQSumSq(df_ankle) # accel data in SVM\n",
    "    df_wrist = df[wrist]\n",
    "    # df_wrist = SQSumSq(df_wrist) # accel data in SVM\n",
    "\n",
    "    df = pd.concat([time, df_ankle, df_wrist], axis=1)\n",
    "    #         df = df.rename(columns={0:'ankle',1:'wrist'})\n",
    "\n",
    "    if withActivities == True:\n",
    "        # load the predicted activities\n",
    "        predAR_path = '/dataDir/origEEm_predActivities/'\n",
    "        predAR = pd.read_csv(predAR_path+data+\"_ARpred_origEEm.csv\", \n",
    "                             header = 0, index_col = None, low_memory=False)\n",
    "        predAR = predAR[['time','label']]\n",
    "\n",
    "        # concat with activities\n",
    "        df = df.merge(predAR, left_on='time', right_on='time', how='left')\n",
    "        df['label'] = df['label'].fillna(method='ffill')\n",
    "        df['label'] = df['label'].fillna(method='bfill')\n",
    "    #eIF\n",
    "\n",
    "    inv_yhat = np.empty((df.shape[0], 2))\n",
    "    inv_yhat.fill(np.nan)\n",
    "    inv_yhat[:df_cosmed.shape[0]] = df_cosmed\n",
    "    df_cosmed = pd.DataFrame(inv_yhat, columns=['EEm','time_cos'])\n",
    "\n",
    "    df = pd.concat([df, df_cosmed], axis=1)\n",
    "    df['participant'] = data\n",
    "    \n",
    "    return df\n",
    "#eDEF"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Feature construction\n",
    "def featConstr(df, windowSize):\n",
    "#     INPUT:\n",
    "#             - df: the dataframe with accel and EEm data\n",
    "#             - windowSize: the aggregation window in seconds to build the features\n",
    "#     OUTPUT: Returns the set of features in a dataframe\n",
    "\n",
    "    print('Building features ...')\n",
    "    # Prepare predictors/target sets\n",
    "    cols = ['time_cos', 'EEm']\n",
    "    df_geneA = df.drop(cols, axis=1)\n",
    "    df_cosmed = df[cols]\n",
    "\n",
    "    df_geneA.set_index('time', inplace=True)\n",
    "    df_geneA.index = pd.to_datetime(df_geneA.index, unit = \"ms\")\n",
    "\n",
    "    df_cosmed.set_index('time_cos', inplace=True)\n",
    "    df_cosmed.index = pd.to_datetime(df_cosmed.index, unit = \"ms\")\n",
    "\n",
    "    preds = df_geneA.drop(['participant'], axis=1).sort_index()\n",
    "    idx = df_cosmed['EEm'].notnull()\n",
    "    target = df_cosmed['EEm'][idx].sort_index()\n",
    "    \n",
    "    # compute features\n",
    "    key = str(windowSize) + 'S' # sampling rate for downsampling\n",
    "    print('=== Mean ...')\n",
    "    preds_mean = preds.resample(key).mean()\n",
    "    aggFun = 'mean'\n",
    "    preds_mean.columns = ['ankle_x_'+aggFun,'ankle_y_'+aggFun,'ankle_z_'+aggFun,\n",
    "                         'wrist_x_'+aggFun,'wrist_y_'+aggFun,'wrist_z_'+aggFun]\n",
    "    print('=== SD ...')\n",
    "    preds_std = preds.resample(key).std()\n",
    "    aggFun = 'std'\n",
    "    preds_std.columns = ['ankle_x_'+aggFun,'ankle_y_'+aggFun,'ankle_z_'+aggFun,\n",
    "                         'wrist_x_'+aggFun,'wrist_y_'+aggFun,'wrist_z_'+aggFun]\n",
    "    print('=== Min ...')\n",
    "    preds_min = preds.resample(key).min()\n",
    "    aggFun = 'min'\n",
    "    preds_min.columns = ['ankle_x_'+aggFun,'ankle_y_'+aggFun,'ankle_z_'+aggFun,\n",
    "                         'wrist_x_'+aggFun,'wrist_y_'+aggFun,'wrist_z_'+aggFun]\n",
    "    print('=== Max ...')\n",
    "    preds_max = preds.resample(key).max()\n",
    "    aggFun = 'max'\n",
    "    preds_max.columns = ['ankle_x_'+aggFun,'ankle_y_'+aggFun,'ankle_z_'+aggFun,\n",
    "                         'wrist_x_'+aggFun,'wrist_y_'+aggFun,'wrist_z_'+aggFun]\n",
    "    print('=== 10th percentile ...')\n",
    "    preds_q10 = preds.resample(key).quantile(.1)\n",
    "    aggFun = 'q10'\n",
    "    preds_q10.columns = ['ankle_x_'+aggFun,'ankle_y_'+aggFun,'ankle_z_'+aggFun,\n",
    "                         'wrist_x_'+aggFun,'wrist_y_'+aggFun,'wrist_z_'+aggFun]\n",
    "    print('=== 25th percentile ...')\n",
    "    preds_q25 = preds.resample(key).quantile(.25)\n",
    "    aggFun = 'q25'\n",
    "    preds_q25.columns = ['ankle_x_'+aggFun,'ankle_y_'+aggFun,'ankle_z_'+aggFun,\n",
    "                         'wrist_x_'+aggFun,'wrist_y_'+aggFun,'wrist_z_'+aggFun]\n",
    "    print('=== 50th percentile ...')\n",
    "    preds_q50 = preds.resample(key).quantile(.5)\n",
    "    aggFun = 'q50'\n",
    "    preds_q50.columns = ['ankle_x_'+aggFun,'ankle_y_'+aggFun,'ankle_z_'+aggFun,\n",
    "                         'wrist_x_'+aggFun,'wrist_y_'+aggFun,'wrist_z_'+aggFun]\n",
    "    print('=== 75th percentile ...')\n",
    "    preds_q75 = preds.resample(key).quantile(.75)\n",
    "    aggFun = 'q75'\n",
    "    preds_q75.columns = ['ankle_x_'+aggFun,'ankle_y_'+aggFun,'ankle_z_'+aggFun,\n",
    "                         'wrist_x_'+aggFun,'wrist_y_'+aggFun,'wrist_z_'+aggFun]\n",
    "    print('=== 90th percentile ...')\n",
    "    preds_q90 = preds.resample(key).quantile(.9)\n",
    "    aggFun = 'q90'\n",
    "    preds_q90.columns = ['ankle_x_'+aggFun,'ankle_y_'+aggFun,'ankle_z_'+aggFun,\n",
    "                         'wrist_x_'+aggFun,'wrist_y_'+aggFun,'wrist_z_'+aggFun]\n",
    "\n",
    "    features = pd.concat([preds_mean, preds_std, \n",
    "                          preds_min, preds_max,\n",
    "                          preds_q10, preds_q25,\n",
    "                          preds_q50, preds_q75,\n",
    "                          preds_q90\n",
    "                         ], axis=1)\n",
    "\n",
    "    print('=== Covariance of adjacent window ...')\n",
    "    ankle_x_cov = []\n",
    "    ankle_y_cov = []\n",
    "    ankle_z_cov = []\n",
    "    wrist_x_cov = []\n",
    "    wrist_y_cov = []\n",
    "    wrist_z_cov = []\n",
    "\n",
    "    for i in range(len(features.index)):\n",
    "        a= features.index[i]\n",
    "        b= a + datetime.timedelta(seconds=30)\n",
    "        c= b + datetime.timedelta(seconds=30)\n",
    "\n",
    "        x= preds.ankle_x.loc[a:b].reset_index().ankle_x\n",
    "        y= preds.ankle_x.loc[b:c].reset_index().ankle_x\n",
    "        ankle_x_cov.append(x.cov(y))\n",
    "\n",
    "        x= preds.ankle_y.loc[a:b].reset_index().ankle_y\n",
    "        y= preds.ankle_y.loc[b:c].reset_index().ankle_y\n",
    "        ankle_y_cov.append(x.cov(y))\n",
    "\n",
    "        x= preds.ankle_z.loc[a:b].reset_index().ankle_z\n",
    "        y= preds.ankle_z.loc[b:c].reset_index().ankle_z\n",
    "        ankle_z_cov.append(x.cov(y))\n",
    "\n",
    "        x= preds.wrist_x.loc[a:b].reset_index().wrist_x\n",
    "        y= preds.wrist_x.loc[b:c].reset_index().wrist_x\n",
    "        wrist_x_cov.append(x.cov(y))\n",
    "\n",
    "        x= preds.wrist_y.loc[a:b].reset_index().wrist_y\n",
    "        y= preds.wrist_y.loc[b:c].reset_index().wrist_y\n",
    "        wrist_y_cov.append(x.cov(y))\n",
    "\n",
    "        x= preds.wrist_z.loc[a:b].reset_index().wrist_z\n",
    "        y= preds.wrist_z.loc[b:c].reset_index().wrist_z\n",
    "        wrist_z_cov.append(x.cov(y))\n",
    "    #eFOR\n",
    "    \n",
    "    features['ankle_x_cov'] = ankle_x_cov\n",
    "    features['ankle_y_cov'] = ankle_y_cov\n",
    "    features['ankle_z_cov'] = ankle_z_cov\n",
    "    features['wrist_x_cov'] = wrist_x_cov\n",
    "    features['wrist_y_cov'] = wrist_y_cov\n",
    "    features['wrist_z_cov'] = wrist_z_cov\n",
    "    \n",
    "    # add target\n",
    "    features['EEm'] = target.resample(key).mean()\n",
    "\n",
    "    # delete nan\n",
    "    features = features.dropna()\n",
    "    \n",
    "    return features\n",
    "#eDEF"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def select_Data(name, randomSeed, numberForVal, df_gene):\n",
    "# selects participant's data and builds the test*, val and training sets\n",
    "# *test set is always the data from the selected (name) participant\n",
    "# INPUT:\n",
    "#     - name: The id-name of test participant (e.g. GOTOV00)\n",
    "#     - randomSeed: the given random seed\n",
    "#     - numberForVal: number of participants used as validation set\n",
    "#     - df_gene: the data frame including predictors (accel,activity data), and target data (EE)\n",
    "# OUTPUT: Returns train,test and validation df for the given name    \n",
    "\n",
    "    df_train, df_val, df_test = pd.DataFrame(), pd.DataFrame(), pd.DataFrame()\n",
    "#     data_to_select_random_ = []\n",
    "    print(name)\n",
    "    random.seed(randomSeed)\n",
    "    print('Random seed is:',randomSeed)\n",
    "\n",
    "    data_to_select_random_ = list(df_gene['ID'].unique())\n",
    "\n",
    "    if name in data_to_select_random_: data_to_select_random_.remove(name) #eIF\n",
    "    all_act_data = ['GOTOV08', 'GOTOV10', 'GOTOV11', 'GOTOV12', #'GOTOV16,'\n",
    "                    'GOTOV17', 'GOTOV20', 'GOTOV21', 'GOTOV28',\n",
    "                    'GOTOV29', 'GOTOV31', 'GOTOV33', 'GOTOV35'] \n",
    "    \n",
    "    if name in all_act_data: all_act_data.remove(name)#eIF\n",
    "    indoors_act_data = ['GOTOV22', 'GOTOV23', 'GOTOV13','GOTOV14', \n",
    "                        'GOTOV24', 'GOTOV25', 'GOTOV26', 'GOTOV27', \n",
    "                        'GOTOV30', 'GOTOV32', 'GOTOV34', 'GOTOV36']\n",
    "    \n",
    "    if name in indoors_act_data: indoors_act_data.remove(name) #fi\n",
    "    validation_data = [random.sample(all_act_data, numberForVal)[0], \n",
    "                       random.sample(indoors_act_data, numberForVal)[0]]\n",
    "\n",
    "    print('Getting val and train data.....')\n",
    "    for data in data_to_select_random_: \n",
    "        if data in validation_data:\n",
    "            df_val = df_val.append(df_gene.query('ID == \"'+data+'\"'))\n",
    "            print('val_data:', data)\n",
    "        else:\n",
    "            df_train = df_train.append(df_gene.query('ID == \"'+data+'\"'))\n",
    "        #eIF\n",
    "    #eFOR\n",
    "\n",
    "    print('Getting test data.....')\n",
    "    df_test = df_test.append(df_gene.query('ID == \"'+name+'\"'))\n",
    "\n",
    "    print('Done creating all dataframes.....')\n",
    "    return df_train, df_val, df_test\n",
    "#eDEF"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dataDir= '/...'\n",
    "saveDir = 'sequences/features_MontoyePaper_replication/'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Read data and build features"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "dataDir= '/...'\n",
    "saveDir = 'sequences/features_MontoyePaper_replication/'\n",
    "features_bindAll = pd.DataFrame()\n",
    "\n",
    "w = 30 # window seconds\n",
    "for data in glob.glob(dataDir+'*csv'):\n",
    "    data = os.path.splitext(os.path.basename(data))[0]\n",
    "    # skip participants with no cosmoed or some accel are missing\n",
    "    if data in ['GOTOV02', 'GOTOV03', 'GOTOV04', 'GOTOV19']:\n",
    "        continue\n",
    "    else:\n",
    "        df = read_data(data, dataDir, withActivities= False)\n",
    "#         print(df)\n",
    "        features = featConstr(df, windowSize=w)\n",
    "        features['ID'] = data\n",
    "#         print(features)\n",
    "        # save features\n",
    "        print('Saving features.... ')\n",
    "        features.to_csv(saveDir+data+'_features_'+str(w)+'sec.csv')  \n",
    "        features_bindAll = features_bindAll.append(features)\n",
    "\n",
    "    #eIF\n",
    "#eFOR"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "features_bindAll.to_csv(saveDir+'features_bindAll_'+str(w)+'sec.csv')  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Linear Model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "w=30\n",
    "print(saveDir+'features_bindAll_'+str(w)+'sec.csv')\n",
    "features_bindAll = pd.read_csv(saveDir+'features_bindAll_'+str(w)+'sec.csv')\n",
    "# print(features_bindAll.columns)\n",
    "\n",
    "features_bindAll.set_index('time', inplace=True)\n",
    "features_bindAll.index = pd.to_datetime(features_bindAll.index)\n",
    "features_bindAll"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ankle_cols = []\n",
    "ankle_cols = [col for col in features_bindAll.columns if 'ankle' in col]\n",
    "ankle_cols.append('EEm')\n",
    "ankle_cols.append('ID')\n",
    "\n",
    "wrist_cols =[]\n",
    "wrist_cols = [col for col in features_bindAll.columns if 'wrist' in col]\n",
    "wrist_cols.append('EEm')\n",
    "wrist_cols.append('ID')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "devices = 'aw'\n",
    "\n",
    "if devices == 'ankle':\n",
    "    features_bindAll = features_bindAll[ankle_cols]\n",
    "elif devices == 'wrist':\n",
    "    features_bindAll = features_bindAll[wrist_cols]\n",
    "features_bindAll"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "features_bindAll.columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "EEm_details = pd.read_csv('EEm_details_per_partic.csv')  \n",
    "\n",
    "\n",
    "numberForVal = 1 # the number of participants you want per gtoup as validation set (random seed is the number of participants' id)\n",
    "# randomSeed = i\n",
    "\n",
    "data = features_bindAll\n",
    "\n",
    "model_name = 'LM_'+devices+'_'+str(w)+'sec'\n",
    "figures_folder = 'images/' + model_name + '/predPlot/'\n",
    "csv_folder = figures_folder\n",
    "\n",
    "if not os.path.exists(figures_folder):\n",
    "    os.makedirs(figures_folder)\n",
    "\n",
    "get_all_patients_rsquared = pd.DataFrame(columns=['participant', 'rsquared', 'inRsquared', 'outRsquared'])\n",
    "get_all_patients_rms = pd.DataFrame(columns=['participant', 'rms', 'inRms', 'outRms'])\n",
    "stats_EEm_results = pd.DataFrame(columns=['participant', \n",
    "                                          'true_tot_mean', 'true_in_mean', 'true_out_mean', \n",
    "                                          'pred_tot_mean', 'pred_in_mean', 'pred_out_mean'])\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "df_allPred = pd.DataFrame()\n",
    "df_allPredin = pd.DataFrame()\n",
    "df_allPredout = pd.DataFrame()\n",
    "\n",
    "for i in range(5, 37):\n",
    "    randomSeed = i\n",
    "    if i in [2,3,4,19]: #[12, 19, 23, 16]: \n",
    "        continue\n",
    "#         print('Skip ID')\n",
    "    else: \n",
    "        if len(str(i)) == 1:\n",
    "            name = 'GOTOV0'+str(i)\n",
    "        else: \n",
    "            name = 'GOTOV'+str(i)\n",
    "        #eIF\n",
    "    #eIF\n",
    "    \n",
    "    # select data for train, val, test sets           \n",
    "    df_train, df_val, df_test = select_Data(name, randomSeed, numberForVal, data)        \n",
    "\n",
    "    time = EEm_details['timeToSplit'][EEm_details['participant']==name].values[0] #time to split in and outdoors\n",
    "    \n",
    "    cols = ['EEm', 'ID']\n",
    "    X_train = df_train.drop(cols, axis=1)\n",
    "    y_train = df_train.EEm\n",
    "\n",
    "    X_test = df_test.drop(cols, axis=1)\n",
    "    y_test = df_test.EEm\n",
    "    \n",
    "    print('train model....')\n",
    "    lm = linear_model.LinearRegression()\n",
    "    model = lm.fit(X_train,y_train)\n",
    "\n",
    "    predictions = model.predict(X_test)\n",
    "    \n",
    "    df_test = pd.DataFrame(data=y_test.values, index=y_test.index, columns=['True'])\n",
    "    df_test['Predicted'] = predictions\n",
    "    df_test['participant'] = name        \n",
    "    df_allPred = df_allPred.append(df_test)\n",
    "\n",
    "    rsquared = r2_score(y_test, predictions)\n",
    "    print('rsquared...', rsquared)\n",
    "\n",
    "    rms = sqrt(mean_squared_error(y_test, predictions))\n",
    "    print('rms...', rms)\n",
    "    \n",
    "    print('Create figures ....')\n",
    "\n",
    "    true_tot_mean = np.mean(df_test['True'])\n",
    "    pred_tot_mean = np.mean(df_test['Predicted'])        \n",
    "    true_tot_std  = np.std(df_test['True'])\n",
    "    pred_tot_std  = np.std(df_test['Predicted'])        \n",
    "\n",
    "    z = np.polyfit(df_test['True'], df_test['Predicted'], 1)\n",
    "    p = np.poly1d(z)\n",
    "    x_max = int(np.max(df_test['True'])+2)\n",
    "    y_max = int(np.max(df_test['Predicted'])+2)\n",
    "\n",
    "    plt.figure(figsize=(15,8))\n",
    "    plt.plot(df_test['True'], label='True_EE')\n",
    "    plt.plot(df_test['Predicted'], label='Predicted_EE')\n",
    "    plt.legend(loc='upper left')\n",
    "    plt.title(name)\n",
    "    plt.savefig(figures_folder+name+'.pdf')\n",
    "    # plt.show()\n",
    "    plt.close()\n",
    "    \n",
    "    if np.isnan(EEm_details['outEEm'][EEm_details['participant']==name].values) == False:\n",
    "        # indoors\n",
    "        indoors = df_test[df_test.index <= time]#'2016-02-19 11:12:25.004000']\n",
    "        indoors['participant'] = name        \n",
    "        df_allPredin = df_allPredin.append(indoors)\n",
    "\n",
    "        in_rsquared = r2_score(indoors['True'], indoors['Predicted'])\n",
    "        in_rms = sqrt(mean_squared_error(indoors['True'], indoors['Predicted']))\n",
    "        print('in :', in_rsquared)\n",
    "        true_in_mean = np.mean(indoors['True'])\n",
    "        pred_in_mean = np.mean(indoors['Predicted'])\n",
    "        true_in_std  = np.std(indoors['True'])\n",
    "        pred_in_std  = np.std(indoors['Predicted'])        \n",
    "\n",
    "        # plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.plot(indoors['True'], label='True_EE')\n",
    "        plt.plot(indoors['Predicted'], label='Predicted_EE')\n",
    "        plt.title(name)\n",
    "        plt.legend(loc='upper left')\n",
    "        plt.savefig(figures_folder+name+'_in.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        # outdoors\n",
    "        outdoors = df_test[df_test.index > time]#'2016-02-19 11:12:25.004000']\n",
    "        outdoors['participant'] = name        \n",
    "        df_allPredout = df_allPredout.append(outdoors)\n",
    "\n",
    "        out_rsquared = r2_score(outdoors['True'], outdoors['Predicted'])\n",
    "        out_rms = sqrt(mean_squared_error(outdoors['True'], outdoors['Predicted']))            \n",
    "        print('out:', out_rsquared)\n",
    "        true_out_mean = np.mean(outdoors['True'])\n",
    "        pred_out_mean = np.mean(outdoors['Predicted'])            \n",
    "        true_out_std  = np.std(outdoors['True'])\n",
    "        pred_out_std  = np.std(outdoors['Predicted'])        \n",
    "\n",
    "        # plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.title(name)\n",
    "        plt.plot(outdoors['True'], label='True_EE')\n",
    "        plt.plot(outdoors['Predicted'], label='Predicted_EE')\n",
    "        plt.legend(loc='upper left')\n",
    "        plt.savefig(figures_folder+name+'_out.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        # scatter plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.xlim([-1,x_max])\n",
    "        plt.ylim([-1,y_max])            \n",
    "        plt.scatter(indoors['True'], indoors['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "        plt.scatter(outdoors['True'], outdoors['Predicted'], c ='orange', alpha=0.5, label='Outdoors')\n",
    "        plt.plot(df_test['True'],p(df_test['True']),\"b\")\n",
    "        plt.legend(loc='upper right')\n",
    "        plt.title(name)\n",
    "        plt.xlabel('True')\n",
    "        plt.ylabel('Predicted')\n",
    "        plt.savefig(figures_folder+name+'_scatter.pdf')            \n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "    else:\n",
    "        in_rsquared   = None\n",
    "        in_rms        = None\n",
    "        out_rsquared  = None\n",
    "        out_rms       = None\n",
    "        true_in_mean  = None\n",
    "        true_out_mean = None\n",
    "        true_in_std   = None\n",
    "        true_out_std  = None\n",
    "        pred_in_mean  = None\n",
    "        pred_out_mean = None\n",
    "        pred_in_std   = None\n",
    "        pred_out_std  = None\n",
    "\n",
    "        # scatter plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.xlim([-1,x_max])\n",
    "        plt.ylim([-1,y_max])\n",
    "        plt.scatter(df_test['True'], df_test['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "        plt.plot(df_test['True'],p(df_test['True']),\"b\")            \n",
    "        plt.legend(loc='upper right')\n",
    "        plt.title(name)\n",
    "        plt.xlabel('True')\n",
    "        plt.ylabel('Predicted')\n",
    "        plt.savefig(figures_folder+name+'_scatter.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        df_allPredin = df_allPredin.append(df_test)\n",
    "    #fi\n",
    "\n",
    "    get_all_patients_rsquared = get_all_patients_rsquared.append({'participant': name,\n",
    "                                                                  'rsquared'   : rsquared,\n",
    "                                                                  'inRsquared' : in_rsquared,\n",
    "                                                                  'outRsquared': out_rsquared},\n",
    "                                                         ignore_index=True) \n",
    "\n",
    "    get_all_patients_rms = get_all_patients_rms.append({'participant': name,\n",
    "                                                        'rms'        : rms,\n",
    "                                                        'inRms'      : in_rms,\n",
    "                                                        'outRms'     : out_rms},\n",
    "                                                         ignore_index=True) \n",
    "\n",
    "    stats_EEm_results = stats_EEm_results.append({'participant'  : name,\n",
    "                                                  'true_tot_mean': true_tot_mean,\n",
    "                                                  'true_in_mean' : true_in_mean,\n",
    "                                                  'true_out_mean': true_out_mean,\n",
    "                                                  'true_tot_std' : true_tot_std,\n",
    "                                                  'true_in_std'  : true_in_std,\n",
    "                                                  'true_out_std' : true_out_std,\n",
    "                                                  'pred_tot_mean': pred_tot_mean,\n",
    "                                                  'pred_in_mean' : pred_in_mean,\n",
    "\n",
    "                                                  'pred_out_mean': pred_out_mean,\n",
    "                                                  'pred_tot_std' : pred_tot_std,\n",
    "                                                  'pred_in_std'  : pred_in_std,\n",
    "                                                  'pred_out_std' : pred_out_std},\n",
    "                                                 ignore_index=True)\n",
    "#efor\n",
    "\n",
    "get_all_patients_rsquared.to_csv(csv_folder+'rsquared.csv',index=False)\n",
    "get_all_patients_rms.to_csv(csv_folder+'rms.csv',index=False)\n",
    "stats_EEm_results.to_csv(csv_folder+'predStats.csv',index=False)\n",
    "\n",
    "# scatter\n",
    "true = stats_EEm_results['true_tot_mean']\n",
    "pred = stats_EEm_results['pred_tot_mean']\n",
    "# std_true  = stats_EEm_results['true_tot_std']\n",
    "# std_pred  = stats_EEm_results['pred_tot_std']\n",
    "# scatter\n",
    "z = np.polyfit(true, pred, 1)\n",
    "p = np.poly1d(z)\n",
    "x_max = int(np.max(true)+1)\n",
    "y_max = int(np.max(pred)+1)\n",
    "plt.figure(figsize=(15,8))\n",
    "plt.xlim([1,x_max])\n",
    "plt.ylim([1,y_max])\n",
    "plt.scatter(true, pred, c ='blue', alpha=0.5)#, label='Indoors')\n",
    "plt.plot(true,p(true),\"b\")     \n",
    "# plt.legend(loc='upper right')\n",
    "plt.title('True Vs Predicted Average EEm per participant')\n",
    "plt.xlabel('True')\n",
    "plt.ylabel('Predicted')\n",
    "plt.savefig(figures_folder+'all_scatter.pdf') \n",
    "\n",
    "print('Rsquared: '+str(get_all_patients_rsquared.rsquared.mean())+' /rms:'+ str(get_all_patients_rms.rms.mean()))\n",
    "print('inRsquared: '+str(get_all_patients_rsquared.inRsquared.mean())+' /inRms:'+ str(get_all_patients_rms.inRms.mean()))\n",
    "print('outRsquared: '+str(get_all_patients_rsquared.outRsquared.mean())+' /outRms:'+ str(get_all_patients_rms.outRms.mean()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Rsquared: '+str(get_all_patients_rsquared.rsquared.mean())+' /rms:'+ str(get_all_patients_rms.rms.mean()))\n",
    "print('inRsquared: '+str(get_all_patients_rsquared.inRsquared.mean())+' /inRms:'+ str(get_all_patients_rms.inRms.mean()))\n",
    "print('outRsquared: '+str(get_all_patients_rsquared.outRsquared.mean())+' /outRms:'+ str(get_all_patients_rms.outRms.mean()))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Linear Mixed Model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "w=30\n",
    "print(saveDir+'features_bindAll_'+str(w)+'sec.csv')\n",
    "features_bindAll = pd.read_csv(saveDir+'features_bindAll_'+str(w)+'sec.csv')\n",
    "# print(features_bindAll.columns)\n",
    "\n",
    "features_bindAll.set_index('time', inplace=True)\n",
    "features_bindAll.index = pd.to_datetime(features_bindAll.index)\n",
    "features_bindAll"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ankle_cols = []\n",
    "ankle_cols = [col for col in features_bindAll.columns if 'ankle' in col]\n",
    "ankle_cols.append('EEm')\n",
    "ankle_cols.append('ID')\n",
    "\n",
    "wrist_cols =[]\n",
    "wrist_cols = [col for col in features_bindAll.columns if 'wrist' in col]\n",
    "wrist_cols.append('EEm')\n",
    "wrist_cols.append('ID')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "devices = 'wrist' # aw or ankle or wrist\n",
    "\n",
    "if devices == 'ankle':\n",
    "    features_bindAll = features_bindAll[ankle_cols]\n",
    "elif devices == 'wrist':\n",
    "    features_bindAll = features_bindAll[wrist_cols]\n",
    "features_bindAll"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "EEm_details = pd.read_csv('EEm_details_per_partic.csv')  \n",
    "\n",
    "\n",
    "numberForVal = 1 # the number of participants you want per gtoup as validation set (random seed is the number of participants' id)\n",
    "# randomSeed = i\n",
    "\n",
    "data = features_bindAll\n",
    "\n",
    "model_name = 'LMM_'+devices+'_'+str(w)+'sec'\n",
    "figures_folder = 'images/' + model_name + '/predPlot/'\n",
    "csv_folder = figures_folder\n",
    "\n",
    "if not os.path.exists(figures_folder):\n",
    "    os.makedirs(figures_folder)\n",
    "\n",
    "get_all_patients_rsquared = pd.DataFrame(columns=['participant', 'rsquared', 'inRsquared', 'outRsquared'])\n",
    "get_all_patients_rms = pd.DataFrame(columns=['participant', 'rms', 'inRms', 'outRms'])\n",
    "stats_EEm_results = pd.DataFrame(columns=['participant', \n",
    "                                          'true_tot_mean', 'true_in_mean', 'true_out_mean', \n",
    "                                          'pred_tot_mean', 'pred_in_mean', 'pred_out_mean'])\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# build formula\n",
    "if devices == 'ankle':\n",
    "    fmla ='EEm ~ ankle_x_mean' \n",
    "    for vname in data.columns[0:30].values:\n",
    "        fmla = fmla + ' + ' +str(vname)\n",
    "elif devices == 'wrist':\n",
    "    fmla ='EEm ~ wrist_x_mean' \n",
    "    for vname in data.columns[0:30].values:\n",
    "        fmla = fmla + ' + ' +str(vname)\n",
    "elif devices == 'aw':\n",
    "    fmla ='EEm ~ ankle_x_mean' \n",
    "    for vname in data.columns[1:60].values:\n",
    "        fmla = fmla + ' + ' +str(vname)\n",
    "   \n",
    "fmla"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "df_allPred = pd.DataFrame()\n",
    "df_allPredin = pd.DataFrame()\n",
    "df_allPredout = pd.DataFrame()\n",
    "\n",
    "for i in range(5, 37):\n",
    "    randomSeed = i\n",
    "    if i in [2,3,4,19]: #[12, 19, 23, 16]: \n",
    "        continue\n",
    "#         print('Skip ID')\n",
    "    else: \n",
    "        if len(str(i)) == 1:\n",
    "            name = 'GOTOV0'+str(i)\n",
    "        else: \n",
    "            name = 'GOTOV'+str(i)\n",
    "        #eIF\n",
    "    #eIF\n",
    "    \n",
    "    # select data for train, val, test sets           \n",
    "    df_train, df_val, df_test = select_Data(name, randomSeed, numberForVal, data)        \n",
    "\n",
    "    time = EEm_details['timeToSplit'][EEm_details['participant']==name].values[0] #time to split in and outdoors\n",
    "    \n",
    "    cols = ['EEm', 'ID']\n",
    "    X_train = df_train.drop(cols, axis=1)\n",
    "    y_train = df_train.EEm\n",
    "\n",
    "    X_test = df_test.drop(cols, axis=1)\n",
    "    y_test = df_test.EEm\n",
    "    \n",
    "    print('train model....')\n",
    "    llm = smf.mixedlm(str(fmla), df_train, groups=df_train[\"ID\"])\n",
    "    model = llm.fit()\n",
    "#     print(model.summary())\n",
    "\n",
    "#     lm = linear_model.LinearRegression()\n",
    "#     model = lm.fit(X_train,y_train)\n",
    "\n",
    "    predictions = model.predict(X_test)\n",
    "    \n",
    "    df_test = pd.DataFrame(data=y_test.values, index=y_test.index, columns=['True'])\n",
    "    df_test['Predicted'] = predictions\n",
    "    df_test['participant'] = name        \n",
    "    df_allPred = df_allPred.append(df_test)\n",
    "\n",
    "    rsquared = r2_score(y_test, predictions)\n",
    "    print('rsquared...', rsquared)\n",
    "\n",
    "    rms = sqrt(mean_squared_error(y_test, predictions))\n",
    "    print('rms...', rms)\n",
    "    \n",
    "    print('Create figures ....')\n",
    "\n",
    "    true_tot_mean = np.mean(df_test['True'])\n",
    "    pred_tot_mean = np.mean(df_test['Predicted'])        \n",
    "    true_tot_std  = np.std(df_test['True'])\n",
    "    pred_tot_std  = np.std(df_test['Predicted'])        \n",
    "\n",
    "    z = np.polyfit(df_test['True'], df_test['Predicted'], 1)\n",
    "    p = np.poly1d(z)\n",
    "    x_max = int(np.max(df_test['True'])+2)\n",
    "    y_max = int(np.max(df_test['Predicted'])+2)\n",
    "\n",
    "    plt.figure(figsize=(15,8))\n",
    "    plt.plot(df_test['True'], label='True_EE')\n",
    "    plt.plot(df_test['Predicted'], label='Predicted_EE')\n",
    "    plt.legend(loc='upper left')\n",
    "    plt.title(name)\n",
    "    plt.savefig(figures_folder+name+'.pdf')\n",
    "    # plt.show()\n",
    "    plt.close()\n",
    "    \n",
    "    if np.isnan(EEm_details['outEEm'][EEm_details['participant']==name].values) == False:\n",
    "        # indoors\n",
    "        indoors = df_test[df_test.index <= time]#'2016-02-19 11:12:25.004000']\n",
    "        indoors['participant'] = name        \n",
    "        df_allPredin = df_allPredin.append(indoors)\n",
    "\n",
    "        in_rsquared = r2_score(indoors['True'], indoors['Predicted'])\n",
    "        in_rms = sqrt(mean_squared_error(indoors['True'], indoors['Predicted']))\n",
    "        print('in :', in_rsquared)\n",
    "        true_in_mean = np.mean(indoors['True'])\n",
    "        pred_in_mean = np.mean(indoors['Predicted'])\n",
    "        true_in_std  = np.std(indoors['True'])\n",
    "        pred_in_std  = np.std(indoors['Predicted'])        \n",
    "\n",
    "        # plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.plot(indoors['True'], label='True_EE')\n",
    "        plt.plot(indoors['Predicted'], label='Predicted_EE')\n",
    "        plt.title(name)\n",
    "        plt.legend(loc='upper left')\n",
    "        plt.savefig(figures_folder+name+'_in.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        # outdoors\n",
    "        outdoors = df_test[df_test.index > time]#'2016-02-19 11:12:25.004000']\n",
    "        outdoors['participant'] = name        \n",
    "        df_allPredout = df_allPredout.append(outdoors)\n",
    "\n",
    "        out_rsquared = r2_score(outdoors['True'], outdoors['Predicted'])\n",
    "        out_rms = sqrt(mean_squared_error(outdoors['True'], outdoors['Predicted']))            \n",
    "        print('out:', out_rsquared)\n",
    "        true_out_mean = np.mean(outdoors['True'])\n",
    "        pred_out_mean = np.mean(outdoors['Predicted'])            \n",
    "        true_out_std  = np.std(outdoors['True'])\n",
    "        pred_out_std  = np.std(outdoors['Predicted'])        \n",
    "\n",
    "        # plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.title(name)\n",
    "        plt.plot(outdoors['True'], label='True_EE')\n",
    "        plt.plot(outdoors['Predicted'], label='Predicted_EE')\n",
    "        plt.legend(loc='upper left')\n",
    "        plt.savefig(figures_folder+name+'_out.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        # scatter plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.xlim([-1,x_max])\n",
    "        plt.ylim([-1,y_max])            \n",
    "        plt.scatter(indoors['True'], indoors['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "        plt.scatter(outdoors['True'], outdoors['Predicted'], c ='orange', alpha=0.5, label='Outdoors')\n",
    "        plt.plot(df_test['True'],p(df_test['True']),\"b\")\n",
    "        plt.legend(loc='upper right')\n",
    "        plt.title(name)\n",
    "        plt.xlabel('True')\n",
    "        plt.ylabel('Predicted')\n",
    "        plt.savefig(figures_folder+name+'_scatter.pdf')            \n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "    else:\n",
    "        in_rsquared   = None\n",
    "        in_rms        = None\n",
    "        out_rsquared  = None\n",
    "        out_rms       = None\n",
    "        true_in_mean  = None\n",
    "        true_out_mean = None\n",
    "        true_in_std   = None\n",
    "        true_out_std  = None\n",
    "        pred_in_mean  = None\n",
    "        pred_out_mean = None\n",
    "        pred_in_std   = None\n",
    "        pred_out_std  = None\n",
    "\n",
    "        # scatter plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.xlim([-1,x_max])\n",
    "        plt.ylim([-1,y_max])\n",
    "        plt.scatter(df_test['True'], df_test['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "        plt.plot(df_test['True'],p(df_test['True']),\"b\")            \n",
    "        plt.legend(loc='upper right')\n",
    "        plt.title(name)\n",
    "        plt.xlabel('True')\n",
    "        plt.ylabel('Predicted')\n",
    "        plt.savefig(figures_folder+name+'_scatter.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        df_allPredin = df_allPredin.append(df_test)\n",
    "    #fi\n",
    "\n",
    "    get_all_patients_rsquared = get_all_patients_rsquared.append({'participant': name,\n",
    "                                                                  'rsquared'   : rsquared,\n",
    "                                                                  'inRsquared' : in_rsquared,\n",
    "                                                                  'outRsquared': out_rsquared},\n",
    "                                                         ignore_index=True) \n",
    "\n",
    "    get_all_patients_rms = get_all_patients_rms.append({'participant': name,\n",
    "                                                        'rms'        : rms,\n",
    "                                                        'inRms'      : in_rms,\n",
    "                                                        'outRms'     : out_rms},\n",
    "                                                         ignore_index=True) \n",
    "\n",
    "    stats_EEm_results = stats_EEm_results.append({'participant'  : name,\n",
    "                                                  'true_tot_mean': true_tot_mean,\n",
    "                                                  'true_in_mean' : true_in_mean,\n",
    "                                                  'true_out_mean': true_out_mean,\n",
    "                                                  'true_tot_std' : true_tot_std,\n",
    "                                                  'true_in_std'  : true_in_std,\n",
    "                                                  'true_out_std' : true_out_std,\n",
    "                                                  'pred_tot_mean': pred_tot_mean,\n",
    "                                                  'pred_in_mean' : pred_in_mean,\n",
    "\n",
    "                                                  'pred_out_mean': pred_out_mean,\n",
    "                                                  'pred_tot_std' : pred_tot_std,\n",
    "                                                  'pred_in_std'  : pred_in_std,\n",
    "                                                  'pred_out_std' : pred_out_std},\n",
    "                                                 ignore_index=True)\n",
    "#efor\n",
    "\n",
    "get_all_patients_rsquared.to_csv(csv_folder+'rsquared.csv',index=False)\n",
    "get_all_patients_rms.to_csv(csv_folder+'rms.csv',index=False)\n",
    "stats_EEm_results.to_csv(csv_folder+'predStats.csv',index=False)\n",
    "\n",
    "# scatter\n",
    "true = stats_EEm_results['true_tot_mean']\n",
    "pred = stats_EEm_results['pred_tot_mean']\n",
    "# std_true  = stats_EEm_results['true_tot_std']\n",
    "# std_pred  = stats_EEm_results['pred_tot_std']\n",
    "# scatter\n",
    "z = np.polyfit(true, pred, 1)\n",
    "p = np.poly1d(z)\n",
    "x_max = int(np.max(true)+1)\n",
    "y_max = int(np.max(pred)+1)\n",
    "plt.figure(figsize=(15,8))\n",
    "plt.xlim([1,x_max])\n",
    "plt.ylim([1,y_max])\n",
    "plt.scatter(true, pred, c ='blue', alpha=0.5)#, label='Indoors')\n",
    "plt.plot(true,p(true),\"b\")     \n",
    "# plt.legend(loc='upper right')\n",
    "plt.title('True Vs Predicted Average EEm per participant')\n",
    "plt.xlabel('True')\n",
    "plt.ylabel('Predicted')\n",
    "plt.savefig(figures_folder+'all_scatter.pdf') \n",
    "\n",
    "print('Rsquared: '+str(get_all_patients_rsquared.rsquared.mean())+' /rms:'+ str(get_all_patients_rms.rms.mean()))\n",
    "print('inRsquared: '+str(get_all_patients_rsquared.inRsquared.mean())+' /inRms:'+ str(get_all_patients_rms.inRms.mean()))\n",
    "print('outRsquared: '+str(get_all_patients_rsquared.outRsquared.mean())+' /outRms:'+ str(get_all_patients_rms.outRms.mean()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Rsquared: '+str(get_all_patients_rsquared.rsquared.mean())+' /rms:'+ str(get_all_patients_rms.rms.mean()))\n",
    "print('inRsquared: '+str(get_all_patients_rsquared.inRsquared.mean())+' /inRms:'+ str(get_all_patients_rms.inRms.mean()))\n",
    "print('outRsquared: '+str(get_all_patients_rsquared.outRsquared.mean())+' /outRms:'+ str(get_all_patients_rms.outRms.mean()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# RandomForest regressor"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "w=30\n",
    "# print('')\n",
    "features_bindAll = pd.read_csv('sequences/features_MontoyePaper_replication/features_bindAll_30sec.csv')\n",
    "# print(features_bindAll.columns)\n",
    "\n",
    "features_bindAll.set_index('time', inplace=True)\n",
    "features_bindAll.index = pd.to_datetime(features_bindAll.index)\n",
    "features_bindAll"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ankle_cols = []\n",
    "ankle_cols = [col for col in features_bindAll.columns if 'ankle' in col]\n",
    "ankle_cols.append('EEm')\n",
    "ankle_cols.append('ID')\n",
    "\n",
    "wrist_cols =[]\n",
    "wrist_cols = [col for col in features_bindAll.columns if 'wrist' in col]\n",
    "wrist_cols.append('EEm')\n",
    "wrist_cols.append('ID')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "devices = 'aw' # aw or ankle or wrist\n",
    "\n",
    "if devices == 'ankle':\n",
    "    features_bindAll = features_bindAll[ankle_cols]\n",
    "elif devices == 'wrist':\n",
    "    features_bindAll = features_bindAll[wrist_cols]\n",
    "features_bindAll"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Fitting Random Forest Regression to the dataset\n",
    "# import the regressor\n",
    "from sklearn.ensemble import RandomForestRegressor\n",
    "  \n",
    "\n",
    "EEm_details = pd.read_csv('EEm_details_per_partic.csv')  \n",
    "\n",
    "\n",
    "numberForVal = 1 # the number of participants you want per gtoup as validation set (random seed is the number of participants' id)\n",
    "# randomSeed = i\n",
    "\n",
    "data = features_bindAll\n",
    "\n",
    "model_name = 'RFregr_1000trees_'+devices+'_'+str(w)+'sec'\n",
    "figures_folder = 'images/' + model_name + '/predPlot/'\n",
    "csv_folder = figures_folder\n",
    "\n",
    "if not os.path.exists(figures_folder):\n",
    "    os.makedirs(figures_folder)\n",
    "\n",
    "get_all_patients_rsquared = pd.DataFrame(columns=['participant', 'rsquared', 'inRsquared', 'outRsquared'])\n",
    "get_all_patients_rms = pd.DataFrame(columns=['participant', 'rms', 'inRms', 'outRms'])\n",
    "stats_EEm_results = pd.DataFrame(columns=['participant', \n",
    "                                          'true_tot_mean', 'true_in_mean', 'true_out_mean', \n",
    "                                          'pred_tot_mean', 'pred_in_mean', 'pred_out_mean'])\n",
    "\n",
    "print(model_name)\n",
    "print(figures_folder)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "df_allPred = pd.DataFrame()\n",
    "df_allPredin = pd.DataFrame()\n",
    "df_allPredout = pd.DataFrame()\n",
    "\n",
    "names =[5]\n",
    "# for i in names:\n",
    "\n",
    "for i in range(5, 37):\n",
    "    randomSeed = i\n",
    "    if i in [2,3,4,19]: #[12, 19, 23, 16]: \n",
    "        continue\n",
    "#         print('Skip ID')\n",
    "    else: \n",
    "        if len(str(i)) == 1:\n",
    "            name = 'GOTOV0'+str(i)\n",
    "        else: \n",
    "            name = 'GOTOV'+str(i)\n",
    "        #eIF\n",
    "    #eIF\n",
    "    \n",
    "    # select data for train, val, test sets           \n",
    "    df_train, df_val, df_test = select_Data(name, randomSeed, numberForVal, data)        \n",
    "\n",
    "    time = EEm_details['timeToSplit'][EEm_details['participant']==name].values[0] #time to split in and outdoors\n",
    "    \n",
    "    cols = ['EEm', 'ID']\n",
    "    X_train = df_train.drop(cols, axis=1)\n",
    "    y_train = df_train.EEm\n",
    "\n",
    "    X_test = df_test.drop(cols, axis=1)\n",
    "    y_test = df_test.EEm\n",
    "    \n",
    "    print('train model....')\n",
    "#   create regressor object\n",
    "    regressor = RandomForestRegressor(n_estimators = 1000, random_state = 0)\n",
    "#   fit the regressor with x and y data\n",
    "    model = regressor.fit(X_train, y_train)  \n",
    "    pickle.dump(model, open('models/2nd_round_of_models_paper_rebuttal/RFregr_1000trees_aw_30sec/'+str(name), 'wb'))\n",
    "\n",
    "    predictions = model.predict(X_test)\n",
    "    \n",
    "    df_test = pd.DataFrame(data=y_test.values, index=y_test.index, columns=['True'])\n",
    "    df_test['Predicted'] = predictions\n",
    "    df_test['participant'] = name        \n",
    "    df_allPred = df_allPred.append(df_test)\n",
    "\n",
    "    rsquared = r2_score(y_test, predictions)\n",
    "    print('rsquared...', rsquared)\n",
    "\n",
    "    rms = sqrt(mean_squared_error(y_test, predictions))\n",
    "    print('rms...', rms)\n",
    "    \n",
    "    print('Create figures ....')\n",
    "\n",
    "    true_tot_mean = np.mean(df_test['True'])\n",
    "    pred_tot_mean = np.mean(df_test['Predicted'])        \n",
    "    true_tot_std  = np.std(df_test['True'])\n",
    "    pred_tot_std  = np.std(df_test['Predicted'])        \n",
    "\n",
    "    z = np.polyfit(df_test['True'], df_test['Predicted'], 1)\n",
    "    p = np.poly1d(z)\n",
    "    x_max = int(np.max(df_test['True'])+2)\n",
    "    y_max = int(np.max(df_test['Predicted'])+2)\n",
    "\n",
    "    plt.figure(figsize=(15,8))\n",
    "    plt.plot(df_test['True'], label='True_EE')\n",
    "    plt.plot(df_test['Predicted'], label='Predicted_EE')\n",
    "    plt.legend(loc='upper left')\n",
    "    plt.title(name)\n",
    "    plt.savefig(figures_folder+name+'.pdf')\n",
    "    # plt.show()\n",
    "    plt.close()\n",
    "    \n",
    "    if np.isnan(EEm_details['outEEm'][EEm_details['participant']==name].values) == False:\n",
    "        # indoors\n",
    "        indoors = df_test[df_test.index <= time]#'2016-02-19 11:12:25.004000']\n",
    "        indoors['participant'] = name        \n",
    "        df_allPredin = df_allPredin.append(indoors)\n",
    "\n",
    "        in_rsquared = r2_score(indoors['True'], indoors['Predicted'])\n",
    "        in_rms = sqrt(mean_squared_error(indoors['True'], indoors['Predicted']))\n",
    "        print('in :', in_rsquared)\n",
    "        true_in_mean = np.mean(indoors['True'])\n",
    "        pred_in_mean = np.mean(indoors['Predicted'])\n",
    "        true_in_std  = np.std(indoors['True'])\n",
    "        pred_in_std  = np.std(indoors['Predicted'])        \n",
    "\n",
    "        # plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.plot(indoors['True'], label='True_EE')\n",
    "        plt.plot(indoors['Predicted'], label='Predicted_EE')\n",
    "        plt.title(name)\n",
    "        plt.legend(loc='upper left')\n",
    "        plt.savefig(figures_folder+name+'_in.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        # outdoors\n",
    "        outdoors = df_test[df_test.index > time]#'2016-02-19 11:12:25.004000']\n",
    "        outdoors['participant'] = name        \n",
    "        df_allPredout = df_allPredout.append(outdoors)\n",
    "\n",
    "        out_rsquared = r2_score(outdoors['True'], outdoors['Predicted'])\n",
    "        out_rms = sqrt(mean_squared_error(outdoors['True'], outdoors['Predicted']))            \n",
    "        print('out:', out_rsquared)\n",
    "        true_out_mean = np.mean(outdoors['True'])\n",
    "        pred_out_mean = np.mean(outdoors['Predicted'])            \n",
    "        true_out_std  = np.std(outdoors['True'])\n",
    "        pred_out_std  = np.std(outdoors['Predicted'])        \n",
    "\n",
    "        # plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.title(name)\n",
    "        plt.plot(outdoors['True'], label='True_EE')\n",
    "        plt.plot(outdoors['Predicted'], label='Predicted_EE')\n",
    "        plt.legend(loc='upper left')\n",
    "        plt.savefig(figures_folder+name+'_out.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        # scatter plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.xlim([-1,x_max])\n",
    "        plt.ylim([-1,y_max])            \n",
    "        plt.scatter(indoors['True'], indoors['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "        plt.scatter(outdoors['True'], outdoors['Predicted'], c ='orange', alpha=0.5, label='Outdoors')\n",
    "        plt.plot(df_test['True'],p(df_test['True']),\"b\")\n",
    "        plt.legend(loc='upper right')\n",
    "        plt.title(name)\n",
    "        plt.xlabel('True')\n",
    "        plt.ylabel('Predicted')\n",
    "        plt.savefig(figures_folder+name+'_scatter.pdf')            \n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "    else:\n",
    "        in_rsquared   = None\n",
    "        in_rms        = None\n",
    "        out_rsquared  = None\n",
    "        out_rms       = None\n",
    "        true_in_mean  = None\n",
    "        true_out_mean = None\n",
    "        true_in_std   = None\n",
    "        true_out_std  = None\n",
    "        pred_in_mean  = None\n",
    "        pred_out_mean = None\n",
    "        pred_in_std   = None\n",
    "        pred_out_std  = None\n",
    "\n",
    "        # scatter plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.xlim([-1,x_max])\n",
    "        plt.ylim([-1,y_max])\n",
    "        plt.scatter(df_test['True'], df_test['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "        plt.plot(df_test['True'],p(df_test['True']),\"b\")            \n",
    "        plt.legend(loc='upper right')\n",
    "        plt.title(name)\n",
    "        plt.xlabel('True')\n",
    "        plt.ylabel('Predicted')\n",
    "        plt.savefig(figures_folder+name+'_scatter.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        df_allPredin = df_allPredin.append(df_test)\n",
    "    #fi\n",
    "\n",
    "    get_all_patients_rsquared = get_all_patients_rsquared.append({'participant': name,\n",
    "                                                                  'rsquared'   : rsquared,\n",
    "                                                                  'inRsquared' : in_rsquared,\n",
    "                                                                  'outRsquared': out_rsquared},\n",
    "                                                         ignore_index=True) \n",
    "\n",
    "    get_all_patients_rms = get_all_patients_rms.append({'participant': name,\n",
    "                                                        'rms'        : rms,\n",
    "                                                        'inRms'      : in_rms,\n",
    "                                                        'outRms'     : out_rms},\n",
    "                                                         ignore_index=True) \n",
    "\n",
    "    stats_EEm_results = stats_EEm_results.append({'participant'  : name,\n",
    "                                                  'true_tot_mean': true_tot_mean,\n",
    "                                                  'true_in_mean' : true_in_mean,\n",
    "                                                  'true_out_mean': true_out_mean,\n",
    "                                                  'true_tot_std' : true_tot_std,\n",
    "                                                  'true_in_std'  : true_in_std,\n",
    "                                                  'true_out_std' : true_out_std,\n",
    "                                                  'pred_tot_mean': pred_tot_mean,\n",
    "                                                  'pred_in_mean' : pred_in_mean,\n",
    "\n",
    "                                                  'pred_out_mean': pred_out_mean,\n",
    "                                                  'pred_tot_std' : pred_tot_std,\n",
    "                                                  'pred_in_std'  : pred_in_std,\n",
    "                                                  'pred_out_std' : pred_out_std},\n",
    "                                                 ignore_index=True)\n",
    "efor\n",
    "\n",
    "get_all_patients_rsquared.to_csv(csv_folder+'rsquared.csv',index=False)\n",
    "get_all_patients_rms.to_csv(csv_folder+'rms.csv',index=False)\n",
    "stats_EEm_results.to_csv(csv_folder+'predStats.csv',index=False)\n",
    "\n",
    "# scatter\n",
    "true = stats_EEm_results['true_tot_mean']\n",
    "pred = stats_EEm_results['pred_tot_mean']\n",
    "# std_true  = stats_EEm_results['true_tot_std']\n",
    "# std_pred  = stats_EEm_results['pred_tot_std']\n",
    "# scatter\n",
    "z = np.polyfit(true, pred, 1)\n",
    "p = np.poly1d(z)\n",
    "x_max = int(np.max(true)+1)\n",
    "y_max = int(np.max(pred)+1)\n",
    "plt.figure(figsize=(15,8))\n",
    "plt.xlim([1,x_max])\n",
    "plt.ylim([1,y_max])\n",
    "plt.scatter(true, pred, c ='blue', alpha=0.5)#, label='Indoors')\n",
    "plt.plot(true,p(true),\"b\")     \n",
    "# plt.legend(loc='upper right')\n",
    "plt.title('True Vs Predicted Average EEm per participant')\n",
    "plt.xlabel('True')\n",
    "plt.ylabel('Predicted')\n",
    "plt.savefig(figures_folder+'all_scatter.pdf') \n",
    "\n",
    "print('Rsquared: '+str(get_all_patients_rsquared.rsquared.mean())+' /rms:'+ str(get_all_patients_rms.rms.mean()))\n",
    "print('inRsquared: '+str(get_all_patients_rsquared.inRsquared.mean())+' /inRms:'+ str(get_all_patients_rms.inRms.mean()))\n",
    "print('outRsquared: '+str(get_all_patients_rsquared.outRsquared.mean())+' /outRms:'+ str(get_all_patients_rms.outRms.mean()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Rsquared: '+str(get_all_patients_rsquared.rsquared.median())+' /rms:'+ str(get_all_patients_rms.rms.median()))\n",
    "print('inRsquared: '+str(get_all_patients_rsquared.inRsquared.median())+' /inRms:'+ str(get_all_patients_rms.inRms.median()))\n",
    "print('outRsquared: '+str(get_all_patients_rsquared.outRsquared.median())+' /outRms:'+ str(get_all_patients_rms.outRms.median()))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# nnet model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "EEm_details = pd.read_csv('EEm_details_per_partic.csv')  \n",
    "\n",
    "# numberForVal = 1 # the number of participants you want per gtoup as validation set (random seed is the number of participants' id)\n",
    "# randomSeed = i\n",
    "w =30\n",
    "devices = 'w' # aw or a or w\n",
    "# data = features_bindAll\n",
    "dataToUse = devices+'/'+str(w)+'sec'\n",
    "model_name = 'nnet_'+dataToUse\n",
    "figures_folder = 'images/' + model_name + '/predPlot/'\n",
    "csv_folder = figures_folder\n",
    "\n",
    "if not os.path.exists(figures_folder):\n",
    "    os.makedirs(figures_folder)\n",
    "\n",
    "get_all_patients_rsquared = pd.DataFrame(columns=['participant', 'rsquared', 'inRsquared', 'outRsquared'])\n",
    "get_all_patients_rms = pd.DataFrame(columns=['participant', 'rms', 'inRms', 'outRms'])\n",
    "stats_EEm_results = pd.DataFrame(columns=['participant', \n",
    "                                          'true_tot_mean', 'true_in_mean', 'true_out_mean', \n",
    "                                          'pred_tot_mean', 'pred_in_mean', 'pred_out_mean'])\n",
    "print(model_name)\n",
    "print(figures_folder)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "df_allPred = pd.DataFrame()\n",
    "df_allPredin = pd.DataFrame()\n",
    "df_allPredout = pd.DataFrame()\n",
    "\n",
    "for i in range(5, 37):\n",
    "    randomSeed = i\n",
    "    if i in [2,3,4,19]: #[12, 19, 23, 16]: \n",
    "        continue\n",
    "#         print('Skip ID')\n",
    "    else: \n",
    "        if len(str(i)) == 1:\n",
    "            name = 'GOTOV0'+str(i)\n",
    "        else: \n",
    "            name = 'GOTOV'+str(i)\n",
    "        #eIF\n",
    "    #eIF\n",
    "    time = EEm_details['timeToSplit'][EEm_details['participant']==name].values[0] #time to split in and outdoors\n",
    "    \n",
    "    # load predictions of R models (they were trained in Shark cluster using another R-Jupyter notebook)\n",
    "    print('read predictions:', name)\n",
    "    \n",
    "    df_test = pd.read_csv('ANNs_predictions_Montoye_repl/nnet/'+dataToUse+'/'+name+\"_preds.csv\", header = 0, index_col = None, low_memory=False)\n",
    "    df_test.columns = ['time', 'True', 'Predicted']\n",
    "    df_test.set_index('time', inplace=True)\n",
    "    df_test.index = pd.to_datetime(df_test.index)\n",
    "\n",
    "    df_test['participant'] = name        \n",
    "    \n",
    "    df_allPred = df_allPred.append(df_test)\n",
    "\n",
    "    rsquared = r2_score(df_test['True'], df_test['Predicted'])\n",
    "    print('rsquared...', rsquared)\n",
    "\n",
    "    rms = sqrt(mean_squared_error(df_test['True'], df_test['Predicted']))\n",
    "    print('rms...', rms)\n",
    "    \n",
    "    print('Create figures ....')\n",
    "\n",
    "    true_tot_mean = np.mean(df_test['True'])\n",
    "    pred_tot_mean = np.mean(df_test['Predicted'])        \n",
    "    true_tot_std  = np.std(df_test['True'])\n",
    "    pred_tot_std  = np.std(df_test['Predicted'])        \n",
    "\n",
    "    z = np.polyfit(df_test['True'], df_test['Predicted'], 1)\n",
    "    p = np.poly1d(z)\n",
    "    x_max = int(np.max(df_test['True'])+2)\n",
    "    y_max = int(np.max(df_test['Predicted'])+2)\n",
    "\n",
    "    plt.figure(figsize=(15,8))\n",
    "    plt.plot(df_test['True'], label='True_EE')\n",
    "    plt.plot(df_test['Predicted'], label='Predicted_EE')\n",
    "    plt.legend(loc='upper left')\n",
    "    plt.title(name)\n",
    "    plt.savefig(figures_folder+name+'.pdf')\n",
    "    # plt.show()\n",
    "    plt.close()\n",
    "    \n",
    "    if np.isnan(EEm_details['outEEm'][EEm_details['participant']==name].values) == False:\n",
    "        # indoors\n",
    "        indoors = df_test[df_test.index <= time]#'2016-02-19 11:12:25.004000']\n",
    "        indoors['participant'] = name        \n",
    "        df_allPredin = df_allPredin.append(indoors)\n",
    "\n",
    "        in_rsquared = r2_score(indoors['True'], indoors['Predicted'])\n",
    "        in_rms = sqrt(mean_squared_error(indoors['True'], indoors['Predicted']))\n",
    "        print('in :', in_rsquared)\n",
    "        true_in_mean = np.mean(indoors['True'])\n",
    "        pred_in_mean = np.mean(indoors['Predicted'])\n",
    "        true_in_std  = np.std(indoors['True'])\n",
    "        pred_in_std  = np.std(indoors['Predicted'])        \n",
    "\n",
    "        # plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.plot(indoors['True'], label='True_EE')\n",
    "        plt.plot(indoors['Predicted'], label='Predicted_EE')\n",
    "        plt.title(name)\n",
    "        plt.legend(loc='upper left')\n",
    "        plt.savefig(figures_folder+name+'_in.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        # outdoors\n",
    "        outdoors = df_test[df_test.index > time]#'2016-02-19 11:12:25.004000']\n",
    "        outdoors['participant'] = name        \n",
    "        df_allPredout = df_allPredout.append(outdoors)\n",
    "\n",
    "        out_rsquared = r2_score(outdoors['True'], outdoors['Predicted'])\n",
    "        out_rms = sqrt(mean_squared_error(outdoors['True'], outdoors['Predicted']))            \n",
    "        print('out:', out_rsquared)\n",
    "        true_out_mean = np.mean(outdoors['True'])\n",
    "        pred_out_mean = np.mean(outdoors['Predicted'])            \n",
    "        true_out_std  = np.std(outdoors['True'])\n",
    "        pred_out_std  = np.std(outdoors['Predicted'])        \n",
    "\n",
    "        # plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.title(name)\n",
    "        plt.plot(outdoors['True'], label='True_EE')\n",
    "        plt.plot(outdoors['Predicted'], label='Predicted_EE')\n",
    "        plt.legend(loc='upper left')\n",
    "        plt.savefig(figures_folder+name+'_out.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        # scatter plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.xlim([-1,x_max])\n",
    "        plt.ylim([-1,y_max])            \n",
    "        plt.scatter(indoors['True'], indoors['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "        plt.scatter(outdoors['True'], outdoors['Predicted'], c ='orange', alpha=0.5, label='Outdoors')\n",
    "        plt.plot(df_test['True'],p(df_test['True']),\"b\")\n",
    "        plt.legend(loc='upper right')\n",
    "        plt.title(name)\n",
    "        plt.xlabel('True')\n",
    "        plt.ylabel('Predicted')\n",
    "        plt.savefig(figures_folder+name+'_scatter.pdf')            \n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "    else:\n",
    "        in_rsquared   = None\n",
    "        in_rms        = None\n",
    "        out_rsquared  = None\n",
    "        out_rms       = None\n",
    "        true_in_mean  = None\n",
    "        true_out_mean = None\n",
    "        true_in_std   = None\n",
    "        true_out_std  = None\n",
    "        pred_in_mean  = None\n",
    "        pred_out_mean = None\n",
    "        pred_in_std   = None\n",
    "        pred_out_std  = None\n",
    "\n",
    "        # scatter plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.xlim([-1,x_max])\n",
    "        plt.ylim([-1,y_max])\n",
    "        plt.scatter(df_test['True'], df_test['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "        plt.plot(df_test['True'],p(df_test['True']),\"b\")            \n",
    "        plt.legend(loc='upper right')\n",
    "        plt.title(name)\n",
    "        plt.xlabel('True')\n",
    "        plt.ylabel('Predicted')\n",
    "        plt.savefig(figures_folder+name+'_scatter.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        df_allPredin = df_allPredin.append(df_test)\n",
    "    #fi\n",
    "\n",
    "    get_all_patients_rsquared = get_all_patients_rsquared.append({'participant': name,\n",
    "                                                                  'rsquared'   : rsquared,\n",
    "                                                                  'inRsquared' : in_rsquared,\n",
    "                                                                  'outRsquared': out_rsquared},\n",
    "                                                         ignore_index=True) \n",
    "\n",
    "    get_all_patients_rms = get_all_patients_rms.append({'participant': name,\n",
    "                                                        'rms'        : rms,\n",
    "                                                        'inRms'      : in_rms,\n",
    "                                                        'outRms'     : out_rms},\n",
    "                                                         ignore_index=True) \n",
    "\n",
    "    stats_EEm_results = stats_EEm_results.append({'participant'  : name,\n",
    "                                                  'true_tot_mean': true_tot_mean,\n",
    "                                                  'true_in_mean' : true_in_mean,\n",
    "                                                  'true_out_mean': true_out_mean,\n",
    "                                                  'true_tot_std' : true_tot_std,\n",
    "                                                  'true_in_std'  : true_in_std,\n",
    "                                                  'true_out_std' : true_out_std,\n",
    "                                                  'pred_tot_mean': pred_tot_mean,\n",
    "                                                  'pred_in_mean' : pred_in_mean,\n",
    "\n",
    "                                                  'pred_out_mean': pred_out_mean,\n",
    "                                                  'pred_tot_std' : pred_tot_std,\n",
    "                                                  'pred_in_std'  : pred_in_std,\n",
    "                                                  'pred_out_std' : pred_out_std},\n",
    "                                                 ignore_index=True)\n",
    "#efor\n",
    "\n",
    "get_all_patients_rsquared.to_csv(csv_folder+'rsquared.csv',index=False)\n",
    "get_all_patients_rms.to_csv(csv_folder+'rms.csv',index=False)\n",
    "stats_EEm_results.to_csv(csv_folder+'predStats.csv',index=False)\n",
    "\n",
    "# scatter\n",
    "true = stats_EEm_results['true_tot_mean']\n",
    "pred = stats_EEm_results['pred_tot_mean']\n",
    "# std_true  = stats_EEm_results['true_tot_std']\n",
    "# std_pred  = stats_EEm_results['pred_tot_std']\n",
    "# scatter\n",
    "z = np.polyfit(true, pred, 1)\n",
    "p = np.poly1d(z)\n",
    "x_max = int(np.max(true)+1)\n",
    "y_max = int(np.max(pred)+1)\n",
    "plt.figure(figsize=(15,8))\n",
    "plt.xlim([1,x_max])\n",
    "plt.ylim([1,y_max])\n",
    "plt.scatter(true, pred, c ='blue', alpha=0.5)#, label='Indoors')\n",
    "plt.plot(true,p(true),\"b\")     \n",
    "# plt.legend(loc='upper right')\n",
    "plt.title('True Vs Predicted Average EEm per participant')\n",
    "plt.xlabel('True')\n",
    "plt.ylabel('Predicted')\n",
    "plt.savefig(figures_folder+'all_scatter.pdf') \n",
    "\n",
    "print('Rsquared: '+str(get_all_patients_rsquared.rsquared.mean())+' /rms:'+ str(get_all_patients_rms.rms.mean()))\n",
    "print('inRsquared: '+str(get_all_patients_rsquared.inRsquared.mean())+' /inRms:'+ str(get_all_patients_rms.inRms.mean()))\n",
    "print('outRsquared: '+str(get_all_patients_rsquared.outRsquared.mean())+' /outRms:'+ str(get_all_patients_rms.outRms.mean()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Rsquared: '+str(get_all_patients_rsquared.rsquared.mean())+' /rms:'+ str(get_all_patients_rms.rms.mean()))\n",
    "print('inRsquared: '+str(get_all_patients_rsquared.inRsquared.mean())+' /inRms:'+ str(get_all_patients_rms.inRms.mean()))\n",
    "print('outRsquared: '+str(get_all_patients_rsquared.outRsquared.mean())+' /outRms:'+ str(get_all_patients_rms.outRms.mean()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "get_all_patients_rsquared"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "get_all_patients_rms"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# neuralnet model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "EEm_details = pd.read_csv('EEm_details_per_partic.csv')  \n",
    "\n",
    "# numberForVal = 1 # the number of participants you want per gtoup as validation set (random seed is the number of participants' id)\n",
    "# randomSeed = i\n",
    "w =30\n",
    "# data = features_bindAll\n",
    "dataToUse = 'aw_'+str(w)+'sec'\n",
    "model_name = 'neuralnet_'+dataToUse\n",
    "figures_folder = 'images/' + model_name + '/predPlot/'\n",
    "csv_folder = figures_folder\n",
    "\n",
    "if not os.path.exists(figures_folder):\n",
    "    os.makedirs(figures_folder)\n",
    "\n",
    "get_all_patients_rsquared = pd.DataFrame(columns=['participant', 'rsquared', 'inRsquared', 'outRsquared'])\n",
    "get_all_patients_rms = pd.DataFrame(columns=['participant', 'rms', 'inRms', 'outRms'])\n",
    "stats_EEm_results = pd.DataFrame(columns=['participant', \n",
    "                                          'true_tot_mean', 'true_in_mean', 'true_out_mean', \n",
    "                                          'pred_tot_mean', 'pred_in_mean', 'pred_out_mean'])\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "df_allPred = pd.DataFrame()\n",
    "df_allPredin = pd.DataFrame()\n",
    "df_allPredout = pd.DataFrame()\n",
    "\n",
    "for i in range(5, 37):\n",
    "    randomSeed = i\n",
    "    if i in [2,3,4,19,9,26,27,28,29]: #[12, 19, 23, 16]: \n",
    "        continue\n",
    "#         print('Skip ID')\n",
    "    else: \n",
    "        if len(str(i)) == 1:\n",
    "            name = 'GOTOV0'+str(i)\n",
    "        else: \n",
    "            name = 'GOTOV'+str(i)\n",
    "        #eIF\n",
    "    #eIF\n",
    "    time = EEm_details['timeToSplit'][EEm_details['participant']==name].values[0] #time to split in and outdoors\n",
    "    \n",
    "    # load predictions of R models (they were trained in Shark cluster using another R-Jupyter notebook)\n",
    "    print('read predictions:', name)\n",
    "    \n",
    "    df_test = pd.read_csv('ANNs_predictions_Montoye_repl/neuralnet/'+dataToUse+'/'+name+\"_preds.csv\", header = 0, index_col = None, low_memory=False)\n",
    "    df_test.columns = ['time', 'True', 'Predicted']\n",
    "    df_test.set_index('time', inplace=True)\n",
    "    df_test.index = pd.to_datetime(df_test.index)\n",
    "\n",
    "    df_test['participant'] = name        \n",
    "    \n",
    "    df_allPred = df_allPred.append(df_test)\n",
    "\n",
    "    rsquared = r2_score(df_test['True'], df_test['Predicted'])\n",
    "    print('rsquared...', rsquared)\n",
    "\n",
    "    rms = sqrt(mean_squared_error(df_test['True'], df_test['Predicted']))\n",
    "    print('rms...', rms)\n",
    "    \n",
    "    print('Create figures ....')\n",
    "\n",
    "    true_tot_mean = np.mean(df_test['True'])\n",
    "    pred_tot_mean = np.mean(df_test['Predicted'])        \n",
    "    true_tot_std  = np.std(df_test['True'])\n",
    "    pred_tot_std  = np.std(df_test['Predicted'])        \n",
    "\n",
    "    z = np.polyfit(df_test['True'], df_test['Predicted'], 1)\n",
    "    p = np.poly1d(z)\n",
    "    x_max = int(np.max(df_test['True'])+2)\n",
    "    y_max = int(np.max(df_test['Predicted'])+2)\n",
    "\n",
    "    plt.figure(figsize=(15,8))\n",
    "    plt.plot(df_test['True'], label='True_EE')\n",
    "    plt.plot(df_test['Predicted'], label='Predicted_EE')\n",
    "    plt.legend(loc='upper left')\n",
    "    plt.title(name)\n",
    "    plt.savefig(figures_folder+name+'.pdf')\n",
    "    # plt.show()\n",
    "    plt.close()\n",
    "    \n",
    "    if np.isnan(EEm_details['outEEm'][EEm_details['participant']==name].values) == False:\n",
    "        # indoors\n",
    "        indoors = df_test[df_test.index <= time]#'2016-02-19 11:12:25.004000']\n",
    "        indoors['participant'] = name        \n",
    "        df_allPredin = df_allPredin.append(indoors)\n",
    "\n",
    "        in_rsquared = r2_score(indoors['True'], indoors['Predicted'])\n",
    "        in_rms = sqrt(mean_squared_error(indoors['True'], indoors['Predicted']))\n",
    "        print('in :', in_rsquared)\n",
    "        true_in_mean = np.mean(indoors['True'])\n",
    "        pred_in_mean = np.mean(indoors['Predicted'])\n",
    "        true_in_std  = np.std(indoors['True'])\n",
    "        pred_in_std  = np.std(indoors['Predicted'])        \n",
    "\n",
    "        # plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.plot(indoors['True'], label='True_EE')\n",
    "        plt.plot(indoors['Predicted'], label='Predicted_EE')\n",
    "        plt.title(name)\n",
    "        plt.legend(loc='upper left')\n",
    "        plt.savefig(figures_folder+name+'_in.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        # outdoors\n",
    "        outdoors = df_test[df_test.index > time]#'2016-02-19 11:12:25.004000']\n",
    "        outdoors['participant'] = name        \n",
    "        df_allPredout = df_allPredout.append(outdoors)\n",
    "\n",
    "        out_rsquared = r2_score(outdoors['True'], outdoors['Predicted'])\n",
    "        out_rms = sqrt(mean_squared_error(outdoors['True'], outdoors['Predicted']))            \n",
    "        print('out:', out_rsquared)\n",
    "        true_out_mean = np.mean(outdoors['True'])\n",
    "        pred_out_mean = np.mean(outdoors['Predicted'])            \n",
    "        true_out_std  = np.std(outdoors['True'])\n",
    "        pred_out_std  = np.std(outdoors['Predicted'])        \n",
    "\n",
    "        # plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.title(name)\n",
    "        plt.plot(outdoors['True'], label='True_EE')\n",
    "        plt.plot(outdoors['Predicted'], label='Predicted_EE')\n",
    "        plt.legend(loc='upper left')\n",
    "        plt.savefig(figures_folder+name+'_out.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        # scatter plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.xlim([-1,x_max])\n",
    "        plt.ylim([-1,y_max])            \n",
    "        plt.scatter(indoors['True'], indoors['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "        plt.scatter(outdoors['True'], outdoors['Predicted'], c ='orange', alpha=0.5, label='Outdoors')\n",
    "        plt.plot(df_test['True'],p(df_test['True']),\"b\")\n",
    "        plt.legend(loc='upper right')\n",
    "        plt.title(name)\n",
    "        plt.xlabel('True')\n",
    "        plt.ylabel('Predicted')\n",
    "        plt.savefig(figures_folder+name+'_scatter.pdf')            \n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "    else:\n",
    "        in_rsquared   = None\n",
    "        in_rms        = None\n",
    "        out_rsquared  = None\n",
    "        out_rms       = None\n",
    "        true_in_mean  = None\n",
    "        true_out_mean = None\n",
    "        true_in_std   = None\n",
    "        true_out_std  = None\n",
    "        pred_in_mean  = None\n",
    "        pred_out_mean = None\n",
    "        pred_in_std   = None\n",
    "        pred_out_std  = None\n",
    "\n",
    "        # scatter plot\n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.xlim([-1,x_max])\n",
    "        plt.ylim([-1,y_max])\n",
    "        plt.scatter(df_test['True'], df_test['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "        plt.plot(df_test['True'],p(df_test['True']),\"b\")            \n",
    "        plt.legend(loc='upper right')\n",
    "        plt.title(name)\n",
    "        plt.xlabel('True')\n",
    "        plt.ylabel('Predicted')\n",
    "        plt.savefig(figures_folder+name+'_scatter.pdf')\n",
    "    #             plt.show()\n",
    "        plt.close()\n",
    "\n",
    "        df_allPredin = df_allPredin.append(df_test)\n",
    "    #fi\n",
    "\n",
    "    get_all_patients_rsquared = get_all_patients_rsquared.append({'participant': name,\n",
    "                                                                  'rsquared'   : rsquared,\n",
    "                                                                  'inRsquared' : in_rsquared,\n",
    "                                                                  'outRsquared': out_rsquared},\n",
    "                                                         ignore_index=True) \n",
    "\n",
    "    get_all_patients_rms = get_all_patients_rms.append({'participant': name,\n",
    "                                                        'rms'        : rms,\n",
    "                                                        'inRms'      : in_rms,\n",
    "                                                        'outRms'     : out_rms},\n",
    "                                                         ignore_index=True) \n",
    "\n",
    "    stats_EEm_results = stats_EEm_results.append({'participant'  : name,\n",
    "                                                  'true_tot_mean': true_tot_mean,\n",
    "                                                  'true_in_mean' : true_in_mean,\n",
    "                                                  'true_out_mean': true_out_mean,\n",
    "                                                  'true_tot_std' : true_tot_std,\n",
    "                                                  'true_in_std'  : true_in_std,\n",
    "                                                  'true_out_std' : true_out_std,\n",
    "                                                  'pred_tot_mean': pred_tot_mean,\n",
    "                                                  'pred_in_mean' : pred_in_mean,\n",
    "\n",
    "                                                  'pred_out_mean': pred_out_mean,\n",
    "                                                  'pred_tot_std' : pred_tot_std,\n",
    "                                                  'pred_in_std'  : pred_in_std,\n",
    "                                                  'pred_out_std' : pred_out_std},\n",
    "                                                 ignore_index=True)\n",
    "#efor\n",
    "\n",
    "get_all_patients_rsquared.to_csv(csv_folder+'rsquared.csv',index=False)\n",
    "get_all_patients_rms.to_csv(csv_folder+'rms.csv',index=False)\n",
    "stats_EEm_results.to_csv(csv_folder+'predStats.csv',index=False)\n",
    "\n",
    "# scatter\n",
    "true = stats_EEm_results['true_tot_mean']\n",
    "pred = stats_EEm_results['pred_tot_mean']\n",
    "# std_true  = stats_EEm_results['true_tot_std']\n",
    "# std_pred  = stats_EEm_results['pred_tot_std']\n",
    "# scatter\n",
    "z = np.polyfit(true, pred, 1)\n",
    "p = np.poly1d(z)\n",
    "x_max = int(np.max(true)+1)\n",
    "y_max = int(np.max(pred)+1)\n",
    "plt.figure(figsize=(15,8))\n",
    "plt.xlim([1,x_max])\n",
    "plt.ylim([1,y_max])\n",
    "plt.scatter(true, pred, c ='blue', alpha=0.5)#, label='Indoors')\n",
    "plt.plot(true,p(true),\"b\")     \n",
    "# plt.legend(loc='upper right')\n",
    "plt.title('True Vs Predicted Average EEm per participant')\n",
    "plt.xlabel('True')\n",
    "plt.ylabel('Predicted')\n",
    "plt.savefig(figures_folder+'all_scatter.pdf') \n",
    "\n",
    "print('Rsquared: '+str(get_all_patients_rsquared.rsquared.mean())+' /rms:'+ str(get_all_patients_rms.rms.mean()))\n",
    "print('inRsquared: '+str(get_all_patients_rsquared.inRsquared.mean())+' /inRms:'+ str(get_all_patients_rms.inRms.mean()))\n",
    "print('outRsquared: '+str(get_all_patients_rsquared.outRsquared.mean())+' /outRms:'+ str(get_all_patients_rms.outRms.mean()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Rsquared: '+str(get_all_patients_rsquared.rsquared.mean())+' /rms:'+ str(get_all_patients_rms.rms.mean()))\n",
    "print('inRsquared: '+str(get_all_patients_rsquared.inRsquared.mean())+' /inRms:'+ str(get_all_patients_rms.inRms.mean()))\n",
    "print('outRsquared: '+str(get_all_patients_rsquared.outRsquared.mean())+' /outRms:'+ str(get_all_patients_rms.outRms.mean()))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# CNN model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "os.environ[\"CUDA_DEVICE_ORDER\"] = \"PCI_BUS_ID\"\n",
    "os.environ[\"CUDA_VISIBLE_DEVICES\"] = \"6\" # set cuda to only use this/these GPU (check if free with nvidia-smi first)`\n",
    "os.getcwd()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import datetime\n",
    "import pickle\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import tensorflow as tf\n",
    "from tensorflow.keras.layers import Dense, Dropout, Input\n",
    "from tensorflow.keras.layers import Conv1D, AveragePooling1D, Flatten\n",
    "from tensorflow.keras.models import Model\n",
    "from tensorflow.keras.layers import Input, concatenate\n",
    "from tensorflow.keras.callbacks import EarlyStopping, ModelCheckpoint, TensorBoard\n",
    "from tensorflow.keras.utils import plot_model\n",
    "from livelossplot.keras import PlotLossesCallback\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "tf.test.gpu_device_name()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def CNN_modeling(name, X_train, y_train, bmi_train, X_val, y_val, bmi_val, n_batch, n_epochs, moldel_dir, image_dir, BMI=True):\n",
    "# the function is used to train the different models per participant based on the decided architecture\n",
    "# INPUTS:\n",
    "#     - name: The id-name of test participant (e.g. GOTOV00)\n",
    "#     - X_train/val: The training/validation predictors (accel)\n",
    "#     - y_train/val: The training/validation targer (EE)\n",
    "#     - bmi_train/val: The training/validation participant level data\n",
    "#     - n_batch: batch size\n",
    "#     - n_epochs: number of epochs\n",
    "#     - moldel_dir: the directory that the models are saved\n",
    "#     - image_dir: the directory that the images are saved\n",
    "# OUTPUTS: The training history (model, images)\n",
    "\n",
    "    print('Creating functional api model')\n",
    "# Use functional API to build cnn layers\n",
    "#     - 1st conv1D layer with filters = 8 and kernel_size = 5,\n",
    "#     - then average pooling factor 2D\n",
    "#     - 2nd conv1D layer with filters = 4 and kernel_size = 5,\n",
    "#     - then average pooling factor 2D\n",
    "#     - Flatten \n",
    "#     - MLP: Dense layer of size = 400\n",
    "\n",
    "#     input1 accelerometer sensor measurements \n",
    "    input_1 = Input(shape=(X_train.shape[1], X_train.shape[2]))\n",
    "#     1st conv layer: filters = 8, size = 5\n",
    "    conv_1 = Conv1D(filters=8, kernel_size=5, activation='tanh')(input_1)\n",
    "    pool_1 = AveragePooling1D(pool_size=2)(conv_1)\n",
    "#     2nd conv layer: filters = 4, size = 5\n",
    "    conv_2 = Conv1D(filters=4, kernel_size=5, activation='tanh')(pool_1)\n",
    "    pool_2 = AveragePooling1D(pool_size=2)(conv_2)    \n",
    "#     image to vector before connecting to dense layer\n",
    "    flat = Flatten()(pool_2)\n",
    "\n",
    "#     input2 BMI of participants \n",
    "    input_2 = Input(shape=(bmi_train.shape[1],))\n",
    "    \n",
    "    if BMI == True:        \n",
    "    #     merge\n",
    "        con = concatenate([input_2, flat])\n",
    "    #         MLP: the dense layer of size 400\n",
    "        mlp_1 = Dense(400)(con)\n",
    "    else:\n",
    "        mlp_1 = Dense(400)(flat)\n",
    "       \n",
    "    output = Dense(1, activation='linear')(mlp_1)\n",
    "\n",
    "    model = Model(inputs=[input_1, input_2], outputs=output)\n",
    "    # plot_model(model, to_file=image_dir)\n",
    "\n",
    "    model.compile(loss='mean_squared_error', optimizer='sgd', metrics=['mse'])\n",
    "    model.summary()\n",
    "    tbCallBack = TensorBoard(log_dir= moldel_dir + '/' + name, histogram_freq=0, write_graph=True, write_images=True)\n",
    "    # monitor = EarlyStopping(monitor='val_loss', min_delta=1e-3, patience=50, mode='auto', verbose=1)\n",
    "    \n",
    "    checkpointer = ModelCheckpoint(filepath= moldel_dir + '/' + name + '.hdf5', \n",
    "                                  save_best_only=True, verbose=1, mode='auto', monitor='val_loss')\n",
    "    \n",
    "    history = model.fit(x=[X_train, bmi_train], y=y_train, \n",
    "                      batch_size=n_batch, epochs=n_epochs, verbose=1, \n",
    "                      validation_data=([X_val, bmi_val], y_val), shuffle=True, \n",
    "                      callbacks=[checkpointer, tbCallBack])\n",
    "    return history\n",
    "#eDEF"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Train"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # # # # # # # # # # # FIXED VARs # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
    "n_batch = 512\n",
    "epochs = 50\n",
    "EEm_rate = 10\n",
    "window = '5.12sec'\n",
    "\n",
    "BMI = True\n",
    "\n",
    "device = 'wrist'\n",
    "#2min_seqs_1a1val_EEm_60ds_std_480\n",
    "\n",
    "seq_path = 'sequences/5.12secs_seqs_1a1val_EEm_10ds_mean_256/'#'sequences/'+str(minutes)+'min_seqs_1val_EEm_'+str(EEm_rate)+'ds_sd/'\n",
    "orig_seq_path = 'sequences/5.12secs_seqs_1a1val_EEm_0ds_mean_256/'#'sequences/'+str(minutes)+'min_seqs_1val_EEm_'+str(EEm_rate)+'ds_sd/'\n",
    "model_name = 'CNNrepl_EEm_'+str(EEm_rate)+'ds_1Seg_'+str(window)+'_'+str(epochs)+'epochs_'+str(n_batch)+'batch_1val_' + device\n",
    "if BMI==False:\n",
    "    model_name = model_name+'_noBMI'\n",
    "else:\n",
    "    model_name = model_name+'_withBMI'\n",
    "    \n",
    "folder = 'models/' + model_name\n",
    "image = 'images/' + model_name + '/Accur_Loss_figures/'\n",
    "predPlot = 'images/' + model_name + '/predPlot/'\n",
    "\n",
    "print('Seqs path :', seq_path)\n",
    "print('Model name:', model_name)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "# run the modeling functions\n",
    "\n",
    "if not os.path.exists(folder):\n",
    "    os.makedirs(folder) \n",
    "if not os.path.exists(image):\n",
    "    os.makedirs(image)\n",
    "    os.makedirs(predPlot)\n",
    "#eIF\n",
    "  \n",
    "# names = [5]\n",
    "# for i in names:\n",
    "for i in range(5, 37):\n",
    "    if i in [2,3,4,19]: \n",
    "    # if i in [12, 16,19,23]:\n",
    "        continue\n",
    "    else:\n",
    "        if len(str(i)) ==1:\n",
    "            name = 'GOTOV0'+str(i)\n",
    "        else:\n",
    "            name = 'GOTOV'+str(i)\n",
    "    if os.path.exists(seq_path+name+'.pkl') == True:\n",
    "        \n",
    "        part =  model_name+'_'+ name\n",
    "        print('Load sequences of participant', name)\n",
    "        with open(seq_path+name+'.pkl','rb') as f:\n",
    "            X_train, y_train, ytrain_time, bmi_train, X_val, y_val, yval_time, bmi_val, X_test, y_test, ytest_time, bmi_test, scaler = pickle.load(f)\n",
    "        # load the orignal val and test\n",
    "        with open(orig_seq_path+name+'.pkl','rb') as f:\n",
    "            X_train1, y_train1, ytrain_time1, bmi_train1, X_val, y_val, yval_time, bmi_val, X_test, y_test, ytest_time, bmi_test, scaler = pickle.load(f)\n",
    "        \n",
    "        X_train1, y_train1, ytrain_time1 = [],[],[]\n",
    "        \n",
    "        # keep only one device data\n",
    "        if device == 'ankle':\n",
    "            X_train = X_train[:,:,0:3]\n",
    "            X_val = X_val[:,:,0:3]\n",
    "            X_test = X_test[:,:,0:3]\n",
    "        elif device == 'wrist':\n",
    "            X_train = X_train[:,:,3:6]\n",
    "            X_val = X_val[:,:,3:6]\n",
    "            X_test = X_test[:,:,3:6]\n",
    "            \n",
    "        print('Training model... ')\n",
    "        model = CNN_modeling(name, X_train, y_train, bmi_train, \n",
    "                             X_val, y_val, bmi_val, n_batch, epochs, folder, image,\n",
    "                             BMI=BMI)\n",
    "        \n",
    "        plt.plot(model.history['mse'])\n",
    "        plt.plot(model.history['val_mse'])\n",
    "        plt.title('model Mean Square Error')\n",
    "        plt.ylabel('Mean Square Error')\n",
    "        plt.xlabel('epoch')\n",
    "        plt.legend(['train', 'val'], loc='upper left')\n",
    "        plt.savefig(image + name + '_mse.pdf')\n",
    "        plt.show()\n",
    "        # summarize history for loss\n",
    "        plt.plot(model.history['loss'])\n",
    "        plt.plot(model.history['val_loss'])\n",
    "        plt.title('model loss')\n",
    "        plt.ylabel('loss')\n",
    "        plt.xlabel('epoch')\n",
    "        plt.legend(['train', 'val'], loc='upper left')\n",
    "        plt.savefig(image + name + '_loss.pdf')\n",
    "        plt.show()\n",
    "    else:\n",
    "        print('filepath doesnt exist')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Predict"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from tensorflow.keras.models import load_model\n",
    "from sklearn.metrics import r2_score\n",
    "import numpy as np\n",
    "from sklearn.metrics import mean_squared_error\n",
    "from math import sqrt\n",
    "import matplotlib.pyplot as plt\n",
    "import csv\n",
    "import glob\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "# Display figures inline in Jupyter notebook\n",
    "import seaborn as sns\n",
    "# Use seaborn style defaults and set the default figure size\n",
    "sns.set(rc={'figure.figsize':(11, 4)})\n",
    "# confs for plotly\n",
    "def configure_plotly_browser_state():\n",
    "  import IPython\n",
    "  display(IPython.core.display.HTML('''\n",
    "        <script src=\"/static/components/requirejs/require.js\"></script>\n",
    "        <script>\n",
    "          requirejs.config({\n",
    "            paths: {\n",
    "              base: '/static/base',\n",
    "              plotly: 'https://cdn.plot.ly/plotly-1.5.1.min.js?noext',\n",
    "            },\n",
    "          });\n",
    "        </script>\n",
    "        '''))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# n_batch = 512\n",
    "# epochs = 50\n",
    "# EEm_rate = 10\n",
    "# windowInMinutes = 4     \n",
    "\n",
    "EEm_details = pd.read_csv('EEm_details_per_partic.csv')  \n",
    "seq_path = orig_seq_path\n",
    "\n",
    "trained_models_location = folder+'/'\n",
    "patient_folder = orig_seq_path # test with orignal target data\n",
    "figures_folder = 'images/' + model_name + '/predPlot/'\n",
    "csv_folder = figures_folder\n",
    "\n",
    "# get_all_patients_rsquared, get_all_patients_rms = {}, {} \n",
    "get_all_patients_rsquared = pd.DataFrame(columns=['participant', 'rsquared', 'inRsquared', 'outRsquared'])\n",
    "get_all_patients_rms = pd.DataFrame(columns=['participant', 'rms', 'inRms', 'outRms'])\n",
    "stats_EEm_results = pd.DataFrame(columns=['participant', \n",
    "                                          'true_tot_mean', 'true_in_mean', 'true_out_mean', \n",
    "                                          'pred_tot_mean', 'pred_in_mean', 'pred_out_mean'])\n",
    "\n",
    "print('Model:', trained_models_location)\n",
    "print('Seqs:', patient_folder)\n",
    "print('Results:',figures_folder)\n",
    "\n",
    "import glob\n",
    "dire= list(glob.glob(trained_models_location +'*.hdf5'))\n",
    "dire.sort()\n",
    "print('N =',len(dire))\n",
    "# dire"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "df_allPred = pd.DataFrame()\n",
    "df_allPredin = pd.DataFrame()\n",
    "df_allPredout = pd.DataFrame()\n",
    "\n",
    "# names = [36]\n",
    "# for i in names:\n",
    "for i in range(5, 37):\n",
    "    if i in [2,3,4,19]: \n",
    "    # if i in [12, 16,19,23]:\n",
    "        continue\n",
    "    else:\n",
    "        if len(str(i)) ==1: \n",
    "            name = 'GOTOV0'+str(i)\n",
    "        else:\n",
    "            name = 'GOTOV'+str(i)\n",
    "        #fi\n",
    "    #fi    \n",
    "    time = EEm_details['timeToSplit'][EEm_details['participant']==name].values[0]    \n",
    "    if os.path.exists(trained_models_location + name +'.hdf5') == True: \n",
    "        print('Predicting test for patient', name)\n",
    "        with open(patient_folder+name+'.pkl','rb') as f:\n",
    "            X_train, y_train, ytrain_time, bmi_train,X_val, y_val, yval_time, bmi_val, X_test, y_test, ytest_time, bmi_test, scaler = pickle.load(f)\n",
    "        \n",
    "        # keep only one device data\n",
    "        if device == 'ankle':\n",
    "            X_train = X_train[:,:,0:3]\n",
    "            X_val = X_val[:,:,0:3]\n",
    "            X_test1 = X_test[:,:,0:3]\n",
    "        elif device == 'wrist':\n",
    "            X_train = X_train[:,:,3:6]\n",
    "            X_val = X_val[:,:,3:6]\n",
    "            X_test1 = X_test[:,:,3:6]\n",
    "        else:\n",
    "            X_test1 = X_test\n",
    "            \n",
    "#         X_train, y_train, ytrain_time, bmi_train, X_val, y_val, yval_time, bmi_val, X_test, y_test, ytest_time, bmi_test, scaler = dataConcatenate(name, mn_seq_path, std_seq_path)\n",
    "        # part = model_name+'_'+name\n",
    "        print('loading model....')\n",
    "        model = load_model(trained_models_location+name+ '.hdf5')\n",
    "        yhat = model.predict([X_test1, bmi_test])\n",
    "        print('predicting....')\n",
    "        # reverse data from standardization to compute the correct MSE\n",
    "#         print('Inversing data.....') \n",
    "        X_test = X_test.reshape(-1, X_test.shape[2])\n",
    "#         d_test = pd.DataFrame(data=X_test)\n",
    "#         labels = d_test[6]\n",
    "#         d_test = d_test.drop([6], axis=1)\n",
    "#         X_test = d_test.to_numpy()\n",
    "\n",
    "        inv_yhat = np.empty((X_test.shape[0], 1))\n",
    "        inv_yhat.fill(np.nan)\n",
    "        inv_yhat[:yhat.shape[0]] = yhat\n",
    "\n",
    "        inv_yhat = np.concatenate((X_test, inv_yhat), axis=1)\n",
    "        inv_yhat = np.ma.array(inv_yhat, mask=np.isnan(inv_yhat))\n",
    "        inv_yhat = scaler.inverse_transform(inv_yhat)\n",
    "        inv_yhat = inv_yhat[:, 6]\n",
    "\n",
    "        y_test = y_test.reshape(len(y_test), 1)\n",
    "        inv_ytest = np.empty((X_test.shape[0], 1))\n",
    "        inv_ytest.fill(np.nan)\n",
    "        inv_ytest[:y_test.shape[0]] = y_test\n",
    "\n",
    "        inv_ytest = np.concatenate((X_test, inv_ytest), axis=1)\n",
    "        inv_ytest = np.ma.array(inv_ytest, mask=np.isnan(inv_ytest))\n",
    "        inv_ytest = scaler.inverse_transform(inv_ytest)\n",
    "        inv_ytest = inv_ytest[:, 6]\n",
    "\n",
    "        inv_yhat = inv_yhat[~np.isnan(inv_yhat)]\n",
    "        inv_ytest = inv_ytest[~np.isnan(inv_ytest)]\n",
    "\n",
    "        df_test = pd.DataFrame(data=inv_ytest, index=ytest_time, columns=['True'])\n",
    "        df_test['Predicted'] = inv_yhat\n",
    "\n",
    "        df_test['participant'] = name        \n",
    "        df_allPred = df_allPred.append(df_test)\n",
    "\n",
    "        # downsample\n",
    "        \n",
    "        rsquared = r2_score(df_test['True'], df_test['Predicted'])\n",
    "        print('rsquared...', rsquared)\n",
    "\n",
    "        rms = sqrt(mean_squared_error(df_test['True'], df_test['Predicted']))\n",
    "        print('rms...', rms)\n",
    "\n",
    "        print('Create figures ....')\n",
    "\n",
    "        true_tot_mean = np.mean(df_test['True'])\n",
    "        pred_tot_mean = np.mean(df_test['Predicted'])        \n",
    "        true_tot_std  = np.std(df_test['True'])\n",
    "        pred_tot_std  = np.std(df_test['Predicted'])        \n",
    "                        \n",
    "        z = np.polyfit(df_test['True'], df_test['Predicted'], 1)\n",
    "        p = np.poly1d(z)\n",
    "        x_max = int(np.max(df_test['True'])+2)\n",
    "        y_max = int(np.max(df_test['Predicted'])+2)\n",
    "        \n",
    "        plt.figure(figsize=(15,8))\n",
    "        plt.plot(df_test['True'], label='True_EE')\n",
    "        plt.plot(df_test['Predicted'], label='Predicted_EE')\n",
    "        plt.legend(loc='upper left')\n",
    "        plt.title(name)\n",
    "        plt.savefig(figures_folder+name+'.pdf')\n",
    "        # plt.show()\n",
    "        plt.close()\n",
    "        \n",
    "        if np.isnan(EEm_details['outEEm'][EEm_details['participant']==name].values) == False:\n",
    "            # indoors\n",
    "            indoors = df_test[df_test.index <= time]#'2016-02-19 11:12:25.004000']\n",
    "            indoors['participant'] = name        \n",
    "            df_allPredin = df_allPredin.append(indoors)\n",
    "            \n",
    "            in_rsquared = r2_score(indoors['True'], indoors['Predicted'])\n",
    "            in_rms = sqrt(mean_squared_error(indoors['True'], indoors['Predicted']))\n",
    "            print('in :', in_rsquared)\n",
    "            true_in_mean = np.mean(indoors['True'])\n",
    "            pred_in_mean = np.mean(indoors['Predicted'])\n",
    "            true_in_std  = np.std(indoors['True'])\n",
    "            pred_in_std  = np.std(indoors['Predicted'])        \n",
    "    \n",
    "            # plot\n",
    "            plt.figure(figsize=(15,8))\n",
    "            plt.plot(indoors['True'], label='True_EE')\n",
    "            plt.plot(indoors['Predicted'], label='Predicted_EE')\n",
    "            plt.title(name)\n",
    "            plt.legend(loc='upper left')\n",
    "            plt.savefig(figures_folder+name+'_in.pdf')\n",
    "#             plt.show()\n",
    "            plt.close()\n",
    "            \n",
    "            # outdoors\n",
    "            outdoors = df_test[df_test.index > time]#'2016-02-19 11:12:25.004000']\n",
    "            outdoors['participant'] = name        \n",
    "            df_allPredout = df_allPredout.append(outdoors)\n",
    "            \n",
    "            out_rsquared = r2_score(outdoors['True'], outdoors['Predicted'])\n",
    "            out_rms = sqrt(mean_squared_error(outdoors['True'], outdoors['Predicted']))            \n",
    "            print('out:', out_rsquared)\n",
    "            true_out_mean = np.mean(outdoors['True'])\n",
    "            pred_out_mean = np.mean(outdoors['Predicted'])            \n",
    "            true_out_std  = np.std(outdoors['True'])\n",
    "            pred_out_std  = np.std(outdoors['Predicted'])        \n",
    "\n",
    "            # plot\n",
    "            plt.figure(figsize=(15,8))\n",
    "            plt.title(name)\n",
    "            plt.plot(outdoors['True'], label='True_EE')\n",
    "            plt.plot(outdoors['Predicted'], label='Predicted_EE')\n",
    "            plt.legend(loc='upper left')\n",
    "            plt.savefig(figures_folder+name+'_out.pdf')\n",
    "#             plt.show()\n",
    "            plt.close()\n",
    "        \n",
    "            # scatter plot\n",
    "            plt.figure(figsize=(15,8))\n",
    "            plt.xlim([-1,x_max])\n",
    "            plt.ylim([-1,y_max])            \n",
    "            plt.scatter(indoors['True'], indoors['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "            plt.scatter(outdoors['True'], outdoors['Predicted'], c ='orange', alpha=0.5, label='Outdoors')\n",
    "            plt.plot(df_test['True'],p(df_test['True']),\"b\")\n",
    "            plt.legend(loc='upper right')\n",
    "            plt.title(name)\n",
    "            plt.xlabel('True')\n",
    "            plt.ylabel('Predicted')\n",
    "            plt.savefig(figures_folder+name+'_scatter.pdf')            \n",
    "#             plt.show()\n",
    "            plt.close()\n",
    "    \n",
    "        else:\n",
    "            in_rsquared   = None\n",
    "            in_rms        = None\n",
    "            out_rsquared  = None\n",
    "            out_rms       = None\n",
    "            true_in_mean  = None\n",
    "            true_out_mean = None\n",
    "            true_in_std   = None\n",
    "            true_out_std  = None\n",
    "            pred_in_mean  = None\n",
    "            pred_out_mean = None\n",
    "            pred_in_std   = None\n",
    "            pred_out_std  = None\n",
    "            \n",
    "            # scatter plot\n",
    "            plt.figure(figsize=(15,8))\n",
    "            plt.xlim([-1,x_max])\n",
    "            plt.ylim([-1,y_max])\n",
    "            plt.scatter(df_test['True'], df_test['Predicted'], c ='green', alpha=0.5, label='Indoors')\n",
    "            plt.plot(df_test['True'],p(df_test['True']),\"b\")            \n",
    "            plt.legend(loc='upper right')\n",
    "            plt.title(name)\n",
    "            plt.xlabel('True')\n",
    "            plt.ylabel('Predicted')\n",
    "            plt.savefig(figures_folder+name+'_scatter.pdf')\n",
    "#             plt.show()\n",
    "            plt.close()\n",
    "    \n",
    "            df_allPredin = df_allPredin.append(df_test)\n",
    "        #fi\n",
    "    \n",
    "        get_all_patients_rsquared = get_all_patients_rsquared.append({'participant': name,\n",
    "                                                                      'rsquared'   : rsquared,\n",
    "                                                                      'inRsquared' : in_rsquared,\n",
    "                                                                      'outRsquared': out_rsquared},\n",
    "                                                             ignore_index=True) \n",
    "        \n",
    "        get_all_patients_rms = get_all_patients_rms.append({'participant': name,\n",
    "                                                            'rms'        : rms,\n",
    "                                                            'inRms'      : in_rms,\n",
    "                                                            'outRms'     : out_rms},\n",
    "                                                             ignore_index=True) \n",
    "        \n",
    "        stats_EEm_results = stats_EEm_results.append({'participant'  : name,\n",
    "                                                      'true_tot_mean': true_tot_mean,\n",
    "                                                      'true_in_mean' : true_in_mean,\n",
    "                                                      'true_out_mean': true_out_mean,\n",
    "                                                      'true_tot_std' : true_tot_std,\n",
    "                                                      'true_in_std'  : true_in_std,\n",
    "                                                      'true_out_std' : true_out_std,\n",
    "                                                      'pred_tot_mean': pred_tot_mean,\n",
    "                                                      'pred_in_mean' : pred_in_mean,\n",
    "     \n",
    "                                                      'pred_out_mean': pred_out_mean,\n",
    "                                                      'pred_tot_std' : pred_tot_std,\n",
    "                                                      'pred_in_std'  : pred_in_std,\n",
    "                                                      'pred_out_std' : pred_out_std},\n",
    "                                                     ignore_index=True)\n",
    "\n",
    "\n",
    "    else: \n",
    "        print(\"Model not found\")\n",
    "    #fi\n",
    "#efor\n",
    "\n",
    "get_all_patients_rsquared.to_csv(csv_folder+'rsquared.csv',index=False)\n",
    "get_all_patients_rms.to_csv(csv_folder+'rms.csv',index=False)\n",
    "stats_EEm_results.to_csv(csv_folder+'predStats.csv',index=False)\n",
    "\n",
    "true = stats_EEm_results['true_tot_mean']\n",
    "pred = stats_EEm_results['pred_tot_mean']\n",
    "# std_true  = stats_EEm_results['true_tot_std']\n",
    "# std_pred  = stats_EEm_results['pred_tot_std']\n",
    "# scatter\n",
    "z = np.polyfit(true, pred, 1)\n",
    "p = np.poly1d(z)\n",
    "x_max = int(np.max(true)+1)\n",
    "y_max = int(np.max(pred)+1)\n",
    "plt.figure(figsize=(15,8))\n",
    "plt.xlim([1,x_max])\n",
    "plt.ylim([1,y_max])\n",
    "plt.scatter(true, pred, c ='blue', alpha=0.5)#, label='Indoors')\n",
    "plt.plot(true,p(true),\"b\")     \n",
    "# plt.legend(loc='upper right')\n",
    "plt.title('True Vs Predicted Average EEm per participant')\n",
    "plt.xlabel('True')\n",
    "plt.ylabel('Predicted')\n",
    "plt.savefig(figures_folder+'all_scatter.pdf')            \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Rsquared MD: '+str(get_all_patients_rsquared.rsquared.median())+'/Rsquared SD: '+str(get_all_patients_rsquared.rsquared.std()))\n",
    "print('Rsquared MN: '+str(get_all_patients_rsquared.rsquared.mean())+' /rms MN:'+ str(get_all_patients_rms.rms.mean()))\n",
    "print('inRsquared: '+str(get_all_patients_rsquared.inRsquared.mean())+' /inRms:'+ str(get_all_patients_rms.inRms.mean()))\n",
    "print('outRsquared: '+str(get_all_patients_rsquared.outRsquared.mean())+' /outRms:'+ str(get_all_patients_rms.outRms.mean()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}