3DSlicerによるJacobianとICEの求め方

OS: Ubuntu 20.04, 3DSlicer Version: 4.11.20210226

3D Slicerを用いてJacobian Determinant (jacobian)画像 およびInverse Consistency (ICE) 画像を計算する方法を記載します.

それぞれの画像から値を取る方法はまた別に記載します.

また, DIRの評価ですので, DIRがすでに終えている状態からの処理になります.

前準備

JacobianやICEを計算するために, RegistrationQA Extensionを使用します. 事前にExtension Managerからインストールしておいてください. 画面をキャプチャした時がたまたまなのかわかりませんが, アイコン画像が表示されませんでした(下のような感じです).

RegistrationQA Extension

前準備2つ目は, DIRで得られたDVF (3D SlicerでDIRした場合TransformNodeでも可)データをRegistrationQAの計算に使えるよう, VectorVolumeデータに変換します.

まずはDIR後のデータを一覧します. ReferenceをCT1, MovingをCT2としたDIRが順方向とし, deformedCT2が変更後のCT画像, ???_goとなっているものがDVFデータです. ICEを計算するのに逆方向も必要ですので逆方向で変形したdeformedCT1が変形画像, ???_returnとついているものがDVFデータです.

DIR後のデータ

順方向のDVFデータ, Transform_goと DVF_goどちらでもいいのですが, 今回名前の短く, モジュール名と重ならないDVF_goを使います.逆方向も同様です.

DVF_go, およびDVF_returnをTransformデータからVectorVolumeデータに変換します. 変換はTransformモジュールを使用します.

Transformモジュールに切り替えると, 下のようなパネルになります. その中で「Active Transform:」(図中上部赤丸)に変換したいDVFデータを指定します. 下の例ではDVF_goが選択されています. 逆方向の場合はここはDVF_returnとなります.

次に間の項目は無視して, Convertの見出しまで飛びます. Convertの設定項目で, Reference Volumeとoutput displacement fieldを入力します.

Reference Volumeは先に指定したDVFを得たときの(DIRしたときの)Reference画像を指定します. ここでは順方向のReference画像はCT1ですので, CT1を指定しています.

Output Displacement Fieldは変換後のデータ名となります. ここで注意ですが, 出力できるデータ形式は次の3つがあります.

Transformデータ(もとのDVF_goの形式と同じ), Volumeデータ, VectorVolumeデータ

必ず「Create new VectorVolume」もしくは「Create new VecorVolume as…」を選んでVectorVolumeに変換するようにしてください. 下の例では「Create new VecorVolume as…」を選んでDVF_go_scalarと言う名前になっています.

以上を設定したら, 後はApplyボタンを押すだけです. これも変換にしばらく時間がかかりますので気長に待ってください.

DVFデータの変換

Jacobian, ICEの計算

モジュールをRegistrationQAに切り替えます. モジュールの中の「Registration」の中に「Registration Quality Assurance」があるので, 選択してください. 事前にRegistrationQA Extensionを入れていればあると思います.

下がRegistration Quality Assuranceを選択した際のパネルです. 今回使うのは今まだ非アクティブですが, 下の方の「Jacobian Determinant」と「Inverse Consistency」ボタンです. 必要なデータが登録されるとアクティブになります.

RegistrationQAのパネル

次に, 使用するデータを登録していきます. この方法が直感的でなく, わからなくなるところです.

使用するデータを登録するにはパネル上部にある「Node」の一覧表で右クリックして出てくるメニューで「RegistrationQA」→「Assign node to:」と選んでいき, それぞれ適切なデータを適切な種類にassignしなくてはなりません.

下の画像では, CT1をReference Imageであるとして登録しました. その他Image系では, Moving image, Forward warped Image, Backward warped imageを登録します. 今回順方向ではCT2をmovingにしましたので, 同様な方法でCT2をMoving imageに登録しました.

データのassign

Forward warped imageは順方向での変形後の画像,Backward warped imageは逆方向での変形後の画像を指定すると思われます. Reference, Movingは順方向(もしくは逆方向)でのセットなのに, 変形後の画像は順方向, 逆方向両方を登録できることになります.

さらに 混乱を招くところですが, パネル一番上の項目にParameter Setというのがありますが, その下のプルダウン(下の図では「RegistrationQAParameters」を広げると「RegistrationQAParameters_Backword」というのが選択しにあります. おそらく, ReferenceとMovingを逆にして登録すればいいと思いますが, 再度変形後の画像を登録するのか, 登録知るのならどっちがforwardになるのか, よくわかりません.

Vector FieldへDVFを登録する際ですが, 方法はCTと変わらず登録するデータを一覧表から選んで右クリックし, 「Assign node to」で「Vector Field」を選んで登録します. Vector Fieldには「Forward」と「Backward」があり(!?) それぞれに, 前準備で変換したVectorVolumeのデータを登録します. Transformデータも登録できてしまうのですが, 計算はされません(しかもなんのエラーも出ません).

今回のようにJacobianとICEを求める場合は, Vector Fieldをassignすれば計算できます. CT画像も必要ありません. JacobianはForwardかBackwardのどちらかさえ登録すればアクティブになります.

データ登録がされると, パネル下方のTableの項目にAssignしたデータ名が表示されます. この表の続きには, JacobianやICEを計算したあとにそれぞれの値(平均や最大, 最小など)が入る項目もあります.

必要なデータの登録がすんだら「Jacobian Determinant」「Inverse Consistency」ボタンが有効になっているはずですので, それぞれをクリックして計算させてください. 一つ一つ計算させることになりますが, この計算もまた少し時間がかかります. 計算が終わるとそれぞれ計算された画像が表示されます. (下図はJacobianを表示)

また, Tableの中に「Jacobian (for)」と「Jacobian (Back)」があり, その並びに数値が入っていると思います.

さらに, データ一覧(パネル上部のIOの表) に, Vector FeildにassignしたDVFデータの尾末に 「-jacobian」や「-invConsist」とついたデータが出来ていると思います. これでJacobianやICEのデータが画像データとして作成されたことになります.

ここまでで, JacobianとICEを画像データとして作成できました. ただし, まだストラクチャに対して個別にそれぞれの値をとることが出来ていません. 次の投稿ではSegment statisticモジュールを使って数値を読む方法を記載します.

Plastimatch のExtension(SlicerRT)でDaemonsができるよう改良しました

OS: Ubuntu 20.04(Linux版はExtensionがInstall出来ない不具合が見つかりました), Windows10

3DSlicer Version: 4.13

Extension SlicerRT1.0.0

今回は具体的な話はさておき, 結果のみの記載しています.

以前PlastimatchのRegistrationのためのExtensionで, B-Splineしか出来ないと書きました. せっかくのポテンシャルがもったいないと思い, SlicerRT Extension (現在PlastimatchのExtensionはSlicerRTに統合されています) のソースコードを覗いたのですが, 実際やっていることはPlastimatchのコマンドの文字列(オプション含む)を作ってPlastimatchに送っているだけでした.

そのため, じゃあその文字列をDaemonsでDIRできるよう変更してやればDaemonsでDIRができるだろうとなりました.

ただし, B-Splineで有効なパラメータとDaemonsで有効なパラメータが若干違うのでその判別ができるようにする必要もありました. ただ, 主なパラメータはB-Splineのほうが多く, B-SplineのパラメータをDaemonsでは無効にすることで行けそうでした.

最終的に,

  1. 今のExtensionの画面を最大限利用して, Daemonsへの切り替えができるようにして, Daemonsに指定できないパラメータはB-Splineだけであることを明示する(値は入れれるけど, いれても無効になる)ようにパラメータ画面を変更する.

2. Daemons用のコマンド文字列を作るようソースをいじる.

3. いじったExtensionのソースも一緒にビルドしてSlicerRT Extensionをデフォルトで入っている状態にする.

としました. ソースコードはGITHUBから最新版をダウンロードしたのでSlicerは4.13になっています.

出来上がったのが, 下の画像です. Algorithmという選択項目が増えて, B-SplineかDaemonsを選べます.(あ, つづりが間違っている)

パラメータの設定が下のようになっています. (Only for B-Spline)とある項目がB-Splineの際にのみ有効な設定項目です.

Daemonsで実行してもちゃんと結果が出ます. 本当にDaemonsで処理されている? と疑問もあるかとおもいますが, 基本的に文字列作っているだけなので, そのコマンド文字列で動いているということはちゃんとDaemonsでDIRされているのだと思います. B-SplineとDaemonsの違いがどれだけかはまだ確認していませんが, 今後の予定です.

3D SlicerのSegment Statisticsの使い方

3D Slicerである点の画素値などを得る場合, マウスのある点の値はSlicerのウィンドウ左下にあるパネルで確認できます.

ある領域の画素値の平均などは, Segment Statistics モジュールを使うことで得ることができきます. 放射線治療のデータで言うと, ある輪郭への線量の平均値などを得るときにもこの Segment Statisticsモジュールが使えます. 単純にセグメントごとに対象データの統計的な値を出してくれるので何かと重宝するのではと思います.

今回はSegment Statisticsモジュールの使い方について記載します.

前準備 Labelmapを作る

Segment Statisticsで使用する”Segment”はBinary Labelmap のことになります. Segment Editor で作成したデータは基本的にvtkMRMLSegmentationNodeなので, labelmapではありません. そのため, SegmentをLabelmapに変換しておく必要があります.

まずDataモジュールで, 対象としたいSegment(ここではRTStructure)を選びます. 目玉マークが目が開いている(表示されている)ものがLabelmapに変換されますので, 必要な輪郭を表示, いらないものを非表示にしておいてください.

その後, SegmentationNodeで右クリックすると下図のようにメニューが出ますので,「Export visible segments to binary labelmap」を選択してください. 少し時間がかかる場合がありますが, これだけで変換は終わりになります. Dataモジュールの一覧に選択したvtkMRMLSgmentationNodeの名前(下図では「1:RTSTRUCT:ANON」)に「-label」と末尾についたvtkMRMLLabelmapVolumeNode (1:RTSTRUCT:ANON-label) が出来ていることを確認してください. ただし, もとのsegmentのようにpatientやGTVなどには別れていません.

DataモジュールでSegmentationNOでをLabelmapに変換

Segment Statisticsでのデータのとり方

Segment Statisticsで得られるのは統計的な値です. つまり平均値や最大値, 最小値, 標準偏差などになります.

ここで重要なのは, 出力されるデータの単位は対象とするデータ(vtkMRMLScalarVolumeNode)によることになります.

つまり, 線量分布を保存するScalarVolumeNodeを対象とすれば, 値はGyなど線量を示す値で出されます. 長さを表すScalarVolumeNodeであれば, mmやmなどの単位となります. そのため, 結果には単位が表記されません. 単位については対象データの方でしっかり確認しておく必要があります(とくに, mmなのかcmなのか, など同じものだけど桁が違う場合もありえます).

Segment Statisticsで対象にできるデータは, vtkMRMLScalarVolumeNode, vtkMRMLBinaryLabelmapNode, vtkMRMLClosedSurfaceNodeとなります. それぞれいられる出力は選択できるので,試して確認してみてください.

今回は, 治療計画のデータからSegmentごとの線量値の平均値などを取りたいと思います.

まず, Segment Statisticsモジュールですが, モジュールのプルダウンのなかの「Quantification」の中にあります.

Segment Statisticsモジュールを開くと, 下のような選択項目があります.

Segment Statisticsモジュールのパネル

ここで「Segmentation」に指定するのがLabelmapです (上図では名前の末尾に−lableがついてないですが). Labelmapにした輪郭ごとに統計的データが計算されることになります.

Scalar Volumeがまとめるデータが含まれるデータとなります.ここではRTDOSEを指定しています. つまり, 線量のデータを保存したデータなので, 輪郭ごとの線量値の統計データが出力されることになります.

Output tableは出力結果を表示する表データの指定です. 基本的には「Create new table」で問題ないと思います.

基本的には, これらを指定すれば後は「Apply」ボタンをクリックするだけです.

結果は下のように画像の下の表として出力されます. 下図の結果では, 出力する値を絞って選択してあります.

SegmentStatisticsの結果表示

出力される値の種類を選びたいときは2つ上の画像(Segment Statisticsモジュールのパネル)にある「Advanced」から, 「Enabled segment statistics plugins:」で対象Volumeに合わせてチェックを外したり, その横の「Options」ボタンをクリックすると, 下図のように計算する値を選択できます.

下図はSclar VolumeのOptionsですが, Labelmapやclosed surfaceのOptionもそれぞれに応じた出力の種類が列挙されており選択できる.

Scalar Volume statisticsのOptionsボタンをクリックすることで表示されるウィンドウ.

Labelmap statisticsはSegmentationで指定したLabelmapから, Scalar Volume statisticsはScalar Volumeで指定したデータから, Closed surface statisticsはLabelmapで閉じた領域になっているものから計算されると思います.

結果の表のデータは, Slicer本体のSave機能でテキスト形式やCSV形式で保存できます.

3DSlicerによるDice係数, Hausdorff distanceの求め方

Ubuntu 20.04, 3DSlicer Version: 4.11.20210226 (Extension SlicerRT)

3DSlicerでDice係数, HDを求める方法を説明します. ただし, 計算するには3DSlicerのExtensionのSlicerRTが必要となります. インストールの方法は「3DSlicerによるDIR: Install」をご覧ください.

前提として, 2つの輪郭セットがSceneに読みこまれている, もしくはDIR後であればhardeningされた状態であることが必要です.

Dice係数, Housdorff Distanceを求める

両係数の計算にはSlicerRTの「Segment comparison」を使用します. 下図のよう「Radiotherapy」モジュールから「Segment comparison」を選択してください.

「Segment comparison」モジュールの画面は下のようになります. ここで操作するのは「Input」の項目のみとなります.

「Input」のしたの「Reference segment」が基準となり, 「Compare segment」が比較するsegmentになります. ここでは,propagationの説明で作成した「ct1roi」と「ct2roi Copy」をReference Segment, Compare segmentとしました.

ここで指定したSegment(vtkMRMLSegmentNodeとなるもの)は輪郭のセットを表すので, 治療計画データを開いていれば複数の輪郭が含まれています. そのため, さらにどの輪郭を対象とするか選ぶ必要があります. Reference segment, Compare segmentそれぞれの下にさらに「segment」という選択項目がありますので, そこで対象の輪郭を指定します. 今回データは適当に開いたので, 一番似かよってそうなそれぞれのpatient輪郭を計算対象として選びました.

この計算対象輪郭は何も別なvtkMRMLSegmentNodeである必要もなく, 同じvtkMRMLSegmentNode(例えばct1roi)の中から別な輪郭を選んだり, 同じ輪郭を選んだりすることも出来ます(もちろんその際DICE係数は1となります).

計算対象の輪郭が設定し終えたら, 後は「Compute Housdorff distance」もしくは「Compute Dice metrics」ボタンをクリックするだけです. 計算はしばらくかかります. どちらかを計算すれば, もう片方はすぐに終わります.

結果は, それぞれの表の一番上の値となります. 下の図ではHDは62.1898, Diceは0.893319です.

これらの値の正しさは確認していないので, 使う際には各自で確認してください.

3D Slicer でのPropagation

OS: Ubuntu 20.04, 3DSlicer Version: 4.11.20210226

3D Slicerを用いた輪郭のpropagationのやり方です.

すでにDIRがなされている状態を想定しています. 説明の画像ではPlastimatch によるB-SplineによるDIRを実施し, その結果をB-SplineパラメータとDVFををMRMLで保存した状態からスタートしています. それぞれvtkMRMLBSplineTransformNodeとvtkMRMLGridTransformNodeです. 2つともある必要はなく, どちらかあればpropagationは可能です.

下準備: 変形したい輪郭のコピーを作成

Propagationする際, Transformsモジュールを使用します. そこで輪郭を変形すると上書きになるので予めもとの輪郭をコピーして変形用の輪郭を作っておきます.

まず, Data モジュール(下図の上の赤丸)に切り替えてコピーしたい輪郭を探します. 下図の画像ではct2roiが変形したい輪郭になります.

変形したい輪郭を選び右クリックするとメニューが出ますので, 「Clone」(下図の下の赤丸)を選択してください. それでコピーが出来ます.

右クリックメニューからCloneを選ぶ

確認すると下図で「ct2roi Copy」(下図赤丸)というvtkMRMLSegmentationNodeが出来ています.

下準備は完了ですが, もし上書きしても構わないというのであれば, この下準備は飛ばしても構いません.

「ct2roi Copy」が出来ている

Transformを用いた変形

実際に変形していきます. 輪郭を変形するためにはTransformsモジュールを使用しますのでモジュールをTransformsモジュールに変更してください.

次に, DIRした際に得られたDVFなどを保存したデータ(vtkMRMLTransformNodeの一つ)を「Active Transform」の欄で選択します. DIRした際にMRMLとして変形を保存していればプルダウンメニューに出てくると思います. Affine変換など線形な変換であれば下図の「Edit」のところを展開すれば手入力で作成できます. B-Splineも作れるようですが, ???な感じです.

ここでは, Active TransformにPlastimatchでデフォルトの名前のまま出力したDVF「Output vector field (MRML)」を指定しています.

Informationは変形のサマリー, DisplayはDVFを矢印やグリッドで表示したりできますが, これも関係ないので飛ばします.

輪郭の変形にはTramsformモジュールの「Apply transform」の設定で行います. Transformableのリストは, 現在Sceneにある変形可能なNodeの(データの)一覧です. この中から変形したい輪郭を指定します.

変形対象の指定は「Transformable」のリストで選択し, 右向き矢印で「Transformed」のリストに移します. 今回の場合変形対象がCloneで作成した「ct2roi Copy」ですので, 選択し右のTransformedのリストに移します.

移した時点で輪郭は変形されます. 画像を表示している画面で右矢印と左矢印で対象を行き来せれば変形しているのが解ると思います.

この時点で変形されているのですが, まだ変形後の輪郭は保存されていません. そのおかげで何度もやり直したりできるのですが, 変形後に解析・評価などする場合では見れるだけだと困ります. そこで変形後のデータを保存する方法ですが, 「Harden transform」で変形をFixしてしまいます. つまり, 今変形している輪郭を上書き保存することになります(そのためコピーを作成しました).

Harden transformの方法は, 保存したい対象をTransformedリストでクリックして選択(水色の背景になる)し, TransformableのリストとTransformedの間にある矢印の下, マス目がずれたアイコン(下図 下の赤丸) をクリックするだけです. アイコンをクリックするとTransformedにあるデータを現在の変形された状態に保存します. うまく保存されると選択された変形対象だけがTransformedリストから消え, Transformableのリストに戻ります.

以上が輪郭のpropagationの方法になります. 今回輪郭を対象にしましたが, 変形できるデータであれば(Transformableリストに上がるデータであれば)同様な操作で変形が出来ます.

注意は, hardenしても3DSlicer上で保存されている(データとなっている)だけで, ファイルには保存されていませんので, 途中経過や結果などを保存したい場合は, 3DSlicerでScene全体や変形した輪郭だけでもファイルで保存する必要がありますので忘れずに.

3DSlicerでCTと輪郭の組み合わせは, おそらく他の商用DIRソフトと比べゆるいです. DICOM画像を読み込む際には組み合わせが違うと言われたりする場合もありますが, Sceneに保存されている輪郭はすべてどの画像上にでも表示できます(というか表示されてしまう). それだけCT-輪郭の関係の自由度が高く, やりたいことができる(手間はかかるけど)ことになります.

3D Slicerを用いた画像のDOCIM化

OS: Ubuntu 20.04, 3DSlicer Version: 4.11.20210226

3D Slicerを用いると, (たぶん読み込める画像なら)どのような画像でもDICOM形式に変換できます.

ここでは, DICOM化の方法を記載しておきます.

Fileの読み込み

基本的に何でもいいのですが, JPEGなどであると1スライスの画像とみなされます. 試していませんが, ファイルを開く際にフォルダを指定できるので, フォルダを指定すると3Dデータになるかもしれません.

今回は本当に一般的な画像がDICOM化できることを示すためにJPEGでやってみました.

3DSlicerのメニューの「File」から「Add Data」を選ぶか, メニューの下の一番左のDataアイコンをクリックすると下のウィンドウが出てきます. 今回, JPEGファイル1つなので「Choose File(s) to Add」を選んでファイルを選択しました. 「Choose File(s) to Add」をクリックするとさらにファイル選択のウィンドウが出てきます.

Add Dataのウィンドウ

ファイルを選択し, 開くと下の図のような画面です. ファイルは研究室のロゴマークのファイルにしました. 1枚の2D画像ですので, Sag, Corの画像は線です.

ロゴファイルをオープンしたところ

DICOM化

もう早速DIROM化です. まず下図のように, Modules横のプルダウンを「Data」にしてください. そうすると左パネルに「Subject hierarchy」タブになっていると思います(なっていなければ選択してください).

DICOMーRTなどを開くと輪郭などがツリー状に表示されますが, これまでの手順でファイルを開くと開いたファイルのファイル名が1つだけ表示されていると思います. (しっかりCT画像などと同じvtkMRMLVolumeNodeの1種になっているのがわかります)

次に, ファイル名を選択し右クリックします. 出てきたメニューに「Export to DICOM…」という項目があるので選択します.

そうすると下の警告が出ます.

DICOM化の注意

警告に書いているとおりですが, DICOM(化)には親となる患者情報が必要でOKを選ぶとデフォルトでつくられる. もし自分で作りたければ一度Cancelしてsubject hierarchyのところで, 「Export to DICOM…」を選んだメニューから「Create new subject」を選んで作成してくださいとのことです. Subjectを作った後にsubject hierarchy上でファイル名をドラッグ&ドロップすれば,作成したSubjectの配下に移動できます.

DICOMの情報は次のところで編集ができるのでそのままOKで問題ないと思います.

そのままOKすると下図右のようなウィンドウが出ます. ここでDICOMのヘッダに含めるMeta情報を入力していきます. とりあえずDICOM化したいだけであれば, なにも変更しなくてもDICOM化できます.

ついでに, 下図左のsubject hierarchyのパネルには, 自動にNewPatientとNewStudyが出来てその下にロゴファイルの名前が移動していると思います.

DICOM化したデータをそのまま3DSlicerのDICOMモジュールが管理するようにしたければ(3DSlicerのDICOM画像管理画面で出てくるようにしたければ), 「Export」をクリックしてください. DICOM画像のリストに適当につけられた患者名となってあるのが確認できます.

ファイルとして出力したい場合は, 下図右のウィンドウ下部の「Export to folder」にチェックを入れて, その右に出力先のフォルダを指定してください.

DICOM情報の入力

3DSlicerのDICOM管理に出力した場合, 以降は通常のDICOM画像と同様なので, もしファイルとしてExportしたければDICOMの管理画面で可能となります.

確認

最後にちゃんと出来たか確認ですが, ImageJを使って確認してみました. Exportしたファイルは1枚の画像ですので指定したフォルダには「IMG001.dcm」というファイルが出来ていました. そのファイルをImageJにドラッグ&ドロップすると白黒化されたロゴが表示されました. DICOMですのでグレーレベルの画像になったようです. ちゃんとDICOMかどうか確認するため, Image infoを表示させると, ちゃんとタグとともにタグ情報も表示されました.

さらにちゃんとしているのがUIDがついていること. 正解で一つしかない一意な番号がつけられるUIDをどのように決めているのか少し謎です. とはいえ, ちゃんとDICOM化が出来ました.

3D Slicerのデータ構造(MRML)

この内容は, 3D SlicerのDeveloper Guideの MRML Overview (https://slicer.readthedocs.io/en/latest/developer_guide/mrml_overview.html)の内容を主観を交えて解説しています. Developer向けガイドですので, 内部構造的な表現だったり, すこし3D Slicerでの表示と異なる用語になっていたりしますが, ご了承ください.

3D Slicerで画像などを開くと単に画像データを読み込んだというだけでなく, その他諸々のデータを一括して管理するSceneに読み込んだとして扱われます. 読み込んだデータもしくは処理により出来たデータは, Sceneの中にいくつかの種類に別れたNodeに適宜分別され管理されます.治療のデータ(DICOM-RTなど)は治療計画時の輪郭データであったり, 線量評価点など各種Nodeに分かれて読み込まれます. そのため, 何か(プログラミングなど)しようとした際にそれらNodeの理解が必要になります.

Sceneのデータ構造はMedical Reality Modeling Language (MRML) と言われ, Visualization Toolkitや Insight Toolkitなどで開発されたライブラリだそうです. そのため, Sceneの各NodeにつけられるIDなどには接頭語としてvtkMRML などが付きます(SlicerのDataモジュールで確認可能).

MRMLに含まれるNodeには大きく7種類のNodeがあります. 順に説明します.

  • DataNodes
  • DisplayNodes
  • StorageNodes
  • ViewNodes
  • PlotNodes
  • SubjectHierarchyNode
  • SequenceNode

DataNodes

DataNodesはもっとっも目にするNodeとなると思います. 様々な画像データを含むデータを保存,管理します. DataNodesの代表的な派生Nodeとして下のものがあります.

  • Volume (vtkMRMLVolumeおよび派生クラス)
    • 一般的な画像データ, つまり, Volumeデータなどマトリックス状にスカラー値をもつボクセルによって作られるデータを保存, 管理します. CTやMRI, 一般撮影画像などの2D画像はこのクラスで管理されます. segmentationされた領域のデータも, ボリュームタイプの3Dデータ(対象領域の画素にのみ0以外の値になり, その他は0となるボリュームデータ)であればこのタイプになります. そのほか, スカラー値を持つボクセルと言いましたが, ボクセルに複数の数値データをもつ多次元データやTensorなども含まれたりします. それらには, RGBデータをもつデータ, Deformation Vector Fieldなども含まれます. よく似た4次元データとして時系列データがありますが, 時系列のデータは別途SequenceNodeがあるのでそちらで扱われます.
  • Model (vtkMRMLModelNode)
    • 体表面を形どったサーフフェイスメッシュや, 内部もメッシュとして分割されたボリュームメッシュを保存・管理します. 治療の輪郭は次のSegmentationとして扱われます.
  • Segmentation (vtkMRMLSegmentationNode)
    • 名前の通りSegmentationされたデータを扱います. 治療の輪郭やラベリングされた領域などがこのNodeになります. データ形式は複数の形式で保管することができ, メッシュやBinaryLabelmapとしても保管されたりします. プログラムコード的にはクラス内部にVolumeやModelのクラス変数を持っているのではと思います.
  • Markups (vtkMRMLMarkupsNodeとその派生クラス)
    • 幾何学的な位置や直線, 曲線などを管理します. Fiducial markerや計測のための位置や領域を表す座標値, アノテーションを付す位置の情報などとなります.
  • Transform(vtkMRMLTransformNode)
    • 画像を変形させるための情報を保管します. 線形, 非線形問わず管理されます. 一応, Linear, B-spline, grid用の3種類の派生クラスがあります.
  • Text(vtkMRMLTextNode)
    • 設定ファイルや説明テキストなどのテキスト情報
  • Table(vtkMRMLTableNode)
    • 表敬式のデータを管理. 計算結果や計測結果などを表示するためのデータ.

DisplayNodes (vtkMRMLDisplayNode)

DIsplayNodeは表示用設定を保存します. 表示のためのデータ自体はないので, DataNodeと連携してデータを表示します.

表示設定ですので, 一つのデータを様々な設定で表示することもありえます. その際は複数のDisplayNodeをもつことになります.

逆にDataNodeは表示設定がないため, とりあえず表示もDisplayNodeに頼ります. そのため, DataNodeが1つ作られると, 必ずデフォルトでDisplayNodeが作られます.

このNodeは1種類ですので, vtkMRMLDisplayNodeとすでにNode IDがついています.

StorageNodes(vtkMRMLStorageNode)

ファイルの保存のしかたを管理するクラス. 保存する情報としては, ファイル名やファイル形式, 座標系など.

ViewNode(vtkMRMLAbstractViewNode)

画像表示のための周辺の設定を管理する. Abstractとついているようにプログラミングでは派生クラスを使うのが通常になる. 派生クラスにはしたの種類がある

  • vtkMRMLCameraNode
    • 3D viewのカメラ位置・向きなどの設定を管理する
  • vtkMRMLClipModelNode
    • ModelNodeをクリップするスライスの設定を管理する
  • vtkMRMLCrosshairNode
    • カーソル位置を示すクロスヘアラインの設定を管理する
  • vtkMRMLLayoutNode
    • 画面レイアウト(分割方法やどのパネルに何を表示するかなど)を管理する
  • vtkMRMLInteractionNode
    • マウス操作などの設定を管理する
  • vtkMRMLSelectionNode
    • Scene全体の情報を管理する. 全体で使用される単位や現在のSceneの状態など.

PlotNodes

TableNodeのデータをどのようにグラフで表示するかを設定する. 以下の3種類がある.

  • PlotSeriesNode
  • PlotViewNode
  • PlotChartNode

SubjectHierarchyNode(vtkMRMLSubjectHierarchyNode)

DataNodeの各ノードを階層形式で管理する. おそらくData moduleを表示した際に表示されるデータのツリー構造を示すためのnode. Sceneに一つのみ作られる.

SequenceNode

DataNodeのタイムシーケンスのリストや多次元データをリストとして管理する.

SequenceBrowserNodeがどのDataNodeを今表示するか(しているか)を管理する.

3DSlicerによるDIR: Plastimatch v.s. Elastix

OS: Ubuntu 20.04

3DSlicer Version: 4.11.20210226

Extension SlicerElastix, SlicerRT1.0.0

Plastimatch v.s. Elastixとしましたが, あくまで1例試して比較した際の個人の意見です.

それぞれの処理の最適化もなにもしていません. 基本デフォルトの設定での処理です.

まず, 処理時間

Ubuntu 20.04, CPU Ryzen7-3700, メモリ 16MB, GPU GeforceRTX 3080 での処理です.

対象画像はDIR-DBにある頭頸部の#12の症例. CT1(実施時にRenameしています)が180枚, CT2が176枚のCT画像です.

PlastimatchはRigitを含めStage3まですべてEnableにして実施. 各ステージの設定はデフォルトのまま. Cost functionはMSE, HardwareはGPUを選択したのですが, リソースモニターではCPUしか動いていませんでした. CUDAは入れてありますが, GPU使うための設定はなにもしていないので, なんらか設定が必要なのかもしれません.

ElastixはPresetをgeneric(all)を選んだだけです.

処理時間

Plastimatch 約45秒

Elastix 約5分

処理時間は, どのような処理をしているかにも大きく依存するので, それらを統一していない処理で比べてもなんとも言えませんが, すこしPlastimatchのほうがいいのかなという個人的な感想です.

結果画像

画像上段左がFixed Image (CT1), 下段左がPlastimatch, 右がElastixです. 少しわかりにくかったのでCT1とCT2の輪郭を一緒に載せています. CT1のピンクの体輪郭(CT2のものもあるのでCT1の画像を参照してください)にどれだけあっているかを見てもらえばいいですが, Plastimatchは左肩あたりが少しElastixより違いが大きいし, すこしボケたような画像になっています.とはいえ僅かな差かと思います.

ついで, もう少し上のスライスで歯冠のアーティファクトがあるスライス前後になると, Elastixは極端に間違った変形となっています. Plastimatchはうまくアーティファクトを対象外にしているのか変な変形はありません. ただし, 顎の形などを見ると変形前の輪郭に近いので, アーティファクトのある部分は変形しないという処理なのかもしれません.

所見

今回の処理では若干Plastimatchのほうが良い結果でした. ただし, 処理自体を何も調整していない状態での比較なので, 結果画像の比較に大きな意味はあまりないようにも思います. おそらくElastixにも対象とする画素値の範囲を設定したりできるのではと思いますので, 今回アーティファクトが原因で変な変形になったとしても処理を最適化していけば良い結果を得れるのではと思います.

処理時間のところで書いたPlastimatchが少し良さそうとしたのは, おもに処理時間の結果とその将来性によります. 今回どちらもCPUのみでの処理で5倍以上の差が出ており, さらにPlastimatchのほうがGPUを使う設定があります. ちゃんとGPUの設定すればもっと早く処理できるのではと思います. 結果画像も大差ないうえでは, 処理時間とその将来性からすこしPlastimatchのほうがいいのではないかと思います.

本来, PlastimatchもDemons法でのDIRができたりするのですが, 3DSlicerのExtensionになった時点で, その処理が選べなくなります. Elastixも同様で, Extensionになることで, 本来できる機能, 処理の1部しか使用できなくなります. そういった意味では, どちらもコマンドラインで処理すればいいのですが, それだとここで3DSlicerのextensionとして説明してきた意味が… Slicerの中で実施できることですぐに結果画像を確認できるメリットや, 他のモジュールを使って評価や検証ができたりします. その意味ではExtensionであることの意味は大きいと思います.

ElastixのExtensionは初心者か上級者向けと書きました. それは, そのままでもある程度使えるようにしてくれていることと, さらに, 処理の変更や調整が設定を書いているファイルを編集することで可能になるので,上級者はさらに深く使えるところにあります.

Plastimatchは中級者クラスの人が, B-Splineで少しステージ数や解像度の調整など微調整ができるので, そのくらいでいいのであれば使いやすいExtensionだと思います. ただ, それ以上の変更をしたい場合は, たぶんこのExtensionではできないので, Elastixやコマンドラインでの実行にせざる得なくなると思います.

どっち使うかは自分のやりたいことにあった方を選ぶことになると思います.

3D SlicerによるDIR: Elastix

OS: Ubuntu 20.04

3DSlicer Version: 4.11.20210226

Extension SlicerElastix

Elastixも3D Slicerで使えるDIRソフトです. Install編で説明したようにExtensionを入れれば使えるようになります.

3D slicerでの設定は, ごく簡単な設定になり, 細かな設定は設定ファイルに書かれておりそれを選んで実行することになります.

ElastixでDIRをしようとした場合, Moduleは「Registration」->「General Registration(Elastix)」を選びます.

下の図が設定項目になります. Plastimatchと比べここで設定できる項目は本当にシンプルです. ElastixのExtensionはplastimatchと比べ一見簡単そうですが, 処理の調整などをしようとするとファイルの記載を変更しなければならないため, とりあえずDIRさせたい, もしくはある程度DIRやプログラミングの知識がある人がガッツリ調整してDIRしたい, どちらか両極端な方向けかと思います. 順に説明していきます.

Elastixのモジュールインターフェイス
  • Reload&Test

この見出しの項目は, このElastix moduleのインターフェイスを変更する場合に使います. 例えばElastixにある機能/設定でこれ以降の見出しの設定項目にない設定を追加したい場合などを想定した機能かと思います. Editボタンをクリックして編集していきますが, インターフェイスを作るPythonコードが立ち上がります. Pythonを理解しかつ3DSlicerのmoduleのコード記載方法を理解していないと変更できないと思います.

そのため, わからない人は触らないほうが無難な設定項目です.

  • Parameter set

この設定は, いまいちなんのためなのか理解できていません. おそらくですが, 新規作成が Create new scripted module となっているので,モジュールの調整ではと推測します. となると, 上記のReload&TestでEditしたインターフェイスを保存するような役割かと思います. あまりいじる必要はないように思います.

  • Inputs

DIRする画像を選択します. Fixed volume,Moving volumeは基準となるボリュームと変形されるボリュームとそれぞれ指定してください.

PresetがDIRの設定となります. 設定となります書きましたが, 設定を細かく変更はできません. すでに作られている設定を選ぶだけになります. プルダウンで選択項目を見ると, CT-CTだけでなくMRI,PETなどマルチモダリティでDIRする設定も多くあります. 自分がしたい処理にあったものを選んでください.

ここで, これらPresetの設定が何をしているかですが, プルダウンで表示される選択項目はSlicerがインストールされているフォルダの下

[Slicerのフォルダ]/lib/Slicer-4.11/qt-scripted-modules/Resources/RegistrationParameters

というフォルダに, ElastixParameterSetDatabase.xmlというファイルがあります. ここで規定されています. ElastixParameterSetDatabase.xmlの冒頭部分は下のようになっています.

<ParameterSet … > </ParameterSet>までが一つのPresetになっており, modality= に記載されている名前がPresetとしてプルダウンメニューに出てきていると思います. <PrameterFiles>以下にあるのがこのPresetで使用する処理の設定ファイルの名前です. Presetの「general(all)」で実行されるのは Parameters_Rigit.txt と Parameters_BSpline.txt ですので, ファイル名的に, rigitをしてB-SplineのDIRを実施する処理のようです.

<?xml version="1.0" encoding="UTF-8"?>
<ElastixParameterSetDatabase>
  <ParameterSet id="default0"
    modality="generic"
    content="all"
    description="'default' parameter files for starting elastix users"
    publications="">
    <ParameterFiles>
      <File Name="Parameters_Rigid.txt" />
      <File Name="Parameters_BSpline.txt" />
    </ParameterFiles>
  </ParameterSet>
  <ParameterSet id="default0"
    modality="generic rigid"
    content="all"
    description="'default' parameter files for starting elastix users"
    publications="">
    <ParameterFiles>
      <File Name="Parameters_Rigid.txt" />
    </ParameterFiles>
  </ParameterSet>

またParameters_Rigit.txt と Parameters_BSpline.txt はElastixParameterSetDatabase.xmlがあるフォルダにあり, そこには他に様々な名前の.txtファイル(設定ファイル)があります. それぞれの具体的な処理の設定が.txtファイルで規定されています.

  • Masking

DIRをする範囲を指定できます. maskに使用するのはLabelmapになり, DICOM-RTを読み込んで輪郭をMaskにしたい場合は, 輪郭をLabelmapに変換しておく必要があります. Labelmapがあれば, どれをMaskに使用するか選択できます.

  • Outputs

出力先を指定します. 基本的に出力しないことはないと思いますので, Create new volume かCreate new volume as…を選ぶことになると思います. いずれも3D slicer内部で保存されている状態になるので, ファイルへは別途保存する必要があります.

output volumeは変形後の画像, output transformはDVFを出力することになります.

  • Advanced

ここはelastix上級者向けに役立つ機能です. 基本的にDIRをの処理を調整したい場合(新たにPresetを作るとか, 設定txtを作るとかの場合)かと思います.

Show detailed log during registrationは詳しいlogを表示してくれますし, Keep temporary filesで途中経過のファイルを残してくれます. Registration presetsは右のボタンで, 上のPresetのところで説明した設定ファイルがあるフォルダを開いてくれます.

Force grid output transformは未確認ですが, Transformモジュールで出力されたtransformのデータを見ると, チェックが入っていないほうがTransformモジュールのinformationのところにより詳しい情報が表示されます. B-Splineの情報などが出ますので, チェックがないとDVFになっていないパラメータの状態でデータを保管しているのかと思います. 3DSlicerのデータ形式でvtkMRMLBSplineTransformNodeとvtkMRMLGridTransformNodeと別になっていますので, チェックを入れることで, BSplineTransformNodeをGridTransformNodeに変換してくれて出力するのかと思います.

Custom Elastix toolbox locationもよくわかっていません. 右の入力項目にはファイルを指定するようです.

  • Applyボタン

Applyボタンでこれまでの設定でDIRを実行します.

3DSlicerによるDIR: Plastimatch

OS: Ubuntu 20.04

3DSlicer Version: 4.11.20210226

SlicerRT(Plastimatch含む) version: 1.0.0

PlastimatchのDIRのパラメータの説明はすでにコマンドライン版として過去に書いていたのでそちらも参考にしてもらえればと思います. 3DSlicerでできるDIRはPlastimachで設定できる項目の一部のようです. 3DSlicerでのPlastimatchの結果に満足行かない場合は, コマンドラインのPlastimatchを使用するか, 他のDIRソフトを使用する必要があると思います.

今回は主にインターフェイスの説明になります.

まずはDIRをしたい2つの画像が開かれている状態を前提とします.

その後 モジュールをPlastimatchのB-Spline deformable registrationに変更します. 下の画面ではすでに選択されていますが, Module右のプルダウンで, 「Plastimatch」 ->「Registration」->「B-Spline deformable registration」を選択してください. そうすると, 下のような設定項目のインターフェイスになります.

Plastimatch Moduleのインターフェイス

見出し頭の ▶ をクリックすることで内容が表示されます.

DIR関係で一番上となる「Parameter set」は「B-Spline deformable registration」が入っていると思いますが, この選択はこれ以降の設定をひとまとめとして保存しておく機能のようです.

Parameter setをプルダウンで開くと「Create a new parameter set」とparameter setを作成できるので良い設定ができれば保存しておけばいいと思います.

左図は少し入力していますが, はじめは何も入力されていない状態だと思います.

各パラメータの説明 (最低限必要な設定は赤字)

  • Fixed volume, Moving volume

Fixed volume, Moving volumeは説明不要かと思います. それぞれ, 3DSlicerで開いている画像をプルダウンで指定できます.

  • Input B-spline transform(MRML)

B-spline で得られた変形を記録したMRMLを指定. 一度行ったDIRから再度DIRしたいときに使う(かな?). MRMLは3DSlicerのデータ保持の内部形式なので, 後にあるOutput B-spline transform (MRML)で出力を指定した場合(それを保存していた場合)に指定可能になると思います.

MRMLは3DSlicerのDataモジュールで左図のインターフェイスの代わりに表示される開いているデータ一覧で確認できます.

  • Input vector field (MRML)

上の B-spline transformと同様 ここでは3DSlicerに保存されてたDeformation vector fieldを指定. B-spline transformと両方指定する必要があるのかは未確認

  • Input transform (file)

上記MRMLで3DSlicerに記録されているもの以外で, transformを保存したファイルを指定. 3DSlicerのTransformモジュールでrigit registrationは手入力で可能なので, そこで保存したファイルも指定できると思われます.

  • Fixed fiducials

FiducialマーカーとなるPoint of Interestを指定. DICOM-RTファイルを読み込んだ際はアイソセンターや線量評価点として指定したものがリストにあげられるので, プルダウンで指定します. プルダウン最後に「Create new MarkupsFiducial」があるので新たに作ることも可能

  • Moving fiducials

Fixed Fiducialsで指定した点と対応する点を指定.

  • Fixed fiducials (Old), Moving fiducials (Old)

上記と基本的には同じと思われるが, 何が違うか不明です. おそらく互換性のため古いバージョンの際の処理を残している?

  • Output volume

変形後の画像を出力名を指定します. ファイルではなく3DSlicer内部に保管する形ですのでもとのCT画像などと同じ扱いで扱えるようになります. もちろん後で3DSlicerで保存すればファイルに保存できます.

  • Output B-spline transform (MRML)

これも計算されたB-Splineの変形パラメーターと言ったほうがいいのかもしれませんが, 3DSlicerの内部に保管されます(MRML形式) . 後でファイルに書き出すことも可能です.

  • Output vector field (MRML)

こちらはDVF版. これも3DSlicerの内部に保管されます.

  • warped fiducials (MRML)

これも変形後のfiducialマーカー(POI) . プルダウンメニューの最後にあるCreate new MakeupsFiducial で新たな名前をつけて保存したほうがいいと思います.

  • Output B-spline transform (file), Output vector field (file), warped fiducials (file)

これら3つは上記の内容をファイルに保存する指定です. ファイル名を指定します.

  • Cost function

選択しとしてMSE (Mean squared error)かMI(Mutual information)のどちらかの選択となっています. 基本的に同じモダリティ-の画像をDIRしたいのであればMSE, 異なるモダリティであればMIが適切かと思います.

  • Hardware

GPUを使うか, CPUのみの処理を指定します. GPUが使えるならそのほうが早いです. ただし, CUDAが使えるビデオカードを搭載していることや, ライブラリを入れておくなど自身で設定する場合は少しハードルが高いかと思います(Linuxの場合. Windowsはライブラリなどは少し楽かも, 未確認ですが).

よく分からなければCPUを選択.

  • Output type

出力されるデータの数値の精度を選択できます. が, 基本Autoでいいのではと思います.

プログラミング経験がある人には馴染みの変数の型ですが, やらない人には馴染みのない名前ですので. ちなみにDoubleが一番精度がいい小数点数を保存できます. 少数はあとFloatのみ. ほかは整数になるので, 選択すると少数が出た場合切り捨てられます. この設定がどこまで(どの出力にたいして)有効なのかは未確認です. 浮動小数点数があるので, DVFなどもここの指定で制御されるのかと思っています.

Stageの設定

設定項目はまだスクロールすると続きがあります. ここから, レジストレーションの設定になります. スクロールしたパラメータの設定を下の図にします.

設定項目はすべて, パラメータを指定するか選択する(しない)かしなければなりません. デフォルトで設定値などが入力されているのでそれを基準に調整してはと思います.

Plastimatchのインターフェース(続き)

Stage 0 Options

  • Enable Stage 0?

このStage 0を実施するかどうかを指定. チェックボックスにチェックを入れると実施されます.

Stage 0はいわゆるrigit registrationになります.

画像全体的な位置のズレがある場合に実施したほうがいい設定. 位置ずれがなくても臨床データの場合, 回転などや体の歪みなどもあるので, チェック入れておいて損はないと思われます.

  • Image subsampling rate (vox)

対象の画像の解像度を粗くする設定. 右図のように5,5,3と指定するとX,Y軸は5voxelごとにサンプリングして解像度を1/5にします. Z軸は3なので3スライスごと.つまり100,100,60のボリュームデータがあれば, この設定だと20,20,20になります.

このサブサンプルの設定があることで, 以降のDIRのStageが複数設定できる意味がうまれます.

はじめは粗い解像度でregistrationすることで画像ノイズに強くなり, Stageを繰り返す中で解像度を上げることで, より精密に合わせていく処理になります.

  • Max iterations

最適化の最大繰り返し回数を指定します. 画像の移動等と画像の違いの計算(類似度)を繰り返し徐々に最適な位置を計算する際の, 最大繰り返し回数. 適当な数値はデータによるので試すしかないですが, このExtensionのDocumentには200以上はあまり意味がないと記載されています. まずはデフォルトで実施して満足行かなければ変更するような感じになるかと思います.

  • Transformation

Stage 0の特有の設定. Rigitの自由度を設定します. TはTranslation, RはRigit, AはAffineの意味.

Tは平行移動のみの3DOF, Rは平行移動と回転を加えた6DOF, Aは更に拡大縮小などの線形な変形を加えた12DOFとなります. 臨床データで1軸だけ拡大や縮小するようなケースはないと思うので, TかRがデフォルト選択になると思われます.

Stage N Options (N=1,2,3)

ここ以降のStageがDIRの設定となります.

  • Enable Stage N?

Stage 1の場合にはない設定項目. このStage (2か3)を実施するかの設定. Image subsampling rateのところで書いたように, subsampling rateなどを変更しつつDIRしていく場合に, (3DSlicerでは最大3段階まで)複数Stageを使って実施します.

Stage 1は必須処理となるので, この設定項目はありません.

Image subsampling rate(vox)

Stage 0と同様, 対象の画像をsubsamplingするための間隔です.

  • Grid size(mm)

B-splineでDIRする場合, コントロールポイントが等間隔に置かれ, その位置を調整することで変形します. そのコントロールポイントの間隔. mm単位での指定. subsamplingの粗さを考慮して設定が必要. subsamplingで下げた解像度よりも細かいGrid sizeを指定してもあまり意味がないことになります.

  • Regularization

DIRでCost functionで指定した類似度を計算する際, ねじれなどが起こらないように類似度を悪くする計算項を入れる(Regularization/正規化と言われる)が, その項の重み付けの値になります. 大きくすれば変な変形が起こりにくいが最終的な結果が出るまで時間がかかる(繰り返し回数がより多く必要になる)ことになります.

  • Landmark Penalty

おそらくfiducialを指定したときにそれらの位置ずれによって類似度を悪くする項があり, その項の重み付けと思われます. PlastimatchのDocumentにも同じ名前の設定がなく, おそらくlandmark_stiffnessにあたると判断(未確認)

  • Max iterations

Stage0に同じ

  • Output volume_1

Stageごとに結果画像を出力できる設定になります. 結果の名前を指定する. 3DSlicer上に保管されるのでファイルが出力されるわけではありません. 後でファイルに保存することは可能. Noneだと出力しない設定になります.

Options

Stage 3の後にある設定.

  • Interactive registration

PopupでマニュアルでDIRをコントロールしたい際にチェックを入れると出るが, チェックを入れて実行しても実行が途中で止まって入力を待つような感じがなくそのまま終了してしまいました. そのため, このOptionsの設定項目がどう反映されるかは不明です.

Apply

設定項目の一番下部右側に「Apply」ボタンがあるので, 上記設定後ボタンを押せばDIRが実施されます.

デフォルト設定(Stage1のみ)で, 160枚程度のCT画像をDIRするのにGPUを使う設定で数分かからないくらいです.

ちなみにPCのスペックは Ryzen 7-3700, メモリ32GB, GeforceRTX-3080といったところです.