在做3D分割任務(wù)中，多數(shù)的方法多采用整體縮放，或裁剪成一個(gè)個(gè)小的patch操作，這樣做的一個(gè)主要原因是內(nèi)存問題。還有就是有些目標(biāo)太小，比如分割結(jié)節(jié)，用整圖直接輸入網(wǎng)絡(luò)，正負(fù)樣本的不均衡是非常大的。

相較于整體縮放，采用裁剪成patch的方法，對(duì)于小目標(biāo)會(huì)更加的魯棒，這也是大多數(shù)3D分割任務(wù)中常選取的方式。尤其是針對(duì)醫(yī)學(xué)影像的器官分割任務(wù)，CT結(jié)節(jié)診斷等等，對(duì)于細(xì)節(jié)的要求是非常高的。采用縮小的方式，反而會(huì)使得目標(biāo)的像素區(qū)域在輸入階段，就損失較多。

本文，就針對(duì)2D、3D的圖像和MR數(shù)據(jù)進(jìn)行有重疊的crop和merge操作，幫助對(duì)其中的細(xì)節(jié)進(jìn)行理解。通過本文的學(xué)習(xí)，對(duì)于下節(jié)推理階段的理解，有較大的幫助。

一、2D crop and merge

對(duì)于一個(gè)[10, 10]大小的示例圖像，采用patch大小為[3, 3]的進(jìn)行裁剪，每次patch與patch之間，在x和y方向重疊1個(gè)像素，無法構(gòu)成一個(gè)patch的部分，選擇丟棄，如下所示：

這里有一個(gè)共識(shí)，那就是當(dāng)輸入圖像的大小，patch的尺寸，以及overlap的size確定后，怎么裁剪，和能裁剪出多少個(gè)patch，就已經(jīng)確定了。下面我們就分析下，這個(gè)過程是什么樣子的。

1. 首先，決定行列可以有多少個(gè)patch，是左上角第一個(gè)patch中右下角的那個(gè)紅色點(diǎn)，因?yàn)樗堑谝粋€(gè)在水平和豎直方向都會(huì)需要重疊的點(diǎn)；
2. 反應(yīng)在行和列上面，可以移動(dòng)的區(qū)域也就是width - patch_width + 1，和height - patch_height + 1，因?yàn)閷?duì)于左上角第一個(gè)patch，只有右下角的坐標(biāo)，是參與到遍歷里面的，可參與遍歷的區(qū)域就是紅色曲線區(qū)域；
3. 對(duì)于行、列的每一步，都會(huì)重疊overlap_size個(gè)像素區(qū)域，所以可走的步長(zhǎng)，是patch_width - overlap_width，和patch_height - overlap_height；
4. 最后以左上角的坐標(biāo)點(diǎn)為準(zhǔn)，x:x+patch_width和y:y+patch_height就表示了一個(gè)區(qū)塊patch的像素范圍，被裁剪下來。

1.1、crop 操作

下面是實(shí)現(xiàn)上述步驟的代碼，主要就是要理解幾個(gè)點(diǎn)：

1. 在水平和數(shù)值方向，能取到多少個(gè)patch？
2. patch的滑動(dòng)選取，一次可以移動(dòng)多大的步長(zhǎng)？
3. 最后取像素塊，就簡(jiǎn)單了許多。

實(shí)操代碼如下：

import numpy as np

def crop_volume(volume, patch_size=[96, 96], overlap_size=[4, 4]):
    """
    Crop a 2D volume into patches with specified size and overlap.
    Args:
        volume (np.ndarray): the 3D volume to be cropped, with shape [width, height]
        patch_size (tuple or list): the size of patch, with format [patch_width, patch_height]
        overlap_size (tuple or list): the size of overlap between adjacent patches, with format [overlap_width, overlap_height]
    Returns:
        np.ndarray: the cropped patches, with shape [num_patches, patch_width, patch_height]
    """
    width, height = volume.shape
    patch_width, patch_height = patch_size

    overlap_width, overlap_height = overlap_size
    patches = []
    # 不夠一個(gè)patch，就丟棄
    for x in range(0, width - patch_width + 1, patch_width - overlap_width):
        for y in range(0, height - patch_height + 1, patch_height - overlap_height):
            print(x, y)
            patch = volume[x:x+patch_width, y:y+patch_height]
            patches.append(patch)
        print('\n')
    patches = np.asarray(patches)

    return patches


# 生成一個(gè)[10, 10]大小的示例圖像
imgs = np.random.rand(10, 10)
patch_size=[3, 3]
overlap_size = [1, 1]
# print('img shape:', imgs.shape)

patches = crop_volume(imgs, patch_size, overlap_size)
print('patches shape:', patches.shape)

驗(yàn)證了前面我們的猜想，后面我們直接取一個(gè)圖片，來驗(yàn)證下我們的思路。如果上述的思路是對(duì)的，那么，在裁剪后保存的圖像，就該是一個(gè)具體部分重疊區(qū)域，但是，還能夠反映全貌的一個(gè)個(gè)小圖。下面就是：

import os
import itk
import cv2
from matplotlib import pylab as plt
from PIL import Image
path = os.path.join(r'F:\tmp\results2', '10.png')
imgs = cv2.imread(path, 0)
volume_size = imgs.shape
patch_size=[96, 96]
overlap_size = [4, 4]
print('img shape:', imgs.shape)

patches = crop_volume(imgs, patch_size, overlap_size)
print('patches shape:', patches.shape)

for i in range(0, patches.shape[0], 1):
    one_patch = patches[i, :, :]
    print(i, one_patch.shape)
    width_p, height_p = one_patch.shape

    img_Image = Image.fromarray(one_patch)
    img_Image = img_Image.convert("L")

    img_Image.save(r"F:\tmp\results1/" + str(i) + ".png")

如下，是讀取的原始圖像，和crop后的一個(gè)個(gè)散落的小圖。盡管一個(gè)個(gè)小圖在我們展示的時(shí)候，他們之間使用間隙的，但并不影響我們看到他的全貌。

還是一種圖的樣子，區(qū)別在與彩色圖像成了灰度圖，最右側(cè)和最下側(cè)的像素像是少了一些。這是因?yàn)椴蛔阋粋€(gè)patch，被丟棄的原因?；诖耍覀円材芙o猜出來，在后續(xù)merge階段，可能會(huì)還原回去的圖像存在些許的差異。

1.2、merge 操作

merge是crop的一個(gè)逆過程。當(dāng)時(shí)我們?cè)趺床鸬?，現(xiàn)在就原路給拼接回去。拼接回去的圖像尺寸和crop前的圖像尺寸是一致的，當(dāng)時(shí)被忽略的部分，都是0。

1. 移動(dòng)的步長(zhǎng)還是一致的，在行列方向分別還是：patch_width - overlap_width，和patch_height - overlap_height；
2. 還要需要知道在行列方向，分別crop了多少次。以行為例，width - patch_width就是把第一個(gè)patch塊去除掉的行長(zhǎng)度，再除以(patch_width - overlap_width)，得到剩余部分可以走多少步，再+1，把第一個(gè)塊的數(shù)量補(bǔ)上，也就得到了在行方向上，有多少個(gè)patch了；
3. 以左上角的坐標(biāo)點(diǎn)為準(zhǔn)，x:x+patch_width和y:y+patch_height就表示了一個(gè)區(qū)塊patch了，現(xiàn)在就把對(duì)應(yīng)的patch像素，給賦值給原始圖了。
4. 因?yàn)?code>overlap重疊區(qū)域，被多次覆蓋，這部分需要求均值，巧妙的采用了一個(gè)新的像素塊，專門記錄被賦值的次數(shù)，最后除以對(duì)應(yīng)的次數(shù)，就可能實(shí)現(xiàn)求均值的過程了。

為什么先需要把第一個(gè)patch塊的行列方向都先去掉呢？

因?yàn)榈谝粋€(gè)塊是最特殊的，它被重疊的區(qū)域，只有overlap_size_w行和overlap_size_h例，而其余的patch塊，重疊區(qū)域都會(huì)是2行和2列，都遵循步長(zhǎng)的節(jié)奏。

width - patch_width把第一個(gè)patch塊去除掉后的行長(zhǎng)度，還能準(zhǔn)確反映有多少個(gè)patch嗎？

答案是可以的，這是因?yàn)闇p去一個(gè)patch，無非就是少了一個(gè)overlap_size_h的長(zhǎng)度，去掉一個(gè)overlap_size_h的長(zhǎng)度，如果恰好整除，那么加上這個(gè)長(zhǎng)度，也是多余的，無法再次構(gòu)成一個(gè)新的patch；即便有剩余，也是無法組成一個(gè)新的patch的。

下面是上面圖像crop階段裁剪得到的patch，加上merge后的操作，如下：

def merge_patches(patches, volume_size, overlap_size):
    """
    Merge the cropped patches into a complete 2D volume.
    Args:
        patches (np.ndarray): the cropped patches, with shape [num_patches, patch_width, patch_height]
        volume_size (tuple or list): the size of the complete volume, with format [width, height]
        overlap_size (tuple or list): the size of overlap between adjacent patches, with format [overlap_width, overlap_height]
    Returns:
        np.ndarray: the merged volume, with shape [width, height]
    """
    width, height = volume_size
    patch_width, patch_height = patches.shape[1:]
    overlap_width, overlap_height = overlap_size
    num_patches_x = (width - patch_width) // (patch_width - overlap_width) + 1
    num_patches_y = (height - patch_height) // (patch_height - overlap_height) + 1
    print('merge:', num_patches_x, num_patches_y)

    merged_volume = np.zeros(volume_size)
    weight_volume = np.zeros(volume_size)   # weight_volume的目的是用于記錄每個(gè)像素在裁剪過程中被遍歷的次數(shù)，最后用于求平均值
    idx = 0
    for x in range(num_patches_x):
        for y in range(num_patches_y):
            x_start = x * (patch_width - overlap_width)
            y_start = y * (patch_height - overlap_height)
            merged_volume[x_start:x_start+patch_width, y_start:y_start+patch_height] += patches[idx]
            weight_volume[x_start:x_start+patch_width, y_start:y_start+patch_height] += 1
            idx += 1
    merged_volume /= weight_volume
    return merged_volume

path = os.path.join(r'./results2', '10.png')
imgs = cv2.imread(path, 0)
volume_size = imgs.shape
patch_size=[96, 96]
overlap_size = [4, 4]
print('img shape:', imgs.shape)

patches = crop_volume(imgs, patch_size, overlap_size)
print('patches shape:', patches.shape)

merged_volume = merge_patches(patches, volume_size, overlap_size)
print('merged_volume shape:', merged_volume.shape)
merged_volume = Image.fromarray(merged_volume)
merged_volume = merged_volume.convert("L")
merged_volume.save(r"./results2/" + "merged_volume.png")

如下，果然和我們前面猜想的一樣，在merge后的圖像，相比于原圖，在右側(cè)和下側(cè)，都少了部分，這個(gè)問題后面在3D時(shí)候，再細(xì)細(xì)的探討。

二、3D crop and merge

3D部分相比于2D，無非就是多了一個(gè)深度信息，也就是z軸信息。所以在crop階段和merge階段，只需要按照2D行列的方式，增加一個(gè)維度的信息即可。

盡管說僅僅是增加了一個(gè)緯度，但是對(duì)于很多人來說，理解一個(gè)二維的平面是好理解的，但是突然變成了三維，還是有些拗不過彎的。此時(shí)，如果能繪制出一個(gè)三維的圖，幫助理解，就再好不過了。

代碼如下所示：

import numpy as np

def crop_volume(volume, patch_size=[96, 96, 96], overlap_size=[4, 4, 4]):
    """
    Crop a 3D volume into patches with specified size and overlap.
    Args:
        volume (np.ndarray): the 3D volume to be cropped, with shape [width, height, depth]
        patch_size (tuple or list): the size of patch, with format [patch_width, patch_height, patch_depth]
        overlap_size (tuple or list): the size of overlap between adjacent patches, with format [overlap_width, overlap_height, overlap_depth]
    Returns:
        np.ndarray: the cropped patches, with shape [num_patches, patch_width, patch_height, patch_depth]
    """
    width, height, depth = volume.shape
    patch_width, patch_height, patch_depth = patch_size

    overlap_width, overlap_height, overlap_depth = overlap_size
    patches = []
    for z in range(0, depth - patch_depth + 1, patch_depth - overlap_depth):
        for y in range(0, height - patch_height + 1, patch_height - overlap_height):
            for x in range(0, width - patch_width + 1, patch_width - overlap_width):
                patch = volume[x:x+patch_width, y:y+patch_height, z:z+patch_depth]
                patches.append(patch)
    patches = np.asarray(patches)
    return patches

def merge_patches(patches, volume_size, overlap_size):
    """
    Merge the cropped patches into a complete 3D volume.
    Args:
        patches (np.ndarray): the cropped patches, with shape [num_patches, patch_width, patch_height, patch_depth]
        volume_size (tuple or list): the size of the complete volume, with format [width, height, depth]
        overlap_size (tuple or list): the size of overlap between adjacent patches, with format [overlap_width, overlap_height, overlap_depth]
    Returns:
        np.ndarray: the merged volume, with shape [width, height, depth]
    """
    width, height, depth = volume_size
    patch_width, patch_height, patch_depth = patches.shape[1:]
    overlap_width, overlap_height, overlap_depth = overlap_size
    num_patches_x = (width - patch_width) // (patch_width - overlap_width) + 1
    num_patches_y = (height - patch_height) // (patch_height - overlap_height) + 1
    num_patches_z = (depth - patch_depth) // (patch_depth - overlap_depth) + 1
    print('merge:', num_patches_x, num_patches_y, num_patches_z)
    merged_volume = np.zeros(volume_size)
    weight_volume = np.zeros(volume_size)   # weight_volume的目的是用于記錄每個(gè)像素在裁剪過程中被遍歷的次數(shù)，最后用于求平均值
    idx = 0
    for z in range(num_patches_z):
        for y in range(num_patches_y):
            for x in range(num_patches_x):
                x_start = x * (patch_width - overlap_width)
                y_start = y * (patch_height - overlap_height)
                z_start = z * (patch_depth - overlap_depth)
                merged_volume[x_start:x_start+patch_width, y_start:y_start+patch_height, z_start:z_start+patch_depth] += patches[idx]
                weight_volume[x_start:x_start+patch_width, y_start:y_start+patch_height, z_start:z_start+patch_depth] += 1
                idx += 1
    merged_volume /= weight_volume
    return merged_volume


import os
import itk

nii_path = os.path.join(r'F:\tmp\results2', 'brain.nii.gz')
imgs = itk.array_from_image(itk.imread(nii_path))
imgs /= np.max(imgs)
volume_size = imgs.shape
patch_size = [96, 96, 96]
overlap_size = [32, 32, 32]
print('img shape:', imgs.shape)

# crop
patches = crop_volume(imgs, patch_size, overlap_size)
print('patches shape:', patches.shape)

print(patches.shape)
for d in range(0, patches.shape[0], 1):
    one_patch = patches[d, :, :, :]*255
    print(d, one_patch.shape)
    width_p, height_p, depth_p = one_patch.shape

    one_patch = itk.image_from_array(one_patch)

    itk.imwrite(one_patch, os.path.join(r'F:\tmp\results2\patch', str(d) + ".nii.gz"))  # 保存nii文件

# merge
merged_volume = merge_patches(patches, volume_size, overlap_size)
print('merged_volume shape:', merged_volume.shape)
merged_volume = itk.image_from_array(merged_volume)
itk.imwrite(merged_volume, r'F:\tmp\results2\merged_volume.nii.gz')

原始的圖像，如下：

crop一個(gè)塊圖像，如下：

patch_size = [96, 96, 96]，overlap_size = [4, 4, 4] 時(shí)，merge后的圖像，發(fā)現(xiàn)丟失像素區(qū)域較多，如下：

patch_size = [96, 96, 96]，overlap_size = [32, 32, 32] 時(shí)，merge后的圖像，發(fā)現(xiàn)丟失像素區(qū)域相對(duì)較少，如下：

對(duì)于大的patch size ，需要采用較大的overlap size，否則會(huì)導(dǎo)致右側(cè)和下側(cè)丟棄的像素區(qū)域過多，在合并的時(shí)候，會(huì)丟失較多信息。所以說為了避免這個(gè)問題，需要使得兩者之間達(dá)到一個(gè)相對(duì)的平衡。

究竟選擇的patch size有多大，這個(gè)還需要綜合來考慮。包括原始圖像的大小和GPU的內(nèi)存大小。太小了，對(duì)于空間相對(duì)位置信息的獲取，肯定是不利的。太大了也可能內(nèi)存占用太高，不太經(jīng)濟(jì)。

三、總結(jié)

這次就是一次對(duì)3D圖像大尺寸的一次處理的記錄，之前對(duì)這塊內(nèi)容涉獵較少，也沒有好好的理解。比如說對(duì)于腦部MRI的處理中，就選擇把無關(guān)區(qū)域去除，只留下可能存在目標(biāo)區(qū)域，用于后續(xù)的預(yù)測(cè)。這也是用全圖訓(xùn)練的一種取巧的方式。

至于新的patch方法，可以與老方法進(jìn)行相互的補(bǔ)充。尤其是對(duì)于尺寸較大的輸入，patch的方法，就更加的取巧了。

除此之外，推薦下MONAI這個(gè)庫和nnU-Net、nnFormer、unetr plus plus系列論文進(jìn)行閱讀和學(xué)習(xí)，對(duì)你做三維立體圖像的分割有非常大的幫助，尤其是數(shù)據(jù)的前處理。趕緊看完了就去搜索看看，麻溜的。文章來源地址http://www.zghlxwxcb.cn/news/detail-676773.html

到了這里，關(guān)于【醫(yī)學(xué)影像數(shù)據(jù)處理】2D/3D patch的crop和merge操作匯總的文章就介紹完了。如果您還想了解更多內(nèi)容，請(qǐng)?jiān)谟疑辖撬阉鱐OY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關(guān)文章，希望大家以后多多支持TOY模板網(wǎng)！