Search code examples
pythonopencvconvex-hullconvexity-defects

How to count protuberances in dendrite with openCV


I'm trying to count dendritic spines (the tiny protuberances) in mouse dendrites obtained by fluorescent microscopy, using Python and OpenCV.

Here is the original image, from which I'm starting:

Raw picture:

enter image description here

After some preprocessing (code below) I've obtained these contours:

Raw picture with contours (White):

enter image description here

What I need to do is to recognize all protuberances, obtaining something like this:

Raw picture with contours in White and expected counts in red:

enter image description here

What I intended to do, after preprocessing the image (binarizing, thresholding and reducing its noise), was drawing the contours and try to find convex defects in them. The problem arose as some of the "spines" (the technical name of those protuberances) are not recognized as they en up bulged together in the same convexity defect, underestimating the result. Is there any way to be more "precise" when marking convexity defects?

Raw image with contour marked in White. Red dots mark spines that were identified with my code. Green dots mark spines I still can't recognize:

enter image description here

My Python code:

import cv2
import numpy as np
from matplotlib import pyplot as plt

#Image loading and preprocessing:

img = cv2.imread('Prueba.jpg')
img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)
img = cv2.pyrMeanShiftFiltering(img,5,11)
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)

ret,thresh1 = cv2.threshold(gray,5,255,0)

kernel = cv2.getStructuringElement(cv2.MORPH_ELLIPSE,(5,5))
img1 = cv2.morphologyEx(thresh1, cv2.MORPH_OPEN, kernel)
img1 = cv2.morphologyEx(img1, cv2.MORPH_OPEN, kernel)
img1 = cv2.dilate(img1,kernel,iterations = 5)

#Drawing of contours. Some spines were dettached of the main shaft due to
#image bad quality. The main idea of the code below is to identify the shaft
#as the biggest contour, and count any smaller as a spine too.

_, contours,_ = cv2.findContours(img1,cv2.RETR_TREE,cv2.CHAIN_APPROX_SIMPLE)
print("Number of contours detected: "+str(len(contours)))

cv2.drawContours(img,contours,-1,(255,255,255),6)
plt.imshow(img)
plt.show()

lengths = [len(i) for i in contours]
cnt = lengths.index(max(lengths))

#The contour of the main shaft is stored in cnt    

cnt = contours.pop(cnt)

#Finding convexity points with hull:

hull = cv2.convexHull(cnt) 

#The next lines are just for visualization. All centroids of smaller contours
#are marked as spines.

for i in contours:

    M = cv2.moments(i)
    centroid_x = int(M['m10']/M['m00'])
    centroid_y = int(M['m01']/M['m00'])
    centroid = np.array([[[centroid_x, centroid_y]]])

    print(centroid)

    cv2.drawContours(img,centroid,-1,(0,255,0),25)
    cv2.drawContours(img,centroid,-1,(255,0,0),10)

cv2.drawContours(img,hull,-1,(0,255,0),25)
cv2.drawContours(img,hull,-1,(255,0,0),10)

plt.imshow(img)
plt.show()

#Finally, the number of spines is computed as the sum between smaller contours
#and protuberances in the main shaft.

spines = len(contours)+len(hull)
print("Number of identified spines: " + str(spines))

I know my code has many smaller problems to solve yet, but I think the biggest one is the one presented here.


Solution

  • I would approximate the contour to a polygon as Silencer suggests (don't use the convex hull). Maybe you should simplify the contour just a little bit to keep most of the detail of the shape.

    This way, you will have many vertices that you have to filter: looking at the angle of each vertex you can tell if it is concave or convex. Each spine is one or more convex vertices between concave vertices (if you have several consecutive convex vertices, you keep only the sharper one).

    EDIT: in order to compute the angle you can do the following: let's say that a, b and c are three consecutive vertices

    angle1 = arctan((by-ay)/(bx-ax))
    angle2 = arctan((cy-by)/(cx-bx))
    angleDiff=angle2-angle1
    if(angleDiff<-PI) angleDiff=angleDiff+2PI
    
    if(angleDiff>0) concave
    Else convex
    

    Or vice versa, depending if your contour is clockwise or counterclockwise, black or white. If you sum all angleDiff of any polygon, the result should be 2PI. If it is -2PI, then the last "if" should be swapped.