{"cells":[{"metadata":{},"cell_type":"markdown","source":"# APTOS 2019 Blindness Detection - Diabetic Retinopathy"},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"# Miscelaneous.\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport time as t\nimport cv2\nimport warnings\nimport shutil\nimport os\n\n# Sklearn utils.\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.metrics import confusion_matrix\nfrom sklearn.utils import shuffle\n\n# Keras.\nfrom keras.applications.resnet50 import ResNet50\nfrom keras.applications import DenseNet121\nfrom keras.utils import Sequence\nfrom keras.layers import Dense, Dropout, Flatten, Input, ZeroPadding2D, GlobalAveragePooling2D\nfrom keras.models import Sequential\nfrom keras.layers import Conv2D, MaxPooling2D\nfrom keras.layers.normalization import BatchNormalization\nfrom keras import regularizers\nfrom keras.models import Model\nfrom keras.callbacks import EarlyStopping, ReduceLROnPlateau, ModelCheckpoint\nfrom keras import optimizers\nfrom keras.utils import to_categorical\nfrom keras.preprocessing.image import ImageDataGenerator\nfrom keras import backend as K\nimport tensorflow as tf","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Some useful variables"},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":true},"cell_type":"code","source":"PATH = '../input/'\nPATH_APTOS = PATH + 'aptos2019-blindness-detection/'\nPATH_2015 = PATH + 'retinopathy-train-2015/rescaled_train_896/'\nPATH_KERAS = PATH + 'keras-pretrained-models/'\ndensenet_weights_path = PATH + 'densenet121full/densenet121_weights_tf_dim_ordering_tf_kernels_notop.h5'\nPATH_AUGM = '/data_augm/'\nHEIGHT, WIDTH = 256, 256\nRANDOM_STATE = 974\nVERBOSE = True\nBATCH_SIZE = 32\nDATA_AUGM = True\nDATA_AUGM_FACTOR = 1\nwarnings.filterwarnings('ignore')\nprint(os.listdir(PATH_KERAS))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 1) Training"},{"metadata":{},"cell_type":"markdown","source":"### a) Getting the data"},{"metadata":{},"cell_type":"markdown","source":"**Loading APTOS dataset.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"df_aptos = pd.read_csv(PATH_APTOS + 'train.csv')\ndf_aptos.rename(index=str, columns={\"id_code\": \"id_code_aptos\"}, inplace=True)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"**Loading 2015 dataset. Since diagnosis 0 is already too prevalent in the aptos dataset, we will remove 80% of it.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"df_2015 = pd.read_csv(PATH_2015 + 'trainLabels.csv')\ndf_2015.rename(index=str, columns={\"image\": \"id_code_2015\", \"level\": \"diagnosis\"}, inplace=True)\ndf_2015 = df_2015.drop(df_2015[df_2015['diagnosis'] == 0].sample(frac=0.8).index)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"**Now we merge both dataset.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"df = df_aptos.append(df_2015)\ndf[\"data_augm\"] = np.nan\nif VERBOSE:\n    display(df.head())\n    print(df.shape)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### b) Analizing the data"},{"metadata":{"trusted":true},"cell_type":"code","source":"if VERBOSE:\n    plt.hist(df.diagnosis, bins=[-0.2, 0.2, 0.8, 1.2, 1.8, 2.2, 2.8, 3.2, 3.8, 4.2])\n    plt.xlabel(\"Severity of Diabetic Retinopathy\")\n    plt.ylabel(\"Count\")\n    plt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### c) Preprocessing the data"},{"metadata":{},"cell_type":"markdown","source":"**Our image preprocessing function removes black borders, corners and return resized image. **"},{"metadata":{"trusted":true},"cell_type":"code","source":"def imagePreprocessing(image, normalize=True):\n    # Cutting black border.\n    row, col = image.shape[0], image.shape[1]\n    center_row, center_col = row // 2, col // 2\n    x_left, x_right, y_top, y_bot = 0, col, 0, row\n    while image[center_row, x_left:x_left + 10].mean().max() <= 10  and x_left < col:\n        x_left += 1\n    while image[center_row, x_right - 10:x_right].mean().max() <= 10 and x_right > 0:\n        x_right -= 1\n    while image[y_top:y_top + 10, center_col].mean().max() <= 10 and y_top < row:\n        y_top += 1\n    while image[y_bot - 10:y_bot, center_col].mean().max() <= 10 and y_bot > 0:\n        y_bot -= 1\n    if y_top < y_bot and x_left < x_right and y_bot - y_top > 0.6 * row and x_right - x_left > 0.6 * col:\n        image = image[y_top:y_bot, x_left:x_right]\n\n    # Cutting to remove black corner.\n    row, col = image.shape[0], image.shape[1]\n    top_left_x, top_left_y = 0, 0\n    while image[0, top_left_x:top_left_x + 10].mean().max() <= 10  and top_left_x < col:\n        top_left_x += 1\n    while image[top_left_y:top_left_y + 10, 0].mean().max() <= 10  and top_left_y < row:\n        top_left_y += 1\n    if 2 * np.abs(top_left_y - top_left_x) / (top_left_y + top_left_x) > 0.85 :\n        crop_left_right = int(0.5 * top_left_x)\n        if crop_left_right < 0.3 * col:\n            image = image[:, crop_left_right:col - crop_left_right]\n    else:\n        crop = int(0.15 * (top_left_x + top_left_y) / 2)\n        if crop < 0.3 * col and crop < 0.3 * row:\n            image = image[crop:row - crop, crop:col - crop]\n    \n    # Resizing image.\n    image = cv2.resize(image, (WIDTH, HEIGHT))\n    \n    # Applying GaussianBlur.\n    blurred = cv2.blur(image, ksize=(int(WIDTH / 6), int(HEIGHT / 6)))\n    image = cv2.addWeighted(image, 4, blurred, -4, 128)\n    \n    try:\n        if normalize:\n            image = image / 255\n            image -= image.mean()\n            return image\n        else:\n            return image\n    except:\n        return np.zeros((WIDTH, HEIGHT, 3))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"**This function just open the image using a dataframe row.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"def openImage(row, train=True, resize=False):\n    image = None\n    if not train:\n        image = cv2.imread('{}test_images/{}.png'.format(PATH_APTOS, row.id_code))\n    elif not pd.isnull(row.data_augm):\n        image = cv2.imread(row.data_augm)\n    elif not pd.isnull(row.id_code_aptos):\n        image = cv2.imread('{}train_images/{}.png'.format(PATH_APTOS, row.id_code_aptos))\n    elif not pd.isnull(row.id_code_2015):\n        image = cv2.imread('{}rescaled_train_896/{}.png'.format(PATH_2015, row.id_code_2015))\n    else:\n        print(\"[Error] Could not open the image. Log: {}\".format(row))\n    if resize and not image is None:\n        return cv2.resize(image, (WIDTH, HEIGHT))\n    return image","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"**Before training the model, let's display few retinas to check if the preprocessing is working correctly.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"if VERBOSE:\n    nb_row = 4\n    nb_col = 6\n    nb = 1\n    plt.figure(figsize=(25, 15))\n    for row in df.itertuples():\n        if nb > nb_col * nb_row:\n            break\n        plt.subplot(nb_row, nb_col, nb)\n        plt.imshow(cv2.cvtColor(imagePreprocessing(openImage(row), normalize=False), cv2.COLOR_BGR2RGB))\n        plt.title('Diagnosed {}'.format(row.diagnosis))\n        nb += 1\n    plt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"**RetinaGenerator is the generator for both training and predicting.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"class RetinaGenerator():\n    def __init__(self, image_df, batch_size, train=True):\n        self.image_df = image_df\n        self.batch_size = batch_size\n        self.train = train\n        self.step_per_epoch = len(self.image_df) // self.batch_size\n        self.step_per_pred = 1 + (len(self.image_df) - 1) // self.batch_size\n    def getGenerator(self):\n        while True:\n            for idx in range(self.step_per_epoch):\n                batch_x = np.array([imagePreprocessing(openImage(row)) for row in self.image_df[idx * self.batch_size:(idx + 1) * self.batch_size].itertuples()])\n                batch_y_cat= to_categorical([row.diagnosis for row in self.image_df[idx * self.batch_size:(idx + 1) * self.batch_size].itertuples()], num_classes=5)\n                batch_y = np.empty(batch_y_cat.shape, dtype=batch_y_cat.dtype)\n                batch_y[:, 4] = batch_y_cat[:, 4]\n                for i in range(3, -1, -1):\n                    batch_y[:, i] = np.logical_or(batch_y_cat[:, i], batch_y[:, i+1])\n                yield batch_x, batch_y\n    def getInputGenerator(self):\n        for idx in range(self.step_per_pred + 1):\n            yield np.array([imagePreprocessing(openImage(row, self.train)) for row in self.image_df[idx * self.batch_size:min((idx + 1) * self.batch_size, len(self.image_df))].itertuples()])","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"**Before spliting our dataframe, we will run a data augmentation. A data augmentation of 1 will allign to the class with the most sample.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"if DATA_AUGM:\n    shutil.rmtree(PATH_AUGM, ignore_errors=True, onerror=None)\n    shutil.os.mkdir(PATH_AUGM)\n    df = shuffle(df, random_state=RANDOM_STATE)\n    max_size = int(df.diagnosis.value_counts().max() * DATA_AUGM_FACTOR)\n    for diag in range(5):\n        shutil.rmtree('{}diag_{}'.format(PATH_AUGM, diag), ignore_errors=True, onerror=None)\n        shutil.os.mkdir('{}diag_{}'.format(PATH_AUGM, diag))\n        diag_df = df[df.diagnosis == diag]\n        size = len(diag_df)\n        to_create = max_size - size\n        augm_per_img = max(to_create // size, 1)\n        while to_create > 0:\n            for row in diag_df.itertuples():\n                if to_create < 0:\n                    break\n                image = np.expand_dims(cv2.cvtColor(openImage(row, resize=True), cv2.COLOR_BGR2RGB), 0)\n                data_generator = ImageDataGenerator(rotation_range=360, vertical_flip=True, horizontal_flip=True, zoom_range=0.1)\n                data_generator.fit(image)\n                id_code = row.id_code_aptos if not pd.isnull(row.id_code_aptos) else row.id_code_2015\n                for x, val in zip(data_generator.flow(image, save_to_dir='{}diag_{}'.format(PATH_AUGM, diag), \n                                                      save_prefix=id_code, save_format='png'), \n                                  range(augm_per_img - 1)):\n                    pass\n                to_create -= augm_per_img\n    for diag in range(5):\n        images = np.array(os.listdir(\"{}diag_{}\".format(PATH_AUGM, diag)))\n        for image in images:#diagnosis\tid_code_2015\tid_code_aptos\tdata_augm\n            df = df.append(pd.DataFrame([[diag, np.nan, np.nan, \"{}diag_{}/{}\".format(PATH_AUGM, diag, image)]], columns=df.columns))\n    df = shuffle(df, random_state=RANDOM_STATE)\n    display(df.head())","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"**Before training the model, let's display few retinas to check if the data augmenation is working correctly.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"if VERBOSE:\n    nb_row = 4\n    nb_col = 6\n    nb = 1\n    plt.figure(figsize=(25, 15))\n    for row in df.itertuples():\n        if nb > nb_col * nb_row:\n            break\n        plt.subplot(nb_row, nb_col, nb)\n        plt.imshow(cv2.cvtColor(imagePreprocessing(openImage(row), normalize=False), cv2.COLOR_BGR2RGB))\n        plt.title('Diagnosed {}'.format(row.diagnosis))\n        nb += 1\n    plt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"**Displaying the histogram of our dataset to check if the data augmentation worked.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"if VERBOSE:\n    plt.hist(df.diagnosis, bins=[-0.2, 0.2, 0.8, 1.2, 1.8, 2.2, 2.8, 3.2, 3.8, 4.2])\n    plt.xlabel(\"Severity of Diabetic Retinopathy\")\n    plt.ylabel(\"Count\")\n    plt.show()\n    print(df.shape)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### d) Training the model"},{"metadata":{},"cell_type":"markdown","source":"**Here we simply split the main dataframe into three for training, validating and testing.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"train_df, val_df = train_test_split(df, test_size=0.2, random_state=RANDOM_STATE)\ntrain_df, test_df = train_test_split(train_df, test_size=0.2, random_state=RANDOM_STATE)\nif VERBOSE:\n    print(train_df.shape, val_df.shape, test_df.shape)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"LEARNING_RATE = 0.005\nEPOCHS = 16\nPATIENCE = 16\nLR_PATIENCE = 8\nVALIDATION_GENERATOR = RetinaGenerator(val_df, BATCH_SIZE)\nTRAINING_GENERATOR = RetinaGenerator(train_df, BATCH_SIZE)\nALL_GENERATOR = RetinaGenerator(df, BATCH_SIZE)\nSTEPS_PER_EPOCH = TRAINING_GENERATOR.step_per_epoch\nVALIDATION_STEPS = VALIDATION_GENERATOR.step_per_epoch","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"visible = Input(shape=(HEIGHT, WIDTH, 3))\ndensenet = DenseNet121(include_top=False,\n                 weights=None,\n                 input_tensor=visible)\ndensenet.load_weights(densenet_weights_path)\n\nmodel = Sequential()\nmodel.add(densenet)\nmodel.add(GlobalAveragePooling2D())\nmodel.add(Dropout(0.5))\nmodel.add(Dense(5, activation='sigmoid'))\n    \nmodel.compile(\n    loss='binary_crossentropy',\n    optimizer=optimizers.Adam(lr=0.00005),\n    metrics=['accuracy']\n)\nprint(model.summary())\n\nopt = optimizers.SGD(lr=LEARNING_RATE, decay=1e-6, momentum=0.9, nesterov=True)\n\nCallbacks=[EarlyStopping(patience=PATIENCE, restore_best_weights=True), \n           ReduceLROnPlateau(monitor='val_loss', factor=0.5, patience=LR_PATIENCE, verbose=VERBOSE)]\nwith tf.device('/device:GPU:0'):\n    H = model.fit_generator(generator=TRAINING_GENERATOR.getGenerator(),\n                            validation_data=VALIDATION_GENERATOR.getGenerator(),\n                            steps_per_epoch=STEPS_PER_EPOCH,\n                            validation_steps=VALIDATION_STEPS,\n                            epochs=EPOCHS,\n                            callbacks=Callbacks,\n                            verbose=VERBOSE)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### e) Testing the model"},{"metadata":{},"cell_type":"markdown","source":"**Some verification of our model predictions using unseen images.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"from sklearn.metrics import accuracy_score, cohen_kappa_score\n\nif VERBOSE:\n    TEST_GENERATOR = RetinaGenerator(test_df, BATCH_SIZE)\n    STEPS = TEST_GENERATOR.step_per_pred\n    Y_pred = (model.predict_generator(TEST_GENERATOR.getInputGenerator(), steps=STEPS) > 0.5).sum(axis=1) - 1\n    Y_test = np.array(test_df.diagnosis)\n    print(\"Average absolute distance is: {:.2f}\".format(np.abs(Y_pred - Y_test).mean()))\n    display(confusion_matrix(Y_test, Y_pred))\n    print(\"Accuracy Score:\" + str(accuracy_score(Y_test, Y_pred)))\n    print(\"Cohen Kappa Score:\" + str(cohen_kappa_score(Y_test, Y_pred, weights='quadratic')))\n    ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Fit regressor"},{"metadata":{"trusted":true},"cell_type":"code","source":"STEPS = ALL_GENERATOR.step_per_pred\nY_pred = model.predict_generator(ALL_GENERATOR.getInputGenerator(), steps=STEPS)\nY_true = to_categorical(np.array(df.diagnosis), 5)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"Y_true.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"post_model = Sequential()\npost_model.add(Dense(5, input_dim=5, activation='relu'))\npost_model.add(Dense(5, activation='softmax'))\npost_model.compile(loss='categorical_crossentropy', optimizer='adam', metrics=['accuracy'])\n    \nopt = optimizers.SGD(lr=LEARNING_RATE, decay=1e-6, momentum=0.9, nesterov=True)\n\nCallbacks=[EarlyStopping(patience=PATIENCE, restore_best_weights=True), \n           ReduceLROnPlateau(monitor='val_loss', factor=0.5, patience=LR_PATIENCE, verbose=VERBOSE)]\nwith tf.device('/device:GPU:0'):\n    H = post_model.fit(Y_pred, Y_true, batch_size=10, epochs=10, verbose=VERBOSE )\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"After many tests, 10 seems to be the best number of epochs"},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"**Once our model is trained and tested, we can delete the data augmentation folder.**"},{"metadata":{"trusted":true},"cell_type":"code","source":"if DATA_AUGM:\n    shutil.rmtree(PATH_AUGM, ignore_errors=True, onerror=None)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## 2) Predicting"},{"metadata":{},"cell_type":"markdown","source":"### a) Getting the data"},{"metadata":{"trusted":true},"cell_type":"code","source":"sub_df = pd.read_csv(PATH_APTOS + 'test.csv')\nif VERBOSE:\n    print(sub_df.shape)\n    display(sub_df.head())","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### b) Predicting"},{"metadata":{"trusted":true},"cell_type":"code","source":"SUB_GENERATOR = RetinaGenerator(sub_df, BATCH_SIZE, False)\nSTEPS = SUB_GENERATOR.step_per_pred\nY_sub_preview = model.predict_generator(SUB_GENERATOR.getInputGenerator(), steps=STEPS)\nY_sub_cat = post_model.predict(Y_sub_preview)\nY_sub = np.argmax(Y_sub_cat, axis=1)\n#Y_sub = np.argmax(model.predict_generator(SUB_GENERATOR.getInputGenerator(), steps=STEPS), axis=1)\nif VERBOSE:\n    print(Y_sub.shape)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"Y_sub","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub_df['diagnosis'] = Y_sub\nif VERBOSE:\n    print(sub_df.shape)\n    display(sub_df)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### c) Testing predictions"},{"metadata":{"trusted":true},"cell_type":"code","source":"if VERBOSE:\n    nb_row = 4\n    nb_col = 6\n    nb = 1\n    plt.figure(figsize=(25, 15))\n    for row in sub_df.itertuples():\n        if nb > nb_col * nb_row:\n            break\n        plt.subplot(nb_row, nb_col, nb)\n        plt.imshow(cv2.cvtColor(imagePreprocessing(openImage(row, False), normalize=False), cv2.COLOR_BGR2RGB))\n        plt.title('Diagnosed {}'.format(row.diagnosis))\n        nb += 1\n    plt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### c) Writing predictions"},{"metadata":{"trusted":true},"cell_type":"code","source":"sub_df.to_csv('submission.csv', index=False)","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}