一种基于SVM的数字PCR假阳性液滴检测方法与流程

一种基于svm的数字pcr假阳性液滴检测方法
技术领域
1.本发明属于图像处理技术领域，具体涉及一种基于svm的数字pcr假阳性液滴检测方法。


背景技术：

2.聚合酶链式反应(pcr)是一种生物体外实现特定dna片段扩增的技术。它具有操作简便、检测灵敏、特异性强等特点，广泛用于临床诊断、法医学调查、生物技术等领域。
3.数字pcr基于聚合酶链式反应的原理精确的测定基因拷贝数，实现对基因突变定性和定量分析。微滴式数字pcr系统是用油水两相间隔得到的液滴当作微反应单元，将dna或 rna样本稀释并分散成数万甚至数百万个独立的反应单元，每个反应单元中包含(或不包含) 1个或1个以上目标分子(dna或rna模板)。微反应单元单层平铺于芯片内，将芯片放入升降温装置，针对所有的微反应单元的靶序列进行分子模板pcr扩增，扩增过程中微反应单元的荧光强度不断增强，完成扩增后。
4.ccd或者cmos相机在光学成像系统中采集用于特定基因片段标记的荧光探针信号，依据荧光信号强弱判断微反应单元的阳性或阴性，最终基于统计学分析(泊松分布)检测样本中核苷酸的浓度。数字pcr不依赖扩增曲线的循环阈值，不受扩增效率的影响，能够直接读出dna分子的个数，是一种核酸分子绝对定量技术。


技术实现要素：

5.微滴式数字pcr低浓度(阳性点数量少)条件下具有极好的灵敏度和准确度。由于在实验室污染或者操作员操作有误，导致配置的试剂可能混入杂质。由于芯片生成过程中质控不严格，导致芯片表面出现划痕。由于清洗不干净或者静电，导致滤波器等成像系统出现粉尘。上述现象在液滴图片形成类似液滴轮廓的异常区域，这些区域可能误判为阳性液滴，影响数字pcr检测的灵敏度和准确度。本发明采用图像处理和机器学习的方法检测这些异常区域，能够准确识别假阳和真阳。
6.本发明的第一个方面公开了一种基于svm的数字pcr假阳性液滴检测方法，包括以下步骤：
7.s1：构建样本库并基于svm训练样本
8.s11：更换不同的滤光片和光源拍摄不同条件下的微液滴图片。拍摄芯片表面有缺陷的液滴图片。
9.s12：截取完整的液滴组成正样本库，归一化为同一尺寸24*24；负样本是芯片表面缺陷导致的“高亮点”、“划痕”、“凹凸点”等情况，负样本归一化为同一尺寸24*24；采用旋转等方法扩增样本量。
10.s13：计算样本的特征值，基于hog特征算子计算样本的特征值，其中hog检测窗口的大小为24*24，窗口滑动块的大小为8*8，滑动块的步长为2，细胞单元的大小为4*4，直方图的bins的数量为9个，特征值的维数为3888。设si表示第i样本， 0≤i≤num_sample，num_
sample，表示正负样本的总数；表示第i样本第k维的特征值， 0≤k≤3888。对于液滴的正、负样本集(s0，cl0)
…
，(si，cli)，...，(s
num_sample
，cl
num_sample
)，其中si表示第i样本，cli＝-1表示负样本，cli＝1表示正样本。
11.s14：生成训练分类器，样本的特征值导入支持向量机训练，生成分类器svm_class.xml 文件，其中选择二分类的svm，核函数为线性核函数，惩罚因子为1，循环迭代的次数为10000。
12.优选地，对s12中的正样本用高斯滤波处理，对s12中的负样本进行直方图均衡化。
13.s2：液滴图片二值化
14.s21：设图像droplet_image，长度为image_height，图像的宽度为image_width，图片横坐标等分num_height份，纵坐标等分num_width份，则图像droplet_image总共分成num_size＝num_width*num_height张子图片，子图片的表示为sub_imager，其中 0≤r≤num_size，子图片高度为子图片宽度为子图片宽度为20≤num_width≤32，16≤num_height≤24，且num_height＜num_width。
15.s22：任选一张子图片sub_imager，其中0≤r≤num_size。g(sub_x，sub_y)表示像素值在点(sub_x，sub_y)的灰度值，0≤sub_x≤sub_cols，0≤sub_y≤sub_rows。计算子图片阈值，生成此子图片的二值化图像：
16.1)统计图片sub_imager灰度直方图：
[0017][0018][0019]
2)针对直方图进行归一化：
[0020]
forii 0：255
[0021][0022]
3)计算概率密度累加和：
[0023][0024]
4)从前往后遍历数组norm_histomw，找到第一个非零的数据位置设置为first_bin，从后往前遍历数组norm_histomr，找到第一个非零的数据位置设置为last_bin；
[0025]
5)计算阈值的位置点：
[0026][0027][0028]
6)生成sub_imager的二值化图像sub_binaryr：
[0029][0030]
7)依据sub_binaryr连通域sub_img_domainr对二值化图像进行修改：
[0031]
设sub_binaryr连通域的像素面积为sub_area_domainr，连通域的孔洞的数量为sub_img_holer，如果或者sub_img_holer＞10，则sub_binaryr的二值化图片sub_binaryr全部像素值为零。
[0032]
s23：通过上面的操作，num_size张子图片对应num_size张二值化图像，所有的二值化图像按原图片droplet_image的位置组合成的二值化image_binary。
[0033]
s3连通域处理
[0034]
s31：设表示二值化图像image_binary的第m个连通域，其中0≤m≤num
i_domain
，num
i_domain
表示连通域的个数。的像素面积为的像素面积为的像素面积为表示第m个连通域像素单元的个数。的像素周长为计算第m个连通域的圆率如果则置第m个连通域所包含位置的像素值为0。
[0035]
s32：连通均中像素集合表示为如果第m个像素单元属于则表示为其中设集合中像素单元横坐标的最大值表示为dmn_x_maxm，横坐标的最小值表示为dmn_x_minm，纵坐标的最大值表示为dmn_y_maxm，纵坐标的最小值表示为dmn_y_minm。
[0036]
s33：如果dmn_x_max
m-dmn_x_minm＞30，则计算连通域在x方向的投影值，其中最小投影值对应的x的坐标值为x_min_locm。基于坐标值x_min_locm沿y方向分割连通域连通域如图2所示。
[0037]
s34：如果dmn_y_max
m-dmn_y_minm＞30，则计算连通域在y方向的投影值，其中最小投影值对应的y的坐标值为y_min_locm。基于坐标值y_min_locm沿x方向分割连通域连通域
[0038]
s35：计算连通域的质心坐标
[0039][0040][0041]
计算连通域的像素单元到质心的最大距离：
[0042][0043][0044]
计算外接圆的面积：
[0045]
如果且b_curveradiom＜0.8，则在图像learning_image_binary的对应
位置标记连通域
[0046]
s4判断假阳性，其中s42和s43是基于液滴像素值的标准差判断液滴假阳性；s44-s46 基于训练的分类器判断液滴假阳性；s47-s411基于液滴局部的像素对比度判断假阳性。
[0047]
s41：任选液滴dropletsi，其中0≤i≤n
droplets
，n
droplets
表示液滴的总数量。如果液滴is_droplets[i]＝1，液滴dropletsi位置(pos_xi，pos_yi)在二值化图像image_binary中对应的像素值为255，则更新is
droplets[i]
＝-1。如果is_droplets[i]＝1则表示真阳性液滴，is
droplets[i]
＝-1则表示假阳性液滴。
[0048]
s42：设液滴dropletsi位置(pos_xi，pos_yi)在二值化图像learning_image_binary中对应的像素值为255，则在灰度图片droplets_image中，以dropletsi的位置(pos_xi，pos_yi)为中心，以4个像素单元为半径，计算此区域内像素的平均值：
[0049]
其中gray(x_lox，y_loc)表示在图像droplets_image位置点(x_lox，y_loc)的像素值，(x_loc-pos_xi)2 (y_loc
‑ꢀ
pos_yi2≤16，numcirclei表示此圆形区域内的像素个数。
[0050]
s43：计算此区域的最大值计算此区域的最小值计算此区域的像素的标准差：
[0051]
如果则is_droplets[i]＝-1，如果则is_droplets[i]＝-1。
[0052]
s44：设表示二值化图像learning_image_binary的第l个连通域，其中0≤l≤num
c_domain
，num
c_domain
表示连通域的个数。
[0053]
计算的像素面积积表示第1个连通域像素单元的个数。连通域中像素集合表示为第l个像素单元属于表示为其中计算连通域的质心坐标
[0054][0055][0056]
计算连通域的像素单元到质心的最大距离：
[0057][0058]
其中计算此连通区域的圆度：如果液滴dropletsi位置点(pos_xi，pos_yi)属于连通域的集合且则is_droplets[i]＝-1；
[0059]
s45：任选一个液滴dropletsi，其中0≤i≤n
droplets
，n
droplets
表示液滴的总数量。液滴dropletsi与连通域的距离为：
[0060][0061]
如果则以(pos_xi，pos_yi)为中心，24个像素单元为矩形宽度，24个像素单元为矩形高度，截取灰度图片droplets_image，形成新的样本图像sample_imgi，sample_imgi表示第i个液滴的样本，样本图像sample_imgi的高度和宽度都为 24个像素。
[0062]
s46：载入训练生成的svm_class.xml文件，计算sample_imgi的hog特征，
[0063]
特征值导入到线性核函数svm中，如果判断结果为1，则is_droplets[i]＝1，否则为is_droplets[i]＝-1。
[0064]
s47：以连通域的质心坐标为中心，以 domain_radius1为半径，其中统计此圆形区域内对应的液滴像素值的灰度直方图histom
l_gray
，直方图以灰度级为横坐标，像素个数为纵坐标。num_i_gray(g)表示此区域内液滴像素值为g的液滴的个数。统计此圆形区域内所有液滴的灰度值，最大灰度值为gray
max
，最小灰度值为gray
min
， 0≤gray
min
≤gray
max
≤255。0≤l≤num
c_domain
，num
c_domain
表示连通域总数量，con_ratio＝1.1。
[0065]
s48：此圆形区域内液滴最大灰度值和液滴最小灰度值分别为gray
max
和gray
min
，灰
度值 g的概率密度为其中num表示此圆形区域内液滴的总个数，ng表示灰度值为g的像素个数，gray
min
≤g≤gray
max
；计算灰度概率密度直方图ck＝pg，其中gray
min
≤k≤gray
max
&&k＝g；
[0066]
计算灰度直方图一阶矩为：
[0067]
二阶矩为：
[0068]
三阶矩为：
[0069]
s49：计算参考值p_refer＝(z
1-m11)/(z
1-z0)，其中：
[0070][0071][0072]
c1＝(-m3
l
m2
l
*m1
l
)/(m2
l-m1
l
*m1
l
)；
[0073]
c0＝(-m2
l
m3
l
*m1
l
)/(m2
l-m1
l
*m1
l
)；
[0074]
s410：灰度值从gray
min
到gray
max
依次遍历灰度概率密度直方图，如果存在ck＞p_refer，则histom
l_gray
_th＝k，其中gray
min
≤k≤gray
max
。
[0075]
s411：以连通域的质心坐标为中心，以domain_radius
l
为半径，统计此区域内小于或等于histom
l_gray
_th液滴的最大像素值为其中gj是以连通域的质心坐标为中心，以domain_radius
l
为半径的圆形区域对应的液滴的像素值。如果存在液滴且且is_droplets[i]＝-1。
[0076]
本发明的另一个方面公开了一种检测装置，所述检测装置包括上述算法。
[0077]
与现有技术相比，本发明取得了如下的显著优势：
[0078]
1)采用一种新的阈值划分算法，对子区域的图片进行分割，组合为原始图片的二值化图片。依据连通域的条件和投影算法，生成新的二值化图image_binary和learning_image_binary，image_binary用于检测异常区域芯片缺陷或者“杂质”等异常的区域。learning_image_binary用于检测轮廓与真实液滴近似的假阳性液滴。
[0079]
2)为了减少假阳性误判，文中建立了液滴的正负样本库，正样本是真实的液滴，负样本为假阳性液滴，包括芯片缺陷或者“杂质”等异常区域。文中选择hog计算样本的特征，采用线性核函数svm训练样本，提高了假阳性识别率。同时使用液滴像素标准差和局部的像素对比度判断假阳性。
附图说明
[0080]
图1为本发明中所述的基于svm的数字pcr假阳性液滴检测方法的流程图；
[0081]
图2为本发明s33步骤的投影图。
具体实施方式
[0082]
下面将结合附图对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
[0083]
一种基于svm的数字pcr假阳性液滴检测方法，具体包括以下步骤：
[0084]
s1：构建样本库并基于svm训练样本
[0085]
s11：更换不同的滤光片和光源拍摄不同条件下的微液滴图片。拍摄芯片表面有缺陷的液滴图片。
[0086]
s12：截取完整的液滴组成正样本库，归一化为同一尺寸24*24，对正样本用高斯滤波处理；负样本是芯片表面缺陷导致的“高亮点”、“划痕”、“凹凸点”等情况，负样本归一化为24*24，对负样本直方图均衡化；采用旋转、翻转、扭曲、仿射变换、缩放、压缩等方法扩增样本量。其中正样本库的样本数量为5734；负样本库的数量为11760个。
[0087]
s13：计算样本的特征值，基于hog特征算子计算样本的特征值，其中hog检测窗口的大小为24*24，窗口滑动块的大小为8*8，滑动块的步长为2，细胞单元的大小为4*4，直方图的bins的数量为9个，特征值的维数为3888。设si表示第i样本， 0≤i≤num_sample，num_sample，表示正负样本的总数；表示第i样本第k维的特征值， 0≤k≤3888。对于液滴的正、负样本集(s0，cl0)
…
，(si，cli)，...，(s
num_sample
，cl
num_sample
)，其中si表示第i样本，cli＝-1表示负样本，cli＝1表示正样本。
[0088]
s14：生成训练分类器，样本的特征值导入支持向量机训练，生成分类器svm_class.xml 文件，其中选择二分类的svm，核函数为线性核函数，惩罚因子为1，循环迭代的次数为10000。
[0089]
s2：液滴图片二值化
[0090]
s21：设图像droplet_image，长度为image_height，图像的宽度为image_width，图片横坐标等分num_height份，纵坐标等分num_width份，则图像droplet_image总共分成num_size＝num_width*num_height张子图片，子图片的表示为sub_imager，其中 0≤r≤num_size，子图片高度为子图片宽度为子图片宽度为20≤num_width≤32，16≤num_height≤24，且num_height＜num_width。
[0091]
s22：任选一张子图片sub_imager，其中0≤r≤num_size。g(sub_x，sub_y)表示像素值在点(sub_x，sub_y)的灰度值，0≤sub_x≤sub_cols，0≤sub_y≤sub_rows。计算子图片阈值，生成此子图片的二值化图像：
[0092]
1)统计图片sub_imager灰度直方图：
[0093]
[0094][0095]
2)针对直方图进行归一化：
[0096]
forii 0:255
[0097][0098]
3)计算概率密度累加和：
[0099][0100]
4)从前往后遍历数组norm_histomr，找到第一个非零的数据位置设置为first_bin，从后往前遍历数组norm_histomr，找到第一个非零的数据位置设置为last_bin；
[0101]
5)计算阈值的位置点：
[0102]
[0103][0104]
6)生成sub_imager的二值化图像sub_binaryr：
[0105][0106]
7)依据sub_binaryr连通域sub_img_domainr对二值化图像进行修改：
[0107]
设sub_binaryr连通域的像素面积为sub_area_domainr，连通域的孔洞的数量为sub_img_holer，如果或者sub_img_holer＞10，则sub_binaryr的二值化图片sub_binaryr全部像素值为零。
[0108]
s23：通过上面的操作，num_size张子图片对应num_size张二值化图像，所有的二值化图像按原图片droplet_image的位置组合成的二值化image_binary。
[0109]
s3连通域处理
[0110]
s31：设定示二值化图像image_binary的第m个连通域，其中0≤m≤num
i_domain
，num
i_domain
表示连通域的个数。的像素面积为的像素面积为的像素面积为表示第m个连通域像素单元的个数。的像素周长为计算第m个连通域的圆率如果则置第m个连通域所包含位置的像素值为0。
[0111]
s32：连通域中像素集合表示为如果第m个像素单元属于则表示为其中设集合中像素单元横坐标的最大值表示为dmn_x_maxm，横坐标的最小值表示为dmn_x_minm，纵坐标的最大值表示为dmn_y_maxm，纵坐标的最小值表示为dmn_y_minm。
[0112]
s33：如果dmn_x_max
m-dmn_x_minm＞30，则计算连通域在x方向的投影值，其中最小投影值对应的x的坐标值为x_min_locm。基于坐标值x_min_locm沿y方向分割连通域连通域如图2所示。
[0113]
s34：如果dmn_y_maxm—dmn_y_minm＞30，则计算连通域在y方向的投影值，其中最小投影值对应的y的坐标值为y_min_locm。基于坐标值y_min_locm沿x方向分割连通域连通域
[0114]
s35：计算连通域的质心坐标
[0115][0116][0117]
计算连通域的像素单元到质心的最大距离：
[0118][0119][0120]
计算外接圆的面积：
[0121]
如果且b_curveradiom＜0.8，则在图像learning_image_binary的对应位置标记连通域
[0122]
s4判断假阳性，其中s42和s43是基于液滴像素值的标准差判断液滴假阳性；s44-s46 基于训练的分类器判断液滴假阳性；s47-s411基于液滴局部的像素对比度判断假阳性。
[0123]
s41：任选液滴dropletsi，其中0≤i≤n
droplets
，n
droplets
表示液滴的总数量。如果液滴is_droplets[i]＝1，液滴dropletsi位置(pos_xi，pos_yi)在二值化图像image_binary中对应的像素值为255，则更新is
droplets[i]
＝-1。如果is_droplets[i]＝1则表示真阳性液滴，is
droplets[i]
＝-1则表示假阳性液滴。
[0124]
s42：设液滴dropletsi位置(pos_xi，pos_yi)在二值化图像learning_image_binary中对应的像素值为255，则在灰度图片droplets_image中，以dropletsi的位置(pos_xi，pos_yi)为中心，以4个像素单元为半径，计算此区域内像素的平均值：
[0125]
其中gray(x_lox，y_loc)表示在图像droplets_image位置点(x_lox，y_loc)的像素值，(x_loc-pos_xi)2 (y_loc
‑ꢀ
pos_yi2≤16，numcirclei表示此圆形区域内的像素个数。
[0126]
s43：计算此区域的最大值计算此区域的最小值计算此区域的像素的标准差：
[0127]
如果则is_droplets[i]＝-1，如果则is_droplets[i]＝-1。
[0128]
s44：设表示二值化图像learning_image_binary的第1个连通域，其中0≤l≤num
c_domain
，num
c_domain
表示连通域的个数。
[0129]
计算的像素面积积表示第l个连通域像素单元的个数。连通域中像素集合表示为第1个像素单元属于表示为其中计算连通域的质心坐标
[0130][0131][0132]
计算连通域的像素单元到质心的最大距离：
[0133][0134]
其中计算此连通区域的圆度：如果液滴dropletsi位置点(pos_xi，pos_yi)属于连通域的集合且则is_droplets[i]＝-1；
[0135]
s45：任选一个液滴dropletsi，其中0≤i≤n
droplets
，n
droplets
表示液滴的总数量。液滴dropletsi与连通域的距离为：
[0136][0137]
如果则以(pos_xi，pos_yi)为中心，24个像素单元为矩形宽度，24个像素单元为矩形高度，截取灰度图片droplets_image，形成新的样本图像sample_imgi，sample_imgi表示第i个液滴的样本，样本图像sample_imgi的高度和宽度都为 24个像素。
[0138]
s46：载入训练生成的svm_class.xml文件，计算sample_imgi的hog特征，
[0139]
特征值导入到线性核函数svm中，如果判断结果为1，则is_droplets[i]＝1，否则为is_droplets[i]＝-1。
[0140]
s47：以连通域的质心坐标为中心，以domain_radius
l
为半径，其中统计此圆形区域内对应的液滴像素值的灰度直方图histom
l-gray
，直方图以灰度级为横坐标，像素个数为纵坐标。num_i_gray(g)表示此区域内液滴像素值为g的液滴的个数。统计此圆形区域内所有液滴的灰度值，最大灰度值为gray
max
，最小灰度值为gray
min
， 0≤gray
min
≤gray
max
≤255。0≤l≤num
c_domain
，num
c_domain
表示连通域总数量，con_ratio＝1.1。
[0141]
s48：此圆形区域内液滴最大灰度值和液滴最小灰度值分别为gray
max
和gray
min
，灰度值 g的概率密度为其中num表示此圆形区域内液滴的总个数，ng表示灰度值为g的像素个数，gray
min
≤g≤gray
max
；计算灰度概率密度直方图ck＝pg，其中gray
min
≤k≤gray
max
&&k＝g；
[0142]
计算灰度直方图一阶矩为：
[0143]
二阶矩为：
[0144]
三阶矩为：
[0145]
s49：计算参考值p_refer＝(z
1-m1
l
)/(z
1-z0)，其中：
[0146][0147][0148]
c1＝(-m3
l
m2
l
*m1
l
)/(m2
l-m1
l
*m1
l
)；
[0149]
c0＝(-m2
l
m3
l
*m1
l
)/(m2
l-m1
l
*m1
l
)；
[0150]
s410：灰度值从gray
min
到gray
max
依次遍历灰度概率密度直方图，如果存在ck＞p_refer，则histom
l-gray_th＝k，其中gray
min
≤k≤gray
max
。
[0151]
s411：以连通域的质心坐标为中心，以domain_radius
l
为半径，统计此区域内小于或等于histom
l-gray
_th液滴的最大像素值为其中gj是以连通域的质心坐标为中心，以domain_radius
l
为半径的圆形区域对应的液滴的像素值。如果存在液滴且且is_droplets[i]＝-1。
[0152]
以上仅以较佳实施例公布了本发明，然其并非用以限制本发明，凡采取等同替换或等效变换的方案所获得的技术方案，均落在本发明的保护范围内。