CN116523936B - 肝脏图像的脉管分割方法、系统及存储介质 - Google Patents

Abstract

本发明涉及医学图像分割技术领域，提供一种肝脏图像的脉管分割方法、系统及存储介质，其中的方法包括：S1、通过图像采集获得包含所有脉管信息的图像；S2、通过图像预处理将脉管腔和肝实质细胞变为高亮，肝叶外部变为黑色背景；S3、先将图像数据集切割为三维数据块，分批次读入内存中，再采用基于区域增长的自动分割算法对脉管进行语义分割，获得连续完整的脉管；S4、采用脉管自动追踪算法对不同类型的脉管进行追踪，完成对脉管的分类。本发明通过并行计算和动态读取数据的方法解决了TB数据处理过程中存在内存不足和运算速度慢的问题，能够在完整肝叶和亚微米分辨率尺度下，自动分割不同类型的脉管，不会丢失肝血窦信号。

Description

肝脏图像的脉管分割方法、系统及存储介质

技术领域

本发明涉及医学图像分割技术领域，特别涉及一种基于并行计算和动态读取数据的肝脏脉管分割方法。

背景技术

肝脏中包含不同类型的脉管，具体包括肝静脉、门静脉、胆管、淋巴管和肝血窦。通过对肝血窦内皮细胞进行标记，并使用显微光学切片断层成像技术，可以在完整肝叶和亚微米分辨率尺度下获取包含所有脉管信息的图像数据。该图像数据包含6.0cm×7.2cm×9.0cm的完整肝叶数据，其中体素分辨率为0.3×0.3×1μm³。原始采集的图像为无符号十六位格式(uint16)，存储这些图像数据大约需要9TB空间。由于这些海量图像数据不能一次性读入到内存中，因此，现有的图像处理算法需要对图像数据进行降采样，从而读入到内存进行后续的图像重建。使用上述方法进行图像重建能够完整获得肝叶中的肝静脉、门静脉、胆管、淋巴管和肝动脉。由于对图像数据进行了降采样处理，导致上述方法存在如需两个问题：

(1)虽然能够获得完整肝叶内的肝静脉、门静脉、胆管、淋巴管和肝动脉图像，但是由于原始数据远远超过内存容量，上述操作是对降采样后的数据进行分割的。因此不能在亚微米分辨率下获得这些脉管信号，降低了脉管的分辨率；

(2)由于肝血窦的直径为3-5微米，降采样后的图像超过肝血窦的直径，就会导致完整肝叶中肝血窦信号的丢失，因此不能重建完整肝叶内的肝血窦结构。

发明内容

本发明的目的是提出一种肝脏图像的脉管分割方法、系统及存储介质，以解决现有的图像处理算法对图像数据进行降采样处理时无法在亚微米分辨率尺度下获得脉管信号以及丢失肝血窦信号的问题。

为实现上述目的，本发明采用以下具体技术方案：

一方面，本发明提供一种基于并行计算和动态读取数据的肝脏图像的脉管分割方法，包括如下步骤：

S1、图像采集：对肝血窦内皮细胞进行标记，并使用显微光学切片断层成像技术进行成像，获得包含所有脉管信息的图像，构成图像数据集；其中，图像中所有脉管的脉管壁和肝实质细胞的细胞膜为前景信号，脉管腔和肝实质细胞的内部为背景信号；

S2、图像预处理：先对图像数据集中的图像进行强度翻转，使脉管腔和肝实质细胞的内部变为前景信号，脉管壁和肝实质细胞的细胞膜变为背景信号；再使用阈值分割算法对肝叶进行分割，获取肝叶区域图像；对肝叶区域图像和强度反转后的图像进行逻辑与操作，保持脉管腔和肝实质细胞为高亮，肝叶外部为黑色背景；

S3、脉管语义分割：先将图像数据集切割为三维数据块，分批次读入内存中；再采用基于区域增长的自动分割算法对脉管进行分割，从而获得连续完整的脉管；其中，脉管类型包括肝静脉、门静脉、淋巴管、胆管、肝动脉和肝血窦；

S4、脉管追踪：采用脉管自动追踪算法对不同类型的脉管进行追踪，以在完整肝叶和亚微米分辨率尺度下完成对脉管的分类。

优选地，采用基于区域增长的自动分割算法对脉管进行分割的具体过程如下：

S31、在三维数据块中等间隔自动设置种子点；

S32、对每个种子点使用区域增长算法进行三维数据块的分割，获得三维数据块中完整连续的脉管和肝实质细胞；

S33、通过欧几里得方法计算脉管和肝实质细胞外轮廓，将脉管和肝实质细胞的外轮廓均向各自的中心缩小，使肝实质细胞和脉管分离，各自成为独立的连通域，而所有相连的脉管仍为一个独立的连通域；

S34、使用体积阈值法去除独立的肝实质细胞；

S35、再次通过欧几里得方法计算肝血管的外层轮廓，然后将肝血管的外层轮廓叠加到现有血管中，使脉管向外扩展到原始尺寸，并将所有的三维数据块进行组合，实现在完整肝叶在亚分辨率尺度下的脉管分割。

优选地，采用脉管自动追踪算法对不同类型的脉管进行追踪的具体过程如下：

S41、使用半径为R1的圆形结构体对脉管语义分割后的图像进行腐蚀获得结果V1，R1的腐蚀半径大于肝血窦的脉管直径，小于淋巴管、胆管和肝动脉的直径，保证只去除肝血窦，保留完整的胆管、淋巴管和肝动脉；

S42、使用半径为R2对圆形结构体对脉管语义分割后的图像进行腐蚀获得结果V2，R2的腐蚀半径大于淋巴管、肝动脉、门静脉的脉管直径，保证只去除门静脉和肝静脉的静脉末端分支，保留门静脉主干和肝静脉主干；

S43、将结果V1和结果V2进行异或操作，获得结果V3，结果V3包含有淋巴管、胆管、肝动脉、门静脉与肝静脉的静脉末端分支；

S44、对结果V3进行降采样获得结果V3_ds，将结果V3_ds读入内存中；

S45、找到结果V3_ds中独立的连通域，确定胆管、淋巴管和肝动脉的种子点，并对种子点进行升采样；

S46、基于胆管、淋巴管和肝动脉的种子点，通过并行动态算法实现在完整肝叶内和亚微米分辨率尺度下对胆管、淋巴管和肝动脉的追踪；

S47、对结果V2进行降采样获得结果V2_ds。

S48、找到结果V2_ds中独立的连通域，确定肝静脉主干和门静脉主干的种子点，并对种子点进行升采样；

S49、基于肝静脉主干和门静脉主干的种子点，通过并行动态算法实现在完整肝叶内和亚微米分辨率尺度下对肝静脉主干和门静脉主干的追踪；

S410、对肝静脉主干和肝静脉的静脉末端分支进行逻辑或操作，通过并行动态算法实现在完整肝叶内和亚微米分辨率尺度下对肝静脉的追踪；

S411、对门静脉主干和门静脉的静脉末端分支进行逻辑或操作，通过并行动态算法实现在完整肝叶内和亚微米分辨率尺度下对门静脉的追踪。

优选地，并行动态算法的处理过程如下：

S461、根据胆管、淋巴管和肝动脉的种子点的三维坐标，在硬盘中分别寻找对应于胆管、淋巴管和肝动脉的种子点的三维数据块，获得三维数据块集合；

S462、将三维数据块集合读入到内存中，使用三维区域生长算法对三维数据块集合中每个三维数据块的脉管语义分割结果进行追踪，同时判断当前三维数据块是否为终端块，如果是则结束追踪，如果否则记录待追踪的三维数据块的编号和对应的种子点；

S463、判断待追踪的三维数据块是否发生追踪冲突，如果发生追踪冲突，则所有的种子点进行顺次追踪，如果未发生追踪冲突，则所有的种子进行并行追踪；

S464、最终在肝叶内和亚微米分辨率尺度下追踪得到肝动脉、淋巴管、胆管、门静脉和肝静脉的静脉末端分支。

优选地，在步骤S4之后，还包括肝小叶分割的步骤，具体如下：

S51、对肝静脉进行降采样；

S52、通过骨架提取算法计算肝静脉的骨架；

S53、根据肝静脉的骨架利用Kalman Palagyi算法计算肝静脉的管径；

S54、根据肝静脉的管径对肝静脉的静脉末端分支进行筛选，获得降采样的肝静脉的静脉末端分支；

S55、通过骨架提取算法计算肝静脉的静脉末端分支的骨架，并进行升采样还原到原始尺寸；

S56、根据肝静脉的静脉末端分支的骨架坐标寻找对应的完整肝叶的位置，并对肝静脉的静脉末端分支周围的数据进行切割，获得包含原始分辨率下的肝静脉的静脉末端分支和肝血窦的数据块；

S57、采用3D腐蚀方法对步骤S56获得的数据块进行操作，获得原始分辨率下的肝静脉的静脉末端分支和肝血窦；

S58、通过骨架提取算法计算肝血窦的骨架，并将肝血窦的骨架转变为图结构；其中图结构中的节点分成终端节点和连接节点，终端节点为图结构的末端节点，而连接节点则为图结构中的中间节点；

S59、将肝静脉的静脉末端分支和图结构中的终端节点进行逻辑与操作，找到与肝静脉的静脉末端分支直接相连的肝血窦；

S510、将与终端节点相连的中间节点记为向量<a，b>，通过余弦定理计算向量<a，b>与肝静脉的静脉末端分支对应的向量<c，d>之间的角度；

S511、当计算出的角度落在角度阈值范围内时，则认为该肝血窦属于肝静脉的静脉末端分支，从而继续重复步骤S511进行追踪，当计算出的角度未落在角度阈值范围内时，则停止追踪。

优选地，在执行步骤S511之后，还包括对肝血窦进行剪枝操作，具体如下：

计算所追踪到的肝血窦的末端相距肝静脉的静脉末端分支的距离，如果该距离超出距离阈值，则删除该终端节点，同时将与该终端节点相连的中间节点设为终端节点，依次进行判断，直到肝血窦的末端相距肝静脉的静脉末端分支的距离未超出距离阈值。

优选地，在步骤S3中，将图像数据集切割成38000个的三维数据块，每个三维数据块的体素分辨率为500×500×500voxel。

另一方面，本发明提供一种计算机存储介质，计算机存储介质中存储有计算机程序，该计算机程序被处理器执行时实现上述的脉管分割方法。

再一方面，本发明提供一种三维细胞图像分割系统，包括存储器、处理器及存储在该存储器上且在处理器上运动的计算机程序，处理器在执行该计算机程序时实现上述的脉管分割方法。

与现有技术相比，本发明能够在完整肝叶和亚微米分辨率尺度下，自动分割不同类型的脉管(包括肝动脉、肝静脉、淋巴管、胆管、门静脉和肝血窦)，不会丢失肝血窦信号。此外本发明还可以自动分割肝小叶。得益于原始精细的成像结果和所采取的基于区域增长的自动分割算法，所分割的脉管结构在完整肝叶和亚微米分辨率尺寸下是连续和完整的。在完整肝叶和亚微米分辨率尺度下，使用基于区域增长的自动分割算法在一台服务器上(CPU为48核96线程，内存为500G)处理一套原始大小为9TB(20000×26000×9000体素)的数据时，原始结构信息数据的预处理时间为需要7.6小时，平均速度为1.7×10⁸voxel/s；语义分割需要7.4小时，平均速度为1.7×10⁸voxel/s；脉管追踪需要20个小时，平均速度为6.5×10⁷voxel/s。图像预处理、脉管语义分割和脉管追踪采用并行化计算，处理器的核心越多，内存越大则运算时间越少。同时，根据本算法动态读取数据的特点，所处理数据量的大小和计算时间为线性关系，可以按照上述平均速度估算所需处理时间。

附图说明

图1是根据本发明提供实施例的基于并行计算和动态读取数据的肝脏图像的脉管分割方法的流程示意图。

图2是根据本发明提供实施例的原始图像和预处理后的图像；

图2中，(A)示出了原始图像中的一个切片；(B)示出了经过预处理后的局部图像；(C)示出了获取的肝叶边缘。

图3是根据本发明提供实施例的脉管语义分割流程示意图；

图3中，(A)示出了对完整肝叶的切块；(B)示出了对每个三维数据块批量设置的种子点；(C)示出了使用区域增长算法分割得到的结果；(D)示出了相连的肝实质细胞和脉管相连；(E)示出了断开的肝实质细胞和脉管；(F)示出了删除后的肝实质细胞；(G)示出了还原后的肝血窦；(H)示出了语义分割后得到的脉管；(I)示出了对于完整肝叶的语义分割结果。

图4是根据本发明实施例提供的脉管追踪和追踪结果示意图；

图4中，(A)示出了完整肝叶内脉管追踪的流程；(B)示出了并行动态算法的流程；(C)示出了肝动脉、胆管和淋巴管的追踪结果；(D)示出了肝静脉和门静脉的追踪结果。

图5是根据本发明实施例提供的肝小叶分割流程的示意图；

图5中，(A)示出了提取的肝静脉的骨架；(B)根据肝静脉的管径对肝静脉的分类结果；(C)示出了获取的肝静脉的静脉末端分支；(D)示出了根据门静脉的骨架在完整肝叶中找到的对应区域；(E)示出了将对应区域内的肝血窦转变的图结构；(F)示出了两个向量的夹角；(G)示出了筛选出的肝血窦；(H)示出了对肝血窦的剪枝操作；(I)示出了实际分割的肝小叶结果。

具体实施方式

在下文中，将参考附图描述本发明的实施例。在下面的描述中，相同的模块使用相同的附图标记表示。在相同的附图标记的情况下，它们的名称和功能也相同。因此，将不重复其详细描述。

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及具体实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本发明，而不构成对本发明的限制。

本发明第一方面提供了一种基于并行计算和动态读取数据的肝脏图像的脉管分割方法。

图1示出了根据本发明提供实施例的基于并行计算和动态读取数据的肝脏图像的脉管分割方法的流程。

如图1所示，本发明提供实施例提供的基于并行计算和动态读取数据的肝脏图像的脉管分割方法，包括如下步骤：

S1、图像采集：对肝血窦内皮细胞进行标记，并使用显微光学切片断层成像技术进行成像，获得包含所有脉管信息的图像，构成图像数据集；其中，图像中所有脉管的脉管壁和肝实质细胞的细胞膜为前景信号，脉管腔和肝实质细胞的内部为背景信号。

本发明使用CD31抗体标记肝血窦内皮细胞，因此采集的图像数据中所有脉管的脉管壁和肝实质细胞的细胞膜为前景信号，脉管腔和肝实质细胞内部为背景信号，如图2中的(A)所示。

S2、图像预处理：先对图像数据集中的图像进行强度翻转，使脉管腔和肝实质细胞的内部变为前景信号，脉管壁和肝实质细胞的细胞膜变为背景信号；再使用阈值分割算法对肝叶进行分割，获取肝叶区域图像；对肝叶区域图像和强度反转后的图像进行逻辑与操作，保持脉管腔和肝实质细胞为高亮，肝叶外部为黑色背景。

预处理实际包括两步：第一步为图像翻转，将图像中的脉管腔和肝实质细胞变为前景信号，将图像中的脉管壁和肝实质细胞的细胞膜变为背景信号，如图2中的(B)所示。第二步，由于图像翻转会将边缘的背景同样变为高亮，因此使用阈值分割方法对肝叶进行分割获取肝叶区域图像，如图2中的(C)所示。在获取肝叶区域图像后，对肝叶区域图像和强度反转后的图像进行逻辑与操作，保持肝叶外为黑色背景。经过预处理操作后，脉管腔和肝实质细胞变为高亮。这一步的图像处理需要4个小时。

S3、脉管语义分割：先将图像数据集切割为三维数据块，分批次读入内存中；再采用基于区域增长的自动分割算法对脉管进行分割，从而获得连续完整的脉管；其中，脉管类型包括肝静脉、门静脉、淋巴管、胆管、肝动脉和肝血窦。

为了在完整肝叶和亚微米分辨率尺度下获得完整连续的肝脏脉管系统，本发明首先将整个数据集切割为38000块，每个块为500×500×500voxel的三维数据块，从而可以分批次读入内存中，如图3中的(A)所示；再采用区域增长方法进行脉管分割，从而获得连续完整的脉管。区域增长算法需要手动指定一个起始种子点，从而开始进行分割。对于整个数据分割需要手动设定3万多个种子点，这需要极大的工作量。为了实现肝脏脉管自动分割算法，本发明根据肝脏脉管系统结构，研发出基于区域增长的自动分割算法。

该自动分割算法共包括五步：

第一步、在三维数据块中等间隔自动设置种子点，如图3中的(B)所示；

第二步、对每个种子点使用区域增长算法对数据块进行分割，如图3中的(C)所示。

因为肝实质细胞的中间也为高亮，因此经过第一步和第二步，可以获得三维数据块中完整连续的脉管和肝实质细胞。肝实质细胞包含有两类，第一类为独立的肝实质细胞，第二类为与肝血窦相连的肝实质细胞。通过对批量种子点进行追踪，能够获得肝血窦和肝实质细胞，当肝实质细胞通过细胞膜和肝血窦相连时，这种连接关系，可以抽象为一条直线和一个圆形相切，而肝血窦和肝实质细胞之间的接触面积小于肝血窦内部的接触面积，如图3中的(D)所示。

第三步、通过欧几里得方法计算脉管和肝实质细胞外轮廓，将脉管和肝实质细胞的外轮廓均向各自的中心缩小一圈，如图3中的(E)所示。由于肝实质细胞和脉管的连接方式具有相切的特点，这样在相连的两个区域都向中心方向缩小后，实现肝实质细胞和脉管的分离，而脉管之间的连接不会被断开。断开之后肝实质细胞与脉管分离成为独立的连通域，而所有相连的脉管仍然为一个独立的连通域。

第四步、使用体积阈值法去除独立的肝实质细胞，如图3中的(F)所示。

计算所有独立连通域的体积。由于肝实质细胞的体积较小，而相连脉管的体积较大。因此设置一个体积阈值，当连通域的体积小于该体积阈值则认为是肝实质细胞，删除该连通域，从而去除独立的肝实质细胞。

第五步、再次通过欧几里得方法计算肝血管的外层轮廓，然后将肝血管的外层轮廓叠加到现有血管中，使脉管向外扩展，恢复到缩小前的尺寸，如图3中的(G)和(H)所示。由于每个三维数据块之间都有密集的种子点，分块操作并不影响最后的结果。对每个三维数据块进行自动分割后，将所有的三维数据块进行组合，实现在完整肝叶在亚微米分辨率尺度下的脉管分割，如图3中的(I)所示。

由于上述自动分割算法完全自动化，块与块之间相互独立，因此使用并行计算方法加快速度计算。在48核96线程服务器上，完整语义分割结果的运行时间需要7h。

S4、脉管追踪：采用脉管自动追踪算法对不同类型的脉管进行追踪，以在完整肝叶和亚微米分辨率尺度下完成对脉管的分类。

在实现对脉管的语义分割之后，得到肝静脉、门静脉，淋巴管、胆管、肝动脉、肝血窦的分割结果，肝静脉和门静脉的静脉末端尺寸较小，接近于胆管。

在完整肝叶和亚微米分辨率尺度下进行脉管分类的前提是脉管内部必须是连通的，不同脉管之间是独立的、不连通的。然而，淋巴管、胆管和肝动脉常伴随着门静脉，因此淋巴管、胆管和肝动脉很可能接触在一起。此外，门静脉和肝静脉分支可能会通过静脉末端接触。这些接触会导致脉管的分割出现错误。而在TB级的海量数据下，通过手工找到脉管相连的地方，然后手工找到接触位置，从而将其手工断开，这个工作量是十分巨大的。

为了在完整肝叶和亚微米分辨率尺度下区分不同类型的脉管，本发明研发一种脉管自动追踪算法，如图4中的(A)所示。由于肝脏中肝静脉和门静脉的管径较大，淋巴管、胆管和肝动脉的管径较小，肝血窦的直径最小。因此对于肝脏中的脉管追踪本发明采用分而治之的方法，具体如下：

Step1、使用半径为R1的圆形结构体对脉管语义分割后的图像进行腐蚀，获得结果V1，R1的腐蚀半径大于肝血窦的脉管直径，小于淋巴管、胆管和肝动脉的直径，从而保证只腐蚀掉肝血窦，保留完整的胆管、淋巴管和肝动脉。

Step2、使用半径R2对脉管语义分割后的图像进行腐蚀获得结果V2，该半径大于淋巴管、肝动脉、门静脉的脉管直径，从而获得门静脉和肝静脉主干，而门静脉和肝静脉的静脉末端分支则会被腐蚀掉。

Step3、将V1结果和V2结果进行异或操作，实现V1结果与V2结果的相减，获得结果V3，结果V3包含有淋巴管、胆管、肝动脉、门静脉和肝静脉的静脉末端分支。

门静脉和肝静脉是肝脏的主脉管。其中，肝静脉的静脉末端分支叫做中央静脉，门静脉的静脉末端分支没有特定的名称。

Step4、对结果V3进行降采样，获得结果V3_ds，将结果V3_ds读入内存中。

Step5、找到结果V3_ds中独立的连通域，确定胆管、淋巴管和肝动脉的种子点，并对胆管、淋巴管和肝动脉的种子点进行升采样。

Step6、利用升采样后的胆管、淋巴管和肝动脉的种子点，采用并行动态算法实现在完整肝叶内和亚微米分辨率尺度下的脉管追踪，如图4中的(B)所示。

并行动态算法的原理如下：

Step6.1、根据胆管、淋巴管和肝动脉的种子点的三维坐标，在硬盘中分别寻找对应于胆管、淋巴管和肝动脉的种子点的三维数据块，获得三维数据块集合。

例如根据肝动脉种子点寻找与肝动脉所对应的三维数据块，根据淋巴管种子点寻找与淋巴管所对应的三维数据块。

Step6.2、将三维数据块集合读入到内存中，使用三维区域生长算法对三维数据块集合中每个三维数据块的脉管语义分割结果进行追踪，同时判断当前三维数据块是否为终端块，如果是则结束追踪，如果否则记录待追踪的三维数据块的编号和对应的种子点。

Step6.3、当完成该三维数据块的追踪之后，将结果保存在硬盘中。同时判断该三维数据块是否为终端块；如果是终端块则结束追踪，如果不是终端块，则记录待追踪的三维数据块的编号和对应的种子点。因为给定的种子点和块编号之后，使用区域生长方法依然是块与块之间独立的，本发明通过并行化算法进行并行追踪从而加快计算速度。

Step6.4、判断待追踪的三维数据块是否发生追踪冲突，如果发生追踪冲突，则所有的种子点进行顺次追踪，如果未发生追踪冲突，则所有的种子点进行并行追踪。

由于追踪过程是并行的，当多个种子点同时追踪同一个三维数据块时，无法同步不同种子点的追踪结果，因此会发生冲突。因此，首先判断待追踪的三维数据块是否发生冲突，如果发生冲突，先选取其中一个种子点返回步骤S461进行追踪，其他种子点等待下一轮再进行追踪；如果不冲突，则所有的种子点直接则返回步骤S461进行并行追踪；

经过并行动态算法之后，本发明在完整肝叶和亚微米分辨率下追踪得到了肝动脉、淋巴管、肝动脉和静脉末端分支。

Step7、对结果V2进行降采样获得结果V2_ds。

Step8、找到结果V2_ds中独立的连通域，确定肝静脉和门静脉主干的种子点，并将肝静脉和门静脉主干的种子点进行升采样。

Step9、使用并行动态算法在完整肝叶和亚微米分辨率尺度上追踪肝静脉主干HV_main和门静脉主干PV_main。

Step10、将肝静脉主干HV_main和肝静脉的中央静脉HV_CV进行逻辑或操作，然后通过并行动态算法在完整肝叶内和亚微米分辨率尺度下追踪肝静脉HV。

Step11、将门静脉主干PV_main和门静脉的静脉末端分支PV_CV进行逻辑或操作，然后通过并行动态算法在完整肝叶内和亚微米分辨率尺度下追踪门静脉PV。

经过上述步骤，本发明在完整肝叶和亚微米分辨率尺度下，实现了对肝静脉、门静脉，淋巴管、胆管和肝动脉的追踪。

本发明通过Step1、Step2、Step3，解决了在完整肝叶内，当肝动脉、淋巴管和胆管与门静脉相连时，不能找到连接点的难题。本发明通过Step4和Step5对海量数据降采样，通过在降采样数据中寻找脉管中的种子点，在经过升采样，将种子点变为完整肝叶下的坐标，从而实现种子点设定，解决了在完整肝叶和亚微米分辨率尺度下进行寻找特定脉管种子点的难题。通过Step6的并行动态算法，能够在数据不完全读入内存中的情况下对脉管进行追踪，解决了内存不足和运算速度慢的问题。通过Step10和Step11解决了当肝静脉和门静脉在末端分支相连时，不能区分肝静脉和门静脉的难题。经过上述步骤最终实现了在完整肝叶和亚微米分辨率尺度下的对不同脉管的追踪，如图4中的(C)和(D)所示。

肝小叶是肝脏组织的最小结构单元。小鼠肝小叶的典型结构为中央一根中央静脉，围绕其周围具有其他许多不同的静脉分支，根据肝小叶的解剖学结构小鼠中的肝小叶组织结构为肝血窦围绕中央静脉进行发散性分布。本发明在获取精细的肝血窦结构和门静脉之后，可以对肝小叶进行分割。

对肝小叶进行分割的步骤如下：

第一步、对肝静脉进行降采样(2×2×2voxel)。

第二步、通过骨架提取算法计算肝静脉的骨架，如图5中的(A)所示。

骨架提取算法为现有技术，参考文献为《Building skeleton models via 3-Dmedial surface/axis thinning algorithms》。

第三步、根据肝静脉的骨架利用Kalman Palagyi算法计算肝静脉的管径，如图5中的(B)所示。

Kalman Palagyi算法将血管分支看成一个圆柱体，对应的管径计算公式如下：

从而计算出脉管的管径。

第四步、根据肝静脉的管径对肝静脉的中央静脉进行筛选，获得降采样的肝静脉的中央静脉，如图5中的(C)所示。

第五步、通过骨架提取算法计算肝静脉的中央静脉的骨架，并进行升采样还原到原始尺寸。

第六步、根据肝静脉的中央静脉的骨架坐标寻找对应的完整肝叶的位置，并对肝静脉的中央静脉的周围进行切割，获得包含原始分辨率下的肝静脉的静脉末端分支和肝血窦的数据块，如图5中的(D)所示。

第七步、采用3D腐蚀方法对步骤S56获得的数据块进行操作，获得原始分辨率下的肝静脉的静脉末端分支和肝血窦。

切割出来的数据块中包含有原始分辨率下的中央静脉和肝血窦。设置一个半径为R的结构体，R大于肝血窦的直径，小于中央静脉直径，然后采用3D腐蚀操作获得原始分辨率下的肝静脉的中央静脉和肝血窦；

第八步、通过骨架提取算法计算肝血窦的骨架，并将肝血窦的骨架转变为图结构；其中，图结构中的节点分为终端节点和连接节点，如图5中的(E)所示，终端节点为图结构的末端节点，而连接节点则为图结构中的中间节点。此时本发明已经获得中央静脉形状、中央静脉对应的中央静脉主轴和对应的肝血窦图结构。

第九步、将肝静脉的中央静脉和图结构中的终端节点进行逻辑与操作，找到与肝静脉的中央静脉直接相连的肝血窦。

第十步、将与肝血窦的终端节点相连的中间节点记为向量<a，b>，然后通过余弦定理计算该向量<a，b>与中央静脉对应的向量<c，d>之间的夹角，如图5中的(F)所示。

第十一步，当计算出的角度落在角度阈值范围内时，则认为该肝血窦属于该中央静脉，从而继续重复步骤S511进行追踪，当计算出的角度未落在角度阈值范围内时，则停止追踪，如图5中的(G)所示。

在上述步骤中，肝小叶分割结果的好坏与所设定的角度相关，而阈值设定高则肝小叶的分割区域则会变小，而阈值设定较低，则肝小叶分割区域则会变大。为了解决阈值角度设定的严格依赖，然后本发明增加了剪枝操作，通过检测肝血窦的末端相距中央静脉的距离，对超过距离阈值的终端节点进行剪枝操作，如图5中的(H)所示。如果肝血窦的末端相距中央静脉的距离超过了距离阈值，则删除该终端节点，同时将与该终端节点相连的中间节点设为终端节点，依次进行判断，直到所有肝血窦的末端相距中央静脉的距离都保持在一个合理范围内，如图5中的(I)所示。

在上述算法中，由于每个肝小叶的分割依然是各自独立进行的，因此本发明采取并行操作的方法进行加快算法的运行。

本发明通过并行计算和动态读取数据的方法解决了TB数据处理过程中存在内存不足和运算速度慢的问题，能够在完整肝叶和亚微米分辨率尺度下，自动分割不同类型的脉管，包括肝动脉、肝静脉、淋巴管、胆管、门静脉和肝血窦。此外本发明还可以自动分割肝小叶。得益于原始精细的成像结果和所采取的基于区域增长算法，所分割的脉管结构在完整肝叶上和亚分辨率上是连续和完整的。在完整肝叶和亚微米分辨率下，使用上述算法在一台服务器上(CPU为48核96线程，内存为500G)处理一套原始大小为9TB(20000×26000×9000体素)的数据时，原始结构信息数据的预处理时间为需要7.6小时，平均速度为1.7×10^⁸voxel/s；语义分割需要7.4小时,平均速度为1.7×10^⁸voxel/s；脉管追踪需要20个小时，平均速度为6.5×10^⁷voxel/s。预处理、语义分割和脉管追踪是并行化计算的，理论上处理器的核心越多，内存越大则运算时间越少。同时，根据本算法动态读取数据的特点，所处理数据量的大小和计算时间为线性关系，可以按照上述平均速度估算所需处理时间。

本发明第二方面提供了一种计算机存储介质，计算机存储介质中存储有计算机程序，计算机程序被处理器执行时实现如本申请第一方面的方法步骤。

存储介质为存储器，存储器可以为非易失性存储介质，示例性地可以包括但不限于只读存储器(Read-Only Memory，ROM)、可编程只读存储器(Programmable ROM，PROM)、可擦除可编程只读存储器(Erasable PROM，EPROM)、电可擦除可编程只读存储器(Electrically EPROM，EEPROM)或闪存(Flash Memory)，例如可以是以下任一种：嵌入式多媒体卡(Embedded Multi Media Card，EMMC)、Nor Flash、Nand Flash等。

示例性地，存储器还可以包括缓存装置，用于缓存数据，例如信号队列。缓存装置可以为易失性存储介质，示例性地可以包括但不限于随机存取存储器(Random AccessMemory，RAM)，静态随机存取存储器(Static RAM，SRAM)、动态随机存取存储器(DynamicRAM，DRAM)、同步动态随机存取存储器(Synchronous DRAM，SDRAM)、双倍数据速率同步动态随机存取存储器(Double Data Rate SDRAM，DDR SDRAM)、DDR2、DDR3、增强型同步动态随机存取存储器(Enhanced SDRAM，ESDRAM)、同步连接动态随机存取存储器(Synchlink DRAM，SLDRAM)和直接内存总线随机存取存储器(Direct Rambus RAM，DRAM)等。

本发明第三方面提供了一种基于并行计算和动态读取数据的肝脏图像的脉管分割系统，脉管分割系统包括存储器、处理器及存储在存储器上且在处理器上运行的计算机程序，处理器执行计算机程序时实现本申请第一方面方法的步骤。

示例性地，存储器例如可以包括智能电话的存储卡、平板电脑的存储部件、个人计算机的硬盘、只读存储器(ROM)、可擦除可编程只读存储器(EPROM)、便携式紧致盘只读存储器(CD-ROM)、USB存储器、或者上述存储介质的任意组合。计算机可读存储介质可以是一个或多个计算机可读存储介质的任意组合。

示例性地，处理器可以是中央处理单元(CPU)、图像处理单元(GPU)、专用集成电路(Application Specific Integrated Circuit，ASIC)、现场可编程门阵列(FieldProgrammable Gate Array，FPGA)或者具有数据处理能力和/或指令执行能力的其它形式的处理单元，并且可以控制系统中的其它组件以执行期望的功能。例如，处理器可以包括一个或多个嵌入式处理器、处理器核心、微型处理器、逻辑电路、硬件有限状态机(FiniteState Machine，FSM)、数字信号处理器(Digital Signal Processing，DSP)或它们的组合。

在本说明书的描述中，参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不必须针对的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任一个或多个实施例或示例中以合适的方式结合。此外，在不相互矛盾的情况下，本领域的技术人员可以将本说明书中描述的不同实施例或示例以及不同实施例或示例的特征进行结合和组合。

尽管上面已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本发明的限制，本领域的普通技术人员在本发明的范围内可以对上述实施例进行变化、修改、替换和变型。

以上本发明的具体实施方式，并不构成对本发明保护范围的限定。任何根据本发明的技术构思所作出的各种其他相应的改变与变形，均应包含在本发明权利要求的保护范围内。

Claims (8)

1.一种基于并行计算和动态读取数据的肝脏图像的脉管分割方法，其特征在于，包括如下步骤：

S1、图像采集：对肝血窦内皮细胞进行标记，并使用显微光学切片断层成像技术进行成像，获得包含所有脉管信息的图像，构成图像数据集；其中，图像中所有脉管的脉管壁和肝实质细胞的细胞膜为前景信号，脉管腔和肝实质细胞的内部为背景信号；

S2、图像预处理：先对图像数据集中的图像进行强度翻转，使脉管腔和肝实质细胞的内部变为前景信号，脉管壁和肝实质细胞的细胞膜变为背景信号；再使用阈值分割算法对肝叶进行分割，获取肝叶区域图像；对肝叶区域图像和强度反转后的图像进行逻辑与操作，保持脉管腔和肝实质细胞为高亮，肝叶外部为黑色背景；

S3、脉管语义分割：先将图像数据集切割为三维数据块，分批次读入内存中；再采用基于区域增长的自动分割算法对脉管进行分割，从而获得连续完整的脉管；其中，脉管类型包括肝静脉、门静脉、淋巴管、胆管、肝动脉和肝血窦；

S4、脉管追踪：采用脉管自动追踪算法对不同类型的脉管进行追踪，以在完整肝叶和亚微米分辨率尺度下完成对脉管的分类；采用脉管自动追踪算法对不同类型的脉管进行追踪的具体过程如下：

S41、使用半径为R1的圆形结构体对脉管语义分割后的图像进行腐蚀获得结果V1，R1的腐蚀半径大于肝血窦的脉管直径，小于淋巴管、胆管和肝动脉的直径，保证只去除肝血窦，保留完整的胆管、淋巴管和肝动脉；

S42、使用半径为R2对圆形结构体对脉管语义分割后的图像进行腐蚀获得结果V2，R2的腐蚀半径大于淋巴管、肝动脉、门静脉的脉管直径，保证只去除门静脉和肝静脉的静脉末端分支，保留门静脉主干和肝静脉主干；

S43、将结果V1和结果V2进行异或操作，获得结果V3，结果V3包含有淋巴管、胆管、肝动脉、门静脉与肝静脉的静脉末端分支；

S44、对结果V3进行降采样获得结果V3_ds，将结果V3_ds读入内存中；

S45、找到结果V3_ds中独立的连通域，确定胆管、淋巴管和肝动脉的种子点，并对各种子点进行升采样；

S46、基于胆管、淋巴管和肝动脉的种子点，通过并行动态算法实现在完整肝叶内和亚微米分辨率尺度下对胆管、淋巴管和肝动脉的追踪；

S47、对结果V2进行降采样获得结果V2_ds；

S48、找到结果V2_ds中独立的连通域，确定肝静脉主干和门静脉主干的种子点，并对种子点进行升采样；

S49、基于肝静脉主干和门静脉主干的种子点，通过并行动态算法实现在完整肝叶内和亚微米分辨率尺度下对肝静脉主干和门静脉主干的追踪；

S410、对肝静脉主干和肝静脉的静脉末端分支进行逻辑或操作，通过并行动态算法实现在完整肝叶内和亚微米分辨率尺度下对肝静脉的追踪；

S411、对门静脉主干和门静脉的静脉末端分支进行逻辑或操作，通过并行动态算法实现在完整肝叶内和亚微米分辨率尺度下对门静脉的追踪。

2.如权利要求1所述的脉管分割方法，其特征在于，采用基于区域增长的自动分割算法对脉管进行分割的具体过程如下：

S31、在三维数据块中等间隔自动设置种子点；

S32、对每个种子点使用区域增长算法进行三维数据块的分割，获得三维数据块中完整连续的脉管和肝实质细胞；

S33、通过欧几里得方法计算脉管和肝实质细胞外轮廓，将脉管和肝实质细胞的外轮廓均向各自的中心缩小，使肝实质细胞和脉管分离，各自成为独立的连通域，而所有相连的脉管仍为一个独立的连通域；

S34、使用体积阈值法去除独立的肝实质细胞；

S35、再次通过欧几里得方法计算肝血管的外层轮廓，然后将肝血管的外层轮廓叠加到现有血管中，使脉管向外扩展到原始尺寸，并将所有的三维数据块进行组合，实现在完整肝叶在亚分辨率尺度下的脉管分割。

3.如权利要求1所述的脉管分割方法，其特征在于，并行动态算法的处理过程如下：

S461、根据胆管、淋巴管和肝动脉的种子点的三维坐标，在硬盘中分别寻找对应于胆管、淋巴管和肝动脉的种子点的三维数据块，获得三维数据块集合；

S462、将三维数据块集合读入到内存中，使用三维区域生长算法对三维数据块集合中每个三维数据块的脉管语义分割结果进行追踪，同时判断当前三维数据块是否为终端块，如果是则结束追踪，如果否则记录待追踪的三维数据块的编号和对应的种子点；

S463、判断待追踪的三维数据块是否发生追踪冲突，如果发生追踪冲突，则所有的种子点进行顺次追踪，如果未发生追踪冲突，则所有的种子点进行并行追踪；

S464、最终在肝叶内和亚微米分辨率尺度下追踪得到肝动脉、淋巴管、胆管、门静脉和肝静脉的静脉末端分支。

4.如权利要求1~3中任一项所述的脉管分割方法，其特征在于，在步骤S4之后，还包括肝小叶分割的步骤，具体如下：

S51、对肝静脉进行降采样；

S52、通过骨架提取算法计算肝静脉的骨架；

S53、根据肝静脉的骨架利用Kalman Palagyi算法计算肝静脉的管径；

S54、根据肝静脉的管径对肝静脉的静脉末端分支进行筛选，获得降采样的肝静脉的静脉末端分支；

S55、通过骨架提取算法计算肝静脉的静脉末端分支的骨架，并进行升采样还原到原始尺寸；

S56、根据肝静脉的静脉末端分支的骨架坐标寻找对应的完整肝叶的位置，并对肝静脉的静脉末端分支的周围进行切割，获得包含原始分辨率下的肝静脉的静脉末端分支和肝血窦的数据块；

S57、采用3D腐蚀方法对步骤S56获得的数据块进行操作，获得原始分辨率下的肝静脉的静脉末端分支和肝血窦；

S58、通过骨架提取算法计算肝血窦的骨架，并将肝血窦的骨架转变为图结构；其中图结构中的节点分成终端节点和连接节点，终端节点为图结构的末端节点，而连接节点则为图结构中的中间节点；

S59、将肝静脉的静脉末端分支和图结构中的终端节点进行逻辑与操作，找到与肝静脉的静脉末端分支直接相连的肝血窦；

S510、将与终端节点相连的中间节点记为向量<a，b>，通过余弦定理计算向量<a，b>与肝静脉的静脉末端分支对应的向量<c，d>之间的角度；

S511、当计算出的角度落在角度阈值范围内时，则认为该肝血窦属于肝静脉的静脉末端分支，从而继续重复步骤S511进行追踪，当计算出的角度未落在角度阈值范围内时，则停止追踪。

5.如权利要求4所述的脉管分割方法，其特征在于，在执行步骤S511之后，还包括对肝血窦进行剪枝操作，具体如下：

计算所追踪到的肝血窦的末端相距肝静脉的静脉末端分支的距离，如果该距离超出距离阈值，则删除该终端节点，同时将与该终端节点相连的中间节点设为终端节点，依次进行判断，直到肝血窦的末端相距肝静脉的静脉末端分支的距离未超出距离阈值。

6.如权利要求1所述的脉管分割方法，其特征在于，在步骤S3中，将图像数据集切割成38000个的三维数据块，每个三维数据块的体素分辨率为500×500×500voxel。

7.一种计算机存储介质，其特征在于，所述计算机存储介质中存储有计算机程序，该计算机程序被处理器执行时实现权利要求1~6中任一项所述的方法。

8.一种三维细胞图像分割系统，其特征在于，包括存储器、处理器及存储在该存储器上且在处理器上运动的计算机程序，处理器在执行该计算机程序时实现权利要求1~6中任一项所述的方法。