{"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":"code","source":"!pip install ../input/hpacellsegmentatormaster/HPA-Cell-Segmentation-master/\n!pip install ../input/hpapytorchzoozip/pytorch_zoo-master/\n!pip install ../input/mmdetection-v280/src/mmpycocotools-12.0.3/mmpycocotools-12.0.3/\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\n\nimport os\nimport sys\nimport cv2 \nimport math\nfrom tqdm import tqdm\nimport pickle\nfrom itertools import groupby\nimport matplotlib.pyplot as plt\nimport os\nimport base64\nimport typing as t\nimport zlib\nimport networkx as nx\nimport random\nfrom timeit import default_timer as timer\nrandom.seed(0)\nfrom tqdm import tqdm\n\nPRINTINFO=False\nSHOWIMAGES=False\n\nfrom pycocotools import mask as mutils\nfrom pycocotools import _mask as coco_mask\n\nfrom skimage import feature\nfrom joblib import dump, load\n\nexp_name = \"v1\"\nconf_name = \"Dataministic-with-CellSeg\"\nmodel_name = 'CellSegOpenCV'\nROOT = '../input/hpa-single-cell-image-classification/'\n\n###### Change here between test and train\n\ntrain_or_test = 'test'\n\n######\n\ndebug = False\ndebug_show_nucleolus_masks=False\ndebug_show_input_images=False\nif train_or_test == 'test':\n data_df = pd.read_csv(os.path.join(ROOT, 'sample_submission.csv'))\n\n if len(data_df) == 559:\n debug = True\n #one_frame = [data_df.loc[(data_df['ID']== '0173029a-161d-40ef-af28-2342915b22fb') |\\\n # (data_df['ID']=='00c9a1c9-2f06-476f-8b0d-6d01032874a2')]]\n\n data_df = data_df[:3]\n #data_df = pd.concat(one_frame)\n #print(data_df)\n else:\n debug = False\n \nelse:\n debug = True\n train_df = pd.read_csv(os.path.join(ROOT, 'train.csv'))\n\n if train_or_test == 'traintest':\n of_each = 3\n tdf0 = train_df.loc[train_df['Label'] == '0'].tail(of_each)\n tdf1 = train_df.loc[train_df['Label'] == '1'].tail(of_each)\n tdf2 = train_df.loc[train_df['Label'] == '2'].tail(of_each)\n tdf3 = train_df.loc[train_df['Label'] == '3'].tail(of_each)\n tdf4 = train_df.loc[train_df['Label'] == '4'].tail(of_each)\n tdf5 = train_df.loc[train_df['Label'] == '5'].tail(of_each)\n tdf6 = train_df.loc[train_df['Label'] == '6'].tail(of_each)\n tdf7 = train_df.loc[train_df['Label'] == '7'].tail(of_each)\n tdf8 = train_df.loc[train_df['Label'] == '8'].tail(of_each)\n tdf9 = train_df.loc[train_df['Label'] == '9'].tail(of_each)\n tdf10 = train_df.loc[train_df['Label'] == '10'].tail(of_each)\n tdf11 = train_df.loc[train_df['Label'] == '11'].tail(of_each)\n tdf12 = train_df.loc[train_df['Label'] == '12'].tail(of_each)\n tdf13 = train_df.loc[train_df['Label'] == '13'].tail(of_each)\n tdf14 = train_df.loc[train_df['Label'] == '14'].tail(of_each)\n tdf15 = train_df.loc[train_df['Label'] == '15'].tail(of_each)\n tdf16 = train_df.loc[train_df['Label'] == '16'].tail(of_each)\n tdf17 = train_df.loc[train_df['Label'] == '17'].tail(of_each)\n tdf18 = train_df.loc[train_df['Label'] == '18'].tail(of_each)\n else:\n of_each = 20\n tdf0 = train_df.loc[train_df['Label'] == '0'].head(of_each)\n tdf1 = train_df.loc[train_df['Label'] == '1'].head(of_each)\n tdf2 = train_df.loc[train_df['Label'] == '2'].head(of_each)\n tdf3 = train_df.loc[train_df['Label'] == '3'].head(of_each)\n tdf4 = train_df.loc[train_df['Label'] == '4'].head(of_each)\n tdf5 = train_df.loc[train_df['Label'] == '5'].head(of_each)\n tdf6 = train_df.loc[train_df['Label'] == '6'].head(of_each)\n tdf7 = train_df.loc[train_df['Label'] == '7'].head(of_each)\n tdf8 = train_df.loc[train_df['Label'] == '8'].head(of_each)\n tdf9 = train_df.loc[train_df['Label'] == '9'].head(of_each)\n tdf10 = train_df.loc[train_df['Label'] == '10'].head(of_each)\n tdf11 = train_df.loc[train_df['Label'] == '11'].head(of_each)\n tdf12 = train_df.loc[train_df['Label'] == '12'].head(of_each)\n tdf13 = train_df.loc[train_df['Label'] == '13'].head(of_each)\n tdf14 = train_df.loc[train_df['Label'] == '14'].head(of_each)\n tdf15 = train_df.loc[train_df['Label'] == '15'].head(of_each)\n tdf16 = train_df.loc[train_df['Label'] == '16'].head(of_each)\n tdf17 = train_df.loc[train_df['Label'] == '17'].head(of_each)\n tdf18 = train_df.loc[train_df['Label'] == '18'].head(of_each) \n \n frames=[\n tdf0,\n tdf1,\n tdf2,\n tdf3,\n tdf4,\n tdf5,\n tdf6,\n tdf7,\n tdf8,\n tdf9,\n tdf10,\n tdf11,\n tdf12,\n tdf13,\n tdf14,\n tdf15,\n tdf16,\n tdf17,\n tdf18 \n ]\n \n \n first_frame = [tdf16]\n \n #one_frame = [train_df.loc[train_df['ID']== '49b573a7-a887-499e-b96b-28fbe3fae7f7'].head(1)]\n one_frame = [train_df.loc[train_df['ID']== '4870f3ee-bb9b-11e8-b2b9-ac1f6b6435d0']]\n\n #data_df = pd.concat(frames)\n data_df = pd.concat(first_frame)\n #data_df = pd.concat(one_frame)\n \n #tdf0nuc = analyzeNucleusImages(tdf0, train_or_test)\n #tdf1nuc = analyzeNucleusImages(tdf1, train_or_test)\n #tdf2nuc = analyzeNucleusImages(tdf2, train_or_test)\n #tdf3nuc = analyzeNucleusImages(tdf3, train_or_test)\n #tdf4nuc = analyzeNucleusImages(tdf4, train_or_test)\n #tdf5nuc = analyzeNucleusImages(tdf5, train_or_test) \nprint(data_df) \nprint(debug)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def encode_binary_mask(mask: np.ndarray) -> t.Text:\n \"\"\"Converts a binary mask into OID challenge encoding ascii text.\"\"\"\n\n # check input mask --\n if mask.dtype != np.bool:\n raise ValueError(\n \"encode_binary_mask expects a binary mask, received dtype == %s\" %\n mask.dtype)\n\n mask = np.squeeze(mask)\n if len(mask.shape) != 2:\n raise ValueError(\n \"encode_binary_mask expects a 2d mask, received shape == %s\" %\n mask.shape)\n\n # convert input mask to expected COCO API input --\n mask_to_encode = mask.reshape(mask.shape[0], mask.shape[1], 1)\n mask_to_encode = mask_to_encode.astype(np.uint8)\n mask_to_encode = np.asfortranarray(mask_to_encode)\n\n # RLE encode mask --\n encoded_mask = coco_mask.encode(mask_to_encode)[0][\"counts\"]\n\n # compress and base64 encoding --\n binary_str = zlib.compress(encoded_mask, zlib.Z_BEST_COMPRESSION)\n base64_str = base64.b64encode(binary_str)\n return base64_str.decode()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!ls ../input/hpa-single-cell-image-classification/test/49b573a7-a887-499e-b96b-28fbe3fae7f7* \n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def printinfo (infostring):\n if PRINTINFO:\n print(infostring)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def find_nuc_contours (nu):\n \n nugauss = cv2.GaussianBlur(nu, (5, 5), 0)\n nuc_thresh = cv2.adaptiveThreshold(nugauss,255,cv2.ADAPTIVE_THRESH_MEAN_C, cv2.THRESH_BINARY,201,0)\n \n \n # Find the contours on the nuclei image\n contours, hierarchy = cv2.findContours(nuc_thresh, cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE)\n \n # Disregard any contours that are too small to be a nucleus\n i=0\n num_nuclei=0\n nuc_cont=[]\n nuc_midpoint=[]\n nuc_rect=[]\n nuc_mask = np.zeros(nu.shape, dtype='uint8')\n \n for x in contours:\n contlength = cv2.arcLength(contours[i],True)\n nucarea = cv2.contourArea(contours[i])\n approx = cv2.approxPolyDP(x,5,True)\n hull = cv2.convexHull(approx)\n hullArcLength = cv2.arcLength(hull,True)\n\n\n printinfo(f'nucleus {num_nuclei} has arc length {contlength}, hullArcLength {hullArcLength},area {nucarea} ')\n \n if ( hullArcLength > 300 and hierarchy[0][i][3] == -1 ):\n nuc_cont.append(hull)\n cv2.drawContours(nuc_mask, nuc_cont,num_nuclei,num_nuclei+1,cv2.FILLED,cv2.LINE_8)\n num_nuclei+=1\n \n M = cv2.moments(contours[i])\n cx = int(M['m10']/M['m00'])\n cy = int(M['m01']/M['m00'])\n nuc_midpoint.append((cx,cy))\n \n rect = cv2.boundingRect(contours[i])\n nuc_rect.append(rect)\n\n #printinfo('Chosen')\n i=i+1\n #display_image_big(nu)\n #display_image_big(nuc_mask)\n return nuc_cont, nuc_midpoint, nuc_mask, nuc_rect","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def show_mt_er_nu(mt, er, nu):\n \n fig, ax = plt.subplots(1,4, figsize=(20,50))\n ax[0].imshow(mt)\n ax[0].axis('off')\n\n ax[1].imshow(er)\n ax[1].axis('off')\n\n ax[2].imshow(nu)\n ax[2].axis('off')\n\n ax[3].imshow(er)\n ax[3].axis('off')\n\n plt.show()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def generateLayoutGraph(image_data, nuc_cont, nuc_midpoint, nuc_rect):\n cell_graph = nx.Graph()\n\n for j in range(len(nuc_cont)):\n cell_graph.add_node(j)\n \n return cell_graph","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def imageLevelStats( img, image_data):\n maxval=img.max()\n minval=img.min()\n sumval=img.sum()\n hist_im = np.histogram(img, bins=[0,31,63,95,127,159,191,223,255])\n hist256=np.histogram(img, bins=256)\n\n maxfreq=0\n for i in range(256):\n if hist256[0][i]>maxfreq:\n maxfreq=hist256[0][i]\n mostfreq=i\n medianval = 0\n \n return maxval, minval, sumval, maxfreq, mostfreq, medianval\n\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def choose_threshold(img, max_pixels_above):\n hist256=np.histogram(img, bins=256)\n\n numpix=0\n for i in range(255,-1,-1):\n numpix = numpix + hist256[0][i]\n if numpix > max_pixels_above:\n break\n return i","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def find_cells( mt, er, nu, nuc_cont, nuc_mask, nuc_midpoint, cell_graph):\n BOXWIDTH = 5\n RAYLENGTH=500\n BLACK_THRESHOLD=10\n BLUE_THRESHOLD=25\n MIN_BOX_AREA=25\n ryb = cv2.merge([mt, er, nu])\n\n #displaycopy= np.copy(ryb)\n #cv2.drawContours(displaycopy,nuc_cont,-1,(200,200,255),9)\n #do_display_image_big(displaycopy)\n image_new_contours=[]\n imagesize=np.array(nu.shape)\n\n #print('imaggesize', imagesize)\n\n for j in range(len(nuc_cont)):\n printinfo(f'Contour {j}')\n cont = nuc_cont[j]\n #print(nuc_cont[j].shape)\n printinfo(f'midpoint {nuc_midpoint[j]}')\n mid = np.array(nuc_midpoint[j])\n #print('mid',mid)\n #print(type(nuc_midpoint[j]))\n #for k in range(len(cont)-1):\n newcont=[]\n\n \n for k in range(len(cont)):\n printinfo(f'Contour {j} Segment {k}')\n edge_start= nuc_cont[j][k][0]\n if (k == len(cont)-1):\n # end of the last segment is the start point\n edge_end= nuc_cont[j][0][0]\n else:\n edge_end= nuc_cont[j][k+1][0]\n\n adjustedEdgeStart=edge_start \n adjustedEdgeEnd=edge_end \n printinfo(f'edge_start {edge_start} edge_end {edge_end}')\n edgevec = edge_end - edge_start\n if (edgevec[0]**2+edgevec[1]**2 < 25):\n printinfo('Warning: Very small contour segment')\n #print('edgevec',edgevec)\n midvec = edge_start + 0.5*edgevec\n #print('midvec',midvec)\n perpvec = midvec - mid\n #print(perpvec)\n scalefac=1/abs(perpvec).sum()\n perpvecW=perpvec* BOXWIDTH *scalefac\n perprayvec=perpvec*RAYLENGTH*scalefac\n #print('perpvecs')\n #print(perpvecW)\n #print(perprayvec)\n boxContour=np.array([edge_start,edge_end,edge_end,edge_start])\n\n #boxContour=np.array([nuc_cont[j][k][0],nuc_cont[j][k+1][0],nuc_cont[j][k+1][0],nuc_cont[j][k][0]])\n #rayBox=np.array([nuc_cont[j][k][0],nuc_cont[j][k+1][0],nuc_cont[j][k+1][0],nuc_cont[j][k][0]])\n\n #print('boxContour')\n #print(boxContour)\n #print(boxContour.shape)\n boxContour[2]=boxContour[2]+perpvecW\n boxContour[3]=boxContour[3]+perpvecW\n # Keep the floating values in boxContourF so that small increments don't get lost in rounding to pixels\n boxContourF=boxContour\n # Continue with rounded values in boxContour\n boxContour=boxContourF.astype(int)\n \n #rayBox[2]=rayBox[2]+perprayvec\n #rayBox[3]=rayBox[3]+perprayvec\n #print('rayBox')\n #print(rayBox)\n #clippedRayBox=np.fmin(rayBox,imagesize)\n #print('clippedRayBox')\n #print(clippedRayBox)\n\n bluefound = False\n blackfound = False\n boxContourAreaInitial=cv2.contourArea(boxContour)\n #print('boxContourAreaInitial', boxContourAreaInitial)\n if boxContourAreaInitial < MIN_BOX_AREA:\n printinfo('Error/ToDo remove small contour edges') \n \n for n in range(1,RAYLENGTH//BOXWIDTH):\n boxContourClipLow = np.fmax(boxContour,[0,0]) \n boxContourClipped = np.fmin(boxContourClipLow,imagesize)\n boxContourArea=cv2.contourArea(boxContourClipped)\n #print('boxContourArea', boxContourArea)\n \n if boxContourClipped[2][0] < boxContour[2][0] or boxContourClipped[3][0] < boxContour[3][0]:\n printinfo('hit right image boundary')\n if boxContourClipped[2][1] < boxContour[2][1] or boxContourClipped[3][1] < boxContour[3][1]:\n printinfo('hit bottom image boundary')\n if boxContourClipped[2][0] > boxContour[2][0] or boxContourClipped[3][0] > boxContour[3][0]:\n printinfo('hit left image boundary')\n if boxContourClipped[2][1] > boxContour[2][1] or boxContourClipped[3][1] > boxContour[3][1]:\n printinfo('hit top image boundary') \n\n if bluefound or blackfound or boxContourArea < MIN_BOX_AREA:\n break\n # boxContour is now clipped so that it remains in the image\n boxContour=boxContourClipped\n boxContourF=np.copy(boxContourClipped) \n #if j==0 and k == 3 :\n #print('boxContour')\n #print(boxContour) \n adjustedEdgeStart=boxContour[0]\n adjustedEdgeEnd=boxContour[1]\n x,y,w,h = cv2.boundingRect(boxContour)\n #print(x,y,w,h)\n enclosingSlice= np.copy(ryb[y:y+h, x:x+w])\n # sliceBoxContour is the contour translated from image coordinates to \n # a slice of the image containing the bounding box of the contour\n sliceBoxContour= boxContour - [x,y]\n #print('sliceBoxContour Shape',sliceBoxContour.shape )\n #print(sliceBoxContour)\n boxMask = np.zeros(enclosingSlice.shape, dtype='uint8')\n cv2.drawContours(boxMask,[sliceBoxContour],-1,(255,255,255),cv2.FILLED,cv2.LINE_8)\n numMaskPixels = np.sum(boxMask==255)/3\n #print('enclosingSlice',enclosingSlice.shape)\n enclosingSlice[boxMask==(0,0,0)]=0\n #do_display_image_big(boxMask)\n #do_display_image_big(enclosingSlice)\n \n #print('numMaskPixels',numMaskPixels)\n #print('boxMask.shape',boxMask.shape)\n #print('enclosingSlice.shape',enclosingSlice.shape)\n #p.set_printoptions(threshold=200)\n #np.set_printoptions(threshold=sys.maxsize)\n \n cropped_mt_sum, cropped_er_sum, cropped_nu_sum = np.sum(enclosingSlice, axis=(0,1) )\n \n #print('cropped_mt_sum',cropped_mt_sum)\n #print('cropped_er_sum',cropped_er_sum)\n #print('cropped_nu_sum',cropped_nu_sum)\n \n if numMaskPixels > 0:\n mt_erMean = (cropped_mt_sum+cropped_er_sum)/numMaskPixels\n nuMean = cropped_nu_sum/numMaskPixels\n #if j==1 and k == 5 :\n #print('numMaskPixels',numMaskPixels)\n #print('cropped_mt_sum',cropped_mt_sum)\n #print('cropped_er_sum',cropped_er_sum)\n #print('cropped_nu_sum',cropped_nu_sum)\n #print(f'nuMean:{nuMean} mt_erMean:{mt_erMean}')\n if nuMean > BLUE_THRESHOLD:\n\n # Take a slice from the nucleus mask corresponding to the \n # current box where blue was found\n # Look at the nucleus mask values to see which nucleus we hit\n # sum the non-zero values and find the mean\n # \n nuc_maskSlice= np.copy(nuc_mask[y:y+h, x:x+w])\n boxMask0 = boxMask[:,:,1]\n nuc_maskSlice[boxMask0==0]=0\n contours_hit=[]\n cont1=np.amax(nuc_maskSlice)\n if cont1 != j+1 and cont1 != 0 and cont1 != 255:\n contours_hit.append(cont1)\n cell_graph.add_edge(j, cont1-1)\n nuc_maskSlice[nuc_maskSlice==0]=255\n cont2=np.amin(nuc_maskSlice)\n if (cont2 != j+1) and (cont2 != cont1) and cont1 != 0 and cont1 != 255:\n contours_hit.append(cont2)\n cell_graph.add_edge(j, cont2-1)\n if len(contours_hit)>0:\n bluefound=True\n if len(contours_hit)>1:\n printinfo(\"Error/Todo: hit two nuclei\")\n printinfo(f'hit contours {contours_hit}')\n\n else:\n if mt_erMean < BLACK_THRESHOLD:\n blackfound=True\n else:\n printinfo('!!! Error: unexpected empty mask !!!')\n boxContourF=boxContourF+perpvecW\n boxContour=boxContourF.astype(int)\n #print(f'Terminated after {n} iterations')\n if bluefound:\n printinfo('Blue found')\n printinfo(f'edge_start:{edge_end}')\n printinfo(f'adjustedEdgeStart: {adjustedEdgeEnd}')\n newstart=(adjustedEdgeStart+edge_start)/2\n newend=(adjustedEdgeEnd+edge_end)/2\n adjustedEdgeStart = newstart.astype(int)\n adjustedEdgeEnd = newend.astype(int)\n printinfo('f {newstart} {newend}')\n if blackfound:\n printinfo('Black found')\n if boxContourArea < MIN_BOX_AREA:\n printinfo('Hit edge of image')\n printinfo(f'Adjusted edge is {adjustedEdgeStart} , {adjustedEdgeEnd}')\n newcont.append(adjustedEdgeStart)\n newcont.append(adjustedEdgeEnd)\n \n cnt=np.array(newcont)\n image_new_contours.append(cnt)\n #print('cnt.shape', cnt.shape )\n #print(cnt)\n #print('cnt.type', cnt.dtype )\n cv2.drawContours(ryb,[cnt],-1,(200,200,(255-10*(min(j,20)))),3)\n\n\n #print('ryb.shape', ryb.shape )\n #cv2.drawContours(ryb,[cnt],-1,(200,200,255),9) \n #do_display_image_big(ryb)\n \n return image_new_contours","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def do_display_image_big(im_data, image_data=None):\n\n dpi = 62\n \n height = im_data.shape[0]\n width = im_data.shape[1]\n\n # What size does the figure need to be in inches to fit the image?\n figsize = width / float(dpi), height / float(dpi)\n\n # Create a figure of the right size with one axes that takes up the full figure\n fig = plt.figure(figsize=figsize)\n fig.set_dpi(dpi)\n ax = fig.add_axes([0, 0, 1, 1])\n\n # Hide spines, ticks, etc.\n ax.axis('off')\n\n # Display the image.\n ax.imshow(im_data)\n \n if not image_data == None:\n for cell in image_data['cells']:\n cellnum = cell['cellnum']\n x,y,w,h = cell['nuc_bbox']\n ax.text(x,y,cellnum , fontsize=80)\n\n plt.show()\n \ndef display_image_big(im_data):\n if SHOWIMAGES:\n do_display_image_big(im_data)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def newCellDict(image_id, image_label, x, y, width, height ):\n cell= {\n 'image_id': image_id,\n 'image_label': image_label,\n 'cell_results': [],\n 'cellnum': None,\n 'bbox': (x,y,width,height),\n 'nuc_bbox': None,\n 'nuc_midpoint': None,\n 'width': width,\n 'height': height,\n 'nucleoli': [],\n 'mask_count': None,\n 'pr_minval': None, \n 'pr_maxval': None, \n 'pr_sumval': None,\n 'pr_meanval': None,\n 'pr_variance': None,\n 'pr_mt_correlation': None,\n 'pr_er_correlation': None,\n 'mt_minval': None, \n 'mt_maxval': None, \n 'mt_sumval': None,\n 'mt_meanval': None,\n 'mt_variance': None,\n 'er_minval': None, \n 'er_maxval': None, \n 'er_sumval': None,\n 'er_meanval': None,\n 'er_variance': None,\n 'nu_minval': None, \n 'nu_maxval': None, \n 'nu_sumval': None,\n 'nu_meanval': None,\n 'nu_variance': None,\n # values over the full cell\n 'mask_count_fc': None,\n 'pr_minval_fc': None,\n 'pr_maxval_fc': None,\n 'pr_sumval_fc': None,\n 'pr_meanval_fc': None,\n 'pr_variance_fc': None,\n 'mt_minval_fc': None,\n 'mt_maxval_fc': None,\n 'mt_sumval_fc': None,\n 'mt_meanval_fc': None,\n 'mt_variance_fc': None,\n 'pr_mt_correlation_fc': None,\n # values over the nucleus\n 'nuc_mask_count': None,\n 'nuc_pr_minval': None, \n 'nuc_pr_maxval': None, \n 'nuc_pr_sumval': None,\n 'nuc_pr_meanval': None,\n 'nuc_pr_variance': None,\n 'nuc_mt_minval': None, \n 'nuc_mt_maxval': None, \n 'nuc_mt_sumval': None,\n 'nuc_mt_meanval': None,\n 'nuc_mt_variance': None,\n 'nuc_nu_minval': None, \n 'nuc_nu_maxval': None, \n 'nuc_nu_sumval': None,\n 'nuc_nu_meanval': None,\n 'nuc_nu_variance': None,\n 'nuc_nu_pr_correlation': None,\n # values over the nucleus inner rim\n 'inrim_mask_count': None,\n 'inrim_pr_minval': None,\n 'inrim_pr_maxval': None,\n 'inrim_pr_sumval': None,\n 'inrim_pr_meanval': None,\n 'inrim_pr_variance': None,\n 'inrim_nu_minval': None,\n 'inrim_nu_maxval': None,\n 'inrim_nu_sumval': None,\n 'inrim_nu_meanval': None,\n 'inrim_nu_variance': None,\n 'inrim_pr_nu_correlation': None,\n # values over the nucleus outer rim\n 'outrim_mask_count': None,\n 'outrim_pr_minval': None,\n 'outrim_pr_maxval': None,\n 'outrim_pr_sumval': None,\n 'outrim_pr_meanval': None,\n 'outrim_pr_variance': None,\n 'outrim_nu_minval': None,\n 'outrim_nu_maxval': None,\n 'outrim_nu_sumval': None,\n 'outrim_nu_meanval': None,\n 'outrim_nu_variance': None,\n 'outrim_pr_nu_correlation': None,\n 'aggresome': False,\n 'aggresome_size': None,\n #number of nuclear bodies\n 'nucbody': 0,\n 'nucbody_size': 0, \n #number of nuclear Speckles\n 'num_nucspeckles': 0,\n 'nucspeckles_size': 0 \n }\n return cell","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def newNucleolusDict(image_id, image_label):\n\n nucleolus= {\n 'image_id': image_id,\n 'image_label': image_label,\n 'mask_count': None,\n 'cellnum': None,\n 'arclength': None,\n 'area': None,\n 'bbox': None,\n 'width': None,\n 'height': None,\n 'nucleoli': [],\n 'pr_minval': None, \n 'pr_maxval': None, \n 'pr_sumval': None,\n 'pr_meanval': None,\n 'pr_variance': None,\n 'pr_nu_correlation' : None,\n 'nu_minval': None, \n 'nu_maxval': None, \n 'nu_sumval': None,\n 'nu_meanval': None,\n 'nu_variance': None\n }\n return nucleolus\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def newImageDict(image_id, trainlabel, width, height ):\n ann = {\n 'image_id': image_id,\n 'train_label': trainlabel,\n 'width': width,\n 'height': height,\n 'cells': [],\n '141617-pred': None,\n 'pr_minval': None, \n 'pr_maxval': None, \n 'pr_sumval': None,\n 'pr_medianval': None,\n 'pr_maxfreq': None,\n 'pr_mostfreq': None,\n 'mt_minval': None, \n 'mt_maxval': None, \n 'mt_sumval': None,\n 'mt_medianval': None,\n 'mt_maxfreq': None,\n 'mt_mostfreq': None,\n 'nu_minval': None, \n 'nu_maxval': None, \n 'nu_sumval': None,\n 'nu_medianval': None,\n 'nu_maxfreq': None,\n 'nu_mostfreq': None,\n 'er_minval': None, \n 'er_maxval': None, \n 'er_sumval': None,\n 'er_medianval': None,\n 'er_maxfreq': None,\n 'er_mostfreq': None,\n 'min_pr_er_correlation': None,\n 'max_pr_er_correlation' : None, \n 'mean_pr_er_correlation': None,\n 'min_pr_mt_correlation': None,\n 'max_pr_mt_correlation' : None, \n 'mean_pr_mt_correlation': None,\n 'min_pr_mt_correlation_fc': None,\n 'max_pr_mt_correlation_fc': None, \n 'mean_pr_mt_correlation_fc': None\n }\n return ann","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def countTrainingScores(image_data_list,falsePos,falseNeg,truePos,trueNeg ):\n\n \n \n for image in image_data_list:\n image_label=image['train_label']\n expected_labels=image_label.split('|')\n\n \n for cell in image['cells']:\n\n results = cell['cell_results']\n for res in results:\n found=False\n for exp in expected_labels:\n if int(exp) == int(res[0]):\n found = True\n if found:\n truePos[int(res[0])]+=1\n else:\n falsePos[int(res[0])]+=1\n for exp in expected_labels:\n found=False\n for res in results:\n if exp == res[0]:\n found=True\n if not found:\n falseNeg[int(exp)]+=1\n \n print('True positives:')\n print(truePos)\n print('False positives:')\n print(falsePos) \n print('False negatives:')\n print(falseNeg) ","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def calculateScores(image_data):\n\n num_cells = len(image_data['cells'])\n #p = image_data['141617-pred']\n #print(f'Image level pred is {p}')\n\n for cell in image_data['cells']:\n cell_scored = False\n if cell['aggresome']:\n cell['cell_results'].append((15,max(cell['aggresome_size']/150,1)))\n cell_scored = True\n\n \n max_mt_corr = max(cell['pr_mt_correlation_fc'], cell['pr_mt_correlation_fc'] )\n if cell['pr_er_correlation'] > 0.4 or max_mt_corr > 0.4 :\n # Correlation with red or yellow layer found\n if cell['pr_er_correlation'] > max_mt_corr:\n cell['cell_results'].append((6,cell['pr_er_correlation']))\n cell_scored = True\n\n else:\n if cell['nuc_mt_meanval'] > 100:\n cell['cell_results'].append((11,max_mt_corr))\n cell_scored = True\n\n else:\n cell['cell_results'].append((10,max_mt_corr))\n cell_scored = True\n \n \n for cell in image_data['cells']: \n max_mean = 0\n nucleolus_detected = False\n nucleolus_fcent_detected= False\n sum_nucleolus_mask_count=0\n sum_nucleolus_sums=0\n # Nucleolus\n for nucleolus in cell['nucleoli']:\n sum_nucleolus_mask_count += nucleolus['mask_count']\n sum_nucleolus_sums+=nucleolus['pr_sumval']\n if max_mean < nucleolus['pr_meanval']:\n max_mean = nucleolus['pr_meanval']\n if nucleolus['pr_meanval'] > 150 \\\n and nucleolus['pr_meanval'] > 1.5 * cell['pr_meanval'] :\n if nucleolus['pr_variance'] < 1000:\n nucleolus_detected = True\n else:\n nucleolus_fcent_detected= True \n \n nucleus_count_without_nucleoli = cell['nuc_mask_count']-sum_nucleolus_mask_count\n adjusted_nuc_mean = (cell['nuc_pr_sumval']-sum_nucleolus_sums)/nucleus_count_without_nucleoli\n\n npmeanval=cell['nuc_pr_meanval']\n im = cell['image_id']\n #if ( adjusted_nuc_mean < npmeanval):\n #print(f'Warning: {im} adjusted_nuc_mean {adjusted_nuc_mean} is less than cell nuc_pr_meanval {npmeanval}')\n \n if nucleolus_detected:\n cell['cell_results'].append((2,max_mean/255))\n cell_scored = True\n else:\n if nucleolus_fcent_detected and not cell['num_nucspeckles'] > 3:\n cell['cell_results'].append((3,max_mean/255))\n cell_scored = True\n else:\n if cell['num_nucspeckles'] > 0 :\n cell['cell_results'].append((4,min(cell['nucspeckles_size']/60,1)))\n cell_scored = True\n else:\n if cell['nucbody'] > 0 and cell['nucbody'] < 3:\n cell['cell_results'].append((5,min(cell['nucbody_size']/18,1)))\n cell_scored = True\n # note (todo): possible nuc fib centre for nucbody > 3\n else:\n inr=cell['inrim_pr_meanval']\n nucr=cell['nuc_pr_meanval']\n if inr > nucr:\n cell['cell_results'].append((1,(inr-nucr)/inr))\n cell_scored = True\n else:\n if cell['inrim_pr_meanval'] > cell['outrim_pr_meanval']*1.1:\n cell['cell_results'].append((0,adjusted_nuc_mean/255))\n cell_scored = True\n\n \n if image_data['141617-pred'] in ['14','16','17']:\n cell['cell_results'].append((image_data['141617-pred'],cell['pr_meanval']/255))\n cell_scored = True\n \n if not cell_scored:\n cell['cell_results'].append((18,1-(cell['pr_meanval']/255)))\n \n \n ","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def maskStats( image_layer, mask_image, masknum, comp_image, statsdebug=False ): \n # image_layer is a single colour layer (or slice of a layer) of the image # mask_image is an image (or slice) of the same size as image_layer with areas # set to different values indicating the number of the corresponding cell # as found in the ordering of the nucleus detection # The shape thus marked in the mask image can be, for instance: # - the area covered by the nucleus # - the area covered by the whole cell # - the area of the cell excluding the nucleus # - the area of a nucleolus # masknum is the value indicating the area of interest in the mask # comp_image is an image (or slice) of the same size to which image_layer is to be compared # Return the following calculated over all pixels with masknum at the same position in mask: # - minimum value of image_layer # - maximum value of image_layer # - sum of pixel values of image_layer # - mean value of image_layer # - variance of image_layer # - correlation coefficient between the layer and the comp_image\n\n # Work on a copies so as not to change the data external to the routine\n image_copy = np.copy(image_layer)\n comp_copy = np.copy(comp_image)\n \n if statsdebug:\n np.set_printoptions(threshold=sys.maxsize)\n print('image_copy')\n print(image_copy)\n print('comp_copy')\n print(comp_copy)\n print('mask_image')\n print(mask_image)\n \n # Find min (set the pixels not belonging to the area of interest to high values)\n image_copy[mask_image!=masknum]=255\n minval = image_copy.min()\n comp_copy[mask_image!=masknum]=255\n compmin = comp_copy.min() \n\n # Find max (set the pixels not belonging to the area of interest to zero)\n image_copy[mask_image!=masknum]=0\n maxval = image_copy.max()\n comp_copy[mask_image!=masknum]=0\n compmax = comp_copy.max()\n\n # Count the number of relevant pixels\n numMaskPixels = np.sum(mask_image==masknum)\n\n # Calculate sum and mean pixel values for image and comp_image\n sumval = image_copy.sum()\n compsum = comp_copy.sum()\n if numMaskPixels == 0:\n meanval = 0\n compmean = 0\n else:\n meanval = sumval / numMaskPixels\n compmean = compsum / numMaskPixels \n\n # Calculate variance, covariance for image and comp_image\n image_diff_mean = image_copy - meanval\n comp_diff_mean = comp_copy - compmean\n image_diff_mean[mask_image!=masknum]=0\n comp_diff_mean[mask_image!=masknum]=0\n \n image_sqdiffs = image_diff_mean**2\n comp_sqdiffs = comp_diff_mean**2\n\n if numMaskPixels == 0:\n im_variance = 0\n comp_variance = 0\n else: \n im_variance = np.sum(image_sqdiffs)/numMaskPixels\n comp_variance = np.sum(comp_sqdiffs)/numMaskPixels\n\n diff_prod = image_diff_mean*comp_diff_mean\n\n if numMaskPixels == 0:\n covar = 0\n else:\n covar=np.sum(diff_prod)/numMaskPixels\n\n if im_variance > 0 and comp_variance > 0:\n corr = covar/math.sqrt(im_variance*comp_variance)\n else:\n corr = 0\n\n return minval, maxval, sumval, meanval, im_variance, compmin, \\\n compmax, compsum, compmean, comp_variance, corr, numMaskPixels","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def analyseMicrotubules(image_data,contours,label,im,mt,pr, nuc_mask):\n\n min_correlation = 1\n max_correlation = -1\n sum_corr = 0 \n if label == '10':\n print(im)\n display_image_big(mt)\n display_image_big(pr)\n for j, (cell_contour, cell_data) in enumerate(zip(contours,image_data['cells'])):\n x,y,w,h = cell_data['bbox']\n full_image_mask=np.zeros(mt.shape, dtype='uint8')\n cv2.drawContours(full_image_mask, contours,j,255,cv2.FILLED,cv2.LINE_8)\n full_image_mask[nuc_mask==j+1]=0\n cropped_mt = mt[y:y+h, x:x+w]\n cropped_pr = pr[y:y+h, x:x+w]\n cropped_mask = full_image_mask[y:y+h, x:x+w]\n\n minval, maxval, sumval, meanval, im_variance, compmin, compmax,\\\n compsum, compmean, comp_variance, corr, mask_count \\\n = maskStats( cropped_pr, cropped_mask, 255, cropped_mt )\n\n cell_data['mask_count']=mask_count\n cell_data['pr_minval']=minval\n cell_data['pr_maxval']=maxval\n cell_data['pr_sumval']=sumval\n cell_data['pr_meanval']=meanval\n cell_data['pr_variance']=im_variance \n cell_data['mt_minval']=compmin\n cell_data['mt_maxval']=compmax\n cell_data['mt_sumval']=compsum\n cell_data['mt_meanval']=compmean\n cell_data['mt_variance']=comp_variance\n cell_data['pr_mt_correlation']=corr\n\n sum_corr = sum_corr + corr\n if corr < min_correlation:\n min_correlation = corr\n if corr > max_correlation:\n max_correlation = corr\n mean_corr = sum_corr/len(contours) \n image_data['min_pr_mt_correlation']= min_correlation\n image_data['max_pr_mt_correlation']= max_correlation \n image_data['mean_pr_mt_correlation']= mean_corr","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def analyseMicrotubulesFullCell(image_data,contours,label,im,mt,pr, nuc_mask):\n\n min_correlation = 1\n max_correlation = -1\n sum_corr = 0 \n if label == '10':\n printinfo(im)\n display_image_big(mt)\n display_image_big(pr)\n for j, (cell_contour, cell_data) in enumerate(zip(contours,image_data['cells'])):\n x,y,w,h = cell_data['bbox']\n full_image_mask=np.zeros(mt.shape, dtype='uint8')\n cv2.drawContours(full_image_mask, contours,j,255,cv2.FILLED,cv2.LINE_8)\n\n cropped_mt = mt[y:y+h, x:x+w]\n cropped_pr = pr[y:y+h, x:x+w]\n cropped_mask = full_image_mask[y:y+h, x:x+w]\n\n minval, maxval, sumval, meanval, im_variance, compmin, compmax,\\\n compsum, compmean, comp_variance, corr, mask_count \\\n = maskStats( cropped_pr, cropped_mask, 255, cropped_mt )\n\n cell_data['mask_count_fc']=mask_count\n cell_data['pr_minval_fc']=minval\n cell_data['pr_maxval_fc']=maxval\n cell_data['pr_sumval_fc']=sumval\n cell_data['pr_meanval_fc']=meanval\n cell_data['pr_variance_fc']=im_variance \n cell_data['mt_minval_fc']=compmin\n cell_data['mt_maxval_fc']=compmax\n cell_data['mt_sumval_fc']=compsum\n cell_data['mt_meanval_fc']=compmean\n cell_data['mt_variance_fc']=comp_variance\n cell_data['pr_mt_correlation_fc']=corr\n\n sum_corr = sum_corr + corr\n if corr < min_correlation:\n min_correlation = corr\n if corr > max_correlation:\n max_correlation = corr\n mean_corr = sum_corr/len(contours) \n image_data['min_pr_mt_correlation_fc']= min_correlation\n image_data['max_pr_mt_correlation_fc']= max_correlation \n image_data['mean_pr_mt_correlation_fc']= mean_corr","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def analyseER(image_data,contours,label,im,er,pr, nuc_mask):\n min_correlation = 1\n max_correlation = -1\n sum_corr = 0 \n\n for j, (cell_contour, cell_data) in enumerate(zip(contours,image_data['cells'])):\n x,y,w,h = cell_data['bbox']\n full_image_mask=np.zeros(er.shape, dtype='uint8')\n cv2.drawContours(full_image_mask, contours,j,255,cv2.FILLED,cv2.LINE_8)\n full_image_mask[nuc_mask==j+1]=0\n \n cropped_er = er[y:y+h, x:x+w]\n cropped_pr = pr[y:y+h, x:x+w]\n cropped_mask = full_image_mask[y:y+h, x:x+w]\n\n if False and debug:\n display_image_big(cropped_mask)\n display_image_big(cropped_pr)\n display_image_big(cropped_er)\n minval, maxval, sumval, meanval, im_variance, compmin, compmax,\\\n compsum, compmean, comp_variance, corr, mask_count \\\n = maskStats( cropped_pr, cropped_mask, 255, cropped_er )\n\n cell_data['mask_count']=mask_count\n cell_data['pr_minval']=minval\n cell_data['pr_maxval']=maxval\n cell_data['pr_sumval']=sumval\n cell_data['pr_meanval']=meanval\n cell_data['pr_variance']=im_variance \n cell_data['er_minval']=compmin\n cell_data['er_maxval']=compmax\n cell_data['er_sumval']=compsum\n cell_data['er_meanval']=compmean\n cell_data['er_variance']=comp_variance\n cell_data['pr_er_correlation']=corr\n\n sum_corr = sum_corr + corr\n if corr < min_correlation:\n min_correlation = corr\n if corr > max_correlation:\n max_correlation = corr\n mean_corr = sum_corr/len(contours) \n image_data['min_pr_er_correlation']= min_correlation\n image_data['max_pr_er_correlation']= max_correlation \n image_data['mean_pr_er_correlation']= mean_corr","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def analyseNU(image_data,contours,label,im,nu,pr, nuc_mask):\n \n\n for j, (cell_contour, cell_data) in enumerate(zip(contours,image_data['cells'])):\n x,y,w,h = cell_data['bbox']\n full_image_mask=np.zeros(nu.shape, dtype='uint8')\n cv2.drawContours(full_image_mask, contours,j,255,cv2.FILLED,cv2.LINE_8)\n full_image_mask[nuc_mask==j+1]=0\n \n cropped_nu = nu[y:y+h, x:x+w]\n cropped_pr = pr[y:y+h, x:x+w]\n cropped_mask = full_image_mask[y:y+h, x:x+w]\n\n if False and debug:\n display_image_big(cropped_mask)\n display_image_big(cropped_pr)\n display_image_big(cropped_nu)\n minval, maxval, sumval, meanval, im_variance, compmin, compmax,\\\n compsum, compmean, comp_variance, corr, mask_count \\\n = maskStats( cropped_pr, cropped_mask, 255, cropped_nu )\n\n cell_data['nu_minval']=compmin\n cell_data['nu_maxval']=compmax\n cell_data['nu_sumval']=compsum\n cell_data['nu_meanval']=compmean\n cell_data['nu_variance']=comp_variance \n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def analyseMTNU(image_data,nuc_cont,label,im,nu,mt, nuc_mask):\n \n\n for j, (nuc_contour, cell_data) in enumerate(zip(nuc_cont,image_data['cells'])):\n x,y,w,h = cell_data['nuc_bbox']\n n_mask=np.zeros(nu.shape, dtype='uint8')\n cv2.drawContours(n_mask, nuc_cont,j,255,cv2.FILLED,cv2.LINE_8)\n \n cropped_nu = nu[y:y+h, x:x+w]\n cropped_mt = mt[y:y+h, x:x+w]\n cropped_mask = n_mask[y:y+h, x:x+w]\n\n if False and debug:\n display_image_big(cropped_mask)\n display_image_big(cropped_mt)\n display_image_big(cropped_nu)\n minval, maxval, sumval, meanval, im_variance, compmin, compmax,\\\n compsum, compmean, comp_variance, corr, mask_count \\\n = maskStats( cropped_mt, cropped_mask, 255, cropped_nu )\n cell_data['nuc_mt_minval']=compmin\n cell_data['nuc_mt_maxval']=compmax\n cell_data['nuc_mt_sumval']=compsum\n cell_data['nuc_mt_meanval']=compmean\n cell_data['nuc_mt_variance']=comp_variance \n ","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def generateNucRims(nuc_cont, nuc_midpoint, nu, rimwidth):\n # generate contours just rimwidth radius wider than the nucleus contours\n image_rim_contours=[]\n for j ,(cont, mid) in enumerate(zip(nuc_cont,nuc_midpoint)):\n\n #print(f'Contour {j}')\n #print(f'midpoint {mid}')\n\n newcont=[]\n \n imagesize=np.array(nu.shape)\n \n for k in range(len(cont)):\n #print(f'Contour {j} Segment {k}')\n vertex= nuc_cont[j][k][0]\n #print('vertex', vertex)\n rayvec = vertex - mid\n #print('rayvec', rayvec)\n #print(type(mid))\n scalefac=1/abs(rayvec).sum()\n rayvecnew=(rayvec*(abs(rayvec).sum()+rimwidth)*scalefac).astype(int)\n vertexnew=rayvecnew+mid\n clipvertex_low = np.fmax(vertexnew,[0,0]) \n clipvertex = np.fmin(clipvertex_low,imagesize)\n newcont.append(clipvertex)\n #print('newcont', newcont)\n cnt=np.array(newcont)\n #print('nuc_cont',nuc_cont[j] )\n #print(type(cnt))\n #print('cnt', cnt)\n\n image_rim_contours.append(cnt) \n #print(f'len(image_rim_contours), {len(image_rim_contours)}')\n #full_image_mask=np.zeros(nu.shape, dtype='uint8') \n #cv2.drawContours(full_image_mask, image_rim_contours,-1,255,1)\n #cv2.drawContours(full_image_mask, nuc_cont,-1,200,1)\n #do_display_image_big(full_image_mask)\n \n\n return image_rim_contours","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def analyseNucleusInnerRims(image_data, nuc_cont, nuc_midpoint, nu, pr):\n INNER_NUCLEUS_RIMWIDTH = -10\n plus10cont=generateNucRims(nuc_cont, nuc_midpoint, nu, INNER_NUCLEUS_RIMWIDTH)\n\n for j, (nuc_contour, rim_contour, cell_data ) in \\\n enumerate(zip(nuc_cont, plus10cont, image_data['cells'])):\n\n full_image_mask=np.zeros(nu.shape, dtype='uint8')\n cv2.drawContours(full_image_mask, nuc_cont,j,255,cv2.FILLED,cv2.LINE_8)\n cv2.drawContours(full_image_mask, plus10cont,j,0,cv2.FILLED,cv2.LINE_8)\n\n #if debug:\n #do_display_image_big(full_image_mask) \n \n x,y,w,h = cell_data['nuc_bbox']\n \n \n #print(f'x,y,w,h:{x},{y},{w},{h}')\n \n #cv2.rectangle(full_image_mask, (x,y), (x+w,y+h), 200,2)\n #display_image_big(full_image_mask) \n cropped_mask = full_image_mask[y:y+h, x:x+w]\n cropped_nu = nu[y:y+h, x:x+w]\n cropped_pr = pr[y:y+h, x:x+w]\n \n\n #display_image_big(cropped_pr)\n #display_image_big(cropped_nu)\n \n minval, maxval, sumval, meanval, im_variance, compmin, compmax,\\\n compsum, compmean, comp_variance, corr, mask_count \\\n = maskStats( cropped_pr, cropped_mask, 255, cropped_nu )\n cell_data['inrim_mask_count']=mask_count\n cell_data['inrim_pr_minval']=minval\n cell_data['inrim_pr_maxval']=maxval\n cell_data['inrim_pr_sumval']=sumval\n cell_data['inrim_pr_meanval']=meanval\n cell_data['inrim_pr_variance']=im_variance \n cell_data['inrim_nu_minval']=compmin\n cell_data['inrim_nu_maxval']=compmax\n cell_data['inrim_nu_sumval']=compsum\n cell_data['inrim_nu_meanval']=compmean\n cell_data['inrim_nu_variance']=comp_variance\n cell_data['inrim_pr_nu_correlation']=corr","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def overlaps(x1,y1,w1,h1,x2,y2,w2,h2):\n xseparate = x1+w1 < x2 or x2+w2 < x1\n yseparate = y1+h1 < y2 or y2+h2 < y1\n return not xseparate and not yseparate\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def lookForAggresomes (image_data, pr):\n MAXWIDTH=150\n MINWIDTH=30\n AREAMAX=22500\n AREAMIN=400\n \n prcopy = np.copy(pr)\n retval, prthresh= cv2.threshold(prcopy, 100,255,cv2.THRESH_BINARY)\n kernel = np.ones((5, 5), 'uint8')\n e=cv2.erode(prthresh,kernel,iterations=2)\n ed=cv2.dilate(e,kernel,iterations=5)\n ede=cv2.erode(ed,kernel,iterations=2)\n #do_display_image_big(ede)\n \n contours, hierarchy = cv2.findContours(ede, cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE)\n #print(f'number of contours: {len(contours)}')\n for i, cont in enumerate(contours):\n area = cv2.contourArea(cont)\n #print(f'area {area}')\n x,y,w,h = cv2.boundingRect(cont)\n #print('x,y,w,h',x,y,w,h)\n if area < AREAMAX and area > AREAMIN and hierarchy[0][i][3] == -1 \\\n and w < MAXWIDTH and w > MINWIDTH \\\n and h > int(w*0.8) and int(h*0.8) < w:\n # found potential age\n # Check whether it touches the bounding box of a cell\n #print('found contour of potential aggresome size, checking overlap with outrim ...')\n for j, cell_data in enumerate(image_data['cells']):\n cellx, celly, cellw, cellh = cell_data['bbox']\n #print('checking overlap with ',x,y,w,h,cellx, celly, cellw, cellh )\n if overlaps(x,y,w,h,cellx, celly, cellw, cellh):\n #print('boxes overlap, checking masks')\n outrim_mask=cell_data['outrim_mask']\n binary_outrim_mask = outrim_mask.astype('bool')\n #do_display_image_big(binary_outrim_mask)\n\n full_image_mask=np.zeros(pr.shape, dtype='uint8')\n cv2.drawContours(full_image_mask, contours,i,255,cv2.FILLED,cv2.LINE_8)\n #do_display_image_big(full_image_mask)\n contour_cell_slice=full_image_mask[celly:celly+cellh, cellx:cellx+cellw]\n #do_display_image_big(contour_cell_slice)\n contour_cell_binmask=contour_cell_slice.astype(bool)\n #do_display_image_big(contour_cell_binmask)\n show_mask = np.logical_or(binary_outrim_mask,contour_cell_binmask)\n #do_display_image_big(show_mask)\n result_mask = np.logical_and(binary_outrim_mask,contour_cell_binmask)\n overlap_found=result_mask.max()\n if overlap_found:\n cell_data['aggresome']=True\n cell_data['aggresome_size'] = w","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def lookForNuclearBodies (image_data, pr, nuc_mask):\n MAXWIDTH=14\n MINWIDTH=3\n AREAMAX=200\n AREAMIN=10\n \n num_cells=len(image_data['cells'])\n thresh = choose_threshold(pr, 2*num_cells*AREAMAX)\n prcopy = np.copy(pr)\n #do_display_image_big(pr)\n retval, prthresh= cv2.threshold(prcopy, thresh,255,cv2.THRESH_BINARY)\n kernel = np.ones((5, 5), 'uint8')\n e=cv2.erode(prthresh,kernel,iterations=1)\n #do_display_image_big(e)\n ed=cv2.dilate(e,kernel,iterations=2)\n #do_display_image_big(ed)\n ede=cv2.erode(ed,kernel,iterations=1)\n #do_display_image_big(ede)\n \n \n\n num_pixels=int(ede.sum()/255)\n\n # only perform on sparse arrays (dull images), not to confuse with nucleoli\n if num_pixels < 2*num_cells*AREAMAX:\n contours, hierarchy = cv2.findContours(ede, cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE)\n #print(f'number of contours: {len(contours)}')\n for i, cont in enumerate(contours):\n area = cv2.contourArea(cont)\n #print(f'area {area}')\n x,y,w,h = cv2.boundingRect(cont)\n #print('x,y,w,h',x,y,w,h)\n if area < AREAMAX and area > AREAMIN and hierarchy[0][i][3] == -1 \\\n and w < MAXWIDTH and w > MINWIDTH \\\n and h > int(w*0.6) and int(h*0.6) < w:\n # found potential nuclear body\n # Check whether it touches the bounding box of a cell\n #print('found contour of potential aggresome size, checking overlap with outrim ...')\n for j, cell_data in enumerate(image_data['cells']):\n cellx, celly, cellw, cellh = cell_data['nuc_bbox']\n #print('checking overlap with ',x,y,w,h,cellx, celly, cellw, cellh )\n if overlaps(x,y,w,h,cellx, celly, cellw, cellh):\n #print('boxes overlap, checking masks')\n nuc_mask_slice=nuc_mask[celly:celly+cellh, cellx:cellx+cellw]\n binary_nuc_mask = nuc_mask_slice.astype('bool')\n #do_display_image_big(binary_outrim_mask)\n\n full_image_mask=np.zeros(pr.shape, dtype='uint8')\n cv2.drawContours(full_image_mask, contours,i,255,cv2.FILLED,cv2.LINE_8)\n #do_display_image_big(full_image_mask)\n contour_cell_slice=full_image_mask[celly:celly+cellh, cellx:cellx+cellw]\n #do_display_image_big(contour_cell_slice)\n contour_cell_binmask=contour_cell_slice.astype(bool)\n #do_display_image_big(contour_cell_binmask)\n show_mask = np.bitwise_xor(binary_nuc_mask,contour_cell_binmask)\n #do_display_image_big(show_mask)\n result_mask = np.logical_and(binary_nuc_mask,contour_cell_binmask)\n overlap_found=result_mask.max()\n if overlap_found:\n cell_data['nucbody']+=1\n cell_data['nucbody_size'] = w","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def lookForNuclearSpeckles (image_data, pr, nuc_mask):\n MAXWIDTH=50\n MINWIDTH=20\n AREAMAX=1500\n AREAMIN=100\n \n num_cells=len(image_data['cells'])\n thresh = choose_threshold(pr, 8*num_cells*AREAMAX)\n prcopy = np.copy(pr)\n #do_display_image_big(pr)\n retval, prthresh= cv2.threshold(prcopy, thresh,255,cv2.THRESH_BINARY)\n kernel = np.ones((5, 5), 'uint8')\n e=cv2.erode(prthresh,kernel,iterations=2)\n #do_display_image_big(e)\n ed=cv2.dilate(e,kernel,iterations=2)\n #do_display_image_big(ed)\n ede=ed\n #cv2.erode(ed,kernel,iterations=1)\n #do_display_image_big(ede)\n \n \n\n num_pixels=int(ede.sum()/255)\n\n # only perform on sparse arrays (dull images), not to confuse with nucleoli\n if num_pixels < 8*num_cells*AREAMAX:\n contours, hierarchy = cv2.findContours(ede, cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE)\n #print(f'number of contours: {len(contours)}')\n for i, cont in enumerate(contours):\n area = cv2.contourArea(cont)\n #print(f'area {area}')\n x,y,w,h = cv2.boundingRect(cont)\n #print('x,y,w,h',x,y,w,h)\n if area < AREAMAX and area > AREAMIN and hierarchy[0][i][3] == -1 \\\n and w < MAXWIDTH and w > MINWIDTH \\\n and h > int(w*0.3) and int(h*0.3) < w:\n # found potential nuclear body\n # Check whether it touches the bounding box of a cell\n #print('found contour of potential aggresome size, checking overlap with outrim ...')\n for j, cell_data in enumerate(image_data['cells']):\n cellx, celly, cellw, cellh = cell_data['nuc_bbox']\n #print('checking overlap with ',x,y,w,h,cellx, celly, cellw, cellh )\n if overlaps(x,y,w,h,cellx, celly, cellw, cellh):\n #print('boxes overlap, checking masks')\n nuc_mask_slice=nuc_mask[celly:celly+cellh, cellx:cellx+cellw]\n binary_nuc_mask = nuc_mask_slice.astype('bool')\n #do_display_image_big(binary_outrim_mask)\n\n full_image_mask=np.zeros(pr.shape, dtype='uint8')\n cv2.drawContours(full_image_mask, contours,i,255,cv2.FILLED,cv2.LINE_8)\n #do_display_image_big(full_image_mask)\n contour_cell_slice=full_image_mask[celly:celly+cellh, cellx:cellx+cellw]\n #do_display_image_big(contour_cell_slice)\n contour_cell_binmask=contour_cell_slice.astype(bool)\n #do_display_image_big(contour_cell_binmask)\n show_mask = np.bitwise_xor(binary_nuc_mask,contour_cell_binmask)\n #do_display_image_big(show_mask)\n result_mask = np.logical_and(binary_nuc_mask,contour_cell_binmask)\n overlap_found=result_mask.max()\n if overlap_found:\n cell_data['num_nucspeckles']+=1\n cell_data['nucspeckles_size'] = w","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def analyseNucleusOuterRims(image_data, nuc_cont, contours, nuc_midpoint, nu, pr):\n OUTER_NUCLEUS_RIMWIDTH = 40\n plus10cont=generateNucRims(nuc_cont, nuc_midpoint, nu, OUTER_NUCLEUS_RIMWIDTH)\n\n for j, (nuc_contour, rim_contour, cell_data, cell_contour) in \\\n enumerate(zip(nuc_cont, plus10cont, image_data['cells'], contours)):\n\n full_image_mask=np.zeros(nu.shape, dtype='uint8')\n cv2.drawContours(full_image_mask, plus10cont,j,255,cv2.FILLED,cv2.LINE_8)\n cv2.drawContours(full_image_mask, nuc_cont,j,0,cv2.FILLED,cv2.LINE_8)\n \n #clip off anything that is outside the cell boundary\n cell_boundary_mask=np.zeros(nu.shape, dtype='uint8')\n cv2.drawContours(cell_boundary_mask, contours,j,255,cv2.FILLED,cv2.LINE_8)\n full_image_mask[cell_boundary_mask==0]=0\n \n if False and debug and j < 3:\n do_display_image_big(full_image_mask) \n \n x,y,w,h = cell_data['bbox']\n \n \n #print(f'x,y,w,h:{x},{y},{w},{h}')\n \n #cv2.rectangle(full_image_mask, (x,y), (x+w,y+h), 200,2)\n #display_image_big(full_image_mask) \n cropped_mask = full_image_mask[y:y+h, x:x+w]\n cropped_nu = nu[y:y+h, x:x+w]\n cropped_pr = pr[y:y+h, x:x+w]\n \n\n #display_image_big(cropped_pr)\n #display_image_big(cropped_nu)\n \n minval, maxval, sumval, meanval, im_variance, compmin, compmax,\\\n compsum, compmean, comp_variance, corr, mask_count \\\n = maskStats( cropped_pr, cropped_mask, 255, cropped_nu )\n cell_data['outrim_mask_count']=mask_count\n cell_data['outrim_pr_minval']=minval\n cell_data['outrim_pr_maxval']=maxval\n cell_data['outrim_pr_sumval']=sumval\n cell_data['outrim_pr_meanval']=meanval\n cell_data['outrim_pr_variance']=im_variance \n cell_data['outrim_nu_minval']=compmin\n cell_data['outrim_nu_maxval']=compmax\n cell_data['outrim_nu_sumval']=compsum\n cell_data['outrim_nu_meanval']=compmean\n cell_data['outrim_nu_variance']=comp_variance\n cell_data['outrim_pr_nu_correlation']=corr\n cell_data['outrim_mask']=np.copy(cropped_mask)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def analyseNucleus( image_data, label, im, nu, er, mt, pr, nuc_cont, nuc_midpoint, nuc_mask, nuc_rect): \n\n np.set_printoptions(threshold=sys.maxsize)\n\n prnuc = np.zeros(nu.shape, dtype='uint8')\n \n for j, (cont, cell_data) in enumerate(zip(nuc_cont,image_data['cells'])):\n x,y,w,h = nuc_rect[j]\n mtbc = np.copy(mt) \n mtbc[nuc_mask!=j+1]=0\n erbc = np.copy(er) \n erbc[nuc_mask!=j+1]=0\n nubc = np.copy(nu) \n nubc[nuc_mask!=j+1]=0\n nubcCont = np.copy(nubc) \n nubcContBlur = cv2.GaussianBlur(nubcCont, (5, 5), 0)\n prbc = np.copy(pr) \n prbc[nuc_mask!=j+1]=0\n prbcCont = np.copy(prbc) \n prbcContBlur = cv2.GaussianBlur(prbcCont, (5, 5), 0)\n \n cropped_mt = mtbc[y:y+h, x:x+w]\n cropped_er = erbc[y:y+h, x:x+w]\n cropped_nu = nubc[y:y+h, x:x+w]\n cropped_pr = prbc[y:y+h, x:x+w]\n cropped_mask = nuc_mask[y:y+h, x:x+w]\n\n minval, maxval, sumval, meanval, im_variance,\\\n compmin, compmax, compsum, compmean, comp_variance, corr, mask_count \\\n = maskStats( pr, nuc_mask, j+1, nu )\n\n cell_data['nuc_mask_count']=mask_count\n cell_data['nuc_pr_minval']=minval \n cell_data['nuc_pr_maxval']=maxval \n cell_data['nuc_pr_sumval']=sumval\n cell_data['nuc_pr_meanval']=meanval\n cell_data['nuc_pr_variance']=im_variance\n cell_data['nuc_nu_minval']=compmin\n cell_data['nuc_nu_maxval']=compmax\n cell_data['nuc_nu_sumval']=compsum\n cell_data['nuc_nu_meanval']=compmean\n cell_data['nuc_nu_variance']=comp_variance\n cell_data['nuc_nu_pr_correlation']=corr\n \n \n \n nuc_contour_area = cv2.contourArea(nuc_cont[j])\n cv2.drawContours(nubcContBlur, nuc_cont,j,255,1)\n cropped_nubcCont = nubcContBlur[y:y+h, x:x+w] \n \n mtMax = cropped_mt.max()\n erMax = cropped_er.max()\n nuMax = cropped_nu.max()\n prMax = cropped_pr.max()\n \n numMaskPixels = np.sum(nuc_mask==j+1)\n \n cropped_prsum = np.sum(cropped_pr)\n cropped_nusum = np.sum(cropped_nu)\n if numMaskPixels == 0:\n prMean = 0\n nuMean = 0\n nu_variance = 0\n pr_variance = 0 \n pearson_nupr = 0\n else:\n prMean = cropped_prsum / (numMaskPixels)\n nuMean = cropped_nusum / (numMaskPixels) \n \n nuDiffMean = cropped_nu - nuMean\n nuDiffMean[cropped_mask!=j+1]=0\n prDiffMean = cropped_pr - prMean\n \n\n \n prDiffMean[cropped_mask!=j+1]=0\n \n\n nuprDiffProd=nuDiffMean*prDiffMean\n \n nuprcovariance=np.sum(nuprDiffProd) / (numMaskPixels) \n \n nu_sqdiffs = nuDiffMean**2\n nu_sqdiffs[cropped_mask!=j+1]=0\n nu_variance = np.sum(nu_sqdiffs)/numMaskPixels\n\n pr_sqdiffs = prDiffMean**2\n pr_sqdiffs[cropped_mask!=j+1]=0\n pr_variance = np.sum(pr_sqdiffs)/numMaskPixels\n \n if nu_variance*pr_variance > 0.001:\n pearson_nupr = nuprcovariance / math.sqrt(nu_variance*pr_variance)\n else:\n pearson_nupr = 0\n\n #if j == 1:\n #print(f'matrix shape : {cropped_pr.shape}')\n #print(f'nuMean : {nuMean}')\n #print(f'numMaskPixels : {numMaskPixels}')\n #print(f'nuprcovariance : , {nuprcovariance}')\n #print(f'nu_variance : , {nu_variance}')\n #print(f'pr_variance : , {pr_variance}')\n #print(f'pearson_nupr : {pearson_nupr}')\n\n #print(f'nuDiffMean[20][60] : {nuDiffMean[20][60]}')\n #print(f'prDiffMean[20][60] : {prDiffMean[20][60]}')\n #print(f'nuprDiffProd[20][60] : {nuprDiffProd[20][60]}')\n #print(f'nu_sqdiffs[20][60] : {nu_sqdiffs[20][60]}')\n #print(f'pr_sqdiffs[20][60] : {pr_sqdiffs[20][60]}')\n #print(nuprDiffProd)\n #print(nu_sqdiffs)\n #print(pr_sqdiffs)\n nucleoliPr=np.copy(cropped_pr)\n \n #if label == '14':\n #print(f'{label}, {im}, cell:{j}')\n #print('nucleus')\n #do_display_image_big(cropped_nubcCont)\n #print('protein in nucleus')\n #do_display_image_big(cropped_pr)\n\n \n cropped_nuc_thresh = cv2.adaptiveThreshold(cropped_nubcCont,255,cv2.ADAPTIVE_THRESH_MEAN_C, cv2.THRESH_BINARY,201,0)\n display_image_big(cropped_nuc_thresh) \n contours, hierarchy = cv2.findContours(cropped_nuc_thresh, cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE)\n \n cropped_pr_thresh = cv2.adaptiveThreshold(cropped_pr,255,cv2.ADAPTIVE_THRESH_MEAN_C, cv2.THRESH_BINARY,201,0) \n pr_contours, pr_hierarchy = cv2.findContours(cropped_pr_thresh, cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE)\n\n num_nu_contours = len(contours)\n printinfo(f'num nu_contours: {num_nu_contours}' )\n\n sum_nuc_arclength = 0\n sum_nuc_area = 0\n min_nuc_area = 1000000\n max_nuc_area = 0\n for yy in contours: \n nuc_arclength = cv2.arcLength(yy,True)\n nuc_area = cv2.contourArea(yy)\n if min_nuc_area > nuc_area :\n min_nuc_area = nuc_area\n if max_nuc_area < nuc_area :\n max_nuc_area = nuc_area \n sum_nuc_arclength = sum_nuc_arclength + nuc_arclength\n sum_nuc_area = sum_nuc_area + nuc_area\n if num_nu_contours > 0:\n mean_nuc_arclength = sum_nuc_arclength / num_nu_contours\n mean_nuc_area = sum_nuc_area / num_nu_contours\n if train_or_test == 'train':\n with open('nucstats.txt', 'a') as outf:\n print(f'nucstats: num_nu_contours: {num_nu_contours}, sum_nuc_arclength: {sum_nuc_arclength}, sum_nuc_area: {sum_nuc_area} \\\n mean_nuc_arclength: {mean_nuc_arclength}, mean_nuc_area: {mean_nuc_area} \\\n min_nuc_area: {min_nuc_area}, max_nuc_area: {max_nuc_area} ', file=outf)\n printinfo(f'nucstats: num_nu_contours: {num_nu_contours}, sum_nuc_arclength: {sum_nuc_arclength}, sum_nuc_area: {sum_nuc_area} \\\n mean_nuc_arclength: {mean_nuc_arclength}, mean_nuc_area: {mean_nuc_area} \\\n min_nuc_area: {min_nuc_area}, max_nuc_area: {max_nuc_area} ')\n \n \n num_pr_contours = len(pr_contours)\n printinfo(f'num pr_contours: {num_pr_contours}' )\n\n sum_prc_arclength = 0\n sum_prc_area = 0\n min_prc_area = 1000000\n max_prc_area = 0\n for y in pr_contours: \n prc_arclength = cv2.arcLength(y,True)\n prc_area = cv2.contourArea(y)\n if min_prc_area > prc_area :\n min_prc_area = prc_area\n if max_prc_area < prc_area :\n max_prc_area = prc_area \n sum_prc_arclength = sum_prc_arclength + prc_arclength\n sum_prc_area = sum_prc_area + prc_area\n if num_pr_contours > 0:\n mean_prc_arclength = sum_prc_arclength / num_pr_contours\n mean_prc_area = sum_prc_area / num_pr_contours\n if train_or_test == 'train':\n with open('nucstats.txt', 'a') as outf:\n print(f'prcstats: num_pr_contours: {num_pr_contours}, sum_prc_arclength: {sum_prc_arclength}, sum_prc_area: {sum_prc_area} \\\n mean_prc_arclength: {mean_prc_arclength}, mean_prc_area: {mean_prc_area} \\\n min_prc_area: {min_prc_area}, max_prc_area: {max_prc_area} ', file=outf)\n printinfo(f'prcstats: num_pr_contours: {num_pr_contours}, sum_prc_arclength: {sum_prc_arclength}, sum_prc_area: {sum_prc_area} \\\n mean_prc_arclength: {mean_prc_arclength}, mean_prc_area: {mean_prc_area} \\\n min_prc_area: {min_prc_area}, max_prc_area: {max_prc_area} ')\n \n i=0\n num_contours_kept=0\n nuc_inner_cont=[]\n nuc_inner_midpoint=[]\n\n nucleolus_detected= False\n\n nucleolus_fcent_detected= False\n\n nuc_inner_mask = np.zeros(cropped_nu.shape, dtype='uint8')\n nucleoli_mask = np.zeros(cropped_nu.shape, dtype='uint8')\n \n max_nucleolus_pr_mean = 0\n # max_nucleolus_pr_mean tracks the maximum mean value of the protein within a recognised nucleolus in this cell\n # A correctly recognised nucleolus that has the nucleolus protein marking tends to have a high mean\n # therefore this is a first stab at a confidence level for the nucleoli detection\n \n for x in contours:\n nucleolus_found= False\n nucleolus_fcent_found = False\n \n arclength = cv2.arcLength(contours[i],True)\n area = cv2.contourArea(contours[i])\n\n #print(f'arclength:{arclength} hierarchy {hierarchy[0][i][3]}')\n if ( arclength > 20 and hierarchy[0][i][3] != -1 and area/nuc_contour_area < 0.25 ):\n num_contours_kept+=1\n nucleolus_data = newNucleolusDict(im,label)\n nucleolus_data['cellnum'] = cell_data['cellnum']\n nucleolus_data['arclength'] = arclength\n nucleolus_data['area'] = area\n nucleolus_data['nucleolus_num'] = num_contours_kept \n cell_data['nucleoli'].append(nucleolus_data)\n \n nuc_inner_cont.append(x)\n\n\n cv2.drawContours(nuc_inner_mask, nuc_inner_cont,num_contours_kept-1,num_contours_kept,cv2.FILLED,cv2.LINE_8)\n cv2.drawContours(nucleoli_mask, nuc_inner_cont,num_contours_kept-1,255,cv2.FILLED,cv2.LINE_8)\n x,y,w,h = cv2.boundingRect(contours[i])\n nlbc = np.copy(cropped_nubcCont[y:y+h, x:x+w])\n nlbcpr = np.copy(cropped_pr[y:y+h, x:x+w])\n nlbc_mask = np.copy(nuc_inner_mask[y:y+h, x:x+w])\n nlbc[nlbc_mask!=num_contours_kept]=0\n nlbcpr[nlbc_mask!=num_contours_kept]=0 \n\n #if num_contours_kept == 21:\n # lookclose = True\n #else:\n lookclose = False\n \n minval, maxval, sumval, meanval, im_variance,\\\n compmin, compmax, compsum, compmean, comp_variance, corr, mask_count \\\n = maskStats( nlbcpr, nlbc_mask, num_contours_kept, nlbc, statsdebug = lookclose )\n \n if False and debug and debug_show_nucleolus_masks:\n print(area, meanval, i, lookclose)\n displaycopy = np.copy(cropped_pr)\n cv2.drawContours(displaycopy,nuc_inner_cont,num_contours_kept-1,255,1)\n do_display_image_big(displaycopy)\n \n nucleolus_data['mask_count']=mask_count\n nucleolus_data['pr_minval']=minval\n nucleolus_data['pr_maxval']=maxval\n nucleolus_data['pr_sumval']=sumval\n nucleolus_data['pr_meanval']=meanval\n nucleolus_data['pr_variance']=im_variance \n nucleolus_data['nu_minval']=compmin\n nucleolus_data['nu_maxval']=compmax\n nucleolus_data['nu_sumval']=compsum\n nucleolus_data['nu_meanval']=compmean\n nucleolus_data['nu_variance']=comp_variance\n nucleolus_data['pr_nu_correlation']=corr \n nucleolus_data['bbox'] = (x,y,w,h) \n\n i=i+1\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def showHistogram(bgrImage):\n bgr_planes = cv2.split(bgrImage)\n histSize = 256\n histRange = (0, 256) # the upper boundary is exclusive\n accumulate = False\n b_hist = cv2.calcHist(bgr_planes, [0], None, [histSize], histRange, accumulate=accumulate)\n g_hist = cv2.calcHist(bgr_planes, [1], None, [histSize], histRange, accumulate=accumulate)\n r_hist = cv2.calcHist(bgr_planes, [2], None, [histSize], histRange, accumulate=accumulate)\n hist_w = 512\n hist_h = 400\n bin_w = int(round( hist_w/histSize ))\n histImage = np.zeros((hist_h, hist_w, 3), dtype=np.uint8)\n cv2.normalize(b_hist, b_hist, alpha=0, beta=hist_h, norm_type=cv2.NORM_MINMAX)\n cv2.normalize(g_hist, g_hist, alpha=0, beta=hist_h, norm_type=cv2.NORM_MINMAX)\n cv2.normalize(r_hist, r_hist, alpha=0, beta=hist_h, norm_type=cv2.NORM_MINMAX)\n for i in range(1, histSize):\n cv2.line(histImage, ( bin_w*(i-1), hist_h - int(b_hist[i-1]) ),\n ( bin_w*(i), hist_h - int(b_hist[i]) ),\n ( 255, 0, 0), thickness=2)\n cv2.line(histImage, ( bin_w*(i-1), hist_h - int(g_hist[i-1]) ),\n ( bin_w*(i), hist_h - int(g_hist[i]) ),\n ( 0, 255, 0), thickness=2)\n cv2.line(histImage, ( bin_w*(i-1), hist_h - int(r_hist[i-1]) ),\n ( bin_w*(i), hist_h - int(r_hist[i]) ),\n ( 0, 0, 255), thickness=2)\n do_display_image_big(histImage)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def printhist(im):\n hist_im = np.histogram(im, bins=[0,31,63,95,127,159,191,223,255])\n hist256=np.histogram(im, bins=256)\n maxbuck=0\n maxfreq=0\n for i in range(256):\n if hist256[0][i]>maxfreq:\n maxfreq=hist256[0][i]\n maxbuck=i\n print(hist_im)\n print(f'most frequent value is {maxbuck}')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from skimage import feature\nimport numpy as np\nclass LocalBinaryPatterns:\n def __init__(self, numPoints, radius):\n # store the number of points and radius\n self.numPoints = numPoints\n self.radius = radius\n def describe(self, image, eps=1e-7):\n # compute the Local Binary Pattern representation\n # of the image, and then use the LBP representation\n # to build the histogram of patterns\n lbp = feature.local_binary_pattern(image, self.numPoints,\n self.radius, method=\"uniform\")\n (hist, _) = np.histogram(lbp.ravel(),\n bins=np.arange(0, self.numPoints + 3),\n range=(0, self.numPoints + 2))\n # normalize the histogram\n hist = hist.astype(\"float\")\n hist /= (hist.sum() + eps)\n # return the histogram of Local Binary Patterns\n return hist\n ","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.svm import LinearSVC\n\ndef addLbpData (desc, image_data, box_sample_list, pr, data, labels ):\n boxwidth=40\n label=image_data['train_label']\n max_boxes_per_image=3\n \n nn0=int(len(box_sample_list)/3)\n nn=nn0\n for b in box_sample_list:\n nn=nn+1\n if nn > nn0+max_boxes_per_image:\n break\n bllx=b[0]\n blly=b[1]\n image_slice = np.copy(pr[blly:blly+boxwidth, bllx:bllx+boxwidth])\n hist = desc.describe(image_slice)\n data.append(hist)\n print(f'adding with label {label}')\n do_display_image_big(image_slice)\n labels.append(label)\n ","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def predictBoxlist(desc, model, image_data, box_sample_list, pr ):\n boxwidth=40\n label=image_data['train_label']\n \n preds=[]\n \n for b in box_sample_list:\n bllx=b[0]\n blly=b[1]\n image_slice = np.copy(pr[blly:blly+boxwidth, bllx:bllx+boxwidth]) \n hist = desc.describe(image_slice)\n prediction = model.predict(hist.reshape(1, -1))\n #print(f'prediction {prediction[0]} label {label}')\n found=False\n\n for pp, pred in enumerate( preds ):\n if prediction[0]==pred[0]:\n preds[pp]=(prediction[0],pred[1]+1)\n found=True\n if not found:\n preds.append((prediction[0],1))\n #printinfo('List result:')\n #printinfo(preds)\n nummaxpred=0\n maxpred='18'\n for pred in preds:\n if pred[1]>nummaxpred:\n maxpred = pred[0]\n nummaxpred=pred[1]\n image_data['141617-pred']=maxpred\n #printinfo(f'train label is {label}, Result is {maxpred}')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def compassSampleBoxes(image_data, nuc_cont, contours, nuc_mask, nuc_midpoint, image):\n contour_box_ll_list=[]\n rimwidth=70\n boxwidth=40\n drawcopy = np.copy(image)\n\n cv2.drawContours(drawcopy, nuc_cont,-1,255,1)\n cv2.drawContours(drawcopy, contours,-1,255,1)\n \n for j ,(cont, mid, cellcont) in enumerate(zip(nuc_cont,nuc_midpoint,contours)):\n\n #print(f'Contour {j}')\n #print(f'midpoint {mid}')\n test_mask=np.zeros((image.shape[0],image.shape[1]), np.uint8)\n cv2.drawContours(test_mask,contours,j,255,cv2.FILLED,cv2.LINE_8)\n cv2.drawContours(test_mask,nuc_cont,j,0,cv2.FILLED,cv2.LINE_8)\n\n imagesize=np.array(image.shape)\n \n for k in range(len(cont)):\n #print(f'Contour {j} Segment {k}')\n vertex= nuc_cont[j][k][0]\n #print('vertex', vertex)\n rayvec = vertex - mid\n #print('rayvec', rayvec)\n #print(type(mid))\n scalefac=1/abs(rayvec).sum()\n rayvecnew=(rayvec*(abs(rayvec).sum()+rimwidth)*scalefac).astype(int)\n vertexnew=rayvecnew+mid\n clipvertex_low = np.fmax(vertexnew,[0,0]) \n clipvertex = np.fmin(clipvertex_low,imagesize) \n \n min_point_x=min(vertex[0],clipvertex[0])\n min_point_y=min(vertex[1],clipvertex[1])\n \n #if j == 4 : \n #if (test_box_slice.sum()/255 > boxwidth*boxwidth*0.8 ) : \n # cv2.rectangle(test_mask,(min_point_x, min_point_y),(min_point_x+boxwidth, min_point_y+boxwidth),125, 2 )\n #else: \n # cv2.rectangle(test_mask,(min_point_x, min_point_y),(min_point_x+boxwidth, min_point_y+boxwidth),20, 2 )\n \n \n test_box_slice=np.copy(test_mask[min_point_y:min_point_y+boxwidth, min_point_x:min_point_x+boxwidth]) \n tbss=test_box_slice.sum()/255\n\n if ( tbss > boxwidth*boxwidth*0.9 and \\\n min_point_x+boxwidth < imagesize[0] and min_point_y+boxwidth < imagesize[1]):\n image_slice = np.copy(image[min_point_y:min_point_y+boxwidth, min_point_x:min_point_x+boxwidth])\n\n info_content=image_slice.sum()/255\n if (info_content > 20):# remove almost-empty boxes\n contour_box_ll_list.append((min_point_x, min_point_y))\n #cv2.rectangle(drawcopy,(min_point_x, min_point_y),(min_point_x+boxwidth, min_point_y+boxwidth), 255, 2 )\n #cv2.circle(drawcopy, (min_point_x, min_point_y), 5, 255, 2)\n\n #do_display_image_big(drawcopy)\n\n return contour_box_ll_list","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"## percentpercentwritefile faster_hpa_cell_segment.py\n\nimport numpy as np\nimport pandas as pd\nimport os\nfrom tqdm import tqdm\n\nfrom hpacellseg.cellsegmentator import *\n\n\nclass CellSegmentator(object):\n \"\"\"Uses pretrained DPN-Unet models to segment cells from images.\"\"\"\n\n def __init__(\n self,\n nuclei_model=\"../input/hpacellsegmentatormodelweights/dpn_unet_nuclei_v1.pth\",\n cell_model=\"../input/hpacellsegmentatormodelweights/dpn_unet_cell_3ch_v1.pth\",\n scale_factor=1.0,\n device=\"cuda\",\n padding=False,\n multi_channel_model=True,\n ):\n \"\"\"Class for segmenting nuclei and whole cells from confocal microscopy images.\n It takes lists of images and returns the raw output from the\n specified segmentation model. Models can be automatically\n downloaded if they are not already available on the system.\n When working with images from the Huan Protein Cell atlas, the\n outputs from this class' methods are well combined with the\n label functions in the utils module.\n Note that for cell segmentation, there are two possible models\n available. One that works with 2 channeled images and one that\n takes 3 channels.\n Keyword arguments:\n nuclei_model -- A loaded torch nuclei segmentation model or the\n path to a file which contains such a model.\n If the argument is a path that points to a non-existant file,\n a pretrained nuclei_model is going to get downloaded to the\n specified path (default: './nuclei_model.pth').\n cell_model -- A loaded torch cell segmentation model or the\n path to a file which contains such a model.\n The cell_model argument can be None if only nuclei\n are to be segmented (default: './cell_model.pth').\n scale_factor -- How much to scale images before they are fed to\n segmentation models. Segmentations will be scaled back\n up by 1/scale_factor to match the original image\n (default: 0.25).\n device -- The device on which to run the models.\n This should either be 'cpu' or 'cuda' or pointed cuda\n device like 'cuda:0' (default: 'cuda').\n padding -- Whether to add padding to the images before feeding the\n images to the network. (default: False).\n multi_channel_model -- Control whether to use the 3-channel cell model or not.\n If True, use the 3-channel model, otherwise use the\n 2-channel version (default: True).\n \"\"\"\n if device != \"cuda\" and device != \"cpu\" and \"cuda\" not in device:\n raise ValueError(f\"{device} is not a valid device (cuda/cpu)\")\n if device != \"cpu\":\n try:\n assert torch.cuda.is_available()\n except AssertionError:\n print(\"No GPU found, using CPU.\", file=sys.stderr)\n device = \"cpu\"\n self.device = device\n\n if isinstance(nuclei_model, str):\n if not os.path.exists(nuclei_model):\n print(\n f\"Could not find {nuclei_model}. Downloading it now\",\n file=sys.stderr,\n )\n download_with_url(NUCLEI_MODEL_URL, nuclei_model)\n nuclei_model = torch.load(\n nuclei_model, map_location=torch.device(self.device)\n )\n if isinstance(nuclei_model, torch.nn.DataParallel) and device == \"cpu\":\n nuclei_model = nuclei_model.module\n\n self.nuclei_model = nuclei_model.to(self.device).eval()\n\n self.multi_channel_model = multi_channel_model\n if isinstance(cell_model, str):\n if not os.path.exists(cell_model):\n print(\n f\"Could not find {cell_model}. Downloading it now\", file=sys.stderr\n )\n if self.multi_channel_model:\n download_with_url(MULTI_CHANNEL_CELL_MODEL_URL, cell_model)\n else:\n download_with_url(TWO_CHANNEL_CELL_MODEL_URL, cell_model)\n cell_model = torch.load(cell_model, map_location=torch.device(self.device))\n self.cell_model = cell_model.to(self.device).eval()\n self.scale_factor = scale_factor\n self.padding = padding\n\n def _image_conversion(self, images):\n \"\"\"Convert/Format images to RGB image arrays list for cell predictions.\n Intended for internal use only.\n Keyword arguments:\n images -- list of lists of image paths/arrays. It should following the\n pattern if with er channel input,\n [\n [microtubule_path0/image_array0, microtubule_path1/image_array1, ...],\n [er_path0/image_array0, er_path1/image_array1, ...],\n [nuclei_path0/image_array0, nuclei_path1/image_array1, ...]\n ]\n or if without er input,\n [\n [microtubule_path0/image_array0, microtubule_path1/image_array1, ...],\n None,\n [nuclei_path0/image_array0, nuclei_path1/image_array1, ...]\n ]\n \"\"\"\n microtubule_imgs, er_imgs, nuclei_imgs = images\n if self.multi_channel_model:\n if not isinstance(er_imgs, list):\n raise ValueError(\"Please speicify the image path(s) for er channels!\")\n else:\n if not er_imgs is None:\n raise ValueError(\n \"second channel should be None for two channel model predition!\"\n )\n\n if not isinstance(microtubule_imgs, list):\n raise ValueError(\"The microtubule images should be a list\")\n if not isinstance(nuclei_imgs, list):\n raise ValueError(\"The microtubule images should be a list\")\n\n if er_imgs:\n if not len(microtubule_imgs) == len(er_imgs) == len(nuclei_imgs):\n raise ValueError(\"The lists of images needs to be the same length\")\n else:\n if not len(microtubule_imgs) == len(nuclei_imgs):\n raise ValueError(\"The lists of images needs to be the same length\")\n\n if not all(isinstance(item, np.ndarray) for item in microtubule_imgs):\n microtubule_imgs = [\n os.path.expanduser(item) for _, item in enumerate(microtubule_imgs)\n ]\n nuclei_imgs = [\n os.path.expanduser(item) for _, item in enumerate(nuclei_imgs)\n ]\n\n microtubule_imgs = list(\n map(lambda item: imageio.imread(item), microtubule_imgs)\n )\n nuclei_imgs = list(map(lambda item: imageio.imread(item), nuclei_imgs))\n if er_imgs:\n er_imgs = [os.path.expanduser(item) for _, item in enumerate(er_imgs)]\n er_imgs = list(map(lambda item: imageio.imread(item), er_imgs))\n\n if not er_imgs:\n er_imgs = [\n np.zeros(item.shape, dtype=item.dtype)\n for _, item in enumerate(microtubule_imgs)\n ]\n cell_imgs = list(\n map(\n lambda item: np.dstack((item[0], item[1], item[2])),\n list(zip(microtubule_imgs, er_imgs, nuclei_imgs)),\n )\n )\n\n return cell_imgs\n\n def pred_nuclei(self, images):\n \"\"\"Predict the nuclei segmentation.\n Keyword arguments:\n images -- A list of image arrays or a list of paths to images.\n If as a list of image arrays, the images could be 2d images\n of nuclei data array only, or must have the nuclei data in\n the blue channel; If as a list of file paths, the images\n could be RGB image files or gray scale nuclei image file\n paths.\n Returns:\n predictions -- A list of predictions of nuclei segmentation for each nuclei image.\n \"\"\"\n\n def _preprocess(image):\n if isinstance(image, str):\n image = imageio.imread(image)\n self.target_shape = image.shape\n if len(image.shape) == 2:\n image = np.dstack((image, image, image))\n image = transform.rescale(image, self.scale_factor, multichannel=True)\n nuc_image = np.dstack((image[..., 2], image[..., 2], image[..., 2]))\n if self.padding:\n rows, cols = nuc_image.shape[:2]\n self.scaled_shape = rows, cols\n nuc_image = cv2.copyMakeBorder(\n nuc_image,\n 32,\n (32 - rows % 32),\n 32,\n (32 - cols % 32),\n cv2.BORDER_REFLECT,\n )\n nuc_image = nuc_image.transpose([2, 0, 1])\n return nuc_image\n\n def _segment_helper(imgs):\n with torch.no_grad():\n mean = torch.as_tensor(NORMALIZE[\"mean\"], device=self.device)\n std = torch.as_tensor(NORMALIZE[\"std\"], device=self.device)\n imgs = torch.tensor(imgs).float()\n imgs = imgs.to(self.device)\n imgs = imgs.sub_(mean[:, None, None]).div_(std[:, None, None])\n\n imgs = self.nuclei_model(imgs)\n imgs = F.softmax(imgs, dim=1)\n return imgs\n\n preprocessed_imgs = list(map(_preprocess, images))\n bs = 24\n predictions = []\n for i in range(0, len(preprocessed_imgs), bs):\n start = i\n end = min(len(preprocessed_imgs), i+bs)\n x = preprocessed_imgs[start:end]\n pred = _segment_helper(x).cpu().numpy()\n predictions.append(pred)\n predictions = list(np.concatenate(predictions, axis=0))\n predictions = map(util.img_as_ubyte, predictions)\n predictions = list(map(self._restore_scaling_padding, predictions))\n return predictions\n\n def _restore_scaling_padding(self, n_prediction):\n \"\"\"Restore an image from scaling and padding.\n This method is intended for internal use.\n It takes the output from the nuclei model as input.\n \"\"\"\n n_prediction = n_prediction.transpose([1, 2, 0])\n if self.padding:\n n_prediction = n_prediction[\n 32 : 32 + self.scaled_shape[0], 32 : 32 + self.scaled_shape[1], ...\n ]\n if not self.scale_factor == 1:\n n_prediction[..., 0] = 0\n n_prediction = cv2.resize(\n n_prediction,\n (self.target_shape[0], self.target_shape[1]),\n interpolation=cv2.INTER_AREA,\n )\n return n_prediction\n\n def pred_cells(self, images, precombined=False):\n \"\"\"Predict the cell segmentation for a list of images.\n Keyword arguments:\n images -- list of lists of image paths/arrays. It should following the\n pattern if with er channel input,\n [\n [microtubule_path0/image_array0, microtubule_path1/image_array1, ...],\n [er_path0/image_array0, er_path1/image_array1, ...],\n [nuclei_path0/image_array0, nuclei_path1/image_array1, ...]\n ]\n or if without er input,\n [\n [microtubule_path0/image_array0, microtubule_path1/image_array1, ...],\n None,\n [nuclei_path0/image_array0, nuclei_path1/image_array1, ...]\n ]\n The ER channel is required when multichannel is True\n and required to be None when multichannel is False.\n The images needs to be of the same size.\n precombined -- If precombined is True, the list of images is instead supposed to be\n a list of RGB numpy arrays (default: False).\n Returns:\n predictions -- a list of predictions of cell segmentations.\n \"\"\"\n\n def _preprocess(image):\n self.target_shape = image.shape\n if not len(image.shape) == 3:\n raise ValueError(\"image should has 3 channels\")\n cell_image = transform.rescale(image, self.scale_factor, multichannel=True)\n if self.padding:\n rows, cols = cell_image.shape[:2]\n self.scaled_shape = rows, cols\n cell_image = cv2.copyMakeBorder(\n cell_image,\n 32,\n (32 - rows % 32),\n 32,\n (32 - cols % 32),\n cv2.BORDER_REFLECT,\n )\n cell_image = cell_image.transpose([2, 0, 1])\n return cell_image\n\n def _segment_helper(imgs):\n with torch.no_grad():\n mean = torch.as_tensor(NORMALIZE[\"mean\"], device=self.device)\n std = torch.as_tensor(NORMALIZE[\"std\"], device=self.device)\n imgs = torch.tensor(imgs).float()\n imgs = imgs.to(self.device)\n imgs = imgs.sub_(mean[:, None, None]).div_(std[:, None, None])\n\n imgs = self.cell_model(imgs)\n imgs = F.softmax(imgs, dim=1)\n return imgs\n\n if not precombined:\n images = self._image_conversion(images)\n preprocessed_imgs = list(map(_preprocess, images))\n bs = 24\n predictions = []\n for i in range(0, len(preprocessed_imgs), bs):\n start = i\n end = min(len(preprocessed_imgs), i+bs)\n x = preprocessed_imgs[start:end]\n pred = _segment_helper(x).cpu().numpy()\n predictions.append(pred)\n predictions = list(np.concatenate(predictions, axis=0))\n predictions = map(self._restore_scaling_padding, predictions)\n predictions = list(map(util.img_as_ubyte, predictions))\n\n return predictions\n \n \n\n\n \n \n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import os.path\nimport urllib\nimport zipfile\n\nimport numpy as np\nimport scipy.ndimage as ndi\nfrom skimage import filters, measure, segmentation\nfrom skimage.morphology import (binary_erosion, closing, disk,\n remove_small_holes, remove_small_objects)\n\nHIGH_THRESHOLD = 0.4\nLOW_THRESHOLD = HIGH_THRESHOLD - 0.25\n\n\ndef download_with_url(url_string, file_path, unzip=False):\n \"\"\"Download file with a link.\"\"\"\n with urllib.request.urlopen(url_string) as response, open(\n file_path, \"wb\"\n ) as out_file:\n data = response.read() # a `bytes` object\n out_file.write(data)\n\n if unzip:\n with zipfile.ZipFile(file_path, \"r\") as zip_ref:\n zip_ref.extractall(os.path.dirname(file_path))\n\n\ndef __fill_holes(image):\n \"\"\"Fill_holes for labelled image, with a unique number.\"\"\"\n boundaries = segmentation.find_boundaries(image)\n image = np.multiply(image, np.invert(boundaries))\n image = ndi.binary_fill_holes(image > 0)\n image = ndi.label(image)[0]\n return image\n\n\n\n\n\ndef label_cell(nuclei_pred, cell_pred):\n \"\"\"Label the cells and the nuclei.\n Keyword arguments:\n nuclei_pred -- a 3D numpy array of a prediction from a nuclei image.\n cell_pred -- a 3D numpy array of a prediction from a cell image.\n Returns:\n A tuple containing:\n nuclei-label -- A nuclei mask data array.\n cell-label -- A cell mask data array.\n 0's in the data arrays indicate background while a continous\n strech of a specific number indicates the area for a specific\n cell.\n The same value in cell mask and nuclei mask refers to the identical cell.\n NOTE: The nuclei labeling from this function will be sligthly\n different from the values in :func:`label_nuclei` as this version\n will use information from the cell-predictions to make better\n estimates.\n \"\"\"\n def __wsh(\n mask_img,\n threshold,\n border_img,\n seeds,\n threshold_adjustment=0.35,\n small_object_size_cutoff=10,\n ):\n img_copy = np.copy(mask_img)\n m = seeds * border_img # * dt\n img_copy[m <= threshold + threshold_adjustment] = 0\n img_copy[m > threshold + threshold_adjustment] = 1\n img_copy = img_copy.astype(np.bool)\n img_copy = remove_small_objects(img_copy, small_object_size_cutoff).astype(\n np.uint8\n )\n\n mask_img[mask_img <= threshold] = 0\n mask_img[mask_img > threshold] = 1\n mask_img = mask_img.astype(np.bool)\n mask_img = remove_small_holes(mask_img, 63)\n mask_img = remove_small_objects(mask_img, 1).astype(np.uint8)\n markers = ndi.label(img_copy, output=np.uint32)[0]\n labeled_array = segmentation.watershed(\n mask_img, markers, mask=mask_img, watershed_line=True\n )\n return labeled_array\n\n nuclei_label = __wsh(\n nuclei_pred[..., 2] / 255.0,\n 0.4,\n 1 - (nuclei_pred[..., 1] + cell_pred[..., 1]) / 255.0 > 0.05,\n nuclei_pred[..., 2] / 255,\n threshold_adjustment=-0.25,\n small_object_size_cutoff=32,\n )\n\n # for hpa_image, to remove the small pseduo nuclei\n nuclei_label = remove_small_objects(nuclei_label, 157)\n nuclei_label = measure.label(nuclei_label)\n # this is to remove the cell borders' signal from cell mask.\n # could use np.logical_and with some revision, to replace this func.\n # Tuned for segmentation hpa images\n threshold_value = max(0.22, filters.threshold_otsu(cell_pred[..., 2] / 255) * 0.5)\n # exclude the green area first\n cell_region = np.multiply(\n cell_pred[..., 2] / 255 > threshold_value,\n np.invert(np.asarray(cell_pred[..., 1] / 255 > 0.05, dtype=np.int8)),\n )\n sk = np.asarray(cell_region, dtype=np.int8)\n distance = np.clip(cell_pred[..., 2], 255 * threshold_value, cell_pred[..., 2])\n cell_label = segmentation.watershed(-distance, nuclei_label, mask=sk)\n cell_label = remove_small_objects(cell_label, 344).astype(np.uint8)\n selem = disk(2)\n cell_label = closing(cell_label, selem)\n cell_label = __fill_holes(cell_label)\n # this part is to use green channel, and extend cell label to green channel\n # benefit is to exclude cells clear on border but without nucleus\n sk = np.asarray(\n np.add(\n np.asarray(cell_label > 0, dtype=np.int8),\n np.asarray(cell_pred[..., 1] / 255 > 0.05, dtype=np.int8),\n )\n > 0,\n dtype=np.int8,\n )\n cell_label = segmentation.watershed(-distance, cell_label, mask=sk)\n cell_label = __fill_holes(cell_label)\n cell_label = np.asarray(cell_label > 0, dtype=np.uint8)\n cell_label = measure.label(cell_label)\n cell_label = remove_small_objects(cell_label, 344)\n cell_label = measure.label(cell_label)\n cell_label = np.asarray(cell_label, dtype=np.uint16)\n nuclei_label = np.multiply(cell_label > 0, nuclei_label) > 0\n nuclei_label = measure.label(nuclei_label)\n nuclei_label = remove_small_objects(nuclei_label, 157)\n nuclei_label = np.multiply(cell_label, nuclei_label > 0)\n\n return nuclei_label, cell_label\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def newImageDict(image_id, trainlabel, width, height ):\n ann = {\n 'image_id': image_id,\n 'train_label': trainlabel,\n 'width': width,\n 'height': height,\n 'cells': [],\n 'pr_minval': None, \n 'pr_maxval': None, \n 'pr_sumval': None,\n 'pr_medianval': None,\n 'pr_maxfreq': None,\n 'pr_mostfreq': None,\n 'mt_minval': None, \n 'mt_maxval': None, \n 'mt_sumval': None,\n 'mt_medianval': None,\n 'mt_maxfreq': None,\n 'mt_mostfreq': None,\n 'nu_minval': None, \n 'nu_maxval': None, \n 'nu_sumval': None,\n 'nu_medianval': None,\n 'nu_maxfreq': None,\n 'nu_mostfreq': None,\n 'er_minval': None, \n 'er_maxval': None, \n 'er_sumval': None,\n 'er_medianval': None,\n 'er_maxfreq': None,\n 'er_mostfreq': None,\n 'min_pr_er_correlation': None,\n 'max_pr_er_correlation' : None, \n 'mean_pr_er_correlation': None,\n 'min_pr_mt_correlation': None,\n 'max_pr_mt_correlation' : None, \n 'mean_pr_mt_correlation': None,\n 'min_pr_mt_correlation_fc': None,\n 'max_pr_mt_correlation_fc': None, \n 'mean_pr_mt_correlation_fc': None\n }\n return ann","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class CellSegmentatorOpenCV(object):\n \"\"\"Uses openCV functionality to identify nucleus and cell masks\"\"\"\n\n def __init__(\n self,\n nuclei_model=\"../input/hpacellsegmentatormodelweights/dpn_unet_nuclei_v1.pth\",\n cell_model=\"../input/hpacellsegmentatormodelweights/dpn_unet_cell_3ch_v1.pth\",\n scale_factor=1.0,\n device=\"cuda\",\n padding=False,\n multi_channel_model=True,\n ):\n self.scale_factor = scale_factor\n self.padding = padding\n \n def _image_conversion(self, images):\n \"\"\"Convert/Format images to RGB image arrays list for cell predictions.\n Intended for internal use only.\n Keyword arguments:\n images -- list of lists of image paths/arrays. It should following the\n pattern if with er channel input,\n [\n [microtubule_path0/image_array0, microtubule_path1/image_array1, ...],\n [er_path0/image_array0, er_path1/image_array1, ...],\n [nuclei_path0/image_array0, nuclei_path1/image_array1, ...]\n ]\n or if without er input,\n [\n [microtubule_path0/image_array0, microtubule_path1/image_array1, ...],\n None,\n [nuclei_path0/image_array0, nuclei_path1/image_array1, ...]\n ]\n \"\"\"\n microtubule_imgs, er_imgs, nuclei_imgs = images\n if self.multi_channel_model:\n if not isinstance(er_imgs, list):\n raise ValueError(\"Please speicify the image path(s) for er channels!\")\n else:\n if not er_imgs is None:\n raise ValueError(\n \"second channel should be None for two channel model predition!\"\n )\n\n if not isinstance(microtubule_imgs, list):\n raise ValueError(\"The microtubule images should be a list\")\n if not isinstance(nuclei_imgs, list):\n raise ValueError(\"The microtubule images should be a list\")\n\n if er_imgs:\n if not len(microtubule_imgs) == len(er_imgs) == len(nuclei_imgs):\n raise ValueError(\"The lists of images needs to be the same length\")\n else:\n if not len(microtubule_imgs) == len(nuclei_imgs):\n raise ValueError(\"The lists of images needs to be the same length\")\n\n if not all(isinstance(item, np.ndarray) for item in microtubule_imgs):\n microtubule_imgs = [\n os.path.expanduser(item) for _, item in enumerate(microtubule_imgs)\n ]\n nuclei_imgs = [\n os.path.expanduser(item) for _, item in enumerate(nuclei_imgs)\n ]\n\n microtubule_imgs = list(\n map(lambda item: imageio.imread(item), microtubule_imgs)\n )\n nuclei_imgs = list(map(lambda item: imageio.imread(item), nuclei_imgs))\n if er_imgs:\n er_imgs = [os.path.expanduser(item) for _, item in enumerate(er_imgs)]\n er_imgs = list(map(lambda item: imageio.imread(item), er_imgs))\n\n if not er_imgs:\n er_imgs = [\n np.zeros(item.shape, dtype=item.dtype)\n for _, item in enumerate(microtubule_imgs)\n ]\n cell_imgs = list(\n map(\n lambda item: np.dstack((item[0], item[1], item[2])),\n list(zip(microtubule_imgs, er_imgs, nuclei_imgs)),\n )\n )\n\n return cell_imgs\n \n\n\n \n \n def pred_nuclei(self, images):\n \"\"\"Predict the nuclei segmentation.\n Keyword arguments:\n images -- A list of image arrays or a list of paths to images.\n If as a list of image arrays, the images could be 2d images\n of nuclei data array only, or must have the nuclei data in\n the blue channel; If as a list of file paths, the images\n could be RGB image files or gray scale nuclei image file\n paths.\n Returns:\n predictions -- A list of predictions of nuclei segmentation for each nuclei image.\n \"\"\"\n \n def _preprocess(image):\n \n if isinstance(image, str):\n image = imageio.imread(image)\n\n if not self.scale_factor == 1:\n nuc_image = transform.rescale(image, self.scale_factor, multichannel=True)\n else:\n nuc_image = image\n\n if self.padding:\n rows, cols = nuc_image.shape\n self.scaled_shape = rows, cols\n nuc_image = cv2.copyMakeBorder(\n nuc_image,\n 32,\n (32 - rows % 32),\n 32,\n (32 - cols % 32),\n cv2.BORDER_REFLECT,\n )\n return nuc_image\n\n preprocessed_imgs = list(map(_preprocess, images))\n\n predictions = []\n for i in range(len(preprocessed_imgs)):\n nuc_cont, nuc_midpoint, nuc_mask, nuc_rect = find_nuc_contours (preprocessed_imgs[i])\n predictions.append(nuc_mask) \n return predictions\n\n\n\n def pred_cells(self, images, prot, image_data_list, precombined=False):\n \"\"\"Predict the cell segmentation for a list of images.\n Keyword arguments:\n images -- list of lists of image paths/arrays. It should following the\n pattern if with er channel input,\n [\n [microtubule_path0/image_array0, microtubule_path1/image_array1, ...],\n [er_path0/image_array0, er_path1/image_array1, ...],\n [nuclei_path0/image_array0, nuclei_path1/image_array1, ...]\n ]\n or if without er input,\n [\n [microtubule_path0/image_array0, microtubule_path1/image_array1, ...],\n None,\n [nuclei_path0/image_array0, nuclei_path1/image_array1, ...]\n ]\n The ER channel is required when multichannel is True\n and required to be None when multichannel is False.\n The images needs to be of the same size.\n precombined -- If precombined is True, the list of images is instead supposed to be\n a list of RGB numpy arrays (default: False).\n Returns:\n predictions -- a list of predictions of cell segmentations.\n \"\"\"\n\n def _preprocess(image):\n self.target_shape = image.shape\n if not len(image.shape) == 3:\n raise ValueError(\"image should have 3 channels\")\n if not self.scale_factor == 1:\n cell_image = transform.rescale(image, self.scale_factor, multichannel=True)\n else:\n cell_image = image\n\n if self.padding:\n rows, cols = cell_image.shape[:2]\n self.scaled_shape = rows, cols\n cell_image = cv2.copyMakeBorder(\n cell_image,\n 32,\n (32 - rows % 32),\n 32,\n (32 - cols % 32),\n cv2.BORDER_REFLECT,\n )\n cell_image = cell_image.transpose([2, 0, 1])\n return cell_image\n\n if not precombined:\n images = self._image_conversion(images)\n preprocessed_imgs = list(map(_preprocess, images))\n predictions=[]\n \n desc = LocalBinaryPatterns(24, 8)\n if train_or_test == 'train':\n # make binary_patterns descriptor\n\n model = LinearSVC(C=100.0, random_state=42, max_iter=10000)\n model_data_list=[]\n model_label_list=[]\n else:\n ROOT = '../input/hpa-single-cell-image-classification/'\n model = load('../input/141617model/celltexture.joblib')\n \n #for i in range(len(preprocessed_imgs)):\n\n for i in tqdm(range(len(preprocessed_imgs))):\n mt,er,nu = preprocessed_imgs[i]\n pr=prot[i]\n image_data = image_data_list[i]\n im =image_data['image_id']\n label=image_data['train_label']\n \n \n maxval, minval, sumval, maxfreq, mostfreq, medianval = imageLevelStats(nu, image_data)\n image_data['nu_maxval']=maxval\n image_data['nu_minval']=minval\n image_data['nu_sumval']=sumval\n image_data['nu_maxfreq']=maxfreq\n image_data['nu_mostfreq']=mostfreq\n image_data['nu_medianval']=medianval\n \n maxval, minval, sumval, maxfreq, mostfreq, medianval = imageLevelStats(pr, image_data)\n image_data['pr_maxval']=maxval\n image_data['pr_minval']=minval\n image_data['pr_sumval']=sumval\n image_data['pr_maxfreq']=maxfreq\n image_data['pr_mostfreq']=mostfreq\n image_data['pr_medianval']=medianval\n \n maxval, minval, sumval, maxfreq, mostfreq, medianval = imageLevelStats(er, image_data)\n image_data['er_maxval']=maxval\n image_data['er_minval']=minval\n image_data['er_sumval']=sumval\n image_data['er_maxfreq']=maxfreq\n image_data['er_mostfreq']=mostfreq\n image_data['er_medianval']=medianval \n \n maxval, minval, sumval, maxfreq, mostfreq, medianval = imageLevelStats(mt, image_data)\n image_data['mt_maxval']=maxval\n image_data['mt_minval']=minval\n image_data['mt_sumval']=sumval\n image_data['mt_maxfreq']=maxfreq\n image_data['mt_mostfreq']=mostfreq\n image_data['mt_medianval']=medianval \n \n\n nuc_cont, nuc_midpoint, nuc_mask, nuc_rect = find_nuc_contours (nu)\n \n cell_graph = generateLayoutGraph(image_data, nuc_cont, nuc_midpoint, nuc_rect)\n contours = find_cells( mt, er, nu, nuc_cont, nuc_mask,nuc_midpoint,cell_graph)\n\n box_sample_list = compassSampleBoxes(image_data, nuc_cont, contours, nuc_mask, nuc_midpoint, pr)\n \n if train_or_test == 'train':\n # make binary_patterns descriptor\n addLbpData( desc, image_data, box_sample_list, pr, model_data_list, model_label_list)\n else:\n predictBoxlist(desc, model, image_data, box_sample_list, pr )\n \n #print('cell nodes', cell_graph.nodes)\n #print('cell edges', cell_graph.edges)\n #print(f'len(nuc_cont),{len(nuc_cont)}')\n for j in range(len(contours)):\n cont = contours[j]\n nucnt = nuc_cont[j]\n x,y,w,h = cv2.boundingRect(cont)\n cell = newCellDict(im, label, x, y, w, h)\n cell['nuc_bbox'] = cv2.boundingRect(nucnt)\n image_data['cells'].append(cell)\n cell['cellnum']=j\n cell['nuc_midpoint']=nuc_midpoint[j]\n\n\n analyseMicrotubules(image_data,contours,label,im,mt,pr,nuc_mask)\n analyseMicrotubulesFullCell(image_data,contours,label,im,mt,pr,nuc_mask)\n analyseER(image_data,contours,label,im,er,pr,nuc_mask)\n analyseNU(image_data,contours,label,im,nu,pr, nuc_mask)\n analyseMTNU(image_data,nuc_cont,label,im,nu,mt, nuc_mask)\n analyseNucleus(image_data,label,im, nu, er, mt, pr, nuc_cont, nuc_midpoint, nuc_mask, nuc_rect)\n #print(f'len(nuc_cont),{len(nuc_cont)}')\n analyseNucleusInnerRims(image_data, nuc_cont, nuc_midpoint, nu, pr)\n analyseNucleusOuterRims(image_data, nuc_cont, contours, nuc_midpoint, nu, pr)\n lookForAggresomes(image_data, pr)\n lookForNuclearBodies(image_data, pr, nuc_mask)\n lookForNuclearSpeckles(image_data, pr, nuc_mask)\n \n #if label == '2':\n if debug and debug_show_input_images:\n printhist(nu)\n printhist(pr)\n printhist(er)\n printhist(mt)\n displaycopy= np.copy(nu)\n erdisplaycopy= np.copy(er)\n mtdisplaycopy= np.copy(mt)\n cv2.drawContours(displaycopy, nuc_cont,-1,255,1)\n cv2.drawContours(displaycopy, contours,-1,255,1)\n cv2.drawContours(pr, nuc_cont,-1,255,1)\n cv2.drawContours(pr, contours,-1,255,1)\n cv2.drawContours(erdisplaycopy, nuc_cont,-1,255,1)\n cv2.drawContours(erdisplaycopy, contours,-1,255,1) \n cv2.drawContours(mtdisplaycopy, contours,-1,255,1)\n cv2.drawContours(mtdisplaycopy, nuc_cont,-1,255,1)\n print(\"Nucleus\")\n do_display_image_big(displaycopy,image_data)\n print(\"Protein\")\n do_display_image_big(pr,image_data)\n print(\"Endoplasmic reticulum\")\n do_display_image_big(erdisplaycopy,image_data)\n print(\"Microtubules\")\n do_display_image_big(mtdisplaycopy,image_data)\n\n\n\n\n calculateScores(image_data)\n\n\n\n return \n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def writeSubmissionFile (i, image_data_list, masks):\n if i == 0:\n f = open('submission.csv', 'w')\n else:\n f = open('submission.csv', 'a')\n with f as outf:\n if i == 0:\n print('ID,ImageWidth,ImageHeight,PredictionString', file=outf) \n for image_data, cellmasklist in tqdm(zip(image_data_list,masks)): \n image_id = image_data['image_id']\n w = image_data['width']\n h = image_data['height'] \n pred_strs=''\n c = 0\n\n for rle, cell in zip(cellmasklist, image_data['cells']):\n if False and debug:\n c = c + 1\n rle = f'cell({i},{c})'\n if len(cell['cell_results']) == 0:\n class_id = 18\n cnf = 0.0001 \n pred_strs += (f' {class_id} {cnf} {rle}')\n else:\n for cell_result in cell['cell_results']:\n class_id, cnf = cell_result\n cnf = round(cnf,4) \n pred_strs += (f' {class_id} {cnf} {rle}') \n print(f'{image_id},{w},{h},{pred_strs}', file=outf)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def outputImageData(image_data_list):\n image_field_excludes=['cells']\n cell_field_excludes=['nucleoli']\n nuc_field_excludes=[]\n \n if train_or_test == 'train' or debug:\n with open('cellstats.txt', 'a') as outf:\n i = newImageDict('', '',0,0)\n c = newCellDict('', '',0,0,0,0)\n n = newNucleolusDict('', '')\n \n header = '' \n for field in i:\n if field not in image_field_excludes:\n header += f'Img {field};'\n for field in c:\n if field not in cell_field_excludes:\n header += f'Cel {field};'\n for field in n:\n if field not in nuc_field_excludes:\n header += f'Nuc {field};'\n print(header, file = outf)\n \n for image in image_data_list:\n for cell in image['cells']:\n for nuc in cell['nucleoli']:\n row = ''\n for field in i:\n if field not in image_field_excludes:\n row += f'{image[field]};' \n for field in c:\n if field not in cell_field_excludes:\n row += f'{cell[field]};' \n for field in n:\n if field not in nuc_field_excludes:\n row += f'{nuc[field]};'\n print(row, file = outf)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\ncellsegmentor = CellSegmentator()\ncellsegmentorOpenCV = CellSegmentatorOpenCV()\n\ndata_size = len(data_df)\n\nbs = 150\n\ndef load_images(df : pd.DataFrame, root=f'../input/hpa-single-cell-image-classification/{train_or_test}/', resize = True):\n gray = []\n ryb = []\n protein=[]\n gray_orig = []\n ryb_orig = []\n protein_orig=[]\n image_data_list = []\n imWidth = []\n imHeight = [] \n \n if train_or_test == 'train' or train_or_test == 'traintest':\n root=f'../input/hpa-single-cell-image-classification/train/'\n else:\n root=f'../input/hpa-single-cell-image-classification/test/'\n \n for i, row in tqdm(df.iterrows(), total=len(df)):\n r = os.path.join(root, f'{row.ID}_red.png')\n y = os.path.join(root, f'{row.ID}_yellow.png')\n b = os.path.join(root, f'{row.ID}_blue.png')\n g = os.path.join(root, f'{row.ID}_green.png')\n r = cv2.imread(r, 0)\n y = cv2.imread(y, 0)\n b = cv2.imread(b, 0)\n g = cv2.imread(g, 0)\n\n \n #do_display_image_big(b)\n #do_display_image_big(g)\n #h1 = np.histogram(g, bins=25)\n #print(h1)\n #g[g>100]=0\n #h2 = np.histogram(g, bins=25) \n #print(h2)\n #g2=g*2.5\n #do_display_image_big(g2)\n \n \n ryb_orig_image = np.stack((r, y, b), axis=2)\n if resize:\n gray_image = cv2.resize(b, (512, 512))\n\n ryb_image = cv2.resize(ryb_orig_image, (512, 512))\n gray.append(gray_image)\n ryb.append(ryb_image)\n protein_image = cv2.resize(g, (512, 512))\n protein.append(g)\n #else:\n ryb_image_orig = np.stack((r, y, b), axis=2)\n gray_orig.append(b)\n ryb_orig.append(ryb_orig_image) \n protein_orig.append(g)\n \n if train_or_test == 'train' or train_or_test == 'traintest':\n label = row.Label\n else:\n label = '' \n\n ann = newImageDict(row.ID, label, b.shape[1],b.shape[0])\n\n image_data_list.append(ann)\n return gray, ryb, protein, gray_orig, ryb_orig, protein_orig,image_data_list\n\nwith open('cellstats.txt', 'w') as outf:\n print('', file = outf)\n \nfalsePos=np.zeros(19)\nfalseNeg=np.zeros(19)\ntruePos=np.zeros(19)\ntrueNeg=np.zeros(19)\n\nprint(data_df)\n\nfor i in range(0, data_size, bs):\n print('!!!!', i, '!!!!')\n start = i\n end = min(len(data_df), start + bs)\n batch_df = data_df[start:end]\n print(len(batch_df))\n print(f'---- start load images ----{start}, {end}')\n gray, ryb, pr, gray_orig, ryb_orig, pr_orig, image_data_list = load_images(batch_df,resize = True)\n\n\n idn=len(image_data_list)\n print( f'image_data_list {idn}')\n \n print(len(gray))\n print(f'---- finish load images --{start}, {end}')\n print(f'---- start pred nuclei --{start}, {end}')\n nuc_segmentations = cellsegmentor.pred_nuclei(gray)\n nsn=len(nuc_segmentations)\n print( f'nuc_segmentations{nsn}') \n print(nuc_segmentations[0].shape)\n print(f'---- finish pred nuclei --{start}, {end}')\n print(f'---- start pred cells --{start}, {end}')\n cell_segmentations = cellsegmentor.pred_cells(ryb, precombined=True)\n \n nsn=len(nuc_segmentations)\n csn=len(cell_segmentations)\n print( f'nuc_segmentations{nsn} cell_segmentations {csn}')\n \n print(cell_segmentations[0].shape)\n print(f'---- finish pred cells --{start}, {end}')\n print(f'---- start pred cells OpenCV --{start}, {end}')\n cellsegmentorOpenCV.pred_cells(ryb_orig,pr_orig,image_data_list,precombined=True)\n print(f'---- finish pred cells --{start}, {end}')\n print(f'---- start calculating masks ----{i}')\n \n outputImageData(image_data_list)\n kernel = np.ones((5, 5), 'uint8') \n maskRLElists=[]\n zcount=0\n nsn=len(nuc_segmentations)\n csn=len(cell_segmentations)\n \n idn=len(image_data_list)\n print( f'nuc_segmentations{nsn} cell_segmentations {csn} \\\n image_data_list {idn}')\n for image_data, nuc_seg, cell_seg in tqdm (zip(image_data_list, nuc_segmentations, \\\n cell_segmentations )):\n zcount= zcount+1\n maskRLElist = []\n maskRLElists.append(maskRLElist)\n ll=len(maskRLElists)\n imid=image_data['image_id']\n\n #print(f'{zcount}: {imid} mask list length {ll}')\n nuc_small, cell_small = label_cell(nuc_seg, cell_seg)\n \n #orig_size = gray_orig[i].shape[0]\n orig_size = image_data['width']\n cell=cv2.resize(cell_small,(orig_size,orig_size), cv2.INTER_NEAREST_EXACT)\n #print(image_data['image_id'])\n cell_data=image_data['cells']\n for cell_data in image_data['cells']:\n nuc_midpoint = cell_data['nuc_midpoint']\n mx,my=nuc_midpoint\n cell_mask_label=cell[my,mx]\n maskImage = np.zeros((orig_size,orig_size), np.uint8)\n maskImageBool = np.zeros((orig_size,orig_size), np.bool8)\n maskImage[cell==cell_mask_label]=255\n #do_display_image_big(maskImage)\n # erode to remove spurious effects from previous resize\n ed=cv2.erode(maskImage,kernel,iterations=2)\n\n maskImageBool[ed==255]=True\n #do_display_image_big(maskImageBool)\n rle = encode_binary_mask(maskImageBool)\n maskRLElist.append(rle)\n \n \n print(f'---- finish calculating masks ----{i}')\n print(f'---- start writing batch to file ----{i}')\n print('image_data length', len(image_data_list))\n #print('mask lists length', len(maskRLElists))\n writeSubmissionFile ( i, image_data_list, maskRLElists)\n print(f'---- finish writing batch to file ----{i}')\n\n \n if train_or_test == 'train' or train_or_test == 'traintest':\n countTrainingScores(image_data_list,falsePos,falseNeg,truePos,trueNeg)\n\n \n\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!wc -l submission.csv","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(len(cell_data))","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}