{"cells":[{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport cv2\nimport tensorflow as tf\nfrom tqdm import tqdm\nimport random\n\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.utils import class_weight\n\nimport scipy as sp\nfrom functools import partial\nfrom sklearn import metrics\n\n%matplotlib inline","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":true},"cell_type":"code","source":"train_df = pd.read_csv('../input/aptos2019-blindness-detection/train.csv')\ntest_df = pd.read_csv('../input/aptos2019-blindness-detection/test.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"IMG_SIZE = 224","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def load_img(path):\n    image = cv2.imread(path)\n    return image","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def remove_unwanted_space(image, threshold=7):\n    gray_image = cv2.cvtColor(image, cv2.COLOR_BGR2GRAY)\n    mask = gray_image > threshold\n    cropped_image = image[np.ix_(mask.any(1), mask.any(0))]\n    if cropped_image.shape[0] == 0:\n        return image\n    return cropped_image","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def preprocess_img(path):\n    image = load_img(path)\n    image = remove_unwanted_space(image, 5)\n    image = cv2.cvtColor(image, cv2.COLOR_BGR2RGB)\n    image = cv2.resize(image, (IMG_SIZE, IMG_SIZE))\n    image = cv2.addWeighted(image,4, cv2.GaussianBlur(image, (0,0), 35), -4, 128)\n    image = image/255.0\n    return image","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.imshow(preprocess_img('../input/aptos2019-blindness-detection/train_images/'+train_df.loc[1]['id_code']+'.png'))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"N = train_df.shape[0]\ntrain = np.empty((N, IMG_SIZE, IMG_SIZE, 3), dtype=np.float32)\nfor i, image_id in enumerate(tqdm(train_df['id_code'])):\n    train[i,:,:,:] = preprocess_img('../input/aptos2019-blindness-detection/train_images/'+image_id+'.png')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Let us not take y as a one hot vecor describing output classes. Instead let us convert y as follows:\n\n0 : [1, 0, 0, 0, 0]\n\n1 : [1, 1, 0, 0, 0]\n\n2 : [1, 1, 1, 0, 0]\n\n3 : [1, 1, 1, 1, 0]\n\n4 : [1, 1, 1, 1, 1]"},{"metadata":{},"cell_type":"markdown","source":"The intuition is that the output classes are not independent. A particular stage of the disease will also have symptoms of the previous stage."},{"metadata":{"trusted":true},"cell_type":"code","source":"y = np.zeros((train_df.shape[0],5), dtype=np.int8)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"for i, diagnosis in enumerate(train_df['diagnosis']):\n    for k in range(diagnosis+1):\n        y[i,k] = 1","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"y[:5,:]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"X_train, X_val, y_train, y_val = train_test_split(train, y, test_size=0.15, random_state=42)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"BATCH_SIZE = 32","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Let us create generators for training, validation and test,"},{"metadata":{"trusted":true},"cell_type":"code","source":"class TestSequence(tf.keras.utils.Sequence):\n    def __init__(self, x_set, batch_size, aug=True):\n        self.x = x_set\n        self.batch_size = batch_size\n        self.aug = aug\n    def __len__(self):\n        return int(np.ceil(len(self.x) / self.batch_size))\n    def __getitem__(self, idx):\n        batch_x = self.x[idx * self.batch_size:(idx + 1) * self.batch_size]\n       \n        images = np.empty((batch_x.shape[0], IMG_SIZE, IMG_SIZE, 3), dtype=np.float32)\n        for i, image_id in enumerate(batch_x):\n            preprocessed_img = preprocess_img('../input/aptos2019-blindness-detection/test_images/'+image_id+'.png')\n            if self.aug:\n                h_flip = bool(random.getrandbits(1))\n                v_flip = bool(random.getrandbits(1))\n                rotate = random.randint(0,20)\n                transform_parameters={'theta':rotate, 'flip_horizontal': h_flip, 'flip_vertical': v_flip}\n                #transform_parameters={'flip_horizontal': h_flip, 'flip_vertical': v_flip}\n                preprocessed_img = tf.keras.preprocessing.image.ImageDataGenerator().apply_transform(preprocessed_img, transform_parameters)\n            images[i,:,:,:] = preprocessed_img\n        return  images","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_data_generator = tf.keras.preprocessing.image.ImageDataGenerator(\n        fill_mode='constant',\n        cval=0.,  # value used for fill_mode = \"constant\"\n        horizontal_flip=True,  # randomly flip images\n        vertical_flip=True,  # randomly flip images\n        rotation_range=15,\n    )","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"val_data_generator = tf.keras.preprocessing.image.ImageDataGenerator()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train_gen = train_data_generator.flow(X_train, y_train, batch_size=BATCH_SIZE)\nval_gen = val_data_generator.flow(X_val, y_val, batch_size=BATCH_SIZE, shuffle=False)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Add callbacks for calculating quadratic weighted kappa score after each epoch."},{"metadata":{"trusted":true},"cell_type":"code","source":"class Metrics(tf.keras.callbacks.Callback):\n    def __init__(self, val_data, batch_size = 20):\n        super().__init__()\n        self.validation_data = val_data\n        self.batch_size = batch_size\n        \n    def on_train_begin(self, logs={}):\n        self.val_kappas = []\n\n    def on_epoch_end(self, epoch, logs={}):\n        X_val, y_val = self.validation_data[:2]\n        y_val = y_val.sum(axis=1) - 1\n        \n        y_pred = self.model.predict(X_val) > 0.5\n        y_pred = y_pred.astype(int).sum(axis=1) - 1\n\n        _val_kappa = metrics.cohen_kappa_score(\n            y_val,\n            y_pred, \n            weights='quadratic'\n        )\n\n        self.val_kappas.append(_val_kappa)\n\n        print(f\" val_kappa: {_val_kappa:.4f}\")\n        \n        if _val_kappa == max(self.val_kappas):\n            print(\"Validation Kappa has improved. Saving model.\")\n            self.model.save('model.h5')\n\n        return","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"kappa_cb = Metrics(((X_val, y_val)))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"densenet = tf.keras.applications.DenseNet121(include_top=False, weights='../input/densenet121/densenet121_weights_tf_dim_ordering_tf_kernels_notop.h5', input_shape=(IMG_SIZE, IMG_SIZE, 3))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"model = tf.keras.Sequential([\n    densenet,\n    tf.keras.layers.GlobalAveragePooling2D(),\n    tf.keras.layers.Dropout(0.5),\n    tf.keras.layers.Dense(5, activation='sigmoid')\n])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"reduce_lr = tf.keras.callbacks.ReduceLROnPlateau(monitor='val_acc', factor=0.6, patience=6, verbose=1)\nclass_weights = class_weight.compute_class_weight('balanced', np.unique(np.sum(y_train, axis=1)), np.sum(y_train, axis=1))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"optimizer = tf.keras.optimizers.Adam(lr=5e-5)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"steps_per_epoch = int(np.ceil(X_train.shape[0]/BATCH_SIZE))\nval_steps_per_epoch = int(np.ceil(X_val.shape[0]/BATCH_SIZE))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"model.compile(optimizer=optimizer, loss='binary_crossentropy', metrics=['accuracy'])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"history = model.fit_generator(train_gen, validation_data=(X_val, y_val), steps_per_epoch=steps_per_epoch, \n                              callbacks=[kappa_cb], epochs=25, class_weight=class_weights)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"model.load_weights('model.h5')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Let us just take a look at the targets and our model outputs."},{"metadata":{"trusted":true},"cell_type":"code","source":"val_preds = np.squeeze(model.predict_generator(val_gen, steps=val_steps_per_epoch))\nlist(zip(val_preds, y_val))[:10]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"pred = val_preds > 0.5\npred = np.sum(pred.astype(np.int8), axis=1)-1\nkappa, acc = metrics.accuracy_score(pred, y_val.sum(axis=1)-1), metrics.cohen_kappa_score(y_val.sum(axis=1)-1, pred, weights='quadratic')\nprint('kappa score: {:.4f}, accuracy: {:.4f}'.format(kappa, acc))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sns.heatmap(metrics.confusion_matrix(y_val.sum(axis=1)-1, pred), annot=True, cmap='Blues', fmt='g')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Now we have to optmize threshold values which determine to which classes our prediction belongs to from our model outputs. Although we have assumed the threshold as  [0.5, 1.5, 2.5, 3.5], we can further improve our quadratic weighted kappa score by optimizing the thresholds."},{"metadata":{"trusted":true},"cell_type":"code","source":"class OptimizedRounder(object):\n    def __init__(self):\n        self.coef_ = 0\n\n    def _kappa_loss(self, coef, X, y):\n        X_p = np.copy(X)\n        for i, pred in enumerate(X_p):\n            if pred < coef[0]:\n                X_p[i] = 0\n            elif pred >= coef[0] and pred < coef[1]:\n                X_p[i] = 1\n            elif pred >= coef[1] and pred < coef[2]:\n                X_p[i] = 2\n            elif pred >= coef[2] and pred < coef[3]:\n                X_p[i] = 3\n            else:\n                X_p[i] = 4\n\n        ll = metrics.cohen_kappa_score(y, X_p, weights='quadratic')\n        return -ll\n\n    def fit(self, X, y):\n        loss_partial = partial(self._kappa_loss, X=X, y=y)\n        initial_coef = [0.5, 1.5, 2.5, 3.5]\n        self.coef_ = sp.optimize.minimize(loss_partial, initial_coef, method='nelder-mead')\n\n    def predict(self, X, coef):\n        X_p = np.copy(X)\n        for i, pred in enumerate(X_p):\n            if pred < coef[0]:\n                X_p[i] = 0\n            elif pred >= coef[0] and pred < coef[1]:\n                X_p[i] = 1\n            elif pred >= coef[1] and pred < coef[2]:\n                X_p[i] = 2\n            elif pred >= coef[2] and pred < coef[3]:\n                X_p[i] = 3\n            else:\n                X_p[i] = 4\n        return X_p\n\n    def coefficients(self):\n        return self.coef_['x']","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"optR = OptimizedRounder()\noptR.fit(np.sum(val_preds,axis=1), np.sum(y_val, axis=1))\ncoefficients = optR.coefficients()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"coefficients","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sample = pd.read_csv(\"../input/aptos2019-blindness-detection/sample_submission.csv\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_gen = TestSequence(sample['id_code'], BATCH_SIZE, True)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"It is time for test time augmentation! Since our test data generator randomly augments the images, we make multiple predictions on the image and agree upon the most voted class."},{"metadata":{"trusted":true},"cell_type":"code","source":"prediction = np.empty((5,len(sample)))\nfor i in tqdm(range(5)):\n    pred = model.predict_generator(test_gen) > 0.5\n    prediction[i,:] = np.sum(pred.astype(np.int8), axis=1)-1\nprediction = prediction.astype(np.int8)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def vote(x):\n    return np.argmax(np.bincount(x))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"preds = np.apply_along_axis(vote, 0, prediction)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.hist(preds)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sample.diagnosis = preds\nsample.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}