016_Pruning.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Callin Switzer\n",
    "### 12 Dec 2018 -- Initial commit\n",
    "___\n",
    "### - Train Dense, Feedforward Neural Network with Keras\n",
    "### - Use data that was generated in Python"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "___"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from matplotlib import pyplot as plt\n",
    "%matplotlib inline\n",
    "from matplotlib.ticker import ScalarFormatter\n",
    "import os\n",
    "import pandas as pd\n",
    "import time\n",
    "from datetime import datetime\n",
    "import sys\n",
    "import seaborn as sns\n",
    "import scipy.io\n",
    "import subprocess\n",
    "import winsound\n",
    "import pickle\n",
    "import glob\n",
    "\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.preprocessing import StandardScaler, MinMaxScaler\n",
    "from sklearn.metrics import mean_squared_error, r2_score\n",
    "import tensorflow as tf\n",
    "print(\"TensorFlow successfully installed.\")\n",
    "if tf.test.is_built_with_cuda():\n",
    "    print(\"The installed version of TensorFlow includes GPU support.\")\n",
    "print(sys.version, \"\\n\")\n",
    "now = datetime.now()\n",
    "print(\"last run on \" + str(now))\n",
    "\n",
    "# define directories\n",
    "baseDir = os.getcwd()\n",
    "dataDir = r'D:\\MothSimulations\\11c-AggressiveManeuver\\Qstore\\hws_am_con'\n",
    "figDir = r'D:\\Dropbox\\AcademiaDropbox\\mothMachineLearning_dataAndFigs\\Figs'\n",
    "dataOutput = r'D:\\Dropbox\\AcademiaDropbox\\mothMachineLearning_dataAndFigs\\DataOutput'\n",
    "savedModels = r'D:\\Dropbox\\AcademiaDropbox\\mothMachineLearning_dataAndFigs\\savedModels'\n",
    "randomRawData = r'D:\\Dropbox\\AcademiaDropbox\\mothMachineLearning_dataAndFigs\\PythonGeneratedData\\TrainingData'\n",
    "if not os.path.exists(dataOutput):\n",
    "    os.mkdir(dataOutput)\n",
    "if not os.path.exists(savedModels):\n",
    "    os.mkdir(savedModels)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# import\n",
    "from keras.layers import Input, Dense, Dropout, BatchNormalization\n",
    "from keras.models import Model\n",
    "from keras import regularizers\n",
    "\n",
    "# Keras callcacks\n",
    "from keras.callbacks import EarlyStopping, ReduceLROnPlateau\n",
    "import keras.backend as K"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Make training and test set"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "# concatenate all files (only need to do this once)\n",
    "# it takes a few minutes\n",
    "all_files = glob.glob(os.path.join(randomRawData, \"*.csv\"))     \n",
    "df_from_each_file = (pd.read_csv(f) for f in all_files)\n",
    "concatenated_df   = pd.concat(df_from_each_file, ignore_index=True)\n",
    "\n",
    "# check for duplicates\n",
    "concatenated_df.drop_duplicates(inplace=True)\n",
    "concatenated_df.shape\n",
    "\n",
    "print(concatenated_df.shape)\n",
    "concatenated_df.tail()\n",
    "\n",
    "# save to hdf5\n",
    "concatenated_df.to_hdf(os.path.join(dataOutput, \"concatenatedRandomICs.h5\"), key = \"data\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# read in data\n",
    "trainDF = pd.read_hdf(os.path.join(dataOutput, \"concatenatedRandomICs.h5\"), key = \"data\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# refref: check for repeats!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# double check for repeats!\n",
    "np.sum(trainDF.iloc[:, [16,17,18]].duplicated()) # 0 means no repeats"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(trainDF.shape)\n",
    "trainDF.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# rename columns to be consistent with other code\n",
    "trainDF.rename(columns={\"x0\" : \"x_0\", \"y0\" : \"y_0\", \"phi0\" : \"phi_0\", \"theta0\" : \"theta_0\", \n",
    "                        \"xf\" : \"x_99\", \"yf\" : \"y_99\", \"phif\" : \"phi_99\", \"thetaf\" : \"theta_99\", \n",
    "                        \"xd0\" : \"x_dot_0\", \"yd0\" : \"y_dot_0\", \"phid0\" : \"phi_dot_0\", \"thetad0\": \"theta_dot_0\", \n",
    "                        \"xdf\" : \"x_dot_99\", \"ydf\": \"y_dot_99\", \"phidf\": \"phi_dot_99\", \"thetadf\": \"theta_dot_99\", \n",
    "                        \"tau0\" : \"tau\"}, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# convert to fx and fy\n",
    "trainDF[\"Fx\"] = trainDF.F * np.cos(trainDF.alpha)\n",
    "trainDF[\"Fy\"] = trainDF.F * np.sin(trainDF.alpha)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# make dataset\n",
    "X = trainDF.loc[:, [ \"phi_0\", \"theta_0\", \n",
    "                    \"x_99\", \"y_99\", \"phi_99\", \"theta_99\", \n",
    "                   \"x_dot_0\", \"y_dot_0\", \"phi_dot_0\", \"theta_dot_0\"]]\n",
    "\n",
    "Y = trainDF.loc[:, [\"Fx\", \"Fy\", \"tau\", \"x_dot_99\", \"y_dot_99\", \n",
    "                    \"phi_dot_99\", \"theta_dot_99\"] ]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "Y.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# test train split\n",
    "Xtrain, Xtest, Ytrain, Ytest = train_test_split(X, Y, test_size=0.2, random_state = 123)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# scale data \n",
    "scalerX = MinMaxScaler([-0.5, 0.5])  \n",
    "scalerY = MinMaxScaler([-0.5, 0.5])  \n",
    "\n",
    "# Don't cheat - fit only on training data\n",
    "scalerX.fit(Xtrain)  \n",
    "scalerY.fit(Ytrain) \n",
    "\n",
    "Xtrain_scaled = scalerX.transform(Xtrain)  \n",
    "Ytrain_scaled = scalerY.transform(Ytrain)  \n",
    "\n",
    "# apply same transformation to test data\n",
    "Xtest_scaled = scalerX.transform(Xtest)\n",
    "Ytest_scaled = scalerY.transform(Ytest)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pd.DataFrame(Xtrain_scaled, columns = X.columns).head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# save scalers, to be used on test set\n",
    "scalerfileX = 'scalerX.pkl'\n",
    "pickle.dump(scalerX, open(os.path.join(dataOutput, scalerfileX), 'wb'))\n",
    "\n",
    "scalerfileY = 'scalerY.pkl'\n",
    "pickle.dump(scalerY, open(os.path.join(dataOutput, scalerfileY), 'wb'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Keras models"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from keras.wrappers.scikit_learn import KerasRegressor\n",
    "from sklearn.model_selection import RandomizedSearchCV"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "K.clear_session()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# refref: start with small network and then build up\n",
    "# refref: start with large network and prune"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create network\n",
    "def create_network(optimizer = 'rmsprop', \n",
    "                    numUnits = [32, 32, 32, 32], \n",
    "                    weightRegularization = 0.0, \n",
    "                    dropout_rate=0.1):\n",
    "    \n",
    "    '''\n",
    "    Create a feed forward network.  Assumes Xtrain & Ytrain have been created and scaled\n",
    "    \n",
    "    Params: \n",
    "    optimizer (str): choice of optimizer\n",
    "    numUnits (list): number of units in each hidden\n",
    "    weightRegularization (float): between 0 and 1\n",
    "    dropout_rate (float): between 0 and 1\n",
    "    \n",
    "    '''\n",
    "    K.clear_session()\n",
    "    inputs = Input(shape=(Xtrain_scaled.shape[1],))    \n",
    "    \n",
    "    # add layers\n",
    "    for ii in np.arange(0, len(numUnits)):\n",
    "        if ii >= 1: \n",
    "            x = Dense(numUnits[ii], activation='tanh', \n",
    "                      kernel_regularizer=regularizers.l1(weightRegularization))(x)\n",
    "\n",
    "        else: \n",
    "            x = Dense(numUnits[ii], activation='tanh')(inputs)\n",
    "\n",
    "\n",
    "        # add dropout\n",
    "        if dropout_rate > 0: \n",
    "            x = Dropout(dropout_rate)(x)\n",
    "    \n",
    "    predictions = Dense(Ytrain_scaled.shape[1], activation='linear')(x)\n",
    "\n",
    "    # create model\n",
    "    model = Model(inputs=inputs, outputs=predictions)\n",
    "    model.compile(loss=\"mean_squared_error\", optimizer = optimizer, metrics = ['mse'])\n",
    "\n",
    "    return(model)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "modelParams = {\"optimizer\": \"rmsprop\", \n",
    "              \"dropout_rate\" : 0.0, \n",
    "               \"numUnits\": [32, 32, 32, 32],\n",
    "               \"weightRegularization\": 0.0\n",
    "              }\n",
    "\n",
    "\n",
    "model = create_network(**modelParams)\n",
    "\n",
    "modelName = ''.join('{}_{}__'.format(key[0:3].capitalize(), val) for  key, val in modelParams.items()).replace(\"[\", \"\").replace(\"]\", \"\").replace(\", \", \"_\")[0:-2]\n",
    "print(modelName)\n",
    "\n",
    "model.summary()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#model.get_config()\n",
    "earlystop = EarlyStopping(monitor='val_mean_squared_error', patience=15, \n",
    "                          verbose=1, mode='auto', min_delta = 0.000001)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# # fit model without regularization\n",
    "stt = time.time()\n",
    "history = model.fit(Xtrain_scaled, Ytrain_scaled, epochs = 2, verbose = 2, \n",
    "                        batch_size=2**14, callbacks = [earlystop], validation_split = 0.3)\n",
    "winsound.PlaySound(\"*\", winsound.SND_ALIAS)\n",
    "endd = time.time() - stt\n",
    "print(endd)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "history.history['val_loss'][-1]# MSE for final epoch testing\n",
    "# why are val MSE vs. loss different????  answer : I think loss incorporates regularization. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "history.history['val_mean_squared_error'][-1] # it is the same when there is no regularization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "history.history['loss'][-1]# MSE for final epoch training"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def format_e(n):\n",
    "    a = '%E' % n\n",
    "    return a.split('E')[0].rstrip('0').rstrip('.') + 'E' + a.split('E')[1]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# get nonzero weights\n",
    "wts = model.get_weights().copy()\n",
    "nzwts = np.sum([np.nonzero(wts[ii].reshape(-1))[0].shape[0] for ii in range(len(wts))])\n",
    "nzwts\n",
    "\n",
    "def plot_model_history(model_history, saveFig = False):\n",
    "    fig, axs = plt.subplots(1,1,figsize=(10,5))\n",
    "    # summarize history for accuracy\n",
    "    axs.plot(range(1,len(model_history.history['mean_squared_error'])+1),\n",
    "             model_history.history['mean_squared_error'])\n",
    "    axs.plot(range(1,len(model_history.history['val_mean_squared_error'])+1),\n",
    "             model_history.history['val_mean_squared_error'])\n",
    "    axs.set_title('Model MSE = '+ str(format_e(model_history.history['val_mean_squared_error'][-1])) + \"\\n\" +  str(nzwts) + \" non-zero weights\")\n",
    "    axs.set_ylabel('mean_squared_error')\n",
    "    axs.set_xlabel('Epoch')\n",
    "    axs.set_xticks(np.arange(1,len(model_history.history['val_mean_squared_error'])+1),\n",
    "                   len(model_history.history['val_mean_squared_error'])/10)\n",
    "    axs.legend(['train', 'val'], loc='best')\n",
    "    plt.yscale('log') #logarithmic scale for y axis\n",
    "    if saveFig:\n",
    "        fig.savefig(os.path.join(figDir, \"ModelTraining_\" + modelName + \".png\"), dpi = 120, bbox_inches='tight')\n",
    "    \n",
    "    plt.show()\n",
    "    \n",
    "plot_model_history(history, saveFig=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# save model\n",
    "model.save(os.path.join(savedModels,  modelName + '.h5'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# save weights\n",
    "wts = model.get_weights().copy()\n",
    "\n",
    "wtsFile = modelName + '_wts.pkl'\n",
    "pickle.dump(wts, open(os.path.join(dataOutput, wtsFile), 'wb'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# read file back in \n",
    "wt2 = pickle.load(open(os.path.join(dataOutput, wtsFile), 'rb'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# END\n",
    "os.path.join(dataOutput, wtsFile)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "___"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## START NEW ITEM: train and trim weights"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# train and trim weights"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "modelParams = {\"optimizer\": \"rmsprop\", \n",
    "              \"dropout_rate\" : 0.0, \n",
    "               \"numUnits\": [32, 32, 32, 32],\n",
    "               \"weightRegularization\": 0\n",
    "              }\n",
    "\n",
    "\n",
    "model = create_network(**modelParams)\n",
    "wts = model.get_weights().copy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.close(\"all\")\n",
    "fig, axs = plt.subplots(2,5, figsize=(30, 10), facecolor='w', edgecolor='k', )\n",
    "fig.subplots_adjust(hspace = 0.1, wspace=0.1)\n",
    "\n",
    "\n",
    "axs = axs.ravel(\"F\")\n",
    "\n",
    "for jj, ii in enumerate(np.arange(0, len(wts), 1)):\n",
    "    \n",
    "    \n",
    "    \n",
    "    if len(wts[ii].shape) < 2: \n",
    "        wtsTmp = wts[ii].reshape(1,-1).copy()\n",
    "    else:\n",
    "        wtsTmp = wts[ii].copy()\n",
    "    \n",
    "    im = axs[jj].matshow(wtsTmp, cmap = \"PRGn\", norm=colors.SymLogNorm(linthresh=0.01, linscale=0.01,\n",
    "                                              vmin=-3.0, vmax=3.0))\n",
    "    if np.mod(ii+1, 2) == 0:\n",
    "        axs[jj].axes.get_yaxis().set_visible(False)\n",
    "\n",
    "plt.colorbar(im,ax=axs[jj], orientation = \"horizontal\")\n",
    "plt.savefig(os.path.join(figDir, \"RandomWeightMatrices.png\"))\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from keras.models import load_model\n",
    "model = load_model(\"D:\\Dropbox\\AcademiaDropbox\\mothMachineLearning_dataAndFigs\\savedModels\\Opt_rmsprop__Dro_0.0__Num_32_32_32_32__Wei_0_pruned_bias.h5\")\n",
    "\n",
    "modelName = ''.join('{}_{}__'.format(key[0:3].capitalize(), val) for  key, val in modelParams.items()).replace(\"[\", \"\").replace(\"]\", \"\").replace(\", \", \"_\")[0:-2]\n",
    "print(modelName)\n",
    "\n",
    "model.summary()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "wts = model.get_weights().copy()\n",
    "\n",
    "wtsFile = modelName + '_wts.pkl'\n",
    "pickle.dump(wts, open(os.path.join(dataOutput, wtsFile), 'wb'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# read in data\n",
    "wts =  pickle.load(open(os.path.join(dataOutput, wtsFile), 'rb'))\n",
    "\n",
    "# print sizes of each weight matrix\n",
    "wtLengths = []\n",
    "for ii in range(len(wts)):\n",
    "    print(wts[ii].shape)\n",
    "    wtLengths.append(np.prod(wts[ii].shape))\n",
    "\n",
    "print(np.sum(wtLengths), \"total weights\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# start training\n",
    "historyDict = {\"mean_squared_error\": [], \n",
    "               \"val_mean_squared_error\": []}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "numCuts = 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# refref: \n",
    "# train until I have very good results\n",
    "# then prune, and retrain until results are close\n",
    "# change pruning rate as data get closer to 100% removed\n",
    "def prune_percent_updater(x):\n",
    "    logit = np.exp(x*8) / (np.exp(x*8) + 1)\n",
    "    return((logit - 0.5)*2*50)\n",
    "\n",
    "\n",
    "# cuts a smaller portion as the percent gets closer to 100%\n",
    "cutPercent = prune_percent_updater(np.linspace(0, 1, 25))\n",
    "\n",
    "while True:   \n",
    "   \n",
    "    for numEpocs in range(100):\n",
    "        \n",
    "        MSE_tmp = []\n",
    "\n",
    "        history = model.fit(Xtrain_scaled, Ytrain_scaled, validation_split = 0.3, \n",
    "                            verbose = 2, batch_size=2**14, epochs = 1, \n",
    "                            callbacks = [earlystop])\n",
    "        \n",
    "        # refref: earlystop doesn't do anything here\n",
    "        \n",
    "        # save history\n",
    "        historyDict[\"mean_squared_error\"].append(history.history[\"mean_squared_error\"][0])\n",
    "        historyDict[\"val_mean_squared_error\"].append(history.history[\"val_mean_squared_error\"][0])\n",
    "        \n",
    "        # local MSE\n",
    "        MSE_tmp.append(history.history[\"mean_squared_error\"][0])\n",
    "\n",
    "        # set weights that are close to 0 all the way back to 0, and then retrain for one epoch\n",
    "        # get nonzero weights\n",
    "        wts = model.get_weights().copy()\n",
    "\n",
    "        # set weights close to 0 to 0 (but ignore biases)\n",
    "        for ii in np.arange(0, len(wts), 1):\n",
    "            qants = np.percentile(np.reshape(wts[ii], -1), q = (50 - cutPercent[numCuts], 50 + cutPercent[numCuts]), )\n",
    "            wts[ii][(wts[ii] > qants[0]) & (wts[ii] < qants[1])] = 0\n",
    "\n",
    "        # print nonzero weights\n",
    "        # calculate number of nonzero weights\n",
    "        nzwts = np.sum([np.nonzero(wts[ii].reshape(-1))[0].shape[0] for ii in range(len(wts))])\n",
    "        print(nzwts, \"of\", np.sum(wtLengths), \"weights retained\")\n",
    "\n",
    "        # set new weights and calculate new loss\n",
    "        model.set_weights(wts)\n",
    "        \n",
    "        # check the change in mean squared error, and if it's not changing much, then cut out more data\n",
    "        # calculate slope of loss, based on previous 5 data points\n",
    "        if numEpocs > 5:\n",
    "            inputData = historyDict[\"mean_squared_error\"][-5:]\n",
    "\n",
    "            m = np.shape(inputData)\n",
    "            X = np.matrix([np.ones(m), np.arange(0, len(inputData))]).T\n",
    "            y = np.matrix(np.log(inputData)).T\n",
    "\n",
    "            # Solve for projection matrix\n",
    "            intercept, slope = np.array(np.linalg.inv(X.T.dot(X)).dot(X.T).dot(y)).reshape(-1,)\n",
    "            print(\"change in log loss:\", slope)\n",
    "    \n",
    "            # break if slope has stopped changing or if the overall min has been surpassed\n",
    "            # in the first training, it will automatically prune after 5 epochs, because the min will be passed\n",
    "            if (np.abs(slope) < 0.001) or (history.history[\"mean_squared_error\"][0] < np.min(historyDict[\"mean_squared_error\"][:-1])): \n",
    "                print(\"************************************************ PRUNING ********************************************************\")\n",
    "                print(\"************************************************ PRUNING ********************************************************\")\n",
    "                print(\"************************************************ PRUNING ********************************************************\")\n",
    "                break\n",
    "                       \n",
    "                    \n",
    "    ## refref: may want to save weights before each pruning, so I can go back, if I need to\n",
    "    ## refref: should I be pruning the biases too?\n",
    "    \n",
    "    ## keep running tally of min mse, and if we can't get back to the min, then break\n",
    "    print(\"Min MSE for this prune \", np.min(MSE_tmp), \"______overall Min MSE \", np.min(historyDict[\"mean_squared_error\"]))\n",
    "    if np.min(MSE_tmp) > np.min(historyDict[\"mean_squared_error\"]):\n",
    "        print(\"no more gain by pruning:  STOPPING Pruning\")\n",
    "        break\n",
    "    \n",
    "    numCuts += 1\n",
    "    if numCuts > len(cutPercent):\n",
    "        break\n",
    "\n",
    "        \n",
    "        #cutPercent += 0.2\n",
    "winsound.PlaySound(\"*\", winsound.SND_ALIAS)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(numCuts)\n",
    "(50 - cutPercent[numCuts]) * 2 # percent of original network size that is used the pruned version"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    " numCuts = 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_model_history_fromDict(model_history_dictionary, saveFig = False):\n",
    "    fig, axs = plt.subplots(1,1,figsize=(10,5))\n",
    "    # summarize history for accuracy\n",
    "    axs.plot(range(1,len(model_history_dictionary['mean_squared_error'])+1),\n",
    "             model_history_dictionary['mean_squared_error'])\n",
    "    axs.plot(range(1,len(model_history_dictionary['val_mean_squared_error'])+1),\n",
    "             model_history_dictionary['val_mean_squared_error'])\n",
    "    axs.set_title('Model MSE = '+ str(format_e(model_history_dictionary['val_mean_squared_error'][-1])) + \"\\n\" +  str(nzwts) + \" non-zero weights\")\n",
    "    axs.set_ylabel('mean_squared_error')\n",
    "    axs.set_xlabel('Epoch')\n",
    "    axs.set_xticks(np.arange(1,len(model_history_dictionary['val_mean_squared_error'])+1),\n",
    "                   len(model_history_dictionary['val_mean_squared_error'])/10)\n",
    "    axs.legend(['train', 'val'], loc='best')\n",
    "    plt.yscale('log') #logarithmic scale for y axis\n",
    "    if saveFig:\n",
    "        fig.savefig(os.path.join(figDir, \"ModelTraining_\" + modelName + \"_pruned_2.png\"), dpi = 120, bbox_inches='tight')\n",
    "        print(os.path.join(figDir, \"ModelTraining_\" + modelName + \"_pruned.png\"))\n",
    "    plt.show()\n",
    "    \n",
    "plot_model_history_fromDict(historyDict, saveFig = True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# how good can I get without trimming\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot matrices\n",
    "import matplotlib.colors as colors\n",
    "wts = model.get_weights().copy()\n",
    "for ii in np.arange(0, len(wts), 2):\n",
    "    plt.matshow(wts[ii], cmap = \"PRGn\", norm=colors.SymLogNorm(linthresh=0.01, linscale=0.01,\n",
    "                                              vmin=-3.0, vmax=3.0))\n",
    "    plt.colorbar()\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot biases\n",
    "\n",
    "for ii in np.arange(1, len(wts), 2):\n",
    "    plt.matshow(wts[ii].reshape(1, -1), cmap = \"PRGn\", norm=colors.SymLogNorm(linthresh=0.01, linscale=0.01,\n",
    "                                              vmin=-3.0, vmax=3.0))\n",
    "    plt.colorbar()\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# refref: save weights and model\n",
    "model.save(os.path.join(savedModels,  modelName + '_pruned_bias.h5'))\n",
    "\n",
    "# save weights\n",
    "wts = model.get_weights().copy()\n",
    "\n",
    "wtsFile = modelName + '_pruned_wts_bias.pkl'\n",
    "pickle.dump(wts, open(os.path.join(dataOutput, wtsFile), 'wb'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.close(\"all\")\n",
    "fig, axs = plt.subplots(1,5, figsize=(20, 5), facecolor='w', edgecolor='k')\n",
    "fig.subplots_adjust(hspace = 0.1, wspace=0.2)\n",
    "\n",
    "\n",
    "axs = axs.ravel()\n",
    "\n",
    "for jj, ii in enumerate(np.arange(0, len(wts), 2)):\n",
    "    im = axs[jj].matshow(wts[ii], cmap = \"PRGn\", norm=colors.SymLogNorm(linthresh=0.01, linscale=0.01,\n",
    "                                              vmin=-3.0, vmax=3.0))\n",
    "    \n",
    "    \n",
    "plt.colorbar(im,ax=axs[jj], orientation = \"vertical\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.close(\"all\")\n",
    "fig, axs = plt.subplots(5,1, figsize=(30, 10), facecolor='w', edgecolor='k')\n",
    "fig.subplots_adjust(hspace = 0.1, wspace=0.5)\n",
    "\n",
    "\n",
    "axs = axs.ravel()\n",
    "\n",
    "for jj, ii in enumerate(np.arange(1, len(wts), 2)):\n",
    "    im = axs[jj].matshow(wts[ii].reshape(1, -1), cmap = \"PRGn\", norm=colors.SymLogNorm(linthresh=0.01, linscale=0.01,\n",
    "                                              vmin=-3.0, vmax=3.0))\n",
    "    axs[jj].axes.get_yaxis().set_visible(False)\n",
    "\n",
    "plt.colorbar(im,ax=axs[jj], orientation = \"horizontal\")\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.close(\"all\")\n",
    "fig, axs = plt.subplots(2,5, figsize=(30, 10), facecolor='w', edgecolor='k', )\n",
    "fig.subplots_adjust(hspace = 0.1, wspace=0.1)\n",
    "\n",
    "\n",
    "axs = axs.ravel(\"F\")\n",
    "\n",
    "for jj, ii in enumerate(np.arange(0, len(wts), 1)):\n",
    "    \n",
    "    \n",
    "    \n",
    "    if len(wts[ii].shape) < 2: \n",
    "        wtsTmp = wts[ii].reshape(1,-1).copy()\n",
    "    else:\n",
    "        wtsTmp = wts[ii].copy()\n",
    "    \n",
    "    im = axs[jj].matshow(wtsTmp, cmap = \"PRGn\", norm=colors.SymLogNorm(linthresh=0.01, linscale=0.01,\n",
    "                                              vmin=-3.0, vmax=3.0))\n",
    "    if np.mod(ii+1, 2) == 0:\n",
    "        axs[jj].axes.get_yaxis().set_visible(False)\n",
    "\n",
    "plt.colorbar(im,ax=axs[jj], orientation = \"horizontal\")\n",
    "plt.savefig(os.path.join(figDir, \"PrunedWeightMatrices.png\"))\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "os.path.join(dataOutput, wtsFile)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# refref: if whole node is basically 0, then remove the node. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "np.std(wts[2].reshape(-1,))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model = create_network(**modelParams)\n",
    "\n",
    "history = model.fit(Xtrain_scaled, Ytrain_scaled, validation_split = 0.3, \n",
    "                        verbose = 2, batch_size=2**14, epochs = 1, \n",
    "                        callbacks = [earlystop])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# check to see if model saved: \n",
    "K.clear_session()\n",
    "from keras.models import load_model\n",
    "model = load_model(os.path.join(savedModels, 'my_model_400Units_newData.h5'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# get nonzero weights\n",
    "wts = model.get_weights().copy()\n",
    "nzwts = np.sum([np.nonzero(wts[ii].reshape(-1))[0].shape[0] for ii in range(len(wts))])\n",
    "nzwts"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "len(wts[4])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# show weights histograms\n",
    "plt.close(\"all\")\n",
    "fig, axs = plt.subplots(2,8, figsize=np.array((20,3)) , facecolor='w', edgecolor='k')\n",
    "fig.subplots_adjust(hspace = 0.3, wspace=0.3)\n",
    "axs = axs.ravel(order = \"F\")\n",
    "\n",
    "for ii in range(int(len(wts) / 2)):  \n",
    "    jj= int(2*ii)\n",
    "    axs[jj].hist(wts[jj].reshape(-1), bins = 100)\n",
    "    axs[jj].set_xlabel(\"Layer \" + str(int(jj/2)) + \" weights\" + \", shape = \" + str(wts[jj].shape))\n",
    "    axs[jj+1].hist(wts[jj+1], bins = 100)\n",
    "    axs[jj+1].set_xlabel(\"Layer \" + str(int(jj/2)) + \" biases\" + \", shape = \" + str(wts[jj+1].shape))\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "len(model.get_weights())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model.evaluate(Xtest_scaled, Ytest_scaled)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# get nonzero weights\n",
    "wts = model.get_weights().copy()\n",
    "nzwts = np.sum([np.nonzero(wts[ii].reshape(-1))[0].shape[0] for ii in range(len(wts))])\n",
    "nzwts"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# trim weights\n",
    "for ii in np.arange(0, 7):\n",
    "    qants = np.percentile(np.reshape(wts[ii], -1), q = (5, 95), )\n",
    "    wts[ii][(wts[ii] > qants[0]) & (wts[ii] < qants[1])] = 0\n",
    "    \n",
    "# calculate number of nonzero weights\n",
    "nzwts = np.sum([np.nonzero(wts[ii].reshape(-1))[0].shape[0] for ii in range(len(wts))])\n",
    "nzwts"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set new weights and calculate new loss\n",
    "model.set_weights(wts)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ValLoss = model.evaluate(Xtrain_scaled[-int(0.05*Xtrain_scaled.shape[0]):, :], \n",
    "                         Ytrain_scaled[-int(0.05*Xtrain_scaled.shape[0]):, :])\n",
    "print(ValLoss[0])"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "K.clear_session()"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "# this is the original model\n",
    "inputs = Input(shape=(Xtrain_scaled.shape[1],))\n",
    "x = Dense(400, activation='tanh')(inputs)\n",
    "x = Dense(400, activation='tanh')(x)\n",
    "x = Dense(400, activation='tanh')(x)\n",
    "x = Dense(16, activation='tanh')(x)\n",
    "predictions = Dense(Ytrain_scaled.shape[1], activation='linear')(x)\n",
    "\n",
    "# create model\n",
    "model = Model(inputs=inputs, outputs=predictions)\n",
    "model.compile(loss=\"mse\", optimizer=\"adam\", metrics = ['mse'])\n",
    "\n",
    "earlystop = EarlyStopping(monitor='val_mean_squared_error', patience=50, \n",
    "                          verbose=1, mode='auto', min_delta = 0.001)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "wts = model.get_weights().copy()\n",
    "nzwts = np.sum([np.nonzero(wts[ii].reshape(-1))[0].shape[0] for ii in range(len(wts))])\n",
    "print(nzwts)\n",
    "\n",
    "# trim weights\n",
    "for ii in np.arange(0, 7):\n",
    "    qants = np.percentile(np.reshape(wts[ii], -1), q = (2, 98), )\n",
    "    wts[ii][(wts[ii] > qants[0]) & (wts[ii] < qants[1])] = 0\n",
    "    \n",
    "# calculate number of nonzero weights\n",
    "nzwts = np.sum([np.nonzero(wts[ii].reshape(-1))[0].shape[0] for ii in range(len(wts))])\n",
    "print(nzwts)\n",
    "model.set_weights(wts)\n",
    "\n",
    "\n",
    "# start training\n",
    "history = model.fit(Xtrain_scaled, Ytrain_scaled, validation_split = 0.3, \n",
    "                    verbose = 2, batch_size=2**14, epochs = 1, \n",
    "                    callbacks = [earlystop])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "wts = model.get_weights().copy()\n",
    "\n",
    "# trim weights\n",
    "for ii in np.arange(0, 7):\n",
    "    qants = np.percentile(np.reshape(wts[ii], -1), q = (2.5, 97.5), )\n",
    "    wts[ii][(wts[ii] > qants[0]) & (wts[ii] < qants[1])] = 0\n",
    "\n",
    "# calculate number of nonzero weights\n",
    "nzwts = np.sum([np.nonzero(wts[ii].reshape(-1))[0].shape[0] for ii in range(len(wts))])\n",
    "print(nzwts)\n",
    "model.set_weights(wts)\n",
    "\n",
    "model.evaluate(Xtest_scaled, Ytest_scaled)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_model_history(model_history, saveFig = False):\n",
    "    fig, axs = plt.subplots(1,1,figsize=(10,5))\n",
    "    # summarize history for accuracy\n",
    "    axs.plot(range(1,len(model_history.history['mean_squared_error'])+1),\n",
    "             model_history.history['mean_squared_error'])\n",
    "    axs.plot(range(1,len(model_history.history['val_mean_squared_error'])+1),\n",
    "             model_history.history['val_mean_squared_error'])\n",
    "    axs.set_title('Model MSE')\n",
    "    axs.set_ylabel('mean_squared_error')\n",
    "    axs.set_xlabel('Epoch')\n",
    "    axs.set_xticks(np.arange(1,len(model_history.history['val_mean_squared_error'])+1),\n",
    "                   len(model_history.history['val_mean_squared_error'])/10)\n",
    "    axs.legend(['train', 'val'], loc='best')\n",
    "    if saveFig:\n",
    "        fig.savefig(os.path.join(figDir, \"ModelTraining.png\"), dpi = 120, bbox_inches='tight')\n",
    "    plt.show()\n",
    "    \n",
    "plot_model_history(history)\n",
    "print(history.history[\"loss\"][-1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# load model that was trained for much longer\n",
    "from keras.models import load_model\n",
    "model = load_model(os.path.join(savedModels, 'my_model.h5'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model.summary()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# get nonzero weights\n",
    "wts = model.get_weights().copy()\n",
    "nzwts = np.sum([np.nonzero(wts[ii].reshape(-1))[0].shape[0] for ii in range(len(wts))])\n",
    "nzwts"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set new weights and calculate new loss\n",
    "model.set_weights(wts)\n",
    "\n",
    "ValLoss = model.evaluate(Xtrain_scaled[-int(0.05*Xtrain_scaled.shape[0]):, :], \n",
    "                         Ytrain_scaled[-int(0.05*Xtrain_scaled.shape[0]):, :])\n",
    "print(ValLoss[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "wts = model.get_weights().copy()\n",
    "nzwts = [np.nonzero(wts[ii].reshape(-1))[0] for ii in range(len(wts))]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# show weights histograms\n",
    "plt.close(\"all\")\n",
    "fig, axs = plt.subplots(2,8, figsize=np.array((20,5)) , facecolor='w', edgecolor='k')\n",
    "#fig.subplots_adjust(hspace = 0.3, wspace=0.3)\n",
    "axs = axs.ravel(order = \"F\")\n",
    "\n",
    "for ii in range(int(len(wts) / 2)):  \n",
    "    jj= int(2*ii)\n",
    "    axs[jj].hist(nzwts[jj].reshape(-1), bins = 30)\n",
    "    axs[jj].set_xlabel(\"Layer \" + str(int(jj/2)) + \" weights\" + \", shape = \" + str(nzwts[jj].shape))\n",
    "    axs[jj+1].hist(nzwts[jj+1], bins = 30)\n",
    "    axs[jj+1].set_xlabel(\"Layer \" + str(int(jj/2)) + \" biases\" + \", shape = \" + str(nzwts[jj+1].shape))\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# show weights histograms\n",
    "plt.close(\"all\")\n",
    "fig, axs = plt.subplots(2,8, figsize=np.array((10,10)) , facecolor='w', edgecolor='k')\n",
    "#fig.subplots_adjust(hspace = 0.3, wspace=0.3)\n",
    "axs = axs.ravel(order = \"F\")\n",
    "\n",
    "for ii in range(int(len(wts) / 2)):  \n",
    "    jj= int(2*ii)\n",
    "    axs[jj].hist(wts[jj].reshape(-1), bins = 30)\n",
    "    axs[jj].set_xlabel(\"Layer \" + str(int(jj/2)) + \" weights\" + \", shape = \" + str(wts[jj].shape))\n",
    "    axs[jj+1].hist(wts[jj+1], bins = 30)\n",
    "    axs[jj+1].set_xlabel(\"Layer \" + str(int(jj/2)) + \" biases\" + \", shape = \" + str(wts[jj+1].shape))\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# predict on test set\n",
    "nnpreds = model.predict(Xtest_scaled[ :])\n",
    "\n",
    "# rescale\n",
    "nnpreds_unscaled = scalerY.inverse_transform(nnpreds)\n",
    "\n",
    "# show residuals\n",
    "# combine residual and regular plots\n",
    "plt.close(\"all\")\n",
    "fig, axs = plt.subplots(2,7, figsize=np.array((30, 8)) / 1.7, facecolor='w', edgecolor='k')\n",
    "fig.subplots_adjust(hspace = 0.1, wspace=0.5)\n",
    "axs = axs.ravel()\n",
    "\n",
    "# replace lightest colors with white\n",
    "import matplotlib.colors\n",
    "cmap = plt.cm.magma_r\n",
    "cmaplist = np.array([cmap(i) for i in range(cmap.N)])\n",
    "cmaplist[:,0:3] = np.divide(cmaplist[:, 0:3], 1.1)\n",
    "cmaplist[0] = (1,1,1,0.5)\n",
    "cmap = matplotlib.colors.LinearSegmentedColormap.from_list('mcm',cmaplist, cmap.N)\n",
    "\n",
    "for ii in range(len(Y.columns)):\n",
    "    \n",
    "    try:\n",
    "        axs[ii].hexbin(y = Ytest.iloc[:,ii],x = nnpreds_unscaled[:,ii], gridsize = 150, cmap = cmap)\n",
    "        axs[ii].yaxis.set_major_formatter(ScalarFormatter(useMathText=True))\n",
    "        axs[ii].ticklabel_format(style='sci',  axis='y', scilimits=(3,4))\n",
    "        axs[ii].axes.xaxis.set_ticklabels([])\n",
    "        if(ii == 0):\n",
    "            axs[ii].set_ylabel(\"Actual Value\")\n",
    "        axs[ii].set_title(nms2[ii])\n",
    "        axs[ii].plot(Ytest.iloc[0:1000,ii], Ytest.iloc[0:1000,ii], 'grey', linewidth = 1, linestyle  = \"--\")\n",
    "    except:\n",
    "        pass\n",
    "\n",
    "for ii in range(len(Y.columns)):\n",
    "    jj = ii + len(Y.columns)\n",
    "    \n",
    "    try:\n",
    "        axs[jj].hexbin(y = Ytest.iloc[:,ii] - nnpreds_unscaled[:,ii],x = nnpreds_unscaled[:,ii], gridsize = 150, cmap = cmap)\n",
    "        axs[jj].set_xlabel(\"Predicted Value\")\n",
    "        axs[jj].xaxis.set_major_formatter(ScalarFormatter(useMathText=True))\n",
    "        axs[jj].yaxis.set_major_formatter(ScalarFormatter(useMathText=True))\n",
    "        mmin = np.min(Ytest.iloc[:,ii] - nnpreds_unscaled[:,ii])\n",
    "        mmax = np.max(Ytest.iloc[:,ii] - nnpreds_unscaled[:,ii])\n",
    "        \n",
    "        upper = np.max([np.abs(mmin), np.abs(mmax)])\n",
    "        axs[jj].set_ylim(-upper, upper)\n",
    "\n",
    "        if(ii == 0):\n",
    "            axs[jj].set_ylabel(\"Actual - Predicted\")\n",
    "        axs[jj].hlines(y = 0, xmin = np.min(nnpreds_unscaled[:,ii]), \n",
    "                       xmax = np.max(nnpreds_unscaled[:,ii]), linestyle =  \"--\", linewidth = 1)\n",
    "    except:\n",
    "        pass\n",
    "plt.tight_layout()\n",
    "#fig.savefig(os.path.join(figDir, \"SmallModelResids.png\"), dpi = 120, bbox_inches='tight')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# trim distribution of weights -- cut out middle 20%\n",
    "for ii in np.arange(0, 7):\n",
    "    qants = np.percentile(np.reshape(wts[ii], -1), q = (40, 60), )\n",
    "    wts[ii][(wts[ii] > qants[0]) & (wts[ii] < qants[1])] = 0\n",
    "    \n",
    "# calculate number of nonzero weights\n",
    "nzwts = np.sum([np.nonzero(wts[ii].reshape(-1))[0].shape[0] for ii in range(len(wts))])\n",
    "nzwts"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# show new histogram of weights (excluding the 0's)\n",
    "# show weights histograms\n",
    "plt.close(\"all\")\n",
    "fig, axs = plt.subplots(2,5, figsize=np.array((15, 6)) , facecolor='w', edgecolor='k')\n",
    "fig.subplots_adjust(hspace = 0.3, wspace=0.3)\n",
    "axs = axs.ravel(order = \"F\")\n",
    "\n",
    "for ii in range(int(len(wts) / 2)):  \n",
    "    jj= int(2*ii)\n",
    "    \n",
    "    d1 = wts[jj].reshape(-1)\n",
    "    axs[jj].hist(d1[d1!=0], bins = 30, facecolor = '#d6bddb' )\n",
    "    axs[jj].set_xlabel(\"Layer \" + str(int(jj/2)) + \" weights\" + \", shape = \" + str(wts[jj].shape))\n",
    "\n",
    "    d2 = wts[jj+1]\n",
    "    axs[jj+1].hist(d2[d2!=0], bins = 30, facecolor = '#d6bddb')\n",
    "    axs[jj+1].set_xlabel(\"Layer \" + str(int(jj/2)) + \" biases\" + \", shape = \" + str(wts[jj+1].shape))\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# the validation.split is the last X% of the data\n",
    "int(0.3*Xtrain_scaled.shape[0])\n",
    "ValLoss = model.evaluate(Xtrain_scaled[-int(0.3*Xtrain_scaled.shape[0]):, :], \n",
    "                         Ytrain_scaled[-int(0.3*Xtrain_scaled.shape[0]):, :])\n",
    "print(ValLoss[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set new weights and calculate new loss\n",
    "model.set_weights(wts)\n",
    "\n",
    "ValLoss = model.evaluate(Xtrain_scaled[-int(0.3*Xtrain_scaled.shape[0]):, :], \n",
    "                         Ytrain_scaled[-int(0.3*Xtrain_scaled.shape[0]):, :])\n",
    "print(ValLoss[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# predict on test set\n",
    "nnpreds = model.predict(Xtest_scaled[ :])\n",
    "\n",
    "# rescale\n",
    "nnpreds_unscaled = scalerY.inverse_transform(nnpreds)\n",
    "\n",
    "# show residuals\n",
    "# combine residual and regular plots\n",
    "plt.close(\"all\")\n",
    "fig, axs = plt.subplots(2,7, figsize=np.array((30, 8)) / 1.7, facecolor='w', edgecolor='k')\n",
    "fig.subplots_adjust(hspace = 0.1, wspace=0.5)\n",
    "axs = axs.ravel()\n",
    "\n",
    "# replace lightest colors with white\n",
    "import matplotlib.colors\n",
    "cmap = plt.cm.magma_r\n",
    "cmaplist = np.array([cmap(i) for i in range(cmap.N)])\n",
    "cmaplist[:,0:3] = np.divide(cmaplist[:, 0:3], 1.1)\n",
    "cmaplist[0] = (1,1,1,0.5)\n",
    "cmap = matplotlib.colors.LinearSegmentedColormap.from_list('mcm',cmaplist, cmap.N)\n",
    "\n",
    "for ii in range(len(Y.columns)):\n",
    "    \n",
    "    try:\n",
    "        axs[ii].hexbin(y = Ytest.iloc[:,ii],x = nnpreds_unscaled[:,ii], gridsize = 150, cmap = cmap)\n",
    "        axs[ii].yaxis.set_major_formatter(ScalarFormatter(useMathText=True))\n",
    "        axs[ii].ticklabel_format(style='sci',  axis='y', scilimits=(3,4))\n",
    "        axs[ii].axes.xaxis.set_ticklabels([])\n",
    "        if(ii == 0):\n",
    "            axs[ii].set_ylabel(\"Actual Value\")\n",
    "        axs[ii].set_title(nms2[ii])\n",
    "        axs[ii].plot(Ytest.iloc[0:1000,ii], Ytest.iloc[0:1000,ii], 'grey', linewidth = 1, linestyle  = \"--\")\n",
    "    except:\n",
    "        pass\n",
    "\n",
    "for ii in range(len(Y.columns)):\n",
    "    jj = ii + len(Y.columns)\n",
    "    \n",
    "    try:\n",
    "        axs[jj].hexbin(y = Ytest.iloc[:,ii] - nnpreds_unscaled[:,ii],x = nnpreds_unscaled[:,ii], gridsize = 150, cmap = cmap)\n",
    "        axs[jj].set_xlabel(\"Predicted Value\")\n",
    "        axs[jj].xaxis.set_major_formatter(ScalarFormatter(useMathText=True))\n",
    "        axs[jj].yaxis.set_major_formatter(ScalarFormatter(useMathText=True))\n",
    "        mmin = np.min(Ytest.iloc[:,ii] - nnpreds_unscaled[:,ii])\n",
    "        mmax = np.max(Ytest.iloc[:,ii] - nnpreds_unscaled[:,ii])\n",
    "        \n",
    "        upper = np.max([np.abs(mmin), np.abs(mmax)])\n",
    "        axs[jj].set_ylim(-upper, upper)\n",
    "\n",
    "        if(ii == 0):\n",
    "            axs[jj].set_ylabel(\"Actual - Predicted\")\n",
    "        axs[jj].hlines(y = 0, xmin = np.min(nnpreds_unscaled[:,ii]), \n",
    "                       xmax = np.max(nnpreds_unscaled[:,ii]), linestyle =  \"--\", linewidth = 1)\n",
    "    except:\n",
    "        pass\n",
    "plt.tight_layout()\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# create dict of hyperparameters\n",
    "\n",
    "\n",
    "# regularization, num layers, num nodes, learning rate, optimizer, activation function, batch size\n",
    "\n",
    "from scipy.stats import randint\n",
    "from scipy.stats import uniform\n",
    "\n",
    "# Create hyperparameter space\n",
    "NumHiddenLayers = randint(low = 2, high = 20)#[4, 2, 8]\n",
    "numUnits  = [2**4, 2**5, 2**6, 2**7, 2**8, 2**9, 2**10]\n",
    "epochs = [200]\n",
    "batches1 = [2**12, 2**10, 2**8, 2**14] \n",
    "optimizers = ['rmsprop', 'adam']\n",
    "dropout_rate =  uniform(loc = 0, scale = 0.5) #[0.0, 0.2, 0.5]\n",
    "weightRegularization = uniform(loc = 0, scale = 0.001) #[0, 0.0001, 0.001, 0.01]\n",
    "secondToLastUnits = [8, 16, 32, 64]\n",
    "\n",
    "# Create hyperparameter options\n",
    "hyperparameters = dict(optimizer=optimizers, \n",
    "                        epochs=epochs, \n",
    "                        batch_size=batches1,\n",
    "                        dropout_rate = dropout_rate, \n",
    "                        numUnits = numUnits, \n",
    "                        NumHiddenLayers = NumHiddenLayers, \n",
    "                        weightRegularization = weightRegularization, \n",
    "                        secondToLastUnits = secondToLastUnits)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# refref: \n",
    "# train until I have very good results\n",
    "# then prune, and retrain until results are close\n",
    "# change pruning rate as data get closer to 100% removed\n",
    "cutPercent = 49.7\n",
    "for numCuts in range(3):\n",
    "    for numEpocs in range(100):\n",
    "\n",
    "        history = model.fit(Xtrain_scaled, Ytrain_scaled, validation_split = 0.3, \n",
    "                            verbose = 2, batch_size=2**14, epochs = 1, \n",
    "                            callbacks = [earlystop])\n",
    "        # save history\n",
    "        historyDict[\"mean_squared_error\"].append(history.history[\"mean_squared_error\"][0])\n",
    "        historyDict[\"val_mean_squared_error\"].append(history.history[\"val_mean_squared_error\"][0])\n",
    "\n",
    "        # set weights that are close to 0 all the way back to 0, and then retrain for one epoch\n",
    "        # get nonzero weights\n",
    "        wts = model.get_weights().copy()\n",
    "\n",
    "        # set weights close to 0 to 0 (but ignore biases)\n",
    "        for ii in np.arange(0, len(wts), 2):\n",
    "            qants = np.percentile(np.reshape(wts[ii], -1), q = (50 - cutPercent, 50 + cutPercent), )\n",
    "            wts[ii][(wts[ii] > qants[0]) & (wts[ii] < qants[1])] = 0\n",
    "\n",
    "        # print nonzero weights\n",
    "        # calculate number of nonzero weights\n",
    "        nzwts = np.sum([np.nonzero(wts[ii].reshape(-1))[0].shape[0] for ii in range(len(wts))])\n",
    "        print(nzwts, \"of\", np.sum(wtLengths), \"weights retained\")\n",
    "\n",
    "        # set new weights and calculate new loss\n",
    "        model.set_weights(wts)\n",
    "        #cutPercent += 0.2\n",
    "winsound.PlaySound(\"*\", winsound.SND_ALIAS)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "deeplearning",
   "language": "python",
   "name": "deeplearning"
  },
  "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.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}