文章目录

    • 数据整理
      • 整理出所有的dcm文件
      • 整理出所有的xml标注文件
      • 整理数据①——舍弃错误的标注文件
      • 整理数据②——两个标注文件指向同一个目标图片的情况
        • 封装函数,传入xml文件,显示标注效果
      • 整理数据③——将PETCT的三通道图像转成平扫CT的单通道图像格式
        • 关于PETCT三通道转单通道灰度图的一些思考
      • 建立Dataset

本文进一步对数据集进行探索,对数据进行清洗、整理、验证。
最后建立Dataset对象,以便进行后续的训练工作。
本文所用代码: 我的Github

数据整理

我们最终的目的是建立Dataset对象,所以我们的目标应该是整理出 数据+标签 的组合。
在目标检测任务中,数据即我们的CT图片(读取自dcm文件);标签即我们的xml文件,内含 Bounding Box框 + 类别标签。
所以我们的目标是对所有可配对的 dcm文件xml文件 进行一一对应,剩下读取则是一件简单的工作。

整理出所有的dcm文件

dcm文件和xml文件对应的方法在上一节中已经阐释。
现在需要将数据集中所有的 dcm文件 和 对应的SOP Instance UID(后通称UID)整理出来,一共有25万条记录。


import pydicomimport matplotlib.pyplot as pltimport osfrom tqdm import tqdmimport pandas as pdimport numpy as np# 查看整个数据集下所有的dcm文件名dcm_file=[]for root, dirs, files in os.walk('manifest-1608669183333/Lung-PET-CT-Dx'):for file in files:file_path = os.path.join(root, file)if 'dcm' in file_path:dcm_file.append(file_path)print(dcm_file[0])len(dcm_file)


导入DataFrame后,获取对应的UID。
由于有25w张图片,我在固态硬盘上跑这段代码耗时约50min,跑完后保存为csv,下次就不用再跑一遍了。

df_dcm = pd.DataFrame(dcm_file)# 建立新列uid,获取dcm文件对应的 SOPInstanceUID 号,这步耗时50min左右df_dcm['uid'] = df_dcm[0].apply(lambda x: pydicom.dcmread(x).SOPInstanceUID)df_dcm['uid_str'] = df_dcm['uid'].apply(lambda x: str(x))df_dcm.to_csv('dcm_file_uid.csv')# 上面一步耗时较久,如果已经读取过一遍并已保存,这里直接读取df_dcm = pd.read_csv('dcm_file_uid.csv', index_col=0)df_dcm


简单检查一下dcm文件是否独一无二。

整理出所有的xml标注文件

我们再把所有的xml文件整理出来,附上对应的UID。


整理数据①——舍弃错误的标注文件

我们先检查一下3万多个xml文件对应的UID是否都对应在df_dcm中。

# 检查Annotation目录下所有的uid文件名是否都在df_dcm中xml_file_df = xml_file_df.loc[xml_file_df['uid_str'].isin(df_dcm['uid_str'])]xml_file_df# 我们发现 30884 < 31562


结果我们发现并不是所有的xml都有对应原图片,这应该属于这个数据集的标注错误问题,有可能标注文件的对应UID写错了。
总之,我们不能将这些xml文件纳入训练,所以将这些找不到对应图片的标注文件舍弃掉。

整理数据②——两个标注文件指向同一个目标图片的情况

我们查看一下XML标注文件是否有重复的文件名。

# 查看Annotation目录下的XML文件是否有重复的文件名len(xml_file_df[0].unique()) == len(xml_file_df['uid_str'].unique())# 如有False,说明Annotation下有重复的文件


这里我们发现该数据集第二个坑,Annotation目录下存在重复的标注文件指向同一个目标图片的情况(否则上面就显示True)。
我们继续进一步探究是哪些标注文件重复了。


我们查看一下原本的XML文件数量 和 去重后的XML文件数量。

# 显示Annotation目录下去重后的XML文件数量、原本的XML文件数量print(len(xml_file_df['uid_str'].unique()), len(xml_file_df[0].unique()) )


我们获取一下那些重名(同一UID)的xml文件。

# 获取那些重复的xml文件的索引idx = xml_file_df['uid_str'][xml_file_df['uid_str'].duplicated(keep=False)].index# 对xml_file_df文件应用idx索引,获取完整文件名xml_file_name_duplicated = pd.DataFrame(xml_file_df.loc[idx,:])# 导出到csv文件(记录一下)xml_file_name_duplicated.to_csv('duplicated.csv')xml_file_name_duplicated[:6]

xml_file_name_duplicated = pd.read_csv('duplicated.csv', index_col=0)print(len(xml_file_name_duplicated))# 显示重复条目数量xml_file_name_duplicated['0'][:12].values# 显示前12行# 可以看到A0003和A0004中存在相同名称的文件。# 我们知道A0003和A0004是两个不同的患者(Clinical Data的“statistics-clinical-20201221.xlsx”文件),理应不可能存在相同的CT影像,所以应该是xml标注文件本身出错了


以上我们可以看到A0003和A0004中存在相同名称的文件。
我们知道A0003和A0004是两个不同的患者(从Clinical Data的“statistics-clinical-20201221.xlsx”文件看出),理应不可能存在相同的CT影像,所以应该是xml标注文件本身出错了。
那么具体是怎么出错了呢?


# 查看重复条目有多少种,如果为重复条目总数的1/2,那么均为两两重复print('重复条目有', len(xml_file_name_duplicated['uid_str'].unique()), '种,总条目为', len(xml_file_name_duplicated), '个')# 下面我们看看这些重复的标注文件的标注情况


上面我们得知重复标注的情况是两两重复。
我们再进一步查看一下重复标注的情况。


封装函数,传入xml文件,显示标注效果

我们设计一个函数,传入xml文件路径,显示图片和锚框。

# 设计一个函数,传入xml文件路径,显示图片和锚框################### 从XML文件获取bounding box信息import xml.etree.ElementTree as ETdef get_labelFromXml(xml_file): label=[]bbox_list=[]an_file = open(xml_file, encoding='utf-8')tree=ET.parse(an_file)root = tree.getroot()for object in root.findall('object'):cancer_type=object.find('name').text.upper()xmin=object.find('bndbox').find('xmin').textxmax=object.find('bndbox').find('xmax').textymin=object.find('bndbox').find('ymin').textymax=object.find('bndbox').find('ymax').textif int(xmin)==0 or int(xmax)==0 or int(ymin)==0 or int(ymax)==0:passelif int(xmin)==int(xmax) or int(ymin)==int(ymax):passelse:bbox=[int(xmin),int(ymin),int(xmax),int(ymax)]bbox_list.append(bbox)label.append(cancer_type)return bbox_list,label# 画布上添加锚框def bbox_to_rect(bbox, color):# 将边界框(左上x,左上y,右下x,右下y)格式转换成matplotlib格式:# ((左上x,左上y),宽,高)return plt.Rectangle(xy=(bbox[0], bbox[1]), width=bbox[2]-bbox[0], height=bbox[3]-bbox[1],fill=False, edgecolor=color, linewidth=2)def show_img_and_anchors(url):# 传入xml文件路径,展示标注效果,一张图片可显示多个bbox框# 获取uid = xml_file_df.loc[xml_file_df[0]==url, 'uid_str'].values[0]im = df_dcm.loc[df_dcm['uid_str']==uid].values[0][0]im = pydicom.read_file(im)# 获取像素矩阵img_arr = im.pixel_arrayfig = plt.imshow(img_arr,cmap=plt.cm.bone)plt.title("UID:{}".format(uid))bbox, label = get_labelFromXml(url)for i, j in zip(bbox, label):fig.axes.add_patch(bbox_to_rect(i, 'red'))fig.axes.text(i[0]+12, i[1]+12, j,va='center', ha='center', fontsize=12, color='red')def show_imgs_and_anchors(urls, num=5):# 传入xml文件路径组成的列表,展示多张图片发标注效果num_toshow = min(num, len(urls))print('载入', len(urls), '张图片,显示前', num, '张。')# 获取temp_df = xml_file_df.loc[xml_file_df[0].isin(urls)]urls = temp_df[0].values# 这里urls顺序被重排,需要重新赋值以此,以和uids建立一一对应uids = temp_df['uid_str'].valuesims = [df_dcm.loc[df_dcm['uid_str']==x].values[0][0] for x in uids]ims_dcm = [pydicom.read_file(x).pixel_array for x in ims]# 获取像素矩阵List# 取标注框和类别bbox_label = [get_labelFromXml(i) for i in urls]# _表示忽略不适用的变量,返回的fig用不上_, axes = plt.subplots(nrows=1, ncols=num_toshow, figsize=(48,48))# 使用zip方法,在for循环中设置axes中的各个子区域的参数并绘图for ax, uid, img, lbl in zip(axes, uids, ims_dcm, bbox_label):ax.imshow(img, cmap=plt.cm.bone)ax.set_title(uid)bbox, label = lblfor i, j in zip(bbox, label):ax.add_patch(bbox_to_rect(i, 'red'))ax.text(i[0]+12, i[1]+12, j,va='center', ha='center', fontsize=12, color='red')ax.axes.get_xaxis().set_visible(False)ax.axes.get_yaxis().set_visible(False)

展示标注效果(一个xml文件):


展示标注效果(多个xml文件):

我们看到上方第1行第1张图片和第2行第1张图片是标注了同一张CT图片的情况。
我们可以看到标注的质量差不多,若标注得差不多,则可以假设不同的标注框可以作为一种数据增强的形式,因此可以无需特意舍弃重复的对应同UID的xml文件。
下面抽样查看一下所有的的重复标记情况。


我们可以再随机抽20张看看

# 随机抽20张看看slice = np.random.choice(len(xml_file_name_duplicated), size=(4,5), replace=False)# 生成4×5整数矩阵for i in slice:urls = xml_file_name_duplicated['0'].values[i]show_imgs_and_anchors(urls, num=5)

效果见源码文件。

整理数据③——将PETCT的三通道图像转成平扫CT的单通道图像格式

普通CT图像是单通道(512×512),但PETCT图像是三通道(512×512×3),需转成单通道进行训练。
因为日后我们希望使用平扫的CT片子进行推理预测。

我们看看读取PETCT和普通CT片子的图片shape:

# 各取一张ct,一张petct图像的xml文件url_ct = 'Lung-PET-CT-Dx-Annotations-XML-Files-rev12222020/Annotation\\A0213\\1.3.6.1.4.1.14519.5.2.1.6655.2359.211226049141098860991228194230.xml'# 普通CT图像url_petct = 'Lung-PET-CT-Dx-Annotations-XML-Files-rev12222020/Annotation\\A0213\\1.3.6.1.4.1.14519.5.2.1.6655.2359.284461122128650673200074494931.xml'# PETCT图像# 从xml文件获取对应ct图像文件信息def get_img_pixel_array(url):uid = xml_file_df.loc[xml_file_df[0]==url, 'uid_str'].values[0]im = df_dcm.loc[df_dcm['uid_str']==uid].values[0][0]im = pydicom.read_file(im)img_arr = im.pixel_array# 获取像素矩阵return img_arrprint('ct图像矩阵的形状为,', get_img_pixel_array(url_ct).shape)print('petct图像矩阵的形状为,', get_img_pixel_array(url_petct).shape)


关于PETCT三通道转单通道灰度图的一些思考

我这里使用了opencv库的方法将PETCT的三通道转成单通道灰度图。
在opencv库中,三通道转单通道主要是两种转换方法:COLOR_BGR2GRAY or COLOR_RGB2GRAY

我们知道,直接使用cv2读取后图片的通道顺序是BGR,而pydicom读取的pixel_array的通道顺序应该是RGB(调用PhotometricInterpretation属性可以看到)。
所以对pixel_array直接使用cv2的转灰度方法时,可以使用COLOR_RGB2GRAY模式。

但是我们知道PETCT的图像信息是 正常的CT影像 和 代谢摄取增高影(高亮部分)两者重合。
我们发现使用COLOR_BGR2GRAY模式转灰度后,代谢摄取增广的高亮信息较COLOR_RGB2GRAY模式不明显。(下面的代码块可见效果)
这样似乎更适合用于训练,因为后面我们需要常用普通平扫的CT片子进行预测。

因此,转单通道我采用了cv2的COLOR_BGR2GRAY模式。

关于三通道转灰度图的一些其他参考资料:

  1. https://note.nkmk.me/python-opencv-numpy-color-to-gray/
  2. https://stackoverflow.com/questions/62855718/why-would-cv2-color-rgb2gray-and-cv2-color-bgr2gray-give-different-results

两种转单通道模式的对比:

import cv2 as cvimg_array = get_img_pixel_array(url_petct)img_array_gray = cv.cvtColor(img_array, cv.COLOR_BGR2GRAY)img_array_gray2 = cv.cvtColor(img_array, cv.COLOR_RGB2GRAY)# # img_array = np.array(img_array, dtype=np.uint16)_, axes = plt.subplots(1,3,figsize=(24,24))axes[0].imshow(img_array, cmap=plt.cm.bone)axes[0].set_title('PETCT')axes[1].imshow(img_array_gray, cmap=plt.cm.bone)axes[1].set_title('COLOR_BGR2GRAY')axes[2].imshow(img_array_gray2, cmap=plt.cm.bone)axes[2].set_title('COLOR_RGB2GRAY')plt.show()plt.close()

建立Dataset

有了前面的数据整理,我们可以建立Dataset数据集对象了。
目标是传入 [dcm图片集] 和 对应的 [xml标注集] ,即生成一个Dataset对象。

import torchfrom torch.utils.data import Dataset, DataLoaderfrom torchvision import transformsfrom PIL import Imagetransform=transforms.Compose([transforms.ToTensor(),])

但是,这里我们要注意,pydicom从CT读取的pixel_array数据类型(dtype)为uint16,从PETCT读取的pixel_array则为uint8
这会导致Image.open(PIL库)读取并使用transforms.ToTensor()转成Tensor后数据类型不一样。前者变为torch.int16,后者变为torch.float32。(如下演示)

因此,我们需要将pixel_array统一转成float32


创建Dataset数据集对象(Pytorch)。
注意一下Dataset输出的数据格式,img为图片Tensor,label则是一个dict,内含bbox的坐标和类别。

class_to_id={'A':1, 'B':2, 'E':3, 'G':4}# 重写get_labelFromXml,返回label的iddef get_labelidFromXml(xml_file): label=[]bbox_list=[]an_file = open(xml_file, encoding='utf-8')tree=ET.parse(an_file)root = tree.getroot()for object in root.findall('object'):cancer_type=object.find('name').text.upper()xmin=object.find('bndbox').find('xmin').textxmax=object.find('bndbox').find('xmax').textymin=object.find('bndbox').find('ymin').textymax=object.find('bndbox').find('ymax').textif int(xmin)==0 or int(xmax)==0 or int(ymin)==0 or int(ymax)==0:passelif int(xmin)==int(xmax) or int(ymin)==int(ymax):passelse:bbox=[int(xmin),int(ymin),int(xmax),int(ymax)]bbox_list.append(bbox)label.append(class_to_id[cancer_type])return bbox_list,label# 计划传入img和label,以此对应dcm文件路径集和xml文件路径集class LungDetection(Dataset):def __init__(self, img, label, transform=transform,):self.img = imgself.label = labelself.transform = transformdef __getitem__(self, index):img = self.img[index]label = self.label[index]img_open=pydicom.read_file(img)img_array=img_open.pixel_arrayif len(img_array.shape) == 3:img_array = cv.cvtColor(img_array, cv.COLOR_BGR2GRAY)img_array = np.array(img_array, dtype=np.float32)img_pic = Image.fromarray(img_array)img_tensor=self.transform(img_pic)bbox, label = get_labelidFromXml(label)bbox_tensor = torch.as_tensor(bbox,dtype=torch.float32)label_tensor = torch.as_tensor(label,dtype=torch.int64)target={}target['boxes']=bbox_tensortarget['labels']=label_tensorreturn img_tensor, targetdef __len__(self):return len(self.img)

xml_file_df 是我们之前整理出的xml标签文件,我们现在将xml文件和dcm文件对应起来。

# xml_file_df 是我们之前整理出的xml标签文件,我们现在将xml文件和dcm文件对应起来xml_file_dataset = xml_file_df.reset_index(drop=True)xml_file_dataset.columns = ['xml', 'uid_str']xml_file_dataset

根据xml的UID对应上dcm文件:
下面这段代码跑起来需要几分钟的时间,所以跑完后可以保存一下csv文件

xml_file_dataset['dcm'] = xml_file_dataset['uid_str'].apply(lambda x : df_dcm.loc[df_dcm['uid_str']==x].values[0][0])xml_file_dataset.to_csv('xml_file_dataset.csv')# 以上代码需要跑8分钟,所以保存起来xml_file_dataset = pd.read_csv('xml_file_dataset.csv', index_col=0)xml_file_dataset


传入Dataset。

img = xml_file_dataset['dcm'].valueslabel = xml_file_dataset['xml'].valueslungdetection_dataset = LungDetection(img=img, label=label, transform=transform)

顺利完成!