{"cells":[{"metadata":{"_uuid":"b8d1c1593b4892a11e30c5362c7557657492c205","_cell_guid":"67e8fd4b-e795-4385-935f-d3db055e71ac"},"cell_type":"markdown","source":"# Overview\nThe notebook aims to organize the data and hack Keras so that we can train a model in a fairly simple way. The aim here is to get a model working that can reliably segment the images into objects and then we can make a model that handles grouping the objects into categories based on the labels. As you will see the Keras requires a fair bit of hackery to get it to load images from a dataframe and then get it to read the label images correctly (uint16 isn't supported well). Once that is done, training a U-Net model is really easy."},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","collapsed":true,"_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"%matplotlib inline\nimport os\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nfrom glob import glob\nimport matplotlib.pyplot as plt\nUSE_CV2 = True\nif USE_CV2:\n    from cv2 import imread # opencv is much faster, but less accurate\n    MIN_OBJ_VAL = 0\nelse:\n    from skimage.io import imread\n    MIN_OBJ_VAL = 1000\n\nfrom skimage.segmentation import mark_boundaries\nDATA_DIR = os.path.join('..', 'input')","execution_count":13,"outputs":[]},{"metadata":{"scrolled":true,"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"all_paths = pd.DataFrame(dict(path = glob(os.path.join(DATA_DIR, '*', '*.*p*g'))))\nall_paths['split'] = all_paths['path'].map(lambda x: x.split('/')[-2].split('_')[0])\nall_paths['group'] = all_paths['path'].map(lambda x: x.split('/')[-2].split('_')[-1])\nall_paths['group'] = all_paths['group'].map(lambda x: 'color' if x == 'test' else x)\nall_paths['id'] = all_paths['path'].map(lambda x: '_'.join(os.path.splitext(os.path.basename(x))[0].split('_')[:4]))\nall_paths.sample(5)","execution_count":14,"outputs":[]},{"metadata":{"_uuid":"6892e98b6dcec247664d1b7a521ec25beeb84619","_cell_guid":"cb0565b2-ff44-4e25-b401-8794472552d7","trusted":true},"cell_type":"code","source":"group_df = all_paths.pivot_table(values = 'path', columns = 'group', aggfunc = 'first', index = ['id', 'split']).reset_index()\ngroup_df.sample(5)","execution_count":15,"outputs":[]},{"metadata":{"_uuid":"7019801430e6b9ebe234e07544d2d4ce47d31015","_cell_guid":"567d4f91-cb7d-4b04-8798-2d1e1a3fcc14"},"cell_type":"markdown","source":"# Explore the training set\nHere we can show the training data image by image to see what exactly we are supposed to detect with the model"},{"metadata":{"_uuid":"d5c261be25e697941928a396e4d6f5a4b1ca4395","_cell_guid":"2a63d977-02c8-43a9-8f8c-02b3499309c4","trusted":true},"cell_type":"code","source":"train_df = group_df.query('split==\"train\"')\nprint(train_df.shape[0], 'rows')\nsample_rows = 6\nfig, m_axs = plt.subplots(sample_rows, 3, figsize = (20, 6*sample_rows))\n[c_ax.axis('off') for c_ax in m_axs.flatten()]\nfor (ax1, ax2, ax3), (_, c_row) in zip(m_axs, train_df.sample(sample_rows).iterrows()):\n    c_img = imread(c_row['color'])\n    l_img = imread(c_row['label'])\n    if l_img.ndim==3: l_img = l_img[:,:,0]\n    ax1.imshow(c_img)\n    ax1.set_title('Color')\n    # make the labels nicer\n    nice_limg = np.zeros(l_img.shape, dtype = np.uint8)\n    for new_idx, old_idx in enumerate(np.unique(l_img[l_img>MIN_OBJ_VAL]), 1):\n        nice_limg[l_img==old_idx]=new_idx\n    ax2.imshow(nice_limg, cmap = 'nipy_spectral')\n    ax2.set_title('Labels')\n    xd, yd = np.where(l_img>MIN_OBJ_VAL)\n    bound_img = mark_boundaries(image = c_img, label_img = l_img, color = (1,0,0), background_label = 255, mode = 'thick')\n    ax3.imshow(bound_img[xd.min():xd.max(), yd.min():yd.max(),:])\n    ax3.set_title('Cropped Overlay')","execution_count":16,"outputs":[]},{"metadata":{"_uuid":"c13081e6618b920cdc36285fde271983d2c3a218","collapsed":true,"_cell_guid":"41ccd780-75e2-4858-8dd7-26ab54c0290e","trusted":true},"cell_type":"code","source":"from sklearn.model_selection import train_test_split\ntrain_split_df, valid_split_df = train_test_split(train_df, random_state = 2018, test_size = 0.25)","execution_count":17,"outputs":[]},{"metadata":{"_uuid":"2c17e0ab0d9347cdf6af56b97aae9047be40d766","_cell_guid":"ddb1c113-1c53-4bc3-948e-e9d1edadc082","trusted":true,"collapsed":true},"cell_type":"code","source":"from keras.preprocessing.image import ImageDataGenerator\nfrom keras.applications.resnet50 import preprocess_input\nIMG_SIZE = (384, 384) # slightly smaller than vgg16 normally expects\nimg_gen_args = dict(samplewise_center=False, \n                              samplewise_std_normalization=False, \n                              horizontal_flip = True, \n                              vertical_flip = False, \n                              height_shift_range = 0.1, \n                              width_shift_range = 0.1, \n                              rotation_range = 3, \n                              shear_range = 0.01,\n                              fill_mode = 'nearest',\n                              zoom_range = 0.05)\n\nrgb_gen = ImageDataGenerator(preprocessing_function = preprocess_input, **img_gen_args)\nlab_gen = ImageDataGenerator(**img_gen_args)","execution_count":18,"outputs":[]},{"metadata":{"_uuid":"0ce7b0154f7d68add10776c10e9053b83b3e255e","collapsed":true,"_cell_guid":"4b67ae84-2f58-425f-9296-dcc7d2d62c9e","trusted":true},"cell_type":"code","source":"def flow_from_dataframe(img_data_gen, in_df, path_col, y_col, seed = None, **dflow_args):\n    base_dir = os.path.dirname(in_df[path_col].values[0])\n    print('## Ignore next message from keras, values are replaced anyways: seed: {}'.format(seed))\n    df_gen = img_data_gen.flow_from_directory(base_dir, \n                                     class_mode = 'sparse',\n                                              seed = seed,\n                                    **dflow_args)\n    df_gen.filenames = in_df[path_col].values\n    df_gen.classes = np.stack(in_df[y_col].values)\n    df_gen.samples = in_df.shape[0]\n    df_gen.n = in_df.shape[0]\n    df_gen._set_index_array()\n    df_gen.directory = '' # since we have the full path\n    print('Reinserting dataframe: {} images'.format(in_df.shape[0]))\n    return df_gen","execution_count":19,"outputs":[]},{"metadata":{"_uuid":"f8d2dfeb980142e54ba517815984bb5570ef91a2","_cell_guid":"28767464-3a4c-4110-a5a2-f0ecbbcf8f51"},"cell_type":"markdown","source":"## Replace PIL with scikit-image \nThis lets us handle the 16bit numbers well in the instanceIds image. This is incredibly, incredibly hacky, please do not use this code outside of this kernel."},{"metadata":{"_uuid":"5b435e6727af094fb76eb157a67950178fc54485","_cell_guid":"5951c0bc-1df1-4a8a-8e87-c076d88f7e12","trusted":true},"cell_type":"code","source":"import keras.preprocessing.image as KPImage\nfrom PIL import Image\nclass pil_image_awesome():\n    @staticmethod\n    def open(in_path):\n        if 'instanceIds' in in_path:\n            # we only want to keep the positive labels not the background\n            in_img = imread(in_path)\n            if in_img.ndim==3:\n                in_img = in_img[:,:,0]\n            return Image.fromarray((in_img>MIN_OBJ_VAL).astype(np.float32))\n        else:\n            return Image.open(in_path)\n    fromarray = Image.fromarray\nKPImage.pil_image = pil_image_awesome","execution_count":20,"outputs":[]},{"metadata":{"_uuid":"75d08811cbb12786d7cdb8f82ac3b8aa95867e41","_cell_guid":"0e9ec0f2-fe4b-4a6f-a9a3-118b4b4f8112"},"cell_type":"markdown","source":"# Create the generators\nWe want to generate parallel streams of images and labels"},{"metadata":{"_uuid":"bb9b79c1b61f3f7948e6fb6c50ae859df9f10c3b","collapsed":true,"_cell_guid":"6c287bc5-9ec0-4291-8747-4d0d80b51dd5","trusted":true},"cell_type":"code","source":"from skimage.filters.rank import maximum\nfrom scipy.ndimage import zoom\ndef lab_read_func(in_path):\n    bin_img = (imread(in_path)>1000).astype(np.uint8)\n    x_dim, y_dim = bin_img.shape\n    max_label_img = maximum(bin_img, np.ones((x_dim//IMG_SIZE[0], y_dim//IMG_SIZE[1])))\n    return np.expand_dims(zoom(max_label_img, (IMG_SIZE[0]/x_dim, IMG_SIZE[1]/y_dim), order = 3), -1)\n\n\ndef train_and_lab_gen_func(in_df, batch_size = 8, seed = None):\n    if seed is None:\n        seed = np.random.choice(range(1000))\n    train_rgb_gen = flow_from_dataframe(rgb_gen, in_df, \n                             path_col = 'color',\n                            y_col = 'id', \n                            target_size = IMG_SIZE,\n                             color_mode = 'rgb',\n                            batch_size = batch_size,\n                                   seed = seed)\n    train_lab_gen = flow_from_dataframe(lab_gen, in_df, \n                             path_col = 'label',\n                            y_col = 'id', \n                            target_size = IMG_SIZE,\n                             color_mode = 'grayscale',\n                            batch_size = batch_size,\n                                   seed = seed)\n    for (x, _), (y, _) in zip(train_rgb_gen, train_lab_gen):\n        yield x, y\n    \ntrain_and_lab_gen = train_and_lab_gen_func(train_split_df, batch_size = 32)\nvalid_and_lab_gen = train_and_lab_gen_func(valid_split_df, batch_size = 32)","execution_count":21,"outputs":[]},{"metadata":{"_uuid":"abbef517deaa6fed4d2988e6306dc0e388a19851","_cell_guid":"7d2e62b7-de3e-4fd6-9fe7-bde1e28f7c19","trusted":true},"cell_type":"code","source":"(rgb_batch, lab_batch) = next(valid_and_lab_gen)\n\nsample_rows = 4\nfig, m_axs = plt.subplots(sample_rows, 3, figsize = (20, 6*sample_rows))\n[c_ax.axis('off') for c_ax in m_axs.flatten()]\nfor (ax1, ax2, ax3), rgb_img, lab_img in zip(m_axs, rgb_batch, lab_batch):\n    # undoing the vgg correction is tedious\n    r_rgb_img = np.clip(rgb_img+110, 0, 255).astype(np.uint8)\n    ax1.imshow(r_rgb_img)\n    ax1.set_title('Color')\n    ax2.imshow(lab_img[:,:,0], cmap = 'nipy_spectral')\n    ax2.set_title('Labels')\n    if lab_img.max()>0.1:\n        xd, yd = np.where(lab_img[:,:,0]>0)\n        bound_img = mark_boundaries(image = r_rgb_img, label_img = lab_img[:,:,0], \n                                    color = (1,0,0), background_label = 255, mode = 'thick')\n        ax3.imshow(bound_img[xd.min():xd.max(), yd.min():yd.max(),:])\n        ax3.set_title('Cropped Overlay')","execution_count":22,"outputs":[]},{"metadata":{"_uuid":"33b097924b25db4b5863e0abce591c03ccddae0a","_cell_guid":"cd584b6c-c185-4519-96b6-047fa587386b","trusted":true},"cell_type":"code","source":"from keras.models import Model, load_model\nfrom keras.layers import Input, BatchNormalization, Dropout, Lambda\nfrom keras.layers.convolutional import Conv2D, Conv2DTranspose\nfrom keras.layers.pooling import MaxPooling2D\nfrom keras.layers.merge import concatenate\n\n# Build U-Net model\ninputs = Input(IMG_SIZE+(3,))\ns = BatchNormalization()(inputs) # we can learn the normalization step\ns = Dropout(0.5)(s)\n\nc1 = Conv2D(8, (3, 3), activation='relu', padding='same') (s)\nc1 = Conv2D(8, (3, 3), activation='relu', padding='same') (c1)\np1 = MaxPooling2D((2, 2)) (c1)\n\nc2 = Conv2D(16, (3, 3), activation='relu', padding='same') (p1)\nc2 = Conv2D(16, (3, 3), activation='relu', padding='same') (c2)\np2 = MaxPooling2D((2, 2)) (c2)\n\nc3 = Conv2D(32, (3, 3), activation='relu', padding='same') (p2)\nc3 = Conv2D(32, (3, 3), activation='relu', padding='same') (c3)\np3 = MaxPooling2D((2, 2)) (c3)\n\nc4 = Conv2D(64, (3, 3), activation='relu', padding='same') (p3)\nc4 = Conv2D(64, (3, 3), activation='relu', padding='same') (c4)\np4 = MaxPooling2D(pool_size=(2, 2)) (c4)\n\nc5 = Conv2D(128, (3, 3), activation='relu', padding='same') (p4)\nc5 = Conv2D(128, (3, 3), activation='relu', padding='same') (c5)\n\nu6 = Conv2DTranspose(64, (2, 2), strides=(2, 2), padding='same') (c5)\nu6 = concatenate([u6, c4])\nc6 = Conv2D(64, (3, 3), activation='relu', padding='same') (u6)\nc6 = Conv2D(64, (3, 3), activation='relu', padding='same') (c6)\n\nu7 = Conv2DTranspose(32, (2, 2), strides=(2, 2), padding='same') (c6)\nu7 = concatenate([u7, c3])\nc7 = Conv2D(32, (3, 3), activation='relu', padding='same') (u7)\nc7 = Conv2D(32, (3, 3), activation='relu', padding='same') (c7)\n\nu8 = Conv2DTranspose(16, (2, 2), strides=(2, 2), padding='same') (c7)\nu8 = concatenate([u8, c2])\nc8 = Conv2D(16, (3, 3), activation='relu', padding='same') (u8)\nc8 = Conv2D(16, (3, 3), activation='relu', padding='same') (c8)\n\nu9 = Conv2DTranspose(8, (2, 2), strides=(2, 2), padding='same') (c8)\nu9 = concatenate([u9, c1], axis=3)\nc9 = Conv2D(8, (3, 3), activation='relu', padding='same') (u9)\nc9 = Conv2D(8, (3, 3), activation='relu', padding='same') (c9)\n\noutputs = Conv2D(1, (1, 1), activation='sigmoid') (c9)\n\nmodel = Model(inputs=[inputs], outputs=[outputs])\nmodel.summary()","execution_count":23,"outputs":[]},{"metadata":{"_uuid":"99c84b4c88ce8a6319bd5ed8830f1d9d864f5035","collapsed":true,"_cell_guid":"1016dae1-b7a7-439e-9d5d-6ee661a4ff6e","trusted":true},"cell_type":"code","source":"import tensorflow as tf\nfrom keras import backend as K\n# Define IoU metric\ndef mean_iou(y_true, y_pred):\n    prec = []\n    for t in np.arange(0.5, 1.0, 0.05):\n        y_pred_ = tf.to_int32(y_pred > t)\n        score, up_opt = tf.metrics.mean_iou(y_true, y_pred_, 2)\n        K.get_session().run(tf.local_variables_initializer())\n        with tf.control_dependencies([up_opt]):\n            score = tf.identity(score)\n        prec.append(score)\n    return K.mean(K.stack(prec))\n\nsmooth = 1.\ndef dice_coef(y_true, y_pred):\n    y_true_f = K.flatten(y_true)\n    y_pred_f = K.flatten(y_pred)\n    intersection = K.sum(y_true_f * y_pred_f)\n    return (2. * intersection + smooth) / (K.sum(y_true_f) + K.sum(y_pred_f) + smooth)\n\ndef dice_coef_loss(y_true, y_pred):\n    return -dice_coef(y_true, y_pred)\n\nmodel.compile(optimizer = 'adam', \n                   loss = dice_coef_loss, \n                   metrics = [dice_coef, 'binary_accuracy', 'mse'])","execution_count":24,"outputs":[]},{"metadata":{"_uuid":"4695ec7bd66e624aee98971a92c97268a5c68b0f","_cell_guid":"fdf99eeb-1bbb-4e19-a106-238a38261edf","trusted":true},"cell_type":"code","source":"from keras.callbacks import ModelCheckpoint, LearningRateScheduler, EarlyStopping, ReduceLROnPlateau\nweight_path=\"{}_weights.best.hdf5\".format('unet')\n\ncheckpoint = ModelCheckpoint(weight_path, monitor='val_loss', verbose=1, \n                             save_best_only=True, mode='min', save_weights_only = True)\n\nreduceLROnPlat = ReduceLROnPlateau(monitor='loss', factor=0.8, patience=10, verbose=1, mode='auto', epsilon=0.0001, cooldown=5, min_lr=0.0001)\nearly = EarlyStopping(monitor=\"val_loss\", \n                      mode=\"min\", \n                      patience=5) # probably needs to be more patient, but kaggle time is limited\ncallbacks_list = [checkpoint, early, reduceLROnPlat]","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"631939949c70b3674b7b093b558fc0719cce0190","_cell_guid":"775fc5c8-b93a-48fa-8711-368f7880985c","trusted":true},"cell_type":"code","source":"# reset the generators so they all have different seeds when multiprocessing lets loose\nbatch_size = 16\ntrain_and_lab_gen = train_and_lab_gen_func(train_split_df, batch_size = batch_size)\nvalid_and_lab_gen = train_and_lab_gen_func(valid_split_df, batch_size = batch_size)\nmodel.fit_generator(train_and_lab_gen, \n                    steps_per_epoch = 2048//batch_size,\n                    validation_data = valid_and_lab_gen,\n                    validation_steps = 256//batch_size,\n                    epochs = 3, \n                    workers = 4,\n                    use_multiprocessing = True,\n                    callbacks = callbacks_list)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"b8e001a984e00b9aa70d9a418cb43e8913e4c8ea","_cell_guid":"7c4caaf8-3c31-4e0e-8298-8d9cc3f44c12"},"cell_type":"markdown","source":"# Showing the results\nHere we can preview the output of the model on a few examples"},{"metadata":{"_uuid":"57800bbeb10c674020f31bc788518b89551626d9","collapsed":true,"_cell_guid":"40f63c1e-4067-4702-878e-2f8c0474c401","trusted":true},"cell_type":"code","source":"(rgb_batch, lab_batch) = next(valid_and_lab_gen)\nsample_rows = 5\nfig, m_axs = plt.subplots(sample_rows, 5, figsize = (20, 6*sample_rows))\n[c_ax.axis('off') for c_ax in m_axs.flatten()]\nfor (ax1, ax2, ax2_pred, ax3, ax3_pred), rgb_img, lab_img in zip(m_axs, rgb_batch, lab_batch):\n    # undoing the vgg correction is tedious\n    r_rgb_img = np.clip(rgb_img+110, 0, 255).astype(np.uint8)\n    lab_pred = model.predict(np.expand_dims(rgb_img, 0))[0]\n    \n    ax1.imshow(r_rgb_img)\n    ax1.set_title('Color')\n    ax2.imshow(lab_img[:,:,0], cmap = 'jet')\n    ax2.set_title('Labels')\n    ax2_pred.imshow(lab_pred[:,:,0], cmap = 'jet')\n    ax2_pred.set_title('Pred Labels')\n    if lab_img.max()>0.1:\n        xd, yd = np.where(lab_img[:,:,0]>0)\n        bound_img = mark_boundaries(image = r_rgb_img, label_img = lab_img[:,:,0], \n                                    color = (1,0,0), background_label = 255, mode = 'thick')\n        ax3.imshow(bound_img[xd.min():xd.max(), yd.min():yd.max(),:])\n        ax3.set_title('Cropped Overlay')\n        ax3_pred.imshow(lab_pred[xd.min():xd.max(), yd.min():yd.max(),0], cmap = 'jet')\n        ax3_pred.set_title('Cropped Overlay')\nfig.savefig('trained_model.png')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"e49199a14b41cff1cc8ee7fc8908219be633ac9e","collapsed":true,"_cell_guid":"b305f765-d263-4792-85f8-0192a8df2a35","trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.6.4","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"}},"nbformat":4,"nbformat_minor":1}