{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import seaborn as sns\nimport matplotlib.pyplot as plt\nfrom catboost import CatBoostRegressor, Pool\nfrom sklearn.preprocessing import StandardScaler\nfrom sklearn.model_selection import GridSearchCV, KFold, train_test_split\nfrom sklearn.svm import SVR, NuSVR\nfrom sklearn.kernel_ridge import KernelRidge\nimport pandas as pd\nimport numpy as np\nimport os\nimport gc\nimport warnings\nwarnings.filterwarnings(\"ignore\")\n\nDATA_DIR = \"../input\"\nTEST_DIR = r'../input/test'\n\nld = os.listdir(TEST_DIR)\nsizes = np.zeros(len(ld))\n\nfrom scipy.signal import hilbert\nfrom scipy.signal import hann\nfrom scipy.signal import convolve\nfrom scipy.stats import pearsonr\nfrom scipy import stats\nfrom sklearn.kernel_ridge import KernelRidge\n\nimport lightgbm as lgb\nfrom sklearn.metrics import mean_squared_error\nfrom sklearn.metrics import mean_absolute_error\n\nfrom tsfresh.feature_extraction import feature_calculators\nfrom tqdm import tqdm\n\n%matplotlib inline\nsns.set_style('darkgrid')","metadata":{"execution":{"iopub.status.busy":"2023-04-10T16:59:32.394028Z","iopub.execute_input":"2023-04-10T16:59:32.394663Z","iopub.status.idle":"2023-04-10T16:59:34.314402Z","shell.execute_reply.started":"2023-04-10T16:59:32.394609Z","shell.execute_reply":"2023-04-10T16:59:34.313449Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def classic_sta_lta(x, length_sta, length_lta):\n    \n    sta = np.cumsum(x ** 2)\n\n    # Zamiana na float\n    sta = np.require(sta, dtype=np.float)\n\n    # Kopia dla LTA\n    lta = sta.copy()\n\n    # Obliczanie STA i LTA\n    sta[length_sta:] = sta[length_sta:] - sta[:-length_sta]\n    sta /= length_sta\n    lta[length_lta:] = lta[length_lta:] - lta[:-length_lta]\n    lta /= length_lta\n\n    # Uzupełnienie zerami\n    sta[:length_lta - 1] = 0\n\n    # Aby nie dzielić przez 0 ustawiamy 0 na małe liczby typu float\n    dtiny = np.finfo(0.0).tiny\n    idx = lta < dtiny\n    lta[idx] = dtiny\n\n    return sta / lta","metadata":{"execution":{"iopub.status.busy":"2023-04-10T16:59:34.31573Z","iopub.execute_input":"2023-04-10T16:59:34.316035Z","iopub.status.idle":"2023-04-10T16:59:34.323527Z","shell.execute_reply.started":"2023-04-10T16:59:34.315987Z","shell.execute_reply":"2023-04-10T16:59:34.321245Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def calc_change_rate(x):\n    change = (np.diff(x) / x[:-1]).values\n    change = change[np.nonzero(change)[0]]\n    change = change[~np.isnan(change)]\n    change = change[change != -np.inf]\n    change = change[change != np.inf]\n    return np.mean(change)","metadata":{"execution":{"iopub.status.busy":"2023-04-10T16:59:34.324727Z","iopub.execute_input":"2023-04-10T16:59:34.325191Z","iopub.status.idle":"2023-04-10T16:59:34.333498Z","shell.execute_reply.started":"2023-04-10T16:59:34.325125Z","shell.execute_reply":"2023-04-10T16:59:34.332777Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"percentiles = [1, 5, 10, 20, 25, 30, 40, 50, 60, 70, 75, 80, 90, 95, 99]\nhann_windows = [50, 150, 1500, 15000]\nspans = [300, 3000, 30000, 50000]\nwindows = [10, 50, 100, 500, 1000, 10000]\nborders = list(range(-4000, 4001, 1000))\npeaks = [10, 20, 50, 100]\ncoefs = [1, 5, 10, 50, 100]\nlags = [10, 100, 1000, 10000]\nautocorr_lags = [5, 10, 50, 100, 500, 1000, 5000, 10000]","metadata":{"execution":{"iopub.status.busy":"2023-04-10T16:59:34.336215Z","iopub.execute_input":"2023-04-10T16:59:34.336671Z","iopub.status.idle":"2023-04-10T16:59:34.343844Z","shell.execute_reply.started":"2023-04-10T16:59:34.336621Z","shell.execute_reply":"2023-04-10T16:59:34.343104Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def gen_features(x, zero_mean=False):\n    if zero_mean==True:\n        x = x-x.mean()\n    strain = {}\n    strain['mean'] = x.mean()\n    strain['std']=x.std()\n    strain['max']=x.max()\n    strain['kurtosis']=x.kurtosis()\n    strain['skew']=x.skew()\n    zc = np.fft.fft(x)\n    realFFT = np.real(zc)\n    imagFFT = np.imag(zc)\n    strain['min']=x.min()\n    strain['sum']=x.sum()\n    strain['mad']=x.mad()\n    strain['median']=x.median()\n    \n    strain['mean_change_abs'] = np.mean(np.diff(x))\n    strain['mean_change_rate'] = np.mean(np.nonzero((np.diff(x) / x[:-1]))[0])\n    strain['abs_max'] = np.abs(x).max()\n    strain['abs_min'] = np.abs(x).min()\n    \n    strain['avg_first_50000'] = x[:50000].mean()\n    strain['avg_last_50000'] = x[-50000:].mean()\n    strain['avg_first_10000'] = x[:10000].mean()\n    strain['avg_last_10000'] = x[-10000:].mean()\n    \n    strain['min_first_50000'] = x[:50000].min()\n    strain['min_last_50000'] = x[-50000:].min()\n    strain['min_first_10000'] = x[:10000].min()\n    strain['min_last_10000'] = x[-10000:].min()\n    \n    strain['max_first_50000'] = x[:50000].max()\n    strain['max_last_50000'] = x[-50000:].max()\n    strain['max_first_10000'] = x[:10000].max()\n    strain['max_last_10000'] = x[-10000:].max()\n    \n    strain['max_to_min'] = x.max() / np.abs(x.min())\n    strain['max_to_min_diff'] = x.max() - np.abs(x.min())\n    strain['count_big'] = len(x[np.abs(x) > 500])\n           \n    strain['mean_change_rate_first_50000'] = calc_change_rate(x[:50000])\n    strain['mean_change_rate_last_50000'] = calc_change_rate(x[-50000:])\n    strain['mean_change_rate_first_10000'] = calc_change_rate(x[:10000])\n    strain['mean_change_rate_last_10000'] = calc_change_rate(x[-10000:])\n    \n    strain['q95'] = np.quantile(x, 0.95)\n    strain['q99'] = np.quantile(x, 0.99)\n    strain['q05'] = np.quantile(x, 0.05)\n    strain['q01'] = np.quantile(x, 0.01)\n    \n    strain['abs_q95'] = np.quantile(np.abs(x), 0.95)\n    strain['abs_q99'] = np.quantile(np.abs(x), 0.99)\n    strain['abs_q05'] = np.quantile(np.abs(x), 0.05)\n    strain['abs_q01'] = np.quantile(np.abs(x), 0.01)\n    \n    for autocorr_lag in autocorr_lags:\n        strain['autocorrelation_' + str(autocorr_lag)] = feature_calculators.autocorrelation(x, autocorr_lag)\n    \n    # percentiles on original and absolute values\n    for p in percentiles:\n        strain['percentile_'+str(p)] = np.percentile(x, p)\n        strain['abs_percentile_'+str(p)] = np.percentile(np.abs(x), p)\n    \n    strain['abs_mean'] = np.abs(x).mean()\n    strain['abs_std'] = np.abs(x).std()\n    \n    strain['quantile_0.95']=np.quantile(x, 0.95)\n    strain['quantile_0.99']=np.quantile(x, 0.99)\n    strain['quantile_0.05']=np.quantile(x, 0.05)\n    strain['realFFT_mean']=realFFT.mean()\n    strain['realFFT_std']=realFFT.std()\n    strain['realFFT_max']=realFFT.max()\n    strain['realFFT_min']=realFFT.min()\n    strain['imagFFT_mean']=imagFFT.mean()\n    strain['imagFFT_std']=realFFT.std()\n    strain['imagFFT_max']=realFFT.max()\n    strain['imaglFFT_min']=realFFT.min()\n    \n    strain['std_first_50000']=x[:50000].std()\n    strain['std_last_50000']=x[-50000:].std()\n    strain['std_first_25000']=x[:25000].std()\n    strain['std_last_25000']=x[-25000:].std()\n    strain['std_first_10000']=x[:10000].std()\n    strain['std_last_10000']=x[-10000:].std()\n    strain['std_first_5000']=x[:5000].std()\n    strain['std_last_5000']=x[-5000:].std()\n        \n    strain['Hilbert_mean'] = np.abs(hilbert(x)).mean()\n    strain['Hann_window_mean'] = (convolve(x, hann(150), mode='same') / sum(hann(150))).mean()\n    strain['classic_sta_lta1_mean'] = classic_sta_lta(x, 500, 10000).mean()\n    strain['classic_sta_lta2_mean'] = classic_sta_lta(x, 5000, 100000).mean()\n    strain['classic_sta_lta3_mean'] = classic_sta_lta(x, 3333, 6666).mean()\n    strain['classic_sta_lta4_mean'] = classic_sta_lta(x, 10000, 25000).mean()\n    strain['classic_sta_lta6_mean'] = classic_sta_lta(x, 100, 5000).mean()\n    strain['classic_sta_lta8_mean'] = classic_sta_lta(x, 4000, 10000).mean()\n    strain['Moving_average_700_mean'] = x.rolling(window=700).mean().mean(skipna=True)\n    moving_average_700_mean = x.rolling(window=700).mean().mean(skipna=True)\n    ewma = pd.Series.ewm\n    strain['exp_Moving_average_300_mean'] = (ewma(x, span=300).mean()).mean(skipna=True)\n    strain['exp_Moving_average_3000_mean'] = ewma(x, span=3000).mean().mean(skipna=True)\n    strain['exp_Moving_average_30000_mean'] = ewma(x, span=30000).mean().mean(skipna=True)\n    no_of_std = 3\n    strain['MA_700MA_std_mean'] = x.rolling(window=700).std().mean()\n    strain['MA_1000MA_std_mean'] = x.rolling(window=1000).std().mean()\n    \n    strain['iqr'] = np.subtract(*np.percentile(x, [75, 25]))\n    strain['q999'] = np.quantile(x,0.999)\n    strain['q001'] = np.quantile(x,0.001)\n    strain['ave10'] = stats.trim_mean(x, 0.1)\n        \n    for window in windows:\n        x_roll_std = x.rolling(window).std().dropna().values\n        x_roll_mean = x.rolling(window).mean().dropna().values\n        \n        strain['ave_roll_std_' + str(window)] = x_roll_std.mean()\n        strain['std_roll_std_' + str(window)] = x_roll_std.std()\n        strain['max_roll_std_' + str(window)] = x_roll_std.max()\n        strain['min_roll_std_' + str(window)] = x_roll_std.min()\n        strain['q01_roll_std_' + str(window)] = np.quantile(x_roll_std, 0.01)\n        strain['q05_roll_std_' + str(window)] = np.quantile(x_roll_std, 0.05)\n        strain['q95_roll_std_' + str(window)] = np.quantile(x_roll_std, 0.95)\n        strain['q99_roll_std_' + str(window)] = np.quantile(x_roll_std, 0.99)\n        strain['av_change_abs_roll_std_' + str(window)] = np.mean(np.diff(x_roll_std))\n        strain['av_change_rate_roll_std_' + str(window)] = np.mean(np.nonzero((np.diff(x_roll_std) / x_roll_std[:-1]))[0])\n        strain['abs_max_roll_std_' + str(window)] = np.abs(x_roll_std).max()\n        \n        for p in percentiles:\n            strain['percentile_roll_std_' + str(p) + '_window_' + str(window)] = np.percentile(x_roll_std, p)\n            strain['percentile_roll_mean_' + str(p) + '_window_' + str(window)] = np.percentile(x_roll_mean, p)\n        \n        strain['ave_roll_mean_' + str(window)] = x_roll_mean.mean()\n        strain['std_roll_mean_' + str(window)] = x_roll_mean.std()\n        strain['max_roll_mean_' + str(window)] = x_roll_mean.max()\n        strain['min_roll_mean_' + str(window)] = x_roll_mean.min()\n        strain['q01_roll_mean_' + str(window)] = np.quantile(x_roll_mean, 0.01)\n        strain['q05_roll_mean_' + str(window)] = np.quantile(x_roll_mean, 0.05)\n        strain['q95_roll_mean_' + str(window)] = np.quantile(x_roll_mean, 0.95)\n        strain['q99_roll_mean_' + str(window)] = np.quantile(x_roll_mean, 0.99)\n        strain['av_change_abs_roll_mean_' + str(window)] = np.mean(np.diff(x_roll_mean))\n        strain['av_change_rate_roll_mean_' + str(window)] = np.mean(np.nonzero((np.diff(x_roll_mean) / x_roll_mean[:-1]))[0])\n        strain['abs_max_roll_mean_' + str(window)] = np.abs(x_roll_mean).max()\n        \n        \n    return pd.Series(strain)","metadata":{"execution":{"iopub.status.busy":"2023-04-10T16:59:34.346647Z","iopub.execute_input":"2023-04-10T16:59:34.346888Z","iopub.status.idle":"2023-04-10T16:59:34.374104Z","shell.execute_reply.started":"2023-04-10T16:59:34.346842Z","shell.execute_reply":"2023-04-10T16:59:34.373225Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import random\n\n# Define the state and action spaces\nacoustic_data = [1,2,3,4,5,6,7,8,9,10,-1,-2,-3,-4,-5,-6,-7,-8,-9,-10,]\ntime_to_failure = [random.uniform(0, 10) for i in range(len(acoustic_data))] # generate random values for time_to_failure\nseg_id = ['seg1', 'seg2', 'seg3', 'seg4', 'seg5']\ntemp=0\nnew=True\n# Create the Q-learning table\nq_table = {}\nwhile new:\n    for s in acoustic_data:\n        for a in time_to_failure:\n            for i in seg_id:\n                q_table[(s, a, i)] = random.uniform(-1, 1)\n    random_idx = random.randint(0, len(time_to_failure)-1)\n    print(time_to_failure[random_idx])\n    temp=temp+1\n    if temp==100:\n        new=False","metadata":{"execution":{"iopub.status.busy":"2023-04-10T16:59:34.375843Z","iopub.execute_input":"2023-04-10T16:59:34.376065Z","iopub.status.idle":"2023-04-10T16:59:34.51871Z","shell.execute_reply.started":"2023-04-10T16:59:34.376023Z","shell.execute_reply":"2023-04-10T16:59:34.518078Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_df = pd.read_csv(os.path.join(DATA_DIR, 'train.csv'), iterator=True, chunksize=150_000, dtype={'acoustic_data': np.int16, 'time_to_failure': np.float32})\nX_train = pd.DataFrame()\ny_train = pd.Series()\n\nfor df in tqdm(train_df):\n    features = gen_features(df['acoustic_data'])\n    X_train = X_train.append(features, ignore_index=True)\n    y_train = y_train.append(pd.Series(df['time_to_failure'].values[-1]), ignore_index=True)\n\nX_train.head()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"del train_df\nX_test = pd.DataFrame()\n\nfor i, f in tqdm(enumerate(ld)):\n    df = pd.read_csv(os.path.join(TEST_DIR, f))\n    features = gen_features(df['acoustic_data'])\n    X_test = X_test.append(features, ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"corelations = np.abs(X_train.corrwith(y_train)).sort_values(ascending=False)\ncorelations_df = pd.DataFrame(data=corelations, columns=['corr'])\nprint(\"Number of high corelated values: \",corelations_df[corelations_df['corr']>=0.55]['corr'].count())\n\nhigh_corr = corelations_df[corelations_df['corr']>=0.55]\nprint(high_corr)\nhigh_corr_labels = high_corr.reset_index()['index'].values\n#print(high_corr_labels)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_train_high_corr = X_train[high_corr_labels]\nX_test_high_corr = X_test[high_corr_labels]","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.preprocessing import MinMaxScaler\n\nscaler = StandardScaler()\nscaler.fit(X_train_high_corr)\nX_train_scaled = pd.DataFrame(scaler.transform(X_train_high_corr), columns=X_train_high_corr.columns)\nX_test_scaled = pd.DataFrame(scaler.transform(X_test_high_corr), columns=X_test_high_corr.columns)\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"p_columns = []\np_corr = []\np_values = []\n\nfor col in X_train_scaled.columns:\n    p_columns.append(col)\n    p_corr.append(abs(pearsonr(X_train_scaled[col], y_train.values)[0]))\n    p_values.append(abs(pearsonr(X_train_scaled[col], y_train.values)[1]))\n\ndf = pd.DataFrame(data={'column': p_columns, 'corr': p_corr, 'p_value': p_values}, index=range(len(p_columns)))\ndf.sort_values(by=['corr', 'p_value'], inplace=True)\ndf.dropna(inplace=True)\ndf = df.loc[df['p_value'] <= 0.05]\n\ndrop_cols = []\n\nfor col in X_train_scaled.columns:\n    if col not in df['column'].tolist():\n        drop_cols.append(col)\n\nprint(drop_cols)\nprint('--------------------')\nprint(X_train_high_corr.columns.values)\n        \nX_train_scaled = X_train_scaled.drop(labels=drop_cols, axis=1)\nX_test_scaled = X_test_scaled.drop(labels=drop_cols, axis=1)\n\nX_train_scaled_minmax = X_train_scaled_minmax.drop(labels=drop_cols, axis=1)\nX_test_scaled_minmax = X_test_scaled_minmax.drop(labels=drop_cols, axis=1)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## DNN-QLearning Model","metadata":{}},{"cell_type":"code","source":"from keras.callbacks import ModelCheckpoint\nfrom keras.models import Sequential\nfrom keras.layers import Dense, Activation, Flatten, Dropout\nfrom sklearn.model_selection import train_test_split\nNN_model = Sequential()\n\n# The Input Layer :\nNN_model.add(Dense(128, kernel_initializer='RandomUniform',input_dim = X_train_scaled.shape[1], activation='relu'))\nNN_model.add(Dropout(0.5))\n# The Hidden Layers :\nNN_model.add(Dense(256, kernel_initializer='RandomUniform',activation='relu'))\nNN_model.add(Dropout(0.5))\nNN_model.add(Dense(256, kernel_initializer='RandomUniform',activation='relu'))\nNN_model.add(Dropout(0.5))\nNN_model.add(Dense(128, kernel_initializer='RandomUniform',activation='relu'))\n# The Output Layer :\nNN_model.add(Dense(1, kernel_initializer='RandomUniform',activation='linear'))\n\n# Compile the network :\nNN_model.compile(loss='mean_absolute_error', optimizer='adam', metrics=['mean_absolute_error'])\nNN_model.summary()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"checkpoint_name = 'Weights-{epoch:03d}--{val_loss:.5f}.hdf5' \ncheckpoint = ModelCheckpoint(checkpoint_name, monitor='val_loss', verbose = 1, save_best_only = True, mode ='auto')\ncallbacks_list = [checkpoint]\nNN_model.fit(X_train_scaled, y_train, epochs=100, batch_size=32, validation_split = 0.2, callbacks=callbacks_list)\npredictions_DNN = NN_model.predict(X_test_scaled)\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submission_DNN = pd.read_csv(os.path.join(DATA_DIR, 'sample_submission.csv'), dtype={\n    'acoustic_data': np.int16, 'time_to_failure': np.float32})\nsubmission_DNN['time_to_failure'] = predictions_DNN\nsubmission_DNN.to_csv('result_DNN.csv', index=False)","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}