使用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()

最后结果如下:

python特征提取代码 pyradiomics特征提取_灰度值