{"cells":[{"metadata":{"_uuid":"4c08ef44cc9bdc0d491b6c5d53387105186c2ab8"},"cell_type":"markdown","source":"**Sections of this kernel**\n- Project understanding\n- Data understanding\n- Data visualization\n- Baseline model (Fastai v1)\n- Validation and analysis\n    - Metrics\n    - Prediction and activation visualizations\n    - ROC & AUC\n- Submit\n- Deploy (example)\n\nSection: Data visualization | Section: Prediction and activation visualizations\n:--------:|:-------:\n<a href=\"https://ibb.co/5Mpm3ZF\"><img src=\"https://i.ibb.co/gjqN8YV/metastases.png\" alt=\"metastases\" border=\"0\" align=”right”></a> |<a href=\"https://ibb.co/Pc6wbyk\"><img src=\"https://i.ibb.co/sKVyD47/gradcam.png\" alt=\"gradcam\" border=\"0\"></a>\n\n---------------------------------------------------\n# Project understanding\n###  What exactly is the problem?\n\n**Binary image classification problem.** Identify the presence of metastases from 96 x 96px digital histopathology images. One key challenge is that the metastases can be as small as single cells in a large area of tissue.\n\n### How would a solution look like?\n\n**Our evaluation metric is [area under the ROC curve](http://en.wikipedia.org/wiki/Receiver_operating_characteristic).** The ROC curve is a plot of *True positive rate* against *False positive rate* at various thresholds and the area under the curve (AUC) is equal to the probability that a classifier will rank a randomly chosen positive instance higher than a randomly chosen negative one. The best possible solution would yield an AUC of 1 which means we would classify all positive samples correctly without getting any false positives. \n\n![ROC curve example](https://i.ibb.co/mBKh6ZB/roc.png)\n<p style=\"text-align: center;\"> ROC curve from a previous run of this kernel </p>\n\n### What is known about the domain?\n\n**The histopathological images are glass slide microscope images of lymph nodes that are stained with hematoxylin and eosin (H&E).** This staining method is one of the most widely used in medical diagnosis and it produces blue, violet and red colors. Dark blue hematoxylin binds to negatively charged substances such as nucleic acids and pink eosin to positively charged substances like amino-acid side chains (most proteins). Typically nuclei are stained blue, whereas cytoplasm and extracellular parts in various shades of pink.\n\n**Low-resolution**             | **Mid-resolution**            | **High-resolution** \n:-------------------------:|:-------------------------:|:-------------------------:\n![](https://camelyon17.grand-challenge.org/site/CAMELYON17/serve/public_html/example_low_resolution.png) | ![Example of a metastatic region](https://camelyon17.grand-challenge.org/site/CAMELYON17/serve/public_html/example_mid_resolution.png) | ![Example of a metastatic region](https://camelyon17.grand-challenge.org/site/CAMELYON17/serve/public_html/example_high_resolution.png)\n**[<p style=\"text-align: center;\"> Example of a metastatic region in lymph nodes, CHAMELYON17 </p>](https://camelyon17.grand-challenge.org/Background/)**\n\nLymph nodes are small glands that filter the fluid in the lymphatic system and they are the first place a breast cancer is likely to spread. Histological assessment of lymph node metastases is part of determining the stage of breast cancer in TNM classification which is a globally recognized standard for classifying the extent of spread of cancer. The diagnostic procedure for pathologists is tedious and time-consuming as a large area of tissue has to be examined and small metastases can be easily missed.\n\n**Useful links for background knowledge**\n- [Patch Camelyon (PCam)](https://github.com/basveeling/pcam)\n- [Hematoxylin and eosin staining of tissue and cell sections](https://www.ncbi.nlm.nih.gov/pubmed/21356829)\n- [H&E-stained sentinel lymph node sections of breast cancer patients: the CAMELYON dataset](https://academic.oup.com/gigascience/article/7/6/giy065/5026175)\n- [CAMELYON16 - background](https://camelyon16.grand-challenge.org/Background/)\n- [CAMELYON17 - background](https://camelyon17.grand-challenge.org/Background/)\n- [TNM classification](https://www.uicc.org/resources/tnm)"},{"metadata":{"_uuid":"19fdd1381f06e397810c3fef5a31be0c9c9024bc"},"cell_type":"markdown","source":"----------------------------------------------\n# Data understanding\n### What data do we have available?\n\n**220k training images and 57k evaluation images.** The dataset is a subset of the [PCam dataset](https://github.com/basveeling/pcam) and the only difference between these two is that all duplicate images have been removed. The PCam dataset is derived from the [Camelyon16 Challenge dataset](https://camelyon16.grand-challenge.org/Data/) which contains 400 H&E stained whole slide images of sentinel lymph node sections that were acquired and digitized at 2 different centers using a 40x objective. The PCam's dataset including this one uses 10x undersampling to increase the field of view, which gives the resultant pixel resolution of 2.43 microns.\n\nAccording to the data description, there is a 50/50 balance between positive and negative examples in the training and test splits. However, **the training distribution seems to be 60/40 (negatives/positives)**. A positive label means that there is at least one pixel of tumor tissue in the center region (32 x 32px) of the image. **Tumor tissue in the outer region of the patch does not influence the label.** This means that a negatively labeled image could contain metastases in the outer region. Thus, it would be a good idea to crop the images to the center region.\n\n**Image file descriptors**\n\nDescription | \n:--------:|:-------:\nFormat | TIF\nSize | 96 x 96\nChannels | 3\nBits per channel | 8\nData type | Unsigned char\nCompression | Jpeg\n\n### Is the data relevant to the problem?\n\nThis dataset is a combination of two independent datasets collected in Radboud University Medical Center (Nijmegen, the Netherlands), and the University Medical Center Utrecht (Utrecht, the Netherlands). The slides are produced by routine clinical practices and a trained pathologist would examine similar images for identifying metastases. However, some relevant information about the surroundings might be left out with these small-sized image samples.\n\n### Is it valid? Does it reflect our expectations?\n\nAccording to the data description, the dataset has been stripped of duplicates. However, this has not been confirmed by testing.\n\n> For the entire dataset, when the slide-level label was unclear during the inspection of the H&E-stained slide, an additional WSI with a consecutive tissue section, immunohistochemically stained for cytokeratin, was used to confirm the classification.\n- [1399 H&E-stained sentinel lymph node sections of breast cancer patients: the CAMELYON dataset](https://academic.oup.com/gigascience/article/7/6/giy065/5026175)\n\n### Is the data quality, quantity, recency sufficient?\n\n> All glass slides included in the CAMELYON dataset were part of routine clinical care and are thus of diagnostic quality. However, during the acquisition process, scanning can fail or result in out-of-focus images. As a quality-control measure, all slides were inspected manually after scanning. The inspection was performed by an experienced technician (Q.M. and N.S. for UMCU, M.H. or R.vd.L. for the other centers) to assess the quality of the scan; when in doubt, a pathologist was consulted on whether scanning issues might affect diagnosis.\n- [1399 H&E-stained sentinel lymph node sections of breast cancer patients: the CAMELYON dataset](https://academic.oup.com/gigascience/article/7/6/giy065/5026175)"},{"metadata":{"_uuid":"dfc41ba2b70e88b92ff41a94f0b6504336f2360d"},"cell_type":"markdown","source":"-----------------------------------------\n# Data visualization"},{"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 os\nimport cv2\nimport matplotlib.pyplot as plt\nimport matplotlib.patches as patches\nimport random\nfrom sklearn.utils import shuffle\nfrom tqdm import tqdm_notebook\n\ndata = pd.read_csv('/kaggle/input/train_labels.csv')\ntrain_path = '/kaggle/input/train/'\ntest_path = '/kaggle/input/test/'\n# quick look at the label stats\ndata['label'].value_counts()","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","collapsed":true,"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","trusted":false},"cell_type":"markdown","source":"We can see that the negative/positive ratio is not entirely 50/50 as there are 130k negatives and 90k negatives. The ratio is closer to 60/40 meaning that there are 1.5 times more negative images than positives.\n\n### Plot some images with and without cancer tissue for comparison"},{"metadata":{"trusted":true,"_uuid":"90811cbf83e65ba2679ed1044798a71b2dd89818","_kg_hide-input":true},"cell_type":"code","source":"def readImage(path):\n    # OpenCV reads the image in bgr format by default\n    bgr_img = cv2.imread(path)\n    # We flip it to rgb for visualization purposes\n    b,g,r = cv2.split(bgr_img)\n    rgb_img = cv2.merge([r,g,b])\n    return rgb_img","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"20e119981333d59b0ba4a52b46a9aa5b4bb2436f","_kg_hide-input":true},"cell_type":"code","source":"# random sampling\nshuffled_data = shuffle(data)\n\nfig, ax = plt.subplots(2,5, figsize=(20,8))\nfig.suptitle('Histopathologic scans of lymph node sections',fontsize=20)\n# Negatives\nfor i, idx in enumerate(shuffled_data[shuffled_data['label'] == 0]['id'][:5]):\n    path = os.path.join(train_path, idx)\n    ax[0,i].imshow(readImage(path + '.tif'))\n    # Create a Rectangle patch\n    box = patches.Rectangle((32,32),32,32,linewidth=4,edgecolor='b',facecolor='none', linestyle=':', capstyle='round')\n    ax[0,i].add_patch(box)\nax[0,0].set_ylabel('Negative samples', size='large')\n# Positives\nfor i, idx in enumerate(shuffled_data[shuffled_data['label'] == 1]['id'][:5]):\n    path = os.path.join(train_path, idx)\n    ax[1,i].imshow(readImage(path + '.tif'))\n    # Create a Rectangle patch\n    box = patches.Rectangle((32,32),32,32,linewidth=4,edgecolor='r',facecolor='none', linestyle=':', capstyle='round')\n    ax[1,i].add_patch(box)\nax[1,0].set_ylabel('Tumor tissue samples', size='large')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"c63783a335ea6f078c07fcd1024baf86e0255234"},"cell_type":"markdown","source":"Classifying metastases is probably not an easy task for a trained pathologist and extremely difficult for an untrained eye. According to [Libre Pathology](https://librepathology.org/wiki/Lymph_node_metastasis), lymph node metastases can have these features:\n\n> - Foreign cell population - key feature (Classic location: subcapsular sinuses)\n- Cells with cytologic features of malignancy\n    - Nuclear pleomorphism (variation in size, shape and staining).\n    - Nuclear atypia:\n        - **Nuclear enlargement**.\n        - **Irregular nuclear membrane**.\n        - **Irregular chromatin pattern, esp. asymmetry**.\n        - **Large or irregular nucleolus**.\n     - Abundant mitotic figures.\n- Cells in architectural arrangements seen in malignancy; highly variable - dependent on tumour type and differentiation.\n    - Gland formation.\n    - Single cells.\n    - Small clusters of cells.\n  \n**The takeaway from this is probably that irregular nuclear shapes, sizes or staining shades can indicate metastases.**\n\n### How is the data best transformed for modeling?\n\nWe know that the label of the image is influenced only by the center region (32 x 32px) so it would make sense to crop our data to that region only. However, some useful information about the surroundings could be lost if we crop too close.  This hypothesis could be confirmed by training models with varying crop sizes. My initial results with 32 x 32px size showed worse performance than with 48 x 48px but I haven't done a search for optimal size.\n\n### How may we increase the data quality?\n\nWe could inspect if the data contains bad data (too unfocused or corrupted) and remove those to increase the overall quality."},{"metadata":{"_uuid":"74f542bfe934472aa9ff06a4031b9a79738f5db5"},"cell_type":"markdown","source":"### Preprocessing and augmentation\nThere are couple of ways we can use to avoid overfitting; more data, augmentation, regularization and less complex model architectures. Here we will define what image augmentations to use and add them directly to our image loader function. Note that if we apply augmentation here, augmentations will also be applied when we are predicting (inference). This is called test time augmentation (TTA) and it can improve our results if we run inference multiple times for each image and average out the predictions. \n\n**The augmentations we can use for this type of data:**\n- random rotation\n- random crop\n- random flip (horizontal and vertical both)\n- random lighting\n- random zoom (not implemented here)\n- Gaussian blur (not implemented here)\n\nWe will use OpenCV with image operations because in my experience, OpenCV is a lot faster than *PIL* or *scikit-image*."},{"metadata":{"trusted":true,"_uuid":"d08cfe0546459f2d9b5c8fd4c5d4994527677173"},"cell_type":"code","source":"import random\nORIGINAL_SIZE = 96      # original size of the images - do not change\n\n# AUGMENTATION VARIABLES\nCROP_SIZE = 90          # final size after crop\nRANDOM_ROTATION = 3    # range (0-180), 180 allows all rotation variations, 0=no change\nRANDOM_SHIFT = 2        # center crop shift in x and y axes, 0=no change. This cannot be more than (ORIGINAL_SIZE - CROP_SIZE)//2 \nRANDOM_BRIGHTNESS = 7  # range (0-100), 0=no change\nRANDOM_CONTRAST = 5    # range (0-100), 0=no change\nRANDOM_90_DEG_TURN = 1  # 0 or 1= random turn to left or right\n\ndef readCroppedImage(path, augmentations = True):\n    # augmentations parameter is included for counting statistics from images, where we don't want augmentations\n    \n    # OpenCV reads the image in bgr format by default\n    bgr_img = cv2.imread(path)\n    # We flip it to rgb for visualization purposes\n    b,g,r = cv2.split(bgr_img)\n    rgb_img = cv2.merge([r,g,b])\n    \n    if(not augmentations):\n        return rgb_img / 255\n    \n    #random rotation\n    rotation = random.randint(-RANDOM_ROTATION,RANDOM_ROTATION)\n    if(RANDOM_90_DEG_TURN == 1):\n        rotation += random.randint(-1,1) * 90\n    M = cv2.getRotationMatrix2D((48,48),rotation,1)   # the center point is the rotation anchor\n    rgb_img = cv2.warpAffine(rgb_img,M,(96,96))\n    \n    #random x,y-shift\n    x = random.randint(-RANDOM_SHIFT, RANDOM_SHIFT)\n    y = random.randint(-RANDOM_SHIFT, RANDOM_SHIFT)\n    \n    # crop to center and normalize to 0-1 range\n    start_crop = (ORIGINAL_SIZE - CROP_SIZE) // 2\n    end_crop = start_crop + CROP_SIZE\n    rgb_img = rgb_img[(start_crop + x):(end_crop + x), (start_crop + y):(end_crop + y)] / 255\n    \n    # Random flip\n    flip_hor = bool(random.getrandbits(1))\n    flip_ver = bool(random.getrandbits(1))\n    if(flip_hor):\n        rgb_img = rgb_img[:, ::-1]\n    if(flip_ver):\n        rgb_img = rgb_img[::-1, :]\n        \n    # Random brightness\n    br = random.randint(-RANDOM_BRIGHTNESS, RANDOM_BRIGHTNESS) / 100.\n    rgb_img = rgb_img + br\n    \n    # Random contrast\n    cr = 1.0 + random.randint(-RANDOM_CONTRAST, RANDOM_CONTRAST) / 100.\n    rgb_img = rgb_img * cr\n    \n    # clip values to 0-1 range\n    rgb_img = np.clip(rgb_img, 0, 1.0)\n    \n    return rgb_img","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"4e36792e155852f9bfcc9f564da0328ff872bf0d","_kg_hide-input":true},"cell_type":"code","source":"fig, ax = plt.subplots(2,5, figsize=(20,8))\nfig.suptitle('Cropped histopathologic scans of lymph node sections',fontsize=20)\n# Negatives\nfor i, idx in enumerate(shuffled_data[shuffled_data['label'] == 0]['id'][:5]):\n    path = os.path.join(train_path, idx)\n    ax[0,i].imshow(readCroppedImage(path + '.tif'))\nax[0,0].set_ylabel('Negative samples', size='large')\n# Positives\nfor i, idx in enumerate(shuffled_data[shuffled_data['label'] == 1]['id'][:5]):\n    path = os.path.join(train_path, idx)\n    ax[1,i].imshow(readCroppedImage(path + '.tif'))\nax[1,0].set_ylabel('Tumor tissue samples', size='large')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"9568cbee8f077da727aa84fd13f459e5d03c5a44"},"cell_type":"markdown","source":"**To see the effects of our augmentation, we can plot one image multiple times.**"},{"metadata":{"trusted":true,"_uuid":"dd77955fd1742b4d72562113beca5712011122db","_kg_hide-input":true},"cell_type":"code","source":"fig, ax = plt.subplots(1,5, figsize=(20,4))\nfig.suptitle('Random augmentations to the same image',fontsize=20)\n# Negatives\nfor i, idx in enumerate(shuffled_data[shuffled_data['label'] == 0]['id'][:1]):\n    for j in range(5):\n        path = os.path.join(train_path, idx)\n        ax[j].imshow(readCroppedImage(path + '.tif'))","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"f4a06249f6905afeb3e333c2e2e940196683461b"},"cell_type":"markdown","source":"### Compute image statistics\n**Do not use augmentation here!**\n\nCalculating statistics will give channel averages of [0.702447, 0.546243, 0.696453],\nand std's of [0.238893, 0.282094, 0.216251].\n\nWhile we are calculating statistics, we can check if there are images that have a very low maximum pixel intensity (almost totally black) or very high minimum pixel intensity (almost totally white). These kind of images could be caused by bad exposure or cropping to an empty area. In case of an empty area, the image would not be an outlier but equally valid negative sample.\n\nWe find that there is at least one very dark and 6 very bright images."},{"metadata":{"_kg_hide-input":true,"trusted":true,"_uuid":"065187a233372649ff42a734ad7380e4a16288b6"},"cell_type":"code","source":"# As we count the statistics, we can check if there are any completely black or white images\ndark_th = 10 / 255      # If no pixel reaches this threshold, image is considered too dark\nbright_th = 245 / 255   # If no pixel is under this threshold, image is considerd too bright\ntoo_dark_idx = []\ntoo_bright_idx = []\n\nx_tot = np.zeros(3)\nx2_tot = np.zeros(3)\ncounted_ones = 0\nfor i, idx in tqdm_notebook(enumerate(shuffled_data['id']), 'computing statistics...(220025 it total)'):\n    path = os.path.join(train_path, idx)\n    imagearray = readCroppedImage(path + '.tif', augmentations = False).reshape(-1,3)\n    # is this too dark\n    if(imagearray.max() < dark_th):\n        too_dark_idx.append(idx)\n        continue # do not include in statistics\n    # is this too bright\n    if(imagearray.min() > bright_th):\n        too_bright_idx.append(idx)\n        continue # do not include in statistics\n    x_tot += imagearray.mean(axis=0)\n    x2_tot += (imagearray**2).mean(axis=0)\n    counted_ones += 1\n    \nchannel_avr = x_tot/counted_ones\nchannel_std = np.sqrt(x2_tot/counted_ones - channel_avr**2)\nchannel_avr,channel_std","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"20dcd95df7b64e785026c8a9816809e5fd0608f9"},"cell_type":"code","source":"print('There was {0} extremely dark image'.format(len(too_dark_idx)))\nprint('and {0} extremely bright images'.format(len(too_bright_idx)))\nprint('Dark one:')\nprint(too_dark_idx)\nprint('Bright ones:')\nprint(too_bright_idx)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"21aa206de3ff96ee7c976f5478a800232a302684"},"cell_type":"markdown","source":"### Plot some of the very bright or very dark images"},{"metadata":{"_kg_hide-input":true,"trusted":true,"_uuid":"e142e95dd1f32e1cc81523b1e266dcf521def1ef"},"cell_type":"code","source":"fig, ax = plt.subplots(2,6, figsize=(25,9))\nfig.suptitle('Almost completely black or white images',fontsize=20)\n# Too dark\ni = 0\nfor idx in np.asarray(too_dark_idx)[:min(6, len(too_dark_idx))]:\n    lbl = shuffled_data[shuffled_data['id'] == idx]['label'].values[0]\n    path = os.path.join(train_path, idx)\n    ax[0,i].imshow(readCroppedImage(path + '.tif', augmentations = False))\n    ax[0,i].set_title(idx + '\\n label=' + str(lbl), fontsize = 8)\n    i += 1\nax[0,0].set_ylabel('Extremely dark images', size='large')\nfor j in range(min(6, len(too_dark_idx)), 6):\n    ax[0,j].axis('off') # hide axes if there are less than 6\n# Too bright\ni = 0\nfor idx in np.asarray(too_bright_idx)[:min(6, len(too_bright_idx))]:\n    lbl = shuffled_data[shuffled_data['id'] == idx]['label'].values[0]\n    path = os.path.join(train_path, idx)\n    ax[1,i].imshow(readCroppedImage(path + '.tif', augmentations = False))\n    ax[1,i].set_title(idx + '\\n label=' + str(lbl), fontsize = 8)\n    i += 1\nax[1,0].set_ylabel('Extremely bright images', size='large')\nfor j in range(min(6, len(too_bright_idx)), 6):\n    ax[1,j].axis('off') # hide axes if there are less than 6","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"4e437a4abe4ca906cf5b580474fed1e060f26c13"},"cell_type":"markdown","source":"All the dark and bright images are labeled negative. I think the bright ones are just cropped from a non-stained part or they don't have any tissue (plain glass?) so the labels are correct. The samples don't have tumor tissue present. I am not so sure about the dark image, is it an outlier crop from badly exposed area or just some very large cell part filling the whole image. Anyway, removing only a small amount of outliers from this size data set has little or no effect on the prediction performance.\n\n-----------------------------------------"},{"metadata":{"_uuid":"2df0d87a864207ce1260acad6926f6dcf8d5fd58"},"cell_type":"markdown","source":"# Baseline model (Fastai v1)\nIn ML production pipeline, it is a good idea to start with a relatively simple model, sort of a minimum viable product (MVP) or a baseline. With MVP, we can very quickly see if there are some unexpected problems like bad data quality that will make any further investments into the model tuning not worth it."},{"metadata":{"_uuid":"fe849b5a44217327e0451c0797bb981e7d700228"},"cell_type":"markdown","source":"### Prepare the data and split train\nSplit train data to 90% training and 10% validation parts. We want to maintain equal ratios of negative/positive (60/40) in both, training and test splits. This is not so crucial here as both labels are almost equally represented but in case we had a rare class, random split could cause severe underrepresentation or in worst case, leave all rare classes out of one split. "},{"metadata":{"trusted":true,"_uuid":"31fef9db25fe198c5698d90db93ac91ddb60886c","_kg_hide-output":false},"cell_type":"code","source":"from sklearn.model_selection import train_test_split\n\n# we read the csv file earlier to pandas dataframe, now we set index to id so we can perform\ntrain_df = data.set_index('id')\n\n#If removing outliers, uncomment the four lines below\n#print('Before removing outliers we had {0} training samples.'.format(train_df.shape[0]))\n#train_df = train_df.drop(labels=too_dark_idx, axis=0)\n#train_df = train_df.drop(labels=too_bright_idx, axis=0)\n#print('After removing outliers we have {0} training samples.'.format(train_df.shape[0]))\n\ntrain_names = train_df.index.values\ntrain_labels = np.asarray(train_df['label'].values)\n\n# split, this function returns more than we need as we only need the validation indexes for fastai\ntr_n, tr_idx, val_n, val_idx = train_test_split(train_names, range(len(train_names)), test_size=0.1, stratify=train_labels, random_state=123)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"a7a0c283990b578078e80772ed44c9f28b11fa79"},"cell_type":"markdown","source":"We will use [**Fast.ai V1**](https://docs.fast.ai/) software library that is built on [PyTorch](https://pytorch.org/). What I like about Fast.ai is that it includes \"out of the box\" support for many recent advancements in deep learning research. If you want to use the 0.7 version of Fast.ai, see commit version 9 of this kernel.\n\nI highly recommend the [Fastai: practical deep learning course for coders, v3](https://course.fast.ai/).\n\nI want to thank [Martijn](https://www.kaggle.com/etinuz) for sharing his examples of Fastai v1 training and visualization! They were helpful when converting this kernel from Fastai v0.7 to v1.0.\n\n### What kind of model architecture suits the problem best?\n\nHere we will be using a pretrained convnet model and transfer learning to adjust the weights to our data. Going for a deeper model architecture will start overfitting faster.\n\nFor differenet pretrained model architectures, check [Fast.ai vision models](https://docs.fast.ai/vision.models.html) and [torchvision models](https://pytorch.org/docs/stable/torchvision/models.html)."},{"metadata":{"trusted":true,"_uuid":"9d217ce198f2f23fc722f8606a722fb9b3e71d38"},"cell_type":"code","source":"# fastai 1.0\nfrom fastai import *\nfrom fastai.vision import *\nfrom torchvision.models import *    # import *=all the models from torchvision  \n\narch = alexnet                  # specify model architecture\nBATCH_SIZE = 128                    # specify batch size, hardware restrics this one. Large batch sizes may run out of GPU memory\nsz = CROP_SIZE                      # input size is the crop size\nMODEL_PATH = str(arch).split()[1]   # this will extrat the model name as the model file name e.g. 'resnet50'","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"0dbd4bddb400e4f0b5aa7d94d8dda552f2e883a5"},"cell_type":"markdown","source":"We load the images to an ImageDataBunch for the training. This fastai's data object is easily customized for loading images using our own `readCroppedImage` function. We just need to subclass ImageList. \n[Data API docs](https://docs.fast.ai/data_block.html)"},{"metadata":{"trusted":true,"_uuid":"c1b95582bc4bf32c8a3686325f165ad72b5a48ed"},"cell_type":"code","source":"# create dataframe for the fastai loader\ntrain_dict = {'name': train_path + train_names, 'label': train_labels}\ndf = pd.DataFrame(data=train_dict)\n# create test dataframe\ntest_names = []\nfor f in os.listdir(test_path):\n    test_names.append(test_path + f)\ndf_test = pd.DataFrame(np.asarray(test_names), columns=['name'])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"5d25df6f20f350b1e796fa4ed34e134a4f11a826","_kg_hide-input":false},"cell_type":"code","source":"# Subclass ImageList to use our own image opening function\nclass MyImageItemList(ImageList):\n    def open(self, fn:PathOrStr)->Image:\n        img = readCroppedImage(fn.replace('/./','').replace('//','/'))\n        # This ndarray image has to be converted to tensor before passing on as fastai Image, we can use pil2tensor\n        return vision.Image(px=pil2tensor(img, np.float32))\n    ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"e82294a7507e48b0481cfb67e3e7c93e9dfc1821","_kg_hide-input":false},"cell_type":"code","source":"# Create ImageDataBunch using fastai data block API\nimgDataBunch = (MyImageItemList.from_df(path='/', df=df, suffix='.tif')\n        #Where to find the data?\n        .split_by_idx(val_idx)\n        #How to split in train/valid?\n        .label_from_df(cols='label')\n        #Where are the labels?\n        .add_test(MyImageItemList.from_df(path='/', df=df_test))\n        #dataframe pointing to the test set?\n        .transform(tfms=[[],[]], size=sz)\n        # We have our custom transformations implemented in the image loader but we could apply transformations also here\n        # Even though we don't apply transformations here, we set two empty lists to tfms. Train and Validation augmentations\n        .databunch(bs=BATCH_SIZE)\n        # convert to databunch\n        .normalize([tensor([0.702447, 0.546243, 0.696453]), tensor([0.238893, 0.282094, 0.216251])])\n        # Normalize with training set stats. These are means and std's of each three channel and we calculated these previously in the stats step.\n       )","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"35a9788f26ff06c92fb6592e3c8ee67db2f7415a","_kg_hide-input":true},"cell_type":"code","source":"# check that the imgDataBunch is loading our images ok\nimgDataBunch.show_batch(rows=2, figsize=(4,4))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"b03f5cf0fd15fe757e5eca23720c5cf39c8cf431"},"cell_type":"markdown","source":"### Training\nWe define a convnet learner object where we set the model architecture and our data bunch.\n[`create_cnn` docs.fast.ai](https://docs.fast.ai/vision.learner.html#create_cnn)"},{"metadata":{"trusted":true,"_uuid":"ebfd9acf7bec2c429dd8499514b9adcc5a235e44","_kg_hide-output":true},"cell_type":"code","source":"# Next, we create a convnet learner object\n# ps = dropout percentage (0-1) in the final layer\ndef getLearner():\n    return create_cnn(imgDataBunch, arch, pretrained=True, path='.', metrics=accuracy, ps=0.5, callback_fns=ShowGraph)\n\nlearner = getLearner()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"ac6b9d503b9144b04c7edd8fe033721ac7269a54"},"cell_type":"markdown","source":"### 1cycle policy\nWe will use the one cycle policy proposed by Leslie Smith, [arXiv, April 2018](https://arxiv.org/abs/1803.09820). The policy brings more disciplined approach for selecting hyperparameters such as learning rate and weight decay.** This can potentially save us a lot of time from training with suboptimal hyperparameters.** In addititon, Fastai library has implemented a training function for one cycle policy that we can use with only a few lines of code.\n\nSylvian Gugger wrote a very clear [explanation](https://sgugger.github.io/the-1cycle-policy.html) of Leslie's proposal that I recommend for reading.\n\nFirst, we find the optimal `learning rate` and `weight decay` values. The optimal lr is just before the base of the loss and before the start of divergence. It is important that the loss is still descending where we select the learning rate.\n\nAs for the `weight decay` that is the L2 penalty of the optimizer, Leslie proposes to select the largest one that will still let us train at a high learning rate so we do a small grid search with 1e-2, 1e-4 and 1e-6 weight decays."},{"metadata":{"trusted":true,"_uuid":"49c4bb9a9d899c9cfbfb483e486ae3338a8373bb","_kg_hide-input":true},"cell_type":"code","source":"# We can use lr_find with different weight decays and record all losses so that we can plot them on the same graph\n# Number of iterations is by default 100, but at this low number of itrations, there might be too much variance\n# from random sampling that makes it difficult to compare WD's. I recommend using an iteration count of at least 300 for more consistent results.\nlrs = []\nlosses = []\nwds = []\niter_count = 600\n\n# WEIGHT DECAY = 1e-6\nlearner.lr_find(wd=1e-6, num_it=iter_count)\nlrs.append(learner.recorder.lrs)\nlosses.append(learner.recorder.losses)\nwds.append('1e-6')\nlearner = getLearner() #reset learner - this gets more consistent starting conditions\n\n# WEIGHT DECAY = 1e-4\nlearner.lr_find(wd=1e-4, num_it=iter_count)\nlrs.append(learner.recorder.lrs)\nlosses.append(learner.recorder.losses)\nwds.append('1e-4')\nlearner = getLearner() #reset learner - this gets more consistent starting conditions\n\n# WEIGHT DECAY = 1e-2\nlearner.lr_find(wd=1e-2, num_it=iter_count)\nlrs.append(learner.recorder.lrs)\nlosses.append(learner.recorder.losses)\nwds.append('1e-2')\nlearner = getLearner() #reset learner","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true,"_uuid":"1d0e7276a4d7d770e920563cd2ad50118f99e7d8"},"cell_type":"code","source":"# Plot weight decays\n_, ax = plt.subplots(1,1)\nmin_y = 0.5\nmax_y = 0.55\nfor i in range(len(losses)):\n    ax.plot(lrs[i], losses[i])\n    min_y = min(np.asarray(losses[i]).min(), min_y)\nax.set_ylabel(\"Loss\")\nax.set_xlabel(\"Learning Rate\")\nax.set_xscale('log')\n#ax ranges may need some tuning with different model architectures \nax.set_xlim((1e-3,3e-1))\nax.set_ylim((min_y - 0.02,max_y))\nax.legend(wds)\nax.xaxis.set_major_formatter(plt.FormatStrFormatter('%.0e'))","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"0adfcfd595dc62d6cf4b8f77f4ba235e82e18980"},"cell_type":"markdown","source":"**We want to select the largest weight decay that gets to a low loss and has the highest learning rate before shooting up.**\nOut of the tested WD's, 1e-4 seems like the largest WD that allow us to train with maximal learning rate. *Note that the plot above may differ between runs as we use a random subset of data for computations.*\n\n> Smaller datasets and architectures seem to require larger values for weight decay while larger datasets and deeper architectures seem to require smaller values. \n[Leslie  Smith](https://arxiv.org/pdf/1803.09820.pdf)\n\nWe can select the `learning rate` around 2e-2 where it is close to the bottom but still descending.\n\nNext, we train only the heads while keeping the rest of the model frozen. Otherwise, the random initialization of the head weights could harm the relatively well-performing pre-trained weights of the model. After the heads have adjusted and the model somewhat works, we can continue to train all the weights."},{"metadata":{"trusted":true,"_uuid":"2f456a26a2e422ded89a6fdee489bec0c061d77f"},"cell_type":"code","source":"max_lr = 2e-2\nwd = 1e-4\n# 1cycle policy\nlearner.fit_one_cycle(cyc_len=8, max_lr=max_lr, wd=wd)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a25db18361929f1642b369c60468de4d86345ce2"},"cell_type":"code","source":"# plot learning rate of the one cycle\nlearner.recorder.plot_lr()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"47f7e7e9eb03a3b1f17f8474dcf7e214162f9d17"},"cell_type":"markdown","source":"We can see that the learning rate starts from lower and reaches the `max_lr` in the middle. Then it slows back down near the end. The idea is that we start with a low warm-up learning rate and gradually increase it to high. The higher rate is having a regularizing effect as it won't allow the model to settle for sharp and narrow local minima but pushes for wider and more stable one.\n\nIn the middle of our cycle, we start to lower the learning rate as we are hopefully in a good stable area. This means that we start to look for the minima within that area."},{"metadata":{"_uuid":"363bd275dbf6bf7e0c2dc2e5b245fb2eee3ea1eb"},"cell_type":"markdown","source":"## Is the model learning?\nOur model should have already learned something and if it hasn't, there's probably something wrong with our code or with the data. "},{"metadata":{"trusted":true,"_uuid":"ea77e45af04b5770f8f5ceffe54076e262394ed4"},"cell_type":"code","source":"# and plot the losses of the first cycle\nlearner.recorder.plot_losses()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"8f3a975ce2c969ae500b0f2f2f08730423057b61"},"cell_type":"markdown","source":"We can see from the plotted losses that there is a small rise after the initial drop which is caused by the increasing learning rate of the first half cycle. The losses are temporarily rising when `max_lr` drives the model out of local minima but this will pay off in the end when the learning rates are decreased."},{"metadata":{"trusted":true,"_uuid":"7f3a660d0be1146e1783aad4a7d6b0bf99bded33"},"cell_type":"code","source":"# predict the validation set with our model\ninterp = ClassificationInterpretation.from_learner(learner)\ninterp.plot_confusion_matrix(title='Confusion matrix')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"39674cab2127108d2e2ceacbd829177747dba19e"},"cell_type":"markdown","source":"**Confusion matrix** can help us understand the ratio of false negatives and positives and it's a fast way looking at our model's performance. This is a simple table that shows the counts in a way of actual label vs. predicted label. Here we can see that the model has learned to distinguish tumor and negative sample and it's already performing well. However, we will continue training further to improve from here."},{"metadata":{"trusted":true,"_uuid":"6e2326fe8a16469f9b120b6e7d19989f8f3affd9"},"cell_type":"code","source":"# before we continue, lets save the model at this stage\nlearner.save(MODEL_PATH + '_stage1')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"4b46c33f3338636560f1ffaca4fd1845e8a4b49a"},"cell_type":"markdown","source":"# Finetuning the baseline model\nNext, we can unfreeze all the trainable parameters from the model and continue its training. \n\nThe model already performs well and now, as we unfreeze the bottom layers that have been pre-trained with a large number of general images to detect common shapes and patterns, all weights are mostly adjusted. We should now train with much lower learning rates."},{"metadata":{"trusted":true,"_uuid":"68d5cd7dac6ee5fa1a41ef6a1029c18c64b709dd"},"cell_type":"code","source":"# load the baseline model\nlearner.load(MODEL_PATH + '_stage1')\n\n# unfreeze and run learning rate finder again\nlearner.unfreeze()\nlearner.lr_find(wd=wd)\n\n# plot learning rate finder results\nlearner.recorder.plot()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"d0f527d1fe80bfe55f2ba0fb29f35bd7e1c8a6d6"},"cell_type":"code","source":"# Now, smaller learning rates. This time we define the min and max lr of the cycle\nlearner.fit_one_cycle(cyc_len=12, max_lr=slice(4e-5,4e-4))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7cfcec9db08bc9ade154907800e9584153f368f3"},"cell_type":"code","source":"learner.recorder.plot_losses()","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"4595b55099a864b7bda54ee19bfbfcb799136b2a"},"cell_type":"markdown","source":"We can see that the validation performance has separated from the training performance a bit in the end of the cycle. This means that our model has started overfitting during the small learning rates. Now, if we would train further, the model would only memorize features from the training set and the validation set performance would rise. This is a good place to stop.  "},{"metadata":{"trusted":true,"_uuid":"d3dec96dfefbf8bce0653122eeaa4539c5ad7c68"},"cell_type":"code","source":"# lets take a second look at the confusion matrix. See if how much we improved.\ninterp = ClassificationInterpretation.from_learner(learner)\ninterp.plot_confusion_matrix(title='Confusion matrix')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"d218bd64e3c68e9814de4058287bee69fae61acf"},"cell_type":"code","source":"# Save the finetuned model\nlearner.save(MODEL_PATH + '_stage2')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a7ab67d367443802daef153ac452053833112001"},"cell_type":"code","source":"# if the model was better before finetuning, uncomment this to load the previous stage\n#learner.load(MODEL_PATH + '_stage1')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"c4b0274418a67c87f5d9dd538a392d3d520cc70f"},"cell_type":"markdown","source":"-------------------------\n# Validation and analysis\nNow the training is done.\n\n### How good does the model perform technically?\n\nWe can only get metrics from our validation set, and the final test metrics will be most likely a bit different."},{"metadata":{"trusted":true,"_uuid":"4891afd988159436ec176dbcca6ddcf7cfd962fa"},"cell_type":"code","source":"preds,y, loss = learner.get_preds(with_loss=True)\n# get accuracy\nacc = accuracy(preds, y)\nprint('The accuracy is {0} %.'.format(acc))","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"e8b994cb6e711b0e0702ecdcbc86168a5d93258e"},"cell_type":"markdown","source":"### How good is the model in terms of project requirements?\nIt is a good idea to look at examples of images from:\n\n- Random samples\n    - Some random predictions\n- The most incorrectly labeled\n    - What our model predicts wrong with very high probability.\n- The most correctly labeled\n    - What our model is most confident of and gets right.\n    \nThis visualization is a good way of understanding what are the images the model struggles with. It might also reveal something about the dataset such as bad quality data."},{"metadata":{"trusted":true,"_uuid":"0399dfe70ae97a156c84034c8093413af61a6334","_kg_hide-input":true},"cell_type":"code","source":"# I modified this from the fastai's plot_top_losses (https://github.com/fastai/fastai/blob/master/fastai/vision/learner.py#L114)\nfrom random import randint\n\ndef plot_overview(interp:ClassificationInterpretation, classes=['Negative','Tumor']):\n    # top losses will return all validation losses and indexes sorted by the largest first\n    tl_val,tl_idx = interp.top_losses()\n    #classes = interp.data.classes\n    fig, ax = plt.subplots(3,4, figsize=(16,12))\n    fig.suptitle('Predicted / Actual / Loss / Probability',fontsize=20)\n    # Random\n    for i in range(4):\n        random_index = randint(0,len(tl_idx))\n        idx = tl_idx[random_index]\n        im,cl = interp.data.dl(DatasetType.Valid).dataset[idx]\n        im = image2np(im.data)\n        cl = int(cl)\n        ax[0,i].imshow(im)\n        ax[0,i].set_xticks([])\n        ax[0,i].set_yticks([])\n        ax[0,i].set_title(f'{classes[interp.pred_class[idx]]} / {classes[cl]} / {interp.losses[idx]:.2f} / {interp.probs[idx][cl]:.2f}')\n    ax[0,0].set_ylabel('Random samples', fontsize=16, rotation=0, labelpad=80)\n    # Most incorrect or top losses\n    for i in range(4):\n        idx = tl_idx[i]\n        im,cl = interp.data.dl(DatasetType.Valid).dataset[idx]\n        cl = int(cl)\n        im = image2np(im.data)\n        ax[1,i].imshow(im)\n        ax[1,i].set_xticks([])\n        ax[1,i].set_yticks([])\n        ax[1,i].set_title(f'{classes[interp.pred_class[idx]]} / {classes[cl]} / {interp.losses[idx]:.2f} / {interp.probs[idx][cl]:.2f}')\n    ax[1,0].set_ylabel('Most incorrect\\nsamples', fontsize=16, rotation=0, labelpad=80)\n    # Most correct or least losses\n    for i in range(4):\n        idx = tl_idx[len(tl_idx) - i - 1]\n        im,cl = interp.data.dl(DatasetType.Valid).dataset[idx]\n        cl = int(cl)\n        im = image2np(im.data)\n        ax[2,i].imshow(im)\n        ax[2,i].set_xticks([])\n        ax[2,i].set_yticks([])\n        ax[2,i].set_title(f'{classes[interp.pred_class[idx]]} / {classes[cl]} / {interp.losses[idx]:.2f} / {interp.probs[idx][cl]:.2f}')\n    ax[2,0].set_ylabel('Most correct\\nsamples', fontsize=16, rotation=0, labelpad=80)","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":false,"_kg_hide-output":false,"trusted":true,"_uuid":"adcd722a33466b0a40a193775026d4f86d1b1129"},"cell_type":"code","source":"#interp = ClassificationInterpretation.from_learner(learner)\nplot_overview(interp, ['Negative','Tumor'])","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"c270a3545fd39bc6f7e02bfda971e59d97b4bb46"},"cell_type":"markdown","source":"### Gradient-weighted Class Activation Mapping (Grad-CAM)\n[Grad-CAM: Visual Explanations from Deep Networks via Gradient-based Localization](https://arxiv.org/abs/1610.02391)\n\nThis method produces a coarse localization map highlighting the areas that the model considers important for the classification decision. The visual explanation gives transparency to the model making it easier to notice if it has learned the wrong things. For example, if we would train a dog breed classifier and all of our images of a certain dog breed would have been taken in a dog show competition. There is a good chance that the model would learn to recognize the competition surroundings instead of doggy features with that breed. Visualizing the localization map would reveal that, and we could focus on getting more diverse data of that breed.  "},{"metadata":{"_kg_hide-input":true,"trusted":true,"_uuid":"34e57f42859587eb528c0c3d07c6f2ee88f58b52"},"cell_type":"code","source":"from fastai.callbacks.hooks import *\n\n# hook into forward pass\ndef hooked_backward(m, oneBatch, cat):\n    # we hook into the convolutional part = m[0] of the model\n    with hook_output(m[0]) as hook_a: \n        with hook_output(m[0], grad=True) as hook_g:\n            preds = m(oneBatch)\n            preds[0,int(cat)].backward()\n    return hook_a,hook_g","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true,"_uuid":"1612cda873b56721486951e702d3271be5a55663"},"cell_type":"code","source":"# We can create a utility function for getting a validation image with an activation map\ndef getHeatmap(val_index):\n    \"\"\"Returns the validation set image and the activation map\"\"\"\n    # this gets the model\n    m = learner.model.eval()\n    tensorImg,cl = imgDataBunch.valid_ds[val_index]\n    # create a batch from the one image\n    oneBatch,_ = imgDataBunch.one_item(tensorImg)\n    oneBatch_im = vision.Image(imgDataBunch.denorm(oneBatch)[0])\n    # convert batch tensor image to grayscale image with opencv\n    cvIm = cv2.cvtColor(image2np(oneBatch_im.data), cv2.COLOR_RGB2GRAY)\n    # attach hooks\n    hook_a,hook_g = hooked_backward(m, oneBatch, cl)\n    # get convolutional activations and average from channels\n    acts = hook_a.stored[0].cpu()\n    #avg_acts = acts.mean(0)\n\n    # Grad-CAM\n    grad = hook_g.stored[0][0].cpu()\n    grad_chan = grad.mean(1).mean(1)\n    grad.shape,grad_chan.shape\n    mult = (acts*grad_chan[...,None,None]).mean(0)\n    return mult, cvIm","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true,"_uuid":"76c4fa90fa4c61e88b122f5af14183135d63da10"},"cell_type":"code","source":"# Then, modify our plotting func a bit\ndef plot_heatmap_overview(interp:ClassificationInterpretation, classes=['Negative','Tumor']):\n    # top losses will return all validation losses and indexes sorted by the largest first\n    tl_val,tl_idx = interp.top_losses()\n    #classes = interp.data.classes\n    fig, ax = plt.subplots(3,4, figsize=(16,12))\n    fig.suptitle('Grad-CAM\\nPredicted / Actual / Loss / Probability',fontsize=20)\n    # Random\n    for i in range(4):\n        random_index = randint(0,len(tl_idx))\n        idx = tl_idx[random_index]\n        act, im = getHeatmap(idx)\n        H,W = im.shape\n        _,cl = interp.data.dl(DatasetType.Valid).dataset[idx]\n        cl = int(cl)\n        ax[0,i].imshow(im)\n        ax[0,i].imshow(im, cmap=plt.cm.gray)\n        ax[0,i].imshow(act, alpha=0.5, extent=(0,H,W,0),\n              interpolation='bilinear', cmap='inferno')\n        ax[0,i].set_xticks([])\n        ax[0,i].set_yticks([])\n        ax[0,i].set_title(f'{classes[interp.pred_class[idx]]} / {classes[cl]} / {interp.losses[idx]:.2f} / {interp.probs[idx][cl]:.2f}')\n    ax[0,0].set_ylabel('Random samples', fontsize=16, rotation=0, labelpad=80)\n    # Most incorrect or top losses\n    for i in range(4):\n        idx = tl_idx[i]\n        act, im = getHeatmap(idx)\n        H,W = im.shape\n        _,cl = interp.data.dl(DatasetType.Valid).dataset[idx]\n        cl = int(cl)\n        ax[1,i].imshow(im)\n        ax[1,i].imshow(im, cmap=plt.cm.gray)\n        ax[1,i].imshow(act, alpha=0.5, extent=(0,H,W,0),\n              interpolation='bilinear', cmap='inferno')\n        ax[1,i].set_xticks([])\n        ax[1,i].set_yticks([])\n        ax[1,i].set_title(f'{classes[interp.pred_class[idx]]} / {classes[cl]} / {interp.losses[idx]:.2f} / {interp.probs[idx][cl]:.2f}')\n    ax[1,0].set_ylabel('Most incorrect\\nsamples', fontsize=16, rotation=0, labelpad=80)\n    # Most correct or least losses\n    for i in range(4):\n        idx = tl_idx[len(tl_idx) - i - 1]\n        act, im = getHeatmap(idx)\n        H,W = im.shape\n        _,cl = interp.data.dl(DatasetType.Valid).dataset[idx]\n        cl = int(cl)\n        ax[2,i].imshow(im)\n        ax[2,i].imshow(im, cmap=plt.cm.gray)\n        ax[2,i].imshow(act, alpha=0.5, extent=(0,H,W,0),\n              interpolation='bilinear', cmap='inferno')\n        ax[2,i].set_xticks([])\n        ax[2,i].set_yticks([])\n        ax[2,i].set_title(f'{classes[interp.pred_class[idx]]} / {classes[cl]} / {interp.losses[idx]:.2f} / {interp.probs[idx][cl]:.2f}')\n    ax[2,0].set_ylabel('Most correct\\nsamples', fontsize=16, rotation=0, labelpad=80)","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":false,"trusted":true,"_uuid":"26965128db3010f2e43dd0875fed848601411265"},"cell_type":"code","source":"plot_heatmap_overview(interp, ['Negative','Tumor'])","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"e50ab570e1cc77979c082c8da561f93f6669f7aa"},"cell_type":"markdown","source":"This shows the activation maps of the predicted category so if the label is ```tumor```, the visualization shows all the places where the model thinks the tumor patterns are."},{"metadata":{"_uuid":"1afca3014ca68714928c96c578b6bd1e4f820fc8"},"cell_type":"markdown","source":"### ROC curve and AUC\nRemember, AUC is the metric that is used for evaluating submissions. We can calculate it here for ou validation set but it will most likely differ from the final score."},{"metadata":{"trusted":true,"_uuid":"02b156c4b335575ac8b9e64e7c9ea8d2a1dead34"},"cell_type":"code","source":"from sklearn.metrics import roc_curve, auc\n# probs from log preds\nprobs = np.exp(preds[:,1])\n# Compute ROC curve\nfpr, tpr, thresholds = roc_curve(y, probs, pos_label=1)\n\n# Compute ROC area\nroc_auc = auc(fpr, tpr)\nprint('ROC area is {0}'.format(roc_auc))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"c6f6818925d9356185761f12f12594f3decd8ded"},"cell_type":"code","source":"plt.figure()\nplt.plot(fpr, tpr, color='darkorange', label='ROC curve (area = %0.2f)' % roc_auc)\nplt.plot([0, 1], [0, 1], color='navy', linestyle='--')\nplt.xlim([-0.01, 1.0])\nplt.ylim([0.0, 1.01])\nplt.xlabel('False Positive Rate')\nplt.ylabel('True Positive Rate')\nplt.title('Receiver operating characteristic')\nplt.legend(loc=\"lower right\")","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"54efc207cf58612e6c27ba69cf7fbcd987a28ce3"},"cell_type":"markdown","source":"----------------\n\n# Submit predictions\n### TTA\nTo evaluate the model, we run inference on all test images. As we have test time augmentation, our results will probably improve if we do predictions multiple times per image and average out the results. "},{"metadata":{"trusted":true,"_uuid":"5b9d4a6eee767bbed3f0ca889ea1ce4c4fffc607"},"cell_type":"code","source":"# make sure we have the best performing model stage loaded\nlearner.load(MODEL_PATH + '_stage2')\n\n# Fastai has a function for this but we don't want the additional augmentations it does (our image loader has augmentations) so we just use the get_preds\n#preds_test,y_test=learner.TTA(ds_type=DatasetType.Test)\n\n# We do a fair number of iterations to cover different combinations of flips and rotations.\n# The predictions are then averaged.\nn_aug = 12\npreds_n_avg = np.zeros((len(learner.data.test_ds.items),2))\nfor n in tqdm_notebook(range(n_aug), 'Running TTA...'):\n    preds,y = learner.get_preds(ds_type=DatasetType.Test, with_loss=False)\n    preds_n_avg = np.sum([preds_n_avg, preds.numpy()], axis=0)\npreds_n_avg = preds_n_avg / n_aug","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"45a59c0e64fc20d1e3a3362f0479b43bf1cbbb3b"},"cell_type":"code","source":"# Next, we will transform class probabilities to just tumor class probabilities\nprint('Negative and Tumor Probabilities: ' + str(preds_n_avg[0]))\ntumor_preds = preds_n_avg[:, 1]\nprint('Tumor probability: ' + str(tumor_preds[0]))\n# If we wanted to get the predicted class, argmax would get the index of the max\nclass_preds = np.argmax(preds_n_avg, axis=1)\nclasses = ['Negative','Tumor']\nprint('Class prediction: ' + classes[class_preds[0]])","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"eaa1a9651dff68793b4c66431d549e7eac2659a5"},"cell_type":"markdown","source":"### Submit the model for evaluation\nWe need to submit the tumor probability for each test sample in this competition. The probability ranges from 0 to 1."},{"metadata":{"trusted":true,"_uuid":"453af28df5ec8db3bf01b51febe332af3622d37c"},"cell_type":"code","source":"# get test id's from the sample_submission.csv and keep their original order\nSAMPLE_SUB = '/kaggle/input/sample_submission.csv'\nsample_df = pd.read_csv(SAMPLE_SUB)\nsample_list = list(sample_df.id)\n\n# List of tumor preds. \n# These are in the order of our test dataset and not necessarily in the same order as in sample_submission\npred_list = [p for p in tumor_preds]\n\n# To know the id's, we create a dict of id:pred\npred_dic = dict((key, value) for (key, value) in zip(learner.data.test_ds.items, pred_list))\n\n# Now, we can create a new list with the same order as in sample_submission\npred_list_cor = [pred_dic['///kaggle/input/test/' + id + '.tif'] for id in sample_list]\n\n# Next, a Pandas dataframe with id and label columns.\ndf_sub = pd.DataFrame({'id':sample_list,'label':pred_list_cor})\n\n# Export to csv\ndf_sub.to_csv('{0}_submission.csv'.format(MODEL_PATH), header=True, index=False)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"377aecc88035df95f2a86d579110f02d2f084573"},"cell_type":"code","source":"# This is what the first 10 items of submission look like\ndf_sub.head(10)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"3769b99bcc4e7515ed7bd441db8c934354e88883"},"cell_type":"markdown","source":"--------------------------------------------------------------\n\n# Deploy (example)\nNow that we have a working model. We could deploy this for inference to another machine, a web server for example.\nFor this, we need our saved model, and then we need to export the  DataBunch.\n\n**Note**. I haven't tried these so I cannot guarantee these will work straight away. You may need to do some debugging. [Official documentation](https://docs.fast.ai/tutorial.inference.html)."},{"metadata":{"trusted":true,"_uuid":"92c565be11aaabe6e2b0c1049415288f7f60f2f7"},"cell_type":"code","source":"# This will create an export.pkl file that you'll need to copy with your model file if you want to deploy it on another device.\n# This saves the internal information (classes, etc) need for inference in a file named 'export.pkl'. \nimgDataBunch.export(fname='./export.pkl')","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"25799c1f6ee7e9b42112deae3029f6881a7dd28d"},"cell_type":"markdown","source":"**On a new machine**\nWe need to create an empty DataBunch and load it to a learner."},{"metadata":{"trusted":true,"_uuid":"aca091a14f2693afb811dd22972952f1f87418ab"},"cell_type":"code","source":"######## RUN THIS ON A NEW MACHINE ##########\n#from fastai.vision import * # fastai 1.0\n#from fastai import *\n#from torchvision.models import *\n#arch = densenet169       # specify model architecture\n#MODEL_PATH = str(arch).split()[1] + '_stage2'\n#empty_data = ImageDataBunch.load_empty('./') #this will look for a file named export.pkl in the specified path\n#learner = create_cnn(empty_data, arch).load(MODEL_PATH)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"5c371cd36577f754b06a65d57a9375e9a10d611c"},"cell_type":"code","source":"## And then we are ready to do predictions\nimport cv2\nsz = 68\n\n# This function will convert image to the prediction format\ndef imageToTensorImage(path):\n    bgr_img = cv2.imread(path)\n    b,g,r = cv2.split(bgr_img)\n    rgb_img = cv2.merge([r,g,b])\n    # crop to center to the correct size and convert from 0-255 range to 0-1 range\n    H,W,C = rgb_img.shape\n    rgb_img = rgb_img[(H-sz)//2:(sz +(H-sz)//2),(H-sz)//2:(sz +(H-sz)//2),:] / 256\n    return vision.Image(px=pil2tensor(rgb_img, np.float32))\n\nimg = imageToTensorImage('/kaggle/input/test/0eb051700fb6b1bf96188f36c8e4889598c6a157.tif')\n\n## predict and visualize\nimg.show(y=learner.predict(img)[0])\nclasses = ['negative', 'tumor']\nprint('This is a ' + classes[int(learner.predict(img)[0])] + ' tissue image.')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"c6d6eefdbeba2be86625abce612a52e7ac19375d"},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat":4,"nbformat_minor":4}