{"cells":[{"metadata":{"_cell_guid":"0951d179-d004-4bff-beae-9c1784d21be8","collapsed":true,"_uuid":"c816cd20d03de3e064f746f9258afca5fb48c5a5","trusted":true},"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport random\nimport skimage.io\nimport matplotlib.pyplot as plt\nfrom skimage import transform\nimport os\nimport shutil\nfrom tqdm import tqdm\nimport tensorflow as tf\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\n# Any results you write to the current directory are saved as output.","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"11614177-8ca4-4e29-8db6-4d701aaa023e","collapsed":true,"_uuid":"f1dba6dc0cb876dcf97cc7d024a4e12039f1f566","trusted":false},"cell_type":"code","source":"def read_image_labels(image_id):\n    # most of the content in this function is taken from 'Example Metric Implementation' kernel \n    # by 'William Cukierski'\n    image_file = \"../input/stage1_train/{}/images/{}.png\".format(image_id,image_id)\n    mask_file = \"../input/stage1_train/{}/masks/*.png\".format(image_id)\n    image = skimage.io.imread(image_file)\n    masks = skimage.io.imread_collection(mask_file).concatenate()    \n    height, width, _ = image.shape\n    num_masks = masks.shape[0]\n    labels = np.zeros((height, width), np.uint16)\n    for index in range(0, num_masks):\n        labels[masks[index] > 0] = index + 1\n    return image, labels\n\ndef data_aug(image,label,angel=30,resize_rate=0.9):\n    flip = random.randint(0, 1)\n    size = image.shape[0]\n    rsize = random.randint(np.floor(resize_rate*size),size)\n    w_s = random.randint(0,size - rsize)\n    h_s = random.randint(0,size - rsize)\n    sh = random.random()/2-0.25\n    rotate_angel = random.random()/180*np.pi*angel\n    # Create Afine transform\n    afine_tf = transform.AffineTransform(shear=sh,rotation=rotate_angel)\n    # Apply transform to image data\n    image = transform.warp(image, inverse_map=afine_tf,mode='edge')\n    label = transform.warp(label, inverse_map=afine_tf,mode='edge')\n    # Randomly corpping image frame\n    image = image[w_s:w_s+size,h_s:h_s+size,:]\n    label = label[w_s:w_s+size,h_s:h_s+size]\n    # Ramdomly flip frame\n    if flip:\n        image = image[:,::-1,:]\n        label = label[:,::-1]\n    return image, label","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"50bde9d7-53ed-47bc-a7a1-879ea908b396","_uuid":"509df62377028f887d9af96cf43707f85e4b97d2"},"cell_type":"markdown","source":"Here comes an example of randomly rotate and resize an image.\n\nThe augmented sample could be not the same shape as the original  image, they should all be resize before feed in deep learning model."},{"metadata":{"_cell_guid":"b666083c-f2d4-449e-9a0d-cb5079c1dcb6","collapsed":true,"_uuid":"3d8141b5833f61547544a5b56b2229cfe77ec8c6","trusted":false},"cell_type":"code","source":"image_ids = check_output([\"ls\", \"../input/stage1_train/\"]).decode(\"utf8\").split()\nimage_id = image_ids[random.randint(0,len(image_ids))]\nimage, labels = read_image_labels(image_id)\nplt.subplot(221)\nplt.imshow(image)\nplt.subplot(222)\nplt.imshow(labels)\n\nnew_image, new_labels = data_aug(image,labels,angel=5,resize_rate=0.9)\nplt.subplot(223)\nplt.imshow(new_image)\nplt.subplot(224)\nplt.imshow(new_labels)","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"930111ce-f7e7-45d3-a770-2a1c9635f5ee","_uuid":"40b4d053997d9a2bc1c8cbaff35fb173b5b63d87"},"cell_type":"markdown","source":"The following code is to save the augmented data locally.\n It is only useful for Linux operating system"},{"metadata":{"_cell_guid":"5dccbfc0-082f-436e-8c38-6c2fd9aa0242","collapsed":true,"_uuid":"a891d486688f455b13aad0f4e14d0b246d832473","trusted":false},"cell_type":"code","source":"def make_data_augmentation(image_ids,split_num):\n    for ax_index, image_id in tqdm(enumerate(image_ids),total=len(image_ids)):\n        image,labels = read_image_labels(image_id)\n        if not os.path.exists(\"../input/stage1_train/{}/augs/\".format(image_id)):\n            os.makedirs(\"../input/stage1_train/{}/augs/\".format(image_id))\n        if not os.path.exists(\"../input/stage1_train/{}/augs_masks/\".format(image_id)):\n            os.makedirs(\"../input/stage1_train/{}/augs_masks/\".format(image_id))\n            \n        # also save the original image in augmented file \n        plt.imsave(fname=\"../input/stage1_train/{}/augs/{}.png\".format(image_id,image_id), arr = image)\n        plt.imsave(fname=\"../input/stage1_train/{}/augs_masks/{}.png\".format(image_id,image_id),arr = labels)\n\n        for i in range(split_num):\n            new_image, new_labels = data_aug(image,labels,angel=5,resize_rate=0.9)\n            aug_img_dir = \"../input/stage1_train/{}/augs/{}_{}.png\".format(image_id,image_id,i)\n            aug_mask_dir = \"../input/stage1_train/{}/augs_masks/{}_{}.png\".format(image_id,image_id,i)\n            plt.imsave(fname=aug_img_dir, arr = new_image)\n            plt.imsave(fname=aug_mask_dir,arr = new_labels)\n\ndef clean_data_augmentation(image_ids):\n    for ax_index, image_id in tqdm(enumerate(image_ids),total=len(image_ids)):\n        if os.path.exists(\"../input/stage1_train/{}/augs/\".format(image_id)):\n            shutil.rmtree(\"../input/stage1_train/{}/augs/\".format(image_id))\n        if os.path.exists(\"../input/stage1_train/{}/augs_masks/\".format(image_id)):\n            shutil.rmtree(\"../input/stage1_train/{}/augs_masks/\".format(image_id))\n\n\nimage_ids = check_output([\"ls\", \"../input/stage1_train/\"]).decode(\"utf8\").split()\nsplit_num = 10\n#make_data_augmentation(image_ids,split_num)\n#clean_data_augmentation(image_ids)","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"eb1351af-4d61-4b51-acb2-a8083e47c6d6","_uuid":"6a031f3563e1445ef1133d28d59999483ec109ee"},"cell_type":"markdown","source":"Here we also provide our tensorflow version UNet, the visualization on tensorboard is shown as below:\n![graph](http://i67.tinypic.com/105v91z.png)"},{"metadata":{"_cell_guid":"fcf9f8e0-74d1-4973-9b09-076437e5d68b","collapsed":true,"_uuid":"ea956e8c718bde7e3e2acc32037bc43a9fb09fa8","trusted":false},"cell_type":"code","source":"def get_variable(name,shape):\n    return tf.get_variable(name, shape, initializer = tf.contrib.layers.xavier_initializer())\n\ndef UNet(X):\n    ### Unit 1 ###\n    with tf.name_scope('Unit1'):\n        W1_1 =   get_variable(\"W1_1\", [3,3,3,16] )\n        Z1 = tf.nn.conv2d(X,W1_1, strides = [1,1,1,1], padding = 'SAME')\n        A1 = tf.nn.relu(Z1)\n        W1_2 =   get_variable(\"W1_2\", [3,3,16,16] )\n        Z2 = tf.nn.conv2d(A1,W1_2, strides = [1,1,1,1], padding = 'SAME')\n        A2 = tf.nn.relu(Z2) \n        P1 = tf.nn.max_pool(A2, ksize = [1,2,2,1], strides = [1,2,2,1], padding = 'SAME')\n    ### Unit 2 ###\n    with tf.name_scope('Unit2'):\n        W2_1 =   get_variable(\"W2_1\", [3,3,16,32] )\n        Z3 = tf.nn.conv2d(P1,W2_1, strides = [1,1,1,1], padding = 'SAME')\n        A3 = tf.nn.relu(Z3)\n        W2_2 =   get_variable(\"W2_2\", [3,3,32,32] )\n        Z4 = tf.nn.conv2d(A3,W2_2, strides = [1,1,1,1], padding = 'SAME')\n        A4 = tf.nn.relu(Z4) \n        P2 = tf.nn.max_pool(A4, ksize = [1,2,2,1], strides = [1,2,2,1], padding = 'SAME')\n    ### Unit 3 ###\n    with tf.name_scope('Unit3'):\n        W3_1 =   get_variable(\"W3_1\", [3,3,32,64] )\n        Z5 = tf.nn.conv2d(P2,W3_1, strides = [1,1,1,1], padding = 'SAME')\n        A5 = tf.nn.relu(Z5)\n        W3_2 =   get_variable(\"W3_2\", [3,3,64,64] )\n        Z6 = tf.nn.conv2d(A5,W3_2, strides = [1,1,1,1], padding = 'SAME')\n        A6 = tf.nn.relu(Z6) \n        P3 = tf.nn.max_pool(A6, ksize = [1,2,2,1], strides = [1,2,2,1], padding = 'SAME')\n    ### Unit 4 ###\n    with tf.name_scope('Unit4'):\n        W4_1 =   get_variable(\"W4_1\", [3,3,64,128] )\n        Z7 = tf.nn.conv2d(P3,W4_1, strides = [1,1,1,1], padding = 'SAME')\n        A7 = tf.nn.relu(Z7)\n        W4_2 =   get_variable(\"W4_2\", [3,3,128,128] )\n        Z8 = tf.nn.conv2d(A7,W4_2, strides = [1,1,1,1], padding = 'SAME')\n        A8 = tf.nn.relu(Z8) \n        P4 = tf.nn.max_pool(A8, ksize = [1,2,2,1], strides = [1,2,2,1], padding = 'SAME')\n    ### Unit 5 ###\n    with tf.name_scope('Unit5'):\n        W5_1 =   get_variable(\"W5_1\", [3,3,128,256] )\n        Z9 = tf.nn.conv2d(P4,W5_1, strides = [1,1,1,1], padding = 'SAME')\n        A9 = tf.nn.relu(Z9)\n        W5_2 =   get_variable(\"W5_2\", [3,3,256,256] )\n        Z10 = tf.nn.conv2d(A9,W5_2, strides = [1,1,1,1], padding = 'SAME')\n        A10 = tf.nn.relu(Z10) \n    ### Unit 6 ###\n    with tf.name_scope('Unit6'):\n        W6_1 =   get_variable(\"W6_1\", [3,3,256,128] )\n        U1 = tf.layers.conv2d_transpose(A10, filters = 128, kernel_size = 2, strides = 2, padding = 'SAME')\n        U1 = tf.concat([U1, A8],3)\n        W6_2 =   get_variable(\"W6_2\", [3,3,128,128] )\n        Z11 = tf.nn.conv2d(U1,W6_1, strides = [1,1,1,1], padding = 'SAME')\n        A11 = tf.nn.relu(Z11)\n        Z12 = tf.nn.conv2d(A11,W6_2, strides = [1,1,1,1], padding = 'SAME')\n        A12 = tf.nn.relu(Z12)\n    ### Unit 7 ###\n    with tf.name_scope('Unit7'):\n        W7_1 =   get_variable(\"W7_1\", [3,3,128,64] )\n        U2 = tf.layers.conv2d_transpose(A12, filters = 64, kernel_size = 2, strides = 2, padding = 'SAME')\n        U2 = tf.concat([U2, A6],3)\n        Z13 = tf.nn.conv2d(U2,W7_1, strides = [1,1,1,1], padding = 'SAME')\n        A13 = tf.nn.relu(Z13)\n        W7_2 =   get_variable(\"W7_2\", [3,3,64,64] )\n        Z14 = tf.nn.conv2d(A13,W7_2, strides = [1,1,1,1], padding = 'SAME')\n        A14 = tf.nn.relu(Z14)\n    ### Unit 8 ###\n    with tf.name_scope('Unit8'):\n        W8_1 =   get_variable(\"W8_1\", [3,3,64,32] )\n        U3 = tf.layers.conv2d_transpose(A14, filters = 32, kernel_size = 2, strides = 2, padding = 'SAME')\n        U3 = tf.concat([U3, A4],3)\n        Z15 = tf.nn.conv2d(U3,W8_1, strides = [1,1,1,1], padding = 'SAME')\n        A15 = tf.nn.relu(Z15)\n        W8_2 =   get_variable(\"W8_2\", [3,3,32,32] )\n        Z16 = tf.nn.conv2d(A15,W8_2, strides = [1,1,1,1], padding = 'SAME')\n        A16 = tf.nn.relu(Z16)\n    ### Unit 9 ###\n    with tf.name_scope('Unit9'):\n        W9_1 =   get_variable(\"W9_1\", [3,3,32,16] )\n        U4 = tf.layers.conv2d_transpose(A16, filters = 16, kernel_size = 2, strides = 2, padding = 'SAME')\n        U4 = tf.concat([U4, A2],3)\n        Z17 = tf.nn.conv2d(U4,W9_1, strides = [1,1,1,1], padding = 'SAME')\n        A17 = tf.nn.relu(Z17)\n        W9_2 =   get_variable(\"W9_2\", [3,3,16,16] )\n        Z18 = tf.nn.conv2d(A17,W9_2, strides = [1,1,1,1], padding = 'SAME')\n        A18 = tf.nn.relu(Z18)\n    ### Unit 10 ###\n    with tf.name_scope('out_put'):\n        W10 =    get_variable(\"W10\", [1,1,16,1] )\n        Z19 = tf.nn.conv2d(A18,W10, strides = [1,1,1,1], padding = 'SAME')\n        A19 = tf.nn.sigmoid(Z19)\n        Y_pred = A19\n    return Y_pred\n\ndef loss_function(y_pred, y_true):\n    cost = tf.reduce_mean(tf.keras.losses.binary_crossentropy(y_true,y_pred))\n    return cost\n\ndef mean_iou(y_pred,y_true):\n    y_pred_ = tf.to_int64(y_pred > 0.5)\n    y_true_ = tf.to_int64(y_true > 0.5)\n    score, up_opt = tf.metrics.mean_iou(y_true_, y_pred_, 2)\n    with tf.control_dependencies([up_opt]):\n        score = tf.identity(score)\n    return score","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-output":false,"_kg_hide-input":false,"_cell_guid":"5f27a81f-310e-4deb-8932-57f35a99876d","_uuid":"d8d9b8d13fb640e0fc3204ec0e41b0ba6f9fae60"},"cell_type":"markdown","source":"To train the model we can use a graph dictionary like this:"},{"metadata":{"_cell_guid":"9f0adc90-b745-4a83-82f5-92ec1bce1b67","collapsed":true,"_uuid":"2a53d94873fa9b06be579efbd909984357196827","trusted":false},"cell_type":"code","source":"# build the graph as a dictionary\ndef build_graph():\n    with tf.Graph().as_default() as g:\n        with tf.device(\"/gpu:0\"):\n            with tf.name_scope('input'):\n                x_ = tf.placeholder(tf.float32, shape=(None,IMG_HEIGHT, IMG_WIDTH, IMG_CHANNELS))\n                y_ = tf.placeholder(tf.float32, shape=(None,IMG_HEIGHT, IMG_WIDTH, 1))\n            y_pred = UNet(x_)\n            with tf.name_scope('loss'):\n                loss = loss_function(y_pred,y_)\n        with tf.device(\"/cpu:0\"):\n            with tf.name_scope(\"metrics\"):\n                iou = mean_iou(y_pred,y_)\n        model_dict = {'graph': g, 'inputs': [x_, y_],'Iou':iou,'Loss':loss, 'y_pred':y_pred}\n    return model_dict","execution_count":null,"outputs":[]},{"metadata":{"_cell_guid":"1c037766-5313-4860-9344-40b9c3414422","_uuid":"e107c483a8b1b705650d521ddab0142fce2a82f6"},"cell_type":"markdown","source":"**prediction result:**\n\nThe model use one channel label to classify if a pixel is belong to mask or background.\n\nAfter passing forward the network, the result is shown as below:\n![result](http://i64.tinypic.com/ruy7sz.png)\nFig[1],[2] is the ground truth and prediction of a vaildation image.\nFig[3],[4] is the input and the prediction on test image\n\nThe overall OUI on the validation set is near 90%"},{"metadata":{"_cell_guid":"d9c86faa-9c4b-465b-9e9b-6680480d58f7","_uuid":"d9d604bd3859629502c2e1d13a5c4f18702475c3"},"cell_type":"markdown","source":"**Three quesion remain for this network:**\n\n1. The threshold of the prediction pixcel is hard to determined. I think the OSTU adative threshold could be useful for the prediction.\n\n2. The prediction would make some very sparse pixel which should not be consider as a Nuclei. Thus a Low pass filter or Markov Random Field could be apply to expel those pixels\n\n3. How to regenerate a reliable seperated mask from the concatenated mask. the result generated by morphology from skimage is not reliable for cell that has interception.\n \nI would be happy if someone could share their points of view on the question above."}],"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}