Radiomics超入門：Intensity-based statistical feature（信号強度に基づく統計的特徴）

体験で学ぶ医用画像情報学超入門

信号強度に基づく統計的特徴は、関心領域（ROI）内のボクセルの信号強度がどのように分布しているかを説明するための統計量です。

この特徴は離散化を必要とせず、連続的な信号強度分布（離散化していないボクセルセット）から計算できます。

ROI マスクに含まれる $${N_v}$$ 個のボクセルの集合は次のように表すことができます。

$$
X_{gl}=\{{X_{gl,1}, X_{gl,2},...,X_{gl,Nv}\}} 
$$

この特徴は、3D ボリュームを用いて計算することが推奨されています。画像スライスごとに特徴を計算し、その後に平均化などで集約するアプローチは推奨されません。

この特徴は、一般的な統計量ですので、ここでは一気通貫で解説してきます。

信号強度の平均：Mean intensity

ROIに含まれるボクセルの信号強度の平均です。
式は以下の通りです。

$$
F_{stat.mean} = \frac {1} {N_v} \displaystyle\sum_{k=1}^{N_v} X_{gl,k}
$$

信号強度の分散：Intensity variance

ROIに含まれるボクセルの信号強度の分散です。
式は以下の通りです。

$$
F_{stat.var} = \frac {1} {N_v} \displaystyle\sum_{k=1}^{N_v} (X_{gl,k}-\mu)^2
$$

$${\mu}$$は$${F_{stat.mean}}$$です。

信号強度の歪度：Intensity skewness

ROIに含まれるボクセルの信号強度の歪度です。
式は以下の通りです。

$$
F_{stat.skew} = \frac {\frac {1} {N_v} \sum_{k=1}^{N_v} (X_{gl,k}-\mu)^3} {(\frac {1} {N_v} \sum_{k=1}^{N_v} (X_{gl,k}-\mu)^2)^{3/2}}
$$

$${\mu}$$は$${F_{stat.mean}}$$です。$${F_{stat.var}=0}$$のとき、$${F_{stat.skew}=0}$$となります。

信号強度の尖度：(Excess) intensity kurtosis

尖度は、強度分布$${X_{gl}}$$の尖り具合を表す尺度です。

$$
F_{stat.kurt} = \frac {\frac {1} {N_v} \sum_{k=1}^{N_v} (X_{gl,k}-\mu)^4} {(\frac {1} {N_v} \sum_{k=1}^{N_v} (X_{gl,k}-\mu)^2)^{2}}-3
$$

ここで$${\mu=F_{stat.mean}}$$です。Excessな尖度は、正規分布を仮定した場合に0が中心となるため、-3のフィッシャー補正をされます。これがない場合は、Excessでない尖度です（"Excess"かどうかは、式中で "-3"をするかしないかのみ）。

$${F_{stat.var}=0}$$のとき、 $${F_{stat.kurt}=0}$$です。

Median intensity

$${X_{gl}}$$の中央値です。

Minimum intensity

$${X_{gl}}$$の最小値です。

$$
F_{stat.min} = min(X_{gl})
$$

$${10^{th}}$$ intensity percentile

$${P_{10}}$$ は$${X_{gl}}$$の10パーセンタイル値です。$${P_{10}}$$ は$${F_{stat.min}}$$に代わるよりロバストな値として利用されることがあります。$${P_{5}}$$などでもよいのですが、IBSIでは$${P_{10}}$$を定義しています。

$${90^{th}}$$ intensity percentile

$${P_{90}}$$ は$${X_{gl}}$$の90パーセンタイル値です。$${P_{90}}$$ は$${F_{stat.max}}$$に代わるよりロバストな値として利用されることがあります。$${P_{95}}$$などでもよいのですが、IBSIでは$${P_{90}}$$を定義しています。

Maximum intensity

$${X_{gl}}$$の最大値です。

$$
F_{stat.max} = max(X_{gl})
$$

信号強度の四分位範囲：Intensity interquartile range

四分位範囲は次のように求められます。
$${F_{stat.iqr} = P_{75} - P_{25}}$$

$${P_{25}}$$と$${P_{75}}$$は、$${X_{gl}}$$の25パーセンタイルと75パーセンタイルです。

信号強度の範囲：Intensity range

信号強度の範囲です。
$${F_{stat.range}=max(X_{gl}) - min(X_{gl})}$$

平均絶対偏差：Intensity-based mean absolute deviation

平均絶対偏差は、$${X_{gl}}$$の平均からのばらつきを示す尺度です。

$$
F_{stat.mad} = \frac {1} {N_v} \displaystyle\sum_{k=1}^{N_v} |X_{gl,k}-\mu|
$$

$$\mu$$は$${F_{stat.mean}}$$です。

ロバストな平均絶対偏差：Intensity-based robust mean absolute deviation

平均絶対偏差は、外れ値の影響を受ける可能性があります。よりロバスト性を高めるために、外れ値を信号強度のボクセルの集合から除外することで、分布の中心に近いものに制限することができます。例えば、下限値を10パーセンタイル、上限値を90パーセンタイルにした場合を考えていきます。

これらの範囲に制限されたボクセルの集合は次のように表すことができます。

$$
X_{gl,10−90} = \{ x∈ X_{gl} | P_{10}(X_{gl})≤x≤P_{90}(X_{gl})\}
$$

ここで、$${X_{gl,10-90}}$$は、$${X_{gl}}$$の10パーセンタイルと90パーセンタイルで囲まれた区間にある信号強度のみをもつ、$${N_{v,10-90}≦Nv}$$個のボクセルの集合です。

そして、ロバスト平均絶対偏差は次のように式で表されます。

$$
F_{stat.rmad} = \frac {1} {N_{v,10-90}} \displaystyle\sum_{k=1}^{N_{v,10-90}} |X_{gl,10-90,k}- \overline{X}_{gl,10-90}|
$$

平均が$${X_{gl,10-90}}$$から計算されています。

中央絶対偏差：Intensity-based median absolute deviation

中央値絶対偏差は、平均絶対偏差とコンセプトは似ています。違いは、平均ではなく中央値からばらつきを計算することです。

$$
F_{stat.medad} = \frac {1} {N_v} \displaystyle\sum_{k=1}^{N_v} |X_{gl,k}-M|
$$

ここで、$${M}$$は$${F_{stat.median}}$$です。

変動係数：Intensity-based coefficient of variation

変動係数は、$${X_{gl}}$$のばらつきを表します。

$$
F_{stat.cov} = \frac {\sigma} {\mu}
$$

$$
\sigma = (F_{stat.var})^{1/2}
$$

$${\sigma}$$は標準偏差、$${\mu}$$は、$${F_{stat.mean}}$$です。

四分位分散係数：Intensity-based quartile coefficient of dispersion

四分位分散係数は、分散係数に代わるよりロバストな係数です。

$$
F_{stat.qcod} = \frac {P_{75}-P_{25}} {P_{75}+P_{25}}
$$

$${P_{25}}$$と$${P_{75}}$$はそれぞれ、$${X_{gl}}$$の25パーセンタイルと75パーセンタイルです。

エネルギー（エナジー）：Intensity-based energy

エネルギー [Aerts2014] は次の式で表されます。

$$
F_{stat.energy} = \displaystyle\sum_{k=1}^{N_v} X_{gl,k}^2
$$

二乗平均平方根：Root mean square intensity

別名で quadratic mean とも呼ばれます[Aerts2014]。

$$
F_{stat.rms} = \sqrt {\frac {\sum_{k=1}^{N_v} X_{gl,k}^2} {N_v}}
$$

実践

RadiomicsJを用いて、IBSIデジタルファントムからこれらの特徴を算出します。

ImagePlus ds_pair[] = TestDataLoader.digital_phantom1();
ImagePlus imp = ds_pair[0];
ImagePlus mask = ds_pair[1];
	
RadiomicsJ.targetLabel = 1;
		
IntensityBasedStatisticalFeatures ibsf = new IntensityBasedStatisticalFeatures(imp, mask, RadiomicsJ.targetLabel);
System.out.println(IntensityBasedStatisticalFeatureType.Mean+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.Mean.id()));//OK
System.out.println(IntensityBasedStatisticalFeatureType.Variance+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.Variance.id()));//OK
System.out.println(IntensityBasedStatisticalFeatureType.Skewness+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.Skewness.id()));//OK
System.out.println(IntensityBasedStatisticalFeatureType.Kurtosis+":"+(ibsf.calculate(IntensityBasedStatisticalFeatureType.Kurtosis.id())-3));//OK, 3.3.4 (Excess) intensity kurtosis
System.out.println(IntensityBasedStatisticalFeatureType.Median+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.Median.id()));//ok
System.out.println(IntensityBasedStatisticalFeatureType.Minimum+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.Minimum.id()));//OK
System.out.println(IntensityBasedStatisticalFeatureType.Percentile10+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.Percentile10.id()));//OK
System.out.println(IntensityBasedStatisticalFeatureType.Percentile90+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.Percentile90.id()));//OK
System.out.println(IntensityBasedStatisticalFeatureType.Maximum+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.Maximum.id()));//OK
System.out.println(IntensityBasedStatisticalFeatureType.Interquartile+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.Interquartile.id()));//OK
System.out.println(IntensityBasedStatisticalFeatureType.Range+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.Range.id()));//ok
System.out.println(IntensityBasedStatisticalFeatureType.MeanAbsoluteDeviation+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.MeanAbsoluteDeviation.id()));
System.out.println(IntensityBasedStatisticalFeatureType.RobustMeanAbsoluteDeviation+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.RobustMeanAbsoluteDeviation.id()));
System.out.println(IntensityBasedStatisticalFeatureType.MedianAbsoluteDeviation+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.MedianAbsoluteDeviation.id()));
System.out.println(IntensityBasedStatisticalFeatureType.CoefficientOfVariation+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.CoefficientOfVariation.id()));
System.out.println(IntensityBasedStatisticalFeatureType.QuartileCoefficientOfDispersion+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.QuartileCoefficientOfDispersion.id()));
System.out.println(IntensityBasedStatisticalFeatureType.Energy+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.Energy.id()));//ok
System.out.println(IntensityBasedStatisticalFeatureType.RootMeanSquared+":"+ibsf.calculate(IntensityBasedStatisticalFeatureType.RootMeanSquared.id()));//ok
		

//出力
Mean:2.1486486486486487
Variance:3.0454711468225
Skewness:1.0838207225574612
Kurtosis:-0.3546204806878417
Median:1.0
Minimum:1.0
Percentile10:1.0
Percentile90:4.0
Maximum:6.0
Interquartile:3.0
Range:5.0
MeanAbsoluteDeviation:1.552227903579255
RobustMeanAbsoluteDeviation:1.1138338159946533
MedianAbsoluteDeviation:1.1486486486486487
CoefficientOfVariation:0.8121978584917317
QuartileCoefficientOfDispersion:0.6
Energy:567.0
RootMeanSquared:2.768061083531605

---

RadiomicsJの引用はこちら

Kobayashi, T. RadiomicsJ: a library to compute radiomic features. Radiol Phys Technol 15, 255–263 (2022). https://doi.org/10.1007/s12194-022-00664-4

RadiomicsJのリンク

https://github.com/tatsunidas/RadiomicsJ

---

Stay visionary

体験で学ぶ医用画像処理プログラミング超入門（Java編）