影像组学特征提取 — 原始图像为dicom格式,mask图像为nrrd格式

该博客详细介绍了如何处理DICOM格式的原始图像和NRRD格式的mask图像,进行Radiomics特征提取。首先,通过SimpleITK库读取两种格式的图像,然后设置日志和计算参数。接着,使用Radiomics库初始化特征提取器,提取包括FirstOrder在内的特征。对于可能出现的问题,如多个mask、图像尺寸不匹配和坐标差异,提出了相应的解决策略。最后,通过循环遍历Excel中的subject ID,批量处理每个subject的多个mask,将计算出的特征和相关信息写入CSV文件。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

原始图像为dicom格式,mask图像为nrrd格式

1. 读取 dicom 和 nrrd

2. 设置 logger 和 setting

3. 初始化特征提取器,设置图像空间和特征类型

4. 特征提取和保存显示

单个subject影像组学特征提取

from __future__ import print_function
import logging
import SimpleITK as sitk
import radiomics
from radiomics import featureextractor
import six

# The original image is in dicom format
readerC = sitk.ImageSeriesReader()
dicom_names = readerC.GetGDCMSeriesFileNames('/project/patient/000000')
readerC.SetFileNames(dicom_names)
readerC.MetaDataDictionaryArrayUpdateOn()
readerC.LoadPrivateTagsOn()
imageName = readerC.Execute()

# The mask image is in nrrd format
maskName = sitk.ReadImage('/project/mask/000000/mask.nrrd')

# Logger setting
logger = radiomics.logger
logger.setLevel(logging.DEBUG)  
handler = logging.FileHandler(filename='testLog.txt', mode='w')
formatter = logging.Formatter("%(levelname)s:%(name)s: %(message)s")
handler.setFormatter(formatter)
logger.addHandler(handler)

# Define settings for signature calculation
# These are currently set equal to the respective default values
settings = {}
settings['binWidth'] = 25
settings['resampledPixelSpacing'] = None  # [3,3,3] is an example for defining resampling (voxels with size 3x3x3mm)
settings['interpolator'] = sitk.sitkBSpline
settings['correctMask'] = True
settings['geometryTolerance'] = 1


# Initialize feature extractor
extractor = featureextractor.RadiomicsFeatureExtractor(**settings)

# Enable the image type
extractor.enableAllImageTypes()
# extractor.enableImageTypes(Original={}, LoG={}, Wavelet={})

# Enable the feature type
# extractor.disableAllFeatures()
extractor.enableFeatureClassByName('firstorder')
# extractor.enableFeaturesByName(firstorder=['Mean', 'Skewness'])
# extractor.enableAllFeatures()

# Output the result
print("Calculating features")
featureVector = extractor.execute(imageName, maskName)
radiomicsList = []
header = []
for key, val in six.iteritems(featureVector):
    if not key.startswith('diagnostics'):
        header.append(key)
        radiomicsList.append(str(val))
print(header)
print(radiomicsList)

for featureName in featureVector.keys():
    print("Computed %s: %s" % (featureName, featureVector[featureName]))

从excel中读取subject名字,进行批量处理

问题一:由于一个 subject 存在多个 mask,所以 excel 中的 subject ID 一列存在重复项

解决方案:读取到一个新的 subject ID 之后,遍历该文件夹内的所有 mask,进行组学特征提取;若当前读取的subject ID为已处理过的重复项,则跳过

问题二:原始图像为 dicom,mask为 nrrd ,存在 dicom 尺寸和 mask 尺寸不匹配的情况

解决方案:读取 dicom 序列文件,计算其层数和尺寸,读取nrrd文件,计算其层数和尺寸,如果两者相等,则计算组学特征并保存,如果不相等,则跳过

问题三:由于影像组学计算过程中需要 image 和 mask 的严格位置匹配,但是由于保存原因,dicom 和 nrrd 的坐标中心可能存在极微小的差异,比如 2.99999 和 3.0 这样子

解决方案:在 setting 中设置 correct mask 和 geometry tolerence 选项,进行位置自动匹配调整

问题四:需要将 subject ID,mask name,和计算出来的组学特征写入 csv 文件,并在第一行保存标题

解决方案:要考虑到问题二中的情况,不匹配的就只写入subject ID 和 mask name;要设置一个count,当 count = 0 才保存标题;要读取 mask 的名字,写入 mask name

import numpy as np
import xlrd
import os
import nrrd
import six
import csv
import SimpleITK as sitk
from radiomics import featureextractor


def count_file_number(filepath, filetype):
    count = 0
    for root, dirname, filenames in os.walk(filepath):
        for filename in filenames:
            if os.path.splitext(filename)[1] == filetype:
                count += 1
    return count, filenames


def dcmseriesread(dicompath):
    readerC = sitk.ImageSeriesReader()
    dicom_names = readerC.GetGDCMSeriesFileNames(dicompath)
    readerC.SetFileNames(dicom_names)
    readerC.MetaDataDictionaryArrayUpdateOn()
    readerC.LoadPrivateTagsOn()
    dicomImage = readerC.Execute()
    return dicomImage


def radiomics_feature_extractor(image, mask):
    settings = {}
    settings['binWidth'] = 25
    settings['resampledPixelSpacing'] = None  # [3,3,3] is an example for defining resampling (voxels with size 3x3x3mm)
    settings['interpolator'] = sitk.sitkBSpline
    settings['correctMask'] = True
    settings['geometryTolerance'] = 1
    extractor = featureextractor.RadiomicsFeatureExtractor(**settings)
    extractor.enableAllImageTypes()
    # extractor.enableFeatureClassByName('firstorder')
    # extractor.enableFeatureClassByName('shape', 'texture')
    featureVector = extractor.execute(image, mask)
    radiomicsList = []
    header = []
    for key, val in six.iteritems(featureVector):
        if not key.startswith('diagnostics'):
            header.append(key)
            radiomicsList.append(str(val))
    return header, radiomicsList


path = '/project/'
SKindex = xlrd.open_workbook(os.path.join(path, 'index.xlsx')).sheets()[0]
subIDtemp = np.array(SKindex.col_values(0))[12:]
subID = [x[:-2].zfill(6) for x in subIDtemp]
# print(subID)
counter = 0

for i in range(len(subID)):
    if i > 0 and subID[i] == subID[i-1]:
        pass
    else:
        patientPath = os.path.join(path, 'patient', subID[i])
        maskPath = os.path.join(path, 'mask', subID[i])
        dicomSlices, dicomNames = count_file_number(patientPath, '.dcm')
        originalImage = dcmseriesread(patientPath)
        maskNumber, maskNames = count_file_number(maskPath, '.nrrd')

        for maskName in maskNames:
            print([str(subID[i]), maskName, 'processing....'])
            maskMatrix, options = nrrd.read(os.path.join(maskPath, maskName))
            maskSlices = maskMatrix.shape[-1]

            if maskSlices == dicomSlices:
                maskImage = sitk.ReadImage(os.path.join(maskPath, maskName))
                header, radiomicsList = radiomics_feature_extractor(originalImage, maskImage)
                with open(os.path.join('/Desktop/data_files', 'radiomics_feature_all.csv'), 'a', newline='') as outcsv:
                    writer = csv.writer(outcsv)
                    if counter == 0:
                        writer.writerow(['patientID', 'maskName'] + header)
                    writer.writerow([str(subID[i]), maskName] + radiomicsList)
                    counter += 1

            else:
                with open(os.path.join('/Desktop/data_files', 'radiomics_feature_all.csv'), 'a', newline='') as outcsv:
                    writer = csv.writer(outcsv)
                    writer.writerow([str(subID[i]), maskName])

### 影像组学特征提取方法及工具 #### 方法概述 影像组学是一种通过高通量技术从医学影像中提取大量定量特征的方法,这些特征能够反映肿瘤的空间异质性和其他生物学特性[^1]。具体来说,影像组学特征提取过程通常分为以下几个阶段: 1. **影像数据获取** 获取高质量的医学影像数据是整个流程的基础。常见的影像模式包括CT、MRI和PET等[^3]。 2. **感兴趣区域(ROI/VOI)分割** 使用手动、半自动或全自动算法对目标区域进行精确分割。这是后续特征提取的关键步骤之一[^2]。 3. **特征提取** 提取的特征主要包括形状特征、强度统计特征以及纹理特征。其中,灰度共生矩阵(GLCM)、灰度游程矩阵(GLRLM)和小波变换常用于描述图像中的空间关系和频率分布[^2]。 4. **特征降维与选择** 鉴于原始特征数量庞大且可能存在冗余,因此需要采用主成分分析(PCA)、LASSO回归或其他机器学习方法来筛选最具代表性的特征集合[^1]。 5. **建模与验证** 利用选定的特征构建预测模型,并评估其性能指标如AUC值、敏感性、特异性等。常用的分类器有支持向量机(SVM)、随机森林(RF)和深度神经网络(DNN)等。 #### 工具推荐 目前存在多种开源软件包可以帮助研究人员完成上述任务,以下是几个主流选项及其特点说明: - **PyRadiomics**: 这是一个基于Python开发的强大库,专为标准化影像组学特征计算设计。它不仅兼容DICOM标准还支持NIfTI/NRRD等多种文件格式输入输出操作。 下面展示如何利用该框架执行简单的特征提取工作: ```python from radiomics import featureextractor extractor = featureextractor.RadiomicsFeatureExtractor() result = extractor.execute('path_to_image.nrrd', 'path_to_mask.nrrd') for key, value in result.items(): print(f"{key}: {value}") ``` - **ITK-SNAP & 3D Slicer**: 虽然它们主要定位于交互式的三维可视化平台,但也提供了插件扩展功能实现自动化批量处理需求[^2]。 - **MATLAB Image Processing Toolbox**: 对熟悉Matlab环境的人来说非常友好,内置丰富的函数可以直接调用来生成所需的各类统计量或者频谱参数表征对象内部结构信息[^2]。 综上所述,在实际应用过程中可以根据项目具体情况灵活选用合适的解决方案组合起来共同发挥作用最大化研究价值所在! 问题
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值