使用pyradiomics提取影像组学特征【详细】

使用pyradiomics提取影像组学特征【详细】

最近由于项目需求要使用pyradiomics提取影像组学特征,网上阅读了很多别人的博客,学到一些,然后去查看了pyradiomics的官方文档,最后自己实现了特征的提取,写下此文记录,方便日后查看。

首先放上官方文档:https://pyradiomics.readthedocs.io/
和pyradiomics源代码地址:https://github.com/AIM-Harvard/pyradiomics

第一步、安装

在已经安装好的Python环境下,执行:

pip install pyradiomics

第二步、提取特征

通过自定义提取特征可以根据自己的需求来设置并提取。

  • 首先,可以指定用于提取特征的图像类型。在pyradiomics包中为我们提供了许多可以使用的滤波器,所以我们可以使用原始图像及经过各种滤波器之后的图像,如下方表格,具体可以参考官方:图像类型
  • 图像类型 解释
    Original 原始图像
    Wavelet 小波变换。产生在三个维度中每个维度分别使用高通、低通滤波器的所有组合(LLH、LHL、LHH、HLL、HLH、HHL、HHH、LLL)
    LoG Laplacian of Gaussian filter高斯滤波器的拉普拉斯算子,是一种边缘增强滤波器。使用它需要指定参数sigma,低 sigma 强调精细纹理,高 sigma 值强调粗糙纹理
    Square 平方。取原始像素的平方并将它们线性缩放回原始范围
    SquareRoot 平方根。取绝对图像强度的平方根并将它们缩放回原始范围
    Logarithm 对数。取绝对强度 + 1 的对数,值缩放到原始范围
    Exponential 指数。采用e^(绝对强度)获取强度的指数值,值被缩放到原始范围
    Gradient 梯度。返回局部梯度的大小
    LocalBinaryPattern2D 在每一片中进行的本地二进制模式
    LocalBinaryPattern3D 在3d中进行的本地二进制模式

    指定方式:设置特征提取器后可以在下边指定(如何设置特征提取器:extractor
    可以直接使用:enableAllImageTypes()启用所有类型
    也可以使用:enableImageTypeByName( imageType , enabled=True , customArgs=None )启用你想用的类型
    例如:

    # 所有类型
    extractor.enableAllFeatures()
    # 指定使用LoG和Wavelet滤波器
    extractor.enableImageTypeByName('LoG')
    extractor.enableImageTypeByName('Wavelet')
    
  • 然后可以指定所要提取的特征。pyradiomics包也为我们提供了很多种可选的特征,如下方表格,具体特征类型详细资料可参考:radiomics-features
  • 特征类型 特征 解释
    First Order Features(共19个) Energy 能量
    Total Energy 总能量
    Entropy
    Minimum 最小特征值
    10Percentile 特征值的百分之10的值
    90Percentile 特征值得百分之90的值
    Maximum 最大特征值
    Mean 均值
    Median 中位数
    InterquartileRange 四分位距离
    Range 灰度值范围
    MeanAbsoluteDeviation (MAD)平均绝对误差
    RobustMeanAbsoluteDeviation (rMAD) 鲁棒平均绝对偏差
    RootMeanSquared (RMS)均方根误差
    StandardDeviation 标准差。测量平均值的变化或离散量,默认不启用,因为与方差相关
    Skewness 偏度。测量值的分布关于平均值的不对称性
    Kurtosis 峰度。是图像 ROI 中值分布的“峰值”的量度
    Variance 方差。是每个强度值与平均值的平方距离的平均值
    Uniformity 均匀度。是每个强度值的平方和
    特征类型 特征 解释
    Shape Features (3D) (共17个) MeshVolume 网格体积
    VoxelVolume 体素体积
    SurfaceArea 表面积
    SurfaceVolumeRatio 表面积体积比
    Sphericity 球形度
    Compactness1 紧凑度1。衡量肿瘤形状相对于球形的紧密程度,由于与球形度、紧凑度2、球形不对称相关,默认不启用
    Compactness2 紧凑度2。衡量肿瘤形状相对于球形的紧密程度,与1的计算公式不同,默认不启用
    SphericalDisproportion 球形不对称。是肿瘤区域的表面积与与肿瘤区域具有相同体积的球体的表面积之比,是球形度的倒数,默认不启用
    Maximum3DDiameter 最大3D直径
    Maximum2DDiameterSlice 最大二维直径(切片),定义为行列(通常是轴向)平面中肿瘤表面网格顶点之间的最大成对欧几里得距离。
    Maximum2DDiameterColumn 最大2D直径(列)。定义为行切片(通常是冠状)平面中肿瘤表面网格顶点之间的最大成对欧几里得距离。
    Maximum2DDiameterRow 最大二维直径(行)。定义为列切片(通常为矢状)平面中肿瘤表面网格顶点之间的最大成对欧几里得距离。
    MajorAxisLength 长轴长度
    MinorAxisLength 短轴长度
    LeastAxisLength 最小轴长度
    Elongation 伸长率。表示ROI形状中两个最大的主成分之间的关系
    Flatness 平整度。表示ROI形状中最大和最小主成分之间的关系
    特征类型 特征 解释
    Shape Features (2D) MeshSurface 网格面。由三角形网格定义,计算ROI的表面
    PixelSurface 像素面。ROI 中的像素数乘以单个像素的表面积
    Perimeter 周长
    PerimeterSurfaceRatio 周长比。值越低越紧凑
    Sphericity 球形度
    SphericalDisproportion 球形不均匀度。是球形度的倒数 ,与球形度相关,默认不启用
    MaximumDiameter 最大二维直径。定义为肿瘤表面网格顶点之间的最大成对欧几里得距离
    MajorAxisLength 长轴长度
    MinorAxisLength 短轴长度
    Elongation 伸长率。显示了 ROI 形状中两个最大主成分之间的关​​系
    特征类型 特征 解释
    Gray Level Co-occurrence Matrix (GLCM) Features灰度共生矩阵特征(共24个) Autocorrelation 自相关。是纹理细度和粗糙度大小的量度
    JointAverage 联合平均
    ClusterProminence 集群显著性。衡量 GLCM 偏度和不对称性的指标
    ClusterShade 集群阴影。衡量 GLCM 偏度和均匀度的指标
    ClusterTendency 集群趋势
    Contrast 对比度
    Correlation 相关性
    DifferenceAverage 平均差
    DifferenceEntropy 差分熵
    DifferenceVariance 差异方差。是异质性的一种度量
    Dissimilarity 差异性。已弃用
    JointEnergy 联合能量
    JointEntropy 联合熵。是邻域强度值的随机性/可变性的度量
    Imc1 相关性信息度量 (IMC) 1
    Imc2 相关性信息度量 (IMC) 2
    Idm 逆差矩
    MCC 最大相关系数
    Idmn 逆差分矩归一化
    Id 逆差
    Idn 逆差分归一化
    InverseVariance 逆方差
    MaximumProbability 最大概率
    SumAverage 求和平均
    SumEntropy 总熵
    SumSquares 平方和
    特征类型 特征 解释
    Gray Level Size Zone Matrix (GLSZM) Features灰度级大小区域矩阵特征 (共16个) SmallAreaEmphasis 是小尺寸区域分布的量度,数值越大表示区域越小,纹理越精细
    LargeAreaEmphasis 是大面积区域分布的度量,值越大表示区域越大,纹理越粗糙。
    GrayLevelNonUniformity 灰度非均匀性。测量图像中灰度级强度值的可变性,值越低表示强度值越均匀。
    GrayLevelNonUniformityNormalized 灰度非均匀归一化。测量图像中灰度强度值的可变性,值越低表示强度值的相似性越大
    SizeZoneNonUniformity 尺寸区域非均匀性。测量图像中大小区域体积的可变性,值越低表示大小区域体积越均匀
    SizeZoneNonUniformityNormalized 尺寸区域非均匀归一化。测量整个图像中大小区域体积的可变性,值越低表示图像中区域大小体积之间的同质性越高
    ZonePercentage 区域百分比。获取 ROI 中的区域数和体素数的比率来测量纹理的粗糙度
    GrayLevelVariance 灰度方差
    ZoneVariance 区域方差
    ZoneEntropy 区域熵
    LowGrayLevelZoneEmphasis 低灰度区域强调
    HighGrayLevelZoneEmphasis 高灰度区域强调
    SmallAreaLowGrayLevelEmphasis 小面积低灰度区域强调
    SmallAreaHighGrayLevelEmphasis 小区域高灰度强调
    LargeAreaLowGrayLevelEmphasis 大面积低灰度强调
    LargeAreaHighGrayLevelEmphasis 大面积高灰度强调
    特征类型 特征 解释
    Gray Level Run Length Matrix (GLRLM) Features灰度游程矩阵 (GLRLM) 特征 (共16个) ShortRunEmphasis 是衡量短运行长度分布的指标,值越大表示运行长度越短,纹理纹理越精细。
    LongRunEmphasis 是对长行程分布的度量,值越大表示行程越长,结构纹理越粗糙。
    GrayLevelNonUniformity 灰度非均匀性
    GrayLevelNonUniformityNormalized 灰度非均匀归一化
    RunLengthNonUniformity 运行长度非均匀性
    RunLengthNonUniformityNormalized 运行长度非均匀归一化
    RunPercentage 运行百分比。通过获取 ROI 中的运行次数和体素数的比率来测量纹理的粗糙度
    GrayLevelVariance 灰度级方差
    RunVariance 运行方差
    RunEntropy 运行熵
    LowGrayLevelRunEmphasis 测量低灰度值的分布,值越高表示图像中低灰度值的集中度越高
    HighGrayLevelRunEmphasis 测量较高灰度值的分布,值越高表示图像中高灰度值的集中度越高
    ShortRunLowGrayLevelEmphasis 测量具有较低灰度值的较短运行长度的联合分布
    ShortRunHighGrayLevelEmphasis 测量具有较高灰度值的较短运行长度的联合分布
    LongRunLowGrayLevelEmphasis 测量具有较低灰度值的长运行长度的联合分布
    LongRunHighGrayLevelEmphasis 测量具有较高灰度值的长运行长度的联合分布
    特征类型 特征 解释
    Neighbouring Gray Tone Difference Matrix (NGTDM) Features相邻灰度差矩阵(共5个) (NGTDM) 特征 Coarseness 粗糙度
    Contrast 对比度
    Busyness 繁忙度。从像素到其邻居的变化的度量。繁忙度的高值表示“繁忙”图像,像素与其邻域之间的强度快速变化
    Complexity 复杂度
    Strength 强度
    特征类型 特征 解释
    Gray Level Dependence Matrix (GLDM) Features 灰度依赖矩阵 (GLDM) 特征(共14个) SmallDependenceEmphasis 小依赖分布的度量,值越大表示依赖越小,纹理越不均匀
    LargeDependenceEmphasis 大依赖分布的度量,值越大表示依赖越大,纹理越均匀
    GrayLevelNonUniformity 灰度非均匀性。测量图像中灰度级强度值的相似性
    DependenceNonUniformity 衡量整个图像中依赖关系的相似性,值越低表示图像中依赖关系之间的同质性越高
    DependenceNonUniformityNormalized 依赖非均匀归一化。衡量整个图像中依赖关系的相似性,值越低表示图像中依赖关系之间的同质性越高
    GrayLevelVariance 灰度方差
    DependenceVariance 依赖方差
    DependenceEntropy 依赖熵
    LowGrayLevelEmphasis 测量低灰度值的分布,值越大表示图像中低灰度值的集中度越高
    HighGrayLevelEmphasis 测量较高灰度值的分布,值越高表示图像中高灰度值的集中度越高
    SmallDependenceLowGrayLevelEmphasis 测量具有较低灰度值的小相关性的联合分布
    SmallDependenceHighGrayLevelEmphasis 测量具有较高灰度值的小相关性的联合分布
    LargeDependenceLowGrayLevelEmphasis 测量具有较低灰度值的大相关性的联合分布
    LargeDependenceHighGrayLevelEmphasis 测量具有较高灰度值的大相关性的联合分布

    注意:除了形状特征类外,其他特征都可以在原始图像和滤波后的图像上进行计算

    具体的指定方法:
    可以直接使用:enableAllFeatures( )启用所有类型
    也可以使用:enableFeatureClassByName(featureClass, enabled=True)启用你想用的类型
    例如:

    # 设置一阶特征
    extractor.enableFeatureClassByName('firstorder')
    # 设置只提取一阶特征的'Mean'和'Skewness'
    extractor.enableFeaturesByName(firstorder=['Mean', 'Skewness'])
    

    也可以直接通过特征类别来启用具体的特征,示例见:https://github.com/AIM-Harvard/pyradiomics/blob/master/examples/helloFeatureClass.py

  • 第三步可以对特征提取器进行一些设置。是否对原图归一化、是否重采样,原文档见:feature-extractor-level
  • 图像归一化:

  • normalize:默认为false。设置为true时进行图像归一化。
  • normalizeScale:确定图像归一化后的比例。默认为1。
  • removeOutliers:定义要从图像中删除的异常值。默认为0。
  • 图像/mask重采样:

  • resampledPixelSpacing:设置重采样时的体素大小。默认无。
  • interpolator:设置用于重采样的插值器。仅适用于重采样图像,sitkNearestNeighbor始终用于重采样掩码以保留标签值。可选的插值器:
  • padDistance:设置重采样期间裁剪肿瘤体时的体素补充数量。
  • 例如:

    settings = {}
    settings['![binWidth](https://img-blog.csdnimg.cn/c9b0896a5eea4eaf8217d0ed7f23e92b.png)
    '] = 25
    settings['resampledPixelSpacing'] = [3,3,3]  # [3,3,3] is an example for defining resampling (voxels with size 3x3x3mm)
    settings['interpolator'] = sitk.sitkBSpline
    

    还有其他的一些设置,根据自己的需要修改即可。
    最后举一个完整的例子,我要处理的是CT数据,所以我参考了pyradiomics的github源码中ct参数的设置,最后我的代码如下:

    from __future__ import print_function
    import six
    import os  # needed navigate the system to get the input data
    import numpy as np
    import radiomics
    import pandas as pd
    import SimpleITK as sitk
    from radiomics import featureextractor  # This module is used for interaction with pyradiomics
    import argparse
    def catch_features(imagePath,maskPath):
        if imagePath is None or maskPath is None:  # Something went wrong, in this case PyRadiomics will also log an error
            raise Exception('Error getting testcase!')  # Raise exception to prevent cells below from running in case of "run all"
        settings = {}
        settings['binWidth'] = 25  # 5
        settings['sigma'] = [3, 5]
        settings['Interpolator'] = sitk.sitkBSpline
        settings['resampledPixelSpacing'] = [1, 1, 1]  # 3,3,3
        settings['voxelArrayShift'] = 1000  # 300
        settings['normalize'] = True
        settings['normalizeScale'] = 100
        extractor = featureextractor.RadiomicsFeatureExtractor(**settings)
        #extractor = featureextractor.RadiomicsFeatureExtractor()
        print('Extraction parameters:\n\t', extractor.settings)
    
        extractor.enableImageTypeByName('LoG')
        extractor.enableImageTypeByName('Wavelet')
        extractor.enableAllFeatures()
        extractor.enableFeaturesByName(firstorder=['Energy', 'TotalEnergy', 'Entropy', 'Minimum', '10Percentile', '90Percentile', 'Maximum', 'Mean', 'Median', 'InterquartileRange', 'Range', 'MeanAbsoluteDeviation', 'RobustMeanAbsoluteDeviation', 'RootMeanSquared', 'StandardDeviation', 'Skewness', 'Kurtosis', 'Variance', 'Uniformity'])
        extractor.enableFeaturesByName(shape=['VoxelVolume', 'MeshVolume', 'SurfaceArea', 'SurfaceVolumeRatio', 'Compactness1', 'Compactness2', 'Sphericity', 'SphericalDisproportion','Maximum3DDiameter','Maximum2DDiameterSlice','Maximum2DDiameterColumn','Maximum2DDiameterRow', 'MajorAxisLength', 'MinorAxisLength', 'LeastAxisLength', 'Elongation', 'Flatness'])
    # 上边两句我将一阶特征和形状特征中的默认禁用的特征都手动启用,为了之后特征筛选
        print('Enabled filters:\n\t', extractor.enabledImagetypes)
        feature_cur = []
        feature_name = []
        result = extractor.execute(imagePath, maskPath, label=255)
            for key, value in six.iteritems(result):
            print('\t', key, ':', value)
            feature_name.append(key)
            feature_cur.append(value)
        print(len(feature_cur[37:]))
        name = feature_name[37:]
        name = np.array(name)
        '''
        flag=1
        if flag:
            name = np.array(feature_name)
            name_df = pd.DataFrame(name)
            writer = pd.ExcelWriter('key.xlsx')
            name_df.to_excel(writer)
            writer.save()
            flag = 0
        '''
        for i in range(len(feature_cur[37:])):
            #if type(feature_cur[i+22]) != type(feature_cur[30]):
            feature_cur[i+37] = float(feature_cur[i+37])
        return feature_cur[37:],name
    
    image_dir = r'D:\study\extract_radiomics_feature\linshi'
    mask_dir = r'D:\study\extract_radiomics_feature\linshi'
    patient_list = os.listdir(image_dir)
    save_file = np.empty(shape=[1,1051])
    id = []
    for patient in patient_list:
        print(patient)
        for file in os.listdir(os.path.join(image_dir,patient)):
            if file =='image.nii':
                imagePath = os.path.join(image_dir,patient,file)
            if file =='mask_gtv-1.nii':
                maskPath = os.path.join(mask_dir,patient,file)
        print(imagePath)
        print(maskPath)
        save_curdata,name = catch_features(imagePath,maskPath)
        save_curdata = np.array(save_curdata)
        save_curdata = save_curdata.reshape([1, 1051])
        id.append(patient.split('.')[0])
        np.concatenate((patient,save_curdata),axis=1)
        save_file = np.append(save_file,save_curdata,axis=0)
        print(save_file.shape)
    save_file = np.delete(save_file,0,0)
    #save_file = save_file.transpose()
    #print(save_file.shape)
    id_num = len(id)
    id = np.array(id)
    name_df = pd.DataFrame(save_file)
    name_df.index = id
    name_df.columns = name
    writer = pd.ExcelWriter('NSCLC-Radiomics-features.xlsx')
    name_df.to_excel(writer)
    writer.save()
    

    最后结果如下:

    来源:fighting_!

    物联沃分享整理
    物联沃-IOTWORD物联网 » 使用pyradiomics提取影像组学特征【详细】

    发表评论