postprocessing.ipynb

{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Postprocessing\n",
    "\n",
    "After training the model and using it to predict submission cases (210-300), in order to send a submission, we must be able to \"inverse\" the processing in order to generate segmentation nifti file that matches voxel dimension and size of the original image, to do that we must store voxdim and shape of the unprocessed test cases.\n",
    "\n",
    "Let's create a json to keep track of essential data in order to complete the submission"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pickle\n",
    "import numpy as np\n",
    "import skimage.transform\n",
    "from scipy.ndimage import zoom\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_path = \"data/case_\"\n",
    "imaging_fname = \"/imaging.nii.gz\"\n",
    "segmentation_fname = \"/segmentation.nii.gz\"\n",
    "\n",
    "class Patient:\n",
    "  def __init__(self, number):\n",
    "      self.number = number\n",
    "      \n",
    "      self.imaging_path = self.get_im_data_path()\n",
    "      self.imaging = nib.load(self.imaging_path)\n",
    "      \n",
    "      if number < 210: # segmentation file only exists for non-submission cases\n",
    "        self.segmentation_path = self.get_seg_data_path()\n",
    "        self.segmentation = nib.load(self.segmentation_path)\n",
    "\n",
    "  def get_im_data_path(self):\n",
    "      '''\n",
    "      returns path to imaging file\n",
    "      '''\n",
    "      out = \"00000\"\n",
    "\n",
    "      out = out[len(str(self.number)):]\n",
    "      out += str(self.number)\n",
    "      return data_path + out + imaging_fname\n",
    "  \n",
    "  def get_seg_data_path(self):\n",
    "      '''\n",
    "      returns path to segmentation file\n",
    "      '''\n",
    "      out = \"00000\"\n",
    "\n",
    "      out = out[len(str(self.number)):]\n",
    "      out += str(self.number)\n",
    "      return data_path + out + segmentation_fname\n",
    "\n",
    "  def rr(self, target_file, target_resolution=[1,1,1], target_shape=[192,192,192]):\n",
    "    '''\n",
    "    Target file can be either self.imaging or self.segmentation\n",
    "    Function that returns reshaped and resized target block\n",
    "    ''' \n",
    "    target_fig = target_file.get_fdata()\n",
    "\n",
    "    # x,y,z\n",
    "    vox_dim = (target_file.header['pixdim'][2], target_file.header['pixdim'][3], target_file.header['pixdim'][1])\n",
    "\n",
    "    scale_vector = (vox_dim[0]/target_resolution[0],\n",
    "                    vox_dim[1]/target_resolution[1],\n",
    "                    vox_dim[2]/target_resolution[2])\n",
    "\n",
    "    \n",
    "\n",
    "    iso_target_fig = skimage.transform.rescale(target_fig, scale_vector, order=3, preserve_range=True,  mode='constant')\n",
    "    \n",
    "    # print(iso_target_fig.dtype) # float64\n",
    "    \n",
    "    factors = ( target_shape[0]/iso_target_fig.shape[0],\n",
    "                target_shape[1]/iso_target_fig.shape[1],\n",
    "                target_shape[2]/iso_target_fig.shape[2],)\n",
    "\n",
    "    rr_target_cube = zoom(iso_target_fig, factors, order=3, mode='nearest')\n",
    "\n",
    "    # print(rr_target_cube.dtype) # float64\n",
    "\n",
    "    # rr = rescaled and reshaped\n",
    "    # rr_target = nib.Nifti1Image(rr_target_fig, target_file.affine)\n",
    "\n",
    "    # # setto le dimensioni manualmente, sono abbastanza sicuro che il rescaling sia avvenuto correttamente\n",
    "\n",
    "    # for i in range(3):\n",
    "    #   rr_target.header['pixdim'][i+1] = target_resolution[i]\n",
    "\n",
    "    # for i in range(3):\n",
    "    #   rr_target.affine[i][2-i] = -target_resolution[i]\n",
    "\n",
    "    return rr_target_cube\n",
    "\n",
    "  def process_imaging(self, target_resolution=[1,1,1], target_shape=[192,256,256]):\n",
    "    '''\n",
    "    function that returns reshaped, resized, normalized imaging block\n",
    "    '''\n",
    "\n",
    "    rr_imaging = self.rr(self.imaging, target_resolution, target_shape)\n",
    "    output = cv2.normalize(rr_imaging, None, alpha = 0, beta = 255, norm_type = cv2.NORM_MINMAX, dtype = cv2.CV_8UC1)\n",
    "    return output\n",
    "\n",
    "  def process_segmentation(self, target_resolution=[1,1,1], target_shape=[192,256,256]):\n",
    "    '''\n",
    "    function that returns reshaped, resized segmentation block\n",
    "    '''\n",
    "\n",
    "    output = self.rr(self.segmentation, target_resolution, target_shape)\n",
    "\n",
    "    # might break things (it does, kinda)\n",
    "    #output = cv2.normalize(rr_segmentation, None, alpha = 0, beta = 3, norm_type = cv2.NORM_MINMAX, dtype = cv2.CV_8UC1)\n",
    "    \n",
    "    output = np.rint(output)\n",
    "    output = output.astype('uint8')\n",
    "\n",
    "    return output\n",
    "\n",
    "def get_number(n):\n",
    "  '''\n",
    "  returns num padded with zeros\n",
    "  '''\n",
    "  out = \"00000\"\n",
    "\n",
    "  out = out[len(str(n)):]\n",
    "  out += str(n)\n",
    "  return out"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "record = []\n",
    "\n",
    "for n_pat in range(0,300):\n",
    "    pat = Patient(n_pat)\n",
    "\n",
    "    #print(pat.imaging.header)\n",
    "    # x, y, z\n",
    "    shape = [pat.imaging.header['dim'][2], pat.imaging.header['dim'][3], pat.imaging.header['dim'][1]]\n",
    "    vox_dim = [pat.imaging.header['pixdim'][2], pat.imaging.header['pixdim'][3], pat.imaging.header['pixdim'][1]]\n",
    "    \n",
    "    data = {\n",
    "        'case_id': n_pat,\n",
    "        'vox_dim': vox_dim,\n",
    "        'shape': shape,\n",
    "    }\n",
    "\n",
    "    record.append(data)\n",
    "\n",
    "with open('original_dims.pkl', 'wb') as f:\n",
    "    pickle.dump(record, f)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Provo con 1 caso (ID168)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open('original_dims.pkl', 'rb') as f:\n",
    "    original_dims = pickle.load(f)\n",
    "\n",
    "patVolList = []\n",
    "\n",
    "path = 'D:\\\\kits19\\\\processed_data\\\\predict_aug_09_12\\\\'\n",
    "\n",
    "for filename in os.listdir(path):\n",
    "    print(filename)\n",
    "    img = Image.open(str(path+filename))\n",
    "    slice = np.asarray(img)\n",
    "    patVolList.append(slice)\n",
    "\n",
    "pat168vol = np.array(patVolList)\n"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Function that revert effects of preprocessing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def revert(volume, n_pat):\n",
    "    \n",
    "    vox_dim = original_dims[n_pat]['vox_dim']\n",
    "    shape = original_dims[n_pat]['shape']\n",
    "\n",
    "    # print(vox_dim)\n",
    "    # z x y\n",
    "    scale_vector = (1/vox_dim[2],\n",
    "                    1/vox_dim[0],\n",
    "                    1/vox_dim[1])\n",
    "    # print(scale_vector)\n",
    "\n",
    "    iso_target_fig = skimage.transform.rescale(volume, scale_vector, order=3, preserve_range=True,  mode='constant')\n",
    "\n",
    "    # z x y\n",
    "    factors = ( shape[2]/iso_target_fig.shape[0],\n",
    "                shape[0]/iso_target_fig.shape[1],\n",
    "                shape[1]/iso_target_fig.shape[2],)\n",
    "\n",
    "    output = zoom(iso_target_fig, factors, order=3, mode='nearest')\n",
    "    \n",
    "\n",
    "    return output"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "original_form = revert(pat168vol, 168)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Evaluate single case"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from starter_code import evaluation\n",
    "\n",
    "print(evaluation.evaluate(168, original_form))"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Evaluate all test cases"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open('original_dims.pkl', 'rb') as f:\n",
    "    original_dims = pickle.load(f)\n",
    "\n",
    "\n",
    "path = 'D:\\\\kits19\\\\predicts\\\\batch_predict_tumor_11_12\\\\'\n",
    "\n",
    "listOfVolumes = []\n",
    "\n",
    "for num_pat in range(168, 210):\n",
    "    volume = []\n",
    "    case_file = 'case_' + get_number(num_pat) + '_predict_img\\\\'\n",
    "    \n",
    "    for filename in os.listdir(str(path+case_file)):\n",
    "        img = Image.open(str(path+case_file+filename))\n",
    "        npimg = np.asarray(img)\n",
    "        volume.append(npimg)\n",
    "    listOfVolumes.append(volume)\n",
    "\n",
    "volumesToRevert = []  # list of ndarray\n",
    "for volume in listOfVolumes:\n",
    "    volumesToRevert.append(np.array(volume))"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Revert all prediction ndarrays into original shape and vox_dim"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def revertGenerator():\n",
    "    for i in range(168,210):\n",
    "        readyVol = revert(volumesToRevert[i-168], i)\n",
    "        yield (i, readyVol)\n",
    "\n",
    "#generator object\n",
    "evalReadyVolGen = revertGenerator()"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally evaluate all test_imgs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from starter_code import evaluation\n",
    "\n",
    "kidneyScoreList = []\n",
    "\n",
    "for i in range(168,210):\n",
    "    predict_data = next(evalReadyVolGen)\n",
    "    print(predict_data[0])\n",
    "    kidneyScoreList.append(evaluation.evaluate(predict_data[0], predict_data[1])[1])\n",
    "    print(kidneyScoreList[-1])\n",
    "\n",
    "partial_sum = 0\n",
    "for val in kidneyScoreList:\n",
    "    partial_sum += val\n",
    "\n",
    "print('Average score:', partial_sum/len(kidneyScoreList))\n",
    "\n",
    "# 57% kidney\n",
    "# 21% cancer"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "Create folder for evaluation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pickle\n",
    "\n",
    "with open('original_dims.pkl', 'rb') as f:\n",
    "    original_dims = pickle.load(f)\n",
    "\n",
    "\n",
    "path = 'D:\\\\kits19\\\\predicts\\\\sub_img_kidney_1\\\\'\n",
    "\n",
    "\n",
    "for num_pat in range(210, 300):\n",
    "    volume = []\n",
    "    case_file = 'case_' + get_number(num_pat) + '_predict_img\\\\'\n",
    "    \n",
    "    for filename in os.listdir(str(path+case_file)):\n",
    "        img = Image.open(str(path+case_file+filename))\n",
    "        npimg = np.asarray(img)\n",
    "        volume.append(npimg)\n",
    "    \n",
    "    volumeToRevert = np.array(volume)\n",
    "    evalReadyVolume = revert(volumeToRevert, num_pat)\n",
    "\n",
    "    evalReadyVolume[evalReadyVolume<0.7]=0\n",
    "    evalReadyVolume[evalReadyVolume>=0.7]=1\n",
    "\n",
    "    pat = Patient(num_pat)\n",
    "    final_predict = nib.Nifti1Image(evalReadyVolume, affine=pat.imaging.affine)\n",
    "\n",
    "    nib.save(final_predict, os.path.join(\"D:\\\\kits19\\\\submission_ready\\\\predictions\", str(\"prediction_ \"+get_number(num_pat)+\".nii.gz\")))"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "# Plot results\n",
    "- first for kidney model\n",
    "- then for tumor model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "predict_168_65 = np.asanyarray(Image.open('D:\\\\kits19\\\\predicts\\\\batch_predict_09_12_56\\\\case_00168_predict_img\\\\case_00168_kid_00065_predict.png'))\n",
    "actual_168_65 = np.asanyarray(Image.open('D:\\\\kits19\\\\processed_data\\\\test_kidney_labels\\\\case_00168_seg_00065.png'))\n",
    "imaging_168_65 = np.asanyarray(Image.open('D:\\\\kits19\\\\processed_data\\\\test_images\\\\case_00168_test_img\\\\case_00168_im_00065.png'))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "fig = plt.figure(figsize=(8,8))\n",
    "a = fig.add_subplot(1,3,1)\n",
    "imgplot = plt.imshow(imaging_168_65[:,:], cmap = 'gray')\n",
    "plt.title(\"Image\")\n",
    "\n",
    "a = fig.add_subplot(1,3,2)\n",
    "imgplot = plt.imshow(predict_168_65[:,:], cmap='gray')\n",
    "plt.title(\"Predicted Segmentation\")\n",
    "\n",
    "a = fig.add_subplot(1,3,3)\n",
    "imgplot = plt.imshow(actual_168_65[:,:], cmap='gray')\n",
    "plt.title(\"Actual Segmentation\")"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "predict_197_61 = np.asanyarray(Image.open('D:\\\\kits19\\\\predicts\\\\batch_predict_tumor_11_12\\\\case_00197_predict_img\\\\case_00197_kid_00061_predict.png'))\n",
    "actual_197_61 = np.asanyarray(Image.open('D:\\\\kits19\\\\processed_data\\\\test_tumor_labels\\\\case_00197_tum_00061.png'))\n",
    "imaging_197_61 = np.asanyarray(Image.open('D:\\\\kits19\\\\processed_data\\\\test_images\\\\case_00197_test_img\\\\case_00197_im_00061.png'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(8,8))\n",
    "a = fig.add_subplot(1,3,1)\n",
    "imgplot = plt.imshow(imaging_197_61[:,:], cmap = 'gray')\n",
    "plt.title(\"Image\")\n",
    "\n",
    "a = fig.add_subplot(1,3,2)\n",
    "imgplot = plt.imshow(predict_197_61[:,:], cmap='gray')\n",
    "plt.title(\"Predicted Segmentation\")\n",
    "\n",
    "a = fig.add_subplot(1,3,3)\n",
    "imgplot = plt.imshow(actual_197_61[:,:], cmap='gray')\n",
    "plt.title(\"Actual Segmentation\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "predict_196_70 = np.asanyarray(Image.open('D:\\\\kits19\\\\predicts\\\\batch_predict_tumor_11_12\\\\case_00196_predict_img\\\\case_00196_kid_00070_predict.png'))\n",
    "actual_196_70 = np.asanyarray(Image.open('D:\\\\kits19\\\\processed_data\\\\test_tumor_labels\\\\case_00196_tum_00070.png'))\n",
    "imaging_196_70 = np.asanyarray(Image.open('D:\\\\kits19\\\\processed_data\\\\test_images\\\\case_00196_test_img\\\\case_00196_im_00070.png'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize=(8,8))\n",
    "a = fig.add_subplot(1,3,1)\n",
    "imgplot = plt.imshow(imaging_196_70[:,:], cmap = 'gray')\n",
    "plt.title(\"Image\")\n",
    "\n",
    "a = fig.add_subplot(1,3,2)\n",
    "imgplot = plt.imshow(predict_196_70[:,:], cmap='gray')\n",
    "plt.title(\"Predicted Segmentation\")\n",
    "\n",
    "a = fig.add_subplot(1,3,3)\n",
    "imgplot = plt.imshow(actual_196_70[:,:], cmap='gray')\n",
    "plt.title(\"Actual Segmentation\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python",
   "version": "3.10.9 (tags/v3.10.9:1dd9be6, Dec  6 2022, 20:01:21) [MSC v.1934 64 bit (AMD64)]"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "602bbe8e8fe3127dc2d2bb8f57c49cfbb507be2cc0e090f10eb4c4d0af354b8b"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}