3-features_engineering.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "# Features Engineering"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "In this notebook we experiment with some features engineering:\n",
    "\n",
    "* using **TSFresh**\n",
    "\n",
    "* and using **ROCKET** (Random Convolutional Kernel Transform)\n",
    "\n",
    "\n",
    "**NOTE:** Before starting exploring this notebook, I recommend checking `1-EDA.ipynb` notebook first - it contains Exploratory Data Analysis and will help you get some understanding of the datasets. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T14:48:55.626975Z",
     "start_time": "2021-11-07T14:48:54.480037Z"
    }
   },
   "outputs": [],
   "source": [
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "\n",
    "from datetime import datetime\n",
    "\n",
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "First, let's read aircraft engines datasets, we start with \"FD001\" dataset, which has:\n",
    "\n",
    "* 100 engines time series in TRAIN set\n",
    "\n",
    "* 100 engines time series in TEST set\n",
    "\n",
    "* 1 Fault condition\n",
    "\n",
    "* 1 Operating condition\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T14:48:55.785509Z",
     "start_time": "2021-11-07T14:48:55.630049Z"
    },
    "tags": [
     "remove"
    ]
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "train.shape = (20631, 27)\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>unit</th>\n",
       "      <th>time_cycles</th>\n",
       "      <th>op_setting_1</th>\n",
       "      <th>op_setting_2</th>\n",
       "      <th>op_setting_3</th>\n",
       "      <th>sensor_1</th>\n",
       "      <th>sensor_2</th>\n",
       "      <th>sensor_3</th>\n",
       "      <th>sensor_4</th>\n",
       "      <th>sensor_5</th>\n",
       "      <th>...</th>\n",
       "      <th>sensor_13</th>\n",
       "      <th>sensor_14</th>\n",
       "      <th>sensor_15</th>\n",
       "      <th>sensor_16</th>\n",
       "      <th>sensor_17</th>\n",
       "      <th>sensor_18</th>\n",
       "      <th>sensor_19</th>\n",
       "      <th>sensor_20</th>\n",
       "      <th>sensor_21</th>\n",
       "      <th>rul</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>-0.0007</td>\n",
       "      <td>-0.0004</td>\n",
       "      <td>100.0</td>\n",
       "      <td>518.67</td>\n",
       "      <td>641.82</td>\n",
       "      <td>1589.70</td>\n",
       "      <td>1400.60</td>\n",
       "      <td>14.62</td>\n",
       "      <td>...</td>\n",
       "      <td>2388.02</td>\n",
       "      <td>8138.62</td>\n",
       "      <td>8.4195</td>\n",
       "      <td>0.03</td>\n",
       "      <td>392</td>\n",
       "      <td>2388</td>\n",
       "      <td>100.0</td>\n",
       "      <td>39.06</td>\n",
       "      <td>23.4190</td>\n",
       "      <td>135</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>1</td>\n",
       "      <td>2</td>\n",
       "      <td>0.0019</td>\n",
       "      <td>-0.0003</td>\n",
       "      <td>100.0</td>\n",
       "      <td>518.67</td>\n",
       "      <td>642.15</td>\n",
       "      <td>1591.82</td>\n",
       "      <td>1403.14</td>\n",
       "      <td>14.62</td>\n",
       "      <td>...</td>\n",
       "      <td>2388.07</td>\n",
       "      <td>8131.49</td>\n",
       "      <td>8.4318</td>\n",
       "      <td>0.03</td>\n",
       "      <td>392</td>\n",
       "      <td>2388</td>\n",
       "      <td>100.0</td>\n",
       "      <td>39.00</td>\n",
       "      <td>23.4236</td>\n",
       "      <td>135</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>2 rows × 27 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "   unit  time_cycles  op_setting_1  op_setting_2  op_setting_3  sensor_1  \\\n",
       "0     1            1       -0.0007       -0.0004         100.0    518.67   \n",
       "1     1            2        0.0019       -0.0003         100.0    518.67   \n",
       "\n",
       "   sensor_2  sensor_3  sensor_4  sensor_5  ...  sensor_13  sensor_14  \\\n",
       "0    641.82   1589.70   1400.60     14.62  ...    2388.02    8138.62   \n",
       "1    642.15   1591.82   1403.14     14.62  ...    2388.07    8131.49   \n",
       "\n",
       "   sensor_15  sensor_16  sensor_17  sensor_18  sensor_19  sensor_20  \\\n",
       "0     8.4195       0.03        392       2388      100.0      39.06   \n",
       "1     8.4318       0.03        392       2388      100.0      39.00   \n",
       "\n",
       "   sensor_21  rul  \n",
       "0    23.4190  135  \n",
       "1    23.4236  135  \n",
       "\n",
       "[2 rows x 27 columns]"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from utils import read_dataset, calculate_RUL, SENSOR_COLUMNS\n",
    "\n",
    "train, test, test_rul = read_dataset('FD001')\n",
    "\n",
    "train['rul'] = calculate_RUL(train, upper_threshold=135)\n",
    "\n",
    "print(f'train.shape = {train.shape}')\n",
    "train.head(2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "Before starting features engineering - we should remove sensors with constant values that we saw during Exploratory Data Analysis. For that, we create a transformer which drops features with variance lower than a `threshold` - see below:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T14:48:55.924496Z",
     "start_time": "2021-11-07T14:48:55.787766Z"
    }
   },
   "outputs": [],
   "source": [
    "from sklearn.base import BaseEstimator, TransformerMixin\n",
    "from sklearn.feature_selection import VarianceThreshold\n",
    "\n",
    "\n",
    "class LowVarianceFeaturesRemover(BaseEstimator, TransformerMixin):\n",
    "    def __init__(self, threshold=0):\n",
    "        self.threshold = threshold\n",
    "        self.selector = VarianceThreshold(threshold=threshold)\n",
    "\n",
    "    def fit(self, X):\n",
    "        self.selector.fit(X)\n",
    "        return self\n",
    "\n",
    "    def transform(self, X):\n",
    "        X_t = self.selector.transform(X)\n",
    "        droped_features = X.columns[~self.selector.get_support()]\n",
    "        print(f'Droped low variance features: {droped_features.to_list()}')\n",
    "        return pd.DataFrame(X_t, columns=self.selector.get_feature_names_out())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "Scaling per engines from `1-EDA.ipynb` notebook:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T14:48:55.934970Z",
     "start_time": "2021-11-07T14:48:55.927792Z"
    }
   },
   "outputs": [],
   "source": [
    "from utils import SENSOR_COLUMNS\n",
    "\n",
    "\n",
    "class ScalePerEngine(BaseEstimator, TransformerMixin):\n",
    "    '''\n",
    "    Scale individual engines time series with respect to its start.\n",
    "    Substract firts `n_first_cycles` AVG values from time series.\n",
    "    '''\n",
    "    def __init__(self, n_first_cycles=20, sensors_columns=SENSOR_COLUMNS):\n",
    "        self.n_first_cycles = n_first_cycles\n",
    "        self.sensors_columns = sensors_columns\n",
    "\n",
    "    def fit(self, X):\n",
    "        return self\n",
    "\n",
    "    def transform(self, X):\n",
    "        self.sensors_columns = [x for x in X.columns if x in self.sensors_columns]\n",
    "\n",
    "        init_sensors_avg = X[X['time_cycles'] <= self.n_first_cycles] \\\n",
    "            .groupby(by=['unit'])[self.sensors_columns] \\\n",
    "            .mean() \\\n",
    "            .reset_index()\n",
    "\n",
    "        X_t = X[X['time_cycles'] > self.n_first_cycles].merge(\n",
    "            init_sensors_avg,\n",
    "            on=['unit'], how='left', suffixes=('', '_init_v')\n",
    "        )\n",
    "\n",
    "        for SENSOR in self.sensors_columns:\n",
    "            X_t[SENSOR] = X_t[SENSOR] - X_t['{}_init_v'.format(SENSOR)]\n",
    "\n",
    "        drop_columns = X_t.columns.str.endswith('init_v')\n",
    "        return X_t[X_t.columns[~drop_columns]]\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "# Features Engineering with TSfresh"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "**TSFresh** (https://github.com/blue-yonder/tsfresh) extracts a large number of features from time series and has a built-in features filtering procedure.\n",
    "\n",
    "You can check a list of features TSfresh extracts here: \n",
    "* https://tsfresh.readthedocs.io/en/latest/api/tsfresh.feature_extraction.html#module-tsfresh.feature_extraction.feature_calculators.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "**NOTE:**\n",
    "\n",
    "While TSfresh automates a lot of features engineering - it also **tends to extract highly correlated features** - see [1]. In order to eliminate this problem, the authors of the paper and the library propose to perform Principal Component Analysis. PCA can be done right after TSfresh generates a candidate features and before the features selection procedure. Or PCA can be performed after TSFresh does features selection. We stick to the first option and apply features selection on principal components.\n",
    "\n",
    "\n",
    "[1] Christ, M., Kempa-Liehr, A.W. and Feindt, M. (2016). Distributed and parallel time series feature extraction for industrial big data applications. ArXiv e-prints: 1610.07717 URL: http://adsabs.harvard.edu/abs/2016arXiv161007717C\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "<br>\n",
    "Our tsfresh features engineering pipeline is:\n",
    "\n",
    "1. transform the original time series into sliding windows of length 30 (see `roll_time_series` method below)\n",
    "\n",
    "2. calculate TSFresh features from the predefined list of features (see `tsfresh_calc` below)\n",
    "\n",
    "3. perform PCA\n",
    "\n",
    "4. do features selection with TSfresh "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T14:48:57.112062Z",
     "start_time": "2021-11-07T14:48:55.937129Z"
    }
   },
   "outputs": [],
   "source": [
    "from tsfresh.utilities.dataframe_functions import roll_time_series\n",
    "\n",
    "\n",
    "class RollTimeSeries(BaseEstimator, TransformerMixin):\n",
    "    def __init__(self, min_timeshift, max_timeshift, rolling_direction):\n",
    "        self.min_timeshift = min_timeshift\n",
    "        self.max_timeshift = max_timeshift\n",
    "        self.rolling_direction = rolling_direction\n",
    "        \n",
    "    def fit(self, X):\n",
    "        return self\n",
    "    \n",
    "    def transform(self, X):\n",
    "        _start = datetime.now()\n",
    "        print('Start Rolling TS')\n",
    "        X_t = roll_time_series(\n",
    "                    X, column_id='unit', column_sort='time_cycles', \n",
    "                    rolling_direction=self.rolling_direction,\n",
    "                    min_timeshift=self.min_timeshift,\n",
    "                    max_timeshift=self.max_timeshift)\n",
    "        print(f'Done Rolling TS in {datetime.now() - _start}')\n",
    "        return X_t"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "<br>\n",
    "\n",
    "**NOTE:** TSfresh suggests a couple of options to choose a set of features from - e.g. `MinimalFCParameters`, `ComprehensiveFCParameters`, `EfficientFCParameters` - see https://tsfresh.readthedocs.io/en/latest/text/feature_extraction_settings.html.\n",
    "\n",
    "However, after reviewing the proposed features lists, a set of specific features was shortlisted - see below: "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T14:48:57.126190Z",
     "start_time": "2021-11-07T14:48:57.114813Z"
    }
   },
   "outputs": [],
   "source": [
    "tsfresh_calc = {\n",
    "    'mean_change': None,\n",
    "    'mean': None,\n",
    "    'standard_deviation': None,\n",
    "    'root_mean_square': None,\n",
    "    'last_location_of_maximum': None,\n",
    "    'first_location_of_maximum': None,\n",
    "    'last_location_of_minimum': None,\n",
    "    'first_location_of_minimum': None,\n",
    "    'maximum': None,\n",
    "    'minimum': None,\n",
    "    'time_reversal_asymmetry_statistic': [{'lag': 1}, {'lag': 2}, {'lag': 3}],\n",
    "    'c3': [{'lag': 1}, {'lag': 2}, {'lag': 3}],\n",
    "    'cid_ce': [{'normalize': True}, {'normalize': False}],\n",
    "    'autocorrelation': [{'lag': 0},\n",
    "        {'lag': 1},\n",
    "        {'lag': 2},\n",
    "        {'lag': 3},],\n",
    "    'partial_autocorrelation': [{'lag': 0},\n",
    "        {'lag': 1},\n",
    "        {'lag': 2},\n",
    "        {'lag': 3},],\n",
    "    'linear_trend': [\n",
    "        {'attr': 'intercept'},\n",
    "        {'attr': 'slope'},\n",
    "        {'attr': 'stderr'}],\n",
    "    'augmented_dickey_fuller': [{'attr': 'teststat'},\n",
    "        {'attr': 'pvalue'},\n",
    "        {'attr': 'usedlag'}],\n",
    "    'linear_trend_timewise': [\n",
    "        {'attr': 'intercept'},\n",
    "        {'attr': 'slope'}],\n",
    "    'lempel_ziv_complexity': [{'bins': 2},\n",
    "        {'bins': 3},\n",
    "        {'bins': 5},\n",
    "        {'bins': 10},\n",
    "        {'bins': 100}],\n",
    "    'permutation_entropy': [{'tau': 1, 'dimension': 3},\n",
    "        {'tau': 1, 'dimension': 4},\n",
    "        {'tau': 1, 'dimension': 5},\n",
    "        {'tau': 1, 'dimension': 6},\n",
    "        {'tau': 1, 'dimension': 7}],\n",
    "    'fft_coefficient': [\n",
    "        {'coeff': 0, 'attr': 'abs'},\n",
    "        {'coeff': 1, 'attr': 'abs'},\n",
    "        {'coeff': 2, 'attr': 'abs'},\n",
    "        {'coeff': 3, 'attr': 'abs'},\n",
    "        {'coeff': 4, 'attr': 'abs'},\n",
    "        {'coeff': 5, 'attr': 'abs'},\n",
    "        {'coeff': 6, 'attr': 'abs'},\n",
    "        {'coeff': 7, 'attr': 'abs'},\n",
    "        {'coeff': 8, 'attr': 'abs'},\n",
    "        {'coeff': 9, 'attr': 'abs'},\n",
    "        {'coeff': 10, 'attr': 'abs'},],\n",
    "    'fft_aggregated': [{'aggtype': 'centroid'},\n",
    "        {'aggtype': 'variance'},\n",
    "        {'aggtype': 'skew'},\n",
    "        {'aggtype': 'kurtosis'}],\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T14:48:57.131313Z",
     "start_time": "2021-11-07T14:48:57.128300Z"
    }
   },
   "outputs": [],
   "source": [
    "from tsfresh import extract_features, select_features\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "from sklearn.decomposition import PCA\n",
    "\n",
    "from utils import calculate_RUL, SENSOR_COLUMNS"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T14:48:57.145858Z",
     "start_time": "2021-11-07T14:48:57.133630Z"
    }
   },
   "outputs": [],
   "source": [
    "class TSFreshFeaturesExtractor(BaseEstimator, TransformerMixin):\n",
    "    def __init__(self, calc=tsfresh_calc):\n",
    "        self.calc = calc\n",
    "        \n",
    "    def _clean_features(self, X):\n",
    "        old_shape = X.shape\n",
    "        X_t = X.T.drop_duplicates().T\n",
    "        print(f'Droped {old_shape[1] - X_t.shape[1]} duplicate features')\n",
    "\n",
    "        old_shape = X_t.shape\n",
    "        X_t = X_t.dropna(axis=1)\n",
    "        print(f'Droped {old_shape[1] - X_t.shape[1]} features with NA values')\n",
    "        return X_t\n",
    "    \n",
    "    def fit(self, X):\n",
    "        return self\n",
    "    \n",
    "    def transform(self, X):\n",
    "        _start = datetime.now()\n",
    "        print('Start Extracting Features')\n",
    "        X_t = extract_features(\n",
    "            X[['id', 'time_cycles'] \n",
    "               + X.columns[X.columns.str.startswith('sensor')].tolist()], \n",
    "            column_id='id', \n",
    "            column_sort='time_cycles', \n",
    "            default_fc_parameters=self.calc)\n",
    "        print(f'Done Extracting Features in {datetime.now() - _start}')\n",
    "        X_t = self._clean_features(X_t)\n",
    "        return X_t\n",
    "\n",
    "    \n",
    "class CustomPCA(BaseEstimator, TransformerMixin):\n",
    "    def __init__(self, n_components=None, random_state=None):\n",
    "        self.n_components = n_components\n",
    "        self.random_state = random_state\n",
    "        \n",
    "    def fit(self, X):\n",
    "        assert 'unit' not in X.columns, \"columns should be only features\"\n",
    "        self.ftr_columns = X.columns\n",
    "\n",
    "        self.scaler = StandardScaler()\n",
    "        X_sc = self.scaler.fit_transform(X[self.ftr_columns].values)\n",
    "\n",
    "        self.pca = PCA(n_components=self.n_components, random_state=self.random_state)\n",
    "        self.pca.fit_transform(X_sc)\n",
    "        return self\n",
    "\n",
    "    def transform(self, X):\n",
    "        X_sc = self.scaler.transform(X[self.ftr_columns].values)\n",
    "        X_pca = self.pca.transform(X_sc)        \n",
    "        return pd.DataFrame(X_pca, index=X.index)\n",
    "\n",
    "    \n",
    "\n",
    "class TSFreshFeaturesSelector(BaseEstimator, TransformerMixin):\n",
    "    def __init__(self, fdr_level=0.001):\n",
    "        self.fdr_level = fdr_level\n",
    "\n",
    "    def fit(self, X):\n",
    "        rul = calculate_RUL(\n",
    "            X.index.to_frame(name=['unit', 'time_cycles']).reset_index(drop=True),\n",
    "            upper_threshold=135)\n",
    "\n",
    "        X_t = select_features(\n",
    "            X, rul, fdr_level=self.fdr_level\n",
    "        )\n",
    "        self.selected_ftr = X_t.columns\n",
    "        \n",
    "        print(f'Selected {len(self.selected_ftr)} out of {X.shape[1]} features: '\n",
    "              f'{self.selected_ftr.to_list()}')\n",
    "        return self\n",
    "\n",
    "    def transform(self, X):\n",
    "        return X[self.selected_ftr]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T14:48:57.157029Z",
     "start_time": "2021-11-07T14:48:57.149117Z"
    }
   },
   "outputs": [],
   "source": [
    "from sklearn.pipeline import Pipeline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T14:48:57.165757Z",
     "start_time": "2021-11-07T14:48:57.161590Z"
    }
   },
   "outputs": [],
   "source": [
    "tsfresh_pipe = Pipeline([\n",
    "    # Cleaning\n",
    "    ('drop-low-variance', LowVarianceFeaturesRemover(threshold=0)),\n",
    "    \n",
    "    # Scaling and Preprocessing\n",
    "    ('scale-per-engine', ScalePerEngine(n_first_cycles=15, sensors_columns=SENSOR_COLUMNS)),\n",
    "    ('roll-time-series', RollTimeSeries(min_timeshift=29, max_timeshift=29, rolling_direction=1)),\n",
    "    \n",
    "    # TSFresh features engineering\n",
    "    ('extract-tsfresh-features', TSFreshFeaturesExtractor(calc=tsfresh_calc)),\n",
    "    ('PCA', CustomPCA(n_components=40)),\n",
    "    ('features-selection', TSFreshFeaturesSelector(fdr_level=0.001)),\n",
    "])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T15:01:48.938067Z",
     "start_time": "2021-11-07T14:48:57.168093Z"
    },
    "tags": [
     "remove"
    ]
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Droped low variance features: ['op_setting_3', 'sensor_1', 'sensor_5', 'sensor_10', 'sensor_16', 'sensor_18', 'sensor_19']\n",
      "Start Rolling TS\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Rolling: 100%|███████████████████████████████████████████████████████████████████| 20/20 [00:10<00:00,  1.90it/s]\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Done Rolling TS in 0:00:11.210756\n",
      "Start Extracting Features\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Feature Extraction: 100%|████████████████████████████████████████████████████████| 20/20 [11:57<00:00, 35.86s/it]\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Done Extracting Features in 0:12:29.266573\n",
      "Droped 19 duplicate features\n",
      "Droped 11 features with NA values\n",
      "Selected 9 out of 40 features: [0, 2, 3, 4, 5, 1, 10, 31, 19]\n"
     ]
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th>0</th>\n",
       "      <th>2</th>\n",
       "      <th>3</th>\n",
       "      <th>4</th>\n",
       "      <th>5</th>\n",
       "      <th>1</th>\n",
       "      <th>10</th>\n",
       "      <th>31</th>\n",
       "      <th>19</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th rowspan=\"2\" valign=\"top\">1.0</th>\n",
       "      <th>45.0</th>\n",
       "      <td>-10.987398</td>\n",
       "      <td>-0.982110</td>\n",
       "      <td>6.989796</td>\n",
       "      <td>-2.297053</td>\n",
       "      <td>-0.054998</td>\n",
       "      <td>0.065238</td>\n",
       "      <td>1.771415</td>\n",
       "      <td>1.013616</td>\n",
       "      <td>-1.271091</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>46.0</th>\n",
       "      <td>-10.781827</td>\n",
       "      <td>-1.067625</td>\n",
       "      <td>6.895214</td>\n",
       "      <td>-2.323534</td>\n",
       "      <td>-0.224924</td>\n",
       "      <td>0.476905</td>\n",
       "      <td>1.526736</td>\n",
       "      <td>1.240907</td>\n",
       "      <td>-1.646843</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                 0         2         3         4         5         1   \\\n",
       "1.0 45.0 -10.987398 -0.982110  6.989796 -2.297053 -0.054998  0.065238   \n",
       "    46.0 -10.781827 -1.067625  6.895214 -2.323534 -0.224924  0.476905   \n",
       "\n",
       "                10        31        19  \n",
       "1.0 45.0  1.771415  1.013616 -1.271091  \n",
       "    46.0  1.526736  1.240907 -1.646843  "
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "train_tsfresh_ftr = tsfresh_pipe.fit_transform(train)\n",
    "\n",
    "train_tsfresh_ftr.head(2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "# Features Engineering with ROCKET"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "ROCKET (Random Convolutional Kernel Transform) transforms time series into features using random convolutional kernels (by default 10000 kernels) and applying max pooling and proportion of positive values (which produces 2 features per kernel). This results in (by default) 20000 features.  \n",
    "\n",
    "[1] Dempster A., Petitjean F. and Webb G., \"ROCKET: Exceptionally fast and accurate time series\n",
    "classification using random convolutional kernels\": https://arxiv.org/pdf/1910.13051.pdf"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "Note: here we add extra step for scaling the time series before running ROCKET - see below:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T15:01:48.947471Z",
     "start_time": "2021-11-07T15:01:48.940817Z"
    }
   },
   "outputs": [],
   "source": [
    "from sklearn.preprocessing import StandardScaler\n",
    "\n",
    "\n",
    "class CustomStandardScaler(BaseEstimator, TransformerMixin):\n",
    "    def fit(self, X):\n",
    "        self.scale_columns = X.columns[X.columns.str.startswith('sensor')\n",
    "                                       | X.columns.str.startswith('op_setting')]\n",
    "        self.scaler = StandardScaler()\n",
    "        self.scaler.fit(X[self.scale_columns])\n",
    "        return self\n",
    "    \n",
    "    def transform(self, X):\n",
    "        X_t = self.scaler.transform(X[self.scale_columns])\n",
    "        df = pd.concat([\n",
    "                X[[c for c in X.columns if c not in self.scale_columns]],\n",
    "                pd.DataFrame(X_t, columns=self.scale_columns)\n",
    "            ], axis=1).sort_values(['unit', 'time_cycles']).reset_index(drop=True)\n",
    "        return df\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T15:01:48.986316Z",
     "start_time": "2021-11-07T15:01:48.949788Z"
    }
   },
   "outputs": [],
   "source": [
    "from datetime import datetime\n",
    "\n",
    "from sktime.datatypes._panel._convert import from_multi_index_to_nested\n",
    "\n",
    "\n",
    "class TransformTS2Nested(BaseEstimator, TransformerMixin):\n",
    "    '''\n",
    "    See https://www.sktime.org/en/stable/examples/loading_data.html#Using-multi-indexed-pandas-DataFrames\n",
    "    '''\n",
    "    def fit(self, X):\n",
    "        return self\n",
    "    \n",
    "    def transform(self, X):\n",
    "        assert 'id' in X.columns, \"X should have `id` column\"\n",
    "        X_mi = X.copy()\n",
    "        X_mi.index = pd.MultiIndex.from_frame(X[['id', 'time_cycles']])\n",
    "        X_mi = X_mi.drop(columns=['rul', 'id'], errors='ignore')\n",
    "        \n",
    "        _start = datetime.now()\n",
    "        print('Start Converting multi-index DF to sktime nested DF')\n",
    "        X_nested = from_multi_index_to_nested(X_mi, instance_index='id')\n",
    "        \n",
    "        X_nested.index = pd.MultiIndex.from_frame(\n",
    "            pd.DataFrame({'unit': X_nested['unit'].apply(lambda x: x.unique()[0]),\n",
    "                          'time_cycles': X_nested['time_cycles'].apply(lambda x: x.iloc[-1])})\n",
    "        )\n",
    "        X_nested = X_nested.drop(columns=['unit', 'time_cycles'])\n",
    "        print(f'Converted multi-index DF to sktime nested DF in {datetime.now() - _start}')\n",
    "\n",
    "        return X_nested\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T15:01:49.596056Z",
     "start_time": "2021-11-07T15:01:48.989189Z"
    }
   },
   "outputs": [],
   "source": [
    "from sktime.transformations.panel.rocket import Rocket\n",
    "\n",
    "\n",
    "class RocketTransform(BaseEstimator, TransformerMixin):\n",
    "    '''\n",
    "    See https://www.sktime.org/en/stable/examples/rocket.html\n",
    "    '''\n",
    "    def __init__(self, num_kernels=10000, normalise=False, n_jobs=-1, random_state=None):\n",
    "        self.num_kernels = num_kernels\n",
    "        self.normalise = normalise\n",
    "        self.n_jobs = n_jobs\n",
    "        self.random_state = random_state\n",
    "\n",
    "    def fit(self, X):\n",
    "        self.rocket = Rocket(\n",
    "            num_kernels=self.num_kernels, \n",
    "            normalise=self.normalise, \n",
    "            n_jobs=self.n_jobs,\n",
    "            random_state=self.random_state)\n",
    "        self.rocket.fit(X)\n",
    "        return self\n",
    "    \n",
    "    def transform(self, X):\n",
    "        rocket_t = self.rocket.transform(X)\n",
    "        rocket_t.index = X.index\n",
    "        return rocket_t\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T15:01:49.602808Z",
     "start_time": "2021-11-07T15:01:49.598109Z"
    }
   },
   "outputs": [],
   "source": [
    "from sklearn.pipeline import Pipeline\n",
    "\n",
    "from utils import SENSOR_COLUMNS\n",
    "\n",
    "\n",
    "rocket_pipe = Pipeline([\n",
    "    # Cleaning constant features\n",
    "    ('drop-low-variance', LowVarianceFeaturesRemover()),\n",
    "    \n",
    "    # Scaling sensors values\n",
    "    ('scale-per-engine', ScalePerEngine(n_first_cycles=15, sensors_columns=SENSOR_COLUMNS)),\n",
    "    ('scale', CustomStandardScaler()),\n",
    "    \n",
    "    # Preprocessing for ROCKET\n",
    "    ('roll-time-series', RollTimeSeries(min_timeshift=29, max_timeshift=29, rolling_direction=1)),\n",
    "    ('nest-time-series', TransformTS2Nested()),\n",
    "    \n",
    "    # Features Engineering w ROCKET\n",
    "    ('rocket', RocketTransform(num_kernels=10000, normalise=False, random_state=2021))\n",
    "])\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2021-11-07T15:14:32.017227Z",
     "start_time": "2021-11-07T15:01:49.604957Z"
    },
    "tags": [
     "remove"
    ]
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Droped low variance features: ['op_setting_3', 'sensor_1', 'sensor_5', 'sensor_10', 'sensor_16', 'sensor_18', 'sensor_19']\n",
      "Start Rolling TS\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Rolling: 100%|███████████████████████████████████████████████████████████████████| 20/20 [00:08<00:00,  2.27it/s]\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Done Rolling TS in 0:00:09.353206\n",
      "Start Converting multi-index DF to sktime nested DF\n",
      "Converted multi-index DF to sktime nested DF in 0:07:41.130546\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "(16231, 20000)"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "train_rocket_ftr = rocket_pipe.fit_transform(train)\n",
    "\n",
    "train_rocket_ftr.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that here we end up with `p >> n` case."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "tags": [
     "remove"
    ]
   },
   "source": [
    "---"
   ]
  }
 ],
 "metadata": {
  "hide_input": false,
  "kernelspec": {
   "display_name": "venv",
   "language": "python",
   "name": "venv"
  },
  "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.8.5"
  },
  "latex_envs": {
   "LaTeX_envs_menu_present": true,
   "autoclose": false,
   "autocomplete": true,
   "bibliofile": "biblio.bib",
   "cite_by": "apalike",
   "current_citInitial": 1,
   "eqLabelWithNumbers": true,
   "eqNumInitial": 1,
   "hotkeys": {
    "equation": "Ctrl-E",
    "itemize": "Ctrl-I"
   },
   "labels_anchors": false,
   "latex_user_defs": false,
   "report_style_numbering": false,
   "user_envs_cfg": false
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {
    "height": "calc(100% - 180px)",
    "left": "10px",
    "top": "150px",
    "width": "255.608px"
   },
   "toc_section_display": true,
   "toc_window_display": true
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}