{"cells":[{"metadata":{"trusted":true},"cell_type":"code","source":"import numpy  as np \nimport pandas as pd \nimport os\nimport cv2 \nimport gc\nimport psutil\nimport matplotlib.pyplot as plt\nfrom sklearn.metrics import cohen_kappa_score, confusion_matrix\n\nimport keras\nfrom keras.preprocessing import image\nfrom keras.applications  import DenseNet121\nfrom keras.activations   import softmax, relu\nfrom keras.optimizers    import Adam\nfrom keras.callbacks     import Callback, ModelCheckpoint, EarlyStopping, ReduceLROnPlateau\nfrom keras.models        import Sequential, Model\nfrom keras.layers        import Conv2D, MaxPooling2D, GlobalAveragePooling2D, Input\nfrom keras.layers        import Dropout, Flatten, Dense, BatchNormalization\n\n# Global constants\nIMG_DIM       = 256\nBATCH_SIZE    = 32\nCHANNELS      = 3\nNUM_CLASSES   = 5\n\nprint(os.listdir(\".\"))\nprint(os.listdir(\"../\"))\nprint(os.listdir(\"../input/\"))\nprint(os.listdir(\"../input/aptos2019-blindness-detection\"))\nprint(os.listdir(\"../input/densenetmulti\"))\n\nINPUT_FOLDER = '../input/aptos2019-blindness-detection/'","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Preprocessing Functions\nWe need to crop and preprocess all the test images into the working directory"},{"metadata":{"trusted":true},"cell_type":"code","source":"\ndef crop(gray, img, percent_smaller):\n    \n    thresh = 8\n    \n    top    = 0\n    left   = 0\n    bottom = gray.shape[0] - 1\n    right  = gray.shape[1] - 1\n    \n    # work in from the top and bottom along the middle collumn\n    middleCol = gray[:, int(gray.shape[1]/2)] > thresh\n    while middleCol[top] == 0:\n        top += 1\n    while middleCol[bottom] == 0:\n        bottom -= 1\n        \n    # work in from the sides along the middle row\n    middleRow = gray[int(gray.shape[0]/2)] > thresh\n    while middleRow[left] == 0:\n        left += 1\n    while middleRow[right] == 0:\n        right -= 1\n        \n    height = bottom - top\n    width  = right - left\n    \n    bottom -= int(percent_smaller*height)\n    top    += int(percent_smaller*height)\n    right  -= int(percent_smaller*width)\n    left   += int(percent_smaller*width)\n        \n    if height < 100 or width < 100:\n        print(\"Error: squareUp: bottom:\", bottom, \"top:\", top)\n        print(\"Error: squareUp: right:\", right, \"left:\", left)\n        return img\n    \n    return img[top:bottom, left:right]\n\n\ndef bensYCC(bgr, weight=4, gamma=20):\n        \n    # convert to y, cr, cb so that we can modify the image based on just the y (brightness)\n    ycc = cv2.cvtColor(bgr, cv2.COLOR_BGR2YCrCb)\n    y, cr, cb = cv2.split(ycc)\n\n    # perform bens algorithm on the y component\n    y = cv2.addWeighted(y, weight, cv2.GaussianBlur(y, (0,0), gamma), -weight, 128)\n\n    # merge the ycc back together, and recolor it\n    ycc_modified = cv2.merge((y, cr, cb))\n    bens = cv2.cvtColor(ycc_modified, cv2.COLOR_YCrCb2BGR)\n    \n    return bens \n\n\n# def claheYCC(bgr, clipLimit=5, grid=8):\n    \n#     # convert to y, cr, cb so that we can modify the image based on just the y (brightness)\n#     ycc = cv2.cvtColor(bgr, cv2.COLOR_BGR2YCrCb)\n#     y, cr, cb = cv2.split(ycc)\n\n#     # perform the clahe algorithm on the y component\n#     clahe = cv2.createCLAHE(clipLimit=clipLimit, tileGridSize=(grid, grid))\n#     y = clahe.apply(y)\n#     y = adjust_gamma(y, 1+np.log(110)-np.log(np.median(y)))\n\n#     # merge the ycc back together, and recolor it\n#     ycc_modified = cv2.merge((y, cr, cb))\n#     img = cv2.cvtColor(ycc_modified, cv2.COLOR_YCrCb2BGR)\n    \n#     return img \n\ndef bensSimple(bgr, weight=4, gamma=20):\n        \n    img = cv2.addWeighted(bgr, weight, cv2.GaussianBlur(bgr, (0,0), gamma), -weight, 128)\n    \n    return img \n\n\n# def reflectAndSquareUp(img):\n    \n#     height = img.shape[0]\n#     width  = img.shape[1]\n    \n#     # if its portrait mode, it's probably already kind of square. Make it properly square by cutting\n#     # down the height until the dimensions match\n#     if (height > width):\n        \n#         offset = int((height - width)/2)\n#         return img[offset:offset+width]\n    \n#     # otherwise, do the whole reflection thingo\n#     else:\n#         if len(img.shape) == 3:\n#             new_img = np.zeros((width, width, img.shape[2]), np.uint8)\n#         else:\n#             new_img = np.zeros((width, width), np.uint8)\n\n#         #  0  |\n#         #     |\n#         #  h1 |####\n#         #     |####\n#         #     |####\n#         #  h2 |\n#         #     |\n\n#         h1 = int((width - height)/2)\n#         h2 = h1 + height\n\n#         # paste the original into the center\n#         new_img[h1:h2,:] = img\n\n#         # paste in the reflections\n#         for i in range(h1):\n#             new_img[h1-i] = img[i]\n\n#         for i in range(width - h2):\n#             new_img[h2+i] = img[height - i - 1]\n\n#         return new_img\n\n# def circleMask(img):\n    \n#     if (img.shape[0] != img.shape[1]):\n#         print(\"Error: circle mask assumes square image\")\n#         return img\n    \n#     dim = img.shape[0]\n#     half = int(dim/2)\n    \n#     # crop out circle:\n#     circle_mask = np.zeros((dim, dim), np.uint8)\n#     circle_mask = cv2.circle(circle_mask, (half, half), half, 1, thickness=-1)\n\n#     return cv2.bitwise_and(img, img, mask=circle_mask)\n\n\n# build a lookup table mapping the pixel values [0, 255] to\n# their adjusted gamma values\n# def adjust_gamma(image, gamma=1.0):\n#     invGamma = 1.0 / gamma\n#     table = np.array([((i / 255.0) ** invGamma) * 255\n#                      for i in np.arange(0, 256)]).astype(\"uint8\")\n#     return cv2.LUT(image, table)\n\n\ndef process(bgr, model):\n    \n    green = bgr[:,:,1] # use green as a greyscale\n    \n    # ========= Crop into 640 x 480 format ========\n    if bgr.shape != (480, 640, 3):\n        \n        cropped = crop(green, bgr, 0.02)\n        width  = int(cropped.shape[1] * 0.9)\n        height = int(width * 480 / 640)\n        if height > cropped.shape[0]:\n            height = cropped.shape[0] - 2\n        h = int((cropped.shape[0] - height) / 2)\n        w = int((cropped.shape[1] - width) / 2)\n\n        test_crop = cropped[h:height+h,w:width+w,:]\n        \n    else:\n        test_crop = bgr\n    \n    if model == \"normal\":\n        colouring_fn = bensYCC\n    elif model == \"weird\":\n        colouring_fn = bensSimple\n    elif model == \"clahe\":\n        colouring_fn = claheYCC\n    else:\n        print(f\"Error, invalid model type given: {model}\")\n        \n    resized = cv2.resize(test_crop, (IMG_DIM, IMG_DIM), interpolation=cv2.INTER_AREA)\n    img     = colouring_fn(resized)\n    \n    return cv2.cvtColor(img, cv2.COLOR_BGR2RGB)\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# Data generators, which put their data into a random crop generator, which is then fed into\n# the network during training\n\ndef dataGenerator(jitter=0.1):\n    datagen = image.ImageDataGenerator(rescale=1./255,\n                                       horizontal_flip  = True and (jitter > 0.01), \n                                       vertical_flip    = True and (jitter > 0.01),\n                                       zoom_range       = [max(0.8, 1-5*jitter), 1],\n                                       rotation_range   = int(600*jitter),\n                                       brightness_range = [1-jitter/3, 1+jitter/3],\n                                       fill_mode        = \"mirror\",\n                                       channel_shift_range=int(30*jitter),\n                                      )\n    return datagen\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def test_datagen_plot(processing_function, jitter=0.03):\n    \n    images_dir  = f\"{INPUT_FOLDER}test_images/\"\n    df          = pd.read_csv(f\"{INPUT_FOLDER}test.csv\")\n    df.id_code  = df.id_code.apply(lambda x: x + \".png\")\n\n    # process images\n    img_block = np.empty((100, IMG_DIM, IMG_DIM, CHANNELS))\n    j = 1\n    for i, filename in enumerate(df.sample(100).id_code):\n        bgr = cv2.imread(images_dir + filename)\n        img_block[i,:,:,:] = process(bgr, processing_function)\n#         if bgr.shape != (480, 640, 3) and j <= 8:\n#             ax = figure.add_subplot(4,4, j)\n#             plt.imshow(img_block[i,:,:,:]/255.)\n#             j += 1\n#         elif bgr.shape == (480, 640, 3) and j > 8:\n#             ax = figure.add_subplot(4,4, j)\n#             plt.imshow(img_block[i,:,:,:]/255.)\n#             j += 1\n#             if j > 16: \n#                 return\n\n    datagen_sample = dataGenerator(jitter).flow(img_block)\n    for x in datagen_sample:\n        for j in range(16):\n            ax = figure.add_subplot(4,4, j+1)\n            plt.imshow(x[j])\n        break\n\nfigure=plt.figure(figsize=(22,20))\ntest_datagen_plot(\"weird\")\ngc.collect()\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Model Importing Function"},{"metadata":{"trusted":true},"cell_type":"code","source":"def load_network(network_name):\n    \n    weights = f\"../input/densenetmulti/{network_name}.h5\"\n    if network_name == \"normal\":\n        weights = f\"../input/densenetmulti/dense-0.800.h5\"\n    \n    model = Sequential()\n    model.add(DenseNet121(weights=None, include_top=False, input_shape=(IMG_DIM, IMG_DIM, CHANNELS)))\n    model.add(GlobalAveragePooling2D())\n    model.add(Dropout(0.5))\n    model.add(Dense(NUM_CLASSES, activation='sigmoid'))\n    \n    model.load_weights(weights)\n    model.compile(optimizer=Adam(lr=0.00005), loss='binary_crossentropy',  metrics=['accuracy'])\n    \n    return model\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Make jittered predictions on either data set"},{"metadata":{"trusted":true},"cell_type":"code","source":"# def prediction_convert_sum(predictions, thresholds):\n    \n#     thresholded = np.zeros(predictions.shape)\n    \n#     for i in range(NUM_CLASSES):\n#         thresholded[:,i] = predictions[:,i] > thresholds[i]\n#     y_val = thresholded.astype(int).sum(axis=1) - 1\n#     return y_val\n\ndef prediction_convert_highest(predictions, thresholds):\n    thresholded = np.zeros(predictions.shape)\n    for i in range(NUM_CLASSES):\n        thresholded[:,i] = predictions[:,i] > thresholds[i]\n        \n    y_val = np.zeros((predictions.shape[0]), dtype=np.int)\n    for i in range(predictions.shape[0]):\n        for j in range(4, -1, -1):\n            if thresholded[i][j]:\n                y_val[i] = j\n                break\n    return y_val\n\ndef make_predictions(d_set, models):\n\n    images_dir  = f\"{INPUT_FOLDER}{d_set}_images/\"\n    df          = pd.read_csv(f\"{INPUT_FOLDER}{d_set}.csv\")\n    df.id_code  = df.id_code.apply(lambda x: x + \".png\")\n\n    block_size  = 512\n    total       = df.index.size\n\n    jitter_amounts = [0, 0.02, 0.02, 0.02, 0.03, 0.03, 0.03, 0.1]\n#     jitter_amounts = [0.01, 0.02, 0.05]\n    \n    ensemble_predictions = np.zeros((df.index.size, len(jitter_amounts)*len(models), NUM_CLASSES))\n    \n    for m, model in enumerate(models):\n    \n        print(f\"Making predictions with the {model} model on the {d_set} dataset.\")\n        neural_net = load_network(model)\n\n        # have to do everything in blocks of images\n        for start in range(0, total, block_size):\n\n            # get subset\n            end = start + block_size\n            if end > total:\n                end = total\n\n            # process images\n            img_block = np.empty((end-start, IMG_DIM, IMG_DIM, CHANNELS))\n            for i, filename in enumerate(df[start:end].id_code):\n                try:\n                    bgr = cv2.imread(images_dir + filename)\n                    img_block[i,:,:,:] = process(bgr, model)\n                except:\n                    print(\"Error opening or manipulating image\")\n                    img_block[i,:,:,:] = 128.\n\n            # make predictions\n            for i, jit in enumerate(jitter_amounts):\n                datagen = dataGenerator(jit).flow(img_block, shuffle=False)\n                ensemble_predictions[start:end, i + len(models)*m] = neural_net.predict_generator(generator=datagen, \n                                                                      steps=len(datagen), workers=4, verbose=1)\n            \n            print(f\"{start} - {end} finished\")\n            gc.collect()\n        \n    return np.median(ensemble_predictions, axis=1)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Find the best Class Thresholds"},{"metadata":{"trusted":true},"cell_type":"code","source":"# def find_best_thresholds(train_predictions):\n    \n#     # get predictions\n#     print(\"Finding best thresholds...\")\n    \n#     prediction_convert = prediction_convert_highest\n    \n#     # make test predictions\n#     gc.collect()\n    \n#     train_df = pd.read_csv(f\"{INPUT_FOLDER}train.csv\")\n#     y_actual = train_df.diagnosis.astype(int).values\n    \n#     thresholds = [0.5 for i in range(NUM_CLASSES)]\n#     d_thresh   = 0.25\n    \n#     for sweep in range(5):\n\n#         for label in range(5):\n            \n#             currKappa = cohen_kappa_score(y_actual, prediction_convert(train_predictions, thresholds), weights='quadratic')\n            \n#             print(currKappa)\n            \n#             thresholds[label] += d_thresh\n#             kappaUp = cohen_kappa_score(y_actual, prediction_convert(train_predictions, thresholds), weights='quadratic')\n\n#             thresholds[label] -= 2*d_thresh\n#             kappaDown = cohen_kappa_score(y_actual, prediction_convert(train_predictions, thresholds), weights='quadratic')\n            \n#             thresholds[label] += d_thresh \n            \n#             if kappaUp > currKappa:\n#                 thresholds[label] += d_thresh\n#             elif kappaDown > currKappa:\n#                 thresholds[label] -= d_thresh\n\n#         d_thresh /= 2\n    \n#     gc.collect()\n#     return thresholds\n\n# train_predictions = make_predictions(\"train\", [\"normal\", \"weird\"])\n# best_thresholds   = find_best_thresholds(train_predictions)\n# print(best_thresholds)\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"\n# thresholds = [0.5 for i in range(5)]\n\n# train_df = pd.read_csv(f\"{INPUT_FOLDER}train.csv\")\n# y_actual = train_df.diagnosis.astype(int).values\n# train_fracs = train_df.diagnosis.value_counts()/train_df.index.size\n# print(train_fracs)\n\n# for sweep in range(5):\n#     for i in range(4):\n#         _, counts = np.unique(prediction_convert_highest(train_predictions, thresholds), return_counts=True)\n#         pred_fracs = counts / len(train_predictions)\n\n#         print(pred_fracs)\n\n#         # if the train fraction is higher than prediction fraction, we need to increase the thresholds\n#         try:\n#             while train_fracs[i] > pred_fracs[i]:\n#                 for j in range(i, 5):\n#                     thresholds[j] += 0.01\n#                 _, counts = np.unique(prediction_convert_highest(train_predictions, thresholds), return_counts=True)\n#                 pred_fracs = counts / len(train_predictions)\n#         except:\n#             for j in range(i, 5):\n#                 thresholds[j] -= 0.01\n#             _, counts = np.unique(prediction_convert_highest(train_predictions, thresholds), return_counts=True)\n#             pred_fracs = counts / len(train_predictions)\n\n#         # conversely decrease them if lower\n#         try:\n#             while train_fracs[i] < pred_fracs[i]:\n#                 for j in range(i, 5):\n#                     thresholds[j] -= 0.01\n#                 _, counts = np.unique(prediction_convert_highest(train_predictions, thresholds), return_counts=True)\n#                 pred_fracs = counts / len(train_predictions)\n#         except:\n#             for j in range(i, 5):\n#                 thresholds[j] += 0.01\n#             _, counts = np.unique(prediction_convert_highest(train_predictions, thresholds), return_counts=True)\n#             pred_fracs = counts / len(train_predictions)\n\n#         print(pred_fracs, \"\\n\")\n    \n#     print(cohen_kappa_score(y_actual, prediction_convert_highest(train_predictions, thresholds), weights='quadratic'), thresholds, \"\\n\")\n    ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"\n# def label_convert_two_stage(stage_1_preds, stage_2_preds):\n    \n#     thresh_1 = np.zeros((stage_1_preds.shape[0], 3))\n#     thresh_2 = np.zeros((stage_2_preds.shape[0], 3))\n    \n#     for i in range(3):\n#         thresh_1[:,i] = stage_1_preds[:,i] > 0.5\n#         thresh_2[:,i] = stage_2_preds[:,i] > 0.5\n    \n#     y_val   = thresh_1.astype(int).sum(axis=1) - 1\n#     y_val_2 = thresh_2.astype(int).sum(axis=1) + 1\n    \n#     for i in range(stage_1_preds.shape[0]):\n#         if y_val[i] == 2:\n#             y_val[i] = y_val_2[i]\n#     return y_val\n\n# def label_convert_two_stage_top(stage_1_preds, stage_2_preds):\n    \n#     thresh_1 = np.zeros((stage_1_preds.shape[0], 3))\n#     thresh_2 = np.zeros((stage_2_preds.shape[0], 3))\n    \n#     for i in range(3):\n#         thresh_1[:,i] = stage_1_preds[:,i] > 0.5\n#         thresh_2[:,i] = stage_2_preds[:,i] > 0.5\n    \n#     y_val   = np.zeros(stage_1_preds.shape[0])\n#     y_val_2 = np.zeros(stage_1_preds.shape[0])\n#     for i in range(stage_1_preds.shape[0]):\n#         for j in range(2, -1, -1):\n#             if thresh_1:\n#                 y_val[i] = j\n#         for j in range(2, -1, -1):\n#             if thresh_2:\n#                 y_val_2[i] = j + 2\n    \n#     for i in range(stage_1_preds.shape[0]):\n#         if y_val[i] == 2:\n#             y_val[i] = y_val_2[i]\n#     return y_val","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Perform predictions on test set "},{"metadata":{"trusted":true},"cell_type":"code","source":"# thresholds = [0.5 for i in range(NUM_CLASSES)]\nthresholds = [0.5, 0.5, 0.4, 0.4, 0.3]\n\n# make test predictions\npredictions = make_predictions(\"test\", [\"normal\", \"weird\"])\n\nas_classes = prediction_convert_highest(predictions, thresholds)\nprint(as_classes[:10])\n\n# save to csv\ntest_df = pd.read_csv(INPUT_FOLDER + 'test.csv')\ntest_df['diagnosis'] = as_classes\ntest_df.to_csv('submission.csv', index=False)","execution_count":null,"outputs":[]}],"metadata":{"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.8"}},"nbformat":4,"nbformat_minor":1}