JP7018297B2 - スペクトル解析装置及び方法 - Google Patents

スペクトル解析装置及び方法 Download PDF

Info

Publication number
JP7018297B2
JP7018297B2 JP2017220787A JP2017220787A JP7018297B2 JP 7018297 B2 JP7018297 B2 JP 7018297B2 JP 2017220787 A JP2017220787 A JP 2017220787A JP 2017220787 A JP2017220787 A JP 2017220787A JP 7018297 B2 JP7018297 B2 JP 7018297B2
Authority
JP
Japan
Prior art keywords
spectrum
residual
evaluation value
space
calculated
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2017220787A
Other languages
English (en)
Other versions
JP2019090747A (ja
Inventor
朋喜 中尾
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Jeol Ltd
Original Assignee
Jeol Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Jeol Ltd filed Critical Jeol Ltd
Priority to JP2017220787A priority Critical patent/JP7018297B2/ja
Priority to US16/190,533 priority patent/US10712414B2/en
Publication of JP2019090747A publication Critical patent/JP2019090747A/ja
Application granted granted Critical
Publication of JP7018297B2 publication Critical patent/JP7018297B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01RMEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R33/00Arrangements or instruments for measuring magnetic variables
    • G01R33/20Arrangements or instruments for measuring magnetic variables involving magnetic resonance
    • G01R33/44Arrangements or instruments for measuring magnetic variables involving magnetic resonance using nuclear magnetic resonance [NMR]
    • G01R33/46NMR spectroscopy
    • G01R33/4625Processing of acquired signals, e.g. elimination of phase errors, baseline fitting, chemometric analysis

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Signal Processing (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Condensed Matter Physics & Semiconductors (AREA)
  • General Physics & Mathematics (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)
  • Investigating Or Analysing Materials By Optical Means (AREA)

Description

本発明は、スペクトル解析装置及び方法に関し、特に、スペクトルフィッティング技術に関する。
スペクトル解析の対象となるスペクトルとして、NMR(Nuclear Magnetic Resonance)スペクトル、X線スペクトル、分光スペクトル、マススペクトル等が挙げられる。スペクトル解析においては、例えば、観測スペクトル(Observed Spectrum)に計算スペクトル(Computed Spectrum, Simulated Spectrum)が最もフィッティングするように、計算スペクトルを定義する関数に与える最適な係数セットが探索される。最適な係数セットを求めるための手法として、各種の手法が提案されている。
例えば、スペクトルフィッティングに際しては、LM(Levenberg-Marquadts)法が利用される。LM法は、非線形最小二乗フィッティング法の一種であり、以下の評価値Jを最小にする係数セット(係数ベクトル)aを求めるものである(非特許文献1を参照)。
Figure 0007018297000001
上記において、yは観測スペクトル(ベクトル)であり、aは係数セットとしての係数ベクトルであり、y(a)は係数セットaをモデル関数に代入することよって生成される計算スペクトル(ベクトル)である。(1)式中の(y-ym(a))は残差スペクトル(ベクトル)である。(1)式の右辺にある1/2は無くてもよく、評価値Jの実体はL2ノルムである。
なお、一般に、Lpノルムに関しては、以下の(2)式及び(3)式に示す2つの表現が認められるが、本願において“Lpノルム”は基本的に(3)式で表現されるノルムである。例えば、上記L2ノルムは(3)式中のpに2を代入したものである。以下の(2)式及び(3)式に含まれる添え字nはベクトル要素数を示している。
Figure 0007018297000002
Figure 0007018297000003
スペクトルフィッティングに際しては、TLS(Total Line Shape fitting)法も利用される。TLS法は、NMRスペクトルの解析において一般的に用いられている手法であり、それも非線形最小二乗フィッティング法の一種である(非特許文献2-4を参照)。この方法では、計算スペクトルを定義する関数にベースライン関数が含まれる。
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).
観測スペクトルには、一般に、注目波形成分及びベースライン成分が含まれる。注目波形成分は、例えば、1又は複数のピークを含む部分であり、それが本来的な解析対象である。一方、ベースライン成分は、本来的な解析対象ではない成分であって、周波数空間において広帯域にわたって存在する成分である。例えば、NMRスペクトルにおいては、測定段階で生じるリンキングノイズ、デジタルフィルタ処理後のデータ欠損、分子構造由来の信号成分等がベースライン成分を生じさせる。NMRスペクトル中の基底を構成する大きな周期をもった変動、湾曲、傾斜等がベースライン成分であり、あるいは、そのような基底それ自体がベースライン成分である。
上記のベースライン成分を含むスペクトルに対して上記LM法を適用した場合、フィッティング精度が低下してしまう。ベースライン成分も、注目波形成分と同様に、フィッティング対象になってしまうからである。なお、TLS法においては、ベースライン関数に対しても少なからずの係数を与える必要がある。これにより、解析者の負担の増大、解析時間の増大、最適解を求められない可能性が高まる、といった問題が生じる。
本発明の目的は、注目波形成分に対して計算スペクトルを精度良くフィッティングさせることにある。あるいは、本発明の目的は、ベースライン成分の影響を受け難いスペクトルフィッティングを実現することにある。
実施形態に係るスペクトル解析装置は、観測スペクトルを受け入れる手段と、前記観測スペクトルと計算スペクトルとの間の残差スペクトルに対応する所定空間内の残差ベクトルについての評価値が最小化されるように、前記計算スペクトルを定義する関数に与える係数セットを変化させることにより、最適係数セットを探索する探索手段と、を含み、前記計算スペクトルは前記観測スペクトル中の注目波形成分にフィッティングされるものであり、前記所定空間は前記観測スペクトル中のベースライン成分がスパースな信号として表現される空間であり、前記評価値は前記所定空間内の残差ベクトルのLpノルム(但しp≦1)である、ことを特徴とする。
上記構成において、所定空間は、観測スペクトル中のベースライン成分がスパースな信号(多くのゼロを含んだ信号)として表現される空間である。すなわち、ベースライン成分は、所定空間軸上において、広域的に存在するのではなく、局所的に存在するものである。一方、Lpノルム(但しp≦1)を最小化する条件の下での最適解の探索において、当該Lpノルムは、ノルム演算対象のスパース性を増進又は強調する作用を発揮することが知られている。そのようなLpノルムを残差ベクトルの評価値とすれば、残差ベクトルのスパース性を高める方向に最適解の探索が進行する。例えば、残差ベクトル全体が均等に抑圧されるのではなく、残差ベクトル中により多くの0又は0近傍値が生じるように、残差ベクトルが抑圧される。結果として、ベースライン成分に対するフィッティングよりも、注目波形成分に対するフィッティングが優先されることになる。一般に、pが1以下の場合にその傾向が認められ、pを小さくすればするほどその傾向が強まる。よって、上記探索手段によれば、ベースライン成分に対してそれ専用の関数を当てはめることなく、注目波形部分に対する計算スペクトルのフィッティングの精度を高められる。
実施形態において、前記所定空間は時間空間である。例えば、周波数軸上の広帯域信号であるベースライン成分は、時間軸上において、幅の狭いパルス状の部分と、0又はほぼ0からなる平坦な部分と、を含む信号として表現される。すなわち、スパースな信号として表現される。上記構成はその性質を利用するものである。
実施形態において、前記所定空間内の残差ベクトルは前記残差スペクトルに対して変換行列を適用することにより求められる。変換行列は、例えば、周波数空間から時間空間への変換を行う行列である。その例として、IFFT(Inverse Fourier Transform)行列があげられる。
実施形態において、前記変換行列は、周波数空間から前記所定空間への変換を行う第1変換行列と、前記周波数空間内の前記残差スペクトルに対して傾斜補正を行う第2変換行列と、を含み、前記傾斜補正前の残差スペクトルに対して前記第2変換行列及び前記第1変換行列がその順で適用される。観測スペクトルにおいて、ベースライン成分が、大きな周期をもった変動を有しつつ、それ全体として傾斜している場合、傾斜の影響が残差ベクトルに現れ、最適解の探索が旨くいかないおそれが生じる。そこで、残差スペクトルに対して第2変換行列を適用し、傾斜成分を除去し又は減少させた上で、その適用結果に対して第1変換行列が適用される。
実施形態においては、前記pを変更する手段が設けられる。上記のように、一般に、pを小さくすればするほど、ベースライン成分に対するフィッティングよりも注目波形成分に対するフィッティングが優先される傾向が強まる。ベースライン成分の態様、必要なフィッティング精度、演算時間等に応じて、pを可変できるように構成するのが望ましい。
実施形態に係るスペクトル解析方法は、NMR測定で得られたNMRスペクトルを受け入れる工程と、前記NMRスペクトルと計算スペクトルとの間の残差スペクトルに対応する所定空間内の残差ベクトルについての評価値が最小化されるように、前記計算スペクトルを定義する関数に与える係数セットを変化させることにより、最適係数セットを探索する工程と、を含み、前記計算スペクトルは前記NMRスペクトル中の注目波形成分にフィッティングされるものであり、前記所定空間は前記NMRスペクトル中のベースライン成分がスパースな信号として表現される空間であり、前記評価値は前記所定空間内の残差ベクトルのLpノルム(但しp≦1)である。上記のNMRスペクトルを受け入れる工程、及び、最適係数セットを探索する工程は、望ましくは情報処理装置において実行される。
上記方法は、スペクトル解析対象ではないベースライン成分が所定空間内においてスパース性を呈することに着目し、且つ、Lpノルム(但しp≦1)がその最小化に際してノルム演算対象のスパース性を高めるように働くことに着目し、最適係数セットの探索において、所定空間内の残差ベクトルのLpノルムを評価値として用いるものである。すなわち、上記方法は、数学的な法則それ自体を特徴事項とするものではなく、スペクトルフィッティングにおいてフィッティング不要な部分が自然に除外されるようにする具体的な処理又は演算を特徴事項とするものである。上記スペクトル解析方法を実行するプログラムは、ネットワークを介して、あるいは、可搬型記憶媒体を介して、情報処理装置(スペクトル測定装置を含む)にインストールされ得る。
本発明によれば、注目波形成分に対して計算スペクトルを高精度にフィッティングさせることが可能となる。あるいは、本発明によれば、ベースライン成分の影響を受け難いスペクトルフィッティングを実現できる。
実施形態に係るNMRシステムを示す概念図である。 実施形態に係るスペクトル解析装置を示すブロック図である。 第1実施例に係るスペクトル解析方法を示す概念図である。 第1実施例に係るアルゴリズムを示す概念図である。 3つのベースライン成分が加えられた3つのスペクトルを示す図である。 3つのスペクトルに対するフィッティング結果を示す図である。 第2実施例に係るスペクトル解析方法を示す概念図である。 第3実施例に係るスペクトル解析方法を示す概念図である。 第4実施例に係るスペクトル解析方法を示す概念図である。
以下、実施形態を図面に基づいて説明する。
図1には、NMR測定装置10とNMRスペクトル解析装置12とからなるNMRシステムが示されている。NMR測定装置10は、大別して、測定部と分光計とにより構成される。測定部は、静磁場発生器及びプローブを含むものである。静磁場発生器は垂直貫通孔としてのボアを有し、そのボアの中にプローブにおける挿入部が挿入される。挿入部のヘッドにおいて、試料で生じた磁気共鳴が検出される。分光計は、送信部、受信部、制御部、演算部等からなるものである。送信部において送信パルスシーケンスに従う送信信号が生成され、その送信信号がプローブへ送られる。プローブからの受信信号が受信部に送られる。受信部において受信信号に基づいてNMRスペクトル(以下、観測スペクトルという。)14が生成される。
観測スペクトル14は、周波数空間(周波数領域)内に存在し、換言すれば、周波数軸上に波形として表現されるものである。観測スペクトル14をベクトルとして表現したものがベクトルyである。ベクトルyは、周波数軸方向に並ぶ複数の強度値からなる。観測スペクトル14は、一般に、注目波形成分16とベースライン成分18とを含む。
注目波形成分16は、一般に、1又は複数のピークを含んだ部分である。典型的なスペクトル解析では、ピーク単位で、振幅、周波数、半値幅等のパラメータが推定される。一方、ベースライン成分18は、フィッティング対象外の部分であり、一般に、それは広帯域成分である。具体的には、ベースライン成分18は、観測スペクトル14の基底を構成する比較的大きな周期をもった変動、湾曲、傾斜等に相当する。あるいは、そのような成分をもった基底自体をベースライン成分として理解してもよい。スペクトル解析においては、そのようなベースライン成分18に影響されずに、注目波形成分16に対して正確なフィッティングを行うことが求められる。
スペクトル解析装置12は、NMR測定装置10で生成された観測スペクトル14を解析する装置であり、それはパーソナルコンピュータ等の情報処理装置により構成される。スペクトル解析装置12は、評価値Jが最小化されるように係数セット(係数ベクトル)aの内容を段階的に優良化することにより、最適解としての最適係数セットを探索するものである。実施形態においては、最適係数セットがスペクトル解析結果を構成する。
具体的には、係数セットaをスペクトル生成関数(モデル関数)に与えることにより、計算スペクトル20が生成される。係数セットaの初期値がa0である。計算スペクトル20をベクトルとして表現したものがベクトルym(a)である。ベクトルym(a)は周波数軸上に並ぶ複数の強度値からなる。図1においては、計算スペクトル20を生成する部分に符号22が付されている。観測スペクトル14と計算スペクトル20の差が残差スペクトル24である。図1においては、残差スペクトル24を演算する部分に符号26が付されている。
残差スペクトル24は周波数空間内に存在するものであり、それを時間空間に表現したものが残差ベクトル(残差信号)30である。その際には残差スペクトル24に対して逆フーリエ変換(IFFT)等が適用される。残差ベクトル30は、時間軸方向に並ぶ複数の振幅値により構成される。図示の例において、残差ベクトル30には、局所的に生じた短い時間幅をもったピーク部分30aと広域に及ぶ平坦な部分30bとが含まれる。ピーク部分30aが周波数空間内においてベースライン成分18の実体をなす。逆に言えば、ベースライン成分18は、時間空間においてスパースな信号として表現されるものである。平坦な部分30bにはランダムノイズが含まれる。図1においては、残差ベクトル30を演算する部分に符号28が付されている。
実施形態においては、評価値Jとして、時間軸上の残差ベクトル30のLpノルム(但しp≦1)が利用される。残差ベクトル30のLpノルムが最小化されるように、係数ベクトル(係数セット)aの内容が順次更新され(符号34参照)、最終的に、最適解としての最適係数セットが特定される。その最適係数セットが外部に出力され(符号36参照)、あるいは、表示される。その最適係数セットが他の演算に利用されてもよい。図1においては、評価値Jを計算し且つ係数ベクトルaを更新する部分に符号32が付されている。
図1においては、最適係数セットの探索過程が模式的に示されているが、実際には、一連の計算式を繰り返し実行することにより、最適係数セットが特定される。例えば、残差スペクトルの計算(符号26参照)、残差ベクトルの計算(符号28参照)、及び、評価値J(Lpノルム)の計算(符号32参照)は、単一の計算式により表現され得る。
評価値としてのLpノルム(但しp≦1)は、その最小化に際して、ノルム計算対象である残差ベクトル30のスパース性を高めるように働く。図示の例では、ピーク部分30aの抑圧よりも、平坦な部分30bの抑圧及び増大が優先されることになる。もっとも、図示された残差ベクトル30は、最小化条件が満たされた時点でのものである。なお、平坦な部分30aにはランダムノイズが含まれ、通常、それは完全に0にはならない。
以上のように、実施形態では、残差スペクトル24それ自体を評価対象とするのではなく、時間空間内の残差ベクトル30を評価対象としており、かつ、評価値としてLpノルム(但しp≦1)を利用している。これにより、最適解の探索に際して、ベースライン成分18へのフィッティングよりも、注目波形成分16へのフィッティングの方が優先される。換言すれば、スペクトルフィッティングにおいて、ベースライン成分18が無視され易くなる。その結果、注目波形成分16に対して計算スペクトル20を正しくフィッティングさせることが可能となる。なお、実施形態によれば、ベースライン成分18を模擬する関数を用いる必要がないので、パラメータ指定に際してのユーザーの負担が大幅に軽減される。
図2には、スペクトル解析装置12の構成例が示されている。スペクトル解析装置12は、CPU100、メモリ102、入力器106、表示器108等を含む。メモリ102にはスペクトル解析プログラム104が格納されている。スペクトル解析プログラム104がCPU100において実行される。スペクトル解析プログラム104の実行過程において、CPU100が、NMRスペクトルを受け入れる手段、最適解としての最適係数セットを探索する探索手段、等として機能する。より詳しくは、CPU100が、計算スペクトル生成手段、残差スペクトル演算手段、残差ベクトル演算手段、評価値演算手段、係数セット更新手段、p変更手段、等として機能する。
入力器106を利用してユーザーによって、初期値、パラメータ、設定値等が入力される。表示器108には解析結果が表示される。例えば、NMR測定装置からNMRスペクトル解析装置12へ、ネットワークを介して、NMRスペクトルデータが転送される。可搬型記憶媒体を介して、NMR測定装置からNMRスペクトル解析装置12へ、NMRスペクトルデータが送られてもよい。NMRスペクトル解析装置12がNMR測定装置内に組み込まれてもよい。
図3には、第1実施例が示されている。第1実施例では、観測スペクトル(ベクトルy)14及び計算スペクトル(ベクトルym(a))20に基づいて、IRLS(Iterative Reweighted Least Square)アルゴリズム38に従って、評価値Jの演算(符号40参照)及び係数ベクトルaの更新(符号42参照)が繰り返し実行される。これにより、評価値J最小化条件を満たす最適解として最適係数ベクトルが求められる。IRLSアルゴリズム38の詳細については後に図4を用いて説明する。なお、後に説明する第2実施例乃至第4実施例においても、IRLSアルゴリズム38が用いられている。もっとも、各実施例において、pが1以下である条件の下でLpノルム最小化問題を解くことが可能な他のアルゴリズムが利用されてもよい。
計算スペクトル20は、例えば、以下の(4)式として示すローレンツ関数に対して係数ベクトル(係数セット)aを与えることによって生成される。
Figure 0007018297000004
上記(4)式に与えられる係数ベクトルaは、例えば、以下の(5)式に示すとおりである。
Figure 0007018297000005
qはピーク数であり、係数の添え字はピーク番号である。cは振幅であり、dは半値幅であり、ωが周波数である。なお、係数ベクトルaは実際には列ベクトルである。qを係数ベクトルaに含めることも可能であるが、一般に、qはユーザーにより指定され、あるいは、試料等に基づいてqが自動的に選択される。
計算スペクトルの生成に際しては、ローレンツ関数の他、ガウス関数、偽フォークト関数、等も利用され得る。溶液サンプルに対するNMR測定の場合、一般に、ローレンツ関数、偽フォークト関数等が利用される。固体サンプルに対するNMR測定の場合、一般に、ガウス関数、偽フォークト関数等が利用される。
評価値Jは、以下の(6)式によって定義される。
Figure 0007018297000006
第1実施例において、pは1である。つまり、評価値JはL1ノルムである。変換行列Aは空間変換を行うものであり、具体的には、それはIFFT行列である。右辺中の1/2は実際の演算に際して無視し得る。それを考慮したとしても、評価値Jの実体はL1ノルムである。なお、pの値の指定や変換行列の指定には入力器が用いられる。初期係数ベクトルa0の指定にも入力器が用いられる。必要なパラメータがプリセットされていてもよい。スペクトル解析に先立って、ユーザーに対してパラメータの指定又は確認を求める画面を提示するようにしてもよい。
既に説明した(1)式と上記の(6)式とを対比すると、(1)式には変換行列Aが含まれていない。(1)式は周波数空間内で残差スペクトルを評価するものである。これに対し、(6)式は変換行列Aを含み、時間空間内で、残差スペクトルに対応する残差ベクトルを評価するものである。また、(1)式においては、pが2であり、L2ノルムが評価値とされている。なお、L2ノルムは一般的なユーグリッドノルムであり、それにはノルム最小化に際してノルム演算対象のスパース性を強調する働きは認められない。これに対し、(6)においては、pが1であり、ノルム最小化に際してノルム演算対象のスパース性を強調する働きをもったL1ノルムが評価値Jとされている。pが1以下の場合にその作用が認められる。演算簡略化のためにはL1ノルムを評価値Jとするのが望ましい。
なお、IRLSアルゴリズム38の実行に先立って、各係数の初期値、並びに、各係数の上限値及び下限値、等がユーザーにより指定される。最大試行回数もユーザーにより指定される。
図4には、第1実施例に係るアルゴリズムの具体的内容が例示されている。S10では初期設定が行われる。λはマルカート数であり、τは更新係数である。kmaxは最大試行回数である。εは0割防止用の係数である。初期係数ベクトルがa0で指定されており、それを上記(4)式に代入し、上記(6)式を実行することにより、初期評価値J0が求められる。S12では第1終了条件を満たすか否かが判断されている。kが最大値を超えるまで、S14以下の計算が実行される。後述するように第2終了条件を組み込んでもよい。
最適解の探索に際しては、一般に、評価値計算式の一階微分及び二階微分が必要となる。L2ノルム計算式については、それに対する一階微分及び二階微分が可能であり、且つ、容易である。一方、L1ノルムの計算式は微分できないものであり、あるいは、それが可能であったとしても多大なる演算量を要する。IRLS法は、重み行列Wを導入することにより、L1ノルム計算式を微分可能なL2ノルム計算式に書き換えるものである。具体的には、IRLS法によれば、上記(5)式は以下の(7)式のように変形される。
Figure 0007018297000007
IRLS法において、重み行列Wを利用する点が上記のLM法と異なっている。換言すれば、アルゴリズムの基本部分はLM法と共通である。
S14では、IRLS法に従って対角要素wが定義される。S16では、S14で定義された対角要素wに基づいて重み行列Wが定義されている。diag([w1,…,wN])は[w1,…,wN]を対角要素とし、それ以外を0とした行列である。S18に示されている式は、上記(7)式を二階微分することにより導出され、Hはヘッセ行列である。∇ym(ak)は、(4)式を、係数ベクトルaを構成する各要素で偏微分することにより求められる。Xは行列Xの転置行列を示しており、Xは行列Xのエルミート転置である。S20では、ヘッセ行列Hがその対角要素によって補正されている。S22に示されている式は、上記(7)式を一階微分することにより導出され、hは勾配ベクトルである。
S24では係数ベクトルaが更新される。その更新式はLM法における更新式と同じである(ヘッセ行列Hが補正されている点を除いて、LM法の基礎をなすガウス・ニュートン法における更新式と同じである)。Hは、行列Hの疑似逆行列である。
S26では評価値Jが更新される。S28では、今回の評価値Jk+1と前回の評価値Jkとが比較され、今回の評価値Jk+1が前回の評価値Jkよりも小さければ、S30においてλが更新される(通常τは1以上の数値であり(例えば10であり)、λが小さくされる)。今回の評価値Jk+1が前回の評価値Jkよりも小さくなければ、S32においてλが更新される(λが大きくされる)。S34では、更新後の係数ベクトルak+1に前回の係数ベクトルakが代入され、S36では、更新後の評価値Jk+1に前回の評価値Jkが代入される。S38では、kがインクリメントされる。
以上のアルゴリズムを繰り返し実行することにより、評価値Jを最小化する問題の解として、最適解つまり最適係数ベクトルが求まる。具体的には、kが最大値に到達した段階でアルゴリズムの実行が終了し、その時点での係数ベクトルaが最適係数ベクトルと判断される。
なお、第2終了条件を加える場合、例えば、S38の直前において、以下の(8)式により解の変化量として指標eを求め、指標eが一定値以下の場合にアルゴリズムを終了させてもよい。
Figure 0007018297000008
あるいは、以下の(9)式により、評価値Jの変化量として指標eを演算し、指標eが一定値以下の場合にアルゴリズムを終了させてもよい。第2終了条件として更に他のもの採用してもよい。
Figure 0007018297000009
次に、図5及び図6を用いて、第1実施例の効果を説明する。元スペクトル(注目波形成分に相当する)とベースライン成分のエネルギー比SBRが以下の(10)式で定義される。
Figure 0007018297000010
上記(10)式において、xは元スペクトル(ベクトル)であり、bは元スペクトルに対して人工的に付加するベースライン成分(ベクトル)である。
図5には、上記SBRを10とした場合における3種類のスペクトルが示されている。A)に示すスペクトルは、余弦関数によって生成されたベースライン成分を含むスペクトルである。B)に示すスペクトルは、直線関数によって生成されたベースライン成分を含むスペクトルである。C)示すスペクトルは、ブロードなローレンツ信号として構成されたベースライン成分を含むスペクトルである。
図6には、SBRを変更しながら、上記の3種類のスペクトルに対して、以下の3種類のフィッティング法を適用した場合の推定誤差が示されている。ここで、LM法に基づくフィッティングを行った場合における推定誤差を示す曲線には複数の○記号が付されており(比較例1)、TLS法に基づくフィッティングを行った場合における推定誤差を示す曲線には複数の×記号が付されており(比較例2)、第1実施例によるフィッティングを行った場合における推定誤差を示す曲線には複数の+記号が付されている。推定誤差は、以下の(11)式により計算される。
Figure 0007018297000011
上記(11)式の分母は、元スペクトルのL2ノルムである。上記(11)式の分子は、元スペクトルと計算スペクトルの差分のL2ノルムである。なお、図6に示されている推定誤差は、1000回のフィッティングを行って得られた1000個の推定誤差の平均値である。各手法において同じ初期パラメータを使用した。TLS法においては、ベースライン関数を余弦関数及びN次多項式でフィッティングし、各関数に含まれるパラメータを10個とし、それらのパラメータの初期値をいずれも0とした。
第1実施例とLM法とを対比した場合、3つのスペクトルのいずれにおいても、一部を除いて、第1実施例の方が優良となっている。第1実施例とTLS法とを対比した場合、A)においては、第1実施例とTLS法との間で推定精度の大差はない。時間空間において、ベースライン成分のスパース性の度合いが高い場合、第1実施例によれば、かなり良好な結果を得られる。しかも、第1実施例によれば、TLS法と同じ推定精度(あるいはそれ以上の推定精度)を、TLS法よりも少ないパラメータで実現できる。なお、B)及びC)においては、SBRが小さい範囲においてTLS法の方が優っているものの、第1実施例においても、概ね良好な推定精度を得られている。一部においては、第1実施例の方がTLS法よりも推定精度において優っている。
図7には第2実施例が示されている。図8には第3実施例が示されている。図9には第4実施例が示されている。なお、各図において、図3に示した要素と同様の要素には同一符号を付しその説明を省略する。
図7に示す第2実施例では、pとして0.75が与えられている。この場合、Lpノルムがノルム計算対象のスパース性をより強調させる作用を発揮することになる。例えば、時間軸上において、ベースライン成分に明確なスパース性が認められる場合には、第2実施例を採用し得る。pとして0.5を与えるようにしてもよい。pの値を連続的に又は段階的に変更できるように構成してもよい。
図8に示す第3実施例では、pとして1が与えられており、変換行列Aとして、IFFT行列(第1変換行列)と、傾斜補正行列D(第2変換行列)と、が与えられている。残差ベクトルに対し、最初に傾斜補正行列Dが適用され、それに続いてIFFT行列が適用される。傾斜補正行列は、周波数空間においてベースライン成分が傾斜成分を有する場合に、その傾斜成分を除外する作用を発揮するものである。傾斜補正行列として、公知の差分行列を利用し得る。傾斜補正がされた後の残差スペクトルが周波数空間に変換されることになる。なお、数式上は、残差スペクトルに対し、その左側から最初に傾斜補正行列が適用され、更にその左側からIFFT行列が適用されることになる。
図9に示す第4実施例では、pとして1が与えられており、変換行列Aとして、離散ウエーブレット変換行列DWTが与えられている。このように時間空間とは異なる他の空間へ残差スペクトルを変換し、残差ベクトルを得てもよい。但し、変換先の空間において、ベースライン成分がスパースな信号として表現されることが条件とされる。
以上のとおり、実施形態によれば、注目波形成分に対して計算スペクトルを高精度にフィッティングさせることが可能となる。すなわち、ベースライン成分の影響を受け難いスペクトルフィッティングを実現できる。NMRスペクトル以外のスペクトルに対して上記実施形態に係る手法を適用してもよい。
10 NMR測定装置、12 NMRスペクトル解析装置、14 NMRスペクトル(観測スペクトル)、16 注目波形成分、18 ベースライン成分、20 計算スペクトル、24 残差スペクトル、30 残差ベクトル(残差信号)。

Claims (9)

  1. 周波数空間内の観測スペクトルを受け入れる手段と、
    前記観測スペクトルと計算スペクトルとの間の残差スペクトルの変換により求められる時間空間内の残差ベクトルについての評価値が最小化されるように、前記計算スペクトルを定義する関数に与える係数セットを変化させることにより、最適係数セットを探索する探索手段と、
    を含み、
    前記計算スペクトルは前記観測スペクトル中の注目波形成分にフィッティングされるものであり
    記評価値は前記時間空間内の残差ベクトルのLpノルム(但しp≦1)である、
    ことを特徴とするスペクトル解析装置。
  2. 請求項1記載の装置において、
    前記pは、1.0、0.75又は0.5である、
    ことを特徴とするスペクトル解析装置。
  3. 請求項1記載の装置において、
    前記時間空間内の残差ベクトルは前記残差スペクトルに対して変換行列を適用することにより求められる、
    ことを特徴とするスペクトル解析装置。
  4. 請求項3記載の装置において、
    前記変換行列は、前記周波数空間から前記時間空間への変換を行う第1変換行列を含む、
    ことを特徴とするスペクトル解析装置。
  5. 請求項4記載の装置において、
    前記変換行列は、更に、前記周波数空間内の前記残差スペクトルに対して傾斜補正を行う第2変換行列を含み、
    前記傾斜補正前の残差スペクトルに対して前記第2変換行列及び前記第1変換行列がその順で適用される、
    ことを特徴とするスペクトル解析装置。
  6. 請求項1記載の装置において、
    前記pを変更する手段を含む、
    ことを特徴とするスペクトル解析装置。
  7. 請求項1記載の装置において、
    前記観測スペクトルはNMRスペクトルである、
    ことを特徴とするスペクトル解析装置。
  8. NMR測定で得られた周波数空間内のNMRスペクトルを受け入れる工程と、
    前記NMRスペクトルと計算スペクトルとの間の残差スペクトルの変換により求められる時間空間内の残差ベクトルについての評価値が最小化されるように、前記計算スペクトルを定義する関数に与える係数セットを変化させることにより、最適係数セットを探索する工程と、
    を含み、
    前記計算スペクトルは前記NMRスペクトル中の注目波形成分にフィッティングされるものであり
    記評価値は前記時間空間内の残差ベクトルのLpノルム(但しp≦1)である、
    ことを特徴とするスペクトル解析方法。
  9. 情報処理装置においてスペクトル解析方法を実行するプログラムであって、
    前記スペクトル解析方法は、
    NMR測定で得られた周波数空間内のNMRスペクトルを受け入れる工程と、
    前記NMRスペクトルと計算スペクトルとの間の残差スペクトルの変換により求められる時間空間内の残差ベクトルについての評価値が最小化されるように、前記計算スペクトルを定義する関数に与える係数セットを変化させることにより、最適係数セットを探索する工程と、
    を含み、
    前記計算スペクトルは前記NMRスペクトル中の注目波形成分にフィッティングされるものであり
    記評価値は前記時間空間内の残差ベクトルのLpノルム(但しp≦1)である、
    ことを特徴とするプログラム。
JP2017220787A 2017-11-16 2017-11-16 スペクトル解析装置及び方法 Active JP7018297B2 (ja)

Priority Applications (2)

Application Number Priority Date Filing Date Title
JP2017220787A JP7018297B2 (ja) 2017-11-16 2017-11-16 スペクトル解析装置及び方法
US16/190,533 US10712414B2 (en) 2017-11-16 2018-11-14 Apparatus and method for analyzing spectrum

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2017220787A JP7018297B2 (ja) 2017-11-16 2017-11-16 スペクトル解析装置及び方法

Publications (2)

Publication Number Publication Date
JP2019090747A JP2019090747A (ja) 2019-06-13
JP7018297B2 true JP7018297B2 (ja) 2022-02-10

Family

ID=66534267

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2017220787A Active JP7018297B2 (ja) 2017-11-16 2017-11-16 スペクトル解析装置及び方法

Country Status (2)

Country Link
US (1) US10712414B2 (ja)
JP (1) JP7018297B2 (ja)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP6817243B2 (ja) * 2018-02-28 2021-01-20 日本電子株式会社 スペクトル処理装置及び方法
JP7390210B2 (ja) * 2020-02-28 2023-12-01 日本電子株式会社 スペクトル処理装置及び方法
WO2022024263A1 (ja) * 2020-07-29 2022-02-03 中国電力株式会社 評価装置及び評価システム
CN116522269B (zh) * 2023-06-28 2023-09-19 厦门炬研电子科技有限公司 一种基于Lp范数非平稳信号稀疏重建的故障诊断方法

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2001159603A (ja) 1999-11-30 2001-06-12 Semiconductor Leading Edge Technologies Inc 分光スペクトル帰属同定方法および分光スペクトル帰属同定装置
JP2008048878A (ja) 2006-08-24 2008-03-06 Hitachi Medical Corp 磁気共鳴イメージング装置
JP2015128568A (ja) 2013-12-05 2015-07-16 キヤノン株式会社 劣決定系線型方程式の求解方法及び求解装置
JP2016057065A (ja) 2014-09-05 2016-04-21 株式会社分光科学研究所 分光定量方法、分光定量装置及びプログラム
JP2016224031A (ja) 2015-04-22 2016-12-28 ザ・ボーイング・カンパニーThe Boeing Company フォービエイテッド圧縮投影を用いたハイパースペクトルデミキシング
WO2017119086A1 (ja) 2016-01-06 2017-07-13 株式会社島津製作所 クロマトグラムデータ処理方法及び装置
JP2017129389A (ja) 2016-01-19 2017-07-27 国立大学法人島根大学 ラマン散乱光の観測方法、ラマン散乱光の観測装置及びプログラム

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH04364490A (ja) * 1991-06-11 1992-12-16 Siemens Asahi Medeitetsuku Kk Nmrスペクトルの基線補正装置
WO2008004324A1 (fr) * 2006-07-07 2008-01-10 Hitachi Medical Corporation Dispositif d'imagerie par résonance magnétique

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2001159603A (ja) 1999-11-30 2001-06-12 Semiconductor Leading Edge Technologies Inc 分光スペクトル帰属同定方法および分光スペクトル帰属同定装置
JP2008048878A (ja) 2006-08-24 2008-03-06 Hitachi Medical Corp 磁気共鳴イメージング装置
JP2015128568A (ja) 2013-12-05 2015-07-16 キヤノン株式会社 劣決定系線型方程式の求解方法及び求解装置
JP2016057065A (ja) 2014-09-05 2016-04-21 株式会社分光科学研究所 分光定量方法、分光定量装置及びプログラム
JP2016224031A (ja) 2015-04-22 2016-12-28 ザ・ボーイング・カンパニーThe Boeing Company フォービエイテッド圧縮投影を用いたハイパースペクトルデミキシング
WO2017119086A1 (ja) 2016-01-06 2017-07-13 株式会社島津製作所 クロマトグラムデータ処理方法及び装置
JP2017129389A (ja) 2016-01-19 2017-07-27 国立大学法人島根大学 ラマン散乱光の観測方法、ラマン散乱光の観測装置及びプログラム

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
圧縮センシングの数理,IEICE Fundamentals Review,電子情報通信学会,2010年07月,Vol.4 No.1,pp.39-47

Also Published As

Publication number Publication date
US20190154778A1 (en) 2019-05-23
JP2019090747A (ja) 2019-06-13
US10712414B2 (en) 2020-07-14

Similar Documents

Publication Publication Date Title
JP7018297B2 (ja) スペクトル解析装置及び方法
Hamaus et al. Minimizing the stochasticity of halos in large-scale structure surveys
Piskunov et al. Spectroscopy made easy: evolution
EP0535700A2 (en) Method and apparatus for comparing spectra
Haney et al. Causal instrument corrections for short‐period and broadband seismometers
US7429860B2 (en) Noise reduction for spectroscopic signal processing
Liland et al. ROSA—a fast extension of partial least squares regression for multiblock data analysis
CN110929217A (zh) 一种含有毛刺、尖峰干扰频响曲线的矫正方法及装置
Šprlák et al. Spherical integral formulas for upward/downward continuation of gravitational gradients onto gravitational gradients
WO2012039061A1 (ja) データ処理方法及び装置
JP7018321B2 (ja) スペクトル処理装置及び方法
US10094899B2 (en) Magnetic resonance T2 preparatory pulses for magnetic field inhomogeneity robustness and contrast optimization
Jazia et al. Wave propagation in a fractional viscoelastic Andrade medium: diffusive approximation and numerical modeling
Zhang et al. The attenuated Ricker wavelet basis for seismic trace decomposition and attenuation analysis
JP6817243B2 (ja) スペクトル処理装置及び方法
CN111505709A (zh) 一种基于稀疏谱分解的衰减定性分析的方法
Martins et al. ALiCCE: atomic lines calibration using the cross-entropy algorithm
Tan et al. A non-parametric approach for linear system identification using principal component analysis
WO2021210228A1 (ja) 波形情報推定方法及び装置、並びに、ピーク波形処理方法及び装置
JP2006194761A (ja) Nmr測定方法および装置
Lebsack et al. Iterative RF pulse refinement for magnetic resonance imaging
Rutledge et al. PoLiSh—smoothed partial least-squares regression
Jenkins et al. The effect of numerical model error on data assimilation
JP6160366B2 (ja) Nmr信号処理システム
Leles et al. Comparison between basic and toeplitiz ssa applied to non-stationary time-series

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20200717

A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20210520

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20210601

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20210702

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20220104

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20220131

R150 Certificate of patent or registration of utility model

Ref document number: 7018297

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150