社団法人 電子情報通信学会
THE INSTITUTE OF ELECTRONICS,
INFORMATION AND COMMUNICATION ENGINEERS

信学技報
TECHNICAL REPORT OF IEICE.

マルチクラス辞書学習を用いた低線量 CT 画像再構成
鴨下

飛竜†

大地†

北原

平林

晃†

† 立命館大学大学院情報理工学研究科 〒 525–8577 滋賀県草津市野路東 1–1–1
E-mail: † is0254iv@ed.ritsumei.ac.jp, d-kita@fc.ritsumei.ac.jp, akirahrb@media.ritsumei.ac.jp
あらまし

X 線 CT における低線量投影データからの高精度な画像再構成手法を提案する．パッチベースの単一辞書を用いる従来の再構成手法では，辞書によってスパースに表現できないパッチが多く存在してしまい，十分な再構成
結果が得られない．そこで本研究では，パッチをクラスタリングし，クラスごとに異なる辞書を用いて再構成を行う．
実際の CT 画像を用いたシミュレーションにより提案法の有効性を示す．
キーワード 低線量 CT，画像再構成，辞書学習，クラスタリング

Low-Dose CT Image Reconstruction with Multiclass Dictionary Learning
Hiryu KAMOSHITA† , Daichi KITAHARA† , and Akira HIRABAYASHI†† Graduate School of Information Science and Engineering, Ritsumeikan University
1–1–1 Nojihigashi, Kusatsu, Shiga, 525–8577 Japan
E-mail: † is0254iv@ed.ritsumei.ac.jp, d-kita@fc.ritsumei.ac.jp, akirahrb@media.ritsumei.ac.jp
Abstract We propose a high-precision CT image reconstruction method from low-dose X-ray projection data. In
conventional reconstruction methods using a single patch-based dictionary, there exist many patches that cannot
be expressed sparsely by the dictionary, and thus the reconstruction quality would not be suﬃcient. The proposed
method adopts patch clustering and uses a diﬀerent dictionary for each patch cluster. Simulations using real CT
data demonstrate the eﬀectiveness of the proposed method.
Key words Low-Dose CT, Image Reconstruction, Dictionary Learning, Clustering

1. は じ め に

低線量の投影データから高品質な画像を再構成するために，
圧縮センシング技術 [4] が用いられている．圧縮センシングは，

X 線コンピュータ断層撮影法（Computed Tomography: CT）

対象信号のスパース表現の存在を仮定した上で，不十分な観測や磁気共鳴画像法（Magnetic Resonance Imaging: MRI）は，から信号を高精度に復元する技術である．画像処理分野を中心

体内を非侵襲的に撮像できるため，様々な疾患の診断や発見にに，対象信号のスパース表現を獲得する手段として辞書学習を役立っている．CT は骨折や肺野病変などの診断に優れており，利用する [1, 3, 5–8]．辞書は，応用ごとに適切な学習用画像を

MRI は水分の多い脳などの診断に優れている．また，CT には用いて基底（またはフレーム）として作成され，基底ベクトルのMRI と比べて撮像時間が短いという利点がある．例えば，全身線形結合により対象信号をスパースに表現できるようになる．

を撮像する際に，MRI は各種パラメータに依存して 30 分から

低線量 CT の分野においては，統計的反復再構成（Statistical

1 時間程度を要するのに対して，CT は数分しか必要としない．

Iterative Reconstruction: SIR）法 [9] に辞書学習を適用した

一方で，CT は撮像に X 線を用いるため，癌性や遺伝性の疾患

画像再構成手法が Xu らにより提案されている [1]．この手法誘発が懸念されており，なるべく少ない X 線量で撮像すること

は，ストリークアーチファクトの除去に優れているものの，再が望まれている [1]．X 線量を低減する方法には大きく分けて，

構成画像の細部構造が過度な平滑化によって失われてしまう

1）X 線の照射強度を下げる，2）X 線の投影回数を減らす，の 2という問題を抱えている．ノイズやアーチファクトを抑制しな通りがある．しかし，どちらの場合においても，低線量の投影

がら，画像の細部構造を鮮明に再構成するために，画像をよりデータに対して，通常の再構成手法であるフィルタ補正逆投影

スパースに表現することが課題となっている．

（Filtered Back Projection: FBP）法を適用すると，画像内にこの課題を解決するために本論文では，Xu らの手法をマル

大きなノイズや不自然なアーチファクトが発生してしまう [2]．

チクラスに拡張した再構成手法を提案する．Xu らの手法では，結果として，正常組織と病理組織の識別が難しくなる [3]．

学習用画像から抽出された全てのパッチに対して単一の辞書を—1—

R ∈ RI×J ，ξ ∈ RI である．通常の X 線 CT では，I ≧ J かつ

𝑦𝑣

ノイズレベル ∥ξ∥2 も小さいため，角度 θg ごとの投影データ

𝑢X線源

lθg ,n （n = 1, . . . , N ）を事前に設計した周波数フィルタにより補正した後に逆投影する（FBP 法または逆ラドン変換と呼ば

𝜇(𝑥, 𝑦)れる）ことで，画像 µ を高精度に再構成できる．最終的な CT𝜃	  値は，水の X 線減衰係数 µwater を基準として，

𝑥(𝛾,,. 𝑣

CT = 1000𝑧,,-µ − µwaterµwater)(3)のように計算される [10]．式（3）から，水の CT 値は 0 であり，
𝑧,,/
𝑧,,.

空気中では X 線が減衰しないため CT 値は −1, 000 である．

検出器

2. 2 低線量 CT における辞書学習を用いた反復再構成

図 1 X 線 CT における投影の模式図

撮像時の X 線量の低減方法として，1）X 線強度 b を下げる，

2）投影回数 G を減らす，の 2 通りが考えられる．しかしながら，用いている．一方，提案法では各パッチをクラスタリングし，

前者では，測定値 zθg ,n が小さくなり投影データ lθg ,n の信頼性

クラスごとに別々の辞書を作成する．そして，これらの辞書をが弱まる（ノイズレベル ∥ξ∥2 が増大する）という問題が発生し，用いて対象画像のパッチをクラスごとに処理し，画像再構成

後者では，I = GN が減少して I < J となる問題が発生する．

を行う．このとき，辞書によるパッチの近似精度はクラスごとしたがって，どちらの場合においても単に FBP 法を適用するに異なることが想定されるため，辞書を用いた正則化項の重みだけでは，画像内に大きなノイズや不自然なアーチファクト（正則化パラメータ）をクラスごとに設定することで，再構成が現れてしまい良好な結果は得られない．そこで Sauer らは，精度を向上させる．実際の CT 画像を使用して投影回数を減ら

最大事後確率（Maximum a Posteriori: MAP）推定により µした場合のシミュレーションを行い，提案法の有効性を示す．

を再構成する手法を提案した [9]．この手法は，統計的なコスト

2. X 線 CT における画像再構成∥Rµ − l ∥22,w + λΨ(µ) :=2. 1 X 線 CT の基本原理

I∑wi (riT µ − li )2 + λΨ(µ) (4)i=1

図 1 の xy 座標系で X 線 CT における画像再構成を説明する．

を最小にする µ を反復的に再構成するため，SIR と呼ばれる．

被写体内部の X 線減衰係数分布を µ(x, y) ≧ 0 とする．投影は，重み付き ℓ2 ノルムの各重みは統計的に wi := zi = b/ exp(li ) と

xy 座標系を角度 θ だけ回転させた uv 座標系

定まり，Ψ : RJ → R は µ の先験情報に基づいた正則化項，i は{投影のインデックス，λ > 0 は正則化パラメータである．

u = x cos θ + y sin θ

Xu らは，式（4）にパッチベースの辞書学習を組み合わせた
√
√
P× P

v = −x sin θ + y cos θ

手法を提案した [1]．この手法では，画像 µ をサイズ⋆を用いて記述される．投影データ lθ,nは，X 線源から第 n 番目のパッチと呼ばれる小領域 Hs µ ∈ RP （s = 1, . . . , S ）に分割の検出器までの直線経路 γθ,n (v) に沿った減衰係数 µ(x, y) の線

し，「それぞれのパッチが辞書 D ∈ RP ×K と呼ばれる適切な⋆
積分 lθ,n:=µ(γθ,n (v)) dv として定義され，この積分はラドン

基底（またはフレーム）によりスパースに表現できる」という変換と呼ばれる．ラドン変換は数値積分により離散化されて，先験情報に基づいて µ を再構成する．ここで，Hs ∈ {0, 1}P ×J

∫

⋆
T
⋆
lθ,n
= rθ,nµ + ξθ,nはパッチ抽出行列，s はパッチのインデックス，S はパッチ数，(1)と表される．ここで，µ ∈ RJ は画素数 J の離散画像の画素値，⋆
rθ,n ∈ RJ は数値積分の重み係数，ξθ,n∈ R は数値積分の誤差を表している．第 n 番目の X 線検出器における測定値 zθ,n ∈ Z はポアソン分布に従う整数であり，その期待値は X 線の照射強度⋆
b > 0 を用いて E[zθ,n ] = b exp(−lθ,n) と与えられる．この関係⋆
から投影データを lθ,n= ln(b/E[zθ,n ]) と表せるが，E[zθ,n ] は

K は基底数であり，D の列ベクトル dk ∈ RP （k = 1, . . . , K ）は ∥dk ∥2 = 1 を満たしている．各パッチのスパース係数を

cs ∈ RK （s = 1, . . . , S ）とおけば，再構成結果は最適化問題
minimize ∥Rµ − l∥22,w + λS∑µ,(D),C(∥Hs µ − Dcs ∥22 + νs ∥cs ∥0 )

s=1

(5)
の解として与えられる．ここで，C := (c1 , . . . , cS ) ∈ R

K×S

で未知である．したがって，通常は，測定値 zθ,n から計算されるあり，νs > 0 はパッチごとのラグランジュ乗数である．Xu らは

lθ,n := ln(b/zθ,n ) を投影データとして用いる．G 個の角度 θg ∈式（5）の辞書 D に関して，学習用画像を用いて事前に作成する[0, π) に対して N 個の検出器により投影データを取得する場合，手法と，再構成中の画像を用いて適応的に更新する手法を提案

観測モデルは式（1）をまとめてしている．前者は GDSIR（Global Dictionary based SIR），

l = Rµ + ξ(2)のように表現される．ここで，I := GN とおけば，l ∈ RI ，

後者は ADSIR（Adaptive Dictionary based SIR）と呼ばれる．
式（5）の最適化問題を直接解くのは困難であるため，FBP
法を適用して得られる画像を初期値 µ0 として，スパース係数—2—

3. 2 マルチクラス GDSIRC（及び辞書 D）と画像 µ を交互に更新する．C（及び D）をm−1

更新する際には，µ が推定値 µminimize∑(D),Cとする．クラス q のパッチを表現するための辞書 Dq ∈ RP ×Kq(∥Hs µm−1 − Dcs ∥22 + νs ∥cs ∥0 )(6)は，学習用画像を用いて事前に作成されるものとする．ここで，

s=1を考える．式（6）の問題に対して，GDSIR では OMP [12] を用いて係数 C のみを更新する．一方 ADSIR では，文献 [13] の手法によって辞書 D を更新した後に，OMP で係数 C を更新するため，GDSIR に比べて多くの計算時間を必要とする．

µ を更新する際には，C （及び D）が固定された最適化問題
minimize ∥Rµ − l∥22,w + λµS∑∥Hs µ − Dcs ∥22(7)を考える．最適解 µ∗ は，対角行列 W := diag(w) を用いて

(

µ =

T

R WR + λ

S
∑

)−1(

HsT HsT

R Wl+λs=1

S∑

Kq は各辞書の基底数である．再構成する画像のパッチ Hs µ はクラスタリングされており，クラス q のパッチのインデックス∪Q

集合を Sq とすれば，

q=1

Sq = {1, . . . , S} かつ Sq ∩ Sq′ = ∅′（q ̸= q のとき）が成り立つ．本研究では，マルチクラス辞書

Dq を用いて式（5）を一般化した最適化問題
minimize ∥Rµ − l∥22,w +

Q
∑

µ,C

λq

q=1

s=1

∗

クラス数を Q とし，クラスのインデックスを q = 1, . . . , Qに固定された最適化問題

S)(∥Hs µ−Dq cs ∥22 +νs ∥cs ∥0 )

s∈Sq

(9)
を解くことにより µ を再構成する手法を提案する（注 1）．ここで，λq > 0 はクラスごとの正則化パラメータである．

HsT Dcs

s=1∑

パッチのクラスタリング手法には様々なものが考えられるが，計算時間と実装の容易さを考慮して，本研究では文献 [8] とと表されるが，画像サイズが大きい場合には逆行列計算が困難

同じく K-means 法を使用する．まず，学習用画像のパッチをとなってしまう．そこで，現在の推定値 µm−1 において式（7）K-means 法により Q 個のクラスに分類し，クラスごとに辞書と関数値が等しくなるような上界関数を設計し，この上界関数

Dq を作成する．この際，各クラスのクラスタリングセンターをを最小にする µ を新たな推定値 µm として用いる [11]．この値

記憶しておく．次に，式（2）の投影データ l に FBP 法を適用

は画素値 µm
j （j = 1, . . . , J ）ごとに式（8）で与えられる．

して初期画像 µ0 を求め，記憶したクラスタリングセンターを

3. マルチクラス辞書学習を用いた反復再構成

利用して µ0 のパッチを K-means 法により分類する．その後，式（9）の問題の解を反復的に計算することで，提案する再構成

3. 1 複数の辞書を用いた再構成手法が達成される．以下で提案アルゴリズムを詳細に説明する．

近年，高速撮像 MRI や低線量 CT の分野において，圧縮率3. 3 提案する反復再構成アルゴリズム

向上やノイズ除去・アーチファクト抑制のために，マルチクラス

提案するアルゴリズムは，事前準備を 1 回した後に，係数更新

辞書学習を用いる再構成手法が提案されている [3,7, 8]．これらステップと画像更新ステップを交互に M 回繰り返すことで，の手法では，まずパッチを幾何学的方向または輝度値によって投影データ l から画像 µ を再構成する．

複数のクラスに分類する．そして，クラスごとに別々の辞書を事前準備 まず，学習用画像からパッチを抽出し，K-means

用意することで，パッチのスパース表現の近似精度を向上させ，法により Q 個のクラスに分類する．再構成する画像のパッチ

画像を良好に再構成できるようにしている．しかし，文献 [7, 8]

も後にクラスタリングするため，各クラスタリングセンターを

は高速撮像 MRI を対象としており，そのままでは低線量 CT

記憶しておく．次に，クラスごとの辞書 Dq を最適化問題には適用できない．文献 [3] は低線量 CT を対象としているが，minimize

最適化問題のデータ整合性として，FBP 法によって得られる

Dq ,C

∑

2
tr
(∥Hs µtr − Dq ctr
s ∥2 + νs ∥cs ∥0 )(10)s∈Sqtr画像との誤差を用いるため，この画像が大きなアーチファクトを含んでいる場合には，再構成結果が著しく劣化してしまう．

の解として作成する．ここで，µtr は学習用画像，ctr
s は第 s 番

特に，投影回数を減らす方式の低線量 CT においては，FBP 法

目の学習用画像パッチに対するスパース係数，Sqtr はクラス q

によって大きなアーチファクトが発生する（第 5 節の図 3（b）の学習用画像パッチのインデックス集合である．Dq が辞書とと図 4（b）を参照）．このため，文献 [3] の手法は適さない．

して適切であるためには，各パッチを少ない誤差で表現できる

一方，GDSIR と ADSIR には，文献 [3] の手法に比べて計算能力が必要と考え，本研究では式（10）を等価な問題である

量は多いものの，X 線量の低減方法を限定しないという利点が

Dq ,Cある．そこで本研究では，GDSIR をマルチクラスへと一般化した，高精度な再構成手法を提案する．提案法では，FBP 法で得られる画像を利用して各パッチのクラスタリングを 1 回のみ行うため，計算量は GDSIR からほとんど増加しない．また，
クラスごとに別々の正則化パラメータを用意しており，これらを適切に設定することで，高精度な再構成が実現可能である．

∑I
m−1
µm
−
j = µji=1

i=1

rij wi∑Jj ′ =1

∥ctr
s ∥0

2
s.t. ∥Hs µtr − Dq ctr
s ∥2 ≦ εs (11)s∈Sqtrに変形する [4]．ここで，εs > 0 は各スパース表現の近似誤差の許容範囲を示している．式（11）の問題に対する高速解法が文献 [13] に提案されており，本研究でもこの解法を使用する．
（注 1）：式（9）において Q = 1 とすれば，式（5）が直ちに導かれる．

rij wi (riT µm−1 − li ) + λ∑I∑minimize∑S
s=1

rij ′ + λ

∑P

∑S

p=1

s=1

hspj ([Hs µm−1 ]p − [Dcs ]p )

∑P

p=1

hspj

∑J

j ′ =1

hspj ′(8)—3—

式（9）の問題を解くに当たり，初期画像 µ0 と集合 Sq を決定する．まず，投影データ l に FBP 法を適用して得られる画像を初期値 µ0 とする．次に，学習用画像のパッチに対するクラスタ
リングセンターを再度利用して，初期画像の各パッチ Hs µ0 を

𝜆" = 600

K-means 法により Q 個のクラスに分類する．これにより Sq が𝜆& = 𝜆' = 60

定まり，以降，画像が更新されてもクラスタリングは行わない．
以上の事前準備が完了したら，以下の係数更新ステップと𝜆( = 𝜆) = 𝜆* = 𝜆+ = 0.06

画像更新ステップを m = 1, . . . , M に渡って繰り返すことで，最終的な再構成結果 µM を得る．
係数更新ステップ（a）学習用画像（No. 1601）

スパース係数 C の第 m 回目の更新時に（b）クラスタリングセンター

図 2 事前辞書学習に用いた画像と各クラスタリングセンター

は，画像 µ が µm−1 に固定された最適化問題

minimizeS∑C

4. 計算機シミュレーション(∥Hs µm−1 − Dq(s) cs ∥22 + νs ∥cs ∥0 )(12)
4. 1 実 験 設 定

s=1を考える．ここで，q(s) は第 s 番目のパッチが属するクラスを

実際の CT 画像を用いた数値実験によって提案法の有効性を表している．式（12）の問題は，係数 cs ごとの最適化問題

示す．学習用画像と対象画像には同一人物の頭部画像を用い，

minimize ∥Hs µm−1 − Dq(s) cs ∥22 + νs ∥cs ∥0

CT 値が記録されていたため，1 cm あたりの水の減衰係数を

サイズはどちらも 256 × 256（J = 65, 536）である．画像には
csµwater = 0.2059 cm−1 ，1 画素の横幅と縦幅を 0.09 cm として，に帰着され，GDSIR と同様に OMP を用いて cs を更新する．
画像更新ステップ

式（3）から真の X 線減衰係数 µ を作成した．通常の X 線 CT

画像 µ の第 m 回目の更新時には，係数において，1 度刻みで 0 から 179 度までの G = 180 個の角度

C が固定された最適化問題
minimize ∥Rµ − l∥22,w +に対して投影を行うと想定し，MATLAB のラドン変換行列 R
Q
∑

µ

q=1

λq

∑

に基づいて，X 線検出器の数は N = 367 とした．このときは，∥Hs µ − Dq cs ∥22 (13)I = GN = 66, 060 > J が成り立つため，FBP 法により良好な

s∈Sq

再構成結果が得られる．投影回数を減らす方式の低線量 CT をを考える．式（13）の問題の解は閉形式で与えられるものの，再現するために，3 度刻みで 0 から 177 度までの G = 60 個の第 2.2 節で説明したように画像サイズが大きくなると数値計算

角度に対してのみ投影を行う（圧縮率 1/3 の投影と呼ぶ）．測定が困難になる．文献 [11] に基づいて，現在の推定値 µm−1 に値は本来ポアソン分布に従うが，本実験では X 線の照射強度をおいて関数値が等しくなるような上界関数を設計することで，

b = 106 とし，b exp(−riT µ) を四捨五入した値を zi とした．

m

新たな推定値 µが画素値ごとに式（14）のように与えられる．

提案法の比較対象を GDSIR 及び ADSIR とし，それぞれの

3. 4 正則化パラメータ λq の重要性手法で M = 1, 000 回の反復計算により再構成を行った．辞書

マルチクラス辞書学習では，クラスごとにパッチ数が違い，

学習と係数更新に用いるパラメータは，全ての手法で K = 256，

更に学習用画像と対象画像でパッチ数の比率も大きく変化する

εs = 0.2，νs = 0.2 に固定した．高速辞書学習 [13] と OMP [12]可能性がある．そのため，辞書によるスパース表現の各パッチの実装には公開されているプログラム [14] を使用し，高速辞書に対する近似精度はクラスごとに異なっている．したがって，学習の設定項目は更新回数 2, 000，バッチサイズ 40，t0 = 0.001

良好な再構成結果を得るために，基底数 Kq や辞書学習のパラとした．各画像からサイズ 8 × 8（P = 64）のパッチを，1 画素

メータ εs ，及び係数更新のパラメータ νs をクラスごとに設定ずつずらしながら S = 62, 001 個抽出する．提案法と GDSIRすることが望ましい．しかしながら，設定するパラメータ数はで用いる辞書は図 2（a）の学習用画像を用いて作成され，提案クラス数に比例して増大するため，これら全てをクラスごとに法では学習用画像のパッチを K-means 法によって Q = 7 個の設定することは現実的ではない．この問題を解決するために，クラスに分類した．また，各クラスタリングセンターは図 2（b）

提案法では式（9）においてクラスごとの正則化パラメータ λqのようになった．正則化パラメータは，GDSIR と ADSIR ではを導入した．これにより，εs , νs , Kq が全てのクラスに対して

λ = 60 とし，提案法では均一な値 λq = 60（q = 1, . . . , 7）に共通の値であっても，λq がクラスごとに適切に設定されていれした場合と，クラスごとに最適な値 λ1 = 600, λ2 = λ3 = 60,ば高精度な画像再構成を実現できる．具体的には，辞書による

λ4 = λ5 = λ6 = λ7 = 0.06 にした場合で実験を行った．計算のスパース表現の近似精度が高いクラスでは，λq を大きくする

環境は MacBook Pro，OS X Yosemite version 10.10.5，プロ

ことで辞書の影響を強調し，近似精度が低いクラスでは，λq を

セッサ 2.9GHz Inrel Core i5，メモリ 8GB 1867MHz DDR3，小さくすることで辞書の影響を抑えるようにすればよい．

プラットフォーム MATLAB R2015a 64-bit（maci64）である．

∑I
m−1
µm
−
j = µji=1

rij wi (riT µm−1 − li ) +

∑I

i=1

rij wi∑Jj ′ =1∑Q
q=1

rij ′ +

λq

∑

∑Q

q=1

∑P

s∈Sq

λq

∑

p=1

s∈Sq

hspj ([Hs µm−1 ]p − [Dq cs ]p )

∑P

p=1

hspj

∑J

j ′ =1

hspj ′

(14)
—4—

表1

圧縮率 1/3 の投影からの再構成時間 [s]再構成手法 事前辞書

事前準備

提案法では，FBP 法で得られる初期画像のパッチを学習用画像

更新 1 回と同じクラスタリングセンターを用いて 1 回のみ分類するため，
マルチクラスでありながら計算量をほとんど増加させずに画像

GDSIR

4.76

1.51 平均 1.31

ADSIR

——

1.55 平均 5.84

提案法

35.08

1.75 平均 1.58

を再構成できる．また，クラスごとの正則化パラメータを新たに導入することで，辞書学習や係数更新のパラメータが全てのクラスに対して共通であっても，高精度な画像再構成が可能と

4. 2 実 験 結 果

なった．実データを用いたシミュレーションにより，提案法がまず，提案法と GDSIR 及び ADSIR の再構成結果の精度をシングルクラスの辞書学習を用いる GDSIR や ADSIR に比べ比較する．再構成する対象画像を学習用画像と同一人物の頭部て，遥かに高精度に画像を再構成できることを示した．なお，

画像とし，スライス番号は No. 1610 と No. 1615 を選んだ．

提案法は画像再構成に 1, 000 回の更新で約 30 分の時間を要し画像の評価指標として PSNR（Peak Signal-to-Noise Ratio）とているので，高速化が今後の課題である．

SSIM（Structural SIMilarity）の 2 つを用いた．図 3 と図 4 はNo. 1610 と No. 1615 の画像を圧縮率 1/3 の投影から再構成
した結果を示しており，（a）が原画像，（b）が FBP 法により得られる初期画像，（c），（d），（e），（f）がそれぞれ GDSIR，ADSIR，λq を均一に設定した提案法，λq を最適に設定した提案法の再構成結果である．図 3（b）と図 4（b）からも分かるよ
うに，FBP 法により得られる初期画像には多くのアーチファクトが含まれている．GDSIR や ADSIR では，正則化パラメータλ を小さくするとアーチファクトが抑制しきれずに残ってしま
い，大きくするとアーチファクトのみならず画像の細部構造も失われた，過度に平滑化された画像が再構成されてしまう [3]．
一方，提案法では，正則化パラメータ λq を全てのクラスに対して均一に設定した場合でも，GDSIR と ADSIR に比べて高い

PSNR と SSIM を与えており，マルチクラス辞書学習により再構成精度が向上していることが確認できる．更に，正則化パラ
メータ λq をクラスごとに最適に設定した場合では，PSNR とSSIM の値が大幅に上昇しており，実際に再構成された画像を見てみると，アーチファクトを抑制しながらも，複雑な構造は平滑化されていない良好な画像が得られていることが分かる．
次に，提案法と GDSIR 及び ADSIR の計算時間を比較する．
表 1 は各手法における事前辞書学習，事前準備，反復再構成の
更新 1 回の計算時間を示している．提案法では，学習用画像のパッチをクラスタリングして Q = 7 個の辞書を作成するため，
事前辞書学習において GDSIR に比べて 7 倍強の時間を要している．一方，ADSIR では事前に辞書を作成しないため，表中で空欄になっている．事前準備は，FBP 法による初期画像 µ0 の作成，初期画像のパッチ抽出，パッチ Hs µ0 のクラスタリング（提案法のみ），式（8）または式（14）の分母計算のことであり，提案法のみクラスタリングを行うため計算時間が約 0.2 秒長くなった．反復再構成の更新 1 回にかかる時間は，ADSIR のみ
辞書を再構成中の画像を用いて適応的に学習するため，非常に長くなる．提案法では式（14）の計算や OMP による係数更新の際に，プログラム内でクラスを区別しなければならないため，

GDSIR に比べて時間が約 0.3 秒増加した．以上から，提案法は高精度な画像を比較的高速に再構成できることが確認された．

5. お わ り に本研究では，低線量 CT における画像再構成に対して，マルチクラス辞書学習を利用する反復再構成アルゴリズムを提案した．

文

献[1] Q. Xu, H. Yu, X. Mou, L. Zhang, J. Hsieh, and G. Wang,“Low-dose X-ray CT reconstruction via dictionary learning,” IEEE Transactions on Medical Imaging, vol. 31, no. 9,
pp. 1682–1697, 2012.
[2] 田中 敏幸, “X 線 CT の原理・現状とさらなる画像の高品質化,”
計測と制御, vol. 56, no. 11, pp. 874–879, 2017.
[3] Y. Chen, L. Shi, Q. Feng, J. Yang, H. Shu, L. Luo, J. L. Coatrieux, and W. Chen, “Artifact suppressed dictionary learning for low-dose CT image processing,” IEEE Transactions
on Medical Imaging, vol. 33, no. 12, pp. 2271–2292, 2014.
[4] M. Elad, Sparse and Redundant Representations: From Theory to Applications in Signal and Image Processing. New
York, NY: Springer, 2010.
[5] C. Zhang, T. Zhang, M. Li, C. Peng, Z. Liu, and J. Zheng,“Low-dose CT reconstruction via ℓ1 dictionary learning regularization using iteratively reweighted least-squares,” Biomedical Engineering Online, vol. 15, no. 66, pp. 1–21, 2016.
[6] Y. Chen, X. Yin, L. Shi, H. Shu, L. Luo, J. L. Coatrieux, andC. Toumoulin, “Improving abdomen tumor low-dose CT
images using a fast dictionary learning based processing,”
Physics in Medicine and Biology, vol. 58, no. 16, pp. 5803–
5820, 2013.
[7] Z. Zhan, J. F. Cai, D. Guo, Y. Liu, Z. Chen, and X. Qu,“Fast multiclass dictionaries learning with geometrical directions in MRI reconstruction,” IEEE Transactions on
Biomedical Engineering, vol. 63, no. 9, pp. 1850–1861, 2016.
[8] 稲室 憲人, 柴田 基, 湯 晨旻, 井尻 敬, 平林 晃, “学習用画像と観測信号を併用した辞書学習による MR 画像高精度再構成,” 電子情報通信学会技術研究報告, vol. 116, no. 395, pp. 123–128, 2017.
[9] K. Sauer and C. Bouman, “A local update strategy for iterative reconstruction from projections,” IEEE Transactions
on Signal Processing, vol. 41, no. 2, pp. 534–548, 1993.
[10] 小西 達郎, 浅井 義行, 上村 正信, 山田 浩司, 松本 政雄, “被写体物質の違いに起因する X 線 CT 画像コントラストの管電圧依存性の変化,” 医用画像情報学会雑誌, vol. 25, no. 4, pp. 90–95, 2008.
[11] I. A. Elbakri and J. A. Fessler, “Statistical image reconstruction for polyenergetic X-ray computed tomography,” IEEE
Transactions on Medical Imaging, vol. 21, no. 2, pp. 89–99,
2002.
[12] Y. C. Pati, R. Rezaiifar, and P. S. Krishnaprasad, “Orthogonal matching pursuit: Recursive function approximation with applications to wavelet decomposition,” in Conference Record of The Twenty-Seventh Asilomar Conferenceon Signals, Systems and Computers, 1993, pp. 40–44.
[13] J. Mairal, F. Bach, J. Ponce, and G. Sapiro, “Online learning for matrix factorization and sparse coding,” Journal of
Machine Learning Research, vol. 11, pp. 19–60, 2010.
[14] J. Mairal, J. P. Chieze, Y. Chen, G. Durif, R. Jenatton,
F. Bach, J. Ponce, G. Obozinski, B. Yu, G. Sapiro, and
Z. Harchaoui, Spams: a SPArse Modeling Software ver. 2.6.
[Online] http://spams-devel.gforge.inria.fr/doc/html—5—

（a）原画像（No. 1610）（b）FBP 法により得られる初期画像（c）GDSIR [31.26 / 0.817]（d）ADSIR [30.48 / 0.769]（e）提案法（均一な λq ）[32.10 / 0.858]（f）提案法（最適な λq ）[37.04 / 0.980]

図 3 No. 1610 の CT 画像の再構成結果 [PSNR / SSIM]（a）原画像（No. 1615）（b）FBP 法により得られる初期画像（c）GDSIR [31.51 / 0.831]（d）ADSIR [30.71 / 0.786]（e）提案法（均一な λq ）[32.28 / 0.865]（f）提案法（最適な λq ）[37.00 / 0.980]

図 4 No. 1615 の CT 画像の再構成結果 [PSNR / SSIM]—6—