立石科学技術振興財団

助成研究成果集(第29号) 2020

［ 研究助成 (A) ］

圧縮センシングを用いた呼吸性移動のある臓器位置・形状高速推定法の開発
Development of Rapid Position Estimation Method using Compressed Sensing
for Abodminal Organ with Respiratory Motion2181014

研究代表者

神戸大学大学院
システム情報学研究科

助 教

國領大介

[研究の目的][研究の内容，成果]

肝臓をはじめとした呼吸性移動のある臓器に先に挙げた目標を達成するために，本研究で生じた腫瘍などの患部に対し，温熱治療，特には MR 画像のデータ領域である k 空間の特徴

治療対象を動きに合わせて加温位置を追尾するを利用するとともに反復回数を減らし再構成時必要のある集束超音波治療を安全かつ効果的に間を短縮するための圧縮センシング技術の改実現するためには，照射位置を患部の動きに合良・評価を実施した。そして，呼吸性移動のあわせて追尾し，適切な加温を実現する必要がある肝臓に含まれる複数の血管に注目し，事前にる。温熱治療時の患部の監視には，軟部組織の取得した画像と圧縮センシング技術により少な描出に優れ，形態情報のみならず代謝・温度変い計測データから復元した再構成画像を比較し，化・機 能 情 報 も 取 得 で き る 磁 気 共 鳴 画 像

患部などの位置移動を予測する位置・形状追尾(Magnetic Resonance Imaging, MRI) の利用が手法を提案し，その精度を検討した。

有用であり，子宮筋腫などの MR ガイド下集束超音波治療において術中の状態監視に適用さ1．位置推定実現に向けた圧縮センシング技術の改良れている。その一方，MR 画像の撮像時間は長いため，呼吸性移動を捉えるためには撮像技術1. 1取得する k 空間領域の設定の高速化，もしくは計測データ数を減らすことスパース性を持つ未知のベクトルを L1 ノルによる高速化が必要である。本研究ではスパーム最小化により推定する技術である圧縮センシ

ス性を持つ高次元信号を復元する方法である圧ングでは，原信号を等間隔に計測するのではな

縮センシング (Compressed Sensing) を適用し， く，ランダムな間隔で計測することで少ない事前情報として複数の呼吸位相において取得しデータ数で再構成画像を復元することが可能でた MR 画像と，計測データ数を削減し圧縮セある。本研究では圧縮センシングをスパース性ンシングを用いて復元した位置推定用 MR 画の高い MR 画像に適用することで撮像の高速像を用いて，臓器内の注目位置・形状の高速推化を目指した。MR 画像は計測データを格納し定・追尾手法を開発することを目標とした。ている k 空間データを逆フーリエ変換することで得られる。k 空間領域は図 1 に示したように，周波数方向と位相方向が直交した領域であり，― 1 ―立石科学技術振興財団

助成研究成果集(第29号) 2020

12min Q xt= ℱ xt−y 2+γ1 ψxt12+γ2 TV 1
where TV 1yx=   x i, , t−x i1, , t +x i, t−x i, 1, t22

1 i1

目的関数 Qxtの第 1 項は復元する信号と収集した信号のずれを表し，ℱ はアンダーサンプ
リングを行いフーリエ変換を行う関数，xtは反図 1 k 空間データの例①低周波数帯，②中間の周波数帯，③高周波数帯

復 t 回目の再構成画像の画素値ベクトル，yはデータ収集に関しては，周波数方向のデータを収集する周波数エンコードは時間がかからないが，位相方向のデータを収集する位相エンコードには時間がかかる。そこで圧縮センシングを適用する際には位相エンコード回数を減少させることができるようにする。また，低周波数帯である k 空間の中心付近は MR 画像にした際に撮影対象の大まかな形状を表すデータが格納されており，一方高周波数帯である k 空間の辺縁部は撮影対象の細かい部分を表すデータが格納されている。本研究では，低周波数帯・中間の周波数帯・高周波数帯に分け，全体の取得データ数を間引きながら，位置・形状追尾が可能なデータ収集範囲を検討した。位置・形状推定手法において対象臓器内の情報を利用するた原画像の信号値ベクトルである。第 2 項は再構成画像をウェーブレット変換を行うことで再構成画像のスパース性を持つように変換している。
第 3 項は再構成画像のエッジ情報の劣化を小さくする全変動項を表している。本研究では信号適の一致度を最優先とし，次に優先する項の決定のために 1 より十分小さいそれぞれの重み係数γ1, γ2 を設定した。また目的関数の勾配の値に着目し，初期勾配値よりも十分に小さくなれば最適解に近い解が得られたとして計算を打ち切ることで反復回数削減による再構成時間の短縮を行った。最適化を行う手順を以下に示す。
STEP 1 : 信号収集，初期値x0ならびに初期解に対する勾配ベクトルg0, 修正ベクトルd 0を設定

めに，今回，大まかな情報を取得できる低周波数帯のデータが重要と考え，その領域のデータ収集量を多く求めた。
1. 2x0=0, g0=d 0=−∇Q0x0
STEP 2 : バックトラック法によりステップ幅αを決定し，t+1回目の画素値を更新

再構成に用いる最適化関数の改良

xt1=xt+αtd t

圧縮センシングを用いた画像再構成では欠損したデータから復元を行うため反復計算が必要となり，再構成時間が大きくなってしまう。そこで，信号復元を最適化問題とし，復元信号と収集した信号のずれだけでなく，再構成画像のSTEP 3 : 反復 t 回目の最大画素値x  に 1 未満の定数rdを掛け，決定した閾値λtを用いたソフト閾値処理を実施

スパース性やエッジ情報の劣化の評価を同時にλt=x  ×rd行うことにより，反復回数を削減し計算時間の0x t, = x t, −λtx t, x t, 短縮を行った。今回提案した最適化関数を次式に示す。

― 2 ―x t, <λtx t, ≥λt立石科学技術振興財団

助成研究成果集(第29号) 2020

STEP 4 : 評価関数を偏微分し，勾配ベクトルgt1, 修正ベクトルd t1を計算
gt1=−∇Q(xt1)
2d t1=gt1+

 gt1 22

 gt2dt
図2

再構成画像の例

2STEP 5 : 規 定 回 数 l に 到 達 す る か， gt12<2 g02×ε な ら 反 復 計 算 終 了，そ う で な け れ ば
l≔l+1として STEP 2 に戻る。
1. 3

異なる収集領域・データ量を用いた時の再構成画像を評価するため，1.5 テスラ MR 装置で撮像した画像サイズ 512×512 画素の肝臓 MR画像を用いた評価実験を行った。提案手法に関連したパラメータは以下の通りである。収束判定係数 ε : 0.001, 規定回数 l : 100, 重み係数 λ 1,λ 2 : ともに 0.01。予備検討より，全データ領域の収集率は 50% とした。画像評価にはミクロ
領 域 の 誤 差 を 評 価 す る RMSE (Root MeanSquare Error) を用いた。収集領域は k 空間全体をランダム (Random) に計測した場合，並びに表 1 に示した収集率で計測した 3 パターンの計 4 パターンを用いた。なお ①, ②, ③は図1 に示した k 空間における周波数領域の範囲を示した。試行回数は 20 回とした。
表 2 に各収集領域パターンにおける目的関数値Q，肝臓実質部分に関心領域を設定して計算した RMSE 値 (0 に近いほど良い) 並びに計算時間を示した。目的関数値Q，RMSE 値とも

Pattern 1
Pattern 2
Pattern 350% にすることで計測にかかる時間は全領域を収集した場合に比べ半分に抑えられるとともに，低周波領域を多く収集することで，呼吸性移動によるノイズも低減でき，位置・形状追尾に使用可能な画像が取得できると考えられた。
計 算 時 間 に 関 し て は，デ ー タ 収 集 率 は 同 じ50% であったため有意な差はなかった。今回使用した画像が 512×512 画素と画 素数が 多かったため，計算時間が 5 秒以上かかっているが，よく用いられる 256×256 画素の画像であれば，計算時間も短縮でき，高速な位置・形状追尾に利用できると考えられた。
2．臓器内位置・形状高速推定・追尾手法の提案
圧縮センシングを用いることにより，計測時間を短縮することで，呼吸性移動のある臓器に対し，臓器内において注目している位置並びに形状を高速に追尾できる可能性が高まった。次に臓器内の複数の血管位置に注目し，事前画像における血管位置と取得画像における血管位置を比較し，現在の呼吸位相並びに取得画像以降各パターンの収集率①100%90%80%なった。図 2 に Random ならびに Pattern 1 で得られた再構成画像を示す。データ収集率を評価実験

表1に低周波数領域を多く収集するにつれて良く②60%80%80%③0%0%
10%の血管位置ならびに治療位置など注目位置の移動位置を捉える手法を提案し，評価を実施した。
2. 1

提案する位置・形状推定・追尾手法の流れ

提案手法の流れを以下に示す。
表2Random
Pattern 1
Pattern 2
Pattern 3

収集領域変更による各値の変化
目的関数値[×10 6]6.87±7.12
1.93±2.32
3.44±3.14
4.46±2.99

RMSE
0.17±0.01
0.11±0.00
0.13±0.01
0.14±0.01

計算時間
[sec]
5.58±0.35
5.85±0.64
5.22±0.10
5.58±0.18

STEP 0 : 事前画像として，複数の呼吸状態における画像を収集する。この画像はデータの収集率は 100% とする。
STEP 1 : 調整した収集データ率ならびに収集― 3 ―立石科学技術振興財団

助成研究成果集(第29号) 2020

領域において信号を計測し，圧縮センシングを用いて再構成画像を作成する。
STEP 2 : 血管領域を抽出しやすくするためにエッジを強調する画像処理を実施する。
STEP 3 : 事前画像における複数の注目血管位置と取得した画像における注目血管位置を比較し，事前画像における呼吸位相を判定する。
図4

STEP 4 : 同じ呼吸位相と判断した事前画像を移動位置追尾・予測方法の概要

用いて，取得画像以降の注目血管ならびに治療置の予測時には，計測中にも呼吸による位置の位置など注目血管以外の移動位置を推定・予測移動が起きていることを想定し，事前画像におする。

いて現在の呼吸位相と同じであると判断した画

STEP 5 : 取得画像は次の移動位置を予測する像に対し，撮像にかかる時間を考慮した数枚後ための事前画像に取り入れ，STEP 1 に戻る。の画像を基準とした (図 4 参照)。決定した基次節以降において STEP 3 ならびに STEP 4

準画像における注目血管の位置は移動後に取得について詳細を述べる。

した画像の血管位置と同じになると考え，治療2. 2

部位など注目血管以外の移動位置を予測する際呼吸位相判断

呼吸により臓器位置は上下に移動しているたには，移動前後の注目血管の位置関係は変化しめ適切な照射位置の決定には術中画像がどの呼ないという関係を用いて導出する。

吸状態で撮像された画像なのかを判断すること

2. 4

評価実験が重要である。そこで対象とした臓器内に含ま提案手法の位置推定精度を評価するためにれる血管位置に注目し現在の呼吸位相判断を行

1.5 テ ス ラ MR 装 置 で 撮 像 し た 画 像 サ イ ズう。図 3 に呼吸位相判断の概要図を示す。エッ

256×256 画素の肝臓 MR 画像を用いた評価実ジを強調した複数の注目した血管の重心位置を験を行った。事前画像は 3 呼吸周期分の MR取得し，各事前画像における血管の重心位置と画像を使用し，取得画像は 1 章において検討し取得画像における血管の重心位置を比較する。

たデータ収集率 50% ならびに低周波数帯を中そして，差が最小となった事前画像の呼吸位相心とした収集領域の計測データに対し，圧縮セが取得画像における呼吸位相であると判断する。

ンシングを用いて再構成した画像を使用した。

2. 3

位置追尾精度は，呼吸約 1 周期分の 11 枚の再取得画像以降の移動位置追尾・予測方法

2. 2 節で決定した呼吸位相情報と複数の注目構成画像を用いて評価した。なお使用した肝臓血管の位置関係を用いて，現在取得した画像以MR 画像は正常の肝臓であったため，位置推定降の注目した血管位置並びに治療部位など注目に用いない血管を治療位置として推定精度の評血管以外の移動位置を推定・予測する。移動位価に用いた。
図 5 に提案手法を用いて得た移動位置と目視により求めた移動位置との差を示した。使用した 11 画像の平均誤差は 2.26±1.18 mm，最大誤差は 4.59 mm であった。超音波治療などの温熱療法では患部の大きさにもよりが 5 mm 程度のマージンを取って治療することを考えると，今回の結果は実際の治療に用いることが出来る

図3

呼吸位相判断の概要精度であると考えられる。またデータ収集率― 4 ―立石科学技術振興財団

助成研究成果集(第29号) 2020

置・形状追尾手法の実現には，精度・計算時間の改善が必要であることがわかっている。今後，精度・計算時間の更なる改善手法を検討するとともに，MR シミュレーション技術を組合わせることによる新しい手法を開発し，更なる精度向上を目指す予定である。

図5[成果の発表，論文等]

位置推定精度結果[ 1 ] 仲川侑介，國領大介，貝原俊也，藤井信忠，熊50% の画像を用いていることから，これまで

本悦子：高速 MR 温度計測のための圧縮センシンよりも高速に追尾することが可能であると考えグを用いた画像再構成手法の検討，第 61 回自動制ている。一方，画像 3 や 8 等，他の画像に比べて誤差が大きくなった場合も存在した。これは御連合講演会，名古屋，4A2, 2018.11.
[2]

Y. Nakagawa, D. Kokuryo, T. Kaihara, N. Fujii, E.

Kumamoto : Image reconstruction method with

最大の呼気もしくは最大の吸気付近の呼吸位相

compressed sensing for high-speed MR tempera-

であり，呼吸状態による違いが大きく現れやす

ture measurement of abdominal organs, Proc. 41st

い部位であったと考えられる。今後，そのような状態においても精度を向上させるための検討

Annual International Conference of the IEEE Engineering in Medicine and Biology Society, Berlin,
Germany. 2019.7[3]が必要であると考えている。

仲川侑介，國領大介，貝原俊也，藤井信忠，熊本悦子：腹腔臓器を対象とした高速 MR 温度計測のための圧縮センシングを用いた画像再構成手法の[今後の研究の方向，課題]検討，第 63 回システム制御情報学会研究発表講演
会，大阪，TS12-1-1, pp. 1115 - 1119, 2019.5

本研究において，臓器位置・形状推定手法において，少ないデータ収集においても圧縮センシングを用いることにより，一定の精度での位置推定・追尾が可能であることが分かった。一方，課題として，3 次元計測データを用いた位― 5 ―[4]仲川侑介，國領大介，貝原俊也，藤井信忠，熊本悦子：腹腔臓器を対象とした MR ガイド下集束
超音波治療のための高速照射位置同定手法の提案，
生体医工学シンポジウム 2019, 徳島，2A-33, 2019.9