JP2018110867A - 医用画像生成装置及び医用画像生成方法 - Google Patents

Abstract

【課題】ＯＳ加速されたアルゴリズムを用いる逐次再構成法において、収束率を改善することができる医用画像生成装置等を提供すること。
【解決手段】本実施形態に係る医用画像生成装置は、取得部、再構成部を具備する。取得部は、Ｘ線発生装置から照射され被検体を透過した透過Ｘ線を用いて検出された、前記被検体に関する検出データを取得する。再構成部は、前記検出データを複数のサブセットに分けて、前記複数のサブセットのうち２以上のサブセットに基づいたオブジェクト関数における第一の勾配情報を取得し、前記サブセットそれぞれに基づいたオブジェクト関数における第二の勾配情報を前記第一の勾配情報に基づき補正し、当該補正後の第二の勾配情報に基づいたオブジェクト関数を用いて前記サブセット毎に逐次近似再構成処理を行い、再構成画像を取得する。
【選択図】 図２

Description

本実施形態は、医用画像生成装置及び医用画像生成方法に関する。

コンピュータ断層撮影（ＣＴ）システムおよび方法は、特に医用撮像および医用診断の為に、幅広く使用されている。ＣＴシステムは、一般的に被検体の身体に関する一枚または複数の断面的なスライスの画像を作成する。Ｘ線源などの放射源は、一側面から身体にＸ線を照射する。身体と反対側にある少なくとも一つの検出器は、身体を通過した放射線を受け取る。身体を通過してきた放射線の減衰は、検出器から受け取った電気信号を処理することで測定される。

ＣＴサイノグラムは、検出器アレイに沿った位置の関数として、またＸ線源からＸ線検出器までの投影角の関数として、身体を通過した減衰を示す。サイノグラムにおいて、空間的次元は、Ｘ線検出器のアレイに沿った位置を参照する。時間／角度次元は、ＣＴスキャンの間の時間の関数として変化するＸ線の投影角を参照する。画像化された被検体の部位から生じる減衰は、投影角に一致して縦軸に沿って正弦波を描くだろう。回転軸から更に遠い画像化された被検体の部位は、より大きな振幅での正弦波に対応し、当該正弦波の位相は、回転軸周辺の角度位置に一致する。逆ラドン変換―また任意のその他の画像再構成法―の実行は、サイノグラムにおいて表された投影データから画像を再構成することである。

断層撮影法における統計的逐次再構成（ＩＲ）アルゴリズムは、フィルタ補正逆投影（ＦＢＰ）法の様により従来的な再構成法と比べて、線量レベルを抑えて改善された画質を提供することが可能である。しかし、特定の実施形態において、統計的アプローチは相当な計算時間を必要とするので遅い。標準的且つ統計的な再構成アプローチの計算的に集約された（ｃｏｍｐｕｔａｔｉｏｎａｌｌｙ ｉｎｔｅｎｓｉｖｅ）演算の遅さを改善するためには、より少ない回数逐次で、より素早く収束する統計的な再構成のための逐次アルゴリズムを使用する改善された方法が、注目を集めている。

加速器法は、オーダードサブセット（ＯＳ：順序部分集合）、ネステロフの加速法を含む、ＩＲ法と様々な組み合わせが出来る。ＯＳ法は、画像再構成アルゴリズムの逐次毎の測定のサブセットのみを使用することで、計算的コストを効果的に減らす。ネステロフの加速法も、ＩＲ法の計算効率とパフォーマンスを向上するために使用することが出来る。

米国特許出願公開第２０１３／０３４３６７３Ａ１号明細書

しかしながら、特定の実施形態について以下に説明するように、ＯＳ加速されたアルゴリズムは収束しない可能性がある。例えば、ネステロフの加速法と組み合わせられた場合に、ＯＳ加速されたアルゴリズムは発散という課題に悩まされることがある。

本実施形態は、上記課題に鑑み、ＯＳ加速されたアルゴリズムを用いる逐次再構成（ＩＲ）法において、分散を縮小することで、収束率を改善することができる医用画像生成装置及び医用画像生成方法を提供することである。

本実施形態に係る医用画像生成装置は、取得部、再構成部を具備する。取得部は、Ｘ線発生装置から照射され被検体を透過した透過Ｘ線を用いて検出された、前記被検体に関する検出データを取得する。再構成部は、前記検出データを複数のサブセットに分けて、前記複数のサブセットのうち２以上のサブセットに基づいたオブジェクト関数における第一の勾配情報を取得し、前記サブセットそれぞれに基づいたオブジェクト関数における第二の勾配情報を前記第一の勾配情報に基づき補正し、当該補正後の第二の勾配情報に基づいたオブジェクト関数を用いて前記サブセット毎に逐次近似再構成処理を行い、再構成画像を取得する。

一実施形態に係る、小数のオーダードサブセット（ＯＳ）と、多数のＯＳと、更にＯＳおよびネステロフの加速法の組み合わせとを使用する、逐次再構成（ＩＲ）方についての収束率のプロットを図示する。 一実施形態に係る、ＩＲの様々な実施形態の収束率についてのシミュレーション結果のプロットを図示する。 一実施形態に係る、分散が縮小されたＯＳ ＩＲ法のフローチャートを図示する。 一実施形態に係る、分散が縮小されたＯＳ ＩＲ法の疑似コードでのフローチャートを図示する。 一実施形態に係る、ネステロフの加速法を伴う分散が縮小されたＯＳ ＩＲ法の疑似コードでのフローチャートを図示する。 腹部断面の参照再構成画像を図示する。 図５に示された参照再構成画像に対応する、オーダードサブセットと分離可能二次サロゲート（ＯＳ−ＳＱＳ）との組み合わせを使用しサブセット２４個とを使用して、六回の逐次後の再構成画像を図示する。 図５の参照画像と図６Ａの再構成画像との間の差分である差分画像を図示する。 図５に示された参照再構成画像に対応する、オーダードサブセットと分離可能二次サロゲートとの組み合わせとネステロフの加速法（ＯＳ−ＳＱＳ−ｍｏｍ）を使用しサブセット２４個とを使用して、六回の逐次後の再構成画像を図示する。 図５の参照画像と図７Ａの再構成画像との間の差分である差分画像を図示する。 図５に示された参照再構成画像に対応する、オーダードサブセットおよび分離可能二次サロゲートに加えて分散縮小の組み合わせとネステロフの加速法（ＯＳ−ＳＱＳ−ｍｏｍ−ＶＲ）とを使用しサブセット２４個とを使用して、六回の逐次後の再構成画像を図示する。 図５の参照画像と図８Ａの再構成画像との間の差分である差分画像を図示する。 頭部断面の参照再構成画像を図示する。 図９Ａに示された参照再構成画像に対応する、オーダードサブセットと分離可能二次サロゲートとの組み合わせとネステロフの加速法（ＯＳ−ＳＱＳ−ｍｏｍ）を使用しサブセット１２個とを使用して、十回の逐次後の再構成画像を図示する。 図９Ａに示された参照再構成画像に対応する、オーダードサブセットおよび分離可能二次サロゲートに加えて分散縮小の組み合わせとネステロフの加速法（ＯＳ−ＳＱＳ−ｍｏｍ−ＶＲ）とを使用しサブセット２４個とを使用して、四回の逐次後の再構成画像を図示する。 一実施形態に係る、１フルセットの勾配と異なる数のサブセットとをアップデートするための異なる期間におけるＯＳ−ＳＱＳ−ｍｏｍ−ＶＲ ＩＲ法の様々な実行を使用する際の、エラー収束率についてのシミュレーション結果のプロットを図示する。 実施形態に係るＣＴスキャナの実行の概要図を示している。

逐次再構成（ＩＲ）法は、従来的なフィルタ補正逆投影（ＦＢＰ）再構成法を上回る優れた画質のために、Ｘ線ＣＴ撮像においてかなりの注目を得つつある。しかし、ＩＲは、具体的には複雑なシステムの順モデル（ｃｏｍｐｌｉｃａｔｅｄ ｓｙｓｔｅｍ ｆｏｒｗａｒｄ ｍｏｄｅｌｓ）および逐次での処理手順（ｉｔｅｒａｔｅｄ ｐｒｏｃｅｓｓｉｎｇ ｐｒｏｃｅｄｕｒｅｓ）を使用するため、より多くの計算的リソースとＦＢＰのための時間とを必要とする。素早い収束率でのアルゴリズムは、この様にして、画質に悪影響を与えることなく、逐次時間を減らすのに極めて有益である。

様々な逐次アルゴリズムは、ＣＴ画像を逐次的に再構成するために使用することが出来る。これらの様々な方法は、勾配に基づく方法、逐次座標降下法、そして最適化遷移（ｏｐｔｉｍｉｚａｔｉｏｎ ｔｒａｎｓｆｅｒ）に基づく方法などを含むことが可能である。追加での方法は、バーグマンの分割スキーム（Ｂｅｒｇｍａｎ’ｓ ｓｐｌｉｔｔｉｎ ｓｃｈｅｍｅ）と同様に、増大されたラグランジュ乗数フレームワークに基づくものでも良い。

これらのアルゴリズムの大半は、その断層撮影再構成問題の不良設定性質（ｉｌｌ−ｐｏｓｅｄ ｎａｔｕｒｅ）の為に、やはり収束が遅く、そのため加速法は収束率を上げるのに有益である。勾配に基づく方法において、素早い収束アルゴリズムは、例えばいくつか所定の前提条件で共役勾配法を使用して、ニュートンタイプのアップデートを実現するために、使用することが出来る。加えて、オーダードサブセット（ＯＳ）法が別の加速法である。ＯＳ加速の一つの利点としては、ＯＳ加速が追加のアルゴリズムの複雑性を加えず、故に大半の存在するアルゴリズムと容易に組み合わせられるという点である。また、特定の場合において、ネステロフのモメンタム法は、存在する逐次アルゴリズムの中の多くに適用された場合も、二次収束率を提供出来る。二つ以上の加速法の組み合わせは、もし本実施形態で説明される方法の様に適切に行われれば、どちらか片方の手法だけの場合と比較して、一層の収束の加速という結果になることがある。例えば、ネステロフのモメンタム法をＯＳに基づく分離可能二次サロゲート（ＳＱＳ）法と組み合わせる、本実施形態で説明される方法を使用して、再構成において相当な時間短縮が達成出来る。

時間短縮を制限する可能性がある一つの潜在的な課題として、ＯＳ型アルゴリズムは、リミットサイクル問題によく悩まされることが挙げられる。ＯＳ型アルゴリズムがモメンタム法と組み合わせられた場合に、当該アルゴリズムの収束率はより複雑になる。実際に、後ほど説明される通り、特定の組み合わせにおいて、ＯＳ型アルゴリズムは、発散的な解という結果になる可能性がある。この問題への解決策としては、確率的勾配最適化法の理論に基づく分析から提案されるが、この解決策の意味するところは、逐次と最終的な解との間の画像エラーは二次的な割合（ｑｕａｄｒａｔｉｃ ｒａｔｅ）で低減することが出来るものの、実際には確率的なノイズが原因の当該エラーは逐次数の線形関数として増える、ということである。この確率的なノイズは、勾配近似エラーと呼ぶことが出来る、つまり、オリジナルのオブジェクト関数（目的関数）の勾配と単一のサブセットのみを使用して定義された勾配との間のミスマッチである。この様にして、より大きな分散におけるより多数のサブセット結果は蓄積されており、逐次数と共に発散率の上昇を引き起こしている。

この問題を取り巻く一つの解決方法は、逐次数が増えるにつれてステップサイズが減るという弛緩したモメンタム法を使用することである。当該弛緩したモメンタム法は、減少するステップサイズを処理する追加のパラメータを取り込むが、この追加のパラメータは困難且つしばしばデータに依存することがある。更に、減少するステップサイズは、収束スピードが長引くことに繋がり、減少するステップサイズの制御には、再構成ワークフローを複雑化させる追加のパラメータ調整が必要となる。

収束問題を解消するその他の方法は、ステップサイズのやり直しスキーム、従来的な勾配に基づくアルゴリズム、そしてその他の全体的な収束法を含む。ステップサイズのやり直しスキームは、ネステロフの加速法を使用する際に、発散問題を制御するために使用することが出来るが、リミットサイクル問題は回避出来ない。ＯＳ型アルゴリズムは、前提条件が調整された共役勾配法など従来的な勾配に基づくアルゴリズムを使用することで、完全に回避が可能である。しかし、効果的な前提条件を計画する事は、特に広いコーン角画像再構成について、とりわけ困難なことが判明する可能性がある。インクリメンタルな最適化遷移アルゴリズムなど、その他の全体的な収束法を適用することが出来る。

上述の方法よりも更に良いアプローチが、本実施形態に説明される分散縮小法である。この優れたアプローチは、リミットサイクルや分散問題を解消するために、ＯＳ型アルゴリズムとネステロフ加速法との組み合わせで分散縮小法を使用する。更に、本実施形態に説明される方法は、最適な収束率を有利に維持しつつ、逐次再構成手順の間に確率的ノイズを効果的に低減する。その上、本実施形態に説明される方法は、余分な追加パラメータを必要としない。

要約すると、高速逐次アルゴリズムは、画質を低下させることなく、全体的な再構成時間を減らす手助けが出来る。本実施形態に説明される方法において、オーダードサブセット（ＯＳ）とネステロフのモメンタム法とによって加速された高速逐次再構成を改善するための、分散縮小方法が提案される。本実施形態で説明される方法の改善無しでは、ＯＳ加速されたアルゴリズムが、収束を長引かせながら、リミットサイクルの制約を受ける可能性がある。ネステロフ法と更に組み合わせられた場合に、アルゴリズムはノイズ蓄積が原因で頻繁に発散する。本実施形態に説明される方法は、このノイズを低減し、本実施形態に説明される結果によって実証されているように、結果として速い率での収束を可能にしている。

参照画像が数枚の図にわたって同一または対応する部分を指し示す、図面を次に参照する。図１Ａは、逐次数の関数として再構成エラーのプロット図を示す。逐次数は、横軸（Ｘ軸）に沿って描かれている。再構成エラーは、縦軸（Ｙ軸）に沿って描かれており、所定の逐次数での再構成画像と参照画像との間の差を表す。プロットは、ネステロフの加速法無しのＯＳを使用する収束を図示する。更に、これらのプロットは、多数のＯＳを使用すると、ＩＲは初めこそ素早く収束するものの、リミットサイクルが多数のＯＳに対し大きなエラー制約を強いる、ということを実証している。ＯＳ法単独の使用よりも更に悪いことに、所定数の逐次の後、ＯＳとネステロフ加速法との組み合わせることにより、再構成画像が参照画像からの発散を現にもたらしていることである。

そのため、大抵のＩＲアルゴリズムはゆっくりと収束することを観察することが出来る。ＯＳを使用するアルゴリズムは、収束スピードを上げるために適用出来るが、ＯＳで加速されたアルゴリズムは、リミットサイクル問題のために収束しないかもしれず、図１Ａに示されているような発散現象という結果になる。基本的に、多数のＯＳは大きなリミットサイクルを作り出す、つまり多数のＯＳは比較的少数且つ早期の逐次で早く収束するが、この収束は分散によって生み出された比較的大きなエラーにより最終的に制限され、逐次の間に投影データの全体ではなくサブセットを使用することに起因する。ネステロフ加速法が、ＯＳを使用するＩＲ法に追加される場合に、ＯＳ-ＡＣＣアルゴリズムと呼ばれ、勾配不整合によるエラーは、逐次処理の間に蓄積される可能性がある。更に、仕様するサブセット数が大きくなるほど、ＯＳ−ＳＱＳアルゴリズムがより素早く発散する様に、不整合エラーが益々大きくなる。

図１Ｂは、実際のデータを使用して、図１Ａの説明において上で特徴付けられた制限された収束および発散現象を例示して、生成されたシミュレーション結果である。図１Ｂにおける結果の詳細説明は、「ＯＳ−ＳＱＳ」、「ＯＳ−ＳＱＳ−ｍｏｍ」、「ＯＳ−ＳＱＳ−ｍｏｍ―ＶＲ」アルゴリズムについて説明された後に、以下で提供される。

特定の実行において、ＩＲ法は、正則化され重み付けられた最小二乗再構成に対するオブジェクト関数を最小化することで、例えば、以下の式（１）を（ただし、第２式はΦ（ｘ）の定義式）、Ｌ（ｘ）を定義する式（２）で解くことにより、画像を生成するために公式化することが出来る。

ここでｙは投影データ、ｘは再構成されている画像、Ｐは順投影作用素を表す。マトリクスＷは、対角に沿う値が統計的重みである対角マトリクスを取ることが出来る。正則化関数Ｕ（・）は、正則化項の強さを調節する正則化パラメータβによって増加される。最後に、符号Ｔはマトリクスまたはベクトルの転置を意味する。

ＯＳ法が適用された場合に、投影データは、Ｎ（Ｎ＞０）重複しないサブセットへと区切られる。各サブセットに対し、部分的なオブジェクト関数は、式（３）の様な定義が出来る。

なお、Ｌｎ（ｘ）は、式（４）によって定義される。

ここでＡ_ｎはＡのｎ番目のサブセットである。多くの場合、Φ_ｎ（ｘ）は、以下の式（５）を満たす。

合計オブジェクト関数Φ（ｘ）は、構成要素サブセットオブジェクト関数Φ_ｎ（ｘ）の平均であり、つまり、式（６）のように表すことができる。

ＯＳに基づく逐次は、ｘ^ｋ+１を獲得するために、まずΦ_ｎ（ｘ）を最小化することで、実施される（ｘ^ｋ+１は、必ずしも厳密な最小化するもの（ｍｉｎｉｍｉｚｅｒ）でなくて良いかもしれない）。その場合に、ｘ^ｋ＋１は、サブセットｎ＋１について定義された次の関数Φ_ｎ（ｘ）を解くための初期画像として、使用が出来る。この処理は、繰り返され、最後のサブセットの最小化は、各サブセットにわたる逐次の新たな一ラウンドを開始するのに第一のサブセットへと供給することが出来る。

ＯＳスキームは、要素Ｎの時間短縮をもたらす可能性がある。しかし、ＯＳを使用して加速されたアルゴリズムは、リミットサイクル効果により課されたエラーリミットに近づくにつれて、収束しなくなる可能性がある。統計的な勾配法の理論からの一解釈としては、次の式（７）で示されるエラー（または分散）Ｅ_ｎが大きくなる。

従って収束を阻止するものであり、ここでサブセット数ｎはランダム変数として扱われており、Ｅ_ｎ（・）はｎについての期待値作用素である。サブセット数が増加した場合は、変数Ｅ_ｎ（・）も増加し、より緩やかな収束に変化するようになっていく（大きなリミットサイクルと同様）。

この問題を解消するために、分散が縮小された方法が使用される。分散が縮小された方法は、ξ_ｎによって示され、∇Φ_ｎ（ｘ）に対するサロゲートを使用し、アルゴリズムが逐次するにつれて、次の式（８）で示されるエラー（又は分散）は、式（７）のエラーよりもずっと小さくなることが出来る具合である。

例えば、特定の実施形態において、ｎ番目のサブセットに対するサロゲート関数は、以下の式（９）のように表すことが出来る。

ここで、μ＝∇Φ（ｘ^〜）およびｘ^〜は、中間画像とすることが出来る（ここで、左記表記では、^〜はｘの上に載っているのもとする）。

このサロゲートは、式（１０）、式（１１）という事実により観察された様に、上の公式化と整合する。

ここで二番目の等式は、以下の式（１２）の事実による。

この事実が暗示するのは、サロゲート関数が不偏であり、従ってオリジナルオブジェクト関数を変更するようなものではない、ということである。

ξ_ｎ、∇Φ_ｎ ^〜（ｘ）の選択は（左記表記では、^〜はΦ_ｎの上に載っているものとする）、分散を大きく縮小させることが出来る。直観的に、これは、μ→０のようにｘ^〜が推定ｘ^＾（左記表記では、^＾はｘの上に載っているのもとする）へと収束する場合を考えることで、理解されよう。これは、更に次の式（１３）、式（１４）という結果になる。

ｘ^ｋもｘ^＾へと収束する場合、次の式（１５）の様になり、勾配近似によるエラーが縮小することを意味する。

図２は、縮小された分散を伴うＯＳを使用してＩＲを実行する方法２００のフロー概要図を図示する。なお、図２に示す縮小された分散を伴うＯＳを使用するＩＲ法は、医用ワークステーション等に代表される医用画像生成装置、或いは専用プログラムをインストールしたコンピュータによって実現される。分散が縮小された方法２００は、リミットサイクルによって生じたノイズを低減する。この低減は、中間画像に基づいた１フルセットの勾配（「第一の勾配情報」とも言う。全てのビューに基づいた勾配情報が典型であるが、少なくとも２以上のサブセットに基づいた勾配情報であってもよい。）を計算することで、達成される。当該１フルセットの勾配は、その後サブセット勾配（サブセットに基づいた勾配情報であり、「第二の勾配情報」とも言う。）と１フルセットの勾配との間の差を取ることで、補正項を構成するために使用される。画像アップデートは、各サブセットに対する補正された勾配を使用して、実行される。都合の良いことに、分散が縮小された方法２００は、一般的に任意のＯＳ型再構成アルゴリズムへの適用が可能である。更に、方法２００は、オリジナルＯＳ型アルゴリズムの収束率を維持することが可能で、追加のパラメータを必要としない。

なお、１フルセットの勾配（第一の勾配情報）が全てのビューに基づいた勾配情報でなく、少なくとも２以上のサブセットに基づいた勾配情報である場合には、当該第一の勾配情報の取得に用いた前２以上のサブセットのみに対して、当該２以上のサブセットに基づいたオブジェクト関数における前記第二の勾配情報を第一の勾配情報に基づき補正し、当該補正後の第二の勾配情報に基づいたオブジェクト関数を用いて前記２以上のサブセットについて逐次近似再構成処理を行い、再構成画像を取得することになる。

方法２００のステップ２１０において、投影データｙは、収集されＯＳへと配置される。加えて、再構成画像とループ変数とが初期化される。例えば、ＦＢＰ、フェルドカンプ再構成法、または再構成画像を初期化するための任意のその他の方法を使用して、初期再構成画像の取得が出来る。

方法２００のステップ２２０において、外側のループ（ボクセルに関するループ）処理変数がインクリメントされる。

方法２００のステップ２３０において、フル勾配μは計算され（例えば、次の式（１６）参照）内側のループ（サブセットに関するループ）処理変数が当該開始値へと初期化／再初期化される。

方法２００のステップ２４０において、補正項ｃ^ｋが計算される（例えば、次の式（１７）参照。

更に、画像ｘ^ｋは、式（１８）に従って、アップデートされる。

特定の実施形態において、画像アップデート法は、ｆ（・）が、次の式（１９）として表されるように、分離可能二次サロゲートアルゴリズムを使用する。

ここでＤ＝ｄｉａｇ｛Ａ^ＴＷＡ１｝および[・]_＋は、任意の負の値をゼロへと設定する作用素である。加えて、特定の実施形態において、ネステロフの加速法は、サブセットに基づく画像アップデートが、式（２０−１）、（２０−２）、（２０−３）、（２０−４）の計算によって実行することが出来る様に、分散が縮小された方法２００と組み合わせることが出来る。

ここで個別の変数は、ｔ_０＝１およびｚ^０＝ｖ^０＝ｘ^０としてステップ２１０においてもともと初期化されている。更に、ｘ^０はステップ２１０からの初期再構成画像である。最後に、ステップ２４０において、ループ変数ｎおよびｋは、インクリメントされる。

方法２００のステップ２４５において、内側のループ処理がサブセットのそれぞれを通してインクリメントされたかどうかを判定する問い合わせが実行される。ステップ２４０は、最後のサブセットまで繰り返される。いずれかの点において、再構成画像ｘ^〜は最後のオーダードサブセット画像ｘ^ｋ、つまりｘ^〜←ｘ^ｋの値を取得するために、アップデートされる。

方法２００のステップ２５５において、外側のループ処理に対する停止基準が満たされているかどうかを判定するための問い合わせが実行される。例えば、外側のループ処理停止基準とは、再構成画像が十分に収束したかどうか、および／または逐次の最大数に達したかどうかの決定を含むことが出来る。停止基準に達していない場合、方法２００はステップ２５５からステップ２２０へと戻る（ｐｒｏｃｅｅｄ）。若しくは、方法２００は完了する。

図３および４は、疑似コードを用いて表された方法２００の様々な実行を図示する。図３は、分離可能二次サロゲートアルゴリズムが方法２００において適用された場合に、アップデート（２．２．２）ｆ（・）は、次の式（２１）の様になる。

ここでＤ＝ｄｉａｇ｛Ａ^ＴＷＡ１｝であり、且つ[・]_＋は、任意の負の値をゼロへと設定する。図３に示された方法２００の実行は、「ＯＳ−ＳＱＳ−ＶＲ」と呼ばれ、分散が縮小されたＯＳ−ＳＱＳの略である。

これに関連して、図４は、分散が縮小された方法２００にネステロフの加速法が適用された場合を図示しており、サブセットに基づく画像アップデート（２．２．２）は、次の式（２２−１）、（２２−２）、（２２−３）、（２２−４）を用いて表すことができる。

ここでｔ_０＝１、ｚ^０＝ｖ^０＝ｘ^０である。図４に図示された方法２００の実行は、ＯＳ−ＳＱＳ−ｍｏｍ−ＶＲアルゴリズムと呼ばれ、ネステロフの加速法を含む分散が縮小されたＯＳ−ＳＱＳの略である。

方法２００は、組み合わされたＯＳ−ＳＱＳアルゴリズムの適用が出来る。特定の実施形態において、ＯＳ−ＳＱＳアルゴリズムは、サブセットサロゲート関数としての式（２３）を使用する。

直行マトリクスＤは、正の値として定義され、拘束条件としての式（２４）を満たすよう選択される。

分散が縮小された方法を適用するために、二次サロゲート関数はΦ_ｎ ^〜（ｘ）に対して、式（２５）に従って構成される。

当該サロゲート関数の最小化は、次の式（２６）とすることが出来る。

同様に、ネステロフの加速法は、上での説明の様に、組み込むことが出来る。x~に対する選択肢は多数ある。例えば、サブセット逐次の一つのループ後の画像出力は、x~のアップデートとして、使用が出来る。

画像アップデート方程式としての式（２７）は、二つのサブセットに基づいた勾配∇Φ_ｎ（ｘ^ｋ）および∇Φ_ｎ（ｘ^〜）との間の差を使用して、計算される。

しかし、勾配を拡大することで、次の式（２８）を観察することが出来る。

これは一つのサブセット順投影と一回の逆投影とを使用して実行出来ることを示している。更に、１フルセットの勾配は、サブセット数とは無関係である、ｘ^〜での正則化の勾配を含む。従って、∇Ｕ（ｘ^〜）は、各逐次についてμ^Ｌ＝∇Ｌ（ｘ^〜）のみが計算された場合に必要とされず、この場合、式（２９）の様になる。

この様にして、分散縮小があるＯＳ−ＳＱＳアルゴリズムは、分散縮小がないＯＳ−ＳＱＳアルゴリズムの複雑性と等価になる。

分散が縮小された方法２００に関する一つの課題は、１フルセットの勾配μ（μ^Lと同様）の評価が計算的に冗長ということである。μの計算は、一回のフル順投影と一回のフル逆投影とを使用して実行することが出来るので、故に当該計算にかなりの計算的なリソースを使用することがある。特定の実施形態において、一逐次の間に計算される１フルセットの勾配は、次に続く逐次の間に計算される予定のフル勾配に近い近似（ｃｌｏｓｅ ａｐｐｒｏｘｉｍａｔｉｏｎ）のままであり、また更にその後の逐次に対しても近い近似のままであることが想定出来る。しかし、各サブセット逐次と組み合わせられる再構成画像ｘ^〜へと変化するにつれて、より早い段階の逐次の間に計算されたフル勾配は、接近が近い状態で現在の逐次の１フルセットを近似している、またある時点で１フルセットの勾配は、分散があまりにかけ離れる（ｌａｒｇｅ）のを留めるために計算する必要がある。この様にして、１フルセットの勾配が、外側のループ処理の各逐次において、計算されずまたアップデートされなくても（つまり、図３および４における指数ｉによって追跡されたループ）、１フルセットの勾配の補正効果は、１フルセットの勾配が外側のループ処理の各逐次よりも低頻度でアップデートされる場合であっても、より低い程度に対してではあるが、持続することが想定出来る。この様にして、各逐次についてμをアップデートするのではなく、特定の実施形態において、１フルセット勾配アップデートは、Ｔ逐次毎に周期的に実行されるのみである。ここで、Ｔは１よりも大きい所定の整数値として事前に定義することが出来る。

図１Ｂに戻って、ＯＳ−ＳＱＳアルゴリズム（つまり、分離可能二次サロゲートを伴うオーダードサブセット）およびＯＳ−ＳＱＳ−ｍｏｍアルゴリズム（つまり、ネステロフの加速法を伴うＯＳ−ＳＱＳアルゴリズム）に優る、ＯＳ−ＳＱＳ−ｍｏｍ−ＶＲアルゴリズムの利点が見て取れる。図１Ｂは、東芝Ａｑｕｉｌｉｏｎ ＯＮＥスキャナを用いて収集された、実際の臨床円形３２０スライス腹部スキャン使用での、提案された再構成アルゴリズムの効率を実証している。関心領域（ＲＯＩ）における収束率を評価するために、ＨＵに基づく二乗平均平方根差（ＲＭＳＤ）は、所定の逐次の画像と収束した画像ｘ^＾（つまり、参照画像）との間で計算され、ＲＭＳＤは、次の式（２９）によって与えられる。

ここで｜ＲＯＩ｜はＲＯＩの内側のボクセル数を意味する。オブジェクト関数で使用される正則化は、次の式（３１）によって与えられる。

ここでＮ_ｊはボクセルｊの近傍を表す。φ（ｔ）はφ（ｔ）＝（ｔ^２＋δ^２）^１／２として定義される。δは所定の定数である。重みｗ_ｉｊは、ボクセルｊとボクセルｉとの間の逆距離を使用して、計算することが出来る。収束した画像ｘ^＾は、一つのサブセットのみでＳＱＳアルゴリズムの１０，０００回の逐次を実行することで取得され、初期画像はＦＢＰ再構成法を使用して再構成されたものである。

図１Ｂは、異なる数のサブセット（ＯＳ＝１２、２４、３０）が、様々なＩＲ法の収束率を評価するために求められる。多数のサブセットを使用することで初期収束率は増えるものの、時間短縮に制限が掛かることが見て取れる。ネステロフの加速法で、ＯＳ−ＳＱＳ−ｍｏｍは初め、非常に早く収束するが、これは例えばＯＳ−ＳＱＳ（ＯＳ＝１２）とＯＳ−ＳＱＳ−ｍｏｍ（ＯＳ＝１２）とを比較することで見て取れる。リミットサイクルが原因なのだが、しかし、ネステロフの加速法でのＯＳ＝１２に対するＲＭＳＤは、約十回の逐次の後に減少が止まる。ＯＳ＝２４が使用された場合に、ネステロフの加速法でのＯＳは、約三回の逐次の後に発散する。より多くのサブセットはより大きなリミットサイクルを作り出し、また発散がより一層目立つようになっている。

本実施形態で説明された分散縮小スキームに付け加えることにより、方法２００は、ＯＳ−ＳＱＳ−ｍｏｍアルゴリズムよりもずっとパフォーマンスが上がる。ＯＳ＝１２の場合、ＯＳ−ＳＱＳ−ｍｏｍ−ＶＲアルゴリズムは、少数の逐次を使用してＯＳ−ＳＱＳ−ｍｏｍアルゴリズムと同様に実行する。しかし、逐次数が大きくなるにつれて、ＯＳ−ＳＱＳ−ｍｏｍ−ＶＲアルゴリズムは、ＯＳ−ＳＱＳ−ｍｏｍアルゴリズムとは異なり、より大きなエラーリミットに収束するのではなく、ＲＭＳＤにおいて連続的な縮小が発生する。これは、方法２００が分散を縮小することでリミットサイクルを回避し、参照画像によって表された真画像により近い解へと繋がっている、ということを実証している。同様にサブセット数が増える場合に、方法２００はリミットサイクルによって逆に影響を受けるように見えない。図１Ｂに示される通り、多数のサブセットに対してであっても、ＲＭＳＤは発散することなく縮小する。これによって、方法２００は、ＯＳ−ＳＱＳ−ｍｏｍ（ＯＳ＝１２）に対して達成されるエラーリミットよりも、ずっと低いＲＭＳＤへと最終的には収束する。

図５は、ＯＳ＝１でのＳＱＳアルゴリズムの１０，０００回の逐次を使用して生成された参照画像を示している。図６Ａ、７Ａ、８Ａは、図１Ｂに示された個別の方法に対する六回目の逐次での再構成画像を示している、つまり、図６Ａ、７Ａ、８Ａは、方法ＯＳ−ＳＱＳ、方法ＯＳ−ＳＱＳ−ｍｏｍ（ＯＳ＝２４）、方法ＯＳ−ＳＱＳ−ｍｏｍ−ＶＲ（ＯＳ＝２４）に個別に対応する。図６Ｂ、７Ｂ、８Ｂは、それぞれ図５と図６Ａ、７Ａ、８Ａとの間の差分画像を、三枚示している。これらの差分画像は、２４個のサブセットでのＯＳ−ＳＱＳ法が六回の逐次では完全には収束しないことを実証している。明らかに、図６Ｂにおけるかなりの構造と詳細とが、多くの高周波数特徴と同様に目立っている。ＯＳ−ＳＱＳ−ｍｏｍアルゴリズムを使用して再構成された画像も、発散から生じるかなりのノイズを示している。他方で、ＯＳ―ＳＱＳ―ｍｏｍ−ＶＲアルゴリズムを使用する方法２００は、たった六回の逐次後でも参照画像からの差の提示はごく僅かである。

図９Ａは、頭部画像についての参照画像を示している。図９Ｂは、ＯＳ−ＳＱＳ−ｍｏｍアルゴリズムを使用して十回の逐次後の再構成画像を示しており、図９Ｃは、ＯＳ−ＳＱＳ−ｍｏｍ−ＶＲアルゴリズムを使用して四回の逐次後の再構成画像を示している。図９Ａ、９Ｂ、９Ｃにおける画像のそれぞれは、同じＦＢＰ画像から開始したものである。図９ＢにおいてＯＳ−ＳＱＳ−ｍｏｍアルゴリズムを使用して生成された十回の逐次後の画像であっても、リミットサイクルが原因で、まだ収束しておらずやはりノイジーに見える。図９Ｂとは対照的に、たった四回の逐次後に、ＯＳ−ＳＱＳ−ｍｏｍ−ＶＲアルゴリズム（つまり、図９Ｃ）を使用して生成された画像は、図９Ａの参照画像と見た目としては（ｖｉｓｕａｌｌｙ）より近いものとなる。

図１０は、エラー測定を示しており、図１Ｂのプロットと同じく、サブセット数とＯＳとの両方、そして１フルセットの勾配計算Ｔの期間が変えられている。上で説明した通り、計算コストを抑えるために、１フルセットの勾配アップデートは、外側のループ処理の各逐次に対する１フルセットの勾配を計算することとは対照的に、周期的に計算するだけで良い。従って、図１０において、１フルセットの勾配は、毎Ｔ逐次にアップデートされ、ここでＴは１、２、３に等しい。当該対応するＲＭＳＤの結果は、図１０において表示された。所定数のサブセット、例えばＯＳ＝１２、そして２４に対し、１から３まで変化するＴは、ＲＭＳＤにほぼ何も変化がないことが見て取れる。３０のサブセットが使用された場合に、ＲＭＳＤにおける差はわずかであり、このことは図１０において逐次が五、六回辺りで観察することが出来る。しかし、当該差は、逐次を多数重ねると、図１０で見ることが出来ない。従って、計算的な条件は、逐次数の関数として収束率を著しく減らすことなく、Ｔ＞1を使用することで、縮らすことが出来る。

図１１は、ＣＴ装置またはＣＴスキャナに含まれる放射線ガントリの実施例を描いている。なお、図２に示した分散が縮小されたＯＳ ＩＲ法は、上述したように医用画像生成装置によって実現される。図１１に示したＣＴ装置は、例えば、システムコントローラ５１０（の一部の機能）、データ／制御バス５１１、格納部５１２は、再構成デバイス５１４、入力部５１５、表示部５１６等によって、上述したＯＳ ＩＲ法を実施す医用画像生成装置を実現する。図１１に図示されるように、放射線ガントリ５００は側面から見て描かれており、Ｘ線管（Ｘ線発生装置）５０１、環状フレーム５０２、そして多列または２次元アレイ型Ｘ線検出器５０３を更に含む。Ｘ線管５０１及びＸ線検出器５０３は、環状フレーム５０２上に被検体ＯＢＪを横切って正反対に取り付けられ、環状フレーム５０２は回転軸ＲＡの回りに回転可能に支持される。被検体ＯＢＪが図示された頁の奥の方向または手前の方向の軸ＲＡに沿って移動されながら、回転ユニット５０７は環状フレーム５０２を０．４秒／回転もの高速で回転させる。

本発明に係るＸ線コンピュータ断層撮影（ＣＴ）装置の第一の実施形態は、付随する図面を参照しながら以下に説明される。Ｘ線ＣＴ装置は、様々なタイプの装置を含むことに留意されたい。具体的には、Ｘ線管とＸ線検出器とが検査される予定の被検体の周辺を一緒に回る回転／回転型機構と、そして多数の検出器素子がリング状または水平状に配置されており、Ｘ線管のみが検査される予定の被検体の周辺を回る固定／回転型機構とがある。本開示は、いずれのタイプにも適用可能である。今回は、現在の主流である回転／回転型機構が例示される。

マルチスライスＸ線ＣＴ装置は高電圧発生器５０９を更に含み、Ｘ線管５０１がＸ線を生成するように、スリップリング５０８を通して、高電圧発生器５０９はＸ線管５０１に印加される管電圧を生成する。Ｘ線は、被検体ＯＢＪに向かって照射され、被検体ＯＢＪの断面領域が円で表される。例えば、第一のスキャンにわたる平均的なＸ線エネルギーを有するＸ線管５０１は、第二のスキャンにわたる平均的なＸ線エネルギーよりも小さい。このようにして、二回以上のスキャンが異なるＸ線エネルギーに対応して、取得することが出来る。Ｘ線検出器５０３は、被検体ＯＢＪを通り抜けてきた照射Ｘ線を検出するために、被検体ＯＢＪを挟んでＸ線管５０１から反対側の位置にある。Ｘ線検出器５０３は、個々の検出器素子または検出器ユニットを更に含む。

ＣＴ装置は、Ｘ線検出器５０３から検出された信号を処理するための、その他のデバイスを更に含む。データ収集回路またはデータ収集システム（ＤＡＳ）５０４は、各チャンネルに対するＸ線検出器５０３からの出力信号を電圧信号に変換し、その電圧信号を増幅し、更にその電圧信号をデジタル信号へと変換する。Ｘ線検出器５０３およびＤＡＳ５０４は、１回転当たりの所定全投影数（ＴＰＰＲ）を処理するよう構成されている。

上に説明されたデータは、非接触データ送信装置５０５を通して、放射線ガントリ５００外部のコンソール内に収容された、前処理デバイス５０６へと送られる。前処理デバイス５０６は、ローデータに関する感度補正など、特定の補正を実行する。格納部５１２は、再構成処理直前のステージで、投影データとも呼ばれる結果データを格納する。格納部５１２は、再構成デバイス５１４、入力部５１５、表示部５１６と共に、データ／制御バス５１１を通して、システムコントローラ５１０に接続されている。システムコントローラ５１０は、ＣＴシステムを駆動させるのに十分なレベルに達するまで電流を制限する電流調整器５１３を制御する。
検出器は、どんな世代のＣＴスキャナシステムであっても、患者に対して回転および／または固定される。一実施形態において、上に説明されたＣＴシステムは、第三世代ジオメトリシステムと第四世代ジオメトリシステムとが組み合わせられた例であってもよい。第三世代ジオメトリシステムにおいて、Ｘ線管５０１とＸ線検出器５０３とは、環状フレーム５０２上に正反対に取り付けられ、環状フレーム５０２が回転軸ＲＡの周りを回転する時に、被検体ＯＢＪの周りを回転する。第四世代ジオメトリシステムにおいて、検出器は患者の周辺に固定して取り付けられており、Ｘ線管は患者の周辺を回転する。代替的な実施形態において、放射線ガントリ５００は、Ｃアームおよびスタンドによって支持されている、環状フレーム５０２上に配置された多数の検出器を有する。

格納部５１２は、Ｘ線検出器５０３でＸ線照射量を示す測定値を格納することが出来る。更に、格納部５１２は、ＣＴ画像再構成用の方法２００を実行するための専用プログラムを格納することが出来る。

再構成デバイス５１４は、方法２００を実行することが出来る。更に、再構成デバイス５１４は、必要に応じてボリュームレンダリング処理や画像差分処理など、画像処理を実行することが出来る。

前処理デバイス５０６によって実行された投影データの前再構成処理は、例えば検出器キャリブレーション、検出器非直線性、極性効果のための補正を含むことが出来る。

再構成デバイス５１４によって実行される後再構成処理は、画像のフィルタリングや平滑化、ボリュームレンダリング処理、そして画像差分処理を、必要に応じて含むことが出来る。画像再構成処理は、方法２００を実行出来る。再構成デバイス５１４は、例えば投影データ、再構成画像、キャリブレーションデータやパラメータ、そしてコンピュータプログラムを格納するのに格納部を使うことが出来る。

再構成デバイス５１４は、特定用途向け集積回路(ＡＳＩＣ)、フィールドプログラマブルゲートアレイ(ＦＰＧＡ)、または複合プログラマブル論理デバイス(ＣＰＬＤ)など、個々の論理ゲートとして実行可能な、ＣＰＵ（処理回路）を含むことが出来る。ＦＰＧＡまたはＣＰＬＤ実行は、ＶＨＤＬ、ベリログ、またはその他のハードウェア記述言語でコード化されていてもよく、そして当該コードはＦＰＧＡまたはＣＰＬＤにおいて直接電子メモリ内に格納されてもよいし、あるいは個別の電子メモリとして格納されてもよい。更に、格納部５１２は、ＲＯＭ、ＥＰＲＯＭ、ＥＥＰＲＯＭ（登録商標）、またはＦＬＡＳＨメモリなど、不揮発性メモリであってもよい。格納部５１２は、静的または動的ＲＡＭなど揮発性でよく、電子メモリおよびＦＰＧＡまたはＣＰＬＤと格納部との間の相互作用を管理するマイクロコントローラやマイクロプロセッサなどプロセッサ（処理回路）が提供されていてもよい。

代替的に、再構成デバイス５１４におけるＣＰＵは、本実施形態で説明された機能を実行するコンピュータ読み取り可能命令のセットを含むコンピュータプログラムを実行することが出来、当該コンピュータプログラムは、任意の上述の非一時的電子メモリおよび／またはハードディスクドライブ、ＣＤ、ＤＶＤ、ＦＬＡＳＨドライブ、またはその他の任意の既知の格納媒体に格納されている。更に、コンピュータ読み取り可能命令は、ユーティリティアプリケーション、バックグラウンドデーモン、またはオペレーティングシステムの構成要素、またはそれらの組み合わせで提供されてもよく、所定のオペレーティングシステムや、当業者にとっては既知のその他のオペレーティングシステムと一体となって実行する。更に、ＣＰＵは、命令を実行するために並行して協同的に動く、マルチプルプロセッサとして実行されてもよい。

一実行において、再構成画像は、表示部５１６上に映し出されてよい。当該表示部５１６は、ＬＣＤディスプレイ、ＣＲＴディスプレイ、プラズマディスプレイ、ＯＬＥＤ、ＬＥＤ、または当業者にとって既知のその他のディスプレイであってもよい。
格納部５１２は、ハードディスクドライブ、ＣＤ-ＲＯＭドライブ、ＤＶＤドライブ、ＦＬＡＳＨドライブ、ＲＡＭ、ＲＯＭ、または当業者にとって既知のその他の格納メディアであってもよい。

特定の実施形態について説明してきたが、これらの実施形態は、一例として提示したに過ぎず、本発明の範囲を限定することを意図したものではない。実際に、本明細書で説明された新規の方法およびシステムはさまざまなその他の形態で具体化することが出来る。更には、本開示の精神から乖離することなく、本開示に説明された方法、装置やシステムの形式で省略、置き換え、変更が可能である。

５００…放射線ガントリ
５０１…Ｘ線管
５０２…環状フレーム
５０３…Ｘ線検出器
５０４…データ収集システム
５０５…非接触データ送信装置
５０６…前処理デバイス
５０７…回転ユニット
５０８…スリップリング
５０９…高電圧発生器
５１０…システムコントローラ
５１１…データ／制御バス
５１２…格納部
５１３…電流調整器
５１４…再構成デバイス
５１５…入力部
５１６…表示部

Claims (10)

1. Ｘ線発生装置から照射され被検体を透過した透過Ｘ線を用いて検出された、前記被検体に関する検出データを取得する取得部と、
   前記検出データを複数のサブセットに分けて、前記複数のサブセットのうち２以上のサブセットに基づいたオブジェクト関数における第一の勾配情報を取得し、前記サブセットそれぞれに基づいたオブジェクト関数における第二の勾配情報を前記第一の勾配情報に基づき補正し、当該補正後の第二の勾配情報に基づいたオブジェクト関数を用いて前記サブセット毎に逐次近似再構成処理を行い、再構成画像を取得する再構成部と、
   を具備する医用画像生成装置。
2. 前記再構成部は、前記複数のサブセット全てに基づいたオブジェクト関数における前記第一の勾配情報を取得し、前記サブセットそれぞれに基づいたオブジェクト関数における第二の勾配情報を前記第一の勾配情報に基づき補正し、当該補正後の第二の勾配情報に基づいたオブジェクト関数を用いて前記各サブセット毎に前記逐次近似再構成処理を行い、前記再構成画像を取得する請求項１記載の医用画像生成装置。
3. 前記再構成部は、前記第一の勾配情報の取得に用いた前記２以上のサブセットのみに対して、当該２以上のサブセットに基づいたオブジェクト関数における前記第二の勾配情報を前記第一の勾配情報に基づき補正し、当該補正後の第二の勾配情報に基づいたオブジェクト関数を用いて前記２以上のサブセットについて前記逐次近似再構成処理を行い、前記再構成画像を取得する請求項１記載の医用画像生成装置。
4. 前記再構成部は、
   前記第二の勾配情報と前記各サブセットについて前記再構成画像の前記アップデートとを実行することで、ボクセルに関するループを繰り返し、
   前記ボクセルに関するループの各逐次は、過去の逐次からの再構成画像を使用して前記再構成画像の前記アップデートを実行する、
   請求項１記載の医用画像生成装置。
5. 前記再構成部は、
   前記サブセットに関するループを繰り返し、
   前記サブセットに関するループの各逐次において、前記再構成画像から別の再構成画像を生成し、
   前記別の再構成画像を使用して前記第一の勾配情報を取得し、
   前記第二の勾配情報を補正するために、前記別の再構成画像を使用して前記サブセットに関するループを繰り返す、
   請求項４記載の医用画像生成装置。
6. 前記再構成部は、
   前記ボクセルに関するループを繰り返し、
   前記ボクセルに関するループの各逐次において、前記再構成画像と等しい別の再構成画像を生成し、
   前記各サブセットに対応する前記サブセットに関するループの各逐次において、前記サブセットに関するループを繰り返し、
   所定の逐次期間の開始で発生する前記ボクセルに関するループの各逐次において、前記別の再構成画像を使用して前記第一の勾配情報を取得し、
   前記所定の逐次期間の開始で発生しない前記ボクセルに関するループの各逐次において、前記第一の勾配の前記計算を省略する、
   請求項４記載の医用画像生成装置。
7. 前記再構成部は、
   前記再構成画像と等しい別の再構成画像を初期化し、
   前記第一の勾配情報を取得し、
   前記第一の勾配情報を取得するために使用される前記一つ以上のサブセットのオブジェクト関数は前記別の再構成画像を含み、
   前記第二の勾配情報を補正するために使用される前記一つのサブセットの前記オブジェクト関数は、前記別の再構成画像を含む、
   請求項４記載の医用画像生成装置。
8. 前記再構成部は、フィルタ補正逆投影法とフェルドカンプ再構成法とのうちの一方を使用して前記再構成画像を初期化する請求項１乃至７のうちいずれか一項記載の医用画像生成装置。
9. 前記再構成部は、分離可能二次サロゲート法とネステロフの加速法との組み合わせ、又はネステロフの加速法を使用して、前記再構成画像をアップデートする請求項１乃至８のうちいずれか一項記載の医用画像生成装置。
10. Ｘ線発生装置から照射され被検体を透過した透過Ｘ線を用いて検出された、前記被検体に関する検出データを取得し、
    前記検出データを複数のサブセットに分けて、前記複数のサブセットのうち２以上のサブセットに基づいたオブジェクト関数における第一の勾配情報を取得し、
    前記サブセットそれぞれに基づいたオブジェクト関数における第二の勾配情報を前記第一の勾配情報に基づき補正し、
    戦記補正後の第二の勾配情報に基づいたオブジェクト関数を用いて前記サブセット毎に逐次近似再構成処理を行い、再構成画像を取得すること、
    を具備する医用画像生成方法。