一、基于單張圖片的分析
對于這張圖片，想要進行細胞核細胞質面積比的計算，一種常用的手段是通過識別影像的灰度，從而劃分圖片中各個部分，
因此首先要做的就是對影像進行二值化處理，影像二值化在數字影像處理領域是一種非常常見的技術，以設定的閾值為參考，將影像上的像素灰度值設定為0或255，便可以區分出細胞核、細胞質和圖片背景，
使用OpenCV庫或Pillow庫均能完成二值化的作業，使用Pillow庫完成這個作業時，我們只需要先將影像轉換為黑白圖片，然后根據每個像素點的明度判斷它是屬于細胞核還是細胞質，值得注意的是，針對每張圖的灰度閾值是不一樣的，因此需要經過多次調整與試錯才能得到合適的閾值，
或者可以通過PS幫助找出最佳的閾值，分別使用選擇工具選中整個細胞和細胞核，調出直方圖工具，就可以查到區域的灰度分布，
可以看到整個細胞的灰度值分布在225以下，而細胞核的灰度分布在153以下，甚至連像素數量都拿到了，因此可以直接粗略地計算出來它的細胞核和細胞質的面積比，大約是0.18，
接下來以這個資料為一個參考值，就可以方便地驗證我們的演算法是否正確，
from PIL import Image
img = Image.open("../cell.jpg")
img = img.convert('L')
# 獲取圖片大小
width = img.size[0]
height = img.size[1]
pix = img.load()
img.show()
cell_whole = 0
cell_nucleus = 0
pixels = 0
for x in range(width):
for y in range(height):
pixels = pixels + 1
grayScale = pix[x, y]
if grayScale <= 225:
cell_whole = cell_whole + 1
if grayScale <= 150:
cell_nucleus = cell_nucleus + 1
print("細胞質面積：",cell_whole)
print("細胞核面積：",cell_nucleus)
print("核質面積比：{:.3f}".format(cell_nucleus / (cell_whole - cell_nucleus)))
基于Pillow完成的測量效果并不理想，最后的結果大概是這樣：
和預期的0.19差得有點遠了，原因其實很好解釋，因為Pillow是對整個圖片中所有的像素點進行統計，因此細胞質中一些較深的像素點被不可避免地統計在了細胞核的面積中，導致細胞核的面積偏大，
既然純粹使用灰度的統計方法無法完成這項作業，那么單純使用Pillow庫就沒法完成這項任務了，接下來就請出我們的OpenCV庫，
OpenCV在Pillow的基礎上多了很多與計算機視覺有關的庫，比如說，可以在灰度識別、二值化的基礎上，再提取出圖形的輪廓，然后相加，
但是可以注意到，在識別的時候，OpenCV識別出的輪廓不只有細胞核，還有細胞質中的一些深色色塊，這也是Pillow庫識別的時候細胞核面積偏大的原因，解決方法也很簡單，只需要找出最大的輪廓就可以了，剩下的深色色塊不需要算進細胞核面積中，
通過這種方法，我們得到了面積比0.175，和之前得到的0.18參考值差不多，在合理的誤差范圍內，
至此，我們已經完成了關于這張圖的細胞核和細胞質面積比的計算，
二、更一般性的考慮
再又進一步考慮，如果上傳的圖片是任意一張單個細胞的圖片，該怎么獲取細胞的核質面積比呢？
于是我與一位大佬經過一個晚上加一個大半天的討論與研究，成功找到了一套能夠獲取任意單個細胞圖片核質面積比的方法，
基于上面已有的獲取單個細胞的方法，獲取多個細胞的難點在于如何獲得兩個灰度值的閾值，這兩個閾值決定著能否正確從圖片中分出細胞質與細胞核，而對于每一張細胞的圖片，它的閾值都是不同的，
起初我們的思路是通過分析整張細胞圖片的灰度值分布，嘗試找到一些特殊的規律，從而找到那兩個閾值，
確實很難從這種雜亂無章的曲線中找出些什么，但實際上它確實是存在著某種規律，第一個細胞的兩個閾值分別是225、153，第二個細胞的兩個閾值分別是252、192，它們大致都處于一個波峰的兩端，這個波峰恰巧是代表細胞質灰度值的像素點的波峰，而它們的左側都有一個波峰，代表細胞核灰度值像素點的波峰，它們右側代表背景的波峰因為數值實在是過大，會影響整個圖的觀察，因此被我們手工切除了，
那天晚上我們討論到很晚，從曲線的平滑化一直討論到函式的凹凸性，想通過這個灰度值找出合適的閾值，一直到最后都還是找不到合理的方案，因為整個資料的噪聲太大了，我們能想到的演算法都非常容易受到干擾，
第二天早上，我帶著一個全新的方案，根據我們找到的兩個圖片，我決定使用橢圓擬合的方法去找到合適的閾值點，
具體的思路就是使用橢圓去擬合根據灰度閾值二值化后的圖片所生成的輪廓，然后比較擬合出的橢圓與輪廓的匹配程度，匹配度越高則說明擬合出的輪廓越規則，則越有可能是正確的灰度值，
這是一張在0-255的不同灰度值作為閾值的情況下，通過二值化得到的輪廓與擬合出的橢圓的誤差，這張圖就明顯比剛剛那個圖噪聲小了很多，兩個我們想要的閾值剛好在從右到左數的兩個波底上，
然后我們想辦法整了個計算方法，即比較右側曲線的平均斜率和左側曲線的平均斜率，如果它們的比大于資料的平均變化程度，并且擬合的誤差小于某個值，那么就認為這是我們要找的灰度閾值，
并且考慮到這樣的演算法并不一定100%準確，因此我們用tkinter寫了一個簡單的GUI界面，以保證在自動閾值檢測失效的情況下，仍然能通過手動調節得到合適的核質面積比，
最后得到了這樣的一個簡單的小工具，當匯入圖片時，它會根據我們寫的小演算法檢測出可能合適的灰度閾值，并直接應用，之后仍然可以通過左中右三張輔助圖和中間的滑塊來微調閾值或者重設閾值，以得到更加準確的面積比，
我們還在網路上多找了幾個細胞圖片進行測驗，發現成功率確實超出了我們原本的預估，但也沒辦法避免一些失敗的情況，這些失敗的情況普遍是只找到了細胞核或細胞質中的其中一個，但畢竟我們的程式只花了一個晚上加一個大半天，所以還有很多思考和改進的空間，
因為我們之前也沒有接觸過太多有關機器視覺的東西，所以這兩天折騰出來的東西純屬是憑借想象力和肝，以后如果進一步系統地去學習這方面的知識，應該又會有不一樣的思路和想法吧，
三、代碼實作
import tkinter
from tkinter import *
from PIL import Image, ImageTk ###這個是沒有想到的模塊，也不確定能不能省
from tkinter.filedialog import askopenfilename
import time
import cv2 as cv
import numpy as np
import math
from time import sleep
import matplotlib.pyplot as plt
# from scipy.signal import savgol_filter
root = Tk()
root.geometry('500x600') ##這個小了一點，不知道怎么自適應
root.title('細胞胞質比計算工具')
def slopeComparison(i,ls,ls1):
left = abs(ls[i] - ls[i - 1]) + abs(ls[i-1] - ls[i-2])
right = ls[i+2] - ls[i]
mean = np.mean(ls1)
mean2 = np.mean(ls)
# 如果右側平均變化率和左側平均變化率的比值大于平均變化程度，并且左邊不是0，并且誤差的值不能太大
if right > mean * left and left != 0 and (abs(ls[i] - ls[i - 1]) > 0.01) and 0 < ls[i] < mean:
return True
else:
return False
def choosepic():
global scaleX1
global scaleX2
global cell_o
global words
global rec
path_ = askopenfilename()
print(path_)
rec['text'] = "推薦數值："
cell_o = cv.imread(path_, 1)
changePic(cell_o)
ls = []
for i in range(0, 256):
scaleX1.set(i)
scaleX2.set(154)
try:
whole_area, nucleus_area, whole_area_ellipse, nucleus_area_ellipse = runDetect(i, 154)
ls.append(abs(whole_area - whole_area_ellipse))
except Exception:
ls.append(0.1)
# ls = list(savgol_filter(ls, 9, 5, mode= 'nearest'))
flag = False
ls1 = []
amount = 0
for i in range(len(ls)-1):
ls1.append(abs(ls[i+1] - ls[i]))
for i in list(range(2, len(ls) - 2))[::-1]:
if slopeComparison(i,ls,ls1):
if flag == False:
flag = True
print(i-1)
print("{} {:.2f}".format(i, ls[i]))
rec['text'] = rec['text'] + " " + str(i-1)
amount = amount + 1
if amount == 1:
scaleX1.set(i-1)
elif amount == 2:
scaleX2.set(i-1)
runDetect(scaleX1.get(),scaleX2.get())
else:
flag = False
# for i in range(0, len(ls)):
# print("{} {:.2f}".format(i, ls[i]))
# plt.plot(range(0, 256), ls, linewidth=1, color='r', marker='o', markersize=0)
# plt.show()
def changePic(cell_o):
dst = cv.cvtColor(cell_o, cv.COLOR_BGR2RGB)
up_width = 400
up_height = 400
up_points = (up_width, up_height)
dst = cv.resize(dst, up_points, interpolation=cv.INTER_CUBIC)
current_image = Image.fromarray(dst)
imgtk = ImageTk.PhotoImage(image=current_image)
image_label.config(image=imgtk)
image_label.image = imgtk # keep a reference
cell_o = cv.imread("../../cell.jpg", 1)
def runDetect(x1, x2):
global cell_o
global words
cell = cv.cvtColor(cell_o, cv.COLOR_BGR2GRAY)
rows, cols = cell_o.shape[:2] # 獲得圖片的高和寬
mapx = np.zeros(cell_o.shape[:2], np.float32) # mapx是一個圖片大小的變換矩陣
mapy = np.zeros(cell_o.shape[:2], np.float32) # mapy也是一個圖片大小的變換矩陣
for i in range(rows): # 遍歷，將獲得到圖片的寬和稿依次賦予兩個變換矩陣
for j in range(cols):
mapx.itemset((i, j), j)
mapy.itemset((i, j), i)
img = cv.remap(cell_o, mapx, mapy, cv.INTER_LINEAR) # 用仿射函式實作仿射變換
# 將影像轉換為二值圖，第二個值為閾值
ret, cell_2v_whole = cv.threshold(cell, x1, 255, cv.THRESH_BINARY_INV) # 獲取整個細胞的二值圖
ret, cell_2v_nucleus = cv.threshold(cell, x2, 255, cv.THRESH_BINARY_INV) # 獲取細胞核的二值圖
cv.namedWindow("cell_whole", 0)
cv.resizeWindow("cell_whole", 500, 500)
cv.imshow("cell_whole", cell_2v_whole)
cv.namedWindow("cell_nucleus", 0)
cv.resizeWindow("cell_nucleus", 500, 500)
cv.imshow("cell_nucleus", cell_2v_nucleus)
# 基于二值影像用外輪廓的模式，通過全點連接輪廓的方法提取輪廓
contours_whole, hierarchy = cv.findContours(cell_2v_whole, cv.RETR_EXTERNAL, cv.CHAIN_APPROX_NONE)
contours_nucleus, hierarchy = cv.findContours(cell_2v_nucleus, cv.RETR_EXTERNAL, cv.CHAIN_APPROX_NONE)
# 獲取目標影像的最小矩陣，找到最大的輪廓
maxArea = 0
maxContour = {}
for item in contours_whole:
x, y, w, h = cv.boundingRect(item)
if w * h > maxArea:
maxArea = w * h
maxContour = item
whole_area = 0
whole_area += cv.contourArea(maxContour)
img = cv.drawContours(cell_o.copy(), maxContour, -1, (0, 255, 0), 1)
cnt = maxContour
ellipse = cv.fitEllipse(cnt)
(x, y), (a, b), angle = cv.fitEllipse(cnt)
cv.ellipse(img, ellipse, (255, 0, 0), 1)
whole_area_ellipse = math.pi * a * b / 4
# 重復，查找細胞核的最大輪廓
maxArea = 0
maxContour = {}
for item in contours_nucleus:
x, y, w, h = cv.boundingRect(item)
if w * h > maxArea:
maxArea = w * h
maxContour = item
nucleus_area = 0
nucleus_area += cv.contourArea(maxContour)
img = cv.drawContours(img, maxContour, -1, (0, 0, 255), 1)
cnt = maxContour
ellipse = cv.fitEllipse(cnt)
(x, y), (a, b), angle = cv.fitEllipse(cnt)
cv.ellipse(img, ellipse, (255, 0, 0), 1)
nucleus_area_ellipse = math.pi * a * b / 4
# img = cv.rectangle(img,(x,y),(x+w,y+h),(0,255,0),1)
changePic(img)
# 繪制目標框
words['text'] = "細胞質面積：{} 細胞核面積：{} \n 擬合面積：{:.3f} {:.3f} \n 核質面積比：{:.3f}".format(whole_area, nucleus_area,
whole_area_ellipse,
nucleus_area_ellipse,
nucleus_area_ellipse / (
whole_area_ellipse - nucleus_area_ellipse))
return whole_area, nucleus_area, whole_area_ellipse, nucleus_area_ellipse
def handleLeftButtonReleaseEvent(e):
global scaleX1
global scaleX2
runDetect(scaleX1.get(), scaleX2.get())
path = StringVar()
Button(root, text='選擇圖片', command=choosepic).pack()
file_entry = Entry(root, state='readonly', text=path)
# file_entry.pack()
image_label = Label(root)
image_label.pack()
scaleX1 = Scale(root, from_=0, to=255, orient=HORIZONTAL, length=255)
scaleX2 = Scale(root, from_=0, to=255, orient=HORIZONTAL, length=255)
scaleX1.bind("<ButtonRelease-1>", handleLeftButtonReleaseEvent)
scaleX2.bind("<ButtonRelease-1>", handleLeftButtonReleaseEvent)
scaleX1.pack()
scaleX2.pack()
words = Label(root)
words.pack()
rec = Label(root)
rec.pack()
# plt.plot(list(ls),list(range(1,256)),linewidth=1,color='r',marker='o',markersize=0)
# plt.xlabel('Position(2-Theta)')
# plt.ylabel('Intensity')
# plt.legend()
# plt.show()
root.mainloop()
轉載請註明出處，本文鏈接：https://www.uj5u.com/qita/348527.html
標籤：其他