Radiomics超入門：Intensity-volume histogram特徴

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

ROI 信号強度マスク内のボクセルの集合$${X_{gl}}$$の（累積）信号強度-体積ヒストグラム（Intensity-volume histogram, IVH）は、離散化された信号強度$${i}$$に対して、$${i}$$を含むボクセルの体積の割合$${\nu}$$の関係を表します[ElNaqa2009]。

画像モダリティによっては、IVH 特徴を計算するために$${X_{gl}}$$を離散化し、離散化ボクセルセット$${X_{d,gl}}$$を作成する必要があります。さらに、離散化された信号強度の全範囲$${G}$$と離散化間隔$${w_d}$$を知っておく必要がある場合があります。$${G}$$は、離散化の範囲を決定し、離散化の間隔$${w_d}$$は、離散化の信号強度間隔を決定するものです。これらがどのように特徴計算に影響するかについては、特徴ごとに説明していきます（離散化については、"Radiomics特徴計算のスキーム：Discretisation（離散化処理）"を参照ください）。

離散化された信号強度に対する離散化

画像モダリティによっては、デフォルトで、離散的な信号強度として較正されたボクセルを提供するものがあります。例えば、CTです。このようなモダリティの場合、離散化された ROI ボクセルセットは$${X_{d,gl}=X_{gl}}$$と考えることもでき、離散化は必ずしも必要ではありません。

再セグメンテーション範囲が指定されている場合、信号強度の範囲$${G}$$は再セグメンテーション範囲に等しくなります。範囲を指定しない場合、$${G=[min(X_{gl}),max(X_{gl})]}$$となります。このときの離散化間隔は $${w_d=1}$$ です。

連続値で較正された信号強度に対する離散化

PETのような連続値で較正された信号強度を持つ画像の場合は、定量的な信号強度情報を保持しながら IVH を決定するために、特殊な離散化が必要です。連続値の場合、Fixed bin size による離散化が推奨されます。

この方法では、最小信号強度$${X_{gl,min}}$$と最大信号強度$${X_{gl,max}}$$ およびヒストグラムのビン幅 $${w_b}$$を指定して離散化を行います。再セグメンテーション範囲が指定されている場合、$${X_{gl,min}}$$は再セグメンテーション範囲の下限に、$${X_{gl,max}}$$は上限に設定されます。そうでない場合は$${X_{gl,min}=min(X_{gl})}$$、$${X_{gl,max}=max(X_{gl})}$$です（すなわち、離散化前のROI 信号強度マスクにおける最小と最大）。ビン幅$${w_b}$$はモダリティに依存しますが、例えば、18F-FDG-PETでは0.10 SUVなどが設定されます。

次に、Fixed bin sizeによる離散化により、離散化済みボクセル集合$${X_d}$$が生成されますが、IVHでは、元の信号強度との関係を維持するために、ビン中心を元の信号強度に対応するように変換する必要があります。したがって、ビンの値$${X_d}$$をビン中心に対応する信号強度に置き換えます。

$$
X_{d,gl}=X_{gl,min}+(X_d−0.5)w_b
$$

これにより、範囲は、$${G=[X_{gl,min}+0.5{w_b}, {  }X_{gl,max}-0.5{  }w_b]}$$となります。この場合、離散化間隔はビン幅と一致することになり、$${ w_d=w_b}$$となります。

任意単位の信号強度に対する離散化

MRI画像のように、画像モダリティの中には、任意の信号強度を持つ画像を生成するものがあります。このような場合、$${N_g=1000}$$ビンの fixed bin number による離散化法が推奨されます。離散化されるビン幅は、$${w_b=(X_{gl,max}-X_{gl,min})/N_g}$$であり、範囲は$${X_{gl,max}=max(X_{gl})}$$、$${X_{gl,min}=min(X_{gl})}$$となります。

Fixed bin numberによる離散化は、ボクセル集合$${X_d∈\{1,2,...,N_g\}}$$を生成します。信号強度のキャリブレーションがないため、$${X_{d,gl}=X_d}$$となり、その結果、離散化間隔は$${w_d=1}$$となり、範囲は$${G=[1,N_g]}$$となります。

IVヒストグラムの計算

$${X_{d,gl}}$$を使って、体積分率（fractional volume）と強度分率（fractional intensity）を計算します。

画像スタックのボクセルは、一般的にすべて同じ次元を持つので、離散化された信号強度$${i}$$に対して体積分率$${\nu}$$は次のように定義されます。

$$
\nu_i = 1 - \frac {1} {N_v} \displaystyle\sum_{k=1}^{N_v} [X_{d,gl,k} < i]
$$

ここで、[...]は条件が真であれば1、そうでなければ0を返すアイバーソン記法のブラケットです。要するに、$${i}$$より小さい離散化信号強度を含むボクセルを数え、これをボクセルの総数で割り、これを全体割合「1」から差分することで信号強度$${i}$$に対する体積分率$${\nu_i}$$が求められます。

範囲$${G}$$における離散化された強度$${i}$$ の強度分率$${\gamma}$$ は次のように計算されます。

$$
\gamma_i = \frac {i - min(G)} {max(G) - min(G)}
$$

IVH特徴は3Dで計算されることが推奨されています。2Dでの画像スライスごとに計算して、平均化して集約するなどは推奨されません。

Volume at intensity fraction

ある強度分率における体積$${V_x}$$は、少なくとも x ％の強度分率$${\gamma}$$を有するときの最大の体積分率$${\nu}$$を意味します。

これは、概念的に類似している放射線治療計画で使用される線量-体積ヒストグラム（dose-volume histogram）とは異なり、$${V_{10}}$$では、少なくとも10 Gyの計画線量を受ける体積分率を意味します。[ElNaqa2009]は、$${V_{10}}と$${V_{90}}$$の両方を特徴として定義しています。

IBSIでは、これら2つの特徴をそれぞれ$${F_{ivh.V10}}$$および$${F_{ivh.V90}}$$と定義しています。

Intensity at volume fraction

ある体積分率における強度$${I_x}$$は、体積の最大 x %における最小の離散化信号強度$${i}$$です。[ElNaqa2009]は、I10とI90の両方を特徴として定義しています。

IBSIでは、これら2つの特徴量をそれぞれ$${F_{ivh.I10}}$$および$${F_{ivh.I90}}$$と定義しています。

Volume fraction difference between intensity fractions

この特徴は、2つの異なる強度分率における体積分率の差です。例えば$${V_{10}-V_{90}}$$など[ElNaqa2009]。

IBSIでは、この特徴を$${F_{ivh.V10minusV90}}$$と定義しています。

Intensity fraction difference between volume fractions

この特徴は、2つの異なる体積分率における強度分率の差です。強度分率は離散化されていることが前提になります[ElNaqa2009]。

IBSIでは、この特徴を$${F_{ivh.I10minusI90}}$$と定義しています。

Area under the IVH curve

IVH曲線下面積$${F_{ivh.auc}}$$は[VanVelden2011]によって定義されました。IVH曲線下面積は、一般的なAUCと同様に、台形則を用いてリーマン和を計算することで近似できます。なお、ROI 内に離散化された信号強度が 1 つしかない場合は、$${F_{ivh.auc}=0 }$$と定義されています。

2023/9 時点では、IVH曲線化面積の利用は非推奨であり、参考値も公開されていません。

実践

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

まずは、インスタンスを用意します。

ImagePlus ds_pair[] = TestDataLoader.digital_phantom1();
ImagePlus imp = ds_pair[0];
ImagePlus mask = ds_pair[1];
		
RadiomicsJ.targetLabel = 1;
RadiomicsJ.nBins = Utils.getNumOfBinsByMinMaxRange(imp, mask, RadiomicsJ.targetLabel);//digital phantom  of roi, max.
RadiomicsJ.debug = false;
		
IntensityVolumeHistogramFeatures ivhf = new IntensityVolumeHistogramFeatures(imp, mask, 1, 0);

体積分率と強度分率の分布を表示します。

ivhf.toString();//show IVH matrix

//出力
===  IVH result (nBins)  ===
discretised value:1.0, gamma:0.0, nu:1.0
discretised value:2.0, gamma:0.2, nu:0.32432432432432434
discretised value:3.0, gamma:0.4, nu:0.32432432432432434
discretised value:4.0, gamma:0.6, nu:0.31081081081081074
discretised value:5.0, gamma:0.8, nu:0.09459459459459452
discretised value:6.0, gamma:1.0, nu:0.09459459459459452
=== IVH result end ===

特徴量を計算してみます。

System.out.println(IntensityVolumeHistogramFeatureType.VolumeAtIntensityFraction10+":"+ivhf.calculate(IntensityVolumeHistogramFeatureType.VolumeAtIntensityFraction10.id()));//OK
System.out.println(IntensityVolumeHistogramFeatureType.VolumeAtIntensityFraction90+":"+ivhf.calculate(IntensityVolumeHistogramFeatureType.VolumeAtIntensityFraction90.id()));//OK
System.out.println(IntensityVolumeHistogramFeatureType.IntensityAtVolumeFraction10+":"+ivhf.calculate(IntensityVolumeHistogramFeatureType.IntensityAtVolumeFraction10.id()));//OK
System.out.println(IntensityVolumeHistogramFeatureType.IntensityAtVolumeFraction90+":"+ivhf.calculate(IntensityVolumeHistogramFeatureType.IntensityAtVolumeFraction90.id()));//OK
System.out.println(IntensityVolumeHistogramFeatureType.VolumeFractionDifferenceBetweenIntensityFractions+":"+ivhf.calculate(IntensityVolumeHistogramFeatureType.VolumeFractionDifferenceBetweenIntensityFractions.id()));//OK
System.out.println(IntensityVolumeHistogramFeatureType.IntensityFractionDifferenceBetweenVolumeFractions+":"+ivhf.calculate(IntensityVolumeHistogramFeatureType.IntensityFractionDifferenceBetweenVolumeFractions.id()));//OK
System.out.println(IntensityVolumeHistogramFeatureType.AreaUnderTheIVHCurve+":"+ivhf.calculate(IntensityVolumeHistogramFeatureType.AreaUnderTheIVHCurve.id()));//OK

//出力
VolumeAtIntensityFraction10:0.32432432432432434
VolumeAtIntensityFraction90:0.09459459459459452
IntensityAtVolumeFraction10:5.0
IntensityAtVolumeFraction90:2.0
VolumeFractionDifferenceBetweenIntensityFractions:0.22972972972972983
IntensityFractionDifferenceBetweenVolumeFractions:3.0
AreaUnderTheIVHCurve:0.26689189189189183

---

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編）