(57)【要約】【課題】注目波形成分に対して計算スペクトルを高精度にフィッティングさせる。ベースライン成分の影響を受け難いスペクトルフィッティングを実現する。【解決手段】ＮＭＲスペクトル（観測スペクトル）１４中のベースライン成分１８は時間空間においてスパースな信号として表現される。一方、Ｌｐノルム（但しｐ≦１）の最小化を条件とした非線形フィッティングにおいて、Ｌｐノルムはノルム計算対象のスパース性を強調する作用を発揮する。この２つの性質を利用して、残差スペクトル２４に対応する時間空間内の残差ベクトル３０についてのＬｐノルムが最小化されるように、計算スペクトル２０に与える係数ベクトルの最適解が探索される。【選択図】図１
2019090747.tif 000002
【特許請求の範囲】
【請求項１】
観測スペクトルを受け入れる手段と、
前記観測スペクトルと計算スペクトルとの間の残差スペクトルに対応する所定空間内の残差ベクトルについての評価値が最小化されるように、前記計算スペクトルを定義する関数に与える係数セットを変化させることにより、最適係数セットを探索する探索手段と、
を含み、
前記計算スペクトルは前記観測スペクトル中の注目波形成分にフィッティングされるものであり、
前記所定空間は前記観測スペクトル中のベースライン成分がスパースな信号として表現される空間であり、
前記評価値は前記所定空間内の残差ベクトルのＬｐノルム（但しｐ≦１）である、
ことを特徴とするスペクトル解析装置。
【請求項２】
請求項１記載の装置において、
前記所定空間は時間空間である、
ことを特徴とするスペクトル解析装置。
【請求項３】
請求項１記載の装置において、
前記所定空間内の残差ベクトルは前記残差スペクトルに対して変換行列を適用することにより求められる、
ことを特徴とするスペクトル解析装置。
【請求項４】
請求項３記載の装置において、
前記変換行列は、周波数空間から前記所定空間への変換を行う第１変換行列を含む、
ことを特徴とするスペクトル解析装置。
【請求項５】
請求項４記載の装置において、
前記変換行列は、更に、前記周波数空間内の前記残差スペクトルに対して傾斜補正を行う第２変換行列を含み、
前記傾斜補正前の残差スペクトルに対して前記第２変換行列及び前記第１変換行列がその順で適用される、
ことを特徴とするスペクトル解析装置。
【請求項６】
請求項１記載の装置において、
前記ｐを変更する手段を含む、
ことを特徴とするスペクトル解析装置。
【請求項７】
請求項１記載の装置において、
前記観測スペクトルはＮＭＲスペクトルである、
ことを特徴とするスペクトル解析装置。
【請求項８】
ＮＭＲ測定で得られたＮＭＲスペクトルを受け入れる工程と、
前記ＮＭＲスペクトルと計算スペクトルとの間の残差スペクトルに対応する所定空間内の残差ベクトルについての評価値が最小化されるように、前記計算スペクトルを定義する関数に与える係数セットを変化させることにより、最適係数セットを探索する工程と、
を含み、
前記計算スペクトルは前記ＮＭＲスペクトル中の注目波形成分にフィッティングされるものであり、
前記所定空間は前記ＮＭＲスペクトル中のベースライン成分がスパースな信号として表現される空間であり、
前記評価値は前記所定空間内の残差ベクトルのＬｐノルム（但しｐ≦１）である、
ことを特徴とするスペクトル解析方法。
【請求項９】
情報処理装置においてスペクトル解析方法を実行するプログラムであって、
前記スペクトル解析方法は、
ＮＭＲ測定で得られたＮＭＲスペクトルを受け入れる工程と、
前記ＮＭＲスペクトルと計算スペクトルとの間の残差スペクトルに対応する所定空間内の残差ベクトルについての評価値が最小化されるように、前記計算スペクトルを定義する関数に与える係数セットを変化させることにより、最適係数セットを探索する工程と、
を含み、
前記計算スペクトルは前記ＮＭＲスペクトル中の注目波形成分にフィッティングされるものであり、
前記所定空間は前記ＮＭＲスペクトル中のベースライン成分がスパースな信号として表現される空間であり、
前記評価値は前記所定空間内の残差ベクトルのＬｐノルム（但しｐ≦１）である、
ことを特徴とするプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、スペクトル解析装置及び方法に関し、特に、スペクトルフィッティング技術に関する。
【背景技術】
【０００２】
スペクトル解析の対象となるスペクトルとして、ＮＭＲ（Nuclear Magnetic Resonance）スペクトル、Ｘ線スペクトル、分光スペクトル、マススペクトル等が挙げられる。スペクトル解析においては、例えば、観測スペクトル（Observed Spectrum）に計算スペクトル（Computed Spectrum, Simulated Spectrum）が最もフィッティングするように、計算スペクトルを定義する関数に与える最適な係数セットが探索される。最適な係数セットを求めるための手法として、各種の手法が提案されている。
【０００３】
例えば、スペクトルフィッティングに際しては、ＬＭ（Levenberg-Marquadts）法が利用される。ＬＭ法は、非線形最小二乗フィッティング法の一種であり、以下の評価値Ｊを最小にする係数セット（係数ベクトル）ａを求めるものである（非特許文献１を参照）。
【数１】
2019090747.tif 000003
【０００４】
上記において、ｙは観測スペクトル（ベクトル）であり、ａは係数セットとしての係数ベクトルであり、ｙ（ａ）は係数セットａをモデル関数に代入することよって生成される計算スペクトル（ベクトル）である。（１）式中の（ｙ−ｙｍ（ａ））は残差スペクトル（ベクトル）である。（１）式の右辺にある１／２は無くてもよく、評価値Ｊの実体はＬ２ノルムである。
【０００５】
なお、一般に、Ｌｐノルムに関しては、以下の（２）式及び（３）式に示す２つの表現が認められるが、本願において“Ｌｐノルム”は基本的に（３）式で表現されるノルムである。例えば、上記Ｌ２ノルムは（３）式中のｐに２を代入したものである。以下の（２）式及び（３）式に含まれる添え字ｎはベクトル要素数を示している。
【数２】
2019090747.tif 000004
【数３】
2019090747.tif 000005
【０００６】
スペクトルフィッティングに際しては、ＴＬＳ（Total Line Shape fitting）法も利用される。ＴＬＳ法は、ＮＭＲスペクトルの解析において一般的に用いられている手法であり、それも非線形最小二乗フィッティング法の一種である（非特許文献２−４を参照）。この方法では、計算スペクトルを定義する関数にベースライン関数が含まれる。
【先行技術文献】
【非特許文献】
【０００７】
【非特許文献１】Jorge J. More, The Levenberg-Marquardt Algorithm: implementation and theory, Conference of numerical analysis (1977).
【非特許文献２】R. Laatikainen, M. Niemitz., W. J. Malaisse, M. Biesemans and R. Willem, A Computational Strategy for the Decovolution of NMR Spectra with Multiplet Structures and Constraits : Analysis of Overlapping 13C-2H multiplets of 13C Enriched Metabolies from Cell Suspensions Incubated in Deuterated Media., Magn. Reson. Med., 36, 359 (1996).
【非特許文献３】P. Soininen, J. Haarala, J. Vepsalainen, M. Niemitz and R. Laatikainen, Strategies for organic impurity quantification by 1H NMR spectroscopy: Constrained total-line-shape fitting, Analytica Chimica Acta., 29, 542 (2005).
【非特許文献４】Juliane Sigl, Nonlinear Residual Minimization by Iteratively Reweighted Least Squares, Comput. Optim. Appl.64, 755 (2015).
【発明の概要】
【発明が解決しようとする課題】
【０００８】
観測スペクトルには、一般に、注目波形成分及びベースライン成分が含まれる。注目波形成分は、例えば、１又は複数のピークを含む部分であり、それが本来的な解析対象である。一方、ベースライン成分は、本来的な解析対象ではない成分であって、周波数空間において広帯域にわたって存在する成分である。例えば、ＮＭＲスペクトルにおいては、測定段階で生じるリンキングノイズ、デジタルフィルタ処理後のデータ欠損、分子構造由来の信号成分等がベースライン成分を生じさせる。ＮＭＲスペクトル中の基底を構成する大きな周期をもった変動、湾曲、傾斜等がベースライン成分であり、あるいは、そのような基底それ自体がベースライン成分である。
【０００９】
上記のベースライン成分を含むスペクトルに対して上記ＬＭ法を適用した場合、フィッティング精度が低下してしまう。ベースライン成分も、注目波形成分と同様に、フィッティング対象になってしまうからである。なお、ＴＬＳ法においては、ベースライン関数に対しても少なからずの係数を与える必要がある。これにより、解析者の負担の増大、解析時間の増大、最適解を求められない可能性が高まる、といった問題が生じる。
【００１０】
本発明の目的は、注目波形成分に対して計算スペクトルを精度良くフィッティングさせることにある。あるいは、本発明の目的は、ベースライン成分の影響を受け難いスペクトルフィッティングを実現することにある。
【課題を解決するための手段】
【００１１】
実施形態に係るスペクトル解析装置は、観測スペクトルを受け入れる手段と、前記観測スペクトルと計算スペクトルとの間の残差スペクトルに対応する所定空間内の残差ベクトルについての評価値が最小化されるように、前記計算スペクトルを定義する関数に与える係数セットを変化させることにより、最適係数セットを探索する探索手段と、を含み、前記計算スペクトルは前記観測スペクトル中の注目波形成分にフィッティングされるものであり、前記所定空間は前記観測スペクトル中のベースライン成分がスパースな信号として表現される空間であり、前記評価値は前記所定空間内の残差ベクトルのＬｐノルム（但しｐ≦１）である、ことを特徴とする。
【００１２】
上記構成において、所定空間は、観測スペクトル中のベースライン成分がスパースな信号（多くのゼロを含んだ信号）として表現される空間である。すなわち、ベースライン成分は、所定空間軸上において、広域的に存在するのではなく、局所的に存在するものである。一方、Ｌｐノルム（但しｐ≦１）を最小化する条件の下での最適解の探索において、当該Ｌｐノルムは、ノルム演算対象のスパース性を増進又は強調する作用を発揮することが知られている。そのようなＬｐノルムを残差ベクトルの評価値とすれば、残差ベクトルのスパース性を高める方向に最適解の探索が進行する。例えば、残差ベクトル全体が均等に抑圧されるのではなく、残差ベクトル中により多くの０又は０近傍値が生じるように、残差ベクトルが抑圧される。結果として、ベースライン成分に対するフィッティングよりも、注目波形成分に対するフィッティングが優先されることになる。一般に、ｐが１以下の場合にその傾向が認められ、ｐを小さくすればするほどその傾向が強まる。よって、上記探索手段によれば、ベースライン成分に対してそれ専用の関数を当てはめることなく、注目波形部分に対する計算スペクトルのフィッティングの精度を高められる。
【００１３】
実施形態において、前記所定空間は時間空間である。例えば、周波数軸上の広帯域信号であるベースライン成分は、時間軸上において、幅の狭いパルス状の部分と、０又はほぼ０からなる平坦な部分と、を含む信号として表現される。すなわち、スパースな信号として表現される。上記構成はその性質を利用するものである。
【００１４】
実施形態において、前記所定空間内の残差ベクトルは前記残差スペクトルに対して変換行列を適用することにより求められる。変換行列は、例えば、周波数空間から時間空間への変換を行う行列である。その例として、ＩＦＦＴ（Inverse Fourier Transform）行列があげられる。
【００１５】
実施形態において、前記変換行列は、周波数空間から前記所定空間への変換を行う第１変換行列と、前記周波数空間内の前記残差スペクトルに対して傾斜補正を行う第２変換行列と、を含み、前記傾斜補正前の残差スペクトルに対して前記第２変換行列及び前記第１変換行列がその順で適用される。観測スペクトルにおいて、ベースライン成分が、大きな周期をもった変動を有しつつ、それ全体として傾斜している場合、傾斜の影響が残差ベクトルに現れ、最適解の探索が旨くいかないおそれが生じる。そこで、残差スペクトルに対して第２変換行列を適用し、傾斜成分を除去し又は減少させた上で、その適用結果に対して第１変換行列が適用される。
【００１６】
実施形態においては、前記ｐを変更する手段が設けられる。上記のように、一般に、ｐを小さくすればするほど、ベースライン成分に対するフィッティングよりも注目波形成分に対するフィッティングが優先される傾向が強まる。ベースライン成分の態様、必要なフィッティング精度、演算時間等に応じて、ｐを可変できるように構成するのが望ましい。
【００１７】
実施形態に係るスペクトル解析方法は、ＮＭＲ測定で得られたＮＭＲスペクトルを受け入れる工程と、前記ＮＭＲスペクトルと計算スペクトルとの間の残差スペクトルに対応する所定空間内の残差ベクトルについての評価値が最小化されるように、前記計算スペクトルを定義する関数に与える係数セットを変化させることにより、最適係数セットを探索する工程と、を含み、前記計算スペクトルは前記ＮＭＲスペクトル中の注目波形成分にフィッティングされるものであり、前記所定空間は前記ＮＭＲスペクトル中のベースライン成分がスパースな信号として表現される空間であり、前記評価値は前記所定空間内の残差ベクトルのＬｐノルム（但しｐ≦１）である。上記のＮＭＲスペクトルを受け入れる工程、及び、最適係数セットを探索する工程は、望ましくは情報処理装置において実行される。
【００１８】
上記方法は、スペクトル解析対象ではないベースライン成分が所定空間内においてスパース性を呈することに着目し、且つ、Ｌｐノルム（但しｐ≦１）がその最小化に際してノルム演算対象のスパース性を高めるように働くことに着目し、最適係数セットの探索において、所定空間内の残差ベクトルのＬｐノルムを評価値として用いるものである。すなわち、上記方法は、数学的な法則それ自体を特徴事項とするものではなく、スペクトルフィッティングにおいてフィッティング不要な部分が自然に除外されるようにする具体的な処理又は演算を特徴事項とするものである。上記スペクトル解析方法を実行するプログラムは、ネットワークを介して、あるいは、可搬型記憶媒体を介して、情報処理装置（スペクトル測定装置を含む）にインストールされ得る。
【発明の効果】
【００１９】
本発明によれば、注目波形成分に対して計算スペクトルを高精度にフィッティングさせることが可能となる。あるいは、本発明によれば、ベースライン成分の影響を受け難いスペクトルフィッティングを実現できる。
【図面の簡単な説明】
【００２０】
【図１】実施形態に係るＮＭＲシステムを示す概念図である。
【図２】実施形態に係るスペクトル解析装置を示すブロック図である。
【図３】第１実施例に係るスペクトル解析方法を示す概念図である。
【図４】第１実施例に係るアルゴリズムを示す概念図である。
【図５】３つのベースライン成分が加えられた３つのスペクトルを示す図である。
【図６】３つのスペクトルに対するフィッティング結果を示す図である。
【図７】第２実施例に係るスペクトル解析方法を示す概念図である。
【図８】第３実施例に係るスペクトル解析方法を示す概念図である。
【図９】第４実施例に係るスペクトル解析方法を示す概念図である。
【発明を実施するための形態】
【００２１】
以下、実施形態を図面に基づいて説明する。
【００２２】
図１には、ＮＭＲ測定装置１０とＮＭＲスペクトル解析装置１２とからなるＮＭＲシステムが示されている。ＮＭＲ測定装置１０は、大別して、測定部と分光計とにより構成される。測定部は、静磁場発生器及びプローブを含むものである。静磁場発生器は垂直貫通孔としてのボアを有し、そのボアの中にプローブにおける挿入部が挿入される。挿入部のヘッドにおいて、試料で生じた磁気共鳴が検出される。分光計は、送信部、受信部、制御部、演算部等からなるものである。送信部において送信パルスシーケンスに従う送信信号が生成され、その送信信号がプローブへ送られる。プローブからの受信信号が受信部に送られる。受信部において受信信号に基づいてＮＭＲスペクトル（以下、観測スペクトルという。）１４が生成される。
【００２３】
観測スペクトル１４は、周波数空間（周波数領域）内に存在し、換言すれば、周波数軸上に波形として表現されるものである。観測スペクトル１４をベクトルとして表現したものがベクトルｙである。ベクトルｙは、周波数軸方向に並ぶ複数の強度値からなる。観測スペクトル１４は、一般に、注目波形成分１６とベースライン成分１８とを含む。
【００２４】
注目波形成分１６は、一般に、１又は複数のピークを含んだ部分である。典型的なスペクトル解析では、ピーク単位で、振幅、周波数、半値幅等のパラメータが推定される。一方、ベースライン成分１８は、フィッティング対象外の部分であり、一般に、それは広帯域成分である。具体的には、ベースライン成分１８は、観測スペクトル１４の基底を構成する比較的大きな周期をもった変動、湾曲、傾斜等に相当する。あるいは、そのような成分をもった基底自体をベースライン成分として理解してもよい。スペクトル解析においては、そのようなベースライン成分１８に影響されずに、注目波形成分１６に対して正確なフィッティングを行うことが求められる。
【００２５】
スペクトル解析装置１２は、ＮＭＲ測定装置１０で生成された観測スペクトル１４を解析する装置であり、それはパーソナルコンピュータ等の情報処理装置により構成される。スペクトル解析装置１２は、評価値Ｊが最小化されるように係数セット（係数ベクトル）ａの内容を段階的に優良化することにより、最適解としての最適係数セットを探索するものである。実施形態においては、最適係数セットがスペクトル解析結果を構成する。
【００２６】
具体的には、係数セットａをスペクトル生成関数（モデル関数）に与えることにより、計算スペクトル２０が生成される。係数セットａの初期値がａ0である。計算スペクトル２０をベクトルとして表現したものがベクトルｙｍ（ａ）である。ベクトルｙｍ（ａ）は周波数軸上に並ぶ複数の強度値からなる。図１においては、計算スペクトル２０を生成する部分に符号２２が付されている。観測スペクトル１４と計算スペクトル２０の差が残差スペクトル２４である。図１においては、残差スペクトル２４を演算する部分に符号２６が付されている。
【００２７】
残差スペクトル２４は周波数空間内に存在するものであり、それを時間空間に表現したものが残差ベクトル（残差信号）３０である。その際には残差スペクトル２４に対して逆フーリエ変換（ＩＦＦＴ）等が適用される。残差ベクトル３０は、時間軸方向に並ぶ複数の振幅値により構成される。図示の例において、残差ベクトル３０には、局所的に生じた短い時間幅をもったピーク部分３０ａと広域に及ぶ平坦な部分３０ｂとが含まれる。ピーク部分３０ａが周波数空間内においてベースライン成分１８の実体をなす。逆に言えば、ベースライン成分１８は、時間空間においてスパースな信号として表現されるものである。平坦な部分３０ｂにはランダムノイズが含まれる。図１においては、残差ベクトル３０を演算する部分に符号２８が付されている。
【００２８】
実施形態においては、評価値Ｊとして、時間軸上の残差ベクトル３０のＬｐノルム（但しｐ≦１）が利用される。残差ベクトル３０のＬｐノルムが最小化されるように、係数ベクトル（係数セット）ａの内容が順次更新され（符号３４参照）、最終的に、最適解としての最適係数セットが特定される。その最適係数セットが外部に出力され（符号３６参照）、あるいは、表示される。その最適係数セットが他の演算に利用されてもよい。図１においては、評価値Ｊを計算し且つ係数ベクトルａを更新する部分に符号３２が付されている。
【００２９】
図１においては、最適係数セットの探索過程が模式的に示されているが、実際には、一連の計算式を繰り返し実行することにより、最適係数セットが特定される。例えば、残差スペクトルの計算（符号２６参照）、残差ベクトルの計算（符号２８参照）、及び、評価値Ｊ（Ｌｐノルム）の計算（符号３２参照）は、単一の計算式により表現され得る。
【００３０】
評価値としてのＬｐノルム（但しｐ≦１）は、その最小化に際して、ノルム計算対象である残差ベクトル３０のスパース性を高めるように働く。図示の例では、ピーク部分３０ａの抑圧よりも、平坦な部分３０ｂの抑圧及び増大が優先されることになる。もっとも、図示された残差ベクトル３０は、最小化条件が満たされた時点でのものである。なお、平坦な部分３０ａにはランダムノイズが含まれ、通常、それは完全に０にはならない。
【００３１】
以上のように、実施形態では、残差スペクトル２４それ自体を評価対象とするのではなく、時間空間内の残差ベクトル３０を評価対象としており、かつ、評価値としてＬｐノルム（但しｐ≦１）を利用している。これにより、最適解の探索に際して、ベースライン成分１８へのフィッティングよりも、注目波形成分１６へのフィッティングの方が優先される。換言すれば、スペクトルフィッティングにおいて、ベースライン成分１８が無視され易くなる。その結果、注目波形成分１６に対して計算スペクトル２０を正しくフィッティングさせることが可能となる。なお、実施形態によれば、ベースライン成分１８を模擬する関数を用いる必要がないので、パラメータ指定に際してのユーザーの負担が大幅に軽減される。
【００３２】
図２には、スペクトル解析装置１２の構成例が示されている。スペクトル解析装置１２は、ＣＰＵ１００、メモリ１０２、入力器１０６、表示器１０８等を含む。メモリ１０２にはスペクトル解析プログラム１０４が格納されている。スペクトル解析プログラム１０４がＣＰＵ１００において実行される。スペクトル解析プログラム１０４の実行過程において、ＣＰＵ１００が、ＮＭＲスペクトルを受け入れる手段、最適解としての最適係数セットを探索する探索手段、等として機能する。より詳しくは、ＣＰＵ１００が、計算スペクトル生成手段、残差スペクトル演算手段、残差ベクトル演算手段、評価値演算手段、係数セット更新手段、ｐ変更手段、等として機能する。
【００３３】
入力器１０６を利用してユーザーによって、初期値、パラメータ、設定値等が入力される。表示器１０８には解析結果が表示される。例えば、ＮＭＲ測定装置からＮＭＲスペクトル解析装置１２へ、ネットワークを介して、ＮＭＲスペクトルデータが転送される。可搬型記憶媒体を介して、ＮＭＲ測定装置からＮＭＲスペクトル解析装置１２へ、ＮＭＲスペクトルデータが送られてもよい。ＮＭＲスペクトル解析装置１２がＮＭＲ測定装置内に組み込まれてもよい。
【００３４】
図３には、第１実施例が示されている。第１実施例では、観測スペクトル（ベクトルｙ）１４及び計算スペクトル（ベクトルｙｍ（ａ））２０に基づいて、ＩＲＬＳ（Iterative Reweighted Least Square）アルゴリズム３８に従って、評価値Ｊの演算（符号４０参照）及び係数ベクトルａの更新（符号４２参照）が繰り返し実行される。これにより、評価値Ｊ最小化条件を満たす最適解として最適係数ベクトルが求められる。ＩＲＬＳアルゴリズム３８の詳細については後に図４を用いて説明する。なお、後に説明する第２実施例乃至第４実施例においても、ＩＲＬＳアルゴリズム３８が用いられている。もっとも、各実施例において、ｐが１以下である条件の下でＬｐノルム最小化問題を解くことが可能な他のアルゴリズムが利用されてもよい。
【００３５】
計算スペクトル２０は、例えば、以下の（４）式として示すローレンツ関数に対して係数ベクトル（係数セット）ａを与えることによって生成される。
【数４】
2019090747.tif 000006
【００３６】
上記（４）式に与えられる係数ベクトルａは、例えば、以下の（５）式に示すとおりである。
【数５】
2019090747.tif 000007
【００３７】
ｑはピーク数であり、係数の添え字はピーク番号である。ｃは振幅であり、ｄは半値幅であり、ωが周波数である。なお、係数ベクトルａは実際には列ベクトルである。ｑを係数ベクトルａに含めることも可能であるが、一般に、ｑはユーザーにより指定され、あるいは、試料等に基づいてｑが自動的に選択される。
【００３８】
計算スペクトルの生成に際しては、ローレンツ関数の他、ガウス関数、偽フォークト関数、等も利用され得る。溶液サンプルに対するＮＭＲ測定の場合、一般に、ローレンツ関数、偽フォークト関数等が利用される。固体サンプルに対するＮＭＲ測定の場合、一般に、ガウス関数、偽フォークト関数等が利用される。
【００３９】
評価値Ｊは、以下の（６）式によって定義される。
【数６】
2019090747.tif 000008
【００４０】
第１実施例において、ｐは１である。つまり、評価値ＪはＬ１ノルムである。変換行列Ａは空間変換を行うものであり、具体的には、それはＩＦＦＴ行列である。右辺中の１／２は実際の演算に際して無視し得る。それを考慮したとしても、評価値Ｊの実体はＬ１ノルムである。なお、ｐの値の指定や変換行列の指定には入力器が用いられる。初期係数ベクトルａ0の指定にも入力器が用いられる。必要なパラメータがプリセットされていてもよい。スペクトル解析に先立って、ユーザーに対してパラメータの指定又は確認を求める画面を提示するようにしてもよい。
【００４１】
既に説明した（１）式と上記の（６）式とを対比すると、（１）式には変換行列Ａが含まれていない。（１）式は周波数空間内で残差スペクトルを評価するものである。これに対し、（６）式は変換行列Ａを含み、時間空間内で、残差スペクトルに対応する残差ベクトルを評価するものである。また、（１）式においては、ｐが２であり、Ｌ２ノルムが評価値とされている。なお、Ｌ２ノルムは一般的なユーグリッドノルムであり、それにはノルム最小化に際してノルム演算対象のスパース性を強調する働きは認められない。これに対し、（６）においては、ｐが１であり、ノルム最小化に際してノルム演算対象のスパース性を強調する働きをもったＬ１ノルムが評価値Ｊとされている。ｐが１以下の場合にその作用が認められる。演算簡略化のためにはＬ１ノルムを評価値Ｊとするのが望ましい。
【００４２】
なお、ＩＲＬＳアルゴリズム３８の実行に先立って、各係数の初期値、並びに、各係数の上限値及び下限値、等がユーザーにより指定される。最大試行回数もユーザーにより指定される。
【００４３】
図４には、第１実施例に係るアルゴリズムの具体的内容が例示されている。Ｓ１０では初期設定が行われる。λはマルカート数であり、τは更新係数である。ｋmaxは最大試行回数である。εは０割防止用の係数である。初期係数ベクトルがａ0で指定されており、それを上記（４）式に代入し、上記（６）式を実行することにより、初期評価値Ｊ0が求められる。Ｓ１２では第１終了条件を満たすか否かが判断されている。ｋが最大値を超えるまで、Ｓ１４以下の計算が実行される。後述するように第２終了条件を組み込んでもよい。
【００４４】
最適解の探索に際しては、一般に、評価値計算式の一階微分及び二階微分が必要となる。Ｌ２ノルム計算式については、それに対する一階微分及び二階微分が可能であり、且つ、容易である。一方、Ｌ１ノルムの計算式は微分できないものであり、あるいは、それが可能であったとしても多大なる演算量を要する。ＩＲＬＳ法は、重み行列Ｗを導入することにより、Ｌ１ノルム計算式を微分可能なＬ２ノルム計算式に書き換えるものである。具体的には、ＩＲＬＳ法によれば、上記（５）式は以下の（７）式のように変形される。
【数７】
2019090747.tif 000009
【００４５】
ＩＲＬＳ法において、重み行列Ｗを利用する点が上記のＬＭ法と異なっている。換言すれば、アルゴリズムの基本部分はＬＭ法と共通である。
【００４６】
Ｓ１４では、ＩＲＬＳ法に従って対角要素ｗｉが定義される。Ｓ１６では、Ｓ１４で定義された対角要素ｗｉに基づいて重み行列Ｗが定義されている。diag([w1,…,wN])は[w1,…,wN]を対角要素とし、それ以外を0とした行列である。Ｓ１８に示されている式は、上記（７）式を二階微分することにより導出され、Ｈはヘッセ行列である。∇ｙｍ（ａk）は、（４）式を、係数ベクトルａを構成する各要素で偏微分することにより求められる。ＸＴは行列Ｘの転置行列を示しており、ＸＨは行列Ｘのエルミート転置である。Ｓ２０では、ヘッセ行列Ｈがその対角要素によって補正されている。Ｓ２２に示されている式は、上記（７）式を一階微分することにより導出され、ｈは勾配ベクトルである。
【００４７】
Ｓ２４では係数ベクトルａが更新される。その更新式はＬＭ法における更新式と同じである（ヘッセ行列Ｈが補正されている点を除いて、ＬＭ法の基礎をなすガウス・ニュートン法における更新式と同じである）。Ｈ＋は、行列Ｈの疑似逆行列である。
【００４８】
Ｓ２６では評価値Ｊが更新される。Ｓ２８では、今回の評価値Ｊk+1と前回の評価値Ｊkとが比較され、今回の評価値Ｊk+1が前回の評価値Ｊkよりも小さければ、Ｓ３０においてλが更新される（通常τは１以上の数値であり（例えば１０であり）、λが小さくされる）。今回の評価値Ｊk+1が前回の評価値Ｊkよりも小さくなければ、Ｓ３２においてλが更新される（λが大きくされる）。Ｓ３４では、更新後の係数ベクトルａk+1に前回の係数ベクトルａkが代入され、Ｓ３６では、更新後の評価値Ｊk+1に前回の評価値Ｊkが代入される。Ｓ３８では、ｋがインクリメントされる。
【００４９】
以上のアルゴリズムを繰り返し実行することにより、評価値Ｊを最小化する問題の解として、最適解つまり最適係数ベクトルが求まる。具体的には、ｋが最大値に到達した段階でアルゴリズムの実行が終了し、その時点での係数ベクトルａが最適係数ベクトルと判断される。
【００５０】
なお、第２終了条件を加える場合、例えば、Ｓ３８の直前において、以下の（８）式により解の変化量として指標ｅを求め、指標ｅが一定値以下の場合にアルゴリズムを終了させてもよい。
【数８】
2019090747.tif 000010
【００５１】
あるいは、以下の（９）式により、評価値Ｊの変化量として指標ｅを演算し、指標ｅが一定値以下の場合にアルゴリズムを終了させてもよい。第２終了条件として更に他のもの採用してもよい。
【数９】
2019090747.tif 000011
【００５２】
次に、図５及び図６を用いて、第１実施例の効果を説明する。元スペクトル（注目波形成分に相当する）とベースライン成分のエネルギー比ＳＢＲが以下の（１０）式で定義される。
【数１０】
2019090747.tif 000012
【００５３】
上記（１０）式において、ｘは元スペクトル（ベクトル）であり、ｂは元スペクトルに対して人工的に付加するベースライン成分（ベクトル）である。
【００５４】
図５には、上記ＳＢＲを１０とした場合における３種類のスペクトルが示されている。Ａ）に示すスペクトルは、余弦関数によって生成されたベースライン成分を含むスペクトルである。Ｂ）に示すスペクトルは、直線関数によって生成されたベースライン成分を含むスペクトルである。Ｃ）示すスペクトルは、ブロードなローレンツ信号として構成されたベースライン成分を含むスペクトルである。
【００５５】
図６には、ＳＢＲを変更しながら、上記の３種類のスペクトルに対して、以下の３種類のフィッティング法を適用した場合の推定誤差が示されている。ここで、ＬＭ法に基づくフィッティングを行った場合における推定誤差を示す曲線には複数の○記号が付されており（比較例１）、ＴＬＳ法に基づくフィッティングを行った場合における推定誤差を示す曲線には複数の×記号が付されており（比較例２）、第１実施例によるフィッティングを行った場合における推定誤差を示す曲線には複数の＋記号が付されている。推定誤差は、以下の（１１）式により計算される。
【数１１】
2019090747.tif 000013
【００５６】
上記（１１）式の分母は、元スペクトルのＬ２ノルムである。上記（１１）式の分子は、元スペクトルと計算スペクトルの差分のＬ２ノルムである。なお、図６に示されている推定誤差は、１０００回のフィッティングを行って得られた１０００個の推定誤差の平均値である。各手法において同じ初期パラメータを使用した。ＴＬＳ法においては、ベースライン関数を余弦関数及びＮ次多項式でフィッティングし、各関数に含まれるパラメータを１０個とし、それらのパラメータの初期値をいずれも０とした。
【００５７】
第１実施例とＬＭ法とを対比した場合、３つのスペクトルのいずれにおいても、一部を除いて、第１実施例の方が優良となっている。第１実施例とＴＬＳ法とを対比した場合、Ａ）においては、第１実施例とＴＬＳ法との間で推定精度の大差はない。時間空間において、ベースライン成分のスパース性の度合いが高い場合、第１実施例によれば、かなり良好な結果を得られる。しかも、第１実施例によれば、ＴＬＳ法と同じ推定精度（あるいはそれ以上の推定精度）を、ＴＬＳ法よりも少ないパラメータで実現できる。なお、Ｂ）及びＣ）においては、ＳＢＲが小さい範囲においてＴＬＳ法の方が優っているものの、第１実施例においても、概ね良好な推定精度を得られている。一部においては、第１実施例の方がＴＬＳ法よりも推定精度において優っている。
【００５８】
図７には第２実施例が示されている。図８には第３実施例が示されている。図９には第４実施例が示されている。なお、各図において、図３に示した要素と同様の要素には同一符号を付しその説明を省略する。
【００５９】
図７に示す第２実施例では、ｐとして０．７５が与えられている。この場合、Ｌｐノルムがノルム計算対象のスパース性をより強調させる作用を発揮することになる。例えば、時間軸上において、ベースライン成分に明確なスパース性が認められる場合には、第２実施例を採用し得る。ｐとして０．５を与えるようにしてもよい。ｐの値を連続的に又は段階的に変更できるように構成してもよい。
【００６０】
図８に示す第３実施例では、ｐとして１が与えられており、変換行列Ａとして、ＩＦＦＴ行列（第１変換行列）と、傾斜補正行列Ｄ（第２変換行列）と、が与えられている。残差ベクトルに対し、最初に傾斜補正行列Ｄが適用され、それに続いてＩＦＦＴ行列が適用される。傾斜補正行列は、周波数空間においてベースライン成分が傾斜成分を有する場合に、その傾斜成分を除外する作用を発揮するものである。傾斜補正行列として、公知の差分行列を利用し得る。傾斜補正がされた後の残差スペクトルが周波数空間に変換されることになる。なお、数式上は、残差スペクトルに対し、その左側から最初に傾斜補正行列が適用され、更にその左側からＩＦＦＴ行列が適用されることになる。
【００６１】
図９に示す第４実施例では、ｐとして１が与えられており、変換行列Ａとして、離散ウエーブレット変換行列ＤＷＴが与えられている。このように時間空間とは異なる他の空間へ残差スペクトルを変換し、残差ベクトルを得てもよい。但し、変換先の空間において、ベースライン成分がスパースな信号として表現されることが条件とされる。
【００６２】
以上のとおり、実施形態によれば、注目波形成分に対して計算スペクトルを高精度にフィッティングさせることが可能となる。すなわち、ベースライン成分の影響を受け難いスペクトルフィッティングを実現できる。ＮＭＲスペクトル以外のスペクトルに対して上記実施形態に係る手法を適用してもよい。
【符号の説明】
【００６３】
１０  ＮＭＲ測定装置、１２  ＮＭＲスペクトル解析装置、１４  ＮＭＲスペクトル（観測スペクトル）、１６  注目波形成分、１８  ベースライン成分、２０  計算スペクトル、２４  残差スペクトル、３０  残差ベクトル（残差信号）。
【図１】
2019090747.tif 000014
【図２】
2019090747.tif 000015
【図３】
2019090747.tif 000016
【図４】
2019090747.tif 000017
【図５】
2019090747.tif 000018
【図６】
2019090747.tif 000019
【図７】
2019090747.tif 000020
【図８】
2019090747.tif 000021
【図９】
2019090747.tif 000022