知財求人 - 知財ポータルサイト「IP Force」
▶ ユーエムセー ユトレヒト ホールディング ベイフェイの特許一覧
(19)【発行国】日本国特許庁(JP)
(12)【公報種別】特許公報(B2)
(11)【特許番号】
(24)【登録日】2024-05-16
(45)【発行日】2024-05-24
(54)【発明の名称】タイムドメイン磁気共鳴のためのパラメータマップ決定
(51)【国際特許分類】
A61B 5/055 20060101AFI20240517BHJP
G01N 24/00 20060101ALI20240517BHJP
【ＦＩ】
A61B5/055 376
A61B5/055 311
G01N24/00 530G
G01N24/00 530Y
【請求項の数】
17
(21)【出願番号】P 2021559438
(86)(22)【出願日】2020-04-08
(65)【公表番号】
(43)【公表日】2022-06-02
(86)【国際出願番号】 EP2020059997
(87)【国際公開番号】W WO2020208066
(87)【国際公開日】2020-10-15
【審査請求日】2021-12-02
(31)【優先権主張番号】2022890
(32)【優先日】2019-04-08
(33)【優先権主張国・地域又は機関】NL
(73)【特許権者】
【識別番号】520099210
【氏名又は名称】ユーエムセー ユトレヒト ホールディング ベイフェイ
【氏名又は名称原語表記】ＵＭＣ Ｕｔｒｅｃｈｔ Ｈｏｌｄｉｎｇ Ｂ．Ｖ．
(74)【代理人】
【識別番号】100080816
【弁理士】
【氏名又は名称】加藤 朝道
(74)【代理人】
【識別番号】100098648
【弁理士】
【氏名又は名称】内田 潔人
(72)【発明者】
【氏名】ファン デル ハイデ、 オスカー
(72)【発明者】
【氏名】スブリッジ、アレッサンドロ
(72)【発明者】
【氏名】ファン デン ベルグ、コルネリス アントニウス テオドルス
【審査官】永田 浩司
(56)【参考文献】
【文献】Alessandro Sbrizzi, et al.，Fast quantitative MRI as a nonlinear tomography problem，arXiv，2017年，1705.03209v2，1-22
【文献】Pierre Ablin, et al.，Faster Independent Component Analysis by Preconditioning With Hessian Approximations，IEEE TRANSACTIONS ON SIGNAL PROCESSING，2018年，VOL.66, N.15，pp.4040-4049
【文献】ANDREA WALTHER，Computing Sparse Hessians with Automatic Differentiation，ACM Transactions on Mathematical Software，2008年，Vol. 34, No. 1, Article 3，1-15
【文献】Ananya Panda, et al.，Magnetic Resonance Fingerprinting-An Overview，Curr Opin Biomed Eng.，2017年，3，pp.56-66
【文献】Gopal Nataraj, et al.，Dictionary-Free MRI PERK: Parameter Estimation via Regression with Kernels，IEEE TRANSACTIONS ON MEDICAL IMAGING，2018年，VOL.37, NO.9，pp.2103-2114
(58)【調査した分野】(Int.Cl.，ＤＢ名)
Ａ６１Ｂ ５／０５５
Ｇ０１Ｎ ２４／００
(57)【特許請求の範囲】
【請求項1】
印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する方法であって、該方法は、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定すること、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定すること；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定すること；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定すること；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットを用いてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得すること；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出すること
を含むこと
を特徴とする方法。
【請求項2】
請求項１に記載の方法において、前記印加パルスシーケンスは勾配エンコードパターン及び／又は無線周波数励起パターンを含むこと、及び、
前記スパースヘッセは複数の対角バンドを含み、該複数の対角バンドの何れか１つの幅は、前記印加パルスシーケンスの前記勾配エンコードパターン及び／又は前記無線周波数励起パターンによって決定されること
を特徴とする方法。
【請求項3】
請求項１又は２に記載の方法において、前記スパースヘッセを計算することは、ヘッセ行列の無視不可能要素を決定することを含むこと、好ましくは該ヘッセ行列の要素は、各要素の絶対値が閾値を上回る場合、無視不可能として決定されること、より好ましくは該ヘッセ行列の要素は、前記残差ベクトルのヤコビアンと前記残差ベクトルのヤコビアンのエルミート転置の積の絶対値が所定の値より大きい場合、無視不可能として決定されること
を特徴とする方法。
【請求項4】
請求項２又は３に記載の方法において、前記複数の対角バンドの何れか１つの幅は１より大きいこと、好ましくは前記複数の対角バンドの何れか１つの幅は２～５５の範囲、より好ましくは３～８の範囲にあり、最も好ましくは３であること
を特徴とする方法。
【請求項5】
請求項１～４の何れかに記載の方法において、前記ヘッセは、前記残差ベクトルのヤコビアンと前記残差ベクトルのヤコビアンのエルミート転置の積として計算されること
を特徴とする方法。
【請求項6】
請求項１～５の何れかに記載の方法において、前記印加パルスシーケンスは、カーテシアン（Cartesian）収集、ラジアル収集又はスパイラル収集の何れか１つをもたらすよう構成されていること
を特徴とする方法。
【請求項7】
請求項６に記載の方法において、前記印加パルスシーケンスの勾配エンコードパターンは、対応する点広がり関数が位相エンコード方向において専ら広がるように、カーテシアン収集をもたらすよう構成されていること
を特徴とする方法。
【請求項8】
請求項１～７の何れかに記載の方法において、前記印加パルスシーケンスは、変化するフリップ角をもたらすよう構成されていること
を特徴とする方法。
【請求項9】
請求項７に記載の方法において、前記印加パルスシークエンスの無線周波数励起パターンは、対応する点広がり関数が幅方向において空間的に限定されるように、滑らかに変化するフリップ角をもたらすよう構成されていること
を特徴とする方法。
【請求項10】
請求項１～９の何れかに記載の方法において、前記少なくとも１つの組織パラメータは、Ｔ１緩和時間、Ｔ２緩和時間、Ｔ２＊緩和時間及びプロトン密度の何れか１つ又はこれらの組み合わせを含むこと
を特徴とする方法。
【請求項11】
請求項１～１０の何れかに記載の方法において、前記残差ベクトルの勾配を計算することは、該残差ベクトルと該残差ベクトルのヤコビアンのエルミート転置を乗算することを含むこと
を特徴とする方法。
【請求項12】
請求項１～１１の何れかに記載の方法において、前記ＴＤＭＲ信号モデルは、ブロッホ型体積信号モデルであること
を特徴とする方法。
【請求項13】
請求項２～１２の何れかに記載の方法において、前記ヘッセにおける複数の対角バンドは、前記印加パルスシーケンスの勾配エンコードパターンによって決定されること
を特徴とする方法。
【請求項14】
請求項１～１３の何れかに記載の方法において、前記スパースヘッセはＴＤＭＲ信号モデルパラメータの各セットについて一定であること
を特徴とする方法。
【請求項15】
請求項１～１４の何れかに記載の方法において、前記更新のステップvii) は、
スパースヘッセにＴＤＭＲ信号モデルパラメータ変化を乗じたものが残差関数の勾配のマイナスに等しいという方程式の線形システムを解くことによって、ＴＤＭＲ信号モデルパラメータ変化のセットを取得すること、及び、
取得したＴＤＭＲ信号モデルパラメータ変化を、ＴＤＭＲ信号モデルパラメータの更新セットを取得するために更新されるべきＴＤＭＲ信号モデルパラメータのセットに加算すること
を含むこと
を特徴とする方法。
【請求項16】
印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する装置であって、該装置は、下記ステップを実行するよう構成されたプロセッサを含むこと、
即ち、該プロセッサは、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定する、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定する；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定する；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定する；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得する；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出する
よう構成されること、
を特徴とする装置。
【請求項17】
プログラムがコンピュータで実行されるとき、請求項１～１５の何れかに記載の方法を行うためのコンピュータで実行可能な命令を含むコンピュータプログラム製品。
【発明の詳細な説明】
【技術分野】
【0001】
本特許開示は、印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する方法及び装置、及び、該方法を行うためのコンピュータプログラム製品に関する。
【背景技術】
【0002】
磁気共鳴画像法（ＭＲＩ）は、多くの応用に使用され、かつ、調整可能な多くのシーケンスパラメータ及び例えば異なる種類の生物学的情報を抽出することが認められ得る多くのイメージングパラメータを有するイメージングモダリティである。伝統的なＭＲＩ画像再構成はｋ空間信号の収集及び収集したデータに対する逆高速フーリエ変換の実行を伴う。伝統的ＭＲＩイメージングは遅い。なぜなら、測定されるべき全てのパラメータ（例えばＴ１又はＴ２）について、別々のＭＲＩ測定（値）が異なる設定を有するＭＲＩ装置によって収集されなければならないからである。スキャンには３０～４５分間もかかり得る。
【0003】
タイムドメイン磁気共鳴スピントモグラフィ（ＭＲ－ＳＴＡＴ）はタイムドメインデータから直接的にＭＲ画像を取得するための定量的方法である。とりわけ、ＭＲ－ＳＴＡＴは単回（シングル）ショートスキャンを用いてマルチパラメトリック定量的ＭＲマップを取得するためのフレームワークである。
【0004】
ＷＯ ２０１６／１８４７７９ Ａ１に記載されているように、非厳密（inexact）ガウス・ニュートン法を用いてブロッホ（Bloch）型体積信号モデルをタイムドメインデータに直接的にフィッティングさせることによって信号の空間的局在化（localization）と組織パラメータの推定（estimation）が同時に実行される大規模最適化問題は解決される。高解像度スキャンについての大規模最適化問題を解決するために、高度に並列化されたマトリックスフリーの非厳密ガウス・ニュートン再構成アルゴリズムを使用することができる。
【0005】
G. Wuebbeler et al.の"A Large-Scale Optimization Method Using a Sparse Approximation of the Hessian for Magnetic Resonance Fingerprinting", SIAM JOURNAL ON IMAGING SCIENCES, vol. 10, nr. 3, 18 July 2017, pp. 979-1004 には、磁気共鳴フィンガープリンティング（Magnetic Resonance Fingerprinting：ＭＲＦ）に対する最小二乗アプローチが記載されており、この場合、ヘッセ（Hessian）のスパース近似が使用される。
【先行技術文献】
【特許文献】
【0006】
【文献】ＷＯ ２０１６／１８４７７９ Ａ１
【非特許文献】
【0007】
【文献】G. Wuebbeler et al., "A Large-Scale Optimization Method Using a Sparse Approximation of the Hessian for Magnetic Resonance Fingerprinting", SIAM JOURNAL ON IMAGING SCIENCES, vol.10, nr.3, 18 July 2017, pp.979-1004
【発明の概要】
【発明が解決しようとする課題】
【0008】
測定時間は数分のオーダーに大きく減少されるが、信号モデルをタイムドメインデータにフィッティングするための演算時間は、単一の２Ｄスライスについて凡そ１時間もかかる。
【0009】
本特許開示の複数の課題のうちの１つの課題は、タイムドメインＭＲ信号の定量的ＭＲマップ（複数）への変換を改善することである。
【課題を解決するための手段】
【0010】
本発明の第１の視点において、印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する方法が提供される。
該方法は、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定すること、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定すること；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定すること；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定すること；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットを用いてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得すること；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出すること
を含むこと
を特徴とする（形態１）。
本発明の第２の視点により、印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する装置であって、
該装置は、下記ステップを実行するよう構成されたプロセッサを含むこと、
即ち、該プロセッサは、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定する、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定する；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定する；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定する；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得する；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出する
よう構成されること、
を特徴とする（形態１６）。
本発明の第３の視点により、プログラムがコンピュータで実行されるとき、本発明の方法を行うためのコンピュータで実行可能な命令を含むコンピュータプログラム製品が提供される（形態１７）。
該方法は、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定すること、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定すること；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定すること；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定すること；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットを用いてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得すること；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出すること
を含むこと
を特徴とする（形態１）。
本発明の第２の視点により、印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する装置であって、
該装置は、下記ステップを実行するよう構成されたプロセッサを含むこと、
即ち、該プロセッサは、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定する、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定する；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定する；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定する；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得する；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出する
よう構成されること、
を特徴とする（形態１６）。
本発明の第３の視点により、プログラムがコンピュータで実行されるとき、本発明の方法を行うためのコンピュータで実行可能な命令を含むコンピュータプログラム製品が提供される（形態１７）。
【発明を実施するための形態】
【0011】
ここに本発明の好ましい形態を示す。
（形態１）上記本発明の第１の視点参照。
（形態２）形態１に記載の方法において、
前記印加パルスシーケンスは勾配エンコードパターン及び／又は無線周波数励起パターンを含むこと、及び、
前記スパースヘッセは複数の対角バンドを含み、該複数の対角バンドの何れか１つの幅は、前記印加パルスシーケンスの前記勾配エンコードパターン及び／又は前記無線周波数励起パターンによって決定されることが好ましい。
（形態３）形態１又は２に記載の方法において、
前記スパースヘッセを計算することは、ヘッセ行列の無視不可能要素を決定することを含むこと、好ましくは該ヘッセ行列の要素は、各要素の絶対値が閾値を上回る場合、無視不可能として決定されること、より好ましくは該ヘッセ行列の要素は、前記残差ベクトルのヤコビアンと前記残差ベクトルのヤコビアンのエルミート転置の積の絶対値が所定の値より大きい場合、無視不可能として決定されることが好ましい。
（形態４）形態２又は３に記載の方法において、
前記複数の対角バンドの何れか１つの幅は１より大きいこと、好ましくは前記複数の対角バンドの何れか１つの幅は２～５５の範囲、より好ましくは３～８の範囲にあり、最も好ましくは３であることが好ましい。
（形態５）形態１～４の何れかに記載の方法において、
前記ヘッセは、前記残差ベクトルのヤコビアンと前記残差ベクトルのヤコビアンのエルミート転置の積として計算されることが好ましい。
（形態６）形態１～５の何れかに記載の方法において、
前記印加パルスシーケンスは、カーテシアン（Cartesian）収集、ラジアル収集又はスパイラル収集の何れか１つをもたらすよう構成されていることが好ましい。
（形態７）形態６に記載の方法において、
前記印加パルスシーケンスの勾配エンコードパターンは、対応する点広がり関数が位相エンコード方向において専ら広がるように、カーテシアン収集をもたらすよう構成されていることが好ましい。
（形態８）形態１～７の何れかに記載の方法において、
前記印加パルスシーケンスは、変化するフリップ角をもたらすよう構成されていることが好ましい。
（形態９）形態７に記載の方法において、
前記印加パルスシークエンスの無線周波数励起パターンは、対応する点広がり関数が幅方向において空間的に限定されるように、滑らかに変化するフリップ角をもたらすよう構成されていることが好ましい。
（形態１０）形態１～９の何れかに記載の方法において、
前記少なくとも１つの組織パラメータは、Ｔ１緩和時間、Ｔ２緩和時間、Ｔ２ ＊ 緩和時間及びプロトン密度の何れか１つ又はこれらの組み合わせを含むことが好ましい。
（形態１１）形態１～１０の何れかに記載の方法において、
前記残差ベクトルの勾配を計算することは、該残差ベクトルと該残差ベクトルのヤコビアンのエルミート転置を乗算することを含むことが好ましい。
（形態１２）形態１～１１の何れかに記載の方法において、
前記ＴＤＭＲ信号モデルは、ブロッホ型体積信号モデルであることが好ましい。
（形態１３）形態２～１２の何れかに記載の方法において、
前記ヘッセにおける複数の対角バンドは、前記印加パルスシーケンスの勾配エンコードパターンによって決定されることが好ましい。
（形態１４）形態１～１３の何れかに記載の方法において、
前記スパースヘッセはＴＤＭＲ信号モデルパラメータの各セットについて一定であることが好ましい。
（形態１５）形態１～１４の何れかに記載の方法において、
前記更新のステップvii) は、
スパースヘッセにＴＤＭＲ信号モデルパラメータ変化を乗じたものが残差関数の勾配のマイナスに等しいという方程式の線形システムを解くことによって、ＴＤＭＲ信号モデルパラメータ変化のセットを取得すること、及び、
取得したＴＤＭＲ信号モデルパラメータ変化を、ＴＤＭＲ信号モデルパラメータの更新セットを取得するために更新されるべきＴＤＭＲ信号モデルパラメータのセットに加算すること
を含むことが好ましい。
（形態１６）上記本発明の第２の視点参照。
（形態１７）上記本発明の第３の視点参照。
（形態１）上記本発明の第１の視点参照。
（形態２）形態１に記載の方法において、
前記印加パルスシーケンスは勾配エンコードパターン及び／又は無線周波数励起パターンを含むこと、及び、
前記スパースヘッセは複数の対角バンドを含み、該複数の対角バンドの何れか１つの幅は、前記印加パルスシーケンスの前記勾配エンコードパターン及び／又は前記無線周波数励起パターンによって決定されることが好ましい。
（形態３）形態１又は２に記載の方法において、
前記スパースヘッセを計算することは、ヘッセ行列の無視不可能要素を決定することを含むこと、好ましくは該ヘッセ行列の要素は、各要素の絶対値が閾値を上回る場合、無視不可能として決定されること、より好ましくは該ヘッセ行列の要素は、前記残差ベクトルのヤコビアンと前記残差ベクトルのヤコビアンのエルミート転置の積の絶対値が所定の値より大きい場合、無視不可能として決定されることが好ましい。
（形態４）形態２又は３に記載の方法において、
前記複数の対角バンドの何れか１つの幅は１より大きいこと、好ましくは前記複数の対角バンドの何れか１つの幅は２～５５の範囲、より好ましくは３～８の範囲にあり、最も好ましくは３であることが好ましい。
（形態５）形態１～４の何れかに記載の方法において、
前記ヘッセは、前記残差ベクトルのヤコビアンと前記残差ベクトルのヤコビアンのエルミート転置の積として計算されることが好ましい。
（形態６）形態１～５の何れかに記載の方法において、
前記印加パルスシーケンスは、カーテシアン（Cartesian）収集、ラジアル収集又はスパイラル収集の何れか１つをもたらすよう構成されていることが好ましい。
（形態７）形態６に記載の方法において、
前記印加パルスシーケンスの勾配エンコードパターンは、対応する点広がり関数が位相エンコード方向において専ら広がるように、カーテシアン収集をもたらすよう構成されていることが好ましい。
（形態８）形態１～７の何れかに記載の方法において、
前記印加パルスシーケンスは、変化するフリップ角をもたらすよう構成されていることが好ましい。
（形態９）形態７に記載の方法において、
前記印加パルスシークエンスの無線周波数励起パターンは、対応する点広がり関数が幅方向において空間的に限定されるように、滑らかに変化するフリップ角をもたらすよう構成されていることが好ましい。
（形態１０）形態１～９の何れかに記載の方法において、
前記少なくとも１つの組織パラメータは、Ｔ１緩和時間、Ｔ２緩和時間、Ｔ２ ＊ 緩和時間及びプロトン密度の何れか１つ又はこれらの組み合わせを含むことが好ましい。
（形態１１）形態１～１０の何れかに記載の方法において、
前記残差ベクトルの勾配を計算することは、該残差ベクトルと該残差ベクトルのヤコビアンのエルミート転置を乗算することを含むことが好ましい。
（形態１２）形態１～１１の何れかに記載の方法において、
前記ＴＤＭＲ信号モデルは、ブロッホ型体積信号モデルであることが好ましい。
（形態１３）形態２～１２の何れかに記載の方法において、
前記ヘッセにおける複数の対角バンドは、前記印加パルスシーケンスの勾配エンコードパターンによって決定されることが好ましい。
（形態１４）形態１～１３の何れかに記載の方法において、
前記スパースヘッセはＴＤＭＲ信号モデルパラメータの各セットについて一定であることが好ましい。
（形態１５）形態１～１４の何れかに記載の方法において、
前記更新のステップvii) は、
スパースヘッセにＴＤＭＲ信号モデルパラメータ変化を乗じたものが残差関数の勾配のマイナスに等しいという方程式の線形システムを解くことによって、ＴＤＭＲ信号モデルパラメータ変化のセットを取得すること、及び、
取得したＴＤＭＲ信号モデルパラメータ変化を、ＴＤＭＲ信号モデルパラメータの更新セットを取得するために更新されるべきＴＤＭＲ信号モデルパラメータのセットに加算すること
を含むことが好ましい。
（形態１６）上記本発明の第２の視点参照。
（形態１７）上記本発明の第３の視点参照。
【0012】
第１の側面により、印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する方法、好ましくはコンピュータで実行される方法が提供される。該方法は、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定すること、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータ（複数）に依存（従属）するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定すること（estimating）；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定すること；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配（gradient）を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセ（sparse Hessian）を計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定すること；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返し（反復）が完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返し（反復し）、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得すること；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出すること
を含む。
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定すること、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータ（複数）に依存（従属）するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定すること（estimating）；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定すること；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配（gradient）を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセ（sparse Hessian）を計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定すること；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返し（反復）が完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返し（反復し）、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得すること；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出すること
を含む。
【0013】
ステップiii)、iv)、v)、vi) 及びvii) を繰り返す場合、ＴＤＭＲ信号モデルパラメータの初期セットの代わりに、ＴＤＭＲ信号モデルパラメータの更新セットが使用される。この方法で、ＴＤＭＲ信号モデルパラメータの更新セットは、ＴＤＭＲ信号モデルパラメータの最終更新セットへと反復的に収斂する。一般的には、所定の閾値より小さい上記の差及び所定回数の繰り返しが例を表す任意の停止基準が使用され得る。従って、１５０ステップの繰り返しが、停止基準が満たされるまで、行われ得る。
【0014】
「スパース（sparse）」の語は、本開示においては、行列（matrix）であるヘッセ（Hessian）がゼロ又は無視可能である要素（複数）を（例えば要素の総数の半分より多く）含むことを示すために使用される。例えば、スパースヘッセのスパース性は、ヘッセのゼロ値の及び無視可能な要素の数を該ヘッセの要素の総数で除算したものである。ステップvi)において計算されたヘッセのスパース性は印加パルスシーケンスに基づいている。より具体的には、ヘッセは、非ゼロ及び無視不可能（無視できない：non-negligible）要素（複数）の対角（ダイアゴナル）バンド（複数）を含むスパース性パターンのような、スパース性パターンを有する。印加パルスシーケンスは、ヘッセが何らかのスパース性パターンを有するスパースヘッセとなるように、選択される。
【0015】
本方法では、演算時間は、非厳密ガウス・ニュートン法と比べて少なくとも１０分の１に減少される。後者の方法では、ヘッセに近似するために、各パラメータについて信号の偏導関数（偏微分係数）を演算しなければならないが、これは多くの演算時間を必要とする。これに対し、本方法では、ヘッセは、これはスパースヘッセであるが、ＴＤＭＲ信号モデルパラメータの各セットについて予め演算されて記憶されるが、そのため、各パラメータについてのこれらの偏導関数を計算する必要はない。更に、本方法の必要記憶容量は非厳密ガウス・ニュートン法の必要記憶容量の僅か０．０４％程度であり得る。
【0016】
G. Wuebbeler et al.では、スパースヘッセの近似は、印加パルスシーケンスに依存しない（から独立の）近似である。なぜなら、
の近似の結果、スパースヘッセ行列の対角バンド当たり１つの単一パラメータのみを常に有するスパースヘッセが得られるからである。一方、本願では、スパースヘッセは印加パルスシーケンスに基づいて決定され、１より大きいバンド幅が許容される。この結果、所要計算時間は少なくとも一桁減少し、更には二桁さえも減少する。
の近似の結果、スパースヘッセ行列の対角バンド当たり１つの単一パラメータのみを常に有するスパースヘッセが得られるからである。一方、本願では、スパースヘッセは印加パルスシーケンスに基づいて決定され、１より大きいバンド幅が許容される。この結果、所要計算時間は少なくとも一桁減少し、更には二桁さえも減少する。
【0017】
なお、残差ベクトルは決定されるが、この残差ベクトルは残差関数とも称され得る。
【0018】
スパースヘッセとしてのＴＤＭＲ信号モデルのヘッセの計算において、該ヘッセのスパース性パターンは印加パルスシーケンスに基づくことが好ましい。
【0019】
一実施形態では、印加パルスシーケンスは勾配エンコード（gradient encoding）パターン及び／又は無線周波数励起パターンを含み、及び、ヘッセのスパース性パターンは勾配エンコードパターン及び／又は無線周波数励起パターンに基づいて決定される。
【0020】
ヘッセは、残差ベクトルのヤコビアンと残差ベクトルのヤコビアンのエルミート転置の積として計算される。残差ベクトルのヤコビアンの代わりに、ＴＤＭＲ信号モデルのヤコビアン、及び相応のエルミート転置、を採用することは等価である。勾配エンコードパターンは、ｋ空間サンプリング軌道、例えばカーテシアン（直交座標系：Cartesian）、ラジアル、スパイラルを決定する。
【0021】
一実施形態では、印加パルスシーケンスは、カーテシアン（式）収集（法）、ラジアル（式）収集（法）又はスパイラル（式）収集（法）の何れか１つをもたらすよう構成されている。これらのタイプのパルスシーケンスについては、スパース性パターンは、ヘッセの無視不可能要素の複数の対角バンドを含むと好適である。他の要素はゼロ又は無視可能であることが好ましい。
【0022】
一実施形態に応じ、印加パルスシーケンスの勾配エンコードパターンは、対応する点広がり関数が位相エンコード方向において専ら広がるよう、カーテシアン収集をもたらすよう構成されている。３Ｄ収集のためのもののような、２つの点広がり関数があってもよく、この場合、各点広がり関数は２つの位相エンコード方向の対応する一方に関連付けられている。これらの２つの位相エンコード方向の各々は相違する。
【0023】
一実施形態では、印加パルスシーケンスは変化するフリップ角をもたらすよう構成されている。
【0024】
好ましくは、印加パルスシークエンスの無線周波数励起パターンは、対応する点広がり関数が幅方向において空間的に限定されるよう、滑らかに変化するフリップ角をもたらすよう構成されている。滑らかに変化するとは、ＲＦ（無線周波数）励起の振幅（大きさ）が限定された量（大きさ）だけで経時的に変化するシーケンスを示し得る。１つのｋ空間の（又は各ｋ空間の）サンプリング中における２つの連続ＲＦ励起の間での変化量は、所定の量より小さく、好ましくは５度より小さい。
【0025】
より好ましくは又は代替的に、印加パルスシーケンスの無線周波数励起パターンは、対応する点広がり関数が幅方向に空間的に制限されるよう、１つのｋ空間又は複数のｋ空間の各ｋ空間のサンプリング内において滑らかに変化するフリップ角をもたらすよう構成されている。
【0026】
一実施形態では、少なくとも１つの組織パラメータは、Ｔ１緩和時間、Ｔ２緩和時間、Ｔ２＊緩和時間及びプロトン密度の何れか１つ又はこれらの組み合わせを含む。
【0027】
一実施形態では、残差ベクトルの勾配を計算することは、残差ベクトルと残差ベクトルのヤコビアンのエルミート転置を乗算することを含む。残差ベクトルのヤコビアンの代わりに、ＴＤＭＲ信号モデルのヤコビアン、及び対応するエルミート転置、を採用することは等価である。
【0028】
一実施形態では、ＴＤＭＲ信号モデルは、ブロッホ（Bloch）型体積信号モデルである。本方法は代替的な体積信号モデルによって実行可能であることは本開示を読み取れば明らかであろう。
【0029】
一実施形態では、ヘッセにおけるスパース性パターンは、印加パルスシーケンスの勾配エンコードパターンによって決定される。印加パルスシーケンスは、とりわけ好適なスパース性パターンを有する、好適にスパースなヘッセを取得するよう選択され得る。
【0030】
一実施形態では、ヘッセは、複数の対角（ダイアゴナル）バンドを含むスパースヘッセとして計算される。スパース性パターンは複数の対角バンドを含み得る。換言すれば、ヘッセの無視不可能要素（複数）の対角バンド（複数）は１つのタイプのスパース性パターンを表す。
【0031】
一実施形態では、ヘッセにおける複数の対角バンドは、印加パルスシーケンスの勾配エンコードパターンによって決定される。第１の例として、カーテシアン収集は、モデルパラメータ当たり１つの対角バンドの好適なスパース性パターンを提供することが決定される。第２の例として、ラジアル収集及びスパイラル収集はモデルパラメータ当たり４つの対角バンドの好適なスパース性パターンを提供することが決定される。
【0032】
一実施形態では、複数の対角バンドの何れか１つの幅は、印加パルスシーケンスの無線周波数励起パターンによって決定される。
【0033】
一実施形態では、複数の対角バンドの何れか１つの幅は１より大きく、好ましくは複数の対角バンドの何れか１つの幅は２～５５の範囲、より好ましくは３～８の範囲にあり、最も好ましくは３である。幅は、スパースヘッセ行列の行当たり複数の要素／係数を示し得る。
【0034】
一実施形態では、スパースヘッセはＴＤＭＲ信号モデルパラメータの各セットについて一定である。ＴＤＭＲ信号モデルパラメータ変化のセットは、スパースヘッセと残差関数の勾配に基づいて方程式（複数）の線形システム（linear system）を解くことによって得られることが好ましい。
【0035】
上記の更新のステップvii) は、
スパースヘッセにＴＤＭＲ信号モデルパラメータ変化を乗じたものが残差関数の勾配をマイナスにしたものに等しいという（スパースヘッセ×ＴＤＭＲ信号モデルパラメータ変化＝残差関数の勾配のマイナスからなる）方程式（複数）の線形システムを解くことによって、ＴＤＭＲ信号モデルパラメータ変化のセットを取得すること、及び、
取得したＴＤＭＲ信号モデルパラメータ変化を、ＴＤＭＲ信号モデルパラメータの更新セットを取得するために更新されるべきＴＤＭＲ信号モデルパラメータのセットに加算すること
を含む。
スパースヘッセにＴＤＭＲ信号モデルパラメータ変化を乗じたものが残差関数の勾配をマイナスにしたものに等しいという（スパースヘッセ×ＴＤＭＲ信号モデルパラメータ変化＝残差関数の勾配のマイナスからなる）方程式（複数）の線形システムを解くことによって、ＴＤＭＲ信号モデルパラメータ変化のセットを取得すること、及び、
取得したＴＤＭＲ信号モデルパラメータ変化を、ＴＤＭＲ信号モデルパラメータの更新セットを取得するために更新されるべきＴＤＭＲ信号モデルパラメータのセットに加算すること
を含む。
【0036】
本特許開示の第２の側面に応じ、印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する装置が提供される。該装置は、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定する、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存（従属）するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定する（estimate）；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定する；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定する；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返し（反復）が完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返し（反復し）、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得する；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出する
よう構成されたプロセッサを含む。
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定する、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存（従属）するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定する（estimate）；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定する；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定する；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返し（反復）が完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返し（反復し）、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得する；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出する
よう構成されたプロセッサを含む。
【0037】
明らかであろうように、第２の側面による装置は、とりわけ、上記の及び／又は下記の方法ステップの何れか１つ又は２つ以上を適用するために構成されている。更に、（１つ以上の）本方法及び本書に記載した方法ステップ（複数）について言及した利点（複数）の何れかは本装置に当て嵌り、本装置について記載した利点（複数）は対応する（１つ以上の）方法及び方法ステップに当て嵌ることは明らかであろう。
【0038】
更なる一側面に応じ、プログラムがコンピュータで実行されるとき、上記及び／又は下記の実施形態（複数）の何れか１つに係るステップ（複数）の何れか１つに応じて何れか１つの方法に係る方法を実行するためのコンピュータで読み取り可能な命令を含むコンピュータプログラム製品が提供される。
【0039】
更なる一側面に応じ、プログラムがコンピュータで実行されるとき、上記及び／又は下記の実施形態（複数）の何れか１つに係るステップ（複数）の何れか１つに応じた方法を実行するためのコンピュータで読み取り可能な命令を含むコンピュータプログラムが提供される。
【0040】
更なる一側面に応じ、上記及び／又は下記の方法の実施形態（複数）の何れか１つに係る１つ以上のステップを実行するためにプログラムされたコンピュータ装置又は他のハードウエア装置が提供される。
【0041】
他の一側面に応じ、上記及び／又は下記の方法の実施形態（複数）の何れか１つに係る１つ以上のステップを実行するための機械可読かつ機械実行可能な形式のプログラムをコードするデータ記憶装置が提供される。
【0042】
添付の図面は、本開示の装置の現在好ましい非限定的例示的実施形態を説明するために使用されている。添付の図面と組み合わせて読めば以下の詳細な説明から、本開示の特徴及び対象の上記の及び他の利点はより明らかになり、視点ないし側面及び実施形態はより良く理解されるであろう。
【図面の簡単な説明】
【0043】
【図1】空間分布を決定する一例示的方法の一実施形態を説明するフローチャート。
【図3】カーテシアンサンプリング法が使用されている、印加パルスシーケンスのためのスパースヘッセ行列の一例の模式図。
【図4】非カーテシアンサンプリング法が使用されている、印加パルスシーケンスのためのスパースヘッセ行列の一例の模式図。
【図5】印加パルスシーケンスの一例の励起回数に対する、フリップ角の一例（上側）及び対応する位相エンコードステップの一例（下側）を説明するグラフ。
【図6】カラーで表したスパースヘッセ行列の一例の模式的マップ。カラーはスパースヘッセ行列の各要素の大きさのｌｏｇ１０を示す。
【図7】グレースケールで表したスパースヘッセ行列の他の一例の模式的マップ。グレーの濃淡はスパースヘッセ行列の各要素の大きさのｌｏｇ１０を示す。
【図8】関連技術の一方法と本特許開示の一例示的方法についての演算時間の対数に対する目的関数の対数のｌｏｇプロット。
【図9】本開示の一例示的方法によって得られた２つの再構成されたインビボ（in-vivo）Ｔ１及びＴ２マップのセット（左側）、関連技術の一方法についての２つの再構成されたインビボＴ１及びＴ２マップのセット（中央）、及び、対比される両方法についての得られた平均Ｔ１及びＴ２値の対比（右側）。
【図12】本特許開示に応じた空間分布の決定を実行するための一例示的装置のブロック図。
【図13】サンプルのタイムドメイン磁気共鳴測定及び分析を実行する方法の一例のフローチャート。
【図14】サンプルのタイムドメイン磁気共鳴測定を実行するための一例示的装置のブロック図。
【図15】本特許開示の方法（複数）の例についての演算時間の対数に対する目的関数の対数のｌｏｇプロット。スパースヘッセ行列の対角バンドの幅は変化される。
【実施例】
【0044】
図１を参照すると、印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する方法１００が示されている。該方法１００は、以下のステップを含む：
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定すること１１０、但し、該ＴＤＭＲ信号モデルは、
サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含むＴＤＭＲ信号モデルパラメータ（複数）に依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータ（複数）の初期セットを推定する（estimating）こと１２０；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータ（複数）の初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定すること１３０；
iv) ＴＤＭＲ信号モデルパラメータ（複数）のインプットセットに基づいてＴＤＭＲ信号モデルパラメータ（複数）の更新セットを決定すること１４０；ＴＤＭＲ信号モデルパラメータ（複数）のインプットセットは、最初の（１回目の）繰り返し（反復）ではＴＤＭＲ信号モデルパラメータ（複数）の初期セットであり、その後は、ＴＤＭＲ信号モデルパラメータ（複数）の更新セットである；従って、信号モデルパラメータ（複数）はフィッティング結果を得るために繰り返し更新される；該決定すること１４０は
v) 残差ベクトルの勾配を計算すること１４２；
vi) スパースヘッセとしてのＴＤＭＲ信号モデルのヘッセを計算すること１４４；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新し、それによってＴＤＭＲ信号モデルパラメータ（複数）の更新セットを取得すること１４６
を含む；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さいことが決定（確認）される１５２まで又は所定回数の繰り返し（反復）が完了するまで、ＴＤＭＲ信号モデルパラメータ（複数）の更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返す（反復する）こと１５０、それによって、ＴＤＭＲ信号モデルパラメータ（複数）の最終更新セットを取得すること；一般的に、任意の停止基準が使用され得る、従って、該繰り返すこと１５０は（１つの）停止基準が満たされるまで行われ得る。
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定すること１１０、但し、該ＴＤＭＲ信号モデルは、
サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含むＴＤＭＲ信号モデルパラメータ（複数）に依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータ（複数）の初期セットを推定する（estimating）こと１２０；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータ（複数）の初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定すること１３０；
iv) ＴＤＭＲ信号モデルパラメータ（複数）のインプットセットに基づいてＴＤＭＲ信号モデルパラメータ（複数）の更新セットを決定すること１４０；ＴＤＭＲ信号モデルパラメータ（複数）のインプットセットは、最初の（１回目の）繰り返し（反復）ではＴＤＭＲ信号モデルパラメータ（複数）の初期セットであり、その後は、ＴＤＭＲ信号モデルパラメータ（複数）の更新セットである；従って、信号モデルパラメータ（複数）はフィッティング結果を得るために繰り返し更新される；該決定すること１４０は
v) 残差ベクトルの勾配を計算すること１４２；
vi) スパースヘッセとしてのＴＤＭＲ信号モデルのヘッセを計算すること１４４；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新し、それによってＴＤＭＲ信号モデルパラメータ（複数）の更新セットを取得すること１４６
を含む；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さいことが決定（確認）される１５２まで又は所定回数の繰り返し（反復）が完了するまで、ＴＤＭＲ信号モデルパラメータ（複数）の更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返す（反復する）こと１５０、それによって、ＴＤＭＲ信号モデルパラメータ（複数）の最終更新セットを取得すること；一般的に、任意の停止基準が使用され得る、従って、該繰り返すこと１５０は（１つの）停止基準が満たされるまで行われ得る。
【0045】
その後、ＴＤＭＲ信号モデルパラメータ（複数）の最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出するステップix) が実行される。
【0046】
抽出された空間分布は、ＭＲ分布マップ、換言すれば、サンプルのスライス等のようなＭＲ画像の取得を可能にする。
【0047】
ステップiii)、iv)、v)、vi) 及びvii) を繰り返す場合、ＴＤＭＲ信号モデルパラメータの初期セットの代わりに、ＴＤＭＲ信号モデルパラメータの更新セットが使用される。この方法で、ＴＤＭＲ信号モデルパラメータの更新セットはＴＤＭＲ信号モデルパラメータの最終更新セットへと反復的に収斂（収束）する。
【0048】
スパースヘッセとしてＴＤＭＲ信号モデルのヘッセを計算する場合、スパースヘッセは、印加パルスシーケンスに基づくスパース性パターンを有し得る。
【0049】
一実施形態では、ヘッセのスパース性パターンは、勾配エンコードパターン及び／又は無線周波数励起パターンに基づいて決定される。
【0050】
一実施形態では、スパースヘッセ
の計算１４４は、残差ベクトルのヤコビアン（関数行列式）と残差ベクトルのヤコビアンのエルミート転置の積としてヘッセを計算することによって行われる。残差ベクトルのヤコビアンの代わりに、ＴＤＭＲ信号モデルのヤコビアン及び対応するエルミート転置を採用することは等価であろう。
【0051】
一実施形態では、残差ベクトルの勾配の計算１４２は、残差ベクトルと残差ベクトルのヤコビアンのエルミート転置の乗算を含む。残差ベクトルのヤコビアンのエルミート転置の代わりに、ＴＤＭＲ信号モデルのエルミート転置を採用することは等価であろう。
【0052】
一実施形態では、ヘッセの計算１４４は、複数の対角（ダイアゴナル）バンドを含むスパースヘッセとしてヘッセを計算することによって行われる。スパース性パターンは、複数の対角バンドを含み得る。換言すれば、ヘッセの無視不可能（無視できない：non-negligible）要素の対角バンドは、スパース性パターンの１つのタイプを表す。他の要素はゼロか、そうでなければ無視可能である。ある要素は、例えば、その（絶対）値がある閾値より小さい場合、無視可能であると見なされ得る。閾値は予め決定され得る。
【0053】
一実施形態では、ヘッセの複数の対角バンド（対角バンドの個数）は印加パルスシーケンスの勾配エンコードパターンによって決定される。第１の例として、カーテシアン収集はモデルパラメータ当たり１つの対角バンドの好適なスパース性パターンを提供することが決定される。第２の例として、ラジアル及びスパイラル収集はモデルパラメータ当たり４つの対角バンドの好適なスパース性パターンを提供することが決定される。
【0054】
一実施形態では、複数の対角バンドの何れか１つの幅は、印加パルスシーケンスの無線周波数励起パターンによって決定される。
【0055】
一例として、図２に示されているように、スパースヘッセの計算１４４を含むＴＤＭＲ信号モデルパラメータの更新セットの決定１４０には、以下のことが含まれ得る：無視不可能であるヘッセ行列の要素／係数の決定２１０。無視不可能とは、例えば残差ベクトルのヤコビアンと残差ベクトルのヤコビアンのエルミート転置の積の各要素の絶対値が予め決定される値、例えば０．１，０．２，０．３，０．４，０．５，...，１よりも大きいことを意味し得る。予め決定される値は、パラメータマップの十分に正確な再構成を可能にするよう選択され得る。決定２１０のステップは任意的である。なぜなら、これは（１つの）印加パルスシーケンスについて予め決定され得るからである。従って、無視不可能要素のみが評価２２０されかつ記憶２３０される。結果として得られる行列は、真のヘッセ行列に対して、スパース近似ヘッセ（sparse approximation Hessian）
（ＳＡＨ）である。この近似は迅速に演算されかつ記憶されることができる。なぜなら、本方法及び印加パルスシーケンスについては係数（複数）の一部分のみが無視不可能であるからである。スパース性パターン（例えば無視不可能ＳＡＨ係数の指標）と収集法の間には明確な関係があるので、演算される必要がある係数は予め分かる。他の（無視可能な）係数はすべて演算されない。この結果、演算と記憶は劇的に削減される。
【0056】
勾配エンコードパターンはＳＡＨ行列のスパース性パターンを決定する。カーテシアン収集については、位相エンコード方向においてのみ広がる点広がり関数である。これの結果、ＳＡＨ行列にバンド状の対角（ダイアゴナル）ブロック（複数）が生じる。図３は、カーテシアン収集についてのＳＡＨ行列３００の決定されたスパース性パターンの一例を示す。白色の要素は、当該要素が無視可能であることを意味する。従って、これらの要素は、ＳＡＨ行列が使用される演算ステップ（複数）の何れにおいても評価されないであろう（されてはならない）。黒色の要素は、当該要素が無視不可能であることを意味する。軸３０２及び３０４は夫々ＳＡＨ行列３００の行及び列を表す。軸３０２及び３０４は［ボクセル当たりパラメータ］×［ボクセル］と見なすことができる。この例では、破線で区分けされた３×３のブロック３１０がある。ＳＡＨ行列３００は無視不可能要素の対角バンド（複数）３２０を有する。ブロックの個数は一般的にＮ×Ｎであり、Ｎは、再構成されるべきパラメータタイプの個数であり、例えば、Ｔ１、Ｔ２及びプロトン密度が再構成される場合、Ｎ＝３である。
【0057】
図４は、図３の場合と同様の方法で、スパイラル及び／又はラジアル収集のような非カーテシアン収集についてのＳＡＨ行列４００について決定されたスパース性パターンの一例を示す。この場合、点広がり関数は両方向に広がり、従って、その結果として、各ブロックに複数のバンド４２０，４２２及び４２４が生じる。軸４０２及び４０４は夫々、軸３０２及び３０４と同様に、ＳＡＨ行列４００の行及び列を表す。破線によって区分けされた３×３のブロック４１０がある。ＳＡＨ行列４００は、無視不可能要素のメイン対角バンド４２０と、対称的サブバンド（複数）４２２及び４２４を有する。個別にマークされた要素によって形成されるラインの太さないし幅は、当該ＳＡＨ行列要素の絶対値の値を示す。
【0058】
非カーテシアンシーケンスはデータ収集の観点からより効率的であるのに対し、カーテシアンシーケンスはハードウエアの不完全性（例えば渦電流）に対し遥かにより安定的（影響を受けにくい：robust）であるという利点を有する。これはとりわけ、ＭＲ－ＳＴＡＴにおける使用に好適な勾配バランス化（安定化）シーケンス（gradient-balanced sequence）についての場合である。
【0059】
無線周波数励起（フリップ角列ないし曲線（train））とスパース性パターンと滑らかさとの間の関係の例
【0060】
フリップ角（ＲＦ）励起の滑らかさは、例えば上記において図４について説明したような、スパース性パターンにおけるバンド（複数）の幅を決定する。滑らかなフリップ角列については、細いバンドのみを考慮すれば十分であろう。これが意味するところは、ＳＡＨ行列の行当たり幾つかの要素／係数のみが評価されかつ記憶される必要があるということである。上記の図におけるバンドは、従って、フリップ角列が滑らかな場合に、細いことになる。滑らかさの一例は、励起当たりのフリップ角の変化が５度以下であることである。迅速に変化するフリップ角列の場合、反対のことが成り立つ。即ち、より幅の広いバンド（複数）が必要になる。これが意味するところは、ＳＡＨ行列の行当たりより多くの係数が評価されかつ記憶される必要があるということである。従って、例えば、図３及び図４のバンド３２０，４２０，４２２及び４２４は幅がより広くなるであろう。
【0061】
フリップ角列の要求される滑らかさはサンプル、測定を行うために使用される装置、その設定等に依存するが、Ｔ１、Ｔ２及びプロトン密度のようなパラメータを正確に定量できるようにフリップ角にある程度の変化を有することが好ましい。フリップ角列の一例は図５に示されている。この図５は、マルチパラメトリック定量的シーケンスについてのフリップ角列の好適な程度の滑らかさの一例を示している。とりわけ、図５は、変化するフリップ角と線形（直線的）カーテシアンサンプリング法による遷移状態２Ｄバランス化勾配－エコーシーケンスを示す。
【0062】
ここで図１の方法１００に戻り、更には図１０を参照すると、好ましい一実施形態では、更新１４６のステップvii) は、スパースヘッセと残差ベクトルの勾配に基づき方程式（複数）の線形システムを解く５１０ことによって得られるＴＤＭＲ信号モデルパラメータ変化（複数）のセットの取得５２０を含む。
【0063】
ここで図１１を参照すると、より好ましい一実施形態では、更新１４６のステップvii) は、スパースヘッセにＴＤＭＲ信号モデルパラメータ変化を乗じたものが残差ベクトルの勾配をマイナスにしたものに等しいという（スパースヘッセ×ＴＤＭＲ信号モデルパラメータ変化＝残差ベクトルの勾配のマイナスとしての）方程式（複数）の線形システムを解く６３０ことによってＴＤＭＲ信号モデルパラメータ変化のセットを取得６２０すること、及び、取得したＴＤＭＲ信号モデルパラメータ変化を、ＴＤＭＲ信号モデルパラメータの更新セットを取得するために更新されるＴＤＭＲ信号モデルパラメータのセットに加える（加算する）６４０こと、を含む。
【0064】
図１２を参照すると、印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づき該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定するための装置の一実施形態である装置（空間分布決定装置）７００は、上記の方法の何れか１つ以上を実行するよう構成されたプロセッサ７１０を含む。装置７００は、モデル、パラメータ及び／又は方法ステップの実行のために必要な他のデータの何れかを記憶するための記憶媒体７２０を含み得る。記憶媒体７２０は、プロセッサ７１０によって実行される場合、上記のような方法ステップの１つ以上を実行する実行可能コードも記憶し得る。装置７００は、ユーザのインプットを受け入れるための、ネットワークインターフェース７３０及び／又はインプット／アウトプットインターフェース７４０も含み得る。単一の装置として図示されているが、装置７００は、並列コンピュータシステム又はスーパーコンピュータ等のような複数の装置のネットワークとしても構成され得る。
【0065】
図１３を参照すると、サンプルの励起１３１０、信号の受信１３２０及びタイムドメインでの該信号の分析１３３０を含む方法（タイムドメイン分析方法）１３００が示されている。分析１３３０は前記方法１００を含み得る。
【0066】
方法ステップ１３１０及び１３２０は、図１４に示したようなＭＲＩ装置１４００によって実行され得る。装置１４００は、勾配コイル１４０２と、サンプル１４０１を一時的に（過渡的に）励起し、以って該サンプル１４０１に放出信号を放出させる励起装置１４００と、該放出信号を受信するための受信コイル１４０３を含む。ＭＲＩ装置１４００はサンプル１４０１から信号を取得し得る。
【0067】
励起装置１４００及び受信コイル１４０３はＭＲＩの技術分野において周知であり、従って、ここではより詳細には説明しない。装置１４００は、更に、励起装置１４００及び受信コイル１４０３を含むシステムのコンポーネントを制御するプロセッサ１４０５を含む。プロセッサ１４０５は、更に、画像及び／又はステータス情報を出力するためのディスプレイ１４０４と、ユーザからの命令及び付加情報を受信するための、タッチスクリーン、マウス、キーボードのような入力装置１４０６を制御し得る。プロセッサ１４０５は、演算のようなタスクを連携して実行するよう構成されている複数の処理装置を含み得る。代替的に、プロセッサ１４０５は単一の処理装置からなる。中央処理装置（ＣＰＵ）、コントローラ又はＦＰＧＡのような、そのような処理装置は当業界において既知である。説明が不明瞭にならないようにするために、ＭＲＩ装置の幾つかの周知の要素については図示及び説明を省略した。
【0068】
装置１４００は、更に、メモリ１４０７のような記憶手段を含む。メモリ１４０７は、プロセッサ１４０５の制御下で、受信コイル１４０３から受信される信号１４１２を記憶するよう構成され得る。再構成画像データ１４１４は、装置１４００それ自体によって生成され得るが、装置７００の特徴の全てのような、本書で説明されるような方法を実行するために必要な特徴をすべて含む。更に、画像データは装置７００のような外部装置によって再構成されることも可能である。再構成は、本書に記載した方法、例えば方法１００によって受信信号１４１２を処理することによって行われる。メモリ１４０７は、更に、プロセッサ１４０５にそのタスクを実行させるためのコンピュータコードを記憶し得る。例えば、コンピュータコードは、図１について説明したような収集スキームに基づいてデータ収集を実行するための画像収集モジュール１４０８を含み得る。該画像収集モジュール１４０８は、励起装置１４００にサンプル１４０１を一時的に（過渡的に）励起させ、それによって、サンプルに放出信号を放出させ得る。更に、画像収集モジュール１４０８は、プロセッサ１４０５に、サンプル１４０１によって放出される信号を受信する受信コイル１４０３からの放出信号を受信かつ記憶させ得る。
【0069】
具体化の例
【0070】
タイムドメインＭＲスピントモグラフィ（MR Spin Tomography in Time-Domain：ＭＲ－ＳＴＡＴ）は、単回（シングル）ショートスキャンからのデータを用いてマルチパラメトリック定量的ＭＲマップを取得するためのフレームワークである。ブロッホ（Bloch）型体積信号モデルをタイムドメインデータに直接的にフィッティングさせることによって信号の空間的局在化と組織パラメータの推定が同時に実行される大規模最適化問題は解決される。関連技術の方法では、高解像度スキャンについての大規模最適化問題を解決するために使用可能な高度に並列化されたマトリックスフリー（matrix-free）の（不正確）ガウス・ニュートン再構成アルゴリズムが使用される。
【0071】
理論
【0072】
このセクションでは、関連技術の方法に関するＭＲ－ＳＴＡＴフレームワークを示し、再構成プロセスにおける演算上の問題の幾つかが従来どのように取り扱われたかについて概要を示す。その後、用いられた印加パルスシーケンスの一例のために、ヘッセ行列はＭＲ－ＳＴＡＴ再構成を加速するために本方法において利用されるスパース構造を有することを理論的に導き出す。
【0073】
スカラ量（実数及び複素数の両方）は小文字で表され、ベクトル量はボールドの小文字で表され、行列はボールドの大文字で表される。
【0074】
ＭＲ－ＳＴＡＴフレームワーク
【0075】
空間座標
による単一の（１つの）磁化ベクトル（single spin isochromat）
の時間展開（変化）はブロッホ方程式で記述可能である。この時間展開は、使用されるパルスシーケンス（例えば磁化ベクトルに作用するＲＦ励起パルス及び空間勾配）に依存し、そのＭＲ関連生物物理学的組織特性
にも依存する。パルスシーケンスに関連する勾配軌道を
で表すものとする。
【0076】
【0077】
フルサンプリング定常状態シーケンスにおいては、横磁化はその時間依存性を喪失し、ＦＦＴは定性的画像を回復するために使用されることができる。過渡状態シーケンスのより一般的な場合では、ＦＦＴは、測定されたタイムドメイン（又はｋ空間）データと画像空間の間での変換のために最早直接使用されることはできない。従って、過渡状態の場合、以下が実行される。まず、視野
が
ボクセルに離散化（量子化）され、方程式（１）は
となる。
となる。
【0078】
ここで、
はボクセル
における磁化である。勾配バランス化（安定化）シーケンスについて、
は個々の磁化ベクトル（複数）についてブロッホ方程式の数値積分によって演算されることができる。
【0079】
【0080】
ここで、
は成分ごとのベクトル乗算を意味する。信号ベクトル
は上記のＴＤＭＲ信号モデルの一例である。ここで、
は（例えばＴ１、Ｔ２、及びプロトン密度の実数及び虚数部分を含む）ボクセル当たりの別異（distinct）パラメータの数を意味するものとする。
【0081】
従って、
は
の別異（different）パラメータに依存する。パラメータ
がボクセル
に関連するパラメータとなるように、全てのパラメータが結合されて（concatenated）、単一ベクトル
が形成される。
【0082】
上述したようなサンプルから放出されるＴＤＭＲ信号の一例である、測定されたタイムドメインサンプル（複数）のベクトル
を与えられると、残差ベクトル
は
として定義される。
【0083】
残差ベクトル
は方法１００について上述した残差ベクトルの一例である。目的関数
は
として定義可能である。
【0084】
目的関数ｆは、上述したような放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を絶対値の項で（in absolute terms）記述するために使用されることができる。上述したような組織パラメータについての空間分布の一例である、パラメータマップ
は、パラメータについてのブロッホ方程式及び取り得るインターバルによって表される物理的制約が課せられる
を数値（計算）的に解くことによって得られる。
【0085】
非厳密（inexact）ガウス・ニュートン法
【0086】
方程式（８）は非線形最適化問題であることに留意すべきである。そのような問題はニュートンベースの方法（ニュートン法をベースとした方法）を用いて解くことができる。ニュートンベースの方法は、初期推定（initial guess）
で開始し、次いで、線形システム
を解くことによって、パラメータ空間における更新ステップ
を得る。
【0087】
ここで、
は目的関数の勾配であり、
は
として定義されるヘッセ行列である。
【0088】
ＭＲ－ＳＴＡＴへのニュートン法の直接適用に伴う困難性は問題の本来的な大規模性である。２Ｄ問題でさえも、パラメータの数
は１０６のオーダーであり得る。従って、ヘッセ行列又はその逆（行列）を明確に形成することは今日のコンピュータアーキテクチャにおいては実現不可能である。（逆）ヘッセ行列それ自体というよりも寧ろ（逆）ヘッセ行列への近似を用いるニュートン法への多くの適合化が既に開発されている。大規模最適化問題に適するそのような方法の１つは、記憶制限ＢＦＧＳ法（“Ｌ－ＢＦＧＳ”）である。この方法は逆ヘッセ行列に直接的に近似するが、従前の繰り返しからの限定的な数の勾配ベクトルの記憶を必要とするだけである。代替として、最小二乗問題における慣用技術は
としてヘッセ行列に近似することである。ここで、
は
として定義されるヤコビ行列である。
【0089】
【0090】
実用上、マトリックスフリーガウス・ニュートン法は、パラメータ空間
において良好な更新ステップを生成することが観察されるが、この方法は、
の列が行列・ベクトル乗算の各々について再演算されることを必要としない。即ち、インナーループの各繰り返しにおいて、
パラメータの各々について信号の偏導関数が演算される。これらの演算は再構成アルゴリズムの演算のボトルネックをなす。
【0091】
スパースヘッセ近似
【0092】
選択された印加パルスシーケンスに対応するＴＤＭＲ信号について、
は、例えばカーテシアン収集を用いる場合にバンド化される、スパース構造を有する。そして、
は予め（アウター（ループ）繰り返し（反復）毎に１回）演算されかつ記憶されることができる。従って、スパース
によるその後の乗算の演算努力（計算量）は無視可能であり、方程式（９）の線形システムは、
を得るための繰り返し（反復）アルゴリズムによって高速に解かれることができ、このため、より高速なＭＲ－ＳＴＡＴ再構成が得られる。
【0093】
【0094】
過渡状態シーケンスの場合、収集されたデータは、上部に（on top）摂動を有する定常状態シーケンスからのデータとして解釈されることができる。そのようなデータにＦＦＴが適用されると、摂動は、画像空間において点広がり関数との畳み込み（摂動のＦＦＴ）をもたらすことになる。即ち、１つのボクセルに由来する信号は他のボクセル内へ漏れ入り得る。しかしながら、カーテシアンサンプリング法の場合、収集に関連する点広がり関数は（良好な近似で）位相エンコード方向へと空間的に制限される。その上部に滑らかに変化するフリップ角が採用されると（これにより滑らかな信号応答がもたらされる）、点広がり関数は幅も制限される。従って、行列
はよりスパースになる。
【0095】
スパース性パターンの定式的導出（formal derivation）の一例は以下の通りである。まず、ボクセル
における磁化は
（これはＲＦ励起パルスとボクセル
のパラメータ、例えば
に依存する）と、空間エンコード勾配からの位相項を含むベクトル
との成分毎の積として表現されることができる。
【0096】
１．表記の単純化のために、各ボクセルは、ボクセル
とその対応パラメータ
を表すために同じ添え字を使用できるよう、１つの関連パラメータのみを有するものとする。
は組織パラメータから独立している（に依存しない）ので、
となる。
【0097】
一定の（regular）完全にサンプリングされた定常状態シークエンスが使用される場合、
は実際上時間から独立しており（に依存せず）、各要素は
で表されるべき同じ値を有する。
【0098】
従って、
は対角行列である。
【0099】
２．ボクセル
当たり複数のパラメータを有する他の一例については、行列
は
ブロック（複数）の対角行列（複数）から構成されることになる。
【0100】
３．ここで、再び単純化のために
である更なる一例については、過渡状態シーケンスが代わりに使用される。内積は、測定したデータセット（ＴＤＭＲ信号）についての内積の和として以下のように記述することができる：
【0101】
方程式（１４）の項
について、リードアウト（readouts）がカーテシアンであれば、
となる。
となる。
【0102】
即ち、異なるｘ座標を有するボクセル（複数）に関連するパラメータ間の依存性（依存関係）は排除される。
【0103】
４．
とすれば、方程式（１４）の角括弧内の他の項はパーセバル（Parseval）の定理を用いて内積として以下のように記述することができる：
【0104】
印加パルス列のフリップ角列が滑らかであれば、その一例は図５に示されているが、信号応答（サンプルから放出されるＴＤＭＲ信号）は、更にはパラメータに関する信号応答の偏導関数もそうなる。従って、フーリエ変換は制限された（有限）帯域幅を有する。従って、ｙ方向において十分に大きく離れたボクセル（複数）については、内積はゼロになる。そのため、位相エンコード方向において同じライン（同じｘ座標）上にありかつ互いに近接するボクセル（複数）に関連するパラメータ（複数）についての相関（関係）があるのみである。これは、行列
の大きなスパース性をもたらす。
【0105】
５．
の場合について、座標（複数）が
として配列（順序付け）される（ordered）とすれば、行列
は帯行列（ないし多重対角行列：banded matrix）となる。
の場合、
ブロックの帯行列（そのうちの
のみが対称性のために一意的ないし固有（値）（unique）である）から構成されることになる。
【0106】
図６には、小さな数値ファントムについての
の大きさ（magnitude）の対数プロットが示されている（このファントムについて完全な
が演算されることができる）。図６に見出すことができるように、図示の
はブロック構造（６×６ブロック、別異パラメータの各ペア毎に１つ）を有し、各ブロックは、概ね、スパース帯行列である。図７には、より多くのパラメータ（従ってより多くのブロックもある）についての
の大きさの類似の対数プロットが図６と同じ縮尺（スケール）で、但しグレースケールで示されている。
【0107】
疑似アルゴリズム２には、上記の方法１００の一例であるスパースガウス・ニュートン法の一例が示されている。
【0108】
この例では、この方法はインナーループにおける各行列・ベクトル積について
の列（複数）を再演算する必要はないため、ステップ４において、必要となる演算時間は遥かに少ないことが分かる。即ち、インナーループの各繰り返しは線形システム
を解くことを含む。
【0109】
ここで、
は残差ベクトルの一例であり、
は放出されたＴＤＭＲ信号の一例であり、
はＴＤＭＲ信号モデルの一例であり、
は残差関数の計算された勾配の一例であり、
はスパースヘッセの一例であり、
は残差ベクトル（又はＴＤＭＲ信号モデル）のヤコビアンのエルミート転置の一例であり、
は残差ベクトル（又はＴＤＭＲ信号モデル）のヤコビアンの一例であり、
はＴＤＭＲ信号モデルパラメータ変化のセットの一例であり、
はＴＤＭＲ信号モデルパラメータの（更新又は初期）セットの一例である。
【0110】
疑似アルゴリズム２において、ステップ１は方法１００における残差ベクトルの決定１３０のステップiii) の一例である。疑似アルゴリズム２のステップ２は、方法１００の残差ベクトルの勾配の計算１４２のステップv) の一例である。疑似アルゴリズム２のステップ３は、方法１００のスパースヘッセとしてＴＤＭＲ信号モデルのヘッセの計算１４４のステップvi) の一例である。疑似アルゴリズム２のステップ４及び５は方法１００のＴＤＭＲ信号モデルパラメータのインプット（入力）セットの更新１４６の一例である。とりわけ、疑似アルゴリズム２のステップ４は、図１１におけるスパースヘッセにＴＤＭＲ信号モデルパラメータ変化を乗じたものが残差関数の勾配をマイナスにしたもの等しいという（スパースヘッセ×ＴＤＭＲ信号モデルパラメータ変化＝残差関数の勾配のマイナスからなる）方程式の線形システムを解くこと６３０によってＴＤＭＲ信号モデルパラメータ変化（ここでは
）のセットを取得すること６２０の一例であり、疑似アルゴリズム２のステップ５は、取得したＴＤＭＲ信号モデルパラメータ変化を、図１１のＴＤＭＲ信号モデルパラメータの更新セットを取得するために更新されたものであるＴＤＭＲ信号モデルパラメータのセットに加算すること６４０の一例である。
【0111】
実施例
【0112】
パルスシーケンスの例
【0113】
関連技術のマトリックスフリーのガウス・ニュートン法と本開示のスパースガウス・ニュートン（ＧＮ）法を、合成的に生成されたデータについて、ゲルファントムからのデータについて及びインビボ脳データについて対比した。
【0114】
すべてのケースにおいて、２Ｄ勾配（グラジエント）バランス化（安定化）された遷移状態パルスシーケンスが、図５に示されているような、各ＴＲ（繰り返し時間）を滑らかに変化するフリップ角と線形カーテシアンｋ空間充填と共に、用いられた。合成実験のために、シーケンスは、７．８秒の（シミュレート）スキャン時間でデータをシミュレートするために用いられた。ゲルファントム及びインビボ脳実験のために、過渡状態シーケンスは、1.5 Philips-Ingenia ＭＲシステムで実行され、健康なボランティアを１３チャンネル受信ヘッドコイルでスキャンするために用いられた。スキャン時間は、１×１×５ｍｍ３のスキャン解像度で６．８秒であった。
【0115】
シミュレートされた／収集したデータから、
、
、
及び
プロトン密度マップを再構成した。マトリックスフリーＧＮ法については、インナー（ループ）繰り返しの回数は、時間上の制約のために、１０回に限定した。スパースＧＮ法のために、再構成を行ったが、この場合、各パラメータペアについて１５のサブダイアゴナル（subdiagonals）を演算した。その結果、インビボ脳データセットについての完全な行列サイズと比べて必要な記憶（容量）は０．０４％であった。再構成アルゴリズムはオープンソースのJuliaプログラミング言語で記述されており、６４コアを用いるコンピュータクラスタで実行された。
プロトン密度マップを再構成した。マトリックスフリーＧＮ法については、インナー（ループ）繰り返しの回数は、時間上の制約のために、１０回に限定した。スパースＧＮ法のために、再構成を行ったが、この場合、各パラメータペアについて１５のサブダイアゴナル（subdiagonals）を演算した。その結果、インビボ脳データセットについての完全な行列サイズと比べて必要な記憶（容量）は０．０４％であった。再構成アルゴリズムはオープンソースのJuliaプログラミング言語で記述されており、６４コアを用いるコンピュータクラスタで実行された。
【0116】
結果
【0117】
図８に、マトリックスフリーＧＮ法とスパースＧＮ法についてのコンバージェンス（convergence）曲線が示されている。得られた再構成時間はマトリックスフリーＧＮ法と比べて一桁小さいことが分かる。
【0118】
図９は、人間の脳の、スパースＧＮ法についての再構成されたインビボ
及び
マップと、マトリックスフリーＧＮ法との対比を示す。紙面右側に示されているのは、周知のパラメータＧＭ left（脳の左側の灰白質）、ＧＭ right（脳の右側の灰白質）、ＷＭ left（脳の左側の白質）及びＷＭ right（脳の右側の白質）についての平均Ｔ１値（紙面上側）及び平均Ｔ２値（紙面下側）である。スパース（Sparse）ＧＮ法を用いると、再構成時間は５０分間であったのに対し、マトリックスフリー（Matrix-free）ＧＮ法では３２０分間であった。
【0119】
スパースＧＮ法では（それ自体がヘッセ行列に近似する）
の近似が使用されているにも拘らず、最適化アルゴリズムは同じパラメータマップ（複数）に収斂する。実際、該方法は（図８に見られるような）より良好な更新ステップをもたらしさえし得る。なぜなら、インナーループのための繰り返し回数は限定されておらず、そのため、本方法は目的関数の曲率（ないし湾曲：curvature）をより良好に考慮することができるからである。
【0120】
図３及び図４の決定されたスパースヘッセ行列の対角バンドのような、対角バンドの幅が、印加パルスシーケンスに基づいて、又は換言すれば印加パルスシーケンスを用いて、決定されることは、上記から明らかなはずである。その１つの付加的利点は、その幅（又はバンド幅）が１より大きい場合、演算時間は短縮されることである。この効果は図１５に示されているが、１の幅は、５の幅が適用される場合と比べて少なくとも５０倍又は更には１００倍をもより長い演算回数をもたらすことが分かる。図１５の（バンド）幅（ＢＷ）５に対応する図１５のデータはスパースＭＲ－ＳＴＡＴについて図８に示したデータに類似する。図１５に示した演算は、例えば、線形カーテシアンサンプリングと（図５のＲＦ列のような）滑らかなＲＦ列を用いて行い得るが、他のタイプのサンプリングとＲＦ列も使用可能であることは明らかであろう。
【0121】
上記の種々の方法のステップ（複数）はプログラムされたコンピュータによって実行可能であることは当業者であれば容易に認識するであろう。ここで、幾つかの実施形態は、機械又はコンピュータ読み取り可能でありかつ命令（指令）の機械実行可能又はコンピュータ実行可能なプログラムをコードするプログラム記憶装置、例えばデジタルデータ記憶媒体をカバーすることも意図されている。なお、該命令（指令）は上記の方法のステップ（複数）の幾つか又はすべてを実行する。プログラム記憶装置は、例えばデジタルメモリ、磁気ディスク及び磁気テープのような磁気記憶媒体、ハードディスク、又は光学的に読み取り可能なデジタルデータ記憶媒体であり得る。実施形態は、上記の方法（複数）の上記ステップ（複数）を実行するようプログラムされたコンピュータをカバーするようにも意図されている。
【0122】
「ユニット」、「プロセッサ」又は「モジュール」として名称が付された機能ブロックを含む、図示された種々の要素の機能は、専用ハードウエアや適切なソフトウエアに関連してソフトウエアを実行可能なハードウエアの使用によって提供され得る。プロセッサによって提供される場合、その機能は、１つの専用プロセッサによって、１つの共有プロセッサによって、又は複数の個別プロセッサ（これらのうちの幾つかは共有され得る）によって提供され得る。更に、用語「ユニット」、「プロセッサ」又は「コントローラ」の明示的な使用は、ソフトウエアを実行可能なハードウエアのみを指すものと解釈されるべきではなく、デジタル信号プロセッサ（ＤＳＰ）ハードウエア、ネットワークプロセッサ、特定用途向け集積回路（ＡＳＩＣ）、フィールド・プログラマブル・ゲート・アレイ（ＦＰＧＡ）、ソフトウエアを記憶するためのリード・オンリー・メモリ（ＲＯＭ）、ランダム・アクセス・メモリ（ＲＡＭ）及び不揮発性記憶装置を、これらに限定されずに、非明示的に含み得る。他のハードウエアも、従来型及び／又はカスタムである否かに拘わらず、含まれ得る。同様に、図示の何れのスイッチも概念的なものに過ぎない。それらの機能は、プログラムロジックの作動によって、専用ロジックによって、プログラム制御と専用ロジックの相互作用によって、又は手動によってさえも、実行され得るが、文脈からより具体的に理解されるように実行者によって特定の技術を選択することができる。
【0123】
本願に含まれるブロック図は何れも本発明の原理を具現化する説明のための回路の概念図を表すことは当業者であれば分かるはずである。同様に、フローチャート、フローダイアグラム、状態遷移図、疑似コード等は何れも、コンピュータ読み取り可能媒体において実質的に表すことが可能であり、コンピュータ又はプロセッサが明示されているか否かに拘わらず、そのようなコンピュータ又はプロセッサによって実行可能な種々のプロセスを表すことは理解されるはずである。
【0124】
説明した方法及び装置の原理を特定の実施形態（ないし実施例）との関連において上記したが、この説明は専ら例示としてなされたものと理解されるべきであり、添付の（特許）請求の範囲によって決定される保護範囲の限定として理解されるべきではない。
【0125】
本開示は、更に、以下の項目を含む：
１．印加パルスシーケンスに応じたサンプルの励起（励磁）後にサンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいてサンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する方法であって、該方法は、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定すること、但し、ＴＤＭＲ信号モデルは、
サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存（従属）するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定する（estimating）こと；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定すること；
iv) ＴＤＭＲ信号モデルパラメータの更新セットをパラメータのインプット（入力）セットに基づき
v) 残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのヘッセをスパースヘッセとして計算すること、但し、ヘッセのスパース性パターンは印加パルスシーケンスに基づくこと；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定すること；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、ＴＤＭＲ信号モデルパラメータの更新セットをインプットとして用いてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得すること；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出すること
を含む。
１．印加パルスシーケンスに応じたサンプルの励起（励磁）後にサンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいてサンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する方法であって、該方法は、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定すること、但し、ＴＤＭＲ信号モデルは、
サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存（従属）するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定する（estimating）こと；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定すること；
iv) ＴＤＭＲ信号モデルパラメータの更新セットをパラメータのインプット（入力）セットに基づき
v) 残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのヘッセをスパースヘッセとして計算すること、但し、ヘッセのスパース性パターンは印加パルスシーケンスに基づくこと；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定すること；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、ＴＤＭＲ信号モデルパラメータの更新セットをインプットとして用いてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得すること；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出すること
を含む。
【0126】
２．項目１に記載の方法において、印加パルスシーケンスは勾配エンコードパターン及び／又は無線周波数励起パターンを含むこと、及び、ヘッセのスパース性パターンは勾配エンコードパターン及び／又は無線周波数励起パターンに基づいて決定されることを特徴とする。
【0127】
３．項目１又は２に記載の方法において、ヘッセは、残差ベクトルのヤコビアンと残差ベクトルのヤコビアンのエルミート転置の積として計算されることを特徴とする。
【0128】
４．項目１～３の何れかに記載の方法において、印加パルスシーケンスは、カーテシアン（Cartesian）収集（acquisition）、ラジアル収集又はスパイラル収集の何れか１つをもたらすよう構成されていることを特徴とする。
【0129】
５．項目４に記載の方法において、印加パルスシーケンスの勾配エンコードパターンは、対応する点広がり関数が位相エンコード方向において専ら広がるように、カーテシアン収集をもたらすよう構成されていることを特徴とする。
【0130】
６．項目１～５の何れかに記載の方法において、印加パルスシーケンスは、変化するフリップ角をもたらすよう構成されていることを特徴とする。
【0131】
７．項目６に記載の方法において、印加パルスシークエンスの無線周波数励起パターンは、対応する点広がり関数が幅方向において空間的に限定されるように、滑らかに変化するフリップ角をもたらすよう構成されていることを特徴とする。
【0132】
８．項目１～７の何れかに記載の方法において、少なくとも１つの組織パラメータは、Ｔ１緩和時間、Ｔ２緩和時間、Ｔ２＊緩和時間及びプロトン密度の何れか１つ又はこれらの組み合わせを含むことを特徴とする。
【0133】
９．項目１～８の何れかに記載の方法において、残差ベクトルの勾配を計算すること（計算するステップ）は、残差ベクトルと残差ベクトルのヤコビアンのエルミート転置を乗算することを含むことを特徴とする。
【0134】
１０．項目１～９の何れかに記載の方法において、ＴＤＭＲ信号モデルは、ブロッホ（Bloch）型体積信号モデルであることを特徴とする。
【0135】
１１．項目１～１０の何れかに記載の方法において、ヘッセのスパース性パターンは、印加パルスシーケンスの勾配エンコードパターンによって決定されることを特徴とする。
【0136】
１２．項目１～１１の何れかに記載の方法において、ヘッセは、複数の対角（ダイアゴナル）バンドを含むスパースヘッセとして計算されることを特徴とする。
【0137】
１３．項目１２に記載の方法において、ヘッセの複数の対角バンドは、印加パルスシーケンスの勾配エンコードパターンによって決定されることを特徴とする。
【0138】
１４．項目１２又は１３に記載の方法において、複数の対角バンドの何れか１つの幅は、印加パルスシーケンスの無線周波数励起パターンによって決定されることを特徴とする。
【0139】
１５．項目１～１４の何れかに記載の方法において、スパースヘッセは、ＴＤＭＲ信号モデルパラメータの各セットについて一定であることを特徴とする。
【0140】
１６．項目１～１５の何れかに記載の方法において、更新のステップvii) は、
スパースヘッセにＴＤＭＲ信号モデルパラメータ変化を乗じたものが残差関数の勾配をマイナスにしたもの等しいという（スパースヘッセ×ＴＤＭＲ信号モデルパラメータ変化＝残差関数の勾配のマイナスからなる）方程式（複数）の線形システムを解くことによって、ＴＤＭＲ信号モデルパラメータ変化のセットを取得すること、及び、
取得したＴＤＭＲ信号モデルパラメータ変化を、ＴＤＭＲ信号モデルパラメータの更新セットを取得するために更新されるべきＴＤＭＲ信号モデルパラメータのセットに加算すること
を含むことを特徴とする。
スパースヘッセにＴＤＭＲ信号モデルパラメータ変化を乗じたものが残差関数の勾配をマイナスにしたもの等しいという（スパースヘッセ×ＴＤＭＲ信号モデルパラメータ変化＝残差関数の勾配のマイナスからなる）方程式（複数）の線形システムを解くことによって、ＴＤＭＲ信号モデルパラメータ変化のセットを取得すること、及び、
取得したＴＤＭＲ信号モデルパラメータ変化を、ＴＤＭＲ信号モデルパラメータの更新セットを取得するために更新されるべきＴＤＭＲ信号モデルパラメータのセットに加算すること
を含むことを特徴とする。
【0141】
１７．印加パルスシーケンスに応じたサンプルの励起後にサンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいてサンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する装置であって、該装置は、下記ステップを実行するよう構成されたプロセッサを含むこと：
即ち、該プロセッサは、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定する、但し、ＴＤＭＲ信号モデルは、
サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存（従属）するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定する（estimate）；
iii) ＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定する；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのヘッセをスパースヘッセとして計算すること、但し、ヘッセのスパース性パターンは印加パルスシーケンスに基づくこと；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定する；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、ＴＤＭＲ信号モデルパラメータの更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得する；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出する
よう構成される。
即ち、該プロセッサは、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定する、但し、ＴＤＭＲ信号モデルは、
サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存（従属）するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定する（estimate）；
iii) ＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定する；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのヘッセをスパースヘッセとして計算すること、但し、ヘッセのスパース性パターンは印加パルスシーケンスに基づくこと；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定する；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、ＴＤＭＲ信号モデルパラメータの更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得する；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出する
よう構成される。
【0142】
１８．プログラムがコンピュータで実行されるとき、項目１～１６の何れかに記載の方法を行うためのコンピュータで実行可能な命令を含むコンピュータプログラム製品。
【0143】
ここに本発明の可能な態様を付記する。
［付記１］印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する方法。
該方法は、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定すること、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定すること；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定すること；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定すること；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットを用いてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得すること；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出すること
を含む。
［付記２］上記の方法において、
前記印加パルスシーケンスは勾配エンコードパターン及び／又は無線周波数励起パターンを含み、及び、
前記スパースヘッセは複数の対角バンドを含み、該複数の対角バンドの何れか１つの幅は、前記印加パルスシーケンスの前記勾配エンコードパターン及び／又は前記無線周波数励起パターンによって決定される。
［付記３］上記の方法において、
前記スパースヘッセを計算することは、ヘッセ行列の無視不可能要素を決定することを含み、好ましくは該ヘッセ行列の要素は、各要素の絶対値が閾値を上回る場合、無視不可能として決定され、より好ましくは該ヘッセ行列の要素は、前記残差ベクトルのヤコビアンと前記残差ベクトルのヤコビアンのエルミート転置の積の絶対値が所定の値より大きい場合、無視不可能として決定される。
［付記４］上記の方法において、
前記複数の対角バンドの何れか１つの幅は１より大きく、好ましくは前記複数の対角バンドの何れか１つの幅は２～５５の範囲、より好ましくは３～８の範囲にあり、最も好ましくは３である。
［付記５］上記の方法において、
前記ヘッセは、前記残差ベクトルのヤコビアンと前記残差ベクトルのヤコビアンのエルミート転置の積として計算される。
［付記６］上記の方法において、
前記印加パルスシーケンスは、カーテシアン（Cartesian）収集、ラジアル収集又はスパイラル収集の何れか１つをもたらすよう構成されている。
［付記７］上記の方法において、
前記印加パルスシーケンスの勾配エンコードパターンは、対応する点広がり関数が位相エンコード方向において専ら広がるように、カーテシアン収集をもたらすよう構成されている。
［付記８］上記の方法において、
前記印加パルスシーケンスは、変化するフリップ角をもたらすよう構成されている。
［付記９］上記の方法において、
前記印加パルスシークエンスの無線周波数励起パターンは、対応する点広がり関数が幅方向において空間的に限定されるように、滑らかに変化するフリップ角をもたらすよう構成されている。
［付記１０］上記の方法において、
前記少なくとも１つの組織パラメータは、Ｔ１緩和時間、Ｔ２緩和時間、Ｔ２ ＊ 緩和時間及びプロトン密度の何れか１つ又はこれらの組み合わせを含む。
［付記１１］上記の方法において、
前記残差ベクトルの勾配を計算することは、該残差ベクトルと該残差ベクトルのヤコビアンのエルミート転置を乗算することを含む。
［付記１２］上記の方法において、
前記ＴＤＭＲ信号モデルは、ブロッホ型体積信号モデルである。
［付記１３］上記の方法において、
前記ヘッセにおける複数の対角バンドは、前記印加パルスシーケンスの勾配エンコードパターンによって決定される。
［付記１４］上記の方法において、
前記スパースヘッセはＴＤＭＲ信号モデルパラメータの各セットについて一定である。
［付記１５］上記の方法において、
前記更新のステップvii) は、
スパースヘッセにＴＤＭＲ信号モデルパラメータ変化を乗じたものが残差関数の勾配のマイナスに等しいという方程式の線形システムを解くことによって、ＴＤＭＲ信号モデルパラメータ変化のセットを取得すること、及び、
取得したＴＤＭＲ信号モデルパラメータ変化を、ＴＤＭＲ信号モデルパラメータの更新セットを取得するために更新されるべきＴＤＭＲ信号モデルパラメータのセットに加算すること
を含む。
［付記１６］印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する装置。
該装置は、下記ステップを実行するよう構成されたプロセッサを含む；
即ち、該プロセッサは、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定する、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定する；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定する；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定する；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得する；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出する
よう構成される。
［付記１７］プログラムがコンピュータで実行されるとき、上記の方法を行うためのコンピュータで実行可能な命令を含むコンピュータプログラム製品。
［付記１］印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する方法。
該方法は、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定すること、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定すること；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定すること；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定すること；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットを用いてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得すること；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出すること
を含む。
［付記２］上記の方法において、
前記印加パルスシーケンスは勾配エンコードパターン及び／又は無線周波数励起パターンを含み、及び、
前記スパースヘッセは複数の対角バンドを含み、該複数の対角バンドの何れか１つの幅は、前記印加パルスシーケンスの前記勾配エンコードパターン及び／又は前記無線周波数励起パターンによって決定される。
［付記３］上記の方法において、
前記スパースヘッセを計算することは、ヘッセ行列の無視不可能要素を決定することを含み、好ましくは該ヘッセ行列の要素は、各要素の絶対値が閾値を上回る場合、無視不可能として決定され、より好ましくは該ヘッセ行列の要素は、前記残差ベクトルのヤコビアンと前記残差ベクトルのヤコビアンのエルミート転置の積の絶対値が所定の値より大きい場合、無視不可能として決定される。
［付記４］上記の方法において、
前記複数の対角バンドの何れか１つの幅は１より大きく、好ましくは前記複数の対角バンドの何れか１つの幅は２～５５の範囲、より好ましくは３～８の範囲にあり、最も好ましくは３である。
［付記５］上記の方法において、
前記ヘッセは、前記残差ベクトルのヤコビアンと前記残差ベクトルのヤコビアンのエルミート転置の積として計算される。
［付記６］上記の方法において、
前記印加パルスシーケンスは、カーテシアン（Cartesian）収集、ラジアル収集又はスパイラル収集の何れか１つをもたらすよう構成されている。
［付記７］上記の方法において、
前記印加パルスシーケンスの勾配エンコードパターンは、対応する点広がり関数が位相エンコード方向において専ら広がるように、カーテシアン収集をもたらすよう構成されている。
［付記８］上記の方法において、
前記印加パルスシーケンスは、変化するフリップ角をもたらすよう構成されている。
［付記９］上記の方法において、
前記印加パルスシークエンスの無線周波数励起パターンは、対応する点広がり関数が幅方向において空間的に限定されるように、滑らかに変化するフリップ角をもたらすよう構成されている。
［付記１０］上記の方法において、
前記少なくとも１つの組織パラメータは、Ｔ１緩和時間、Ｔ２緩和時間、Ｔ２ ＊ 緩和時間及びプロトン密度の何れか１つ又はこれらの組み合わせを含む。
［付記１１］上記の方法において、
前記残差ベクトルの勾配を計算することは、該残差ベクトルと該残差ベクトルのヤコビアンのエルミート転置を乗算することを含む。
［付記１２］上記の方法において、
前記ＴＤＭＲ信号モデルは、ブロッホ型体積信号モデルである。
［付記１３］上記の方法において、
前記ヘッセにおける複数の対角バンドは、前記印加パルスシーケンスの勾配エンコードパターンによって決定される。
［付記１４］上記の方法において、
前記スパースヘッセはＴＤＭＲ信号モデルパラメータの各セットについて一定である。
［付記１５］上記の方法において、
前記更新のステップvii) は、
スパースヘッセにＴＤＭＲ信号モデルパラメータ変化を乗じたものが残差関数の勾配のマイナスに等しいという方程式の線形システムを解くことによって、ＴＤＭＲ信号モデルパラメータ変化のセットを取得すること、及び、
取得したＴＤＭＲ信号モデルパラメータ変化を、ＴＤＭＲ信号モデルパラメータの更新セットを取得するために更新されるべきＴＤＭＲ信号モデルパラメータのセットに加算すること
を含む。
［付記１６］印加パルスシーケンスに応じたサンプルの励起後に該サンプルから放出されたタイムドメイン磁気共鳴（ＴＤＭＲ）信号に基づいて該サンプルの内部の少なくとも１つの組織パラメータの空間分布を決定する装置。
該装置は、下記ステップを実行するよう構成されたプロセッサを含む；
即ち、該プロセッサは、
i) 放出されたタイムドメイン磁気共鳴信号に近似するＴＤＭＲ信号モデルを決定する、但し、該ＴＤＭＲ信号モデルは、
前記サンプルの内部の少なくとも１つの組織パラメータの空間分布、及び、
印加パルスシーケンス
を含む、ＴＤＭＲ信号モデルパラメータに依存するよう構成されていること；
ii) ＴＤＭＲ信号モデルパラメータの初期セットを推定する；
iii) インプットとしてのＴＤＭＲ信号モデルパラメータの初期セットに基づいて、放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差を表す残差ベクトルを決定する；
iv) ＴＤＭＲ信号モデルパラメータの更新セットを
v) 前記残差ベクトルの勾配を計算すること；
vi) ＴＤＭＲ信号モデルのスパースヘッセを計算すること、但し、計算されたスパースヘッセは印加パルスシーケンスを用いて計算されるスパース性パターンを有する；及び、
vii) 計算された勾配と計算されたヘッセに基づいてＴＤＭＲ信号モデルパラメータの初期セットを更新すること
によって決定する；
viii) 放出されたＴＤＭＲ信号とＴＤＭＲ信号モデルとの間の差が所定の閾値より小さくなるまで又は所定回数の繰り返しが完了するまで、インプットとしてのＴＤＭＲ信号モデルパラメータの更新セットについてステップiii)、iv)、v)、vi) 及びvii) を繰り返し、それによって、ＴＤＭＲ信号モデルパラメータの最終更新セットを取得する；及び、
ix) ＴＤＭＲ信号モデルパラメータの最終更新セットから、少なくとも１つの組織パラメータの空間分布を抽出する
よう構成される。
［付記１７］プログラムがコンピュータで実行されるとき、上記の方法を行うためのコンピュータで実行可能な命令を含むコンピュータプログラム製品。
【図1】
【図2】
【図3】
【図4】
【図5】
【図6】
【図7】
【図8】
【図9】
【図10】
【図11】
【図12】
【図13】
【図14】
【図15】
