title: Supervised Time Series efficiency improvement
body: <p>The data that I have is hourly recorded over the past 4 months. I am building a time series model and I've tried several methods so far: Arima, LSTMs, Prophet but they can be quite slow for my task since I have to run the model on thousands of time series in different locations. So then I thought it might be interesting to transform it into a supervised problem and use regression.</p>
  
  <p>I extracted 4 features from the univariate time series and its time index, namely: dayofweek, hour, daily average and hourly average. So at the moment I am using these 4 predictors but could possibly extract more(like beginning of the day, noon, etc-also if you have any other suggestions here they are very welcomed :) ) </p>
  
  <p>I've used XGBoost for the regression and here are parts of the code: </p>
  
  <pre><code># XGB
  import xgboost as xgb 
  from sklearn.model_selection import train_test_split
  from sklearn.metrics import r2_score
  
  # Functions needed 
  def convert_dates(x):
          x['date'] = pd.to_datetime(x['date'])
          #x['month'] = x['date'].dt.month
          #x['year'] = x['date'].dt.year
          x['dayofweek'] = x['date'].dt.dayofweek
          x['hour'] = x['date'].dt.hour
          #x['week_no'] = pd.to_numeric(x['date'].index.strftime("%V"))
          x.pop('date')
          return(x)
  
  def add_avg(x):
          x['daily_avg']=x.groupby(['dayofweek'])['y'].transform('mean')
          x['hourly_avg'] = x.groupby(['dayofweek','hour'])['y'].transform('mean')
          #x['monthly_avg']=x.groupby(['month'])['y'].transform('mean')
          #x['weekly_avg']=x.groupby(['week_no'])['y'].transform('mean')
          return x
  
  xgb_mape_r2_dict = {}
  </code></pre>
  
  <p>I then run a for loop in which I select a location and build the model for it. Here I split the data into a train and test part. I knew there might be problems due to the Easter holidays in my country last week because those are rare events so that is why I split the training and test data in that manner. So I actually consider the data from the beginning of the year up to two weeks ago as training data and the very next week after that as test data. </p>
  
  <pre><code>for j in range(10,20):
      data = df_all.loc[df_all['Cell_Id']==top_cells[j]]
      data.drop(['Cell_Id', 'WDay'], axis = 1, inplace = True)
      data['date'] = data.index
      period = 168
      data_train = data.iloc[:-2*period,:]
      data_test = data.iloc[-2*period:-period,:]
  
  
      data_train = convert_dates(data_train)
      data_test = convert_dates(data_test)
      data_train.columns = ['y', 'dayofweek', 'hour']
      data_test.columns = ['y', 'dayofweek', 'hour']
  
      data_train = add_avg(data_train)
      daily_avg = data_train.groupby(['dayofweek'])['y'].mean().reset_index()
      hourly_avg = data_train.groupby(['dayofweek', 'hour'])['y'].mean().reset_index()
  </code></pre>
  
  <p>Now, for the test data I add the past averages, namely the 7 daily averages from the past and the 168 hourly averages from the past as well. This is actually the part that takes the longest amount of time to run and I would like to improve its efficiency. </p>
  
  <pre><code>    value_dict ={}
      for k in range(168):
          value_dict[tuple(hourly_avg.iloc[k])[:2]] = tuple(hourly_avg.iloc[k])[2]
  
  
      data_test['daily_avg'] = 0
      data_test['hourly_avg'] = 0
      for i in range(len(data_test)):
          data_test['daily_avg'][i] = daily_avg['y'][data_test['dayofweek'][i]]
          data_test['hourly_avg'][i] = value_dict[(data_test['dayofweek'][i], data_test['hour'][i])]
  </code></pre>
  
  <p>My current run time is of 30 seconds for every iteration in the for loop which is way too slow because of the poor way that I use to add the averages in the test data. I would really appreciate if anyone could point out how could I implement this bit faster.</p>
  
  <p>I will also add the rest of my code and make some other observations as well: </p>
  
  <pre><code>    x_train = data_train.drop('y',axis=1)
      x_test = data_test.drop('y',axis=1)
      y_train = data_train['y']
      y_test = data_test['y']
  
      def XGBmodel(x_train,x_test,y_train,y_test):
          matrix_train = xgb.DMatrix(x_train,label=y_train)
          matrix_test = xgb.DMatrix(x_test,label=y_test)
          model=xgb.train(params={'objective':'reg:linear','eval_metric':'mae'}
                          ,dtrain=matrix_train,num_boost_round=500, 
                          early_stopping_rounds=20,evals=[(matrix_test,'test')],)
          return model
  
      model=XGBmodel(x_train,x_test,y_train,y_test)
  
  
      #submission = pd.DataFrame(x_pred.pop('id'))
      y_pred = model.predict(xgb.DMatrix(x_test), ntree_limit = model.best_ntree_limit)
  
      #submission['sales']= y_pred
  
      y_pred = pd.DataFrame(y_pred)
      y_test = pd.DataFrame(y_test)
      y_test.reset_index(inplace = True, drop = True)
      compare_df = pd.concat([y_test, y_pred], axis = 1)
      compare_df.columns = ['Real', 'Predicted']
      compare_df.plot()
  
      mape = (np.abs((y_test['y'] - y_pred[0])/y_test['y']).mean())*100
      r2 = r2_score(y_test['y'], y_pred[0])
      xgb_mape_r2_dict[top_cells[j]] = [mape,r2]
  </code></pre>
  
  <p>I've used both R-squared and MAPE as accuracy measures although I don't think MAPE is indicated anymore since I've transformed the time series problem into a regression problem. Any thoughts on your part on this subject?</p>
  
  <p>Thank you very much for your time and consideration. Any help is very much appreciated.</p>
  
  <p>Update: I have managed to fix the issue using pandas' merge. I've first created two dataframes containing the daily averges and hourly averages from the training data and then merged these ataframes with the test data:</p>
  
  <pre><code>data_test = merge(data_test, daily_avg,['dayofweek'],'daily_avg')
  data_test = merge(data_test, hourly_av['dayofweek','hour'],'hourly_avg')
  data_test.columns = ['y', 'dayofweek', 'hour', 'daily_avg', 'hourly_avg']
  </code></pre>
  
  <p>where we used the merge function defined as:</p>
  
  <pre><code>def merge(x,y,col,col_name):
      x =pd.merge(x, y, how='left', on=None, left_on=col, right_on=col,
              left_index=False, right_index=False, sort=True,
               copy=True, indicator=False,validate=None)
      x=x.rename(columns={'sales':col_name})
      return x
  </code></pre>
  
  <p>I can now run the model for 2000 locations per hour on a laptop with decent results but I will try to improve it while keeping it fast. Thank you very much once again.</p>