{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "_cell_guid": "a64533a8-cc5d-2784-cb70-e3ad95940bf4"
      },
      "source": [
        "analysis"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "58508802-01c7-8410-404a-151e9e9b1e30"
      },
      "outputs": [],
      "source": [
        "import numpy as np\n",
        "import pandas as pd\n",
        "import matplotlib.pyplot as plt\n",
        "import seaborn as sns\n",
        "\n",
        "import cv2\n",
        "\n",
        "%matplotlib inline\n",
        "import plotly.offline as py\n",
        "py.init_notebook_mode(connected=True)\n",
        "import plotly.graph_objs as go\n",
        "import plotly.tools as tls\n",
        "\n",
        "from subprocess import check_output\n",
        "print(check_output([\"ls\", \"../input/train\"]).decode(\"utf8\"))"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "02f80a3d-d84a-eeee-f5c8-8a7aaec38068"
      },
      "outputs": [],
      "source": [
        "from glob import glob\n",
        "basepath = '../input/train/'\n",
        "\n",
        "all_cervix_images = []\n",
        "\n",
        "for path in sorted(glob(basepath + \"*\")):\n",
        "    cervix_type = path.split(\"/\")[-1]\n",
        "    cervix_images = sorted(glob(basepath + cervix_type + \"/*\"))\n",
        "    all_cervix_images = all_cervix_images + cervix_images\n",
        "\n",
        "all_cervix_images = pd.DataFrame({'imagepath': all_cervix_images})\n",
        "all_cervix_images['filetype'] = all_cervix_images.apply(lambda row: row.imagepath.split(\".\")[-1], axis=1)\n",
        "all_cervix_images['type'] = all_cervix_images.apply(lambda row: row.imagepath.split(\"/\")[-2], axis=1)\n",
        "all_cervix_images.head()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "_cell_guid": "1e8fc7e3-9fd7-ed36-0e8e-ba42032f012a"
      },
      "source": [
        "Image types\n",
        "\n",
        "Now that we have the data in a handy dataframe we can do a few aggregations on the data. Let us first see how many images there are for each cervix type and which file types they have.\n",
        "\n",
        "All files are in JPG format and Type 2 is the most common one with a little bit more than 50% in the training data in total, Type 1 on the other hand has a little bit less than 20% in the training data.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "fb388ed1-6709-ed17-02e0-189c6c86fb85"
      },
      "outputs": [],
      "source": [
        "print('We have a total of {} images in the whole dataset'.format(all_cervix_images.shape[0]))\n",
        "type_aggregation = all_cervix_images.groupby(['type', 'filetype']).agg('count')\n",
        "type_aggregation_p = type_aggregation.apply(lambda row: 1.0*row['imagepath']/all_cervix_images.shape[0], axis=1)\n",
        "\n",
        "fig, axes = plt.subplots(nrows=2, ncols=1, figsize=(10, 8))\n",
        "\n",
        "type_aggregation.plot.barh(ax=axes[0])\n",
        "axes[0].set_xlabel(\"image count\")\n",
        "type_aggregation_p.plot.barh(ax=axes[1])\n",
        "axes[1].set_xlabel(\"training size fraction\")\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "_cell_guid": "1d2a27f8-470f-f432-0fd3-1a3e7ad3a7be"
      },
      "source": [
        "Now, lets read the files for each type to get an idea about how the images look like."
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "ae499bbb-131d-a7e4-c35b-0a8581811ca5"
      },
      "outputs": [],
      "source": [
        "\n",
        "\n",
        "fig = plt.figure(figsize=(12,8))\n",
        "\n",
        "i = 1\n",
        "for t in all_cervix_images['type'].unique():\n",
        "    ax = fig.add_subplot(1,3,i)\n",
        "    i+=1\n",
        "    f = all_cervix_images[all_cervix_images['type'] == t]['imagepath'].values[0]\n",
        "    plt.imshow(plt.imread(f))\n",
        "    plt.title('sample for cervix {}'.format(t))\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "_cell_guid": "73a7974e-4948-67e8-08da-833b6dc8038e"
      },
      "source": [
        "Image dimension\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "a8facbc8-f854-d9d7-855a-c8301ab2f268"
      },
      "outputs": [],
      "source": [
        "from collections import defaultdict\n",
        "from skimage.io import imread, imshow\n",
        "\n",
        "images = defaultdict(list)\n",
        "\n",
        "for t in all_cervix_images['type'].unique():\n",
        "    sample_counter = 0\n",
        "    for _, row in all_cervix_images[all_cervix_images['type'] == t].iterrows():\n",
        "        try:\n",
        "            img = imread(row.imagepath)\n",
        "            sample_counter +=1\n",
        "            images[t].append(img)\n",
        "        except:\n",
        "             print('image read failed for {}'.format(row.imagepath))\n",
        "        if sample_counter > 35:\n",
        "             break\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "fe6cf431-ce54-1605-feb8-6da1ca967ade"
      },
      "outputs": [],
      "source": [
        "\n",
        "            \n",
        "dfs = []\n",
        "for t in all_cervix_images['type'].unique():\n",
        "    t_ = pd.DataFrame(\n",
        "        {\n",
        "            'nrows': list(map(lambda i: i.shape[0], images[t])),\n",
        "            'ncols': list(map(lambda i: i.shape[1], images[t])),\n",
        "            'nchans': list(map(lambda i: i.shape[2], images[t])),\n",
        "            'type': t\n",
        "        }\n",
        "    )\n",
        "    dfs.append(t_)\n",
        "\n",
        "shapes_df = pd.concat(dfs, axis=0)\n",
        "shapes_df_grouped = shapes_df.groupby(by=['nchans', 'ncols', 'nrows', 'type']).size().reset_index().sort_values(['type', 0], ascending=False)\n",
        "shapes_df_grouped"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "99cbb14f-ea81-1ae6-0d15-1734a52ff8c2"
      },
      "outputs": [],
      "source": [
        "\n",
        "\n",
        "shapes_df_grouped['size_with_type'] = shapes_df_grouped.apply(lambda row: '{}-{}-{}'.format(row.ncols, row.nrows, row.type), axis=1)\n",
        "shapes_df_grouped = shapes_df_grouped.set_index(shapes_df_grouped['size_with_type'].values)\n",
        "shapes_df_grouped['count'] = shapes_df_grouped[[0]]\n",
        "\n",
        "plt.figure(figsize=(10,8))\n",
        "#shapes_df_grouped['count'].plot.barh(figsize=(10,8))\n",
        "sns.barplot(x=\"count\", y=\"size_with_type\", data=shapes_df_grouped)\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "51c3daf0-5dec-4d53-d2b1-42452a64de6b"
      },
      "outputs": [],
      "source": [
        "def transform_image(img, rescaled_dim, to_gray=False):\n",
        "    resized = cv2.resize(img, (rescaled_dim, rescaled_dim), cv2.INTER_LINEAR)\n",
        "\n",
        "    if to_gray:\n",
        "        resized = cv2.cvtColor(resized, cv2.COLOR_RGB2GRAY).astype('float')\n",
        "    else:\n",
        "        resized = resized.astype('float')\n",
        "\n",
        "    normalized = cv2.normalize(resized, None, 0.0, 1.0, cv2.NORM_MINMAX)\n",
        "    timg = normalized.reshape(1, np.prod(normalized.shape))\n",
        "\n",
        "    return timg/np.linalg.norm(timg)\n",
        "\n",
        "rescaled_dim = 100\n",
        "\n",
        "all_images = []\n",
        "all_image_types = []\n",
        "\n",
        "for t in all_cervix_images['type'].unique():\n",
        "    all_images = all_images + images[t]\n",
        "    all_image_types = all_image_types + len(images[t])*[t]\n",
        "\n",
        "# - normalize each uint8 image to the value interval [0, 1] as float image\n",
        "# - rgb to gray\n",
        "# - downsample image to rescaled_dim X rescaled_dim\n",
        "# - L2 norm of each sample = 1\n",
        "gray_all_images_as_vecs = [transform_image(img, rescaled_dim) for img in all_images]\n",
        "\n",
        "gray_imgs_mat = np.array(gray_all_images_as_vecs).squeeze()\n",
        "all_image_types = np.array(all_image_types)\n",
        "gray_imgs_mat.shape, all_image_types.shape"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "9b477383-40ea-3e5d-3298-95fae59e892c"
      },
      "outputs": [],
      "source": [
        "\n",
        "\n",
        "from sklearn.manifold import TSNE\n",
        "tsne = TSNE(\n",
        "    n_components=3,\n",
        "    init='random', # pca\n",
        "    random_state=101,\n",
        "    method='barnes_hut',\n",
        "    n_iter=500,\n",
        "    verbose=2\n",
        ").fit_transform(gray_imgs_mat)\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "2465358f-6a87-ae88-ae15-a4789db87d38"
      },
      "outputs": [],
      "source": [
        "from sklearn import preprocessing\n",
        "\n",
        "trace1 = go.Scatter3d(\n",
        "    x=tsne[:,0],\n",
        "    y=tsne[:,1],\n",
        "    z=tsne[:,2],\n",
        "    mode='markers',\n",
        "    marker=dict(\n",
        "        sizemode='diameter',\n",
        "        color = preprocessing.LabelEncoder().fit_transform(all_image_types),\n",
        "        colorscale = 'Portland',\n",
        "        colorbar = dict(title = 'cervix types'),\n",
        "        line=dict(color='rgb(255, 255, 255)'),\n",
        "        opacity=0.9\n",
        "    )\n",
        ")\n",
        "\n",
        "data=[trace1]\n",
        "layout=dict(height=800, width=800, title='3D embedding of images')\n",
        "fig=dict(data=data, layout=layout)\n",
        "py.iplot(fig, filename='3DBubble')\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "ca0605c0-e976-348d-685e-6f093f3ecebc"
      },
      "outputs": [],
      "source": [
        "for t in all_cervix_images['type'].unique():\n",
        "    tsne_t = tsne[np.where(all_image_types == t), :][0]\n",
        "    plt.scatter(tsne_t[:, 0], tsne_t[:, 1])\n",
        "plt.legend(all_cervix_images['type'].unique())"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "00071c3d-ae8b-5e0a-4c42-3fcd3da2eb65"
      },
      "outputs": [],
      "source": [
        "from matplotlib.offsetbox import OffsetImage, AnnotationBbox\n",
        "def imscatter(x, y, images, ax=None, zoom=0.01):\n",
        "    ax = plt.gca()\n",
        "    images = [OffsetImage(image, zoom=zoom) for image in images]\n",
        "    artists = []\n",
        "    for x0, y0, im0 in zip(x, y, images):\n",
        "        ab = AnnotationBbox(im0, (x0, y0), xycoords='data', frameon=False)\n",
        "        artists.append(ax.add_artist(ab))\n",
        "    ax.update_datalim(np.column_stack([x, y]))\n",
        "    ax.autoscale()\n",
        "    #return artists\n",
        "\n",
        "nimgs = 60\n",
        "plt.figure(figsize=(10,8))\n",
        "imscatter(tsne[0:nimgs,0], tsne[0:nimgs,1], all_images[0:nimgs])\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "_cell_guid": "8a28e3c9-bdc1-4482-cb83-79f6b17f260f"
      },
      "source": [
        "\n",
        "Clustering of pairwise image distances\n",
        "\n",
        "To get a different view of how images relate to each other from a purely numerical point of view, lets now look at pairwise distances. For that we'll use scipy's pdist method.\n",
        "\n",
        "The yellow somewhat clustered area tells us there are a few images that have relatively high distance to all other images in the sample of our training images we read. On the left and top of the clustermap we find one of three colors for each row and column, this color indicates the type of cervix.\n",
        "\n",
        "    Type 1: Red\n",
        "    Type 2: Green\n",
        "    Type 3: Blue\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "71ed3b49-3564-d071-0713-c10f75e6d1b6"
      },
      "outputs": [],
      "source": [
        "pal = sns.color_palette(\"hls\", 3)\n",
        "sns.palplot(pal)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "588234f0-b0f1-6221-1f42-fa9310ceb670"
      },
      "outputs": [],
      "source": [
        "\n",
        "\n",
        "from scipy.spatial.distance import pdist, squareform\n",
        "\n",
        "sq_dists = squareform(pdist(gray_imgs_mat))\n",
        "\n",
        "all_image_types = list(all_image_types)\n",
        "\n",
        "d = {\n",
        "    'Type_1': pal[0],\n",
        "    'Type_2': pal[1],\n",
        "    'Type_3': pal[2]\n",
        "}\n",
        "\n",
        "# translate each sample to its color\n",
        "colors = list(map(lambda t: d[t], all_image_types))\n",
        "\n",
        "sns.clustermap(\n",
        "    sq_dists,\n",
        "    figsize=(12,12),\n",
        "    row_colors=colors, col_colors=colors,\n",
        "    cmap=plt.get_cmap('viridis')\n",
        ")\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "143ba258-04c0-80c2-1b02-bd379d282a3c"
      },
      "outputs": [],
      "source": [
        "mask = np.zeros_like(sq_dists, dtype=np.bool)\n",
        "mask[np.triu_indices_from(mask)] = True\n",
        "\n",
        "plt.figure(figsize=(12,12))\n",
        "sns.heatmap(sq_dists, cmap=plt.get_cmap('viridis'), square=True, mask=mask)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "_cell_guid": "e4c1a172-99fa-49dc-2aea-4417bc9547e6"
      },
      "source": [
        "\n",
        "Model Selection\n",
        "\n",
        "Now that we've established a basic idea about the data, let's do the most straightforward approach, where we take the resized color images and labels and train a, most likely quite heavily regularized, linear model like logistic regression on it.\n",
        "\n",
        "It is quite important to understand that we only have read a few training instances, 108, and have thousands of dimensions. To be able to cope with that we'll most likely end up using L1 regularization.\n",
        "\n",
        "For the multi-class problem we are faced with here, we'll use standard approach of OVR (one vs rest), meaning we will train three models where each of them is designed to distinguish class 1, 2 and 3 from the others respectively.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "2573e775-e0a9-70a0-5616-bbe7f748eacf"
      },
      "outputs": [],
      "source": [
        "from sklearn.preprocessing import LabelEncoder\n",
        "from sklearn.preprocessing import Normalizer\n",
        "y = LabelEncoder().fit_transform(all_image_types).reshape(-1)\n",
        "X = gray_imgs_mat # no need for normalizing, we already did this earlier Normalizer().fit_transform(gray_imgs_mat)\n",
        "X.shape, y.shape"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "dc9de723-ed5a-210c-fac9-c0d6cdda9820"
      },
      "outputs": [],
      "source": [
        "from sklearn.linear_model import LogisticRegression\n",
        "from sklearn.model_selection import GridSearchCV, train_test_split\n",
        "\n",
        "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.33, random_state=42)\n",
        "\n",
        "X_train.shape, X_test.shape, y_train.shape, y_test.shape"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "f559e124-4052-7fe3-f890-993cee846d15"
      },
      "outputs": [],
      "source": [
        "clf = LogisticRegression()\n",
        "grid = {\n",
        "    'C': [1e-9, 1e-6, 1e-3, 1e0],\n",
        "    'penalty': ['l1', 'l2']\n",
        "}\n",
        "cv = GridSearchCV(clf, grid, scoring='neg_log_loss', n_jobs=-1, verbose=1)\n",
        "cv.fit(X_train, y_train)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "59041a04-44a2-7630-c85a-1ae4478af8e4"
      },
      "outputs": [],
      "source": [
        "\n",
        "\n",
        "for i in range(1, len(cv.cv_results_['params'])+1):\n",
        "    rank = cv.cv_results_['rank_test_score'][i-1]\n",
        "    s = cv.cv_results_['mean_test_score'][i-1]\n",
        "    sd = cv.cv_results_['std_test_score'][i-1]\n",
        "    params = cv.cv_results_['params'][i-1]\n",
        "    print(\"{0}. Mean validation neg log loss: {1:.6f} (std: {2:.6f}) - {3}\".format(\n",
        "        rank,\n",
        "        s,\n",
        "        sd,\n",
        "        params\n",
        "    ))\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "0e77a54a-e0ec-12eb-8ec8-6803da921b75"
      },
      "outputs": [],
      "source": [
        "y_test_hat_p = cv.predict_proba(X_test)\n",
        "y_test_hat_p"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "b89b4ac8-d2ff-55ee-9cd9-3123c5cc52bd"
      },
      "outputs": [],
      "source": [
        "\n",
        "plt.figure(figsize=(12,8))\n",
        "sns.distplot(y_test_hat_p[:,0], color='red')\n",
        "sns.distplot(y_test_hat_p[:,1], color='blue')\n",
        "sns.distplot(y_test_hat_p[:,2], color='green')"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "7690ed96-284d-6ff3-b5a3-18aef10bb49e"
      },
      "outputs": [],
      "source": [
        "dfy = pd.DataFrame({'0': y_test_hat_p[:,0], '1': y_test_hat_p[:,1], '2': y_test_hat_p[:,2]})\n",
        "sns.pairplot(dfy)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "_cell_guid": "1ef590d6-4886-c49d-6f65-6233ee9d2767"
      },
      "outputs": [],
      "source": [
        "from sklearn.metrics import confusion_matrix\n",
        "\n",
        "y_test_hat = cv.predict(X_test)\n",
        "\n",
        "data = [\n",
        "    go.Heatmap(\n",
        "        z=confusion_matrix(y_test, y_test_hat),\n",
        "        x=[0, 1, 2],\n",
        "        y=[0, 1, 2],\n",
        "        colorscale='Viridis',\n",
        "        text = True ,\n",
        "        opacity = 1.0\n",
        "    )\n",
        "]\n",
        "\n",
        "layout = go.Layout(\n",
        "    title='Test Confusion matrix',\n",
        "    xaxis = dict(ticks='', nticks=36),\n",
        "    yaxis = dict(ticks='' ),\n",
        "    width = 900, height = 700,\n",
        "    \n",
        ")\n",
        "\n",
        "\n",
        "fig = go.Figure(data=data, layout=layout)\n",
        "py.iplot(fig, filename='labelled-heatmap')"
      ]
    }
  ],
  "metadata": {
    "_change_revision": 0,
    "_is_fork": false,
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "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.0"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}