capstone.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "## General initialization\n",
    "\n",
    "import warnings\n",
    "with warnings.catch_warnings():\n",
    "    warnings.simplefilter('ignore')\n",
    "import numpy as np\n",
    "seed = 42\n",
    "if seed is not None:\n",
    "   np.random.seed(seed)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Read Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Size of train data:1604\n",
      "size of image in band_1: 5625\n",
      "Shape of image: 75x75\n",
      "In train data: percentage of iceberg: 46.94513715710723%\n"
     ]
    }
   ],
   "source": [
    "## Reading input training Data\n",
    "\n",
    "#configuration path Data\n",
    "path_train_data = \"./data/processed/train.json\"\n",
    "\n",
    "#Read Data\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "train_data = pd.read_json(path_train_data)\n",
    "print(\"Size of train data:{}\".format(len(train_data)))\n",
    "\n",
    "#Size of images\n",
    "size_image = len(train_data.iloc[0]['band_1'])\n",
    "print(\"size of image in band_1: {}\".format(size_image))\n",
    "width_image = int(np.sqrt(size_image))\n",
    "print(\"Shape of image: {}x{}\".format(width_image, width_image))\n",
    "\n",
    "#Percentage Ships/icebergs\n",
    "idx_iceberg_td = train_data[\"is_iceberg\"] == 1\n",
    "idx_ship_td = train_data[\"is_iceberg\"] == 0\n",
    "num_iceberg_in_train_data = train_data[idx_iceberg_td].shape[0]\n",
    "print(\"In train data: percentage of iceberg: {}%\".format(100 * num_iceberg_in_train_data / train_data.shape[0]))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "=> The training data set is quite balanced between the two classes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "## Reading input Test Competition Kaggle Data\n",
    "\n",
    "#configuration path Data\n",
    "path_test_data = \"./data/processed/test.json\"\n",
    "\n",
    "#Read Data\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "test_data = pd.read_json(path_test_data)\n",
    "print(\"size of test data:{}\".format(len(test_data)))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "=> The test dataset of the competition is much larger than the training set. \n",
    "It may be interesting to learn by unsupervised learning from this data set\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Let's have a look at the dirst datas\n",
    "from IPython.display import display # Allows the use of display() for DataFrames\n",
    "display(train_data.head())\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of missing value for inc angle in training: 133\n",
      "Mean value for missing value: 39.26870747790618, std=3.8397444116664374, median:39.5015\n",
      "Mean Squared error for predicting always mean value: 14.733614280066448\n",
      "Mean absolute  error for predicting always mean value: 3.288992522093814\n",
      "Median absolute  error for predicting always mean value: 3.288992522093814\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.figure.Figure at 0x116773da0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.figure.Figure at 0x11f4b7358>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "##################################################\n",
    "##  pre processing inc_angle NaN values\n",
    "\n",
    "#Some of the incident angle from the satellite are unknown and marked as \"na\". \n",
    "#Replace these na with 0 and find the indices where the incident angle is >0 \n",
    "#This way we can use a truncated set or the full set of training data.\n",
    "\n",
    "# Find indexes where there is no angle\n",
    "idx_inc_angle_nan_td = np.where(train_data[\"inc_angle\"] == 'na')\n",
    "idx_inc_angle_not_nan_td = np.where(train_data[\"inc_angle\"] != 'na')\n",
    "print(\"Number of missing value for inc angle in training: {}\".format(idx_inc_angle_nan_td[0].shape[0]))\n",
    "\n",
    "\n",
    "\n",
    "#transform in numeric the inc_angle in new column inc_angle_numeric\n",
    "train_data[\"inc_angle_numeric\"] = pd.to_numeric(train_data.inc_angle, errors='coerce')\n",
    "\n",
    "#replace the nan value for the mean in new columns inc_angle_corrected_mean\n",
    "train_data[\"inc_angle_corrected_mean\"] = train_data[\"inc_angle_numeric\"]\n",
    "mean_angle = train_data.inc_angle_corrected_mean.mean(skipna=True)\n",
    "std_angle = train_data.inc_angle_corrected_mean.std(skipna=True)\n",
    "median_angle = train_data.inc_angle_corrected_mean.median(skipna=True)\n",
    "train_data[\"inc_angle_corrected_mean\"].fillna(mean_angle, inplace=True)\n",
    "print(\"Mean value for missing value: {}, std={}, median:{}\".format(mean_angle, std_angle, median_angle))\n",
    "\n",
    "#Check Mean squared error if we had a regressor which always returns the mean\n",
    "from sklearn.metrics import mean_squared_error, mean_absolute_error, median_absolute_error\n",
    "pred = [mean_angle]*len(idx_inc_angle_not_nan_td[0])\n",
    "groundTruth = train_data['inc_angle_corrected_mean'][idx_inc_angle_not_nan_td[0]]\n",
    "mse = mean_squared_error(groundTruth, pred)\n",
    "\n",
    "meanae =  median_absolute_error(groundTruth, pred) \n",
    "#mae is more robust to outliers\n",
    "medianae =  median_absolute_error(groundTruth, pred) \n",
    "print(\"Mean Squared error for predicting always mean value: {}\".format(mse))\n",
    "print(\"Mean absolute  error for predicting always mean value: {}\".format(meanae))\n",
    "print(\"Median absolute  error for predicting always mean value: {}\".format(medianae))\n",
    "\n",
    "\n",
    "\n",
    "#for backward compatibility\n",
    "train_data['inc_angle_corrected'] = train_data['inc_angle_corrected_mean']\n",
    "\n",
    "# Plot the histogram of inc_angle fo original and corrected with mean\n",
    "import matplotlib.pyplot as plt\n",
    "plt.title(\"histogram inc_angle corected\")\n",
    "train_data.inc_angle_corrected_mean.hist()\n",
    "plt.show()\n",
    "\n",
    "plt.title(\"histogram inc_angle original\")\n",
    "train_data.inc_angle_numeric.hist()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "=> It seems that if  we replace the Nan value by the mean, we did not make a big mistake "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Visual data exploration"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Basic Ships and Icebergs Images"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Visualization of examples of Band_1 band_2 as images\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "num_show = 10\n",
    "def plot_images(idx,num,label):\n",
    "    fig = plt.figure(figsize=(50,50))\n",
    "    plt.suptitle('Data with '+label, fontsize=40)\n",
    "    for i in range(num_show):\n",
    "        plt.subplot(num_show,2,2*i + 1)\n",
    "        plt.title(\"band_1\")\n",
    "        img = np.array(train_data[idx].iloc[i]['band_1']).reshape(75,75)\n",
    "        plt.imshow(img)\n",
    "    \n",
    "        plt.subplot(num_show,2,2*i + 2)\n",
    "        plt.title(\"band_2\")\n",
    "        img = np.array(train_data[idx].iloc[i]['band_2']).reshape(75,75)\n",
    "        plt.imshow(img)\n",
    "    plt.show()\n",
    "    \n",
    "# visualization of the first 10 instances classified as iceberg\n",
    "plot_images(idx_iceberg_td,10,'Iceberg')\n",
    "\n",
    "# visualization of the first 10 instances classified as Ships\n",
    "plot_images(idx_ship_td,10,'Ships')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "=> It seems that the iceberg images are more blurry, and the ships images have more strong light pattern."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Color Composites visualization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "########################################\n",
    "## Definition of utils function to create images to tensors\n",
    "\n",
    "\n",
    "# Definition of pre processing function\n",
    "\n",
    "def no_process(t):\n",
    "    return t\n",
    "\n",
    "def min_max_scaler(t, scale=1):\n",
    "    tmin = t.min()\n",
    "    return  scale*(t - tmin) / (t.max() - tmin)\n",
    "\n",
    "def min_max_mean_scaler(t):\n",
    "    return  (t - t.mean()) / (t.max() - t.min())\n",
    "\n",
    "\n",
    "\n",
    "# Definition of 4 functions to create the composite channel with the band1 and band2\n",
    "\n",
    "def create_band3_1(band1, band2):\n",
    "    return band1/band2\n",
    "\n",
    "def create_band3_2(band1, band2):\n",
    "    return (band1 + band2) /2\n",
    "\n",
    "def create_band3_3(band1, band2):\n",
    "    return band1 - band2\n",
    "\n",
    "def create_band3_4(band1, band2):\n",
    "    return np.ones(band1.shape)\n",
    "\n",
    "\n",
    "\n",
    "# Definition of function to create tensor from image data according pre processing filter\n",
    "\n",
    "def to_tensor(data, \n",
    "              fn_process_band = no_process,\n",
    "              fn_create_band_3 = create_band3_1):\n",
    "    \"\"\"\n",
    "    Create a tensor of the data\n",
    "    \"\"\"\n",
    "    tensor = []\n",
    "    for i, row in data.iterrows():\n",
    "        band_1 = np.array(row['band_1']).reshape(75, 75)\n",
    "        band_2 = np.array(row['band_2']).reshape(75, 75)\n",
    "        if fn_create_band_3:\n",
    "            band_3 = fn_create_band_3(band_1,band_2)\n",
    "        \n",
    "        band_1 = fn_process_band(band_1)\n",
    "        band_2 = fn_process_band(band_2)\n",
    "        if fn_create_band_3:\n",
    "            band_3 = fn_process_band(band_3)\n",
    "        \n",
    "        if fn_create_band_3:\n",
    "            t = np.dstack((band_1, band_2, band_3))\n",
    "        else:\n",
    "            t = np.dstack((band_1, band_2))\n",
    "        tensor.append(t)\n",
    "    return np.array(tensor)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "def plot_stat_train_tensor(label_stat,\n",
    "                     stat,\n",
    "                     label_process = \"None to Band - minmaxscaler to tensor\",\n",
    "                     fn_process_band = no_process,\n",
    "                     fn_process_tensor = min_max_scaler,\n",
    "                           \n",
    "                     label_band3 = \"HH/HV\",\n",
    "                     fn_create_band_3 = create_band3_1):\n",
    "    tensor_data = to_tensor(train_data,\n",
    "                            fn_process_band = fn_process_band,\n",
    "                            fn_create_band_3 = fn_create_band_3)\n",
    "    tensor_iceberg_data = to_tensor(train_data[idx_iceberg_td],\n",
    "                                    fn_process_band = fn_process_band,\n",
    "                                    fn_create_band_3 = fn_create_band_3)\n",
    "    tensor_ship_data = to_tensor(train_data[idx_ship_td],\n",
    "                                 fn_process_band = fn_process_band,\n",
    "                                 fn_create_band_3 = fn_create_band_3)\n",
    "\n",
    "    img_stat_data = stat(tensor_data, axis=0) \n",
    "    img_stat_iceberg_data = stat(tensor_iceberg_data,axis=0) \n",
    "    img_stat_ship_data = stat(tensor_ship_data,axis=0) \n",
    "    \n",
    "    img_stat_data = fn_process_tensor(img_stat_data)\n",
    "    img_stat_iceberg_data = fn_process_tensor(img_stat_iceberg_data)\n",
    "    img_stat_ship_data = fn_process_tensor(img_stat_ship_data)\n",
    "    \n",
    "    \n",
    "    plt.figure()\n",
    "    plt.suptitle(label_stat+\" images processing:\"+label_process+', band3='+label_band3)\n",
    "\n",
    "    plt.subplot(1,3,1)\n",
    "    plt.title(\"Iceberg\")\n",
    "    plt.imshow(img_stat_iceberg_data)\n",
    "\n",
    "    plt.subplot(1,3,2)\n",
    "    plt.title(\"Ship\")\n",
    "    plt.imshow(img_stat_ship_data)\n",
    "\n",
    "    plt.subplot(1,3,3)\n",
    "    plt.title(\"Ship+Iceberg\")\n",
    "    plt.imshow(img_stat_data)\n",
    "\n",
    "    plt.show()\n",
    "\n",
    "\n",
    "    \n",
    "#No Process, band3 = HH/HV\n",
    "#tensor_data = to_tensor(train_data, fn_process_band=no_process, fn_create_band_3=create_band3_1)\n",
    "#tensor_iceberg_data = to_tensor(train_data[idx_iceberg_td], fn_process_band=no_process, fn_create_band_3=create_band3_1)\n",
    "#tensor_ship_data = to_tensor(train_data[idx_ship_td], fn_process_band=no_process, fn_create_band_3=create_band3_1)\n",
    "\n",
    "plot_stat_train_tensor(\"Training Mean\",np.mean)\n",
    "plot_stat_train_tensor(\"Training Std\",np.std)\n",
    "\n",
    "#No Process, band3 = (HH + HV)/2\n",
    "plot_stat_train_tensor(\"Training Mean\",np.mean, \n",
    "                       label_band3=\"(HH + HV)/2\",\n",
    "                       fn_create_band_3=create_band3_2)\n",
    "plot_stat_train_tensor(\"Training Std\",np.std, \n",
    "                       label_band3=\"(HH + HV)/2\",\n",
    "                       fn_create_band_3=create_band3_2)\n",
    "\n",
    "\n",
    "#No Process, band3 = HH - HV\n",
    "plot_stat_train_tensor(\"Training Mean\",np.mean, \n",
    "                       label_band3=\"(HH - HV)\", \n",
    "                       fn_create_band_3=create_band3_3)\n",
    "plot_stat_train_tensor(\"Training Std\",np.std,\n",
    "                       label_band3=\"(HH - HV)\",\n",
    "                       fn_create_band_3=create_band3_3)\n",
    "\n",
    "#No Process, band3 = ones\n",
    "plot_stat_train_tensor(\"Training Mean\",np.mean, \n",
    "                       label_band3=\"ones\", \n",
    "                       fn_create_band_3=create_band3_4)\n",
    "plot_stat_train_tensor(\"Training Std\",np.std, \n",
    "                       label_band3=\"ones\",\n",
    "                       fn_create_band_3=create_band3_4)\n",
    "\n",
    "\n",
    "\n",
    "#Min Max Sacler Process, band3 = default\n",
    "plot_stat_train_tensor(\"Training Mean\",np.mean, \n",
    "                       label_process =  \"min_max_scaler to Band - None to tensor \",\n",
    "                       fn_process_band = min_max_scaler,\n",
    "                       fn_process_tensor = no_process)\n",
    "plot_stat_train_tensor(\"Training Std\",np.std, \n",
    "                       label_process =  \"min_max_scaler to Band - None to tensor \",\n",
    "                       fn_process_band = min_max_scaler,\n",
    "                       fn_process_tensor = no_process)\n",
    "\n",
    "\n",
    "\n",
    "#Min Max Sacler Process, band3 = (HH + HV)/2\n",
    "plot_stat_train_tensor(\"Training Mean\",np.mean, label_band3=\"(HH + HV)/2\",\n",
    "                       label_process =  \"min_max_scaler to Band - None to tensor \",\n",
    "                       fn_process_band = min_max_scaler,\n",
    "                       fn_process_tensor = no_process,\n",
    "                       fn_create_band_3=create_band3_2)\n",
    "plot_stat_train_tensor(\"Training Std\",np.std, \n",
    "                       label_band3=\"(HH + HV)/2\",\n",
    "                       label_process =  \"min_max_scaler to Band - None to tensor \",\n",
    "                       fn_process_band = min_max_scaler,\n",
    "                       fn_process_tensor = no_process,\n",
    "                       fn_create_band_3=create_band3_2)\n",
    "\n",
    "\n",
    "\n",
    "#Min Max Sacler Process, band3 = ones\n",
    "plot_stat_train_tensor(\"Training Mean\",np.mean, label_band3=\"ones\",\n",
    "                       label_process =  \"min_max_scaler to Band - None to tensor \",\n",
    "                       fn_process_band = min_max_scaler,\n",
    "                       fn_process_tensor = no_process,\n",
    "                       fn_create_band_3=create_band3_4)\n",
    "plot_stat_train_tensor(\"Training Std\",np.std, \n",
    "                       label_band3=\"ones\",\n",
    "                       label_process =  \"min_max_scaler to Band - None to tensor \",\n",
    "                       fn_process_band = min_max_scaler,\n",
    "                       fn_process_tensor = no_process,\n",
    "                       fn_create_band_3=create_band3_4)\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "def plot_stat_test_tensor(label_stat,\n",
    "                     stat,\n",
    "                     label_process=\"None to band, min_max_scaler to tensor\",\n",
    "                     fn_process_band = no_process,\n",
    "                     fn_process_tensor = min_max_scaler,\n",
    "                     label_band3=\"HH/HV\",\n",
    "                     fn_create_band_3=create_band3_1):\n",
    "    tensor_data = to_tensor(test_data, fn_process_band=fn_process_band, fn_create_band_3=fn_create_band_3)\n",
    "    img_mean_data = stat(tensor_data,axis=0) \n",
    "    \n",
    "    img_mean_data = fn_process_tensor(img_mean_data)\n",
    "    \n",
    "    plt.figure()\n",
    "    plt.suptitle(label_stat+\" images process:\"+label_process+' band3='+label_band3)\n",
    "\n",
    "    plt.subplot(1,1,1)\n",
    "    plt.title(\"All\")\n",
    "    plt.imshow(img_mean_data)\n",
    "\n",
    "    plt.show()\n",
    "\n",
    "if  'test_data' not in locals():\n",
    "    print(\"Skip this Kernel because test_data not defined\")\n",
    "else:\n",
    "    plot_stat_test_tensor(\"Test Mean\",np.mean)\n",
    "    plot_stat_test_tensor(\"Test Std\",np.std)\n",
    "\n",
    "\n",
    "    #Min Max Sacler Process, band3 = default\n",
    "    plot_stat_test_tensor(\"Test Mean\",np.mean, \n",
    "                           label_process =  \"min_max_scaler to Band - None to tensor \",\n",
    "                           fn_process_band = min_max_scaler,\n",
    "                           fn_process_tensor = no_process)\n",
    "    plot_stat_test_tensor(\"Test Std\",np.std, \n",
    "                           label_process =  \"min_max_scaler to Band - None to tensor \",\n",
    "                           fn_process_band = min_max_scaler,\n",
    "                           fn_process_tensor = no_process)\n",
    "\n",
    "    plot_stat_test_tensor(\"Test Mean\",np.mean, label_band3=\"(HH + HV)/2\",\n",
    "                           label_process =  \"min_max_scaler to Band - None to tensor \",\n",
    "                           fn_process_band = min_max_scaler,\n",
    "                           fn_process_tensor = no_process,\n",
    "                           fn_create_band_3=create_band3_2)\n",
    "    plot_stat_test_tensor(\"Test Std\",np.std, \n",
    "                           label_band3=\"(HH + HV)/2\",\n",
    "                           label_process =  \"min_max_scaler to Band - None to tensor \",\n",
    "                           fn_process_band = min_max_scaler,\n",
    "                           fn_process_tensor = no_process,\n",
    "                           fn_create_band_3=create_band3_2)\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Plot Mean images for Iceberg and Ship and For Each bands\n",
    "def extract_band(data,\n",
    "                 pre_process = no_process):\n",
    "    tensor_band_1 = []\n",
    "    tensor_band_2 = []\n",
    "    \n",
    "    for i, row in data.iterrows():\n",
    "        band_1 = np.array(row['band_1']).reshape(75, 75)\n",
    "        band_2 = np.array(row['band_2']).reshape(75, 75)\n",
    "        \n",
    "        band_1 = pre_process(band_1)\n",
    "        band_2 = pre_process(band_2)\n",
    "        \n",
    "        \n",
    "        tensor_band_1.append(band_1)\n",
    "        tensor_band_2.append(band_2)\n",
    "        \n",
    "    return (np.array(tensor_band_1),np.array(tensor_band_2))\n",
    "\n",
    "def plot_mean_by_band(label_process=\"default\",\n",
    "                      pre_process = no_process):\n",
    "    tensor_data = extract_band(train_data, pre_process=pre_process)\n",
    "    tensor_iceberg_data = extract_band(train_data[idx_iceberg_td], pre_process=pre_process)\n",
    "    tensor_ship_data = extract_band(train_data[idx_ship_td], pre_process=pre_process)\n",
    "\n",
    "    img_mean_data_1, img_mean_data_2 = [np.mean(td,axis=0) for td in tensor_data]\n",
    "    img_mean_iceberg_data_1, img_mean_iceberg_data_2 = [np.mean(td,axis=0) for td in tensor_iceberg_data]\n",
    "    img_mean_ship_data_1, img_mean_ship_data_2 = [np.mean(td,axis=0) for td in tensor_ship_data]\n",
    "\n",
    "\n",
    "    plt.figure()\n",
    "    plt.suptitle(\"Mean images (Band_1,Band_2) process:\"+label_process)\n",
    "\n",
    "    plt.subplot(2,3,1)\n",
    "    plt.title(\"Iceberg Band_1\")\n",
    "    plt.imshow(img_mean_iceberg_data_1)\n",
    "\n",
    "    plt.subplot(2,3,2)\n",
    "    plt.title(\"Ship Band_1\")\n",
    "    plt.imshow(img_mean_ship_data_1)\n",
    "\n",
    "    plt.subplot(2,3,3)\n",
    "    plt.title(\"Ship+Iceberg Band 1\")\n",
    "    plt.imshow(img_mean_data_1)\n",
    "\n",
    "\n",
    "    plt.subplot(2,3,4)\n",
    "    plt.title(\"Iceberg Band_2\")\n",
    "    plt.imshow(img_mean_iceberg_data_2)\n",
    "\n",
    "    plt.subplot(2,3,5)\n",
    "    plt.title(\"Ship Band_2\")\n",
    "    plt.imshow(img_mean_ship_data_2)\n",
    "\n",
    "    plt.subplot(2,3,6)\n",
    "    plt.title(\"Ship+Iceberg Band 2\")\n",
    "    plt.imshow(img_mean_data_2)\n",
    "\n",
    "    plt.show()\n",
    "    \n",
    "plot_mean_by_band()\n",
    "#plot_mean_by_band(label_process='Max Scaling', pre_process= lambda t : min_max_scaler(t,1))\n",
    "plot_mean_by_band(label_process='Max Scaling and recenter', pre_process= lambda t : min_max_mean_scaler(t))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Extraction highest pixel "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def extract2(v):\n",
    "    v_arr = np.array(v)\n",
    "    return v_arr*(v_arr > (v_arr.mean() + 2*v_arr.std()))\n",
    "\n",
    "from scipy.ndimage import gaussian_filter\n",
    "from skimage.morphology import opening, closing\n",
    "from skimage.morphology import erosion, dilation\n",
    "from skimage.morphology import square\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "\n",
    "def show_image(data, title):\n",
    "    img_original = np.array(data).reshape((75,75))\n",
    "    img_blured = gaussian_filter(img_original, 2)\n",
    "    \n",
    "\n",
    "    img_original_high = extract2(img_original.reshape(75*75)).reshape((75,75))\n",
    "    img_blured_high = extract2(img_blured.reshape(75*75)).reshape((75,75))\n",
    "\n",
    "    img_orig_high_opening = opening(img_original_high, square(3))\n",
    "    img_blured_high_opening = opening(img_blured_high, square(3))\n",
    "\n",
    "    img_orig_high_closing = closing(img_original_high, square(3))\n",
    "    img_blured_high_closing = closing(img_blured_high, square(3))\n",
    "\n",
    "\n",
    "    \n",
    "    img_orig_high_erosion = erosion(img_original_high)\n",
    "    img_blured_high_erosion = erosion(img_blured_high)\n",
    "\n",
    "    img_orig_high_dilation = dilation(img_original_high)\n",
    "    img_blured_high_dilation = dilation(img_blured_high)\n",
    "\n",
    "\n",
    "    nRow = 6\n",
    "    fig = plt.figure(figsize=(30,30))\n",
    "    plt.suptitle(title, fontsize=40)\n",
    "    plt.subplot(nRow,2,1)\n",
    "    plt.title(\"Original\")\n",
    "    plt.imshow(img_original)\n",
    "    plt.subplot(nRow,2,2)\n",
    "    plt.title(\"Blured\")\n",
    "    plt.imshow(img_blured)\n",
    "\n",
    "    plt.subplot(nRow,2,3)\n",
    "    plt.title(\"Original: High pixel extracted\")\n",
    "    plt.imshow(img_original_high)\n",
    "    plt.subplot(nRow,2,4)\n",
    "    plt.title(\"Blurred: High pixel extracted\")\n",
    "    plt.imshow(img_blured_high)\n",
    "\n",
    "    plt.subplot(nRow,2,5)\n",
    "    plt.title(\"Original + opening: High pixel extracted\")\n",
    "    plt.imshow(img_orig_high_opening)\n",
    "    plt.subplot(nRow,2,6)\n",
    "    plt.title(\"blured + opening: High pixel extracted\")\n",
    "    plt.imshow(img_blured_high_opening)\n",
    "\n",
    "    plt.subplot(nRow,2,7)\n",
    "    plt.title(\"Original + closing: High pixel extracted\")\n",
    "    plt.imshow(img_orig_high_closing)\n",
    "    plt.subplot(nRow,2,8)\n",
    "    plt.title(\"blured + closing: High pixel extracted\")\n",
    "    plt.imshow(img_blured_high_closing)\n",
    "    \n",
    "    plt.subplot(nRow,2,9)\n",
    "    plt.imshow(img_orig_high_erosion)\n",
    "    plt.title(\"Original + Erosion: High pixel extracted\")\n",
    "    plt.subplot(nRow,2,10)\n",
    "    plt.imshow(img_blured_high_erosion)\n",
    "    plt.title(\"Blured + erosion: High pixel extracted\")\n",
    "    \n",
    "\n",
    "    plt.subplot(nRow,2,11)\n",
    "    plt.title(\"Original + dilation: High pixel extracted\")\n",
    "    plt.imshow(img_orig_high_dilation)\n",
    "    plt.subplot(nRow,2,12)\n",
    "    plt.title(\"Blured + Dilation: High pixel extracted\")\n",
    "    plt.imshow(img_blured_high_dilation)\n",
    "    \n",
    "icebergs_samples = train_data[train_data['is_iceberg'] == True].sample(5)\n",
    "ships_samples = train_data[train_data['is_iceberg'] == False].sample(5)\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "for i, sample in icebergs_samples.iterrows():\n",
    "    show_image(sample['band_1'], \"Iceberg: Band_1 for sample \" + str(i))\n",
    "for i, sample in ships_samples.iterrows():\n",
    "    show_image(sample['band_1'], \"Ships: Band_1 for sample \" + str(i))\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Shape Exploration"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from skimage import measure\n",
    "\n",
    "def show_contour(data):\n",
    "    img_original = np.array(data).reshape((75,75))\n",
    "    img_original_high = extract2(img_original.reshape(75*75)).reshape((75,75))\n",
    "    img_orig_high_dilation = dilation(img_original_high)\n",
    "    img_orig_closed = closing(img_original_high)\n",
    "    \n",
    "    img_blured = gaussian_filter(img_original, 2)\n",
    "    img_blured_high = extract2(img_blured.reshape(75*75)).reshape((75,75))\n",
    "    img_blured_high_dilation = dilation(img_blured_high)\n",
    "    img_blured_closed = closing(img_blured_high) \n",
    "\n",
    "    \n",
    "    r = img_orig_closed\n",
    "    #r = img_blured_high_dilation\n",
    "    # Find contours at a constant value of 0.8\n",
    "    contours = measure.find_contours(r, 1.0)\n",
    "    label_img = measure.label(r)\n",
    "    props = measure.regionprops(label_img)\n",
    "    \n",
    "\n",
    "    # Display the image and plot all contours found\n",
    "    fig, ax = plt.subplots()\n",
    "    ax.imshow(r, interpolation='nearest')\n",
    "    \n",
    "    print('num contours: {}'.format(len(contours)))\n",
    "    for n, contour in enumerate(contours):\n",
    "        ax.plot(contour[:, 1], contour[:, 0], linewidth=2)\n",
    "    return props\n",
    "\n",
    "for i, sample in icebergs_samples.iterrows():\n",
    "    show_contour(sample['band_1'])\n",
    "    \n",
    "for i, sample in ships_samples.iterrows():\n",
    "    show_contour(sample['band_1'])\n",
    " \n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Statistical analisis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# #############################################################################\n",
    "# Extract Statistical features\n",
    "\n",
    "from skimage.morphology import erosion, dilation, closing\n",
    "from scipy.ndimage import gaussian_filter\n",
    "\n",
    "# create statistical dataframe object\n",
    "stat_feature = pd.DataFrame(train_data['is_iceberg'])\n",
    "stat_feature['inc_angle'] = train_data['inc_angle_corrected']\n",
    "\n",
    "#some usefull functions to extract number of most relevant pixels\n",
    "def pass_dilation(v):\n",
    "    v_arr = np.array(v).reshape((75,75))\n",
    "    return dilation(v_arr)\n",
    "def pass_closing(v):\n",
    "    v_arr = np.array(v).reshape((75,75))\n",
    "    return closing(v_arr)\n",
    "\n",
    "def pass_guassian_filter(v):\n",
    "    v_arr = np.array(v).reshape((75,75))\n",
    "    return gaussian_filter(v_arr, 2)\n",
    "def extract_high_value_pixel(v):\n",
    "    v_arr = np.array(v)\n",
    "    return v_arr*(v_arr > (v_arr.mean() + 2*v_arr.std()))\n",
    "def count_pixels_positive(v):\n",
    "    v_arr = np.array(v)\n",
    "    return (v_arr > 0).sum() / v_arr.shape[0] \n",
    "\n",
    "# extend train_data with band_e_1 and band_e_2 with high values pixels \n",
    "train_data['band_blurred_1'] = train_data.iloc[:,0].apply(pass_guassian_filter)\n",
    "train_data['band_blurred_2'] = train_data.iloc[:,1].apply(pass_guassian_filter)\n",
    "\n",
    "#train_data['band_dilation_1'] = train_data.iloc[:,0].apply(pass_dilation)\n",
    "#train_data['band_dilation_2'] = train_data.iloc[:,1].apply(pass_dilation)\n",
    "\n",
    "#train_data['band_highvalue_1'] = train_data.iloc[:,0].apply(extract_high_value_pixel)\n",
    "#train_data['band_highvalue_2'] = train_data.iloc[:,1].apply(extract_high_value_pixel)\n",
    "\n",
    "train_data['band_highvalue_blurred_1'] = train_data['band_blurred_1'].apply(extract_high_value_pixel)\n",
    "train_data['band_highvalue_blurred_2'] = train_data['band_blurred_2'].apply(extract_high_value_pixel)\n",
    "\n",
    "train_data['band_highvalue_original_1'] = train_data['band_1'].apply(extract_high_value_pixel)\n",
    "train_data['band_highvalue_original_2'] = train_data['band_2'].apply(extract_high_value_pixel)\n",
    "\n",
    "train_data['band_closed_highvalue_original_1'] = train_data['band_highvalue_original_1'].apply(pass_closing)\n",
    "train_data['band_closed_highvalue_original_2'] = train_data['band_highvalue_original_2'].apply(pass_closing)\n",
    "\n",
    "def create_stat_features(data):\n",
    "    #fill state_features\n",
    "    stat = pd.DataFrame()\n",
    "    for i in [1,2]:\n",
    "        print(\"calculating fmax_...\")\n",
    "        stat['fmax_'+str(i)] = data['band_'+str(i)].map(lambda x: np.array(x).max())\n",
    "        print(\"calculating fmin_...\")\n",
    "        stat['fmin_'+str(i)] = data['band_'+str(i)].map(lambda x: np.array(x).min())\n",
    "        print(\"calculating fmean_...\")\n",
    "        stat['fmean_'+str(i)] = data['band_'+str(i)].map(lambda x: np.array(x).mean())\n",
    "        print(\"calculating fstd_...\")\n",
    "        stat['fstd_'+str(i)] = data['band_'+str(i)].map(lambda x: np.array(x).std())\n",
    "        print(\"calculating p25_...\")\n",
    "        stat['p25_'+str(i)] = [np.sort(np.array(x))[int(0.25*75*75)] for x in data['band_'+str(i)] ]\n",
    "        print(\"calculating fmedian_...\")\n",
    "        stat['fmedian_'+str(i)] = [np.sort(np.array(x))[int(0.50*75*75)] for x in data['band_'+str(i)] ]\n",
    "        print(\"calculating p75_...\")\n",
    "        stat['p75_'+str(i)] = [np.sort(np.array(x))[int(0.75*75*75)] for x in data['band_'+str(i)] ]\n",
    "        print(\"calculating fargmax_...\")\n",
    "        stat['fargmax_'+str(i)] = data['band_'+str(i)].map(lambda x: np.argmax(np.array(x)))\n",
    "        print(\"calculating fhighpixel...\")\n",
    "        stat['fhighpixel_highvalue_original'+str(i)] = data['band_highvalue_original_'+str(i)].map(count_pixels_positive)\n",
    "        print(\"calculating fhighpixel on blurred image...\")\n",
    "        stat['fhighpixel_highvalue_blurred_'+str(i)] = data['band_highvalue_blurred_'+str(i)].map(count_pixels_positive)\n",
    "        print(\"calculating fhighpixel on closed image...\")\n",
    "        stat['fhighpixel_closed_highvalue_original_'+str(i)] = data['band_closed_highvalue_original_'+str(i)].map(count_pixels_positive)\n",
    "        \n",
    "    return stat\n",
    "\n",
    "\n",
    "stat_features = create_stat_features(train_data)\n",
    "\n",
    "stat_features = pd.concat([pd.DataFrame(train_data['is_iceberg']),\n",
    "                           pd.DataFrame(train_data['inc_angle_corrected']), stat_features], axis=1)\n",
    "    \n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# #########################################################################\n",
    "# Study of correlation of statistical features\n",
    "import seaborn as sns\n",
    "sns.heatmap(stat_features.corr())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We observe a relative strong negative correlation between the target is_iceberg and the features max value of pixel in band1 (f_max1) and in band2 (f_max2), and also the standart desviacion in band2. The ship may reflects more the signal of the radar since they presents more straight angle in the architecture of the ship which favorize the reflectance. It is algo logic to observe the same correlation with the features fhighpixel in both band which sets the percentage of high value pixel. \n",
    "\n",
    "It is also interesting to note the negative correlation between the inc_angle, and the following variables: min and mean value of band1 and band2, and the different percentil of band1 and band2. The strong correlation shows that we must take into account this variable in our future model.\n",
    "\n",
    "The value of argmax1 and argmax2 does not show any interesting correlated inormation. All the image must have been previously centered.\n",
    "\n",
    "The different percentils variable are highly correlated in the same band, also between bands."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# #########################################################################\n",
    "# Study of hisogram of features and scatter plot with each pair of parameters\n",
    "import seaborn as sns\n",
    "sns.pairplot(stat_features, hue=\"is_iceberg\", plot_kws={\"s\": 3})\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The histograms pictured on the diagonal show some variable have different distribution depending \n",
    "of the target class.\n",
    "The study of pair plot in 2D showed also that in some cases the 2 classes seems separable even in 2D. (for example std_2/fmax_2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Benchmark based on SVN and Random Forest"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "##############################################################################\n",
    "# Split into a training set and a test set using a stratified k fold\n",
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "print(\"Splitting data into training and test data set...\")\n",
    "X_train, X_test, y_train, y_test = train_test_split(train_data, train_data['is_iceberg'], \n",
    "                                                    test_size=0.33, random_state=seed)\n",
    "print(\"shape X_train {}\".format(X_train.shape))\n",
    "print(\"shape X_test {}\".format(X_test.shape))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Random Forest on statistical features and feature importance"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "##############################################################################\n",
    "# Random Forest\n",
    "##############################################################################\n",
    "\n",
    "##############################################################################\n",
    "# Creating features for training\n",
    "\n",
    "print('Creating X_train_stat_feature. and X_test_stat_feature..')\n",
    "X_train_rforest_feature= pd.concat([create_stat_features(pd.DataFrame(X_train)),\n",
    "                                 pd.DataFrame(X_train['inc_angle_corrected'])], axis=1)\n",
    "\n",
    "X_test_rforest_feature= pd.concat([create_stat_features(pd.DataFrame(X_test)),\n",
    "                                 pd.DataFrame(X_test['inc_angle_corrected'])], axis=1)\n",
    "\n",
    "print(\"shape X_train_rforest_feature {}\".format(X_train_rforest_feature.shape))\n",
    "print(\"shape X_test_rforest_feature {}\".format(X_test_rforest_feature.shape))\n",
    "\n",
    "\n",
    "##############################################################################\n",
    "# fitting Random Forest\n",
    "from sklearn.ensemble import RandomForestClassifier\n",
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "clf = RandomForestClassifier(max_depth=3, random_state=0)\n",
    "clf.fit(X_train_rforest_feature, y_train.values.ravel())\n",
    "\n",
    "\n",
    "##############################################################################\n",
    "# Evaluating importance\n",
    "importances = clf.feature_importances_\n",
    "indices = np.argsort(importances)[::-1]\n",
    "std = np.std([tree.feature_importances_ for tree in clf.estimators_],\n",
    "             axis=0)\n",
    "\n",
    "# Print the feature ranking\n",
    "print(\"Feature ranking:\")\n",
    "for f in range(X_train_rforest_feature.shape[1]):\n",
    "    print(\"%d. feature %d (%f) %s\" % (f + 1, indices[f], importances[indices[f]], X_train_rforest_feature.columns[indices[f]]))\n",
    "\n",
    "# Plot the feature importances of the forest\n",
    "plt.figure()\n",
    "plt.title(\"Feature importances\")\n",
    "plt.bar(range(X_train_rforest_feature.shape[1]), importances[indices],\n",
    "       color=\"r\", yerr=std[indices], align=\"center\")\n",
    "plt.xticks(range(X_train_rforest_feature.shape[1]), indices)\n",
    "plt.xlim([-1, X_train_rforest_feature.shape[1]])\n",
    "plt.show()\n",
    "\n",
    "\n",
    "\n",
    "# #############################################################################\n",
    "# Analisis Results on stat features\n",
    "\n",
    "\n",
    "from sklearn.metrics import average_precision_score\n",
    "from sklearn.metrics import log_loss\n",
    "from sklearn.metrics import classification_report\n",
    "from sklearn.metrics import confusion_matrix\n",
    "\n",
    "y_pred = clf.predict(X_test_rforest_feature)\n",
    "\n",
    "def print_result(X_test, y_test, y_pred):\n",
    "    print('log loss score: {0:0.2f}'.format(log_loss(y_test, y_pred)))\n",
    "    print('Mean Accuracy score: {0:0.2f}'.format(clf.score(X_test, y_test)))\n",
    "    print(classification_report(y_test, y_pred, target_names=['Ship','Iceberg']))\n",
    "    print('Confusion Matrix: {}'.format( confusion_matrix(y_test, y_pred)))\n",
    "    \n",
    "print_result(X_test_rforest_feature, y_test, y_pred)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The feature importance vector extracted from Random Classifier confirms that inc_angle, fmax_1, fmin_1\n",
    "f_std1 and f_mean1 are the most discriminative features.\n",
    "\n",
    "Interestingly enough, no feature from band2 is selected\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## SVN based on Statistical features"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "##############################################################################\n",
    "# SVN on Statistical features\n",
    "##############################################################################\n",
    "\n",
    "##############################################################################\n",
    "# Creating features for training\n",
    "\n",
    "print('Creating X_train_stat_feature. and X_test_stat_feature..')\n",
    "X_train_stat_feature= pd.concat([create_stat_features(pd.DataFrame(X_train)),\n",
    "                                 pd.DataFrame(X_train['inc_angle_corrected'])], axis=1)\n",
    "\n",
    "X_test_stat_feature= pd.concat([create_stat_features(pd.DataFrame(X_test)),\n",
    "                                 pd.DataFrame(X_test['inc_angle_corrected'])], axis=1)\n",
    "print(\"shape X_train_stat_feature {}\".format(X_train_stat_feature.shape))\n",
    "print(\"shape X_test_stat_feature {}\".format(X_test_stat_feature.shape))\n",
    "\n",
    "drop_columns = True\n",
    "if drop_columns:\n",
    "    drop_cols = ['fargmax_2', 'fhighpixel_highvalue_original2','fmedian_2',\n",
    "                 'fargmax_1','p75_2','fmin_2','fhighpixel_closed_highvalue_original_2', 'p25_2', 'fmean_2', 'fhighpixel_highvalue_blurred_2']\n",
    "    X_train_stat_feature = X_train_stat_feature.drop(drop_cols, axis=1)\n",
    "    X_test_stat_feature = X_test_stat_feature.drop(drop_cols, axis=1)\n",
    "    print(\"shape X_train_stat_feature {} After droping useless features\".format(X_train_stat_feature.shape))\n",
    "    print(\"shape X_test_stat_feature {} After droping useless features\".format(X_test_stat_feature.shape))\n",
    "\n",
    "\n",
    "\n",
    "###########################################################################\n",
    "## Scaling feature\n",
    "from sklearn.preprocessing import MaxAbsScaler, StandardScaler\n",
    "\n",
    "print(\"Scaling Features...\")\n",
    "scaler = MaxAbsScaler()\n",
    "scaler = StandardScaler()\n",
    "X_train_stat_feature = scaler.fit_transform(X_train_stat_feature)\n",
    "\n",
    "# #############################################################################\n",
    "# Train a SVM classification model on stat features\n",
    "from sklearn.svm import SVC\n",
    "from sklearn.model_selection import GridSearchCV\n",
    "\n",
    "print(\"Fitting a SVM...\")\n",
    "param_grid = {'C': [1e2, 1e3, 5e3, 1e4, 5e4, 1e5],\n",
    "              'gamma': [0.0001, 0.0005, 0.001, 0.005, 0.01, 0.1], }\n",
    "clf = GridSearchCV(SVC(kernel='rbf', class_weight='balanced'), param_grid)\n",
    "clf = clf.fit(X_train_stat_feature, y_train.values.ravel())\n",
    "print(\"Best estimator found by grid search:\")\n",
    "print(clf.best_estimator_)\n",
    "\n",
    "\n",
    "# #############################################################################\n",
    "# Analisis Results on stat features\n",
    "\n",
    "\n",
    "from sklearn.metrics import average_precision_score\n",
    "from sklearn.metrics import log_loss\n",
    "from sklearn.metrics import classification_report\n",
    "from sklearn.metrics import confusion_matrix\n",
    "\n",
    "\n",
    "X_test_stat_feature = scaler.transform(X_test_stat_feature)\n",
    "y_pred = clf.predict(X_test_stat_feature)\n",
    "\n",
    "def print_result(X_test, y_test, y_pred):\n",
    "    print('log loss score: {0:0.2f}'.format(log_loss(y_test, y_pred)))\n",
    "    print('Mean Accuracy score: {0:0.2f}'.format(clf.score(X_test, y_test)))\n",
    "    print(classification_report(y_test, y_pred, target_names=['Ship','Iceberg']))\n",
    "    print('Confusion Matrix: {}'.format( confusion_matrix(y_test, y_pred)))\n",
    "    \n",
    "print_result(X_test_stat_feature, y_test, y_pred)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Result when using all features:\n",
    "  \n",
    "log loss score: 6.45\n",
    "Mean Accuracy score: 0.81\n",
    "             precision    recall  f1-score   support\n",
    "\n",
    "       Ship       0.84      0.80      0.82       281\n",
    "    Iceberg       0.79      0.83      0.81       249\n",
    "\n",
    "avg / total       0.81      0.81      0.81       530\n",
    "\n",
    "Confusion Matrix: [[225  56]\n",
    " [ 43 206]]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## SVM based on PCA "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# #############################################################################\n",
    "# utils function\n",
    "\n",
    "\n",
    "def prepare_pca_for_band(band_id, X_train, X_test):\n",
    "    '''Extract band_1 and band_2 series and prepare to pass a PCA'''\n",
    "    train_band_data_pca = pd.DataFrame()\n",
    "    test_band_data_pca = pd.DataFrame()\n",
    "    for index, row in X_train.iterrows():\n",
    "        train_band_data_pca[index] = pd.Series(np.array(row['band_'+str(band_id)]))\n",
    "    for index, row in X_test.iterrows():\n",
    "        test_band_data_pca[index] = pd.Series(np.array(row['band_'+str(band_id)]))\n",
    "    train_band_data_pca = train_band_data_pca.T\n",
    "    test_band_data_pca = test_band_data_pca.T\n",
    "    return (train_band_data_pca, test_band_data_pca)\n",
    "\n",
    "from sklearn.decomposition import PCA\n",
    "def process_pca(X_train, X_test, n_comp):\n",
    "    '''Process PCA on X_train. \n",
    "    return a transform version of X_train and X_test'''\n",
    "    pca = PCA(n_components=n_comp, svd_solver='randomized', whiten=True).fit(X_train)\n",
    "    X_train_pca = pca.transform(X_train)\n",
    "    X_test_pca = pca.transform(X_test)\n",
    "    return (pca, X_train_pca, X_test_pca)\n",
    "\n",
    "\n",
    "def create_pca_features(X_pca_band_1,X_pca_band_2, data):\n",
    "    '''Creating dataframe from info of band1, band2 and inc_angle'''\n",
    "    train_pca_features = pd.concat([pd.DataFrame(X_pca_band_1),\n",
    "                                    pd.DataFrame(X_pca_band_2),\n",
    "                                    pd.DataFrame(data['inc_angle_corrected'].values)\n",
    "                                   ], axis=1)\n",
    "    return train_pca_features\n",
    "\n",
    "\n",
    "from sklearn.svm import SVC\n",
    "from sklearn.model_selection import GridSearchCV\n",
    "\n",
    "def train_svm(X, y):\n",
    "    param_grid = {'C': [1e3, 5e3, 1e4, 5e4, 1e5],\n",
    "              'gamma': [0.0001, 0.0005, 0.001, 0.005, 0.01, 0.1]}\n",
    "    clf = GridSearchCV(SVC(kernel='rbf', class_weight='balanced', verbose= True), param_grid)\n",
    "    clf = clf.fit(X, y)\n",
    "    print(\"Best estimator found by grid search:\")\n",
    "    print(clf.best_estimator_)\n",
    "    return clf\n",
    "\n",
    "\n",
    "\n",
    "from sklearn.metrics import average_precision_score\n",
    "from sklearn.metrics import log_loss\n",
    "from sklearn.metrics import classification_report\n",
    "from sklearn.metrics import confusion_matrix\n",
    "\n",
    "def eval_clf(clf, X_test, y_test):\n",
    "    '''Print evaluation of classifier'''\n",
    "    y_pred = clf.predict(X_test)\n",
    "    print('log loss score: {0:0.2f}'.format(log_loss(y_test, y_pred)))\n",
    "    print('Mean Accuracy score: {0:0.2f}'.format(clf.score(X_test, y_test)))\n",
    "    print(classification_report(y_test, y_pred, target_names=['Ship','Iceberg']))\n",
    "    print('Confusion Matrix: {}'.format( confusion_matrix(y_test, y_pred)))\n",
    "\n",
    "\n",
    "# #############################################################################\n",
    "# Train a SVM classification model on PCA features\n",
    "\n",
    "print(\"Preparing info of band_1 for pca...\")\n",
    "b1 = prepare_pca_for_band(1, X_train, X_test)\n",
    "print(\"Preparing info of band_2 for pca...\")\n",
    "b2 = prepare_pca_for_band(2, X_train, X_test)\n",
    "\n",
    "n_components = 10\n",
    "print(\"Processing PCA with {} components for band_1...\".format(n_components))\n",
    "pca1 = process_pca(b1[0], b1[1], 10)\n",
    "print(\"Processing PCA with {} components for band_2...\".format(n_components))\n",
    "pca2 = process_pca(b2[0], b2[1], 10)\n",
    "\n",
    "print(\"Aggregating features for training...\")\n",
    "pca_features_train = create_pca_features(pca1[1], pca2[1], X_train)\n",
    "print(\"Aggregating features for testing...\")\n",
    "pca_features_test  = create_pca_features(pca1[2], pca2[2], X_test)\n",
    "\n",
    "\n",
    "## Scaling feature\n",
    "from sklearn.preprocessing import MaxAbsScaler, StandardScaler\n",
    "print(\"Scaling Features...\")\n",
    "scaler = MaxAbsScaler()\n",
    "scaler = StandardScaler()\n",
    "pca_features_train = scaler.fit_transform(pca_features_train)\n",
    "pca_features_test = scaler.transform(pca_features_test)\n",
    "\n",
    "print(\"Fitting a SVM...\")\n",
    "svm_pca = train_svm(pca_features_train, y_train.values.ravel())\n",
    "\n",
    "print(\"Evaluating a SVM...\")\n",
    "eval_clf(svm_pca, pca_features_test, y_test)    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Learning basic CNN"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Using TensorFlow backend.\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "keras version 2.1.2\n"
     ]
    }
   ],
   "source": [
    "###############################################\n",
    "## Importing CNN stuff\n",
    "\n",
    "import mxnet\n",
    "mxnet.__version__\n",
    "import keras\n",
    "print('keras version {}'.format(keras.__version__))\n",
    "\n",
    "from keras.models import Sequential\n",
    "from keras.models import Model\n",
    "from keras.layers import Conv2D, MaxPooling2D, Dense, Dropout, Input, Flatten, Activation\n",
    "\n",
    "from keras.layers.normalization import BatchNormalization\n",
    "from keras.optimizers import Adam, SGD, Adagrad, RMSprop, Nadam\n",
    "from keras.callbacks import ModelCheckpoint, Callback, ReduceLROnPlateau, EarlyStopping, TensorBoard, LearningRateScheduler\n",
    "from keras.callbacks import CSVLogger\n",
    "import matplotlib.pyplot as plt\n",
    "from IPython.display import clear_output\n",
    "\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.metrics import log_loss\n",
    "\n",
    "\n",
    "from keras import backend as K\n",
    "#K.set_image_dim_ordering('th')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/bicho/miniconda3/envs/py3_6_keras2/lib/python3.6/site-packages/sklearn/model_selection/_split.py:2026: FutureWarning: From version 0.21, test_size will always complement train_size unless both are specified.\n",
      "  FutureWarning)\n"
     ]
    }
   ],
   "source": [
    "##################################\n",
    "## Split data into train, validation and test set\n",
    "from sklearn.model_selection import train_test_split\n",
    "def indexes_train_test_split(num, random_state, train_size=0.80):\n",
    "    idx_train,idx_test= train_test_split(np.arange(num), random_state=random_state, train_size=train_size, shuffle=True)\n",
    "    return idx_train, idx_test\n",
    "    \n",
    "    \n",
    "\n",
    "#first split the training_validation and test set\n",
    "idx_train_valid, idx_test = indexes_train_test_split(train_data.shape[0],\n",
    "                                               random_state=seed,\n",
    "                                               train_size=0.80)\n",
    "\n",
    "\n",
    "#then split the training_validation and into training and validation set\n",
    "_idx_train, _idx_valid = indexes_train_test_split(len(idx_train_valid),\n",
    "                                               random_state=seed,\n",
    "                                               train_size=0.80)\n",
    "idx_train, idx_valid = np.array(idx_train_valid)[_idx_train], np.array(idx_train_valid)[_idx_valid]\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "###############################################\n",
    "## Defining callbacks\n",
    "import csv\n",
    "\n",
    "class PlotLearning(Callback):\n",
    "    def on_train_begin(self, logs={}):\n",
    "        self.i = 0\n",
    "        self.x = []\n",
    "        self.losses = []\n",
    "        self.val_losses = []\n",
    "        self.acc = []\n",
    "        self.val_acc = []\n",
    "        self.fig = plt.figure()\n",
    "        \n",
    "        self.logs = []\n",
    "\n",
    "    def on_epoch_end(self, epoch, logs={}):\n",
    "        \n",
    "        self.logs.append(logs)\n",
    "        self.x.append(self.i)\n",
    "        self.losses.append(logs.get('loss'))\n",
    "        self.val_losses.append(logs.get('val_loss'))\n",
    "        self.acc.append(logs.get('acc'))\n",
    "        self.val_acc.append(logs.get('val_acc'))\n",
    "        self.i += 1\n",
    "        f, (ax1, ax2) = plt.subplots(1, 2, sharex=True)\n",
    "        \n",
    "        clear_output(wait=True)\n",
    "        \n",
    "        ax1.set_yscale('log')\n",
    "        ax1.plot(self.x, self.losses, label=\"loss\")\n",
    "        ax1.plot(self.x, self.val_losses, label=\"val_loss\")\n",
    "        ax1.legend()\n",
    "        \n",
    "        ax2.plot(self.x, self.acc, label=\"accuracy\")\n",
    "        ax2.plot(self.x, self.val_acc, label=\"validation accuracy\")\n",
    "        ax2.legend()\n",
    "        \n",
    "        plt.show();\n",
    "\n",
    "def lr_schedule(epoch):\n",
    "    return lr * (0.1 ** int(epoch / 10))\n",
    "\n",
    "def create_standart_callbacks(filepath):\n",
    "    #es = EarlyStopping('val_loss', patience=patience, mode=\"min\")\n",
    "    #lrr = LearningRateScheduler(lr_schedule)\n",
    "    msave = ModelCheckpoint(filepath, save_best_only=True, save_weights_only=True)\n",
    "    csv_logger = CSVLogger(filepath+'.csv', append=True, separator=';')\n",
    "    cs = [msave, csv_logger]\n",
    "    if K._BACKEND == 'tensorflow':\n",
    "        #tb= TensorBoard(log_dir='./logs_tensorboard', histogram_freq=0, batch_size=32, write_graph=True, write_grads=False, write_images=False, embeddings_freq=0, embeddings_layer_names=None, embeddings_metadata=None)\n",
    "        pl = PlotLearning()\n",
    "        cs.append(pl)\n",
    "    elif K._BACKEND == 'mxnet':\n",
    "        pl = PlotLearning()\n",
    "        cs.append(pl)\n",
    "    return cs\n",
    "\n",
    "\n",
    "    \n",
    "###############################################\n",
    "## Utils functions for plotting\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "def subplot_hist(history, title, num_rows, idx_row):\n",
    "    plt.subplot(num_rows,2,2*idx_row +1)\n",
    "    plt.plot(history['acc'])\n",
    "    plt.plot(history['val_acc'])\n",
    "    plt.title(title + ' model accuracy')\n",
    "    plt.ylabel('accuracy')\n",
    "    plt.xlabel('epoch')\n",
    "    plt.legend(['train', 'test'], loc='upper left')\n",
    "\n",
    "    plt.subplot(num_rows ,2,2*idx_row + 2)\n",
    "    # summarize history for loss\n",
    "    plt.plot(history['loss'])\n",
    "    plt.plot(history['val_loss'])\n",
    "    plt.title(title + ' model loss')\n",
    "    plt.ylabel('loss')\n",
    "    plt.xlabel('epoch')\n",
    "    plt.legend(['train', 'test'], loc='upper left')\n",
    "\n",
    "def plot_models(models):\n",
    "    #fig = plt.figure()\n",
    "    fig = plt.figure(figsize=(50,50))\n",
    "    \n",
    "    num_plots = len(models)\n",
    "    idx = 0\n",
    "    for m in models:\n",
    "        subplot_hist(m.history, m.name, num_plots, idx)\n",
    "        idx = idx + 1\n",
    "\n",
    "def history_from_csv(path_file):\n",
    "    headers_ordered = [\"epoch\", \"acc\", \"loss\", \"val_acc\", \"val_loss\"]\n",
    "    history = {\"loss\":[], \"val_loss\":[], \"acc\":[], \"val_acc\":[],}\n",
    "    with open(path_file) as f:\n",
    "        reader = csv.reader(f, delimiter=';', quotechar='|')\n",
    "        i = 0;\n",
    "        for row in reader:\n",
    "            i = i + 1\n",
    "            if i > 2:    \n",
    "                for title,idx in zip(headers_ordered[1:],range(1,len(headers_ordered))):\n",
    "                    history[title].append(float(row[idx]))\n",
    "    return history\n",
    "\n",
    "\n",
    "def plot_csv_result(name_files, path_data):\n",
    "    #fig = plt.figure()\n",
    "    fig = plt.figure(figsize=(50,50))\n",
    "    \n",
    "    num_plots = len(name_files)\n",
    "    idx = 0\n",
    "    for name_file in name_files:\n",
    "        path_file = path_data + K._BACKEND + \"/\" +  name_file\n",
    "        history = history_from_csv(path_file)\n",
    "        subplot_hist(history, name_file, num_plots, idx)\n",
    "        idx = idx + 1\n",
    "\n",
    "###############################################\n",
    "## Defining Wrapper on keras Datagenerator model\n",
    "from os import path, makedirs\n",
    "import shutil\n",
    "class MyDataGenerator:\n",
    "    def __init__(self,\n",
    "                 datagen,\n",
    "                 need_fit,\n",
    "                 **kwargs):\n",
    "        self.datagen = datagen\n",
    "        self.need_fit = need_fit\n",
    "        self.kwargs = kwargs\n",
    "        if self.kwargs[\"save_to_dir\"]:\n",
    "            p = self.kwargs[\"save_to_dir\"]\n",
    "            if path.exists(p):\n",
    "                print(\"Cleaning existing directory {}\".format(p))\n",
    "                shutil.rmtree(p)\n",
    "            print(\"creating new directory {}\".format(p))\n",
    "            makedirs(p, exist_ok=True)\n",
    "\n",
    "\n",
    "###############################################\n",
    "## Defining Wrapper on keras model\n",
    "\n",
    "class MyModel:\n",
    "    def __init__(self,\n",
    "                 name,\n",
    "                 model,\n",
    "                 tensor_generator,\n",
    "                 callbacks,\n",
    "                 nb_epoch=20,\n",
    "                 batch_size=32):\n",
    "        self.name = name\n",
    "        self.model = model\n",
    "        self.tensor_generator = tensor_generator\n",
    "        self.callbacks = callbacks\n",
    "        self.nb_epoch=nb_epoch\n",
    "        self.batch_size=batch_size\n",
    "        self.hist = None\n",
    "        \n",
    "    def load_trained_model(self, weights_path):\n",
    "        self.model.load_weights(weights_path)\n",
    "        \n",
    "    def load_training_history(self, hist_csv_path):\n",
    "        self.history = history_from_csv(hist_csv_path)\n",
    "        \n",
    "    def reload(self):\n",
    "        print(\"Model {}: \".format(m.name))\n",
    "        for c in self.callbacks:\n",
    "            if isinstance(c, ModelCheckpoint):\n",
    "                print(\"   - loading weights from file {}\".format(c.filepath))\n",
    "                self.load_trained_model(c.filepath)\n",
    "            if isinstance(c, CSVLogger):\n",
    "                print(\"   - loading history from file {}\".format(c.filename))\n",
    "                self.load_training_history(c.filename)\n",
    "                \n",
    "        \n",
    "    def plot_results(self):\n",
    "        if self.history:\n",
    "            plot_models([self.history])\n",
    "    \n",
    "    def test(self, data, target, idx_test):\n",
    "        print(\"Preparing data for testing model:\"+self.name)\n",
    "        X_test_data, y_test = data.iloc[idx_test], target.iloc[idx_test]\n",
    "        X_test = self.tensor_generator(X_test_data)\n",
    "        \n",
    "        print(X_test.shape)\n",
    "        print(y_test.shape)\n",
    "        \n",
    "        print(\"Reloading best weights\")\n",
    "        self.reload()\n",
    "        \n",
    "        probas_res =  self.model.predict_proba(X_test)\n",
    "        lloss = log_loss(y_test, probas_res)\n",
    "        return probas_res,lloss\n",
    "        \n",
    "    def train2(self, data, target, validation_split=0.25,\n",
    "              verbose=1):\n",
    "        #fit\n",
    "        print(\"Preparing data for fitting model:\"+self.name)\n",
    "        X_train_data, y_train = data, target\n",
    "        \n",
    "        X_train = self.tensor_generator(X_train_data)\n",
    "        \n",
    "        \n",
    "        print(\"Callbacks:\"+str(self.callbacks))\n",
    "        \n",
    "        print(\"Fitting Model\"+self.name)\n",
    "        \n",
    "        print(\"   params:batch_size:{}, nb_epoch:{}\".format(self.batch_size,self.nb_epoch))\n",
    "        self.hist = self.model.fit(X_train, y_train,\n",
    "                  batch_size=self.batch_size,\n",
    "                  nb_epoch=self.nb_epoch,\n",
    "                  verbose=verbose,\n",
    "                  validation_split=validation_split,\n",
    "                  callbacks=self.callbacks)\n",
    "        self.history = self.hist.history\n",
    "    \n",
    "    def train(self, data, target, idx_train, idx_valid, \n",
    "              train_img_generator = None, \n",
    "              valid_img_generator = None,\n",
    "              verbose=1):\n",
    "        #fit\n",
    "        print(\"Preparing data for fitting model:\"+self.name)\n",
    "        X_train_data, y_train = data.iloc[idx_train], target.iloc[idx_train]\n",
    "        X_valid_data, y_valid = data.iloc[idx_valid], target.iloc[idx_valid]\n",
    "        \n",
    "        \n",
    "        X_train = self.tensor_generator(X_train_data)\n",
    "        X_valid = self.tensor_generator(X_valid_data)\n",
    "        \n",
    "        \n",
    "        #print(X_train.shape)\n",
    "        #print(X_valid.shape)\n",
    "        #print(y_train.shape)\n",
    "        #print(y_valid.shape)\n",
    "        #X_train = np.swapaxes(X_train, 1,3)\n",
    "        #X_valid = np.swapaxes(X_valid, 1,3)\n",
    "        \n",
    "        print(\"Callbacks:\"+str(self.callbacks))\n",
    "        \n",
    "        print(\"Fitting Model\"+self.name)\n",
    "        \n",
    "        if train_img_generator == None:\n",
    "            print(\"   params:batch_size:{}, nb_epoch:{}\".format(self.batch_size,self.nb_epoch))\n",
    "            self.hist = self.model.fit(X_train, y_train,\n",
    "                  batch_size=self.batch_size,\n",
    "                  nb_epoch=self.nb_epoch,\n",
    "                  verbose=verbose,\n",
    "                  validation_data=(X_valid, y_valid),\n",
    "                  callbacks=self.callbacks)\n",
    "        else:\n",
    "            if train_img_generator.need_fit:\n",
    "                #Only required if featurewise_center or featurewise_std_normalization or zca_whitening\n",
    "                print(\"    Fitting Training Image generator...\")\n",
    "                train_img_generator.datagen.fit(X_train)\n",
    "            steps_per_epoch = int(X_train.shape[0]/self.batch_size  )                                   \n",
    "            \n",
    "            if valid_img_generator:\n",
    "                if valid_img_generator.need_fit:\n",
    "                    #Only required if featurewise_center or featurewise_std_normalization or zca_whitening\n",
    "                    print(\"    Fitting validation Image generator...\")\n",
    "                    valid_img_generator.datagen.fit(X_valid)\n",
    "                validation_steps = int(X_valid.shape[0]/self.batch_size)\n",
    "                print(\"   Training with train and valid generator, params: batch_size:{} steps_per_epoch:{}, nb_epoch:{} validation_steps:{}\".format(self.batch_size,\n",
    "                                                                                    steps_per_epoch,\n",
    "                                                                                    self.nb_epoch,\n",
    "                                                                                    validation_steps))\n",
    "                self.hist = self.model.fit_generator(train_img_generator.datagen.flow(X_train,\n",
    "                                                                      y_train,\n",
    "                                                                      batch_size=self.batch_size,\n",
    "                                                                      **train_img_generator.kwargs),\n",
    "                                                 steps_per_epoch = steps_per_epoch,\n",
    "                                                 nb_epoch=self.nb_epoch,\n",
    "                                                 verbose=verbose,\n",
    "                                                 validation_data=valid_img_generator.datagen,\n",
    "                                                 validation_steps=validation_steps,\n",
    "                                                 callbacks=self.callbacks,\n",
    "                                                 )\n",
    "            else:\n",
    "                print(\"   Training with train generator, params: batch_size:{} steps_per_epoch:{}, nb_epoch:{}\".format(self.batch_size,\n",
    "                                                                                    steps_per_epoch,\n",
    "                                                                                    self.nb_epoch))\n",
    "            \n",
    "                self.hist = self.model.fit_generator(train_img_generator.datagen.flow(X_train,\n",
    "                                                                      y_train,\n",
    "                                                                      batch_size=self.batch_size,\n",
    "                                                                      **train_img_generator.kwargs),\n",
    "                                                 steps_per_epoch = steps_per_epoch,\n",
    "                                                 nb_epoch=self.nb_epoch,\n",
    "                                                 verbose=verbose,\n",
    "                                                 validation_data=(X_valid, y_valid),\n",
    "                                                 callbacks=self.callbacks)\n",
    "            \n",
    "        self.history = self.hist.history\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "############################################\n",
    "## define our Basic Model 1, without regularization techniques\n",
    "\n",
    "def getModel_1(with_compile=True,\n",
    "               with_bn=False,\n",
    "               with_dropout=False,\n",
    "               optimizer=Adam(lr=0.001, beta_1=0.9, beta_2=0.999, epsilon=1e-08, decay=0.0)):\n",
    "    #Building the model\n",
    "    gmodel=Sequential()\n",
    "    #Conv Layer 1\n",
    "    gmodel.add(Conv2D(64, (3, 3),activation='relu', input_shape=(75, 75,3)))\n",
    "    if with_bn:\n",
    "        gmodel.add(BatchNormalization())\n",
    "    gmodel.add(MaxPooling2D(pool_size=(3, 3), strides=(2, 2)))\n",
    "    if with_dropout:\n",
    "        gmodel.add(Dropout(0.2))\n",
    "    \n",
    "    #Conv Layer 2\n",
    "    gmodel.add(Conv2D(128, (3, 3), activation='relu' ))\n",
    "    if with_bn:\n",
    "        gmodel.add(BatchNormalization())\n",
    "    gmodel.add(MaxPooling2D(pool_size=(2, 2), strides=(2, 2)))\n",
    "    if with_dropout:\n",
    "        gmodel.add(Dropout(0.2))\n",
    "    \n",
    "    #Conv Layer 3\n",
    "    gmodel.add(Conv2D(128, (3, 3), activation='relu'))\n",
    "    if with_bn:\n",
    "        gmodel.add(BatchNormalization())\n",
    "    gmodel.add(MaxPooling2D(pool_size=(2, 2), strides=(2, 2)))\n",
    "    if with_dropout:\n",
    "        gmodel.add(Dropout(0.2))\n",
    "    \n",
    "    #Conv Layer 4\n",
    "    gmodel.add(Conv2D(64, (3, 3), activation='relu'))\n",
    "    if with_bn:\n",
    "        gmodel.add(BatchNormalization())\n",
    "    gmodel.add(MaxPooling2D(pool_size=(2, 2), strides=(2, 2)))\n",
    "    if with_dropout:\n",
    "        gmodel.add(Dropout(0.2))\n",
    "    \n",
    "    #Flatten the data for upcoming dense layers\n",
    "    gmodel.add(Flatten())\n",
    "\n",
    "    #Dense Layers\n",
    "    gmodel.add(Dense(512))\n",
    "    if with_bn:\n",
    "        gmodel.add(BatchNormalization())\n",
    "    gmodel.add(Activation('relu'))\n",
    "    if with_dropout:\n",
    "        gmodel.add(Dropout(0.2))\n",
    "    \n",
    "    #Dense Layer 2\n",
    "    gmodel.add(Dense(256))\n",
    "    if with_bn:\n",
    "        gmodel.add(BatchNormalization())\n",
    "    gmodel.add(Activation('relu'))\n",
    "    if with_dropout:\n",
    "        gmodel.add(Dropout(0.2))\n",
    "    \n",
    "    #Sigmoid Layer\n",
    "    gmodel.add(Dense(1))\n",
    "    gmodel.add(Activation('sigmoid'))\n",
    "\n",
    "\n",
    "    if with_compile:\n",
    "        gmodel.compile(loss='binary_crossentropy',\n",
    "                  optimizer=optimizer,\n",
    "                  metrics=['accuracy'])\n",
    "    return gmodel\n",
    "\n",
    "\n",
    "\n",
    "def get_callback_1(file_weight, path_data = \"./data/models/\"):\n",
    "    path_weight = path_data + K._BACKEND + \"/\" +  file_weight\n",
    "    callbacks = create_standart_callbacks(path_weight)\n",
    "    callbacks.append(ReduceLROnPlateau(monitor='val_loss',\n",
    "                                   factor=0.5, \n",
    "                                   patience=3, \n",
    "                                   verbose=1, \n",
    "                                   mode='auto', epsilon=0.0001, cooldown=0, min_lr=0))\n",
    "    return callbacks\n",
    "\n",
    "def get_callback_simple(file_weight, path_data = \"./data/models/\"):\n",
    "    path_weight = path_data + K._BACKEND + \"/\" +  file_weight\n",
    "    callbacks = create_standart_callbacks(path_weight)\n",
    "    return callbacks\n",
    "\n",
    "\n",
    "         "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Study on pre-process functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#######################################################\n",
    "##  creations of models according preprocessing functions\n",
    "## and third band\n",
    "\n",
    "m1_0_a = MyModel(\"Model 1 - (HH/HV) - minmaxscaler\",\n",
    "             getModel_1(), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_scaler, fn_create_band_3=create_band3_1),\n",
    "             get_callback_1(\"models1/model1_band3_1_min_max_scaler_weights.hdf5\"))\n",
    "m1_1_a = MyModel(\"Model 1 - (HH+HV)/2 - minmaxscaler\",\n",
    "             getModel_1(), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_scaler, fn_create_band_3=create_band3_2),\n",
    "             get_callback_1(\"models1/model1_band3_2_min_max_scaler_weights.hdf5\"))\n",
    "m1_2_a = MyModel(\"Model 1 - (HH-HV) - minmaxscaler\",\n",
    "             getModel_1(), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_scaler, fn_create_band_3=create_band3_3),\n",
    "             get_callback_1(\"models1/model1_band3_3_min_max_scaler_weights.hdf5\"))\n",
    "\n",
    "m1_0_b = MyModel(\"Model 1 - (HH/HV) - no process\",\n",
    "             getModel_1(), \n",
    "             lambda data: to_tensor(data, fn_process_band=no_process, fn_create_band_3=create_band3_1),\n",
    "             get_callback_1(\"models1/model1_band3_1_no_process_weights.hdf5\"))\n",
    "m1_1_b = MyModel(\"Model 1 - (HH+HV)/2 - no process\",\n",
    "             getModel_1(), \n",
    "             lambda data: to_tensor(data, fn_process_band=no_process, fn_create_band_3=create_band3_2),\n",
    "             get_callback_1(\"models1/model1_band3_2_no_process_weights.hdf5\"))\n",
    "m1_2_b = MyModel(\"Model 1 - (HH-HV) - no process\",\n",
    "             getModel_1(), \n",
    "             lambda data: to_tensor(data, fn_process_band=no_process, fn_create_band_3=create_band3_3),\n",
    "             get_callback_1(\"models1/model1_band3_3_no_process_weights.hdf5\"))\n",
    "\n",
    "m1_0_c = MyModel(\"Model 1 - (HH/HV) - minmaxmean scaler\",\n",
    "             getModel_1(), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_mean_scaler, fn_create_band_3=create_band3_1),\n",
    "             get_callback_1(\"models1/model1_band3_1_minmaxmean_weights.hdf5\"))\n",
    "m1_1_c = MyModel(\"Model 1 - (HH+HV)/2 -  minmaxmean scaler\",\n",
    "             getModel_1(), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_mean_scaler, fn_create_band_3=create_band3_2),\n",
    "             get_callback_1(\"models1/model1_band3_2_minmaxmean_weights.hdf5\"))\n",
    "m1_2_c = MyModel(\"Model 1 - (HH-HV) -  minmaxmean scaler\",\n",
    "             getModel_1(), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_mean_scaler, fn_create_band_3=create_band3_3),\n",
    "             get_callback_1(\"models1/model1_band3_3_minmaxmean_weights.hdf5\"))\n",
    "\n",
    "\n",
    "\n",
    "models_1 = [m1_0_a, m1_1_a, m1_2_a, m1_0_b, m1_1_b, m1_2_b, m1_0_c, m1_1_c, m1_2_c]\n",
    "models_1 = [m1_1_c]\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "##################################\n",
    "## train  from training all models  1 y 2\n",
    "\n",
    "for m in models_1:\n",
    "    m.train(train_data, train_data['is_iceberg'], idx_train, idx_valid)\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "##################################\n",
    "## reload  from training all models\n",
    "\n",
    "for m in models_1:\n",
    "    m.reload()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "##################################\n",
    "## plot acc and val training data\n",
    "\n",
    "plot_models(models_1)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "=> EL model with HH/HV as third band, and with a pre processing of min max scaling seems the most promising.\n",
    "We can also observe that 20 epochs seem sufficient to reach the highest validation loss. After that the model seems to overfitt since the training process still increase but the validation loss starts to decrease.\n",
    "From then, we will continue to investigate with the HH/HV and min max scaling. We will also go no longer than 20 epochs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "###################################\n",
    "## Test\n",
    "#probas_res = models_1[1].test(train_data, train_data['is_iceberg'], idx_test)\n",
    "for m in models_1:\n",
    "    _, lloss = m.test(train_data, train_data['is_iceberg'], idx_test)\n",
    "    print(\"loggloss of model {} = {}\".format(m.name, lloss))\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Optimizers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from keras.optimizers import Adam, SGD, Adagrad, RMSprop, Nadam, Adadelta\n",
    "\n",
    "# Info of default value\n",
    "#keras.optimizers.Adagrad(lr=0.01, epsilon=1e-08, decay=0.0)\n",
    "#keras.optimizers.Adadelta(lr=1.0, rho=0.95, epsilon=1e-08, decay=0.0)\n",
    "#keras.optimizers.RMSprop(lr=0.001, rho=0.9, epsilon=1e-08, decay=0.0)\n",
    "#keras.optimizers.Adam(lr=0.001, beta_1=0.9, beta_2=0.999, epsilon=1e-08, decay=0.0)\n",
    "\n",
    "optimizers = [(SGD(), \"SGD\"),\n",
    "                    (RMSprop(), \"RMSprop\"),\n",
    "                    (Adagrad(),\"Adagrad\"),\n",
    "                    (Adadelta(),\"Adadelta\"),\n",
    "                    (Adam(),\"Adam\"),\n",
    "                    (Nadam(),\"Nadam\")]\n",
    "\n",
    "models_optimizer = []\n",
    "for o,name in optimizers:\n",
    "    models_optimizer.append(MyModel(\"Model 1 - (HH+HV)/2 -  minmaxmean scaler\" + name,\n",
    "             getModel_1(optimizer=o), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_mean_scaler, fn_create_band_3=create_band3_3),\n",
    "             get_callback_1(\"models1/model1_band3_3_minmaxmean_opt_\"+name+\"_weights.hdf5\")))\n",
    "    \n",
    "    \n",
    "#### ESpecial SGD\n",
    "\n",
    "import math\n",
    "def step_decay(epoch):\n",
    "   initial_lrate = 0.01\n",
    "   drop = 0.5\n",
    "   epochs_drop = 10.0\n",
    "   lrate = initial_lrate * math.pow(drop,  \n",
    "           math.floor((1+epoch)/epochs_drop))\n",
    "   return lrate\n",
    "\n",
    "from keras.callbacks import  LearningRateScheduler\n",
    "lrate = LearningRateScheduler(step_decay)\n",
    "\n",
    "def getCallbackLearningRate(path):\n",
    "    cb = get_callback_simple(path)\n",
    "    cb.append(lrate)\n",
    "    return cb\n",
    "\n",
    "\n",
    "optimizer_SGD = [\n",
    "                        (SGD(lr=0.01, momentum=0.9, decay=0.0, nesterov=False),\"SGD_2\"),\n",
    "                        (SGD(lr=0.01, momentum=0.9, decay=0.01/50, nesterov=False),\"SGD_3\"),\n",
    "                        (SGD(lr=0.01, momentum=0.9, decay=0.01/50, nesterov=True),\"SGD_4\")]\n",
    "models_optimizer_SGD = [] \n",
    "models_optimizer_SGD.append(MyModel(\"Model 1 - (HH+HV)/2 -  minmaxmean scaler SGD step decay\",\n",
    "             getModel_1(optimizer=SGD()), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_mean_scaler, fn_create_band_3=create_band3_3),\n",
    "             get_callback_1(\"models1/model1_band3_3_minmaxmean_opt_SGD_step_decay_weights.hdf5\")))\n",
    "\n",
    "\n",
    "for o,name in optimizer_SGD:\n",
    "    models_optimizer_SGD.append(MyModel(\"Model 1 - (HH+HV)/2 -  minmaxmean scaler\",\n",
    "             getModel_1(optimizer=o), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_mean_scaler, fn_create_band_3=create_band3_3),\n",
    "             get_callback_1(\"models1/model1_band3_3_minmaxmean_opt_\"+name+\"_weights.hdf5\")))\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "#SGD optimizer also has an argument called nesterov which is set to false by default.\n",
    "#Nesterov momentum is a different version of the momentum method which has stronger theoretical converge \n",
    "#guarantees for convex functions. In practice, it works slightly better than standard momentum.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for m in models_optimizer:\n",
    "    m.train(train_data, train_data['is_iceberg'], idx_train, idx_valid)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for m in models_optimizer_SGD:\n",
    "    m.train(train_data, train_data['is_iceberg'], idx_train, idx_valid)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_models(models_optimizer)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "link: https://towardsdatascience.com/learning-rate-schedules-and-adaptive-learning-rate-methods-for-deep-learning-2c8f433990d1\n",
    "\n",
    "learning rate schedules is hard to tune in advance and they depend heavily on the type of model and problem. Another problem is that the same learning rate is applied to all parameter updates. \n",
    "Adaptive gradient descent algorithms such as Adagrad, Adadelta, RMSprop, Adam, provide a solution to these problems providing heuristics to automatically set learning rate parameters. \n",
    "\n",
    "Adagrad performs larger updates for more sparse parameters and smaller updates for less sparse parameter. It has good performance with sparse data and training large-scale neural network. However, its monotonic learning rate usually proves too aggressive and stops learning too early when training deep neural networks. Adadelta is an extension of Adagrad that seeks to reduce its aggressive, monotonically decreasing learning rate. RMSprop adjusts the Adagrad method in a very simple way in an attempt to reduce its aggressive, monotonically decreasing learning rate. Adam is an update to the RMSProp optimizer which is like RMSprop with momentum.\n",
    "        "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Regularization techniques"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Dropout"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m_dropout = MyModel(\"Model 1 - (HH+HV)/2 -  minmaxmean scaler Adam with dropout\",\n",
    "             getModel_1(with_dropout=True, \n",
    "                        optimizer=SGD(lr=0.01, momentum=0.9, decay=0.0, nesterov=False)), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_mean_scaler, fn_create_band_3=create_band3_3),\n",
    "             get_callback_1(\"models1/model1_band3_3_minmaxmean_opt_SGD_with_dropout_weights2.hdf5\"))\n",
    "\n",
    "m_dropout.train(train_data, train_data['is_iceberg'], idx_train, idx_valid)\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "=> Dropout is a technique where randomly selected neurons are ignored during training.\n",
    "They are “dropped-out” randomly. This means that their contribution to the activation\n",
    "of downstream neurons is temporally removed on the forward pass and any weight updates \n",
    "are not applied to the neuron on the backward pass.\n",
    "\n",
    "=> recomendation:\n",
    "Increase your learning rate by a factor of 10 to 100 and use a high momentum value of 0.9 or 0.99.\n",
    "Increase complexity of network\n",
    "Constrain the size of network weights.\n",
    "\n",
    "(Srivastava et al., 2014) brought Dropout as a simple way to prevent neural networks from \n",
    "overfitting. It has been proved to be significantly effective over a large range of machine\n",
    "learning areas, such as image classification (Szegedy et al., 2015), \n",
    "speech recognition (Hannun et al., 2014) and even natural language processing \n",
    "(Kim et al., 2016). Before the birth of Batch Normalization, it became a necessity of \n",
    "almost all the state-of-the-art networks and successfully boosted their performances against \n",
    "overfitting risks, despite its amazing simplicity.\n",
    "\n",
    "(Ioffe & Szegedy, 2015) demonstrated Batch Normalization (BN), a powerful skill \n",
    "that not only speeded up all the modern architectures but also improved upon their strong \n",
    "baselines by acting as regularizers. \n",
    "Therefore, BN has been implemented in nearly all the recent network structures \n",
    "(Szegedy et al., 2016; 2017; Howard et al., 2017; Zhang et al., 2017) and demonstrates \n",
    "its great practicability and effectiveness."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Batch normalization"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Training Deep Neural Networks is complicated by the fact that the distribution of each layer's inputs changes during training, as the parameters of the previous layers change. This slows down the training by requiring lower learning rates and careful parameter initialization, and makes it notoriously hard to train models with saturating nonlinearities. We refer to this phenomenon as internal covariate shift, and address the problem by normalizing layer inputs. Our method draws its strength from making normalization a part of the model architecture and performing the normalization for each training mini-batch. Batch Normalization allows us to use much higher learning rates and be less careful about initialization. It also acts as a regularizer, in some cases eliminating the need for Dropout. Applied to a state-of-the-art image classification model, Batch Normalization achieves the same accuracy with 14 times fewer training steps, and beats the original model by a s"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m_bn = MyModel(\"Model 1 - (HH+HV)/2 -  minmaxmean scaler Adam with bn\",\n",
    "             getModel_1(with_bn=True), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_mean_scaler, fn_create_band_3=create_band3_3),\n",
    "             get_callback_1(\"models1/model1_band3_3_minmaxmean_opt_Adam_with_bn_weights.hdf5\"))\n",
    "\n",
    "m_bn.train(train_data, train_data['is_iceberg'], idx_train, idx_valid)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Multi input CNN "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX0AAAD8CAYAAACb4nSYAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4xLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvAOZPmwAAIABJREFUeJztnXmYFOW1/z9nuns2GHBg2NdBUfZ1QNwQXNEroLjhDnFJXDDGn8YlUYnLjVdv1JvEmIu5iN6gQDAYUINXDIomYAADyC4CyoDCMMDAbEz3zPv7oxd6Zqq7q2frnqnzeZ556K6qt+p0U/3t0+c97zlijEFRFEVxBimJNkBRFEVpOlT0FUVRHISKvqIoioNQ0VcURXEQKvqKoigOQkVfURTFQajoK4qiOAgVfUVRFAehoq8oiuIg3Ik2ACAnJ8f07t070WYoLZS1a9ceNMZ0SMS19d5WGpO63NtJIfq9e/dmzZo1iTZDaaGIyDeJurbe20pjUpd7W8M7iqIoDkJFX1EUxUGo6CuKojiIBo/pi0gf4GdAW2PMVQ19fkVp6Xi9XvLz8ykvL0+0KUqSkJ6eTvfu3fF4PPU+ly3RF5HZwGXAAWPMoLDtE4D/AlzAH4wxzxpjdgK3isjCelunKA4kPz+frKwsevfujYgk2hwlwRhjKCwsJD8/n9zc3Hqfz254Zw4wIXyDiLiAl4FLgAHAdSIyoN4WKYrDKS8vp3379ir4CgAiQvv27Rvsl58t0TfGrAAO1dg8GthhjNlpjKkA5gGTG8QqRXE4KvhKOA15P9Qnpt8N2BP2PB84XUTaA88Aw0XkEWPML60Gi8gdwB0APXv2rIcZdaD4ABzeDRUl/j9vafV/K71Na4/SMIx9ENypibYiNke/gzWzYfBV0OG0+MZ6y6HKC2lZjWOb0uKpj+hbffUYY0wh8KNYg40xs4BZAHl5eU3bqPf3Z0Px/hgHqafV7Dj7PqAZiH5pIax4DjoPil/0i/dDRTF0Gtg4tiktnvqIfj7QI+x5d2Bf/cxpAozxf3CGXAsjp4EnE1JbQ2pm4HErcNV/hlxRIuLJ8P/rrUOM1lRCVWXD2pNAfD4fbndSFAZwDPXJ018N9BWRXBFJBaYCixvGrEbEF/igdegHvc6ErsMg5xRo0xUyTlLBVxqfoOj7yuIfa4z/rwm4/PLLGTlyJAMHDmTWrFkALF26lBEjRjB06FDOP/98AIqLi5k+fTqDBw9myJAhvP322wC0bt06dK6FCxcybdo0AKZNm8b999/P+PHjeeihh/jnP//JmWeeyfDhwznzzDPZtm0bAJWVlTzwwAOh8/7mN7/ho48+4oorrgid98MPP2TKlClN8Xa0GOymbL4FjANyRCQfeMIY8z8icg/wAf6UzdnGmE2NZmlD4Q180DyZibVDcS4hTz+26P9iySY27zt6YoO3zO/tp66s8+UHdG3DExNjh4dmz55Nu3btKCsrY9SoUUyePJnbb7+dFStWkJuby6FD/tyOp556irZt2/Lll18CcPjw4Zjn3r59O8uWLcPlcnH06FFWrFiB2+1m2bJlPProo7z99tvMmjWLXbt28a9//Qu3282hQ4fIzs7m7rvvpqCggA4dOvDaa68xffr0Or8XTsSW6Btjrouw/X3g/Qa1qLHxlvr/DX7wFKWpcdsX/do03fTXr3/9axYtWgTAnj17mDVrFmPHjg3lirdr1w6AZcuWMW/evNC47OzsmOe++uqrcblcABQVFXHLLbfw1VdfISJ4vd7QeX/0ox+Fwj/B691000388Y9/ZPr06axcuZI33nijgV6xM3BeMC0YR1VPX0kU7jRAbIl+LY/8wBZ/iLLTYHA13sf3448/ZtmyZaxcuZLMzEzGjRvH0KFDQ6GXcIwxlimF4dtq5pi3atUq9Pixxx5j/PjxLFq0iN27dzNu3Lio550+fToTJ04kPT2dq6++WucE4sR5tXfU01cSjYj//qtrTB/AVDWsTTUoKioiOzubzMxMtm7dyqpVqzh+/DiffPIJu3btAgiFdy666CJ++9vfhsYGwzudOnViy5YtVFVVhX4xRLpWt27dAJgzZ05o+0UXXcTvf/97fD5ftet17dqVrl278vTTT4fmCRT7JFT0RWSiiMwqKipquouGYvoq+koCcafXLbwTEvvGFf0JEybg8/kYMmQIjz32GGPGjKFDhw7MmjWLKVOmMHToUK699loAfv7zn3P48GEGDRrE0KFDWb58OQDPPvssl112Geeddx5dunSJeK2f/vSnPPLII5x11llUVp7ITLrtttvo2bMnQ4YMYejQobz55puhfTfccAM9evRgwAAtAhAvYpooEyAaeXl5pskaTXy9HP73cpi+FHqd0TTXVBKKiKw1xuQl4toR7+0XBkKfcXD5y7V2bdmyhf79+1uf8LsN/oncDqc5OkR5zz33MHz4cG699dZEm9JkWN0Xdbm3nRcMC3n66Ym1Q3E2nvQ6hncCHn4SOGuJYuTIkbRq1Ypf/epXiTalWeLgmL5zvSQnIyITRGSbiOwQkYcjHHONiGwWkU0i8mbY9ltE5KvA3y31MsSTEX94xxhC2TuNHNNPZtauXcuKFStIS0tLtCnNEgd7+hrTdxphlWEvxL+ifLWILDbGbA47pi/wCHCWMeawiHQMbG8HPAHk4VfetYGxsZPSrXDXRfSrrB8rShw4z9P3acqmg7FTGfZ24OWgmBtjDgS2Xwx8aIw5FNj3ITXKjcdFnT394GMVfaVuOE/0NWXTyVhVhu1W45hTgVNF5O8isirQKMjuWMBfQVZE1ojImoKCAmtL6pKyqZ6+0gA4UPQDHzS3ir4DsawMW+O5G+iLv+zIdcAfROQkm2P9G42ZZYzJM8bkdejQwdqSuqRsVhN9507kKvXDgXn6pf4PXIrzvu8UW5Vh84G/GGO8xphdwDb8XwINW1XWkxl/lc0k9vSDxdX27dvHVVdZt8YeN24csVKzX3rpJUpLS0PPL730Uo4cOdJwhiqJncg1xiwBluTl5d3eZBf1lvlFvw4cOFrOn9bmU+FLrg+c4ufu8aeQ6o76ZR6qDAvsxV8Z9voax7yD38OfIyI5+MM9O4GvgX8XkWBhmYvwT/jWjTqlbCZ/TL9r164sXFj39tgvvfQSN954I5mZ/jm3999vXqW9jDEYY0hJYqfSgdk7pXWaxN17pIzrZq3i20OlsQ9WEsIPz+1DapQfr8YYn1VlWBF5ElhjjFkc2HeRiGwGKoEHA42BEJGn8H9xADxpjKnZQtQ+nsykzd556KGH6NWrF3fddRcAM2fOJCsrix/+8IdMnjyZw4cP4/V6efrpp5k8ufo8+O7du7nsssvYuHEjZWVlTJ8+nc2bN9O/f3/Kyk683jvvvJPVq1dTVlbGVVddxS9+8Qt+/etfs2/fPsaPH09OTg7Lly+nd+/erFmzhpycHF544QVmz54N+Ffr3nfffezevZtLLrmEs88+m3/84x9069aNv/zlL2RkVA/fLlmyhKeffpqKigrat2/P3Llz6dSpE8XFxcyYMYM1a9YgIjzxxBNceeWVLF26lEcffZTKykpycnL46KOPmDlzJq1bt+aBBx4AYNCgQbz77rsAXHLJJYwfP56VK1fyzjvv8Oyzz9Z6fQCrV6/mxz/+MSUlJaSlpfHRRx9x6aWX8pvf/IZhw4YBcNZZZ/HKK68wZMiQRvjfdaTol8U9ibvnUCnXvbqKojIvi+46k+E9Y1cRVJITq8qwxpjHwx4b4P7AX82xs4HZDWJIMKZvjL8WTyT++jB87y9ZTJXvxK8DVyq46pin3nkwXPJsxN1Tp07lvvvuC4n+ggULWLp0Kenp6SxatIg2bdpw8OBBxowZw6RJkyL2b33llVfIzMxkw4YNbNiwgREjRoT2PfPMM7Rr147KykrOP/98NmzYwL333ssLL7zA8uXLycnJqXautWvX8tprr/H5559jjOH000/n3HPPJTs7m6+++oq33nqLV199lWuuuYa3336bG2+8sdr4s88+m1WrViEi/OEPf+C5557jV7/6lWVZ6IKCAssS0tHYtm0br732Gr/73e8ivr5+/fpx7bXXMn/+fEaNGsXRo0fJyMjgtttuY86cObz00kts376d48ePN5rggyMncsvj8vS/KSzh2v9eybFyH2/eNkYFX2kYPBn+cgpx9WMOD+803kTu8OHDOXDgAPv27WP9+vVkZ2fTs2dPjDE8+uijDBkyhAsuuIC9e/eyf3/ktqMrVqwIie+QIUOqCdmCBQsYMWIEw4cPZ9OmTWzevDnSaQD47LPPuOKKK2jVqhWtW7dmypQpfPrppwDk5uaGvOSRI0eye/fuWuPz8/O5+OKLGTx4MM8//zybNvlbfyxbtoy77747dFx2djarVq2yLCEdjV69ejFmzJior2/btm106dKFUaNGAdCmTRvcbjdXX3017777Ll6vl9mzZzd6ETkHevqltj39nQXFXP/q5xz3VfLm7aczsGvbRjZOcQzh3bOiNXMP98hLC+HIt/7HGe0gu1ejmXfVVVexcOFCvv/+e6ZOnQrA3LlzKSgoYO3atXg8Hnr37l2rZHJNrH4F7Nq1i//8z/9k9erVZGdnM23atJjniVYjLHxlrsvlqhZGCjJjxgzuv/9+Jk2axMcff8zMmTND561pY6SSzm63m6qqE2G1cJvDS0VHen2RzpuZmcmFF17IX/7yFxYsWBBzsru+ONDTtxfe2XGgmKmzVuGtrOKtO8ao4CsNSzCZIJ64fjCOL65Gn8idOnUq8+bNY+HChaFsnKKiIjp27IjH42H58uV88803Uc8xduxY5s6dC8DGjRvZsGEDAEePHqVVq1a0bduW/fv389e//jU0Jisri2PHjlme65133qG0tJSSkhIWLVrEOeecY/v1hJdvfv3110PbrcpCn3HGGZYlpHv37s0XX3wBwBdffBHaX5NIr69fv37s27eP1av900LHjh0LlY2+7bbbuPfeexk1apStXxb1wYGiH3sid/v+Y0ydtZIqA/PuGEO/zm2ayDjFMQTvwbhEP+DtpjS+6A8cOJBjx47RrVu3UFnkG264gTVr1pCXl8fcuXPp169f1HPceeedFBcXM2TIEJ577jlGjx4NwNChQxk+fDgDBw7kBz/4AWeddVZozB133BGaFA1nxIgRTJs2jdGjR3P66adz2223MXz4cNuvZ+bMmVx99dWcc8451eYLrMpCRyohfeWVV3Lo0CGGDRvGK6+8wqmnnmp5rUivLzU1lfnz5zNjxgyGDh3KhRdeGPq1MHLkSNq0adMkrR+dV1r5N3nQeRBcPcdy97bvj3Hdq6twpwhv3j6GUzq2tjxOaT4kZWnlTYvgT9PgzpXQqXpN+IillY99D8e+839hSArk9G0co5UmZ9++fYwbN46tW7dGTPdsqNLKDvT0y6J6+v/10XaqjGH+D89QwVcaj+A9GE+ufii8k6IrclsQb7zxBqeffjrPPPNMk+T3J3QiV0QmAhNPOeWUprtojIncg8cqOK1TFrk5rSIeoyj1JhTTj2NVrqnyC76kQFU8WT9KMnPzzTdz8803N9n1EurpG2OWGGPuaNu2CSdJfeVRRb+ozEvbDE/T2aM4kxgxfcuwa7joJ+mKXKVxaMgwvLPCO8bEnMhV0VeahGDnNovwTnp6OoWFhbU/6Mao6DsQYwyFhYWkpzdMtz9n5emHaumrp68kmGCVVwtPv3v37uTn51OrLHPJQX9Yx10EFaVw2NUEhirJQHp6Ot27d2+QczlL9GOUVa7wVVHmrVTRVxofT2TR93g8odWg1Zh7DRR/D73PgTWz4WffNbKRSkvEWeGdGA1Uisr8k2NtM1X0lUYmiuhHxFfmd1iCXbc0g0epAw4T/WB/XOuYfkj01dNXGpvwMgx2Ca4m92QABnzHG8U0pWXjMNG35+m3UdFXGpu6pGx6A5ln7jp8YShKAIeJfvSJ3KPq6StNhYhfvL1x9GcIrjGpS2hIUQI4TPSDnr6Gd5QkwJN+IqPMDr7yEzF9UNFX6oTDRD8Y048xkauirzQFcXv6ZdU9/Xi+MBQlgMNE32b2joq+0hR4MuKM6Zf5fx1EyfFXlFgkVPRFZKKIzCoqKmqaC9rw9FuluvC4nPVdqCSIYOqlHYwJS9msQy1+RQngrNo7NlI21ctXmgxPhv0MnPDV5HWpxa8oAZzl0vpie/qarqk0Ge50++Gd8F+p7sh1exQlFs4S/RhlGNTTV5oUT6b9idxw0dfsHaUeOEz0S/1eUoRGBUdV9JWmJJ6UzeBxmrKp1BOHiX70pujq6StNSjwpm1aevqZsKnXAYaJfGjG0Ayr6ShMTT8pmtZh+0NOPI8dfUQI4TPQje/reyipKK7SsstKExJOyGZy0daeDOw2Q+HL8FSWAA0U/RgkGLausNBXxpGyGpxuLBL4w1NNX4seBoq+rcZUkwZ0BVT6otNHkPCT6gXRNd5x1exQlgIp+AC2rrDQ58WTh1FxN7snU7B2lTjhM9CM3RVdPX2lyQs3RbXjsvhprTDzpKvpKnXCY6Ef29LWWvtLkxJOFU7MXhCdDwztKnXBewbUIon+kVEVfaWJC4R0b4l2zQmy8ZZkVJYDDCq6V6kSukjx44vD0feWAgCs1MDaOuj2KEoaGdwIUlXnJ1LLKLR4RmSAi20Rkh4g8bLF/mogUiMi6wN9tYfsqw7Yvrrcx8aysDaYbiwTGxlG3R1HCcCfagCYjWI9cyyo7FhFxAS8DFwL5wGoRWWyM2Vzj0PnGmHssTlFmjBnWYAbFFdMvOzHxC5qyqdQZ57i1vuhN0VX0HcFoYIcxZqcxpgKYB0xOmDVxxfRrOCyasqnUEeeIvo0GKpqj3+LpBuwJe54f2FaTK0Vkg4gsFJEeYdvTRWSNiKwSkcvrbU084R1f2Yk6+qApm0qdcZDoR++Pq2WVHYFYbDM1ni8BehtjhgDLgNfD9vU0xuQB1wMvicjJlhcRuSPw5bCmoKAgsjVBEbebshl+73oyNbyj1AkHib42UFHIB8I99+7AvvADjDGFxpjjgaevAiPD9u0L/LsT+BgYbnURY8wsY0yeMSavQ4cOka0JtT20mbIZLvrudP82U/M7S1Gi4yDRj+7pq+g7gtVAXxHJFZFUYCpQLQtHRLqEPZ0EbAlszxaRtMDjHOAsoOYEcHx44vD0feW1wzumyl7dHkUJwznZOzVrl4Tv0rLKjsAY4xORe4APABcw2xizSUSeBNYYYxYD94rIJMAHHAKmBYb3B/5bRKrwO0vPWmT9xIc7npTNUshod+J56FdCKbhT62WG4iwcKPq1J3J1YZZzMMa8D7xfY9vjYY8fAR6xGPcPYHCDGpOSAq60OGL6NVI2QeP6Stw4KLwT2dNX0VcSht3uWVYpm6ALtJS4cZDoB2P66ukrSYTdRipWKZugpRiUuHGQ6Mf29DVPX2ly3Dbz7a1SNsF+5y1FCaCij5ZVVhKInZW1xlinbIIu0FLixkGiHzllU8M7SsKws7K2sgIwNcI7cZRwUJQwHCT6kRdnFWktfSVR2FlZazUfFU9ZZkUJwzlNVIITYSm1X3JRmZcMj4tUt3O+A5UkIbiyNhqhrlnhKZtx5PgrShjOaaISo5a+evlKQrCTsqmevtKAOMe1jdEU/aRMFX0lAdhJ2Qx68xrTVxoAB4l+dE9f0zWVhGAnZTMU3rHw9DVlU4kTZ4m+VthUkg1PZhzhHYsyDJqyqcSJg0Q/clN0raWvJAyPjYncUHgn7P4Vsb+wS1HCcJDol+tErpJ8eDKhyguVvsjHRFpj4slQ0VfixkGibz2R662sokTLKiuJIlQtM4p4eyP0d3bbrNsDVFUZ/vODbew4cKwORiotCQeJvvVErpZgUBKKnSycBvD0vy4o5rfLd/DE4k11MFJpSThM9LXCppJk2Mm3t0rZDI61mbK5bs8RAP6+o5DPdxbGa6XSgnCQ6JdWz34IoKKvJBQ7zVAiNQCyW5YZv+i3TnPTMSuNF5dtr4OhSkvBQaJvHd7RsspKQgk1Q4kW0y8DBNxp1bfHkb2zPv8IQ3u05a5xJ7Nq5yH+8fXButmrNHucIfrG+D0iDe8oyYbHRr59sG6USI2x9mL65d5Ktn53jKHdT2Lq6J50bpPOix9uxxhTD8OV5oozRN8XIfsBnchVEoydZiiRVpPbFP1N+47iqzIM7XES6R4Xd48/mdW7D/P3HRrbdyLOEH1tiq4kK3ZW1kZaY2IzZXN9YBJ3WI+TALhmVA+6tk3nhQ+3qbfvQBwi+pEbqBwp1bLKSgKxFdOPsJrcpqe/bs8RurRNp1Mb/xdMmtvF3eedwhffHmHFVxrbdxrOULoYnr56+UrCsBXTL7euG2UzZXN9/hGGdj+p2rarR/ag20kZvKCxfcfhENEPePo185xR0VcSjJ1mKNFi+jHCO4dLKvimsJShPaqLfqo7hRnnncL6PUf4eFtBvFYrzRiHiH7kpugq+kpCCS3OijWRW9thwZ0BVT6o9EYcuj7fH88f2qN2o6IrR3anR7sMXlym3r6TcIjoW3QeCqC19JWEYkf0fRHKgtsIDa3fU4QIDO5WW/Q9rhRmnNeXDflFfLTlQDxWK80Yh4h+9JRN9fSVhJHiAldq3VM2g/sjsG7PYfp2bE1WuvU9PmV4N3q1z1Rv30G4E3lxEZkITDzllFMa90IxPH0VfSWhuGNk4URL2YSIXxjGGNbnF3F+v46RL+1K4d7z+vL//rSeP3y6i/5d2sRjudIEtG+d2qD/LwkVfWPMEmBJXl7e7Y16oQgxfS2rrCQFsVIvo6VsQsSx+YfLOFRSUWsStyaTh3Xl5Y938Mz7W+xarDQh/za4Cy/fMKLBzpdQ0W8yIqRsnliN64y3QUlSPDFq6ERL2YSIY9fVWJQVCbcrhQU/PINdB0tsmas0LdmZqQ16PmeonVWPUcJW42aqp680LbsOlpCb08r/JNrKWmNix/QjpHuu33OENHcKp3XOimlPTus0clqnxTxOaf44ZCI38IGq4S1pCQYlEazefYgLXviEWSu+9m+Itsiq0gumMnLKJkSsxb9uzxEGdWuLx+WMj7liD+d4+u50SKl+86voK4lgWI+TmDCwM//+/lYqq+DOaDF9n7XDAoSlbNb+wvBWVrFxXxHXj+7VQFYrLQVniL7POvtBRV9JBB5XCv81dRgpKcJ/LN3KpV0MvTIjePpRFhZGq9uzff8xyr1VlouyFGfjjN99EZqiH9UGKo5DRCaIyDYR2SEiD1vsnyYiBSKyLvB3W9i+W0Tkq8DfLfWxw+1K4cVrhjJ5WFe2HPRy8PAR6wOjiX6Upurr9xQBsSdxFefhENGP3jVLPX1nICIu4GXgEmAAcJ2IDLA4dL4xZljg7w+Bse2AJ4DTgdHAEyKSXR973K4UXrhmGB2yT6K4uJiXrNoY1tHTX7fnMNmZHnq2q+3sKM7G8aKf7kkhze1KgFFKAhgN7DDG7DTGVADzgMk2x14MfGiMOWSMOQx8CEyor0GuFGH4yV3ITvXx0rKvale9tBXTt/b0h/Y4CanZbUtxPA4R/VLLD42uxnUc3YA9Yc/zA9tqcqWIbBCRhSLSI86xcZOSmkkbt49r8rrz64++4lf/Fyb8UcM71imbxcd9bD9wrFY5ZUUBx4h+ZE9fRd9RWLm9NQvOLAF6G2OGAMuA1+MY6z9Q5A4RWSMiawoKbJQtdqcj3jKenTKE60b34LfLd/Da33f790WpG0VKCrjSaqVsbtxbhDEaz1escZDoa90dhXygR9jz7sC+8AOMMYXGmOOBp68CI+2ODTvHLGNMnjEmr0OHDrGt8mRAZQUpVPHM5YM58+T2vPrpTiqrTFh4xyJPHwKreat7+sH2iEO6a+aOUhsHib6Vp+9T0XcWq4G+IpIrIqnAVGBx+AEi0iXs6SQgWJDmA+AiEckOTOBeFNhWf8JW1qakCDef0Yvvisr5ZPuBqF3fQttrePrr9hyhZ7tM2usKW8UCB4m+dcqmpms6B2OMD7gHv1hvARYYYzaJyJMiMilw2L0isklE1gP3AtMCYw8BT+H/4lgNPBnYVn/c1WvonN+/Ezmt03jz8z1hoh/B03en14rpr99zJGaRNcW5OGNxVoQqhUVlXk7KaNhiRkpyY4x5H3i/xrbHwx4/AjwSYexsYHaDG1WjcJrHlcJVI7vz6qc7OdrrKG0ghqd/InvnwNFy9hWV8wMN7SgRcJCnX130fZVVFB/X8I6SBFhUy5w6qgeVVYaN33zv3xA1pn9i3Pp8XZSlRKfli35VlX8yrIboHy33AVpWWUkCQjH9E+LdO6cVZ57cnq17Am0MI4p+ZrXwzvo9R3ClCIMs2iMqCjhB9H3WKW9aVllJGtzWi6ymju5JRXkpla60WsUCq40Nm8hdt+cI/Tpnke7RBYeKNQ4S/eoxUS3BoCQNEZqhXDywEye5fZSbKPNOYSmbVVWG9fk6iatEp+WLfqiBitbSV5KUCM1Q0twu+uW4OVbpoeDYcYuBVEvZ3FVYwrFyH8N0Ja4SBQeIvnWes4q+kjREaYZycraLMuPh7S/yI4w9kbL5P5/tAmBEr3rVgVNaOA4Q/eievubpKwknFN6pXVM/K8WHpGYy75/fVi/EFhrrT9n835W7efPzb/nhuX04pWPrxrVXadY4QPStC1YdVU9fSRY8kT19fGW0ycpid2EpK3cWWoxNp8pbxswlmzm/X0d+enG/xrVVafY4QPQDHySL/rhaVllJCqI1OPeW0bZNW9qku5n3zz21dh/2ukmp8tI3J52Xpg7DlaKllJXoOED0rT39olIttqYkCW7r7J3gNldqBlNGdGfpxu85VFIR2lVU5mX+On8Vz1evG0RWut7PSmwcJPrVJ3KPlFWo6CvJgcsNKW5r0feVgzudqaN7UFFZxZ8DE7q+yipmvPUvvi/xe/Y9stTDV+zhINGvHd5R0VeShhora0ME+jv369yGYT1OYt7qPRhj+OVft7JiewEThueeOE5RbOAg0a+Zsql1d5QkosbK2hDe8lCFzetH92THgWIeXbSR//lsF9PO7M2YUwPNu6y+MBTFAgeIvnXKppZVVpIKT4ZlymZ4WfDLhnahdZqbt/75Lef0zeHn/9Y/rDm6evqKPRIq+iLf1IUCAAAXIUlEQVQyUURmFRUVNd5FvNadhzS8oyQVnoyIKZvBezcz1c0Pzs5laPe2/Pa6EbhdKWHN0dXTV+yRUNE3xiwxxtzRtm0jVgT0lvo/NGEFq7SsspJ0eDJqh2gqvVDlqxaavP/CU3nn7rNOFAoM7vNZTAIrigUOCO9EK6usoq8kCe6M2tk7EbpmiYRl6kSo0KkokWj5ou+r3SpR6+4oSUeNZijACc8/Ui19iFihU1Ei0fJF38LTV9FXko4abQ+BsCSECK0SQUVfiRsVfUVJBtzptePywcnZSE3R4cRqXk3ZVGziANEv1fCOkvxYpWzG5elryqZiDweIfplluiao6CtJhFXKZlwxffX0FXs4QPRre/pHtZa+kmxYpWza8fRTXOBK1ZRNxTYOEH3rmH6aO0WbRyvJgzsg+lVVJ7bZiekHx+pErmITB4h+ee2YvpZVVpKNoLCHe/uh8E5G7eNrjlXRV2ziANEv1QqbSvITqqETJt4R6kbVHquevmIfB4i+dXhHRV9JKoKTteGx+VB4J4bouzM0pq/YpmWLflVVxBW5KvpKUhHy9MPCO+rpK41AyxZ9n/VEmIq+knSEqmWGpW3ajulHKMusKBa0bNGP2EBFa+krSYZVc3RvKbjSqlWIjThWwzuKTVq26Ptqt0o8XFJB8XEf3bNjeE+K0pS4LVbWhnXNij5Ws3cU+7Rs0bfw9HceLAEgN6dVIixSFGusmqH4ymKHdkBj+kpctHDRrz0RtrOgGIA+HVonwiJFscaq7aFF5pn1WBV9xT4tXPRrh3d2HSzBnSL00PCOkky4LRZn2RV9TdlU4qCFi37t2iU7C0ro2T7T319UcRwiMkFEtonIDhF5OMpxV4mIEZG8wPPeIlImIusCf79vUMMsF2epp680PO5EG9CoWDRF33WwhD45GtpxIiLiAl4GLgTygdUistgYs7nGcVnAvcDnNU7xtTFmWKMY57Foe+grtx/Tr6yAqkp/ATZFiULLdndrTORWVhl2FZbQp4NO4jqU0cAOY8xOY0wFMA+YbHHcU8BzQNMlv1s1OLcoIWI9VhupKPZxiOj7PxT7jpRR4auij2buOJVuwJ6w5/mBbSFEZDjQwxjzrsX4XBH5l4h8IiLnNKhlLg+Iq0Z4x27KprZMVOzjjPBOwIvSdE3HIxbbTGinSArwIjDN4rjvgJ7GmEIRGQm8IyIDjTFHa11E5A7gDoCePXvat67mylrbKZsWoSFFiUAL9/Srp2xquqbjyQd6hD3vDuwLe54FDAI+FpHdwBhgsYjkGWOOG2MKAYwxa4GvgVOtLmKMmWWMyTPG5HXo0MG+dTW7Z9meyLWYBFaUCLRw0a8+kbvrYAlZaW5yWqcm0CglgawG+opIroikAlOBxcGdxpgiY0yOMaa3MaY3sAqYZIxZIyIdAhPBiEgfoC+ws0GtCzZSCeItt5myaVGhU1Ei0MLDO6X+D1KgdsnOAv8krojVr3ylpWOM8YnIPcAHgAuYbYzZJCJPAmuMMYujDB8LPCkiPqAS+JEx5lCDGlgz9TLeiVz19BUbtHDRL6s2EbbrYAmjc9sl0CAl0Rhj3gfer7Ht8QjHjgt7/DbwdqMaF94Bq9IHVV77KZugoq/YouWHdwLxzrKKSvYeKdNJXCV58WSeCNFYFAuMPE5TNhX7tGzR952YCNsVyNzRHH0laQmvlmm3axZYV+hUlAi0bNH31hZ99fSVpCU8ZdNXezV55HEWFToVJQItXPRLT+ToB9I1VfSVpCU8ZdOiWGDkcRYVOhUlAi1c9Kt7+l3appOZ2rLnrpVmTHjKZjyib1WhU1Ei0MJF/4Sn//VBrbmjJDl19vQ1e0exTwsX/TJwp2OMYVdBsVbXVJIbT7pFTN+G6Ls8kOJW0Vds0fJF35NJYUkFR8t9Gs9XkptgyqYx8Xn6obEa3lFi4wDRz9B0TaV5EB6bj1f03ek6kavYwhGiHyq0puEdJZkJj80HvXY7KZtQPTSkKFFouaJfVRVYnJXJzoISUl0pdNO+uEoyEy76NcqCxx6bqZ6+YouWK/q+Eysadx4soVf7TFwpWmhNSWLcYeUUQqJv09N3p2tMX7FFyxX9sJjozoJijecryY8nrJyCN47sHQh4+pq9o8QmqUX/64JiNu4tqtvgwE/dSlc63x4qJVfj+UqyExL9cn9oMsUDLpuLCcMrdCpKFJJa9H/4v2t59q9b6zY48AE4VOHCW2nU01eSn1C1zLJqFWLtjdWUTcUeSS36F/TvxKqdhRwt98Y/OLC45ftS/0s8WUVfSXbcNSZy7cbzQVM2FdsktehfOKAjvirDJ9sK4h8c8PTz/dmaGt5Rkp/wBue+cvvpmsGxmrKp2CCpRX9Yj2zat0pl2Zb98Q8OeD3fFBtOyvTQrpX2xVWSnGopm6Xxh3fU01ds0OAlJ0WkFfA7oAL42Bgzt67ncqUI5/XryAebvsdbWYXHFcd3VMDT/6bIaPkFpXngDo/pl8cf3tGYvmIDWyoqIrNF5ICIbKyxfYKIbBORHSLycGDzFGChMeZ2YFJ9DbxgQCeOlvtYvTvOHtQB0f/6iE9X4irNg5qLs+oykVtV1Ti2KS0Gu67zHGBC+AYRcQEvA5cAA4DrRGQA0B3YEzissr4GntM3h1R3Css2H4hvYOCnbv4x0cwdpXlQrQxDWfwxfVBvX4mJLdE3xqwAarrao4EdxpidxpgKYB4wGcjHL/y2zx+NzFQ3Z5+Sw4dbvscYY39gYFKrjFT6aHhHaQ64UkFSTqzItVtsDU78KlDRV2JQn5h+N0549OAX+9OBXwO/FZF/A5ZEGiwidwB3APTs2TPqhS7o34m/bT3AVweKObVTlj3rAp5+GWn06dC8wzter5f8/HzKy/UDHY309HS6d++Ox+NJtCl1Q8Qf1w+Fd+IQ/eCvAm8p0K5RzFNaBvURfatCNsYYUwJMjzXYGDMLmAWQl5cX1YU/v39HWAQfbt4fh+j7Y/oV4qFX+zhio0lIfn4+WVlZ9O7dGxGtH2SFMYbCwkLy8/PJzc1NtDl1J7iyNu6UzbDVvIoShfqEX/KBHmHPuwP76meONZ3apDO0e9v4Uje9pVRIGl1PakW6x9UYZjUZ5eXltG/fXgU/CiJC+/btm/+voWANnbhTNsPq9ihKFOoj+quBviKSKyKpwFRgccOYVZsL+ndi3Z4jHDhm80PtLaOctBaTrqmCH5sW8R650+uYshlWoVNRomA3ZfMtYCVwmojki8itxhgfcA/wAbAFWGCM2dRYhl4woBPGwPKt9rJ4jLeUkioPJzfzeH6y0Lq1vo9NgicDKkqg8rh6+kqjYCumb4y5LsL294H3G9SiCPTrnEW3kzL4cPN+rh0VfeIX4HhZMaUmtcV4+opD8GRA2WH/47qkbGpMX4lBUpdhCEdEuHBAJz796iBlFbHT/0tLSignVXP0GxhjDA8++CCDBg1i8ODBzJ8/H4DvvvuOsWPHMmzYMAYNGsSnn35KZWUl06ZNCx374osvJtj6ZoAnA0oD2dHxLs6CUKFBRYlEg5dhaEwu6N+JOf/YzWc7DnLhgE5Rjz1eVkxZC4rpB/nFkk1s3ne0Qc85oGsbnpg40Naxf/7zn1m3bh3r16/n4MGDjBo1irFjx/Lmm29y8cUX87Of/YzKykpKS0tZt24de/fuZeNG/0LuI0eONKjdLRJ3mKcfbxkG0Jr6SkwS6umLyEQRmVVUZK9RyujcdmSluVm2OXYWj6+8hOOSRte22he3Ifnss8+47rrrcLlcdOrUiXPPPZfVq1czatQoXnvtNWbOnMmXX35JVlYWffr0YefOncyYMYOlS5fSpk2bRJuf/HjSw8I78SzOClvNqyhRSKinb4xZAizJy8u73c7xqe4Uzj2tAx9t3U9VlSElSs/bKm8pKZ6cqMc0R+x65I1FpFXRY8eOZcWKFbz33nvcdNNNPPjgg9x8882sX7+eDz74gJdffpkFCxYwe/bsJra4meHJBALvcVwrclX0FXs0m5h+kAsHdOJgcQXr8qOHCsRbhju9eS/KSkbGjh3L/PnzqayspKCggBUrVjB69Gi++eYbOnbsyO23386tt97KF198wcGDB6mqquLKK6/kqaee4osvvki0+clP+ORtnVI2VfSV6DSrmD7AuFM74k4Rlm3ez4ie2ZbHVPiq8FSVk56haYYNzRVXXMHKlSsZOnQoIsJzzz1H586def3113n++efxeDy0bt2aN954g7179zJ9+nSqApUff/nLXybY+mZAuHcfz0SuywPiUk9fiUmzE/22mR5G57Zj2Zb9/HRCP8tj/r7jIMOoILOVin5DUVzsb0EmIjz//PM8//zz1fbfcsst3HLLLbXGJZt3LyITgP8CXMAfjDHPRjjuKuBPwChjzJrAtkeAW/FXj73XGPNBgxsYLvrxpGyK+MdqyqYSg2YX3gF/Fs/2/cV8U1gS2uarrOK9Dd8x5Xd/Z/qc1aRLBR3bW/8SUJxJlHLgNY/LAu4FPg/bNgD/qvOB+MuM/y5wvoalrp5+cKyGd5QYNFvRB1i25QBHy728umIn5z7/MXe/+QUHiyt44rJ+ZFBBVmvNFlGqEakceE2eAp4Dwt3mycA8Y8xxY8wuYEfgfA1LeMZOPDH94FgN7ygxaHbhHYCe7TM5rVMW//3J17zwf9soqajk9Nx2PD5xABf074TLVwbLiC/7QXECkcqBhxCR4UAPY8y7IvJAjbGraozt1uAWhgt9PCmbwbEq+koMmqXoA0wa1pWXlm1n4pCu/ODsXAZ1a3tiZ/DGj/fnsdLSsSwHHtopkgK8CEyLd2y1A+PoFVGL8Hs2XqfFo56+EpuEir6ITAQmnnLKKXGPvfPck/nBWblkpFqEVYNFp9TTV6oTqxx4FjAI+DhQsbMzsFhEJtkYGyKeXhG1qJayGef969aYvhKbhMb0jTFLjDF3tG3bNvbBNUhJEWvBhzBPX0VfqUbUcuDGmCJjTI4xprcxpjf+cM6kQPbOYmCqiKSJSC7QF/hng1sY9PRT3P40zLjGqqevxKbZhneiop6+YoExxiciwXLgLmC2MWaTiDwJrDHGROwHEThuAbAZ8AF3G2NiV/6Ll2BMP954Pvjv95KDDWuP0uJoeaJvDOz6xP9YY/oJoXXr1qG8/prs3r2byy67LFSEramxKgdujHk8wrHjajx/Bnim0YyDE/dsXRwWTdlUbNCyRL/8KLz7E9i4EHLPhV5nJtoiRYmPYEw/3nRN0JRNxRYtR/T3/Qv+NB2OfAvn/RzOvh9SmndvXEv++jB8/2XDnrPzYLjEcmEqAA899BC9evXirrvuAmDmzJmICCtWrODw4cN4vV6efvppJk+2SnmPTHl5OXfeeSdr1qzB7XbzwgsvMH78eDZt2sT06dOpqKigqqqKt99+m65du3LNNdeQn59PZWUljz32GNdee229XnZSEvTw6xTe0ZRNJTbNX/SNgVWvwIePQ+tOMO096HVGoq1qUUydOpX77rsvJPoLFixg6dKl/OQnP6FNmzYcPHiQMWPGMGnSpLj61L788ssAfPnll2zdupWLLrqI7du38/vf/54f//jH3HDDDVRUVFBZWcn7779P165dee+99wCwW4672REU/bqGd1T0lRg0b9EvPQTv3Anbl8Jp/waTfwuZ7RJtVeMSxSNvLIYPH86BAwfYt28fBQUFZGdn06VLF37yk5+wYsUKUlJS2Lt3L/v376dz5862z/vZZ58xY8YMAPr160evXr3Yvn07Z5xxBs888wz5+flMmTKFvn37MnjwYB544AEeeughLrvsMs4555zGermJJRTeqYPoB1M2jfHX4lEUC5plGQYqSmH9fHjlLPj6b3DJczB1bssX/ARy1VVXsXDhQubPn8/UqVOZO3cuBQUFrF27lnXr1tGpUyfKy+Mr9hWpNv/111/P4sWLycjI4OKLL+Zvf/sbp556KmvXrmXw4ME88sgjPPnkkw3xspKP+k7kAvi06JoSmebj6RsDe/4J6/4IGxdBxTHo0A+unwddhibauhbP1KlTuf322zl48CCffPIJCxYsoGPHjng8HpYvX84333wT9znHjh3L3LlzOe+889i+fTvffvstp512Gjt37qRPnz7ce++97Ny5kw0bNtCvXz/atWvHjTfeSOvWrZkzZ07Dv8hkwJ0GSN1TNsEf4tF0ZSUCyb8it2gvbJgH696Ewh3gaQUDL4dh10PPMyGlef5YaW4MHDiQY8eO0a1bN7p06cINN9zAxIkTycvLY9iwYfTrZ13mOhp33XUXP/rRjxg8eDBut5s5c+aQlpbG/Pnz+eMf/4jH46Fz5848/vjjrF69mgcffJCUlBQ8Hg+vvPJKI7zKJCBYIlk9faWRkEg/sZuSvLw8s2bNmuobjYEFN8HW98BUQa+zYNgNMGAypDmrTv6WLVvo379/os1oFli9VyKy1hiTlwh7LO/tWPxHLvS7FCa/HN+49fNh0R3Qrg+4UuMbqyQvJ58PE/7dcldd7u3kDe+I+G/ecx6AYdf5HyuKEzj/cehUh17IuWNhyLXq6bc02nRt0NMlr+gDXNhCJ+scwJdffslNN91UbVtaWhqff/55hBFKiLzpdRvXpgtMmdWwtigtjuQWfaXZMnjwYNatW5doMxRFqYHOgjYTkmHuJdnR90hRYqOi3wxIT0+nsLBQRS0KxhgKCwtJT69DzRpFcRAa3mkGdO/enfz8fAoKChJtSlKTnp5O9+7dE22GoiQ1KvrNAI/HQ25ubqLNUBSlBaDhHUVRFAehoq8oiuIgEir6IjJRRGa12DK5iqIoSUZSlGEQkQIgUsWuHCCZGn8mmz2QfDYlmz29gJ8ZY5p85ZLe2/Ui2eyB5LPpNGNMVjwDkkL0oyEiaxJVN8WKZLMHks+mZLMH1CY7qD2xSTab6mKPxvQVRVEchIq+oiiKg2gOop9sFaSSzR5IPpuSzR5Qm+yg9sQm2WyK256kj+kriqIoDUdz8PQVRVGUBiJpRV9EJojINhHZISIPJ9oeABHZLSJfisg6EYmzHVKD2TBbRA6IyMawbe1E5EMR+Srwb3aC7ZkpInsD79M6Ebm0Ce3pISLLRWSLiGwSkR8HtifsPbKwUe/t2tdPqvs6ik3N/t5OStEXERfwMnAJMAC4TkQGJNaqEOONMcMSmLY1B5hQY9vDwEfGmL7AR4HnibQH4MXA+zTMGPN+E9rjA/6fMaY/MAa4O3DvJPI9CqH3dkTmkFz3dSSboJnf20kp+sBoYIcxZqcxpgKYB0xOsE1JgTFmBXCoxubJwOuBx68DlyfYnoRhjPnOGPNF4PExYAvQjQS+RzXQe9uCZLuvo9iUMBrq3k5W0e8G7Al7nh/YlmgM8H8islZE7ki0MWF0MsZ8B/4bA+iYYHsA7hGRDYGfyAkJpYhIb2A48DnJ8x7pvW2fZPk/q0mzvreTVfTFYlsypBmdZYwZgf+n+d0iMjbRBiUprwAnA8OA74BfNbUBItIaeBu4zxhztKmvHwW9t5s3zf7eTlbRzwd6hD3vDuxLkC0hjDH7Av8eABbh/6meDOwXkS4AgX8PJNIYY8x+Y0ylMaYKeJUmfp9ExIP/QzHXGPPnwOZkeY/03rZPsvyfhWgJ93ayiv5qoK+I5IpIKjAVWJxIg0SklYhkBR8DFwEbo49qMhYDtwQe3wL8JYG2BG+8IFfQhO+TiAjwP8AWY8wLYbuS5T3Se9s+yfJ/FqJF3NvGmKT8Ay4FtgNf46+QmGh7+gDrA3+bEmUT8Bb+n5Ve/F7jrUB7/LP2XwX+bZdge/4X+BLYELghuzShPWfjD5dsANYF/i5N5HtkYaPe2/buo4T+n7XUe1tX5CqKojiIZA3vKIqiKI2Air6iKIqDUNFXFEVxECr6iqIoDkJFX1EUxUGo6CuKojgIFX1FURQHoaKvKIriIP4/F8mPThQYZfQAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x13ab082e8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Epoch 00020: ReduceLROnPlateau reducing learning rate to 0.00015624999650754035.\n",
      "\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\b\r",
      "1074/1074 [==============================] - 80s 75ms/step - loss: 7.7931 - acc: 0.5112 - val_loss: 9.8061 - val_acc: 0.3849\n"
     ]
    }
   ],
   "source": [
    "############################################\n",
    "## define our Basic Model 1, without regularization techniques\n",
    "\n",
    "def getModel_2_with_angle(with_compile=True,\n",
    "               with_bn=False,\n",
    "               with_dropout=False,\n",
    "               optimizer=Adam(lr=0.001, beta_1=0.9, beta_2=0.999, epsilon=1e-08, decay=0.0)):\n",
    "    #Building the model\n",
    "    gmodel=Sequential()\n",
    "    #Conv Layer 1\n",
    "    gmodel.add(Conv2D(64, (3, 3),activation='relu', input_shape=(75, 75,3)))\n",
    "    if with_bn:\n",
    "        gmodel.add(BatchNormalization())\n",
    "    gmodel.add(MaxPooling2D(pool_size=(3, 3), strides=(2, 2)))\n",
    "    if with_dropout:\n",
    "        gmodel.add(Dropout(0.2))\n",
    "    \n",
    "    #Conv Layer 2\n",
    "    gmodel.add(Conv2D(128, (3, 3), activation='relu' ))\n",
    "    if with_bn:\n",
    "        gmodel.add(BatchNormalization())\n",
    "    gmodel.add(MaxPooling2D(pool_size=(2, 2), strides=(2, 2)))\n",
    "    if with_dropout:\n",
    "        gmodel.add(Dropout(0.2))\n",
    "    \n",
    "    #Conv Layer 3\n",
    "    gmodel.add(Conv2D(128, (3, 3), activation='relu'))\n",
    "    if with_bn:\n",
    "        gmodel.add(BatchNormalization())\n",
    "    gmodel.add(MaxPooling2D(pool_size=(2, 2), strides=(2, 2)))\n",
    "    if with_dropout:\n",
    "        gmodel.add(Dropout(0.2))\n",
    "    \n",
    "    #Conv Layer 4\n",
    "    gmodel.add(Conv2D(64, (3, 3), activation='relu'))\n",
    "    if with_bn:\n",
    "        gmodel.add(BatchNormalization())\n",
    "    gmodel.add(MaxPooling2D(pool_size=(2, 2), strides=(2, 2)))\n",
    "    if with_dropout:\n",
    "        gmodel.add(Dropout(0.2))\n",
    "    \n",
    "    #Flatten the data for upcoming dense layers\n",
    "    gmodel.add(Flatten())\n",
    "\n",
    "    \n",
    "    \n",
    "    \n",
    "    \n",
    "    # Image input encoding\n",
    "    image_input = Input(shape=(75,75,3))\n",
    "    encoded_image = gmodel(image_input)\n",
    "\n",
    "    # Inc angle input\n",
    "    inc_angle_input = Input(shape=(1,))\n",
    "\n",
    "    # Combine image and inc angle\n",
    "    combined= keras.layers.concatenate([encoded_image,inc_angle_input])\n",
    "    \n",
    "    dense_model = Sequential()\n",
    "    #Dense Layers\n",
    "    dense_model.add(Dense(512, input_shape=(257,)))\n",
    "    if with_bn:\n",
    "        dense_model.add(BatchNormalization())\n",
    "    dense_model.add(Activation('relu'))\n",
    "    if with_dropout:\n",
    "        dense_model.add(Dropout(0.2))\n",
    "    \n",
    "    #Dense Layer 2\n",
    "    dense_model.add(Dense(256))\n",
    "    if with_bn:\n",
    "        dense_model.add(BatchNormalization())\n",
    "    dense_model.add(Activation('relu'))\n",
    "    if with_dropout:\n",
    "        dense_model.add(Dropout(0.2))\n",
    "    \n",
    "    #Sigmoid Layer\n",
    "    dense_model.add(Dense(1))\n",
    "    dense_model.add(Activation('sigmoid'))\n",
    "\n",
    "    \n",
    "    output = dense_model(combined)\n",
    "\n",
    "    # Final model\n",
    "    combined_model= Model(inputs=[image_input,inc_angle_input],outputs= output)\n",
    "\n",
    "    if with_compile:\n",
    "        combined_model.compile(loss='binary_crossentropy',\n",
    "                  optimizer=optimizer,\n",
    "                  metrics=['accuracy'])\n",
    "    return combined_model\n",
    "\n",
    "\n",
    "def to_tensor_multi_input(data, fn_process_band=None, fn_create_band_3=None):\n",
    "    print(\"to_tensor_multi_input\")\n",
    "    print(data.shape)\n",
    "    a = to_tensor(data, fn_process_band=fn_process_band, fn_create_band_3=fn_create_band_3)\n",
    "    print(a.shape)\n",
    "    b = data['inc_angle_corrected']\n",
    "    print(b.shape)\n",
    "    \n",
    "    return [a,b]\n",
    "\n",
    "m2_bn = MyModel(\"Model 2 - (HH+HV)/2 -  minmaxmean scaler Adam with bn\",\n",
    "             getModel_2_with_angle(with_dropout=True, optimizer=SGD(0.01)), \n",
    "             lambda data: to_tensor_multi_input(data, \n",
    "                                                fn_process_band=min_max_mean_scaler, \n",
    "                                                fn_create_band_3=create_band3_3),\n",
    "             get_callback_1(\"models2/model2_band3_3_minmaxmean_opt_Adam_with_bn_weights.hdf5\"))\n",
    "\n",
    "m2_bn.model.summary()\n",
    "\n",
    "m2_bn.train2(train_data,\n",
    "            train_data['is_iceberg'],\n",
    "            validation_split=0.33)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Data Augmentation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "###########################################################\n",
    "## Prepare models and datagen\n",
    "\n",
    "m1_datagen_0 = MyModel(\"Model 1 - (HH/HV) - minmaxscaler\",\n",
    "             getModel_1(), \n",
    "             lambda data: to_tensor(data, fn_process_band=min_max_scaler, fn_create_band_3=create_band3_1),\n",
    "             get_callback_1(\"model1_band3_1_min_max_scaler_weights_datagen.hdf5\"))\n",
    "\n",
    "from keras.preprocessing.image import ImageDataGenerator\n",
    "datagen1 = ImageDataGenerator(featurewise_center=False,\n",
    "                             featurewise_std_normalization=False,\n",
    "                             samplewise_center=True,\n",
    "                             samplewise_std_normalization=True,\n",
    "                             \n",
    "                             width_shift_range=0.1,\n",
    "                             height_shift_range=0.1,\n",
    "                             #zoom_range=0.1,\n",
    "                             #shear_range=0.1,\n",
    "                             rotation_range=30.,\n",
    "                             horizontal_flip=True,\n",
    "                             vertical_flip=True,\n",
    "                             fill_mode='reflect')\n",
    "\n",
    "datagen2 = ImageDataGenerator(featurewise_center=False,\n",
    "                             featurewise_std_normalization=False,\n",
    "                             samplewise_center=True,\n",
    "                             samplewise_std_normalization=True,\n",
    "                             \n",
    "                             width_shift_range=0.1,\n",
    "                             height_shift_range=0.1,\n",
    "                             #zoom_range=0.1,\n",
    "                             #shear_range=0.1,\n",
    "                             rotation_range=30.,\n",
    "                             horizontal_flip=True,\n",
    "                             vertical_flip=True,\n",
    "                             fill_mode='reflect')\n",
    "\n",
    "train_generator = MyDataGenerator(datagen1, True, \n",
    "                                  save_to_dir=\"/Users/bicho/Documents/GitHub/udacity-capstone/data/generator/model1_band3_1_min_max_scaler_weights_datagen/train/\")\n",
    "                            \n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "###########################################################\n",
    "## train withd datagen\n",
    "\n",
    "m1_datagen_0.train(train_data, train_data['is_iceberg'], idx_train, idx_valid, \n",
    "        train_img_generator=train_generator)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "###########################################################\n",
    "## a manual way\n",
    "\n",
    "_to_process = False\n",
    "\n",
    "if _to_process:\n",
    "    for e in range(1):\n",
    "        print('Epoch', e)\n",
    "        batches = 0\n",
    "        for batch in datagen_original.flow(X_train,\n",
    "                 y_train,\n",
    "                 batch_size=32, \n",
    "                 save_to_dir='/Users/bicho/Documents/GitHub/udacity-capstone/data/original',\n",
    "                                      seed =  seed):\n",
    "            break\n",
    "            if batches > 0:\n",
    "                break\n",
    "            batches += 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using Already existing CNN"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Transfer learning"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from keras.applications import vgg16\n",
    "from keras.applications import resnet50\n",
    "from keras.applications import vgg19\n",
    "from keras.applications import inception_v3\n",
    "\n",
    "apps = [vgg16.VGG16(weights='imagenet', include_top=False, input_shape=(75,75,3))]\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "inn this kernel we will use pretrained VGG-16 network which performs very well on small size images.\n",
    "\n",
    "5. VGG architecture has proved to worked well on small sized images(CIFAR-10) I expected it to work well for this dataset as well.\n",
    "\n",
    "Keras provide the implementation of pretrained VGG, it in it's library so we don't have to build the net by ourselves. Here we are removing the last layer of VGG and putting our sigmoid layer for binary predictions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "apps[0].output_shape\n",
    "keras.__version__"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def getTopModel(base_model, with_compile=True):\n",
    "    print(base_model.output_shape[1:])\n",
    "    top_model = Sequential()\n",
    "    top_model.add(Flatten(input_shape=base_model.output_shape[1:]))\n",
    "    top_model.add(Dense(256, activation='relu'))\n",
    "    top_model.add(Dropout(0.5))\n",
    "    top_model.add(Dense(1, activation='sigmoid'))\n",
    "    if with_compile:\n",
    "        top_model.compile(loss='binary_crossentropy',\n",
    "                  optimizer=Adam(lr=0.001, beta_1=0.9, beta_2=0.999, epsilon=1e-08, decay=0.0),\n",
    "                  metrics=['accuracy'])\n",
    "    return top_model\n",
    "\n",
    "\n",
    "\"\"\"\n",
    "train_data -> X_Train, X_Valid (pre_process)\n",
    "generator(train_data), generator(valid_data)\n",
    "fit base_model\n",
    "save output\n",
    "\n",
    "\"\"\"\n",
    "def train_topmodel_with_inputof_basemodel(base_model ,\n",
    "                                          idxs, \n",
    "                                          data, target,\n",
    "                                          path_save_weights,\n",
    "                                          train_gen=None, valid_gen=None, \n",
    "                                           fit_base=True,\n",
    "                                         epochs=50):\n",
    "    idx_train, idx_valid = idxs\n",
    "    X_train, y_train= data[idx_train], target[idx_train] \n",
    "    X_valid, y_valid= data[idx_valid], target[idx_valid] \n",
    "    \n",
    "    batch_size = 32\n",
    "    if fit_base:\n",
    "        print(\"Predicting base model for training...\")\n",
    "        if train_gen:\n",
    "            train_gen.fit(X_train)\n",
    "            train_generator = train_gen.flow(X_train,\n",
    "                                     batch_size=batch_size)\n",
    "            bottleneck_features_train = base_model.predict_generator(train_generator, X_train.shape[0]/batch_size)\n",
    "        else:\n",
    "            bottleneck_features_train = base_model.predict(X_train)\n",
    "        np.save(open('bottleneck_features_train.npy', 'wb'), bottleneck_features_train)\n",
    "\n",
    "        print(\"Predicting base model for validating...\")\n",
    "        #bottleneck_features_validation = base_model.predict_generator(valid_generator,  X_valid.shape[0]/batch_size)\n",
    "        if valid_gen:\n",
    "            valid_gen.fit(X_valid)\n",
    "            valid_generator = valid_gen.flow(X_valid,\n",
    "                                     batch_size=batch_size)\n",
    "            bottleneck_features_validation = base_model.predict_generator(valid_generator,  X_valid.shape[0]/batch_size)\n",
    "        else:\n",
    "            bottleneck_features_validation = base_model.predict(X_valid)\n",
    "        np.save(open('bottleneck_features_validation.npy', 'wb'), bottleneck_features_validation)\n",
    "    \n",
    "    \n",
    "    ####\n",
    "    print(\"Loading input for top model for training...\")\n",
    "    train_data = np.load(open('bottleneck_features_train.npy', 'rb'))\n",
    "    # the features were saved in order, so recreating the labels is easy\n",
    "    train_labels = y_train\n",
    "\n",
    "    print(\"Loading input for top model for validating...\")\n",
    "    validation_data = np.load(open('bottleneck_features_validation.npy', 'rb'))\n",
    "    validation_labels = y_valid\n",
    "\n",
    "    print(\"Fitting Top model...\")\n",
    "    top_model = getTopModel(base_model)\n",
    "    top_model.fit(train_data, train_labels,\n",
    "              epochs=epochs,\n",
    "              batch_size=batch_size,\n",
    "              validation_data=(validation_data, validation_labels))\n",
    "    print(\"Saving weights Top model...\")\n",
    "    top_model.save_weights(path_save_weights)\n",
    "    return top_model\n",
    "\n",
    "base_model = apps[0]\n",
    "_datas = to_tensor(train_data, fn_process_band=min_max_scaler, fn_create_band_3=create_band3_1)\n",
    "_target = train_data['is_iceberg']\n",
    "\n",
    "\n",
    "train_gen = ImageDataGenerator(featurewise_center=False,\n",
    "                             featurewise_std_normalization=False,\n",
    "                             samplewise_center=True,\n",
    "                             samplewise_std_normalization=True,                            \n",
    "                             width_shift_range=0.1,\n",
    "                             height_shift_range=0.1,\n",
    "                             #zoom_range=0.1,\n",
    "                             #shear_range=0.1,\n",
    "                             rotation_range=30.,\n",
    "                             horizontal_flip=True,\n",
    "                             vertical_flip=True,\n",
    "                             fill_mode='reflect')\n",
    "valid_gen = ImageDataGenerator(featurewise_center=False,\n",
    "                             featurewise_std_normalization=False,\n",
    "                             samplewise_center=True,\n",
    "                             samplewise_std_normalization=True,                            \n",
    "                             width_shift_range=0.1,\n",
    "                             height_shift_range=0.1,\n",
    "                             #zoom_range=0.1,\n",
    "                             #shear_range=0.1,\n",
    "                             rotation_range=30.,\n",
    "                             horizontal_flip=True,\n",
    "                             vertical_flip=True,\n",
    "                             fill_mode='reflect')\n",
    "\n",
    "\n",
    "\n",
    "top_model = train_topmodel_with_inputof_basemodel(base_model,  \n",
    "                                          (idx_train[:5], idx_valid[:5]), \n",
    "                                          _datas, _target,\n",
    "                                          \"weight_top_model.hf5\",\n",
    "                                           train_gen=None, valid_gen=None,  \n",
    "                                           fit_base=True)\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Fine Tuning"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "base_model.summary()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "len(base_model.layers)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "def train_topmodel_and_basemodel(base_model ,top_model_weights_path,\n",
    "                                          idxs, \n",
    "                                          data, target,\n",
    "                                          train_gen=None, valid_gen=None, \n",
    "                                           fit_base=True, epochs=50):\n",
    "    batch_size = 32\n",
    "    idx_train, idx_valid = idxs\n",
    "    X_train, y_train= data[idx_train], target[idx_train] \n",
    "    X_valid, y_valid= data[idx_valid], target[idx_valid] \n",
    "    \n",
    "    top_model = getTopModel(base_model)\n",
    "    top_model.load_weights(top_model_weights_path)\n",
    "\n",
    "    # CREATE AN \"REAL\" MODEL FROM VGG16\n",
    "    # BY COPYING ALL THE LAYERS OF VGG16\n",
    "    new_model = Sequential()\n",
    "    for l in base_model.layers:\n",
    "        new_model.add(l)\n",
    "    \n",
    "    # set the first 15 layers (up to the last conv block)\n",
    "    # to non-trainable (weights will not be updated)\n",
    "    for layer in new_model.layers[:15]:\n",
    "        layer.trainable = False\n",
    "        \n",
    "\n",
    "    # CONCATENATE THE TWO MODELS\n",
    "    new_model.add(top_model)\n",
    "\n",
    "    # compile the model with a SGD/momentum optimizer\n",
    "    # and a very slow learning rate.\n",
    "    new_model.compile(loss='binary_crossentropy',\n",
    "              optimizer=SGD(lr=1e-4, momentum=0.9),\n",
    "              metrics=['accuracy'])\n",
    "\n",
    "    # fine-tune the model\n",
    "    if train_gen and valid_gen:\n",
    "        train_gen.fit(X_train)\n",
    "        train_generator = train_gen.flow(X_train,\n",
    "                                     batch_size=32)\n",
    "        \n",
    "        valid_gen.fit(X_valid)\n",
    "        valid_generator = valid_gen.flow(X_valid,\n",
    "                                     batch_size=32)\n",
    "            \n",
    "        hist = new_model.fit_generator(\n",
    "            train_generator,\n",
    "            samples_per_epoch=len(X_train),\n",
    "            epochs=epochs,\n",
    "            validation_data=valid_generator,\n",
    "            nb_val_samples=len(X_valid))\n",
    "    else:\n",
    "        hist = new_model.fit(X_train, y_train,\n",
    "                  batch_size=batch_size,\n",
    "                  nb_epoch=epochs,\n",
    "                  verbose=True,\n",
    "                  validation_data=(X_valid, y_valid))\n",
    "    return new_model\n",
    "\n",
    "base_model = apps[0]\n",
    "_datas = to_tensor(train_data, fn_process_band=min_max_scaler, fn_create_band_3=create_band3_1)\n",
    "_target = train_data['is_iceberg']\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "train_gen = ImageDataGenerator(featurewise_center=False,\n",
    "                             featurewise_std_normalization=False,\n",
    "                             samplewise_center=True,\n",
    "                             samplewise_std_normalization=True,                            \n",
    "                             width_shift_range=0.1,\n",
    "                             height_shift_range=0.1,\n",
    "                             #zoom_range=0.1,\n",
    "                             #shear_range=0.1,\n",
    "                             rotation_range=30.,\n",
    "                             horizontal_flip=True,\n",
    "                             vertical_flip=True,\n",
    "                             fill_mode='reflect')\n",
    "valid_gen = ImageDataGenerator(featurewise_center=False,\n",
    "                             featurewise_std_normalization=False,\n",
    "                             samplewise_center=True,\n",
    "                             samplewise_std_normalization=True,                            \n",
    "                             width_shift_range=0.1,\n",
    "                             height_shift_range=0.1,\n",
    "                             #zoom_range=0.1,\n",
    "                             #shear_range=0.1,\n",
    "                             rotation_range=30.,\n",
    "                             horizontal_flip=True,\n",
    "                             vertical_flip=True,\n",
    "                             fill_mode='reflect')\n",
    "\n",
    "\n",
    "\n",
    "top_model = train_topmodel_and_basemodel(base_model,  \"weight_top_model.hf5\",\n",
    "                                          (idx_train, idx_valid), \n",
    "                                          _datas, _target,\n",
    "                                         train_gen=None, valid_gen=None)\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "help(base_model)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# Improvement\n",
    "\n",
    "- Use of Incident angle\n",
    "- Batch normalization\n",
    "- Use Data augmentation\n",
    "- Use spekkle filtering"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# Submission"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def predict(data, gmodel, preprocess):\n",
    "    X_test = preprocess(data)\n",
    "    return gmodel.predict_proba(X_test)\n",
    "\n",
    "def create_csv(path, predicted_test):\n",
    "    result = pd.DataFrame()\n",
    "    result['id']=test_data['id']\n",
    "    result['is_iceberg']=predicted_test.reshape((predicted_test.shape[0]))\n",
    "    result.to_csv(path, index=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_data[:10]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_test = to_tensor(test_data, fn_process_band=no_process, fn_create_band_3=create_band3_2)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_test.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "predicted_test = m1_0_b.model.predict_proba(X_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "predicted_test.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "create_csv('res.csv', predicted_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "py3_6_keras2",
   "language": "python",
   "name": "py3_6_keras2"
  },
  "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
}