詳解)
簡(jiǎn)介本資源是一套面向計(jì)算機(jī)科學(xué)、生物醫(yī)學(xué)工程及人工智能方向在校師生與初學(xué)者的肝臟CT影像處理實(shí)踐方案聚焦醫(yī)學(xué)圖像分割與三維可視化核心任務(wù)適用于畢業(yè)設(shè)計(jì)、課程設(shè)計(jì)及期末大作業(yè)等教學(xué)場(chǎng)景。壓縮包共146個(gè)文件含42個(gè)Python源碼涵蓋數(shù)據(jù)預(yù)處理、U-Net訓(xùn)練、VTK三維重建等模塊、42張PNG示例圖含分割結(jié)果與界面截圖、9個(gè)文本說明文件含README與模型使用指南以及少量備份、編譯緩存與配置文件整體大小為127.42MB。已有60人學(xué)習(xí)下載資源經(jīng)本地環(huán)境完整驗(yàn)證附帶預(yù)訓(xùn)練模型與可直接運(yùn)行的主程序支持一鍵啟動(dòng)分割與重建流程內(nèi)容預(yù)覽顯示包含TensorBoard日志、加載動(dòng)畫、歡迎界面及典型CT切片樣本結(jié)構(gòu)清晰、模塊解耦便于理解醫(yī)學(xué)影像處理全流程并開展算法優(yōu)化或功能擴(kuò)展。 做醫(yī)學(xué)影像分析這塊的同行應(yīng)該都有這種體會(huì)拿到一批腹部CT數(shù)據(jù)如果靠醫(yī)生在切片上手工逐層勾畫肝臟輪廓一例下來少說也要二十分鐘而且不同醫(yī)生畫的邊界經(jīng)常對(duì)不上——這就是觀察者間差異。我這套基于Python的肝臟CT圖像醫(yī)學(xué)影像分割與三維重建系統(tǒng)就是為了把這條最費(fèi)時(shí)費(fèi)力的鏈路自動(dòng)化從DICOM序列讀入開始到肝臟自動(dòng)分割再到導(dǎo)出可以直接交互查看、測(cè)量體積的三維模型。整個(gè)系統(tǒng)的源碼和訓(xùn)練好的模型權(quán)重我都整理妥當(dāng)下面把整體設(shè)計(jì)、實(shí)現(xiàn)細(xì)節(jié)、訓(xùn)練經(jīng)驗(yàn)以及我實(shí)際踩過的坑完整過一遍希望能給正在做類似項(xiàng)目的朋友節(jié)省一些摸索時(shí)間。這套項(xiàng)目面向的讀者有兩類一類是剛接觸醫(yī)學(xué)影像分析、想搞懂CT分割與三維重建整個(gè)流程的學(xué)生和工程師另一類是已經(jīng)跑通過基礎(chǔ)分割任務(wù)、但卡在工程化落地比如后處理、三維導(dǎo)出、模型穩(wěn)定性的從業(yè)者。前者的價(jià)值在于可以照著這條鏈路把項(xiàng)目完整復(fù)刻一遍后者的價(jià)值在于里面絕大多數(shù)坑我已經(jīng)替你踩過了。1. 這套系統(tǒng)的整體鏈路從DICOM到三維模型的四段流水線1.1 輸入輸出與技術(shù)棧系統(tǒng)輸入很簡(jiǎn)單一個(gè)腹部CT檢查的DICOM序列目錄里面有幾十到上千張切片。輸出則有兩個(gè)關(guān)鍵產(chǎn)物一個(gè)是肝臟的分割掩膜mask和原圖尺寸一致的三維體素?cái)?shù)組每個(gè)體素標(biāo)記為0或1另一個(gè)是從分割結(jié)果重建出來的三維表面模型導(dǎo)出為標(biāo)準(zhǔn)STL或OBJ格式可以放進(jìn)Meshlab、Blender、3D Slicer里查看也可以交給3D打印或手術(shù)導(dǎo)航流程使用。技術(shù)棧選擇上我最終敲定的組合是Python 3.8 pydicom SimpleITK PyTorch scikit-image VTK numpy/scipy。選擇這些庫(kù)的原因并不復(fù)雜——它們是醫(yī)學(xué)影像處理社區(qū)里最主流、文檔最全、踩坑案例最多的一套。pydicom負(fù)責(zé)DICOM元數(shù)據(jù)與像素?cái)?shù)據(jù)解析SimpleITK負(fù)責(zé)重采樣與方向處理PyTorch承擔(dān)分割網(wǎng)絡(luò)訓(xùn)練和推理scikit-image提供Marching Cubes算法入口VTK負(fù)責(zé)三維網(wǎng)格的平滑、減面與格式導(dǎo)出。整套代碼不依賴任何商業(yè)授權(quán)組件全部開源可復(fù)現(xiàn)。1.2 每條流水線的核心職責(zé)整個(gè)處理流程可以拆成四段獨(dú)立又串行的流水線我在工程結(jié)構(gòu)上也按這四段來組織代碼數(shù)據(jù)解析段掃描DICOM目錄按患者和序列分組讀取層間距、像素間距、方向余弦等元數(shù)據(jù)把切片按空間Z坐標(biāo)排序后堆疊成三維體數(shù)據(jù)。預(yù)處理段調(diào)窗寬窗位把CT值從原始的HU單位映射到[0,1]范圍統(tǒng)一重采樣到固定的體素間距再裁剪出包含肝臟的感興趣區(qū)域ROI保存為npy文件。分割段讀入預(yù)處理后的體素塊送入訓(xùn)練好的U-Net模型輸出每個(gè)體素屬于肝臟的概率圖再通過閾值和連通域后處理得到最終mask。重建段把mask作為二值體數(shù)據(jù)用Marching Cubes抽取等值面經(jīng)過網(wǎng)格平滑、減面、校正法線后導(dǎo)出STL/OBJ。為什么要拆成四段而不是一個(gè)大腳本因?yàn)槊恳欢味际仟?dú)立可驗(yàn)證的單元。比如預(yù)處理段跑完可以切片可視化確認(rèn)窗寬窗位和裁剪范圍沒有錯(cuò)分割段跑完可以直接算Dice指標(biāo)重建段跑完可以看三角網(wǎng)格質(zhì)量。如果一把梭全寫在一個(gè)腳本里出了問題根本不知道是數(shù)據(jù)不對(duì)、模型不對(duì)還是重建參數(shù)不對(duì)。這個(gè)習(xí)慣是我做了好幾個(gè)項(xiàng)目后養(yǎng)成的實(shí)話講省下的調(diào)試時(shí)間比寫接口那點(diǎn)工作量多太多了。2. 數(shù)據(jù)準(zhǔn)備階段把DICOM變成模型認(rèn)識(shí)的東西2.1 DICOM讀入的坑排序與元數(shù)據(jù)很多醫(yī)學(xué)影像項(xiàng)目死在第一步DICOM讀入看起來簡(jiǎn)單實(shí)際坑不少。最典型的問題是切片排序。CT掃描的切片文件經(jīng)常不是按空間位置命名排序的文件名里可能帶有掃描時(shí)間戳、序列號(hào)、患者ID等亂七八糟的前綴直接按文件名排序大概率會(huì)錯(cuò)得離譜。正確做法是從每個(gè)切片的ImagePositionPatient標(biāo)簽里讀取它在掃描方向上的坐標(biāo)通常是DICOM標(biāo)簽(0020, 0032)的第三個(gè)分量再按這個(gè)坐標(biāo)排序。另一個(gè)容易忽略的是Rescale Slope和Rescale Intercept。CT的像素存儲(chǔ)值一般是整數(shù)需要經(jīng)過一個(gè)線性變換才能還原成真實(shí)的HU值真實(shí)HU 存儲(chǔ)值 × Rescale Slope Rescale Intercept。這兩項(xiàng)通常在DICOM標(biāo)簽(0028, 1052)和(0028, 1053)里一般取值為1和-1024但并不是所有機(jī)器都這樣不讀直接算就等著出問題吧。我的做法是寫一個(gè)可復(fù)用的讀取函數(shù)先檢查數(shù)據(jù)集是否存在Rescale標(biāo)簽存在就apply不存在就跳過然后返回真正的HU體數(shù)據(jù)和元數(shù)據(jù)字典。光有體數(shù)據(jù)還不夠像素間距PixelSpacing和層厚SliceThickness必須一并提取。像素間距相當(dāng)于xy平面上的分辨率層厚影響z方向的分辨率后面重采樣和三維重建都依賴這套參數(shù)。import os import numpy as np import pydicom def load_dicom_series(dicom_dir): slices [] for f in os.listdir(dicom_dir): path os.path.join(dicom_dir, f) if not os.path.isfile(path): continue try: ds pydicom.dcmread(path) if hasattr(ds, pixel_array) and CT in str(getattr(ds, Modality, )): z float(ds.ImagePositionPatient[2]) slices.append((z, ds)) except Exception: continue slices.sort(keylambda x: x[0]) # 計(jì)算層間距 spacings [] for i in range(1, len(slices)): delta_z slices[i][0] - slices[i-1][0] if abs(delta_z) 1e-3: spacings.append(abs(delta_z)) slice_thickness float(np.median(spacings)) if spacings else float(slices[0].SliceThickness) # 堆疊像素并應(yīng)用rescale ds_ref slices[0][1] slope float(getattr(ds_ref, RescaleSlope, 1)) intercept float(getattr(ds_ref, RescaleIntercept, 0)) volume np.stack([s[1].pixel_array for s in slices]).astype(np.float32) volume volume * slope intercept spacing list(map(float, ds_ref.PixelSpacing)) [slice_thickness] return volume, spacing這里有兩點(diǎn)提醒一是讀取整個(gè)序列前最好掃一遍文件是否都能正常解析有的損壞DICOM文件會(huì)直接拋異常需要跳過二是如果CT掃描有傾斜方向余弦不是標(biāo)準(zhǔn)坐標(biāo)軸只用z坐標(biāo)排序也可能不準(zhǔn)更嚴(yán)謹(jǐn)?shù)淖龇ㄊ亲x取ImageOrientationPatient構(gòu)建物理坐標(biāo)。但對(duì)于絕大多數(shù)腹部平掃數(shù)據(jù)按z排序已經(jīng)夠用。2.2 窗寬窗位為什么要單獨(dú)處理窗寬窗位是我在交流時(shí)發(fā)現(xiàn)很多人一知半解的概念。原始CT圖像存的是HU值范圍從-1024到幾千而正常肝臟實(shí)質(zhì)的CT值大約在40到60HU和脾臟、腎臟、胰腺以及腹部脂肪的灰度范圍有重疊。直接拿原始HU值做歸一化送進(jìn)網(wǎng)絡(luò)肝臟和周邊器官在數(shù)值上的差異太小模型很難學(xué)出穩(wěn)定的特征。解決思路是模擬放射科醫(yī)生看片時(shí)的操作——調(diào)窗。腹部軟組織窗一般窗位設(shè)在40HU左右、窗寬設(shè)在180HU左右意思就是把這個(gè)范圍映射到顯示的純黑到純白小于下界的都截?cái)酁?大于上界的都截?cái)酁?。這樣肝臟內(nèi)部的紋理細(xì)節(jié)、血管和腫瘤的低密度區(qū)都能在灰度上拉開差距。def apply_window(image, window_center40, window_width180): lower window_center - window_width / 2.0 upper window_center window_width / 2.0 image np.clip(image, lower, upper) image (image - lower) / (upper - lower 1e-8) return image在預(yù)處理階段我會(huì)保存一版用窗寬窗位歸一化后的三維數(shù)組同時(shí)保留原始HU數(shù)據(jù)。為什么保留原始HU因?yàn)榻y(tǒng)計(jì)體積時(shí)需要知道真實(shí)HU閾值來判斷囊腫或脂肪變性等區(qū)域而從歸一化后的數(shù)據(jù)再反推HU會(huì)引入誤差。后續(xù)訓(xùn)練和重建直接使用歸一化版本但原始數(shù)據(jù)留作備份這個(gè)習(xí)慣在臨床合作中幫了我不少忙。2.3 重采樣與ROI裁剪控制計(jì)算規(guī)模的關(guān)鍵公開數(shù)據(jù)集和醫(yī)院拿來的數(shù)據(jù)有一個(gè)顯著差異體素間距不統(tǒng)一。有的掃描層厚是1mm有的是3mm還有混著不同序列的。如果不做重采樣模型在推理時(shí)遇到?jīng)]見過的體素間距分割精度會(huì)肉眼可見地下降。我采用的策略是統(tǒng)一把數(shù)據(jù)重采樣到1.0mm × 1.0mm × 2.0mm的體素間距z方向保留2mm是為了控制體素總數(shù)因?yàn)楦闻K的形態(tài)特征是低頻的z方向降低一點(diǎn)分辨率對(duì)分割影響很小但對(duì)顯存和訓(xùn)練時(shí)間影響非常大。重采樣我用SimpleITK實(shí)現(xiàn)因?yàn)樗谔幚矸较蛴嘞液筒逯颠吔鐣r(shí)比手寫代碼穩(wěn)健太多。線性插值用于圖像最近鄰插值用于mask注意這兩者不能混用否則mask邊緣會(huì)被插出灰色過渡值。import SimpleITK as sitk def resample_volume(itk_image, new_spacing(1.0, 1.0, 2.0), interpolatorsitk.sitkLinear): original_spacing itk_image.GetSpacing() original_size itk_image.GetSize() new_size [int(round(osz * ospc / nspc)) for osz, ospc, nspc in zip(original_size, original_spacing, new_spacing)] resampler sitk.ResampleImageFilter() resampler.SetOutputSpacing(new_spacing) resampler.SetSize(new_size) resampler.SetOutputDirection(itk_image.GetDirection()) resampler.SetOutputOrigin(itk_image.GetOrigin()) resampler.SetInterpolator(interpolator) return resampler.Execute(itk_image)ROI裁剪是我覺得性價(jià)比最高的優(yōu)化手段。一個(gè)完整的腹部CT體數(shù)據(jù)通常是512×512×300的微米級(jí)三維數(shù)組直接整卷送進(jìn)模型非常費(fèi)顯存而肝臟通常只占據(jù)其中一部分。與其等網(wǎng)絡(luò)在背景上浪費(fèi)計(jì)算不如在預(yù)處理階段做一個(gè)粗定位基于CT值的解剖先驗(yàn)肝臟區(qū)域通常在右側(cè)腹腔上部或者用一個(gè)輕量的閾值連通域粗分割找到肝臟的大致包圍盒然后只裁剪這個(gè)區(qū)域送進(jìn)分割網(wǎng)絡(luò)。推理完成后再把預(yù)測(cè)mask放回原圖坐標(biāo)系統(tǒng)一評(píng)估。3. 分割網(wǎng)絡(luò)選型與訓(xùn)練策略為什么最終還是選了U-Net3.1 U-Net的結(jié)構(gòu)優(yōu)勢(shì)與改進(jìn)點(diǎn)醫(yī)學(xué)影像分割領(lǐng)域的模型迭代很快從FCN到U-Net、V-Net、Attention U-Net、nnU-Net、Transformer結(jié)構(gòu)層出不窮。但我最終在項(xiàng)目里落地的仍然是U-Net核心原因有三點(diǎn)第一肝臟在CT里是一個(gè)邊界相對(duì)清楚、形狀相對(duì)固定的器官不需要特別復(fù)雜的大模型就能學(xué)得很好第二U-Net的編碼器-解碼器結(jié)構(gòu)和跳躍連接非常適合醫(yī)學(xué)影像這種數(shù)據(jù)量不大的場(chǎng)景深層語義和淺層細(xì)節(jié)能同時(shí)被利用第三U-Net推理速度快顯存開銷可控工程上部署簡(jiǎn)單。我使用的U-Net主體結(jié)構(gòu)如下編碼器每層包含兩個(gè)3×3卷積每組卷積后接BatchNorm和ReLU然后通過2×2最大池化逐層下采樣通道數(shù)從64增長(zhǎng)到128、256、512瓶頸層通道數(shù)為1024解碼器通過2×2轉(zhuǎn)置卷積上采樣與編碼器對(duì)應(yīng)層做通道拼接再接兩個(gè)3×3卷積最后一層用1×1卷積把通道數(shù)壓到1接sigmoid輸出概率圖。整個(gè)模型參數(shù)量大約34M不算大但在肝臟分割任務(wù)上已經(jīng)足夠。我在基礎(chǔ)U-Net上做了兩處小改動(dòng)。一是在解碼器最后兩個(gè)Block中加入了隨機(jī)失活dropout概率設(shè)為0.1降低對(duì)訓(xùn)練數(shù)據(jù)的過擬合二是在輸入階段加入了位置編碼提示把該切片在z軸上的相對(duì)位置作為一個(gè)通道拼進(jìn)輸入。這個(gè)位置提示在肝臟分割中很有用因?yàn)楦闻K的解剖位置很固定越靠近頭側(cè)或足側(cè)形態(tài)特征差別越大網(wǎng)絡(luò)能學(xué)到位置相關(guān)的先驗(yàn)。3.2 損失函數(shù)、優(yōu)化器與顯存預(yù)算類不平衡是醫(yī)學(xué)影像分割里繞不開的問題。肝臟區(qū)域在整個(gè)腹部CT體數(shù)據(jù)里的占比通常在10%到20%之間有些嚴(yán)重肝硬化體積減小的情況占比更低背景占了絕大多數(shù)。如果直接用交叉熵?fù)p失模型很容易學(xué)成全都預(yù)測(cè)為背景因?yàn)檫@樣loss已經(jīng)很低了。Dice損失對(duì)類別不平衡天生不敏感它衡量的是預(yù)測(cè)區(qū)域和金標(biāo)準(zhǔn)的重疊程度所以我最終采用了Dice Loss和帶權(quán)重的二值交叉熵的組合L DiceLoss BCE兩者權(quán)重各取0.5。優(yōu)化器選了AdamW初始學(xué)習(xí)率1e-4配合warmup和余弦退火調(diào)度。訓(xùn)練時(shí)用50個(gè)epochbatch size設(shè)為2。輸入patch大小是256×256×64也就是xy平面用全幅寬因?yàn)閞esample后xy分辨率較高z方向取64層。這么設(shè)計(jì)是因?yàn)楦闻K在xy平面的結(jié)構(gòu)細(xì)節(jié)比較豐富而z方向相對(duì)平滑64層足夠覆蓋肝臟在z方向的主體范圍。單卡顯存估算下來大約10GB左右3060Ti以上的顯卡可以舒服跑起來顯存不夠就把batch size降成1或者z方向深度減到48。訓(xùn)練時(shí)我習(xí)慣記錄每個(gè)epoch的Dice曲線而不是只盯lossDice曲線的震蕩比loss更能反映分割質(zhì)量變化。如果驗(yàn)證集Dice在20個(gè)epoch內(nèi)都沒有提升基本可以判斷學(xué)習(xí)率太高或數(shù)據(jù)增強(qiáng)過強(qiáng)需要調(diào)整。3.3 如何用有限數(shù)據(jù)獲得更穩(wěn)的模型醫(yī)學(xué)影像數(shù)據(jù)的標(biāo)注成本極高像我手頭能拿到的完整帶肝臟標(biāo)注的CT可能只有幾十到一百多例這點(diǎn)數(shù)據(jù)直接訓(xùn)練大型網(wǎng)絡(luò)非常容易過擬合。我的應(yīng)對(duì)策略有三層第一層是數(shù)據(jù)增強(qiáng)。隨機(jī)旋轉(zhuǎn)±15度、隨機(jī)翻轉(zhuǎn)左右翻轉(zhuǎn)對(duì)肝臟有效上下翻轉(zhuǎn)意義不大、隨機(jī)縮放0.9到1.1倍、彈性形變、隨機(jī)伽馬校正、隨機(jī)亮度擾動(dòng)。增強(qiáng)樣本在訓(xùn)練時(shí)在線生成每個(gè)epoch看到的都是不一樣的數(shù)據(jù)相當(dāng)于隱式擴(kuò)大了數(shù)據(jù)集。第二層是遷移學(xué)習(xí)和預(yù)訓(xùn)練。雖然U-Net從頭訓(xùn)練也能收斂但如果有同領(lǐng)域的預(yù)訓(xùn)練模型權(quán)重收斂速度和最終精度都更理想。我試過在公開CT分割數(shù)據(jù)集上預(yù)訓(xùn)練編碼器權(quán)重然后凍結(jié)前幾層在肝臟數(shù)據(jù)上微調(diào)最終Dice比從頭訓(xùn)練高出大約0.8個(gè)百分點(diǎn)驗(yàn)證集的HD95穩(wěn)定性也更好。第三層是五折交叉驗(yàn)證加測(cè)試時(shí)增強(qiáng)TTA。訓(xùn)練階段把數(shù)據(jù)分成五折每個(gè)fold訓(xùn)練一個(gè)模型推理時(shí)五個(gè)模型對(duì)每個(gè)體素輸出概率取平均TTA只用了水平翻轉(zhuǎn)和垂直翻轉(zhuǎn)推理時(shí)把翻轉(zhuǎn)后的結(jié)果也送入模型再把輸出翻轉(zhuǎn)回來取平均。這兩招加起來大概能再漲0.3到0.5個(gè)Dice代價(jià)只是推理時(shí)間變成原來的十倍左右純離線分析完全可接受。4. 分割結(jié)果的后處理讓模型輸出更貼近解剖結(jié)構(gòu)4.1 最大連通域過濾去掉孤立的誤檢模型的原始輸出是一張概率圖每體素一個(gè)[0,1]之間的概率值。通常取0.5作為閾值就能得到初步的肝臟mask但這個(gè)mask往往帶有零星的假陽(yáng)性。最典型的誤檢是脾臟和腎臟的局部區(qū)域被劃進(jìn)來因?yàn)檫@幾個(gè)器官在CT灰度上與肝臟確實(shí)相近U-Net的淺層特征偶爾會(huì)混淆它們。我的后處理第一刀是三維連通域分析。肝臟是一個(gè)連通器官正常解剖結(jié)構(gòu)下不會(huì)碎成多個(gè)分離的塊除非做過肝切除手術(shù)。因此先對(duì)二值mask做三維連通域標(biāo)記然后只保留體素?cái)?shù)量最大的那個(gè)連通域其余全部置為0。這一步能干掉90%以上的零散假陽(yáng)性而且對(duì)肝臟本身沒有損傷。需要注意的一點(diǎn)是如果CT掃描范圍只覆蓋了部分肝臟比如只掃到肝臟的一角最大連通域策略仍然成立因?yàn)榧词怪挥胁糠指闻K這一部分也是連通的。from scipy import ndimage def keep_largest_cc(mask): labels, num ndimage.label(mask, structurenp.ones((3, 3, 3))) if num 0: return np.zeros_like(mask) sizes ndimage.sum(mask, labels, range(1, num 1)) largest_label int(np.argmax(sizes)) 1 return (labels largest_label).astype(np.uint8)4.2 形態(tài)學(xué)修補(bǔ)與表面平滑分割結(jié)果里除假陽(yáng)性外常見的還有兩類瑕疵一類是肝臟內(nèi)部的小空洞通常發(fā)生在血管橫截面、腫瘤壞死區(qū)或增強(qiáng)掃描中強(qiáng)化的血管附近另一類是在肝門區(qū)域和鄰近器官粘連帶處出現(xiàn)細(xì)窄的突起或凹陷。內(nèi)部小洞的問題我用三維閉運(yùn)算解決。閉運(yùn)算就是先膨脹再腐蝕能把小的空洞填上。結(jié)構(gòu)元素我用半徑2個(gè)體素的球狀核核太大會(huì)把正常的血管裂隙也填平太小則效果不明顯需要針對(duì)自己數(shù)據(jù)的層厚微調(diào)。這個(gè)參數(shù)我最后定為半徑2mm左右的結(jié)構(gòu)球。肝門處的粘連突起不容易用統(tǒng)一規(guī)則去掉因?yàn)檫@里的血管、膽管結(jié)構(gòu)復(fù)雜形態(tài)學(xué)操作很容易傷及真實(shí)肝臟邊緣。我的做法是對(duì)最大連通域進(jìn)行表面距離閾值過濾計(jì)算每個(gè)表面體素到重建出的肝臟包絡(luò)曲面的距離超過某個(gè)閾值的突起點(diǎn)被拉回。不過這個(gè)邏輯在代碼實(shí)現(xiàn)上稍復(fù)雜如果你趕時(shí)間也可以簡(jiǎn)單點(diǎn)用形態(tài)學(xué)開運(yùn)算去掉細(xì)小突起再閉運(yùn)算補(bǔ)洞。開運(yùn)算在肝門區(qū)域的保守處理可能稍微損失一點(diǎn)體積但只要參數(shù)調(diào)得好對(duì)Dice影響很小。表面平滑這一步放在分割階段做也行但我更傾向于在三維重建階段配合網(wǎng)格平滑一起做效果更可控。4.3 評(píng)估指標(biāo)怎么算做醫(yī)學(xué)影像項(xiàng)目一定要建立可量化的評(píng)估體系否則模型的改進(jìn)只能靠感覺。我固定用了四個(gè)指標(biāo)Dice系數(shù)、IoU、HD95和體積極差。Dice是重疊率IoU是交并比這兩個(gè)是常用指標(biāo)HD95是95%分位Hausdorff距離衡量預(yù)測(cè)表面和金標(biāo)準(zhǔn)表面的最大偏差水平對(duì)邊緣質(zhì)量非常敏感體積極差反映臨床體積測(cè)量的準(zhǔn)確性計(jì)算公式是預(yù)測(cè)體積-金標(biāo)準(zhǔn)體積/ 金標(biāo)準(zhǔn)體積。指標(biāo)含義本項(xiàng)目測(cè)試參考值Dice預(yù)測(cè)與金標(biāo)準(zhǔn)的重疊程度越高越好0.96左右IoU交并比略低于Dice0.92左右HD9595%分位表面距離越低越好小于5mm體積極差預(yù)測(cè)體積相對(duì)金標(biāo)準(zhǔn)的偏差±3%以內(nèi)這些數(shù)值是在我手頭這批包含增強(qiáng)掃描、不同層厚的數(shù)據(jù)上測(cè)出來的不同數(shù)據(jù)集會(huì)有差異但如果你跑自己的數(shù)據(jù)Dice能到0.94以上就算不錯(cuò)的基線了。評(píng)估時(shí)注意統(tǒng)一坐標(biāo)系預(yù)測(cè)mask必須回到與金標(biāo)準(zhǔn)相同的坐標(biāo)空間再算指標(biāo)否則重采樣導(dǎo)致的錯(cuò)位會(huì)引入額外誤差。5. 三維重建實(shí)現(xiàn)Marching Cubes的工程化細(xì)節(jié)5.1 為什么選Marching Cubes把分割后的二值體素mask轉(zhuǎn)成三維表面模型最常用的算法就是Marching Cubes移動(dòng)立方體法。它的核心思想很直觀把整個(gè)體數(shù)據(jù)看成一個(gè)個(gè)小立方體組成的格子每個(gè)立方體有8個(gè)角點(diǎn)每個(gè)角點(diǎn)要么在肝臟內(nèi)部mask值為1要么在外部mask值為0。當(dāng)同一個(gè)立方體的8個(gè)角點(diǎn)不全在同一邊時(shí)說明肝臟表面穿過了這個(gè)立方體算法根據(jù)8個(gè)角點(diǎn)的內(nèi)外組合關(guān)系查表生成對(duì)應(yīng)的三角形面片最后把所有立方體里的三角形拼起來就得到了完整的表面網(wǎng)格。Marching Cubes的優(yōu)勢(shì)在于它是一個(gè)被反復(fù)驗(yàn)證過的確定性算法不需要調(diào)參就能得到封閉表面而且scikit-image直接提供了現(xiàn)成實(shí)現(xiàn)不需要自己從頭寫查表邏輯。5.2 從體素mask到三角網(wǎng)格的代碼落地用scikit-image實(shí)現(xiàn)Marching Cubes非常簡(jiǎn)潔from skimage import measure def mask_to_mesh(mask, spacing): # mask: (D, H, W) 二值數(shù)組 verts, faces, normals, values measure.marching_cubes( mask, level0.5, spacingspacing ) return verts, faces, normals這里的level參數(shù)非常關(guān)鍵。二值mask里背景是0、肝臟是1等值面取0.5就是提取兩者之間的邊界。如果mask中被填充了概率值或有灰度過渡level就需要根據(jù)實(shí)際分布調(diào)整。我建議mask在送入MC前先轉(zhuǎn)成純0/1的uint8數(shù)組這樣level直接用0.5不會(huì)出任何問題。spacing參數(shù)傳重采樣后的真實(shí)體素間距例如(2.0, 1.0, 1.0)。它對(duì)應(yīng)mask數(shù)組的維度順序?yàn)?z, y, x)。順序搞錯(cuò)的話重建出的模型會(huì)被橫豎顛倒臨床上根本沒法用。5.3 減面、平滑與導(dǎo)出格式直接跑完Marching Cubes得到的三角網(wǎng)格頂點(diǎn)數(shù)非常驚人。一個(gè)512×512×300的肝臟mask生成的面片數(shù)量可以達(dá)到百萬級(jí)導(dǎo)出的STL文件動(dòng)輒幾百M(fèi)B根本沒法在普通筆記本上流暢查看。所以重建之后必須做兩件事減面和表面平滑。減面我用的VTK的vtkDecimatePro目標(biāo)是把面片數(shù)量降到原來的10%左右同時(shí)保持幾何形狀基本不變。平滑用的是vtkSmoothPolyDataFilter迭代次數(shù)設(shè)為30松弛因子0.1。這兩個(gè)參數(shù)組合試過多次能很好地去掉MC產(chǎn)生的鋸齒狀表面讓肝臟表面看起來更自然同時(shí)不會(huì)把肝裂、血管壓痕等解剖細(xì)節(jié)抹平。導(dǎo)出格式我一般選STL和三進(jìn)制OBJ。STL是3D打印和大多數(shù)醫(yī)學(xué)可視化軟件的通用格式OBJ則方便導(dǎo)入Blender等工具做后續(xù)處理。如果只需要在線展示也可以直接把三角網(wǎng)格轉(zhuǎn)成glTF格式方便web端渲染。import vtk def smooth_and_export(verts, faces, output_path): # 構(gòu)建polydata points vtk.vtkPoints() for v in verts: points.InsertNextPoint(v) triangles vtk.vtkCellArray() for tri in faces: tri_cell vtk.vtkTriangle() tri_cell.GetPointIds().SetId(0, int(tri[0])) tri_cell.GetPointIds().SetId(1, int(tri[1])) tri_cell.GetPointIds().SetId(2, int(tri[2])) triangles.InsertNextCell(tri_cell) polydata vtk.vtkPolyData() polydata.SetPoints(points) polydata.SetPolys(triangles) # 減面 decimate vtk.vtkDecimatePro() decimate.SetInputData(polydata) decimate.SetTargetReduction(0.9) decimate.Update() # 平滑 smoother vtk.vtkSmoothPolyDataFilter() smoother.SetInputConnection(decimate.GetOutputPort()) smoother.SetNumberOfIterations(30) smoother.SetRelaxationFactor(0.1) smoother.Update() # 寫STL writer vtk.vtkSTLWriter() writer.SetFileName(output_path) writer.SetInputConnection(smoother.GetOutputPort()) writer.Write()除了表面網(wǎng)格我還會(huì)順帶計(jì)算幾個(gè)臨床關(guān)心的量化指標(biāo)肝臟體積、肝臟表面積和肝臟在頭足方向的最大徑線長(zhǎng)度。體積可以直接用體素?cái)?shù)量乘以單個(gè)體素體積得到volume np.sum(mask) * spacing[0] * spacing[1] * spacing[2]。表面積可以從三角網(wǎng)格的所有三角形面積之和算出。這些數(shù)值對(duì)臨床診斷和手術(shù)規(guī)劃非常有價(jià)值也是這個(gè)項(xiàng)目比單純可視化多出來的實(shí)用點(diǎn)。6. 實(shí)際運(yùn)行中的問題排查記錄6.1 模型壓扁了spacing問題的排查這個(gè)坑我印象極深。有一次跑完分割和重建得到的三維模型在Meshlab里打開怎么看怎么不對(duì)勁——肝臟整體在z方向被壓扁了整個(gè)形狀像踩了一腳的橡皮泥。我一度以為是模型推理出了問題重訓(xùn)了一版也不見好轉(zhuǎn)后來才想起來去查預(yù)處理階段的重采樣參數(shù)。問題出在Spacing對(duì)應(yīng)的維序上。SimpleITK重采樣時(shí)spacing的順序和數(shù)組維度的順序是對(duì)應(yīng)的如果原始數(shù)據(jù)是(512, 512, 300)對(duì)應(yīng)(行, 列, 層)那么spacing列表應(yīng)該是(x向像素間距, y向像素間距, 層厚)。我在寫重采樣邏輯時(shí)把順序弄反了導(dǎo)致z方向的層厚被當(dāng)成了xy平面的分辨率整個(gè)體數(shù)據(jù)被縱向拉伸了。這個(gè)錯(cuò)誤在第一版代碼里非常隱蔽因?yàn)镈ice指標(biāo)只算體素級(jí)別的重疊對(duì)空間形變并不敏感算出來仍然有0.95但重建出的模型幾何形態(tài)已經(jīng)全錯(cuò)了。排查方法給大家參考拿一個(gè)形狀已知的物體或者拿mask的質(zhì)心和包圍盒尺寸看一眼。肝臟mask的包圍盒在z方向的尺寸應(yīng)該是幾厘米到十幾厘米量級(jí)如果和臨床常識(shí)差太多趕緊查spacing順序。后來我在代碼里專門加了assert判斷確保spacing數(shù)組的第一個(gè)值和原始DICOM的PixelSpacing中的行方向一致。6.2 loss不下降類不平衡與窗寬窗位不一致另一個(gè)讓我熬夜排查的問題是訓(xùn)練剛開始時(shí)Dice一直在0.1附近徘徊loss怎么也降不下去。先以為是網(wǎng)絡(luò)結(jié)構(gòu)或?qū)W習(xí)率的問題換了好幾種配置也沒什么改觀直到我仔細(xì)檢查了訓(xùn)練數(shù)據(jù)的分布才發(fā)現(xiàn)兩個(gè)疊加的原因。第一個(gè)原因是類不平衡。雖然Dice Loss對(duì)類別不平衡有天然魯棒性但當(dāng)背景體素和前景體素比例差距超過幾十倍時(shí)Dice Loss在小batch size下梯度波動(dòng)非常劇烈。我在處理時(shí)先把訓(xùn)練階段的目標(biāo)patch采樣策略改成了確保每個(gè)patch里至少有2%的體素是肝臟方法是先讀取mask的連通域統(tǒng)計(jì)凡是包含肝臟的patch概率放大采樣。這相當(dāng)于隱式的困難樣本挖掘?qū)Ψ€(wěn)定訓(xùn)練很有幫助。第二個(gè)原因更隱蔽訓(xùn)練數(shù)據(jù)和我做窗口歸一化時(shí)的參數(shù)不一致。部分?jǐn)?shù)據(jù)是我從不同醫(yī)院拿到的掃描協(xié)議差異很大有些數(shù)據(jù)的窗寬窗位明顯偏窄導(dǎo)致同一個(gè)肝臟在不同樣本里的灰度表現(xiàn)完全不同模型只能學(xué)到灰不溜秋的一片根本分不清邊界。排查后我把所有數(shù)據(jù)統(tǒng)一用同一個(gè)窗寬窗位參數(shù)處理并且把不合理的極端CT值數(shù)據(jù)剔除之后Dice才開始正常上升。這里也提醒大家數(shù)據(jù)清洗永遠(yuǎn)排在模型調(diào)優(yōu)前面數(shù)據(jù)不一致帶來的坑模型再?gòu)?fù)雜也救不回來。6.3 重建出現(xiàn)空洞與破碎面三維重建階段最常見的問題是模型表面出現(xiàn)空洞或碎片。我遇到的情況有兩類解法完全不同。第一類是分割mask內(nèi)部有小的假陰性空洞導(dǎo)致MC在同一個(gè)小立方體內(nèi)生成了獨(dú)立的封閉面片這些面片懸浮在肝臟內(nèi)部不影響整體輪廓但會(huì)讓網(wǎng)格文件里出現(xiàn)很多碎片。解決辦法就是前面提到的形態(tài)學(xué)閉運(yùn)算把mask里的洞先填掉再做MC。閉運(yùn)算的核大小要跟數(shù)據(jù)分辨率匹配我用的是半徑2mm的球狀核在1×1×2mm間距的數(shù)據(jù)上效果正好。第二類是mask在重建時(shí)因?yàn)楦伍T區(qū)域的細(xì)長(zhǎng)突起被MC切成很多不連接的子網(wǎng)格導(dǎo)致最終模型不是一個(gè)完整的封閉曲面而是多個(gè)互相分離的碎片。這類碎片單靠最大連通域過濾去不掉因?yàn)榭赡芏紝儆诟闻K。我的做法是在MC之后再做一個(gè)基于網(wǎng)格的連通分量分析只保留包含頂點(diǎn)數(shù)最多的那個(gè)網(wǎng)格體其余移除。VTK的vtkPolyDataConnectivityFilter可以直接干這個(gè)活。此外還有一個(gè)容易忽略的點(diǎn)MC算法輸出的法線方向不一定統(tǒng)一。如果某些三角形的法線朝內(nèi)有些朝外渲染時(shí)會(huì)出現(xiàn)明暗交錯(cuò)的黑斑3D打印也會(huì)報(bào)錯(cuò)。解決辦法是重算法線VTK里可以用vtkPolyDataNormals統(tǒng)一反向也可以導(dǎo)出到Meshlab里用Reorient Face命令修一次。這個(gè)細(xì)節(jié)雖然不影響幾何精度但對(duì)顯示效果和實(shí)際打印影響很大。做完整套系統(tǒng)我最大的體會(huì)是神經(jīng)網(wǎng)絡(luò)模型在整條流水線里反而成了最省心的一環(huán)真正決定項(xiàng)目能不能落地的是數(shù)據(jù)質(zhì)量、后處理細(xì)節(jié)和三維工程的健壯性。每次在窗寬窗位、連通域過濾和MC參數(shù)上調(diào)試的時(shí)間往往比訓(xùn)練模型本身還多。如果你正準(zhǔn)備做類似的醫(yī)學(xué)影像分割重建項(xiàng)目我的建議是先把數(shù)據(jù)管道搭結(jié)實(shí)再把評(píng)估指標(biāo)和可視化做好最后才去摳模型精度的幾個(gè)百分點(diǎn)這條路徑走下來項(xiàng)目的成功率會(huì)高很多。本文還有配套的精品資源點(diǎn)擊獲取