{"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"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":22990,"databundleVersionId":2048213,"sourceType":"competition"},{"sourceId":1139701,"sourceType":"datasetVersion","datasetId":642364},{"sourceId":1383977,"sourceType":"datasetVersion","datasetId":807680},{"sourceId":1765169,"sourceType":"datasetVersion","datasetId":1049257},{"sourceId":49232793,"sourceType":"kernelVersion"},{"sourceId":56405004,"sourceType":"kernelVersion"}],"dockerImageVersionId":30066,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"

HuBMAP - Hacking the Kidney

","metadata":{}},{"cell_type":"markdown","source":"

Let’s understand the problem statement first

\n\nThe aim of this competition is to develop a segmentation problem to identify ‘Glomerulus’ in the kidney. For this we have to train a segmentation model that takes PAS kidney image input and identify the segments of glomeruli FTU in the PAS-stained microscopy data.\nFor this task, we are given historical images of kidney and annotation information representing the glomerular segmentation.\nI believe I used some terminologies, that are confusing. So, lets understand these terms and get some Domain Knowledge.\n","metadata":{}},{"cell_type":"markdown","source":"

Domain Knowledge

","metadata":{}},{"cell_type":"markdown","source":"

1. What is Glomerulus?

\nFirst, we understand what ‘Nephron’ is. The nephron is the microscopic structural and functional unit of the kidney. It is composed of a renal corpuscle and a renal tubule. The renal corpuscle consists of a tuft of capillaries called a glomerulus and an encompassing Bowman's capsule. \nGlomerulus a cluster of nerve endings, spores, or small blood vessels, in particular a cluster of capillaries around the end of a kidney tubule, where waste products are filtered from the blood.\n\n\n\nIn short, each nephron in your kidneys has a microscopic filter, called a glomerulus that is constantly filtering your blood.\n","metadata":{}},{"cell_type":"markdown","source":"

2. PAS (Periodic acid-Schiff) Stain Microscopy

\nPAS is a histology stain that detects complex sugars in tissue sections. Periodic acid is used to break specific bonds within these sugars. The resulting aldehydes react with the Schiff reagent to produce the purple-magenta color exhibited by these images. Glomeruli can be observed as the circular areas of dark stain.","metadata":{}},{"cell_type":"markdown","source":"

3.\tFTU (Functional Tissue Unit)

\nFTU or functional tissue unit, is a three-dimensional, maximally connected, block of cells centered around a capillary, such that each cell in this block is within diffusion distance from any other cell in the same block.","metadata":{}},{"cell_type":"markdown","source":"

Let’s understand the DATA now,

\nThe Dataset is comprised of very large TIFF files (TIFF is described later in the notebook).\n\n•\tThe training set has 8 files.\n\n•\tThe public test set has 5 files.\n\n•\tThe private test set is larger than the public test set. I suppose there will be 7 files.\n\nThe train set includes annotations in both RLE-encoded (RLE is explained later in the notebook) and unencoded (JSON) forms. The annotations denote segmentations of glomeruli.\nBoth training and public test sets include anatomical structure segmentations. I suppose this can be used for pretraining.\nJSON files are structured as follows\n\n•\tA type (Feature) and object type id (PathAnnotationObject). Note that these fields are the same between all files and do not offer signal.\n\n•\tA geometry containing a Polygon with coordinates for the feature's enclosing volume\n\n•\tAdditional properties, including the name and color of the feature in the image.\n\n•\tThe IsLocked field is the same across file types (locked in glomerulus, unlocked for anatomical structure) and is not signal bearing.","metadata":{}},{"cell_type":"markdown","source":"

TIFF File Format

\nTag Image File Format, abbreviated TIFF or TIF, is a computer file format for storing raster graphics images, popular among graphic artists, the publishing industry, and photographers. TIFF is a flexible, adaptable file format for handling images and data within a single file, by including the header tags (size, definition, image-data arrangement, applied image compression) defining the image's geometry.","metadata":{}},{"cell_type":"markdown","source":"

RLE

\nThe masks provided in the train.csv is in Running Length Encoding format. This encoding comes in pairs of pixel values as follows:\n\n1.\tThe starting pixel.\n2.\tNumber of pixels from the starting pixel.\n\nSo, to specify 10 pixels starting from pixel number 200 would be written as: 200 10","metadata":{}},{"cell_type":"markdown","source":"

Evaluation Metric

","metadata":{}},{"cell_type":"markdown","source":"

Dice Coefficient

\nDice coefficient is a statistical tool which measures the similarity between two sets of data. This index has become arguably the most broadly used tool in the validation of image segmentation algorithms created with AI, but it is a much more general concept which can be applied sets of data for a variety of applications including NLP.\nThe Dice coefficient can be used to compare the pixel-wise agreement between a predicted segmentation and its corresponding ground truth.","metadata":{}},{"cell_type":"markdown","source":"

Jaccard Score

\nThe Jaccard similarity index (sometimes called the Jaccard similarity coefficient) compares members for two sets to see which members are shared and which are distinct. It’s a measure of similarity for the two sets of data, with a range from 0% to 100%. The higher the percentage, the more similar the two populations. Although it’s easy to interpret, it is extremely sensitive to small samples sizes and may give erroneous results, especially with very small samples or data sets with missing observations.\nThe formula to find the Index is:\nJaccard Index = (the number in both sets) / (the number in either set) * 100","metadata":{}},{"cell_type":"markdown","source":"

Relation between Dice Coefficient and Jaccard Score

\n

J=D/2-D

\n","metadata":{}},{"cell_type":"markdown","source":"

Let’s do the EDA now,

","metadata":{}},{"cell_type":"code","source":"# Libraries\nimport cv2\nimport datetime\nimport gc\nimport glob\nimport math\nimport matplotlib.pyplot as plt\nimport numpy as np\nimport os\nimport pandas as pd\nimport skimage.morphology\nimport sys\nimport tensorflow as tf\nimport tifffile","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Parameters\nbase_path = '../input/hubmap-kidney-segmentation'\n\nplot_full_image = True\n\n# Number of glomeruli to display for each image\nnum_glom_display = 5\n\n# Number of glomberuli to save as tiff files.\nnum_glom_save = 5\n\nglob_scale = 0.25","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Utility Functions\ndef rle_to_image(rle_mask, image_shape):\n \"\"\"\n Converts an rle string to an image represented as a numpy array.\n Reference: https://www.kaggle.com/paulorzp/rle-functions-run-lenght-encode-decode\n\n :param rle_mask: string with rle mask.\n :param image_shape: (width, height) of array to return\n :return: Image as a numpy array. 1 = mask, 0 = background.\n \"\"\"\n\n # Processing\n s = rle_mask.split()\n starts, lengths = [np.asarray(x, dtype=int) for x in (s[0:][::2], s[1:][::2])]\n starts -= 1\n ends = starts + lengths\n image = np.zeros(image_shape[0] * image_shape[1], dtype=np.uint8)\n for lo, hi in zip(starts, ends):\n image[lo:hi] = 1\n\n return image.reshape(image_shape).T","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

File Structure

\nThe files in the root of the dataset are shown below. The dataset consists of 2 directories that contain training and test images and 3 csv-files with additional information about the images.","metadata":{}},{"cell_type":"code","source":"#Directory Contents\nprint('\\n'.join(os.listdir(base_path)))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Training Images\ntrain_files = sorted(glob.glob(os.path.join(base_path, 'train/*.tiff')))\nprint(f'Number of training images: {len(train_files)}')\nprint('\\n'.join(train_files))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Test Images\ntest_files = sorted(glob.glob(os.path.join(base_path, 'test/*.tiff')))\nprint(f'Number of test images: {len(test_files)}')\nprint('\\n'.join(test_files))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Train CSV\n#The masks indicating a glomeruli FTUs are stored in rle format in the train.csv for each training image id.\ndf_train = pd.read_csv(os.path.join(base_path, 'train.csv'))\ndisplay(df_train)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Sample_Submission.csv\n#The sample_submission.csv files shows the format of the submissions files consisting of the test image id and an rle encoded masks.\ndf_submission = pd.read_csv(os.path.join(base_path,'sample_submission.csv'))\ndisplay(df_submission)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Paitient Data\n#HuBMAP-20-dataset_information.csv contains additional information about each image such as image size and anonymized patient data.\ndf_info = pd.read_csv(os.path.join(base_path,'HuBMAP-20-dataset_information.csv'))\ndisplay(df_info)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

Training Image Analysis

\n

Width and Height Distribution

\n\nThe training images do not have consistent dimensions. This has to be corrected when loading the images. They have on of the following shapes:\n\n[height, width, channel]\n \n[channel, height, width]\n\n[1, 1, channel, height, width]","metadata":{}},{"cell_type":"code","source":"for f in train_files + test_files:\n image = tifffile.imread(f)\n print(f'Image {f} shape: {image.shape}', flush=True)\n del image\n gc.collect()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#The size of the images varies greatly as well.\nplt.scatter(df_info['width_pixels'], df_info['height_pixels'])\nplt.title('Image Height and Width')\nplt.xlabel('Width')\nplt.ylabel('Height')\nplt.xlim(0, df_info['width_pixels'].max() * 1.1)\nplt.ylim(0, df_info['height_pixels'].max() * 1.1)\nplt.grid()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

Image Utilitity Functions

","metadata":{}},{"cell_type":"code","source":"def overlay_image_mask(image, mask, mask_color=(0,255,0), alpha=1.0):\n im_f= image.astype(np.float32)\n# if mask.ndim == 2:\n# mask = np.expand_dims(mask,-1) \n mask_col = np.expand_dims(np.array(mask_color)/255.0, axis=(0,1))\n return (im_f + alpha * mask * (np.mean(0.8 * im_f + 0.2 * 255, axis=2, keepdims=True) * mask_col - im_f)).astype(np.uint8)\n\n\ndef overlay_image_mask_original(image, mask, mask_color=(0,255,0), alpha=1.0):\n return np.concatenate((image, overlay_image_mask(image, mask)), axis=1)\n\ndef get_image_id(image_file):\n return os.path.splitext(os.path.split(image_file)[1])[0]\n\n\ndef read_image(image_file, scale=1.0):\n image = tifffile.imread(image_file).squeeze()\n if image.shape[0] == 3:\n image = np.transpose(image, (1,2,0))\n \n orig_shape = image.shape\n if scale != 1.0:\n image = cv2.resize(image, (0,0), fx=scale, fy=scale)\n return image, orig_shape\n\n\ndef read_mask(image_file, image_shape, scale=1.0):\n image_id = get_image_id(image_file)\n train_info = df_train.loc[df_train['id'] == image_id]\n rle = train_info['encoding'].values[0] if len(train_info) > 0 else None\n if rle is not None:\n mask = rle_to_image(rle, (image_shape[1], image_shape[0]))\n if scale != 1.0:\n mask = cv2.resize(mask, (0,0), fx=scale, fy=scale)\n return np.expand_dims(mask,-1)\n else:\n return None \n\n \ndef read_image_mask(image_file, scale=1.0):\n image, image_shape = read_image(image_file, scale)\n mask = read_mask(image_file, image_shape, scale)\n return image, mask\n\n\ndef get_tile(image, mask, x, y, tile_size, scale=1.0):\n x = round(x * scale)\n y = round(y * scale)\n size = int(round(tile_size / 2 * scale))\n image_s = image[y-size:y+size, x-size:x+size, :] \n mask_s = mask[y-size:y+size, x-size:x+size, :]\n return image_s, mask_s\n\n\ndef get_particles(mask, scale=1.0):\n num, labels, stats, centroids = cv2.connectedComponentsWithStats(mask)\n df_particles = pd.DataFrame(dict(zip(['x','y','left','top','width','height','area'],\n [(centroids[1:,0]) / scale,\n (centroids[1:,1]) / scale,\n (stats[1:,cv2.CC_STAT_LEFT]) / scale,\n (stats[1:,cv2.CC_STAT_TOP]) / scale,\n (stats[1:,cv2.CC_STAT_WIDTH]) / scale,\n (stats[1:,cv2.CC_STAT_HEIGHT]) / scale,\n (stats[1:,cv2.CC_STAT_AREA]) / (scale * scale)])))\n df_particles.sort_values(['x','y'], inplace=True, ignore_index=True)\n df_particles['no'] = range(len(df_particles))\n return df_particles\n\n\ndef analyze_image(image_file):\n image_id = get_image_id(image_file)\n image, image_shape = read_image(image_file, glob_scale)\n mask = read_mask(image_file, image_shape, glob_scale)\n\n mask_full = read_mask(image_file, image_shape, scale=1.0)\n df_glom = get_particles(mask_full, scale=1.0)\n df_glom['id'] = image_id\n del mask_full\n gc.collect()\n \n info = df_info[df_info['image_file'] == f'{image_id}.tiff']\n print(f'Image ID: {image_id:}')\n print(f'Image Size: {info[\"width_pixels\"].values[0]} x {info[\"height_pixels\"].values[0]}')\n print(f'Patient No: {info[\"patient_number\"].values[0]}')\n print(f'Sex: {info[\"sex\"].values[0]}')\n print(f'Age: {info[\"age\"].values[0]}')\n print(f'Race: {info[\"race\"].values[0]}')\n print(f'Height: {info[\"height_centimeters\"].values[0]} cm')\n print(f'Weight: {info[\"weight_kilograms\"].values[0]} kg')\n print(f'BMI: {info[\"bmi_kg/m^2\"].values[0]} kg/m^2')\n print(f'Laterality: {info[\"laterality\"].values[0]}')\n print(f'Percent Cortex: {info[\"percent_cortex\"].values[0]} %')\n print(f'Percent Medulla: {info[\"percent_medulla\"].values[0]} %')\n \n # Plot full image\n if plot_full_image:\n scale = 0.1\n image_small = cv2.resize(image, (0,0), fx=scale, fy=scale)\n mask_small = cv2.resize(mask, (0,0), fx=scale, fy=scale)\n mask_small = np.expand_dims(mask_small,-1) \n \n plt.figure(figsize=(16, 16))\n plt.imshow(overlay_image_mask(image_small, mask_small))\n plt.axis('off')\n\n # Plot glomeruli images\n fig_cols = 5\n fig_rows = int(math.ceil(num_glom_display/fig_cols))\n plt.figure(figsize=(4 * fig_cols, 4 * fig_rows))\n if num_glom_save > 0 and not os.path.exists(image_id):\n os.mkdir(image_id)\n for i in range(min(max(num_glom_display, num_glom_save), len(df_glom))):\n image_s, mask_s = get_tile(image,mask, df_glom['x'][i], df_glom['y'][i], 1000, scale=glob_scale)\n ovl = overlay_image_mask(image_s, mask_s)\n if i < num_glom_display:\n plt.subplot(fig_rows, fig_cols, i+1)\n plt.imshow(ovl)\n plt.axis('off')\n if i < num_glom_save:\n cv2.imwrite(f'{image_id}_{i:03}.png', cv2.cvtColor(ovl, cv2.COLOR_RGB2BGR)) \n \n del image, mask\n gc.collect()\n return df_glom\n\n\ndef plot_glom(df, image_id, glom_no):\n image, mask = read_image_mask(os.path.join(base_path, f'train/{image_id}.tiff'), scale=glob_scale)\n glom = df.loc[(df['id'] == image_id) & (df['no'] == glom_no)]\n im, ma = get_tile(image, mask, glom['x'].iloc[0], glom['y'].iloc[0], 1000, scale=glob_scale)\n del image, mask\n gc.collect()\n plt.figure(figsize=(16,8))\n plt.imshow(overlay_image_mask_original(im, ma))\n plt.title(f'Image: {image_id}, Glomeruli No: {glom_no}, Area: {glom[\"area\"].iloc[0]}')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

Training Images With Glomerulis

","metadata":{}},{"cell_type":"code","source":"df_glom = pd.DataFrame()\ndf_glom = df_glom.append(analyze_image(train_files[0]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[1]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[2]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[3]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[4]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[5]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[6]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[7]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[8]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[9]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[10]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[11]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[12]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[13]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom = df_glom.append(analyze_image(train_files[14]), ignore_index=True)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

Glomerulis

","metadata":{}},{"cell_type":"markdown","source":"

Basic Statistics

","metadata":{}},{"cell_type":"code","source":"df_glom.to_csv('glomeruli.csv')\ndisplay(df_glom)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_glom.describe()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

Glomerulis Per Image

","metadata":{}},{"cell_type":"code","source":"g = df_glom.groupby('id')\nplt.bar(g.size().index, g.size().values)\nplt.title('Number of Glomerulis in Image')\nplt.xticks(rotation=90)\nplt.grid()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

Glomeruli Width, Height and Area Distribution

","metadata":{}},{"cell_type":"code","source":"plt.figure(figsize=(20,5))\nplt.subplot(1,3,1)\nplt.hist(df_glom['width'], bins=40, density=True)\nplt.title('Width Distribution')\nplt.grid()\nplt.subplot(1,3,2)\nplt.hist(df_glom['height'], bins=40, density=True)\nplt.title('Height Distribution')\nplt.grid()\nplt.subplot(1,3,3)\nplt.hist(df_glom['area'], bins=40, density=True)\nplt.title('Area Distribution')\nplt.grid()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

Glomerulis by Size

","metadata":{}},{"cell_type":"code","source":"df_glom.sort_values('area', inplace=True)\ndf_glom","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

Five Smallest Glomerulis

","metadata":{}},{"cell_type":"code","source":"for i in range(5):\n plot_glom(df_glom, df_glom['id'].iloc[i], df_glom['no'].iloc[i])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

Five Largest Glomerulis

","metadata":{}},{"cell_type":"code","source":"for i in range(len(df_glom)-5, len(df_glom)):\n plot_glom(df_glom, df_glom['id'].iloc[i], df_glom['no'].iloc[i])","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

SUB EffUNet5 + TPU EfficientUNet 512x512

","metadata":{}},{"cell_type":"markdown","source":"In the code below I will try to use TPU EfficientUNet 512x512 with freeze-pretrained SUB EffUNet5. This will help to improve the existing algorithms used to detect functional tissue units (FTUs) across different tissue preparation pipelines.","metadata":{}},{"cell_type":"code","source":"mod_path = '/kaggle/input/hubmap-tf-with-tpu-efficientunet-512x512-train/'\nimport yaml\nimport pprint\nwith open(mod_path+'params.yaml') as file:\n P = yaml.load(file, Loader=yaml.FullLoader)\n pprint.pprint(P)\n \nTHRESHOLD = 0.4\nWINDOW = 1024\nMIN_OVERLAP = 300\nNEW_SIZE = P['DIM']\n\nSUBMISSION_MODE = 'PUBLIC_TFREC' \n# 'PUBLIC_TFREC' = use created tfrecords for public test set with MIN_OVERLAP = 300 tiling 1024-512, ignore other (private test) data\n# 'FULL' do not use tfrecords, just full submission \n\nCHECKSUM = True # compute mask sum for each image","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# METRICS\n\nimport json\n\nwith open(mod_path + 'metrics.json') as json_file:\n M = json.load(json_file)\nprint('Model run datetime: '+M['datetime'])\nprint('OOF val_dice_coe: ' + str(M['oof_dice_coe']))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"! pip install ../input/kerasapplications/keras-team-keras-applications-3b180cb -f ./ --no-index -q\n! pip install ../input/efficientnet/efficientnet-1.1.0/ -f ./ --no-index -q\nimport numpy as np\nimport pandas as pd\nimport os\nimport glob\nimport gc\n\nimport rasterio\nfrom rasterio.windows import Window\n\nimport pathlib\nfrom tqdm.notebook import tqdm\nimport cv2\n\nimport tensorflow as tf\nimport efficientnet as efn\nimport efficientnet.tfkeras\n\nimport os, glob, gc\nimport json\n\nosj = os.path.join","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def rle_encode_less_memory(img):\n pixels = img.T.flatten()\n pixels[0] = 0\n pixels[-1] = 0\n runs = np.where(pixels[1:] != pixels[:-1])[0] + 2\n runs[1::2] -= runs[::2]\n return ' '.join(str(x) for x in runs)\n\ndef make_grid(shape, window=256, min_overlap=32):\n \"\"\"\n Return Array of size (N,4), where N - number of tiles,\n 2nd axis represente slices: x1,x2,y1,y2 \n \"\"\"\n x, y = shape\n nx = x // (window - min_overlap) + 1\n x1 = np.linspace(0, x, num=nx, endpoint=False, dtype=np.int64)\n x1[-1] = x - window\n x2 = (x1 + window).clip(0, x)\n ny = y // (window - min_overlap) + 1\n y1 = np.linspace(0, y, num=ny, endpoint=False, dtype=np.int64)\n y1[-1] = y - window\n y2 = (y1 + window).clip(0, y)\n slices = np.zeros((nx,ny, 4), dtype=np.int64)\n \n for i in range(nx):\n for j in range(ny):\n slices[i,j] = x1[i], x2[i], y1[j], y2[j] \n return slices.reshape(nx*ny,4)\n\ndef global_shift_mask(maskpred1, y_shift, x_shift):\n \"\"\"\n applies a global shift to a mask by \n padding one side and cropping from the other\n \"\"\"\n if y_shift < 0 and x_shift >=0:\n maskpred2 = np.pad(maskpred1, \n [(0,abs(y_shift)), (abs(x_shift), 0)], \n mode='constant', constant_values=0)\n maskpred3 = maskpred2[abs(y_shift):, :maskpred1.shape[1]]\n elif y_shift >=0 and x_shift <0:\n maskpred2 = np.pad(maskpred1, \n [(abs(y_shift),0), (0, abs(x_shift))], \n mode='constant', constant_values=0)\n maskpred3 = maskpred2[:maskpred1.shape[0], abs(x_shift):]\n elif y_shift >=0 and x_shift >=0:\n maskpred2 = np.pad(maskpred1,\n [(abs(y_shift),0), (abs(x_shift), 0)], \n mode='constant', constant_values=0)\n maskpred3 = maskpred2[:maskpred1.shape[0], :maskpred1.shape[1]]\n elif y_shift < 0 and x_shift < 0:\n maskpred2 = np.pad(maskpred1, \n [(0, abs(y_shift)), (0, abs(x_shift))], \n mode='constant', constant_values=0)\n maskpred3 = maskpred2[abs(y_shift):, abs(x_shift):]\n return maskpred3","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"##MODEL\nidentity = rasterio.Affine(1, 0, 0, 0, 1, 0)\nfold_models_1 = []\nfor fold_model_path in glob.glob(mod_path+'*.h5'):\n fold_models_1.append(tf.keras.models.load_model(fold_model_path,compile = False))\nprint(len(fold_models_1))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"AUTO = tf.data.experimental.AUTOTUNE\nimage_feature = {\n 'image': tf.io.FixedLenFeature([], tf.string),\n 'x1': tf.io.FixedLenFeature([], tf.int64),\n 'y1': tf.io.FixedLenFeature([], tf.int64)\n}\ndef _parse_image(example_proto):\n example = tf.io.parse_single_example(example_proto, image_feature)\n image = tf.reshape( tf.io.decode_raw(example['image'],out_type=np.dtype('uint8')), (P['DIM'],P['DIM'], 3))\n return image, example['x1'], example['y1']\n\ndef load_dataset(filenames, ordered=True):\n ignore_order = tf.data.Options()\n if not ordered:\n ignore_order.experimental_deterministic = False\n dataset = tf.data.TFRecordDataset(filenames)\n dataset = dataset.with_options(ignore_order)\n dataset = dataset.map(_parse_image)\n return dataset\n\ndef get_dataset(FILENAME):\n dataset = load_dataset(FILENAME)\n dataset = dataset.batch(64)\n dataset = dataset.prefetch(AUTO)\n return dataset","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"debug = True # True False\nn_debug_images = 1 if debug else 1000000000\nn_debug_slices = 20 if debug else 1000000000\n\n# whether to run prediction when committing. WILL RUN predictions during submission in any case\ndo_predict = False if not debug else True\n\nmodels_dir = '../input/hubmap-models-cv-08848-pl-0847'\nmodel_filepaths = [ os.path.join(models_dir, f\"model-fold-{i}.h5\") for i in range(4)]\n\nassert len(model_filepaths)==len(np.unique(model_filepaths))\n#folds_to_predict = [i for (i, fn) in enumerate(model_filepaths) if os.path.isfile(fn)]\nmodel_dirnames = [os.path.dirname(filepath) for filepath in model_filepaths]\n\n#check_order = [fn.split('.')[-2].split('-')[-1] == i for (i,fn) in enumerate(model_filepaths) if fn.strip()!='']\n#assert np.sum(check_order)==0, 'models should be in folds order or empty string'\n\nimport yaml\nimport pprint\nwith open(osj(model_dirnames[0],'params.yaml')) as file:\n P = yaml.load(file, Loader=yaml.FullLoader)\n pprint.pprint(P)\n\nTHRESHOLD = 0.3\nWINDOW = 1024\nMIN_OVERLAP = 32\nNEW_SIZE = P['DIM']\n\nassert sum([not os.path.isfile(path_) for path_ in model_filepaths]) == 0\nprint(\"\\n Number of models:: {}\".format(len(model_filepaths)))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"ave_score = 0\nfor i, m_path in enumerate(model_filepaths):\n fold_ = int(m_path.split('.')[-2].split('-')[-1])\n with open(osj(model_dirnames[i],'metrics.json')) as json_file:\n M = json.load(json_file)\n print(f\"\\n ----------- \\nModel {model_dirnames[i].split('/')[-1]}\" +\n '\\nval_dice_coe: '+ str(round(M['val_dice_coe'][fold_], 5)) +\n '\\tval_loss: ' + str(round(M['val_loss'][fold_], 5)) +\n '\\tval_accuracy: '+ str(round(M['val_accuracy'][fold_], 5))\n )\n\n\n\nfor model_group in np.unique(model_dirnames):\n with open(osj(model_group,'metrics.json')) as json_file:\n M = json.load(json_file)\n ave_dice = np.mean(M['val_dice_coe']) \n ave_loss = np.mean(M['val_loss']) # /len(folds_to_predict)\n ave_accuracy = np.mean(M['val_accuracy'])\n print(f\"\\n ============ MODEL GROUP {model_group} ==============\")\n print(\" ------------ \\nAVERAGE DICE SCORE = {}\".format(round(ave_dice, 5)))\n print(\" ------------ \\nAVERAGE VALIDATION LOSS = {}\".format(round(ave_loss, 5)))\n print(\" ------------ \\nAVERAGE VALIDATION ACCURACY = {}\".format(round(ave_accuracy, 5)))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nif do_predict:\n identity = rasterio.Affine(1, 0, 0, 0, 1, 0)\n fold_models_2 = []\n \n for fold_model_path in model_filepaths:\n fold_models_2.append(tf.keras.models.load_model(fold_model_path,compile = False))\n print(len(fold_models_2))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"

References

\n\nhttps://www.kaggle.com/kwk100/hubmap-exploratory-data-analysis-eda\n\nhttps://www.kaggle.com/roydatascience/hubmap-sub-effunet5-tpu-efficientunet-512x512","metadata":{}},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}