{"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"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Crop single cell images in the HPA dataset\n\nThis notebook shows an approach to crop images to obtain single-cell level images from the Human Protein Atlas (HPA)-Single Cell Classification competition. Functions and other ideas were used from other Kaggle notebooks listed below. \n\n## 1 Transform images in the training set to single-cell images\nFor a single image and using the cell segmentations mask from CellSegmentator we can generate a RGB image for each cell in the image as follows:\n\n1. Segment image using CellSegmentator.\n2. Add bounding boxes using cell mask from the first step. \n3. Crop using the bounding boxes.\n\nRepeating this process for all images generates a new dataset containing single-cell level images for each original image. Additionally we should generate a new .csv file containing the new single-cell identifiers and their labels. \n\n## 2 General comments\nIn this notebook only the first 10 images in the train set were used to demonstrate how the cropping process works. However the train set is so large it will require more computing power to complete the whole task in a reasonable time. Some image quality evaluation could be implemented, for example removing images that are too small or that only contain a small area of a single cell could be removed. Maybe this is not the most efficient way to complete this task, so any suggestions or comments are always welcome.  \n\n\nReference Notebooks:\n\n(1) [HPA_segmentation_and_BBOXES](https://www.kaggle.com/philipjamessullivan/hpa-segmentation-and-bboxes) by philipjamessullivan and acqua.\n\n(2) [Crop images using bounding box](https://www.kaggle.com/whizzkid/crop-images-using-bounding-box) by Ashutosh Chandra.\n\n(3) [pretrained ResNet34 with RGBY (0.460 public LB)](https://www.kaggle.com/iafoss/pretrained-resnet34-with-rgby-0-460-public-lb) by Iafoss.\n\n\n","metadata":{}},{"cell_type":"code","source":"# libraries\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport glob\nimport warnings \nimport os,gc,cv2\nimport shutil\nimport random\nfrom tqdm.notebook import tqdm\nfrom PIL import Image, ImageDraw\n\n%matplotlib inline\nwarnings.filterwarnings('ignore')\n\n# HPA single-cell image segmentation \n!pip install https://github.com/CellProfiling/HPA-Cell-Segmentation/archive/master.zip\nimport hpacellseg.cellsegmentator as cellsegmentator\nfrom hpacellseg.utils import label_cell, label_nuclei","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# directories \nDIR = '../input/hpa-single-cell-image-classification/'\nos.listdir(DIR)\n# csv files:\nCSV_PATH=os.path.join(DIR,'train.csv') \n# image files:\nIMG_FOLDER_PATH=os.path.join(DIR, 'train/')\n# create new directories in the  output directories for masks and cropped images\nCROPS = './sc_cropped_img_HPA/'\nif not os.path.exists(CROPS):\n       os.makedirs(CROPS)        \n# directory to save dataframe containing new single cell image identifiers and labels\nNEW_DF='./train_sc_df.csv'","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# read id/label csv to array\nid_labels_array=pd.read_csv(CSV_PATH)\n# fix labels (convert to arrays)\nid_labels_array[\"Label\"]=id_labels_array[\"Label\"].apply(lambda x:list(map(int, x.split(\"|\"))))     \n# create list of all ids from the id_labels_array\nids_images=(id_labels_array[\"ID\"]).tolist()\n# create dictionary of all unique labels from the id_labels_array\nlabels=id_labels_array.set_index('ID').T.to_dict('list')\n# fix the dictionary format (2d arrays to 1D arrays)\nlabels = {num: labels[0] for num, labels in labels.items()}","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Define CellSegmentator class\n# [source: https://github.com/CellProfiling/HPA-Cell-Segmentation]\n#-----------------------------------------------------------\n\n# [1] path to the nuclei model weights:\nNUC_MODEL = './nuclei-model.pth'\n#-----------------------------------------------------------\n# [2] path to the cell model weights:\nCELL_MODEL = './cell-model.pth'\n#-----------------------------------------------------------\n# [3] scale_factor: determines how much the images should be \n# scaled before being fed to the models. For HPA Cell images, \n# a value of 0.25 (default) is good.\n#-----------------------------------------------------------\n# [4] device: Inform Torch which device to put the model on. \n#Valid  values are ‘cpu’ or ‘cuda’ or pointed cuda device  \n# like 'cuda:0’. Defaults to cuda.\n#-----------------------------------------------------------\n# [5] padding: If True, add some padding before feeding the \n# images to the neural networks. This is not required but \n# can make segmentations, especially cell segmentations,\n# more accurate. Defaults to False. Note: If you have issues \n# running the segmentation due to image dimensions, setting \n# padding to True may help.\n#-----------------------------------------------------------\n# [6] multi_channel_model: If True, use the pretrained \n# three-channel version of the model. Having this set to \n# True gives you better cell segmentations but requires \n# you to give the model endoplasmic reticulum images as \n# part of the cell segmentation. Otherwise, the version \n# trained with only two channels, microtubules and nuclei, \n# will be used. Defaults to True\n#-----------------------------------------------------------\nsegmentator = cellsegmentator.CellSegmentator(\n    NUC_MODEL,\n    CELL_MODEL,\n    scale_factor=0.25,\n    device=\"cpu\",\n    padding=True,\n    multi_channel_model=True,\n)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# functions:\n\ndef open_rgby(img_id, img_dir): \n    '''\n    Reads individual images\n    of different filters (R, G, B, Y)\n    and stack them.\n    ---------------------------------\n    Arguments:\n    img_id -- image file identifier\n    img_dir -- the path containing all image files\n    Returns:\n    stacked (RGBY) image '''\n    colors=['red','green','blue','yellow']\n    flags=cv2.IMREAD_UNCHANGED\n    img=[cv2.imread(os.path.join(img_dir, img_id+'_'+color+'.png'), flags).astype(np.float32)/255 for color in colors]\n    return np.stack(img, axis=-1)\n\n\ndef generate_single_mask(img_id, img_dir):\n    '''\n    Computes cell mask for a single image\n    using HPA CellSegmentator tool\n    ---------------------------------\n    Arguments:\n    img_id -- image file identifier\n    img_dir -- the path containing all image files\n    Returns:\n    cell_mask as ndarray'''\n    path=os.path.join( img_dir + img_id )\n    ch_r = Image.open(path+\"_red.png\")\n    ch_y = Image.open(path+\"_yellow.png\")\n    ch_b = Image.open(path+\"_blue.png\")\n    nuc_segmentations = segmentator.pred_nuclei([np.asarray(ch_b)])\n    cell_segmentations = segmentator.pred_cells([\n                [np.asarray(ch_r)],\n                [np.asarray(ch_y)],\n                [np.asarray(ch_b)] ])\n    cell_nuclei_mask, cell_mask = label_cell(nuc_segmentations[0], cell_segmentations[0])\n    cell_mask = np.uint8(cell_mask)\n    return cell_mask\n\n\ndef bbox(mask,cell_id):\n    '''\n    Finds the bounding box for each \n    single cell in an image\n    ---------------------------------\n    Arguments:\n    mask -- cell mask for single image as ndarray\n    cell_id -- cell identifier (integer) from HPA CellSegmentator\n    Returns:\n    bounding-box coordinates as list: xmin, ymin, xmax, ymax'''\n    a = np.where(mask == cell_id)\n    xmin, ymin, xmax, ymax = np.min(a[1]), np.min(a[0]), np.max(a[1]), np.max(a[0])\n    return xmin, ymin, xmax, ymax\n\n\ndef crop_sc_imgs(imgid_list, dtrain, img_dir, out_dir, csv_dir):\n    '''\n    Crops images to obtain \n    single-cell images\n    ---------------------------------\n    Arguments:\n    imgid_list -- image file identifier\n    dtrain -- train set dataframe containing image indentifiers and labels\n    img_dir -- the path containing all image files\n    out_dir -- path for saving single-cell cropped images\n    csv_dir -- path to save new train set dataframe\n    Returns:\n    stacked (RGBY) images for cells in each image\n    new dataframe containing cell image identifiers and their labels\n    '''\n    new_ids=[]\n    new_labels=[]\n    for imgid in tqdm(imgid_list):\n        img=open_rgby(img_id=imgid, img_dir=img_dir)\n        # read image level multi-labels\n        df_img=dtrain[ dtrain[\"ID\"]==imgid ]\n        cell_label=df_img[\"Label\"].iloc[0]\n        # image cell mask\n        cell_mask=generate_single_mask(img_id=imgid, img_dir=img_dir)\n        cell_mask_flattened=np.ravel(cell_mask)\n        cell_ids=set(cell_mask_flattened); cell_ids.remove(0)\n        for cell_id in cell_ids:\n            cell_img_id=imgid+'_sc_'+str(cell_id)\n            new_ids.append(cell_img_id)\n            new_labels.append(cell_label)\n            xmin, ymin, xmax, ymax = bbox(cell_mask,cell_id)\n            new_img=img[ymin:ymax, xmin:xmax]\n            new_img=cv2.convertScaleAbs(new_img, alpha=(255.0))\n            cv2.imwrite( os.path.join(out_dir+ cell_img_id +'.png'), new_img )\n    new_df=pd.DataFrame({'ID': new_ids, 'Label': new_labels })\n    new_df.to_csv(csv_dir, index=False)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# reads train set dataframe again\ntrain_df=pd.read_csv(CSV_PATH)\n# select only 10 images and use the function\nimages=ids_images[0:10]\ncrop_sc_imgs(imgid_list=images, dtrain=train_df, img_dir=IMG_FOLDER_PATH, out_dir=CROPS, csv_dir=NEW_DF)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# checks old and new dataframes\nnew_train_df=pd.read_csv(NEW_DF)\nprint( train_df.head(), '\\n' )\nprint( new_train_df.head(), '\\n' )","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# plot a randomly selected single-cell cropped image\nsc_imgs=glob.glob( os.path.join(CROPS + '*.png') )\nrandom_img=random.sample(sc_imgs, 1)[0]\nex_sc_img=cv2.imread(random_img)\nplt.imshow(ex_sc_img)\nplt.axis('off')\nplt.title('single-cell image example')   \nplt.show()  ","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}