(57)【特許請求の範囲】
【請求項１】
測定されたｎ／２点のスペクトルデータに対して平滑化フィルタを適用して平滑化スペクトルを得る平滑化工程と、
前記平滑化スペクトルと測定スペクトルとの差スペクトルを演算する差スペクトル演算工程と、
差スペクトルデータを対称に折返してｎ点のデータとするデータ調製工程と、
フーリエ変換によりｎ点の対称インターフェログラムデータを得るフーリエ変換工程と、
差スペクトルデータに基づくｎ点のインターフェログラムから大きさｎ×ｐ（ｐ＜ｎ）の巡回行列Ｘを作成する巡回行列作成工程と、
前記巡回行列Ｘに特異値分解を施し、所定値以上の値を持つＬ（＜ｐ）個の特異値に相当する行列成分Ｓを抽出する特異値分解工程と、
前記行列Ｓの各列ベクトルを巡回的にシフトさせた行列の各行の平均値から得られるｎ点のインターフェログラムを逆フーリエ変換することにより、ノイズが除去された演算差スペクトルを得る演算差スペクトル再生工程と、
前記平滑化スペクトルと演算差スペクトルの和合成を行う合成工程と、
を備えることを特徴とするスペクトルデータからのノイズ除去方法。
【請求項２】
請求項１記載の方法において、
平滑化スペクトルと差スペクトルの各々をＭ分割する分割工程を有し、
Ｍ分割された各々の部分差スペクトルに対して等しい値のｐパラメータ値で巡回行列作成工程、特異値分解工程及び演算差スペクトル再生工程、合成工程を行い、得られた部分再生スペクトルを結ぶことを特徴とするスペクトルデータからのノイズ除去方法。
【請求項３】
請求項１に記載の方法において、ｐがｎ／８≦ｐ≦ｎ／３であることを特徴とするスペクトルデータからのノイズ除去方法。
【請求項４】
請求項２に記載の方法において、ｐがｎ／８Ｍ≦ｐ≦ｎ／３Ｍであることを特徴とするスペクトルデータからのノイズ除去方法。
【発明の詳細な説明】
【０００１】
【産業上の利用分野】
本発明はノイズ除去方法、特に連続スペクトルデータからの効率的なノイズ除去方法に関する。
【０００２】
【従来の技術】
スペクトルデータの信号−ノイズ比を改善する一般的な方法として、サビツキー（Savitzky）及びゴレイ（Golay）により提案されたものがあり、これは最小二乗法に基づく多項式曲線近似法である。
【０００３】
【発明が解決しようとする課題】
しかしながら、このノイズ除去法は、信号及びノイズの統計学的特性を考慮しなければ、ノイズを抑制するというよりもむしろ信号の波形を劣化を引起こしてしまう場合がある。このような主な問題点としては、スペクトル線幅の拡大と位相シフトが挙げられる。大きな信号ピークの近傍の弱い信号ピークを検出する場合には、これらの問題点を確実に除去しなければならない。
本発明は前記従来技術の課題に鑑みなされたものであり、その目的は信号波形の劣化を伴うことなく、ノイズを除去することにある。
【０００４】
【課題を解決するための手段】
前記目的を達成するために本発明にかかるスペクトルデータのノイズ除去方法は、平滑化工程と、差スペクトル演算工程と、データ調整工程と、フーリエ変換工程と、巡回行列作成工程と、特異値分解工程と、演算差スペクトル再生工程と、合成工程とを備える。
前記平滑化工程は、測定されたｎ／２点のスペクトルデータに対して平滑化フィルタを適用して平滑化スペクトルを得る。
前記差スペクトル演算工程は、平滑化スペクトルと測定スペクトルとの差スペクトルを演算する。
【０００５】
前記データ調整工程では、差スペクトルデータを対称に折返してｎ点のデータとする。
前記フーリエ変換工程では、フーリエ変換によりｎ点の対称インターフェログラムデータを得る。
そして、巡回行列作成工程では、スペクトルデータに基づくｎ点のインターフェログラムから大きさｎ×ｐ（ｐ＜ｎ）の巡回行列Ｘを作成する。
特異値分解工程では、前記巡回行列Ｘに特異値分解を施し、所定値以上の値を持つＬ（＜ｐ）個の特異値に相当する行列成分Ｓを抽出する。
【０００６】
演算差スペクトル再生工程では、前記行列Ｓの各列ベクトルを巡回的にシフトさせた行列の各行の平均値から得られるｎ点のインターフェログラムを逆フーリエ変換することにより、ノイズが除去された演算差スペクトルを得る。
合成工程では、前記平滑化スペクトルと演算差スペクトルの和合成を行う。
以上の請求項１記載の方法によれば、差スペクトルに対して特異値分解法を適用するため、輝線スペクトルではない連続スペクトルに対しても短時間で良好なノイズ除去を行うことができる。
【０００７】
すなわち、差スペクトル演算工程と、スペクトルと合成工程を含まないで、巡回行列作成工程と特異値分解工程のみを用いたノイズ除去方法は、原理的には連続スペクトルにも適用できるが、連続スペクトルに適用する場合にはｐとＬのパラメータ値を非常に大きくする必要がある。さらに、Ｌを大きくしても、信号とノイズの区別がつきにくい場合が多い。このため、計算時間が長くなり、しかもあまり良好な結果が得られない場合がある。
【０００８】
本方法によれば、連続スペクトルより輝線スペクトル的な差スペクトル部分を抽出し、該差スペクトルに特異値分解を適用するため、Ｌの値が小さくでき、したがって計算時間が短くしかも良好な結果を得ることができる。
又、請求項２記載の方法は、さらに平滑化スペクトルと差スペクトルの各々をＭ分割する分割工程を有し、
Ｍ分割された各々の部分差スペクトルに対して等しい値のｐパラメータ値で巡回行列作成工程、特異値分解工程及び演算差スペクトル再生工程及び合成工程を行い、得られた部分再生スペクトルを結ぶことを特徴とする。
【０００９】
この請求項２の方法によれば、さらに計算時間の短縮を行うことができる。すなわち、特にデータ点数が多い場合、一般に特異値分解に要する時間は著しく増加する。しかしながら、請求項２記載の方法によれば、Ｍ分割した場合には、ｐの値が小さくできるため、分割していない場合と比較して計算時間は約１／Ｍとなる。
また、ｐの値は測定データを分割しない場合にはｎ／８≦ｐ＜ｎ／３であることが好適である。測定データをＭ分割した場合にはｎ／８Ｍ≦ｐ＜ｎ／３Ｍであることが好適である。
【００１０】
【実施例】
以下、図面に基づき本発明の好適な実施例を説明する。
［第１実施例］
本実施例においては、たとえば赤外吸収スペクトルなどの連続スペクトルからノイズを除去するため、次のステップを要件とする。
▲１▼平滑化工程：測定されたｎ／２点のスペクトルデータに対して幅Ｗ点の平滑化フィルタを適用して平滑化スペクトルを得る。
▲２▼差スペクトル演算工程：前記平滑化スペクトルと測定スペクトルとの差スペクトルを演算する。
【００１１】
▲３▼巡回行列作成工程：差スペクトルデータを対称に折り返して、そのデータのフーリエ変換により得られた対称インタフェログラムデータを第１列のデータとする巡回行列の作成。
▲４▼特異値分解工程：上記巡回行列に対して特異値分解（ＳＶＤ）を適用。この工程によりノイズが除去された行列を得る。
▲５▼演算差スペクトル再生工程：上記ノイズが除去された行列からの前記差スペクトルに相当するノイズが除去されたスペクトルの再生。
【００１２】
▲６▼合成工程：平滑化スペクトルと演算差スペクトルの和合成を行い、ノイズ除去された連続スペクトルを得る。
以下、各工程を説明する。
▲１▼平滑化工程
この平滑化工程においては、たとえば赤外吸収スペクトルとして採取されたｎ／２点のスペクトルデータに対して幅Ｗ点の平滑化フィルタを適用する。平滑化フィルタの適用に関しては、例えばサビツキ−ゴレイフィルタを適用することができる。なお、この平滑化フィルタを適用することにより、元連続スペクトル中の意味のある急峻なピークも除かれてしまうことがある。
【００１３】
▲２▼差スペクトル演算工程
この差スペクトル演算工程においては、平滑化工程において得られた平滑化スペクトルと測定スペクトルとの差スペクトルを演算する。この結果、測定スペクトルに存在し、かつ平滑化により消失した比較的急峻なピークがノイズとともに差スペクトルとして採取されることとなる。
【００１４】
一方、以下の巡回行列作成工程〜演算差スペクトル再生工程のみからなる特異値分解によるノイズ除去は、急峻なピークを有するスペクトルにおいて特に好適であるが、前記赤外スペクトルなどの連続スペクトルに適用することは困難であった。そこで、前記連続スペクトルから急峻なピーク部分（差スペクトル）のみを抽出し、これに対して前記特異値分解によるノイズ除去を適用することにより、特異値分解による効率的かつ適正なノイズ除去を行うことができるのである。
【００１５】
▲３▼巡回行列作成工程
差スペクトルのサンプルポイントの数（巡回行列の行）は、次段階で要求される特異値分解（ＳＶＤ）工程の計算時間を削減するため、できるだけ小さくすべきである。このため、全差スペクトルの内の小さい特定部分を抽出する（図１（ａ））。この部分スペクトル（ｎ／２点）は、対称に折返され、ｎ点のサンプルのインターフェログラムを得るためにフーリエ変換される（図１（ｂ））。なお、スペクトルの折返しを行なっているため、フーリエ変換の虚数部分は零となる対称な虚数部分のインターフェログラムデータのみを考慮すれば良い。
このインターフェログラムデータを用いて、ｎ×ｐの行列Ｘを構成する（図１（ｃ））。この行列Ｘは、各列ベクトルが一つ前の列に対して一つづつ下方にずれて巡回するように構成されている。図において｛ｘk｝（ｋ＝１，２，…ｎ）がインターフェログラムのサンプルバージョンである。行列Ｘの列の数ｐ（＜ｎ））は適宜選択される第１パラメータである。
【００１６】
▲４▼ＳＶＤ工程
図２に示されるように、ｎ×ｐの巡回行列をＸ＝ＵＤＶt＝ＴＶtに再構成する。ここで、Ｕ及びＶはそれぞれｎ×ｎ及びｐ×ｐのユニタリー行列であり、Ｄはｎ×ｐの対角行列である。この対角行列の対角成分は行列Ｘの特異値である。行列Ｔ及びＶは多変量解析の分野でそれぞれスコア行列及びローディング行列と呼ばれている。これらの特異値に対応させて、行列Ｘを、階数１のＬ個の行列Ｘi（ｉ＝１，２，…，Ｌ）及びＸi（ｉ＝Ｌ＋１，…，ｐ）に分解する。そして、最初の行列を統合して行列Ｓと残りを行列Ｅとする。最初の差スペクトルの信号成分は行列Ｓに濃縮され、ノイズ成分は行列Ｅに濃縮される。従って、当初の行列Ｘの代りに行列Ｓを以後の計算に用いれば、ノイズ成分を効率的に除去することができる。Ｌの値は、行列Ｘの階数と同じである。最初の差スペクトルがＭ本の急峻なピークからなるときには、Ｌ≧２Ｍとする。これは１本のピークは、行列Ｘにおいて階数２（インターフェログラム中のサイン、コサイン項に相当する）となることによる。本実施例においては、Ｌの値は適当なスレッショルドレベルを設けて特異値分析を行なうことにより決定している。このＬは適宜選択される第２のパラメータである。
【００１７】
本実施例においては、非線形反復部分最小二乗法（ＮＩＰＡＬＳ）法を、ＳＶＤの実施のために、採用している。これはコンピュータによる計算に適合しやすいためである。このＮＩＰＡＬＳは、すべての特異値を一度に計算するのではなく、最も大きい値から最も小さい値へ、順次計算を行なうものである。このため、対応する特異値を選択することにより、すべての差スペクトルから特定のピークのみを選択することも可能となる。
【００１８】
▲５▼演算差スペクトル再生工程（差スペクトルの再構成）
図１に示す工程とちょうど逆の工程をたどって、ｎ×ｐの行列Ｓから、ノイズが除去されたインターフェログラムを得る。行列Ｓのｉ番目（２≦ｉ≦ｐ）の列ベクトルは、上方にｉ個、巡回的にシフトする。そして、行列の各行ベクトルの平均の計算により、ｎサンプル点のインターフェログラムを得ることができる。このようにして得た、インターフェログラムデータに対して逆フーリエ変換を施すことによりノイズが除去された差スペクトル（演算差スペクトル）を再生することができる。
【００１９】
［第２実施例］
本発明において平滑化スペクトルと差スペクトルの各々を等波長又は等波数間隔になるようにＭ分割し、Ｍ分割された各差スペクトル部に対して等しい値のｐパラメータ値で前記特異値分解法を適用することが好適である。
すなわち、図３に示すように、前記第１実施例と同じステップ１００，１０２，１０４を介して差スペクトルを得た後、該差スペクトルをＭ分割し（ステップ１０６）、各部分差スペクトルに対して特異値分解及び各部分差スペクトルの再生、部分再生スペクトルの合成を行うのである（ステップ１０８，１１０）。そして、再生部分差スペクトルはそれぞれ対応する部分平滑化スペクトルと合成され、さらにスペクトルの再結合が行われる（ステップ１１２）。
【００２０】
ここで、分割された各領域のＬパラメータは図４に示すように選択されることが好適である。
すなわち、特異値曲線▲１▼〜▲４▼は、以下で示す図５（Ｂ）のスペクトルデータに対してＭ＝４の場合の各々にＳＶＤを適用した場合の、特異値ナンバー対特異値プロットを示す。
【００２１】
図４（Ａ）において、
▲３▼はＬ＝Ｐ（ノイズありのスペクトル、つまり平滑化をする前のスペクトル）▲４▼はＬ＝０（平滑化スペクトルそのもの）
▲１▼、▲２▼はＬ＝Ｌ1，Ｌ＝Ｌ2
また、図４（Ｂ）においては、▲１▼、▲２▼、▲３▼及び▲４▼に対し、Ｌ＝Ｌ1，Ｌ2，Ｌ3，Ｌ4
なお、上記図４において、閾値レベルはＬｉ≧２ｌ（ここでｌはスペクトルを特徴づける本数）のように定めることができる。
【００２２】
以上のようにして得られたＭ個の部分差スペクトルと、Ｍ個の部分平滑化スペクトルの各々を加算し、得られたＭ個の和スペクトルを再度接続して一本の連続スペクトルとする。この分割処理によって計算時間が１／Ｍ程度に短縮可能であり、しかも正確にノイズ除去が行われたスペクトルを得ることができる。
次に、図５〜図９に基づき、本発明にかかる方法を、ＦＴ−ＩＲ法により得られたポリスチレンの実際のスペクトルデータに適用した場合について説明する。
【００２３】
まず、図５（ａ）は、ポリスチレンのＦＴ−ＩＲスペクトル（分解能４cm-1）の１７０６．７〜１２１４．９cm-1領域を拡大したものである。データ点数は２５７点で、これを真のスペクトルとする。
そして、この図５（ａ）の真のスペクトルに対し、５％のガウスノイズを重畳させたものを図５（ｂ）に示しており、これを処理対象のスペクトルとする。
この処理対象スペクトルに対し、１３点の窓幅でサビツキ−ゴレイの平滑化を行ったものを図６（ａ）に示す。この図６（ａ）と図５（ｂ）を比較すると、平滑化処理によりポリスチレンが本来有するピークも消失することが理解される。
【００２４】
一方、この図５（ｂ）の処理対象スペクトルから図６（ａ）の平滑化スペクトルを減算したものを図６（ｂ）に示す。この図６（ｂ）に示す差スペクトルは、ＦＴ−ＩＲ特有の連続スペクトルの特性が消失し、輝線スペクトルに類似した特性を有する。
この図６（ｂ）に示す差スペクトルに対し、特異値分解法を利用したノイズ除去をｐ＝１００，Ｌ＝２０というパラメータ設定で行った結果得られるノイズ除去（演算）差スペクトルを、図６（ｃ）に示す。
【００２５】
この演算差スペクトルを前記平滑化スペクトルに加えた結果（再生スペクトル）を図６（ｄ）に示す。該再生スペクトルを前記図５（ａ）に示す真スペクトル及び図６（ａ）に示す平滑化スペクトルと比較すると、平滑化スペクトルよりもはるかに真スペクトルに近いことが示される。なお、図７には特異値成分数（Ｌ）と特異値の関係と、特異値成分数（Ｌ）と真のスペクトル（図５（ａ））とノイズ除去された結果との残差二乗和をプロットしたものが示されている。同図より明らかなように残差二乗和はＬ＝２０で極小値をとる。
【００２６】
これに対し、図８（ａ）は図５（ｂ）を１００％ラインから引いた後、一次のベースラインを取り除いたものに従来の特異値分解法を利用したノイズ除去法をｐ＝１００，Ｌ＝７０というパラメータ設定で適用した結果を示している。同図にはスパイクノイズの残存が認められ、さらに演算時間もＬ＝７０と大きな値のため長時間にわたる。
また、図８（ｂ）には、図５（ｂ）に局所分散量を考慮した適応化平滑化処理を局所分散窓幅１３点で適用した結果が示されている。このときのノイズの分散はスペクトル中１２２０cm-1付近の２１点の分散値を使用した。この場合には演算時間の短縮ははかられるが、前記図８（ａ）の場合と同じくスパイクノイズの残存が認められる。
【００２７】
一方、図９（ａ）には図５（ｂ）を４分割し、その第一エリアに提案手法をｐ＝２５，Ｌ＝６というパラメータ設定で適用したものを示す。平滑化点数は１３点、データ点数は６４点である。
同じく図９（ｂ）は第二エリアに提案手法をＰ＝２５、Ｌ＝１０というパラメータ設定で適用したものである。平滑化点数は１３点、データ点数は６４点である。
図９（ｃ）は第三エリアに提案手法をＰ＝２５，Ｌ＝４というパラメータ設定で適用したものである。平滑化点数は１３点、データ点数は６５点である。
【００２８】
図９（ｄ）は第四エリアに提案手法をＰ＝２５，Ｌ＝０というパラメータ設定で適用したものである。平滑化点数は１３点、データ点数は６５点である。
以上の図９（ａ）〜図９（ｄ）をつなぎ合わせたものを図９（ｅ）に示す。前記図６（ｄ）の再生スペクトルと全く同一の結果が得られることが理解される。そして、本実施例のように４分割してノイズ除去を行った場合、演算に要する時間は約１／４になる。
図９（ｆ）には、図９（ａ）〜（ｄ）のそれぞれの特異値成分数（Ｌ）と特異値の関係をプロットしたものを示す。なお、図中の点線は図８（ｄ）の結果を得たときに採用した特異値である。
【００２９】
【発明の効果】
以上説明したように本発明にかかる方法によれば、測定スペクトルと平滑化スペクトルの差スペクトルのインタフェログラムデータから成る巡回行列にＳＶＤ法を適用することにより、効率的に且つデジタル的に連続スペクトルデータからのノイズ除去を行なうことが可能となる。
【図面の簡単な説明】
【図１】本発明の巡回行列作成工程の説明図である。
【図２】本発明の特異値分解工程の説明図である。
【図３】本発明の第２実施例にかかるノイズ除去方法の説明図である。
【図４】本発明の第２実施例にかかる特異値と特異指数の関係の説明図である。
【図５】同図（ａ）はポリスチレンのＦＴ−ＩＲスペクトルの説明図、同図（ｂ）は該スペクトルにガウスノイズを重畳させて処理対象スペクトルとしたものである。
【図６】図５に示した処理対象スペクトルに本発明の第１実施例にかかるノイズ除去方法を適用した過程を示す説明図である。
【図７】本発明の第１実施例にかかる特異値成分数と特異値の関係と、特異値成分数と真スペクトルとノイズ除去された結果との残差二乗和をプロットした説明図である。
【図８】前記図５に示す処理対象スペクトルに従来の特異値分解法を利用したノイズ除去方法を適用した場合の説明図（ａ）と従来の適応化平滑化処理によるノイズ除去法を適用した場合の説明図（ｂ）である。
【図９】前記図５に示す処理対象スペクトルに本発明の第２実施例にかかるノイズ除去方法を適用した過程を示す説明図である。
【図１】
0003541083.tif 000002
【図２】
0003541083.tif 000003
【図３】
0003541083.tif 000004
【図４】
0003541083.tif 000005
【図５】
0003541083.tif 000006
【図６】
0003541083.tif 000007
【図７】
0003541083.tif 000008
【図８】
0003541083.tif 000009
【図９】
0003541083.tif 000010