{"cells": [{"outputs": [], "cell_type": "markdown", "source": "This script shows the full training and prediction pipeline for a pixel-based classifier: we create a mask, train logistic regression on one-pixel patches, make prediction for all pixels, create and smooth polygons from pixels.", "metadata": {"_cell_guid": "9e340962-ed4a-8e0f-6877-f504cffc25ab", "_uuid": "141bfbfc4b1a43804c6417da796b15db2e3e08ff"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "from collections import defaultdict\nimport csv\nimport sys\n\nimport cv2\nfrom shapely.geometry import MultiPolygon, Polygon\nimport shapely.wkt\nimport shapely.affinity\nimport numpy as np\nimport tifffile as tiff\n\ncsv.field_size_limit(sys.maxsize);", "metadata": {"trusted": false, "_cell_guid": "e388c17e-209b-97bb-a4d1-5ee5b8a9dbdc", "_uuid": "5903ca2e3119d2d773b2421faac1717680200880", "_execution_state": "busy"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "We'll work on buildings (class 1) from image 6120_2_2. Fist load grid sizes and polygons.", "metadata": {"_cell_guid": "f419e71d-3ae0-9cdf-b10c-88049fafc541", "_uuid": "e35f2e95601cd6bbc1681dca5475bbb89c8707d6"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "IM_ID = '6120_2_2'\nPOLY_TYPE = '1'  # buildings\n\n# Load grid size\nx_max = y_min = None\nfor _im_id, _x, _y in csv.reader(open('../input/grid_sizes.csv')):\n    if _im_id == IM_ID:\n        x_max, y_min = float(_x), float(_y)\n        break\nprint (x_max)\nprint (y_min)\n# Load train poly with shapely\ntrain_polygons = None\nfor _im_id, _poly_type, _poly in csv.reader(open('../input/train_wkt_v4.csv')):\n    if _im_id == IM_ID and _poly_type == POLY_TYPE:\n        train_polygons = shapely.wkt.loads(_poly)\n        break\n\n# Read image with tiff\nim_rgb = tiff.imread('../input/three_band/{}.tif'.format(IM_ID)).transpose([1, 2, 0])\nim_size = im_rgb.shape[:2]\nprint (im_size)", "metadata": {"trusted": false, "_cell_guid": "8ba03b67-e2a0-d855-2fe6-ae2f8e4ad897", "_uuid": "fc1188d2e5990f048b00f6c911b2951652b960db", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "Scale polygons to match image:", "metadata": {"_cell_guid": "738a815a-a1ee-aed4-19e5-a6f745458bbf", "_uuid": "57f6b84a793cd0d168d988870cef6e149cf2298e"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "def get_scalers():\n    h, w = im_size  # they are flipped so that mask_for_polygons works correctly\n    w_ = w * (w / (w + 1))\n    h_ = h * (h / (h + 1))\n    return w_ / x_max, h_ / y_min\n\nx_scaler, y_scaler = get_scalers()\n\ntrain_polygons_scaled = shapely.affinity.scale(\n    train_polygons, xfact=x_scaler, yfact=y_scaler, origin=(0, 0, 0))", "metadata": {"trusted": false, "_cell_guid": "9352582b-4da8-52bd-d951-96ae2e157416", "_uuid": "452aa03ac5a4d924439e5dfba80c87eb58c477a7", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "Create a mask from polygons:", "metadata": {"_cell_guid": "8ce835a0-2e65-4685-6001-e3aa211c05fa", "_uuid": "e8074e2cc21bdd8d3b795445e443115edb6b00f3"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "def mask_for_polygons(polygons):\n    #print (polygons)\n    img_mask = np.zeros(im_size, np.uint8)\n    if not polygons:\n        return img_mask\n    int_coords = lambda x: np.array(x).round().astype(np.int32)\n    exteriors = [int_coords(poly.exterior.coords) for poly in polygons]\n    interiors = [int_coords(pi.coords) for poly in polygons\n                 for pi in poly.interiors]\n    cv2.fillPoly(img_mask, exteriors, 1)\n    cv2.fillPoly(img_mask, interiors, 0)\n    return img_mask\n\ntrain_mask = mask_for_polygons(train_polygons_scaled)", "metadata": {"trusted": false, "_cell_guid": "49b3b539-00da-2adb-4d23-eeedc896a7eb", "_uuid": "652ae0ae2bd567c56c0e0284964137b761e46c82", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "A helper for nicer display", "metadata": {"_cell_guid": "ed5d76c8-47a3-2a65-f6b4-46ef0e53af65", "_uuid": "dbd6cefeac14f3c8bf32222f72ee91279cc265b0"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "def scale_percentile(matrix):\n    w, h, d = matrix.shape\n    matrix = np.reshape(matrix, [w * h, d]).astype(np.float64)\n    # Get 2nd and 98th percentile\n    mins = np.percentile(matrix, 1, axis=0)\n    maxs = np.percentile(matrix, 99, axis=0) - mins\n    matrix = (matrix - mins[None, :]) / maxs[None, :]\n    matrix = np.reshape(matrix, [w, h, d])\n    matrix = matrix.clip(0, 1)\n    return matrix", "metadata": {"trusted": false, "_cell_guid": "78357c93-f936-a843-110a-009455a25830", "_uuid": "d7977b2bbe06759b30093f0198365ff53161f0f2", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "Check that image and mask are aligned.\nImage:", "metadata": {"_cell_guid": "be829690-6fcb-8bf9-6b75-6ca4808bd6c3", "_uuid": "5db1f9428a2c1033895c2aa7c69931fb30c046ea"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "tiff.imshow(255 * scale_percentile(im_rgb[2900:3200,2000:2300]));", "metadata": {"trusted": false, "_cell_guid": "92545b40-b831-8165-1dc7-fe395f900741", "_uuid": "7cd727dd3b982096f58674ca87dcbf42d0c11c2e", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "And mask:", "metadata": {"_cell_guid": "e8c166c3-eb24-a9a3-998e-0b60cbb7dd6b", "_uuid": "ac7fe03201a7ab3748412726842954b16598e5c3"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "def show_mask(m):\n    # hack for nice display\n    tiff.imshow(255 * np.stack([m, m, m]));\nshow_mask(train_mask[2900:3200,2000:2300])", "metadata": {"trusted": false, "_cell_guid": "74b5d9ba-fdde-671c-757e-2608cf2d356f", "_uuid": "5af4ba8da828f10043a2e14948aa0bb189cf2b64", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "Now, let's train a very simple logistic regression classifier, just to get some noisy prediction to show how output mask is processed.", "metadata": {"_cell_guid": "d7d5b159-0deb-6154-5e42-eceb73290c9a", "_uuid": "45a8e7a6ffa0352ce763804417d79805d77fa0ba"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "from sklearn.linear_model import SGDClassifier\nfrom sklearn.preprocessing import StandardScaler\nfrom sklearn.pipeline import make_pipeline\nfrom sklearn.model_selection import cross_val_predict\nfrom sklearn.metrics import average_precision_score\n\nxs = im_rgb.reshape(-1, 3).astype(np.float32)\nys = train_mask.reshape(-1)\npipeline = make_pipeline(StandardScaler(), SGDClassifier(loss='log'))\n\nprint('training...')\n# do not care about overfitting here\npipeline.fit(xs, ys)\npred_ys = pipeline.predict_proba(xs)[:, 1]\nprint('average precision', average_precision_score(ys, pred_ys))\npred_mask = pred_ys.reshape(train_mask.shape)", "metadata": {"trusted": false, "_cell_guid": "36724807-e62f-6d9e-d074-2b10bdeb5f2b", "_uuid": "5e40501a417390a40c01a3155c6ea2755d129d7c", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "Now check predictions:", "metadata": {"_cell_guid": "82444eb3-6619-c5e0-8b75-e84861fc02b8", "_uuid": "15c15fefaca5c0adb963112affecd88b10bdea64"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "show_mask(pred_mask[2900:3200,2000:2300])", "metadata": {"trusted": false, "_cell_guid": "62b059c5-3b85-ff83-844f-28b0eb27f113", "_uuid": "6b7844843e92a30754bd555ec7771619bca4b080", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "We must choose a threshold to turn it into a binary mask:", "metadata": {"_cell_guid": "1e45d1b3-ff2a-f922-9d9f-8685043aef26", "_uuid": "f0f2daeccbd0f523677095da63f237ce67163f08"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "threshold = 0.3\npred_binary_mask = pred_mask >= threshold\nshow_mask(pred_binary_mask[2900:3200,2000:2300])", "metadata": {"trusted": false, "_cell_guid": "bd6bce17-ea6a-b5b3-5e39-6973717585cf", "_uuid": "7d5cbbbffea659ff414d05f2dd203b62f138b418", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "Now it's possible to check Jaccard on the pixel level:", "metadata": {"_cell_guid": "37a76e2e-d88e-13e1-8959-9b6abc3bca09", "_uuid": "97e410571925c36e38e7480e1c65c999286e25cd"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "# check jaccard on the pixel level\ntp, fp, fn = (( pred_binary_mask &  train_mask).sum(),\n              ( pred_binary_mask & ~train_mask).sum(),\n              (~pred_binary_mask &  train_mask).sum())\nprint('Pixel jaccard', tp / (tp + fp + fn))", "metadata": {"trusted": false, "_cell_guid": "da81d13e-9555-67b3-31fe-d336c01af0d3", "_uuid": "c961ea884ab58c9cdc6980d711a850303f0c09e3", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "Next is the most interesting bit, creating polygons from bit masks. Please see inline comments:", "metadata": {"_cell_guid": "7dc1b380-3e44-9169-e6da-ae65982cd071", "_uuid": "9b05c0e92ea7939a6283c9f2c08b1fde1568b31d"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "def mask_to_polygons(mask, epsilon=10., min_area=10.):\n    # first, find contours with cv2: it's much faster than shapely\n    image, contours, hierarchy = cv2.findContours(\n        ((mask == 1) * 255).astype(np.uint8),\n        cv2.RETR_CCOMP, cv2.CHAIN_APPROX_TC89_KCOS)\n    # create approximate contours to have reasonable submission size\n    approx_contours = [cv2.approxPolyDP(cnt, epsilon, True)\n                       for cnt in contours]\n    if not contours:\n        return MultiPolygon()\n    # now messy stuff to associate parent and child contours\n    cnt_children = defaultdict(list)\n    child_contours = set()\n    assert hierarchy.shape[0] == 1\n    # http://docs.opencv.org/3.1.0/d9/d8b/tutorial_py_contours_hierarchy.html\n    for idx, (_, _, _, parent_idx) in enumerate(hierarchy[0]):\n        if parent_idx != -1:\n            child_contours.add(idx)\n            cnt_children[parent_idx].append(approx_contours[idx])\n    # create actual polygons filtering by area (removes artifacts)\n    all_polygons = []\n    for idx, cnt in enumerate(approx_contours):\n        if idx not in child_contours and cv2.contourArea(cnt) >= min_area:\n            assert cnt.shape[1] == 1\n            poly = Polygon(\n                shell=cnt[:, 0, :],\n                holes=[c[:, 0, :] for c in cnt_children.get(idx, [])\n                       if cv2.contourArea(c) >= min_area])\n            all_polygons.append(poly)\n    # approximating polygons might have created invalid ones, fix them\n    all_polygons = MultiPolygon(all_polygons)\n    if not all_polygons.is_valid:\n        all_polygons = all_polygons.buffer(0)\n        # Sometimes buffer() converts a simple Multipolygon to just a Polygon,\n        # need to keep it a Multi throughout\n        if all_polygons.type == 'Polygon':\n            all_polygons = MultiPolygon([all_polygons])\n    return all_polygons", "metadata": {"trusted": false, "_cell_guid": "ca144939-f1bd-c3c1-83c4-56d0febf7695", "_uuid": "8baeb4f16d6d60d3d7256e72951cf5870a3de46f", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "Turn our prediction to polygons, and then turn back into a mask to check what it looks like:", "metadata": {"_cell_guid": "d26dc686-a0f6-979a-4cf6-1dd00541766f", "_uuid": "23a009ee5f1aba1a073799bf4a3440b02a83ba68"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "pred_polygons = mask_to_polygons(pred_binary_mask)\npred_poly_mask = mask_for_polygons(pred_polygons)\nshow_mask(pred_poly_mask[2900:3200,2000:2300])", "metadata": {"trusted": false, "_cell_guid": "f8e31a8f-b2e3-7f84-f733-b0b8602520df", "_uuid": "c0474df9455f128b278b930b12a9ee1b6f9402b3", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "Now to create a submission we just scale back to original coordinates", "metadata": {"_cell_guid": "d18792ef-7138-4fc0-6ff9-1dc3f860f104", "_uuid": "7d7ae3b4b242ab8218f148e3894b756019efe062"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "scaled_pred_polygons = shapely.affinity.scale(\n    pred_polygons, xfact=1 / x_scaler, yfact=1 / y_scaler, origin=(0, 0, 0))", "metadata": {"trusted": false, "_cell_guid": "6f9f1e0b-465b-2e86-f059-1b76e5505c49", "_uuid": "a8b3c8d3548812b38078cb3fccc71100b40230c7", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "Checking submission size:", "metadata": {"_cell_guid": "9d066236-0de6-d722-8247-129a6bf32661", "_uuid": "11576b36a4ad77a3a257201c5dcf0f3dfd7fc446"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "dumped_prediction = shapely.wkt.dumps(scaled_pred_polygons)\nprint('Prediction size: {:,} bytes'.format(len(dumped_prediction)))\nfinal_polygons = shapely.wkt.loads(dumped_prediction)", "metadata": {"trusted": false, "_cell_guid": "7ffd7429-1931-5a7f-acd3-8f5c314b38b6", "_uuid": "b80123756e792a2b6713e88bb3f7b283b420cc56", "_execution_state": "idle"}, "execution_count": null}, {"outputs": [], "cell_type": "markdown", "source": "Now the litmus test: check Jaccard compared to **original** polygons\n", "metadata": {"_cell_guid": "27cbb32c-a33c-820f-9279-0a70b54f6b99", "_uuid": "ab4a9f5edc0595489d4d3b378bee02bc5f844876"}, "execution_count": null}, {"outputs": [], "cell_type": "code", "source": "print('Final jaccard',\n      final_polygons.intersection(train_polygons).area /\n      final_polygons.union(train_polygons).area)", "metadata": {"trusted": false, "_cell_guid": "fa3b414c-aa0c-29d6-50c5-a2bae13e4414", "_uuid": "ff74fade950284353c674bab317a0c2c5eaad7c1", "_execution_state": "idle"}, "execution_count": null}], "metadata": {"language_info": {"nbconvert_exporter": "python", "file_extension": ".py", "version": "3.6.1", "name": "python", "mimetype": "text/x-python", "codemirror_mode": {"name": "ipython", "version": 3}, "pygments_lexer": "ipython3"}, "_is_fork": false, "kernelspec": {"language": "python", "name": "python3", "display_name": "Python 3"}, "_change_revision": 0}, "nbformat": 4, "nbformat_minor": 0}