{"cells":[{"metadata":{"_cell_guid":"c5e1d71f-2c71-3d3a-53f8-e16bbfedc215","_uuid":"e015bc7b3d6c08016c2fcbe93832dd08e4c366d0","trusted":true},"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load in \n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\n# Input data files are available in the \"../input/\" directory.\n# For example, running this (by clicking run or pressing Shift+Enter) will list the files in the input directory\n\nfrom subprocess import check_output\nprint(check_output([\"ls\", \"../input/\"]).decode(\"utf8\"))\n\nimport matplotlib.pyplot as plt # showing and rendering figures\n# io related\nfrom skimage.io import imread\nimport os\nfrom glob import glob\n# not needed in Kaggle, but required in Jupyter\n%matplotlib inline \n\nimport tifffile as tiff\nimport cv2 as cv2\nfrom skimage.segmentation import slic, mark_boundaries\n\n# Any results you write to the current directory are saved as output.","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a8eb3828afc6ffac746767a058b5535b42290de2"},"cell_type":"code","source":"#base_image_dir = os.path.join('..', 'input', 'diabetic-retinopathy-detection')\nbase_image_dir = os.path.join('..', 'input')\nretina_df = pd.read_csv(os.path.join(base_image_dir, 'trainLabels.csv'))\nretina_df['PatientId'] = retina_df['image'].map(lambda x: x.split('_')[0])\nretina_df['path'] = retina_df['image'].map(lambda x: os.path.join(base_image_dir,\n                                                         '{}.jpeg'.format(x)))\nretina_df['exists'] = retina_df['path'].map(os.path.exists)\nprint(retina_df['exists'].sum(), 'images found of', retina_df.shape[0], 'total')\nretina_df['eye'] = retina_df['image'].map(lambda x: 1 if x.split('_')[-1]=='left' else 0)\nfrom keras.utils.np_utils import to_categorical\nretina_df['level_cat'] = retina_df['level'].map(lambda x: to_categorical(x, 1+retina_df['level'].max()))\n\nretina_df.dropna(inplace = True)\nretina_df = retina_df[retina_df['exists']]\nretina_df.sample(3)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"85f43db3a0159c4420233d71bd99893c14815b6b"},"cell_type":"code","source":"from sklearn.model_selection import train_test_split\nrr_df = retina_df[['PatientId', 'level']].drop_duplicates()\ntrain_ids, valid_ids = train_test_split(rr_df['PatientId'], \n                                   test_size = 0.25, \n                                   random_state = 2018,\n                                   stratify = rr_df['level'])\nraw_train_df = retina_df[retina_df['PatientId'].isin(train_ids)]\nvalid_df = retina_df[retina_df['PatientId'].isin(valid_ids)]\nprint('train', raw_train_df.shape[0], 'validation', valid_df.shape[0])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"a743a022a4f8a0dd4c9c03be4430a67f6ed6f10e"},"cell_type":"code","source":"#balance the train set:\ntrain_df = raw_train_df.groupby(['level', 'eye']).apply(lambda x: x.sample(75, replace = True)\n                                                      ).reset_index(drop = True)\nprint('New Data Size:', train_df.shape[0], 'Old Size:', raw_train_df.shape[0])\ntrain_df[['level', 'eye']].hist(figsize = (10, 5))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"3f156c03432c82aca3868a439a32a7fb451ec15e"},"cell_type":"code","source":"import tensorflow as tf\nfrom keras import backend as K\nfrom keras.applications.inception_v3 import preprocess_input\nimport numpy as np\nIMG_SIZE = (512, 512) # slightly smaller than vgg16 normally expects\ndef tf_image_loader(out_size, \n                      horizontal_flip = True, \n                      vertical_flip = False, \n                     random_brightness = True,\n                     random_contrast = True,\n                    random_saturation = True,\n                    random_hue = True,\n                      color_mode = 'rgb',\n                       preproc_func = preprocess_input,\n                       on_batch = False):\n    def _func(X):\n        with tf.name_scope('image_augmentation'):\n            with tf.name_scope('input'):\n                X = tf.image.decode_png(tf.read_file(X), channels = 3 if color_mode == 'rgb' else 0)\n                X = tf.image.resize_images(X, out_size)\n            with tf.name_scope('augmentation'):\n                if horizontal_flip:\n                    X = tf.image.random_flip_left_right(X)\n                if vertical_flip:\n                    X = tf.image.random_flip_up_down(X)\n                if random_brightness:\n                    X = tf.image.random_brightness(X, max_delta = 0.1)\n                if random_saturation:\n                    X = tf.image.random_saturation(X, lower = 0.75, upper = 1.5)\n                if random_hue:\n                    X = tf.image.random_hue(X, max_delta = 0.15)\n                if random_contrast:\n                    X = tf.image.random_contrast(X, lower = 0.75, upper = 1.5)\n                return preproc_func(X)\n    if on_batch: \n        # we are meant to use it on a batch\n        def _batch_func(X, y):\n            return tf.map_fn(_func, X), y\n        return _batch_func\n    else:\n        # we apply it to everything\n        def _all_func(X, y):\n            return _func(X), y         \n        return _all_func\n    \ndef tf_augmentor(out_size,\n                intermediate_size = (640, 640),\n                 intermediate_trans = 'crop',\n                 batch_size = 16,\n                   horizontal_flip = True, \n                  vertical_flip = False, \n                 random_brightness = True,\n                 random_contrast = True,\n                 random_saturation = True,\n                    random_hue = True,\n                  color_mode = 'rgb',\n                   preproc_func = preprocess_input,\n                   min_crop_percent = 0.001,\n                   max_crop_percent = 0.005,\n                   crop_probability = 0.5,\n                   rotation_range = 10):\n    \n    load_ops = tf_image_loader(out_size = intermediate_size, \n                               horizontal_flip=horizontal_flip, \n                               vertical_flip=vertical_flip, \n                               random_brightness = random_brightness,\n                               random_contrast = random_contrast,\n                               random_saturation = random_saturation,\n                               random_hue = random_hue,\n                               color_mode = color_mode,\n                               preproc_func = preproc_func,\n                               on_batch=False)\n    def batch_ops(X, y):\n        batch_size = tf.shape(X)[0]\n        with tf.name_scope('transformation'):\n            # code borrowed from https://becominghuman.ai/data-augmentation-on-gpu-in-tensorflow-13d14ecf2b19\n            # The list of affine transformations that our image will go under.\n            # Every element is Nx8 tensor, where N is a batch size.\n            transforms = []\n            identity = tf.constant([1, 0, 0, 0, 1, 0, 0, 0], dtype=tf.float32)\n            if rotation_range > 0:\n                angle_rad = rotation_range / 180 * np.pi\n                angles = tf.random_uniform([batch_size], -angle_rad, angle_rad)\n                transforms += [tf.contrib.image.angles_to_projective_transforms(angles, intermediate_size[0], intermediate_size[1])]\n\n            if crop_probability > 0:\n                crop_pct = tf.random_uniform([batch_size], min_crop_percent, max_crop_percent)\n                left = tf.random_uniform([batch_size], 0, intermediate_size[0] * (1.0 - crop_pct))\n                top = tf.random_uniform([batch_size], 0, intermediate_size[1] * (1.0 - crop_pct))\n                crop_transform = tf.stack([\n                      crop_pct,\n                      tf.zeros([batch_size]), top,\n                      tf.zeros([batch_size]), crop_pct, left,\n                      tf.zeros([batch_size]),\n                      tf.zeros([batch_size])\n                  ], 1)\n                coin = tf.less(tf.random_uniform([batch_size], 0, 1.0), crop_probability)\n                transforms += [tf.where(coin, crop_transform, tf.tile(tf.expand_dims(identity, 0), [batch_size, 1]))]\n            if len(transforms)>0:\n                X = tf.contrib.image.transform(X,\n                      tf.contrib.image.compose_transforms(*transforms),\n                      interpolation='BILINEAR') # or 'NEAREST'\n            if intermediate_trans=='scale':\n                X = tf.image.resize_images(X, out_size)\n            elif intermediate_trans=='crop':\n                X = tf.image.resize_image_with_crop_or_pad(X, out_size[0], out_size[1])\n            else:\n                raise ValueError('Invalid Operation {}'.format(intermediate_trans))\n            return X, y\n    def _create_pipeline(in_ds):\n        batch_ds = in_ds.map(load_ops, num_parallel_calls=4).batch(batch_size)\n        return batch_ds.map(batch_ops)\n    return _create_pipeline","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"5e4dcf642ac9c985b776e3023c3c633d2e818c22"},"cell_type":"code","source":"def flow_from_dataframe(idg, \n                        in_df, \n                        path_col,\n                        y_col, \n                        shuffle = True, \n                        color_mode = 'rgb'):\n    files_ds = tf.data.Dataset.from_tensor_slices((in_df[path_col].values, \n                                                   np.stack(in_df[y_col].values,0)))\n    in_len = in_df[path_col].values.shape[0]\n    while True:\n        if shuffle:\n            files_ds = files_ds.shuffle(in_len) # shuffle the whole dataset\n        \n        next_batch = idg(files_ds).repeat().make_one_shot_iterator().get_next()\n        for i in range(max(in_len//32,1)):\n            # NOTE: if we loop here it is 'thread-safe-ish' if we loop on the outside it is completely unsafe\n            yield K.get_session().run(next_batch)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"4646d20a78adc00a1a8074c4644245145ef9c5e9"},"cell_type":"code","source":"batch_size = 48\ncore_idg = tf_augmentor(out_size = IMG_SIZE, \n                        color_mode = 'rgb', \n                        vertical_flip = True,\n                        crop_probability=0.0, # crop doesn't work yet\n                        batch_size = batch_size) \nvalid_idg = tf_augmentor(out_size = IMG_SIZE, color_mode = 'rgb', \n                         crop_probability=0.0, \n                         horizontal_flip = False, \n                         vertical_flip = False, \n                         random_brightness = False,\n                         random_contrast = False,\n                         random_saturation = False,\n                         random_hue = False,\n                         rotation_range = 0,\n                        batch_size = batch_size)\n\ntrain_gen = flow_from_dataframe(core_idg, train_df, \n                             path_col = 'path',\n                            y_col = 'level_cat')\n\nvalid_gen = flow_from_dataframe(valid_idg, valid_df, \n                             path_col = 'path',\n                            y_col = 'level_cat') # we can use much larger batches for evaluation","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"1215778863d4b4f16b0c31bfe91143b5c92c2498"},"cell_type":"code","source":"#display images from validation set\nt_x, t_y = next(valid_gen)\nfig, m_axs = plt.subplots(2, 4, figsize = (16, 8))\nfor (c_x, c_y, c_ax) in zip(t_x, t_y, m_axs.flatten()):\n    c_ax.imshow(np.clip(c_x*127+127, 0, 255).astype(np.uint8))\n    c_ax.set_title('Severity {}'.format(np.argmax(c_y, -1)))\n    c_ax.axis('off')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"954fb65cf97e82a854233d9daabe1312f869a607"},"cell_type":"code","source":"#we focus on a healthy eye and display it alone\nprint('t_x shape: ', t_x.shape)\nfig, m_axs = plt.subplots(2, 1, figsize = (16, 8))\nfor (c_x, c_y, c_ax) in zip(t_x, t_y, m_axs.flatten()):\n    c_ax.imshow(np.clip(c_x*127+127, 0, 255).astype(np.uint8))\n    c_ax.set_title('Severity {}'.format(np.argmax(c_y, -1)))\n    c_ax.axis('off')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"7f8f225a72cd914ca5eb7e7daee8b905fb4e355c"},"cell_type":"code","source":"edges = cv2.Canny(np.clip(c_x*127+127, 0, 255).astype(np.uint8),70,130)\nplt.subplot(121),plt.imshow(np.clip(c_x*127+127, 0, 255).astype(np.uint8),cmap = 'gray')\nplt.title('Original Image'), plt.xticks([]), plt.yticks([])\nplt.subplot(122),plt.imshow(edges,cmap = 'gray')\nplt.title('Edge Image'), plt.xticks([]), plt.yticks([])\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"61f8b26fdd36cd8a4409b2899c22224de67e10b5"},"cell_type":"code","source":"#https://www.pyimagesearch.com/2015/04/06/zero-parameter-automatic-canny-edge-detection-with-python-and-opencv/\nimg = cv2.imread(os.path.join(base_image_dir, '15_right.jpeg'), cv2.IMREAD_COLOR)\nprint('img type ', type(img), img.shape)\nprint('img min max ', img.min(), img.max())\ngray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)\nblurred = cv2.GaussianBlur(gray, (3, 3), 0)\nv = np.median(img)\nsigma =  0.33      #0.93 \n# apply automatic Canny edge detection using the computed median\nlower = int(max(0, (1.0 - sigma) * v))\nupper = int(min(255, (1.0 + sigma) * v))\n#imgProd = img*127+127\nimgProd = gray*127      #+127\nimgClip = np.clip(imgProd, 0, 255).astype(np.uint8)\nprint('imgProd min max ', imgProd.min(), imgProd.max())\n\n# apply Canny edge detection using a wide threshold, tight\n# threshold, and automatically determined threshold\n\nwide  = cv2.Canny(blurred, 10, 200)\ntight = cv2.Canny(imgProd, 200, 250)\nauto  = cv2.Canny(imgProd, lower, upper)\n    \n#edges = cv2.Canny(np.clip(imgClip, 0, 255).astype(np.uint8), lower, upper)\n#edges = cv.Canny(img,lower,upper)\nfig, m_axs = plt.subplots(1, 1, figsize = (16, 8))\n#plt.imshow(img)    #,cmap = 'gray')\n#plt.title('Original Image'), plt.xticks([]), plt.yticks([])\n\nplt.imshow(tight,cmap = 'gray')\n#plt.title('Edge Image'), plt.xticks([]), plt.yticks([])\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"5653db613fcda2df212297e658ab1bc5db6244e8"},"cell_type":"code","source":"#another path to segmentation\ndef stretch_8bit(bands, lower_percent=2, higher_percent=98):\n    out = np.zeros_like(bands)\n    for i in range(3):\n        a = 0 #np.min(band)\n        b = 255  #np.max(band)\n        c = np.percentile(bands[:,:,i], lower_percent)\n        d = np.percentile(bands[:,:,i], higher_percent)        \n        t = a + (bands[:,:,i] - c) * (b - a) / (d - c)    \n        t[t<a] = a\n        t[t>b] = b\n        out[:,:,i] =t\n    return out.astype(np.uint8)    \n    \ndef RGB(image_id):\n    filename = os.path.join(base_image_dir,  '{}.jpeg'.format(image_id))\n    #img = tiff.imread(filename)\n    img = cv2.imread(filename, cv2.IMREAD_COLOR)\n    img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)\n    #img = np.rollaxis(img, 0, 3)    \n    return img\n\ndef GRAY(image_id):\n    filename = os.path.join(base_image_dir,  '{}.jpeg'.format(image_id))\n    img = cv2.imread(filename, cv2.IMREAD_COLOR)\n    img = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)\n    blurred = cv2.GaussianBlur(img, (3, 3), 0)\n    # invert the image\n    img = cv2.bitwise_not(blurred)\n    return img\n    \ndef M(image_id):\n    filename = os.path.join('..', 'input', 'sixteen_band', '{}_M.tif'.format(image_id))\n    img = tiff.imread(filename)    \n    img = np.rollaxis(img, 0, 3)\n    return img","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"0c82376c49a971db5103c4eb15038bdaf98f966c"},"cell_type":"code","source":"image_id = '15_right'\nrgb = GRAY(image_id)\nprint('shape rgb: ', rgb.shape, rgb.min(), rgb.max())\n#rgb1 = stretch_8bit(rgb)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"c9178914428282fef0ed2a4b0c2716b800290f9e"},"cell_type":"code","source":" y1,y2,x1,x2 =  1900, 2364, 3000, 4000\nregion = rgb[y1:y2, x1:x2]\nplt.figure()\nplt.imshow(region)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"946d7747361b211d02da469f656b7e99482b332c"},"cell_type":"code","source":"#blue_mask = cv2.inRange(region, np.array([175,70,20]), np.array([189,140,150])) #img minthresh maxthre (BGR)   \nblue_mask = cv2.inRange(region, 19, 113)  #img minthresh maxthre (BGR)   \nprint('blue_mask shape: ', blue_mask.shape, ' region shape: ', region.shape, 'min max bmask', blue_mask.min(), blue_mask.max())\nmask = cv2.bitwise_and(region, region, mask=blue_mask)\n#mask = cv2.cvtColor(mask, cv2.COLOR_RGB2GRAY) \n#print('mask shape ', mask.shape, mask.max(), region[:,:,2].max())\nplt.figure()\nplt.imshow(blue_mask, cmap='gray')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"d8bff86030fbd1c992717ac6311d6bfa0f72dbe9"},"cell_type":"code","source":"'''\nsegments = slic(region, n_segments=100, compactness=20.0, \n                    max_iter=10, sigma=5, spacing=None, multichannel=True, \n                    convert2lab=True, enforce_connectivity=False, \n                    min_size_factor=10, max_size_factor=3, slic_zero=False)\nboundaries = mark_boundaries(region, segments, color=(0,255,0))\nplt.figure()\nplt.imshow(boundaries)\n'''","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"e16ddd21e101d8cfe63006e6aff098adc12b8344"},"cell_type":"code","source":"'''\nout = np.zeros_like(mask)\nfor i in range(np.max(segments)):\n    s = segments == i\n    s_size = np.sum(s)\n    s_count = np.sum([1 for x in mask[s].ravel() if x>0])\n    #print(s_count, s_size)\n    if s_count > 0.1*s_size:\n        out[s] = 255\n        \nplt.figure()\nplt.imshow(out, cmap='gray')\n'''","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"8064caf5c027134609db2f6b5953500ec05ef021"},"cell_type":"code","source":" y1,y2,x1,x2 =  100, 3364, 500, 4000\nregion = rgb[y1:y2, x1:x2]\nplt.figure()\nplt.imshow(region)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"82bab9db7321b7b99874f5379bda38418545e20e"},"cell_type":"code","source":"#https://docs.opencv.org/3.4/d3/db4/tutorial_py_watershed.html\nret, thresh = cv2.threshold(region,125, 255,cv2.THRESH_TOZERO_INV)    #_INV+cv2.THRESH_OTSU)\n'''\nret,thresh1 = cv.threshold(img,127,255,cv.THRESH_BINARY)\nret,thresh2 = cv.threshold(img,127,255,cv.THRESH_BINARY_INV)\nret,thresh3 = cv.threshold(img,127,255,cv.THRESH_TRUNC)\nret,thresh4 = cv.threshold(img,127,255,cv.THRESH_TOZERO)\nret,thresh5 = cv.threshold(img,127,255,cv.THRESH_TOZERO_INV)\n'''\nprint('threshold ret shapes ', thresh.min(), thresh.max())\nplt.figure()\nplt.imshow(thresh)\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"d0d4463bb9a54ff8261e4875af531a24eecbb3d8"},"cell_type":"code","source":"th = cv2.adaptiveThreshold(region, 255, cv2.ADAPTIVE_THRESH_GAUSSIAN_C, cv2.THRESH_BINARY, 115, 1)\nplt.figure()\nplt.imshow(th)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_uuid":"f6637e03ac73af73bd1298fc9deec2ee5c70ec12"},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"metadata":{"_change_revision":0,"_is_fork":false,"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.6.6","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"}},"nbformat":4,"nbformat_minor":1}