(57)【要約】  情報処理装置（２０００）は、センサ（１０）によって出力された時系列データ（１４）を取得し、複数の特徴定数θi と、各特徴定数θi の時系列データ（１４）に対する寄与を表す寄与値ξi を算出する。そして情報処理装置（２０００）は、特徴定数θi の集合Θと、寄与値ξi の集合Ξとを対応づけた情報を、対象ガスの特徴量として出力する。特徴定数θとしては、速度定数βや、速度定数の逆数である時定数τを採用することができる。
2020026327.tif 000002
【特許請求の範囲】
【請求項１】
対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得部と、
前記時系列データに対して寄与する複数の特徴定数と、各前記特徴定数の前記時系列データに対する寄与の大きさを表す寄与値とを算出する算出部と、
複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、前記センサによってセンシングされたガスの特徴量として出力する出力部と、を有し、
前記特徴定数は、前記センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である、情報処理装置。
【請求項２】
前記算出部は、
前記時系列データの測定期間から複数の部分期間を抽出し、
前記部分期間ごとに、その部分期間における前記検出値の時間変化率の対数に基づいて、前記特徴定数を算出し、
前記部分期間は、その中に含まれる検出値の時間変化率の対数が略同一の期間である、請求項１に記載の情報処理装置。
【請求項３】
前記算出部は、
前記時系列データを用いて、各時刻における検出値及びその時刻における検出値の時間変化率を要素とする時系列のベクトルデータを算出し、
各前記算出した時系列のベクトルデータそれぞれについて速度ベクトルを算出し、
前記速度ベクトルの向きに基づいて、前記時系列データの測定期間から、複数の部分期間を抽出し、
前記部分期間ごとに、その部分期間における前記速度ベクトルの向きに基づいて、前記特徴定数を算出し、
前記部分期間は、その中に含まれる前記速度ベクトルの向きが略同一の期間である、請求項１に記載の情報処理装置。
【請求項４】
前記算出部は、複数の特徴定数それぞれの寄与値をパラメータとする前記センサの検出値の予測モデルについて、前記取得した時系列データを用いたパラメータ推定を行うことで、各寄与値を算出する、請求項１乃至３いずれか一項に記載の情報処理装置。
【請求項５】
前記算出部は、前記予測モデルから得られる時系列データと前記取得した時系列データとについて、最小二乗法による最尤推定を行うことで、各前記寄与値を算出する、請求項４に記載の情報処理装置。
【請求項６】
前記最小二乗法における最尤推定において、目的関数に正則化項が含まれている、請求項５に記載の情報処理装置。
【請求項７】
前記算出部は、各前記寄与値の事前分布と前記取得した時系列データとを用いた MAP（Maximum a Posteriori）推定又はベイズ推定により、各前記寄与値を算出する、請求項４に記載の情報処理装置。
【請求項８】
前記事前分布は、多変量正規分布又はガウス過程である、請求項７に記載の情報処理装置。
【請求項９】
前記算出部は、前記予測モデルから得られる時系列データと前記取得した時系列データとの二乗誤差を表す目的関数について、その目的関数の最小値を複数の特徴定数に関して最小化することにより、複数の特徴定数及び複数の寄与値を算出する、請求項４に記載の情報処理装置。
【請求項１０】
前記予測モデルには、バイアスを表すパラメータが含まれており、
前記算出部は、前記予測モデルについて寄与値とバイアスそれぞれを表すパラメータを推定する、請求項４乃至９いずれか一項に記載の情報処理装置。
【請求項１１】
前記時系列データ取得部は、複数の時系列データを取得し、
前記算出部は、複数の時系列データそれぞれについて、特徴定数の集合と寄与値の集合の組みを算出し、
前記出力部は、前記算出された特徴定数の集合と寄与値の集合の組みを複数まとめた情報を、前記対象ガスの特徴量として出力する、請求項１乃至１０いずれか一項に記載の情報処理装置。
【請求項１２】
前記複数の時系列データは、前記センサを前記対象ガスに曝す際に得られる時系列データと、前記センサから前記対象ガスを取り除く際に得られる時系列データとの双方を含む、請求項１１に記載の情報処理装置。
【請求項１３】
前記複数の時系列データは、特性の異なる複数の前記センサそれぞれから得られる時系列データを含む、請求項１１に記載の情報処理装置。
【請求項１４】
コンピュータによって実行させる制御方法であって、
対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得ステップと、
前記時系列データに対して寄与する複数の特徴定数と、各前記特徴定数の前記時系列データに対する寄与の大きさを表す寄与値とを算出する算出ステップと、
複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、前記センサによってセンシングされたガスの特徴量として出力する出力ステップと、を有し、
前記特徴定数は、前記センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である、制御方法。
【請求項１５】
前記算出ステップにおいて、
前記時系列データの測定期間から複数の部分期間を抽出し、
前記部分期間ごとに、その部分期間における前記検出値の時間変化率の対数に基づいて、前記特徴定数を算出し、
前記部分期間は、その中に含まれる検出値の時間変化率の対数が略同一の期間である、請求項１４に記載の制御方法。
【請求項１６】
前記算出ステップにおいて、
前記時系列データを用いて、各時刻における検出値及びその時刻における検出値の時間変化率を要素とする時系列のベクトルデータを算出し、
各前記算出した時系列のベクトルデータそれぞれについて速度ベクトルを算出し、
前記速度ベクトルの向きに基づいて、前記時系列データの測定期間から、複数の部分期間を抽出し、
前記部分期間ごとに、その部分期間における前記速度ベクトルの向きに基づいて、前記特徴定数を算出し、
前記部分期間は、その中に含まれる前記速度ベクトルの向きが略同一の期間である、請求項１４に記載の制御方法。
【請求項１７】
前記算出ステップにおいて、複数の特徴定数それぞれの寄与値をパラメータとする前記センサの検出値の予測モデルについて、前記取得した時系列データを用いたパラメータ推定を行うことで、各寄与値を算出する、請求項１４乃至１６いずれか一項に記載の制御方法。
【請求項１８】
前記算出ステップにおいて、前記予測モデルから得られる時系列データと前記取得した時系列データとについて、最小二乗法による最尤推定を行うことで、各前記寄与値を算出する、請求項１７に記載の制御方法。
【請求項１９】
前記最小二乗法における最尤推定において、目的関数に正則化項が含まれている、請求項１８に記載の制御方法。
【請求項２０】
前記算出ステップにおいて、各前記寄与値の事前分布と前記取得した時系列データとを用いた MAP（Maximum a Posteriori）推定又はベイズ推定により、各前記寄与値を算出する、請求項１７に記載の制御方法。
【請求項２１】
前記事前分布は、多変量正規分布又はガウス過程である、請求項２０に記載の制御方法。
【請求項２２】
前記算出ステップにおいて、前記予測モデルから得られる時系列データと前記取得した時系列データとの二乗誤差を表す目的関数について、その目的関数の最小値を複数の特徴定数に関して最小化することにより、複数の特徴定数及び複数の寄与値を算出する、請求項１７に記載の制御方法。
【請求項２３】
前記予測モデルには、バイアスを表すパラメータが含まれており、
前記算出ステップにおいて、前記予測モデルについて寄与値とバイアスそれぞれを表すパラメータを推定する、請求項１７乃至２２いずれか一項に記載の制御方法。
【請求項２４】
前記時系列データ取得ステップにおいて、複数の時系列データを取得し、
前記算出ステップにおいて、複数の時系列データそれぞれについて、特徴定数の集合と寄与値の集合の組みを算出し、
前記出力ステップにおいて、前記算出された特徴定数の集合と寄与値の集合の組みを複数まとめた情報を、前記対象ガスの特徴量として出力する、請求項１４乃至２３いずれか一項に記載の制御方法。
【請求項２５】
前記複数の時系列データは、前記センサを前記対象ガスに曝す際に得られる時系列データと、前記センサから前記対象ガスを取り除く際に得られる時系列データとの双方を含む、請求項２４に記載の制御方法。
【請求項２６】
前記複数の時系列データは、特性の異なる複数の前記センサそれぞれから得られる時系列データを含む、請求項２４に記載の制御方法。
【請求項２７】
請求項１４乃至２６いずれか一項に記載の制御方法の各ステップをコンピュータに実行させるプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明はガスの特徴の解析に関する。
【背景技術】
【０００２】
ガスをセンサで測定することにより、ガスに関する情報を得る技術が開発されている。特許文献１は、ナノメカニカルセンサで試料ガスを測定することで得られるシグナル（検出値の時系列データ）を利用して、試料ガスの種類を判別する技術を開示している。具体的には、センサの受容体に対する試料ガスの拡散時定数が、受容体の種類と試料ガスの種類の組み合わせによって決まるため、シグナルから得られる拡散時定数と、受容体の種類とに基づいて、試料ガスの種類を判別できることが開示されている。
【先行技術文献】
【特許文献】
【０００３】
【特許文献１】特開２０１７−１５６２５４号公報
【発明の概要】
【発明が解決しようとする課題】
【０００４】
特許文献１では、試料ガスに含まれている分子が１種類であることが前提となっており、複数種類の分子が混合している試料ガスを扱うことが想定されていない。本願発明は上記の課題に鑑みてなされたものであり、複数種類の分子が混合しているガスの特徴を抽出する技術を提供することである。
【課題を解決するための手段】
【０００５】
本発明の情報処理装置は、１）対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得部と、２）時系列データに対して寄与する複数の特徴定数と、各特徴定数の時系列データに対する寄与の大きさを表す寄与値とを算出する算出部と、３）複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、センサによってセンシングされたガスの特徴量として出力する出力部と、を有する。特徴定数は、センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である。【０００６】
本発明の制御方法は、コンピュータによって実行される。当該制御方法は、１）対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得ステップと、２）時系列データに対して寄与する複数の特徴定数と、各特徴定数の時系列データに対する寄与の大きさを表す寄与値とを算出する算出ステップと、３）複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、センサによってセンシングされたガスの特徴量として出力する出力ステップと、を有する。特徴定数は、センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である。【０００７】
本発明のプログラムは、コンピュータに、本発明の制御方法が有する各ステップを実行させる。
【発明の効果】
【０００８】
本発明によれば、複数種類の分子が混合しているガスの特徴を抽出する技術が提供される。
【図面の簡単な説明】
【０００９】
上述した目的、およびその他の目的、特徴および利点は、以下に述べる好適な実施の形態、およびそれに付随する以下の図面によってさらに明らかになる。
【図１】実施形態１の情報処理装置の概要を例示する図である。
【図２】情報処理装置が取得するデータを得るためのセンサを例示する図である。
【図３】実施形態１の情報処理装置の機能構成を例示する図である。
【図４】情報処理装置を実現するための計算機を例示する図である。
【図５】実施形態１の情報処理装置によって実行される処理の流れを例示するフローチャートである。
【図６】センサから得られる複数の時系列データを例示する図である。
【図７】g(t) のグラフを例示する図である。
【図８】(y(t), y'(t)) の速度ベクトルを例示する図である。
【図９】AIC をグラフで表す図である。
【図１０】対象ガスの特徴量をグラフで例示する図である。
【図１１】立ち上がりの時系列データと立ち下がりの時系列データそれぞれから特徴行列を得るケースを例示する図である。
【図１２】複数のセンサそれぞれから時系列データを得ることで、複数の特徴行列を得るケースを例示する図である。
【発明を実施するための形態】
【００１０】
以下、本発明の実施の形態について、図面を用いて説明する。尚、すべての図面において、同様な構成要素には同様の符号を付し、適宜説明を省略する。また、特に説明する場合を除き、各ブロック図において、各ブロックは、ハードウエア単位の構成ではなく、機能単位の構成を表している。
【００１１】
［実施形態１］
＜発明の概要と理論的背景＞
図１は、実施形態１の情報処理装置２０００の概要を例示する図である。また、図２は、情報処理装置２０００が取得するデータを得るためのセンサ１０を例示する図である。センサ１０は、分子が付着する受容体を有し、その受容体における分子の付着と離脱に応じて検出値が変化するセンサである。なお、センサ１０によってセンシングされているガスを、対象ガスと呼ぶ。また、センサ１０から出力される検出値の時系列データを、時系列データ１４と呼ぶ。ここで、必要に応じ、時系列データ１４を Y とも表記し、時刻 t の検出値を y(t) とも表記する。Y は、y(t) が列挙されたベクトルとなる。
【００１２】
例えばセンサ１０は、膜型表面応力（Membrane-type Surface Stress; MSS)センサである。MSS センサは、受容体として、分子が付着する官能膜を有しており、その官能膜に対する分子の付着と離脱によってその官能膜の支持部材に生じる応力が変化する。MSS センサは、この応力の変化に基づく検出値を出力する。なお、センサ１０は、MSS センサには限定されず、受容体に対する分子の付着と離脱に応じて生じる、センサ１０の部材の粘弾性や動力学特性(質量や慣性モーメントなど）に関連する物理量の変化に基づいて検出値を出力するものであればよく、カンチレバー式、膜型、光学式、ピエゾ、振動応答などの様々なタイプのセンサを採用することができる。
【００１３】
ここで、説明のため、センサ１０によるセンシングを以下のようにモデル化する。
（１）センサ１０は、K 種類の分子を含む対象ガスに曝されている。
（２）対象ガスにおける各分子 k の濃度は一定のρkである。
（３）センサ１０には、合計 N 個の分子が吸着可能である。
（４）時刻t においてセンサ１０に付着している分子k の数は nk(t) 個である。
【００１４】
センサ１０に付着している分子 k の数 nk(t) の時間変化は、以下のように定式化できる。
【数１】
2020026327.tif 000003
【００１５】
式（１）の右辺の第１項と第２項はそれぞれ、単位時間当たりの分子 k の増加量（新たにセンサ１０に付着する分子 k の数）と減少量（センサ１０から離脱する分子 k の数）を表している。また、αk とβk はそれぞれ、分子 k がセンサ１０に付着する速度を表す速度定数と、分子 k がセンサ１０から離脱する速度を表す速度定数である。
【００１６】
ここで、濃度ρkが一定であるため、上記式（１）から、時刻t における分子 k の数 nk(t) は、以下のように定式化できる。
【数２】
2020026327.tif 000004
【００１７】
また、時刻 t0（初期状態）でセンサ１０に分子が付着していないと仮定すれば、nk(t) は以下のように表される。
【数３】
2020026327.tif 000005
【００１８】
センサ１０の検出値は、対象ガスに含まれる分子によってセンサ１０に働く応力によって定まる。そして、複数の分子によってセンサ１０に働く応力は、個々の分子に働く応力の線形和で表すことができると考えられる。ただし、分子によって生じる応力は、分子の種類によって異なると考えられる。すなわち、センサ１０の検出値に対する分子の寄与は、その分子の種類によって異なると言える。
【００１９】
そこで、センサ１０の検出値 y（t） は、以下のように定式化できる。
【数４】
2020026327.tif 000006
ここで、γk とξk はいずれも、センサ１０の検出値に対する分子 k の寄与を表す。なお、「立ち上がり」と「立ち下がり」の意味については、後述する。
【００２０】
ここで、対象ガスをセンシングしたセンサ１０から得た時系列データ１４を上述の式（４）のように分解できれば、対象ガスに含まれる分子の種類や、各種類の分子が対象ガスに含まれる割合を把握することができる。すなわち、式（４）に示す分解によって、対象ガスの特徴を表すデータ（すなわち、対象ガスの特徴量）が得られる。
【００２１】
そこで情報処理装置２０００は、センサ１０によって出力された時系列データ１４を取得し、時系列データ１４を以下の式（５）に示すように分解する。
【数５】
2020026327.tif 000007
ここで、ξi は、センサ１０の検出値に対する特徴定数θi の寄与を表す寄与値である。
【００２２】
具体的には、まず情報処理装置２０００は、時系列データ１４を用いて、時系列データ１４に対して寄与する複数の特徴定数の集合Θ={θ1,..., θm} と、各特徴定数θi の時系列データ１４に対する寄与の大きさを表す寄与値ξi を算出する。なお、後述するように、特徴定数θi を算出してから寄与値ξi が算出されるケースと、特徴定数θi と共に寄与値ξi が算出されるケースとがある。
【００２３】
さらに情報処理装置２０００は、特徴定数の集合Θと寄与値の集合Ξとを対応づけた情報を、対象ガスの特徴を表す特徴量として出力する。特徴定数の集合Θと寄与値の集合Ξとの対応づけは、例えば、m 行２列の特徴行列 F で表される（m は特徴定数と寄与値それぞれの数）。例えばこの行列 F は、特徴定数の集合を表す特徴定数ベクトルΘ=(θ1,..., θm) を第１列に有し、なおかつ寄与値の集合を表す寄与ベクトルΞ=(ξ1,..., ξm) を第２列に有する。すなわち、F=(ΘT, ΞT) である。以降の説明では、特に断らない限り、対象ガスの特徴量は、特徴行列 F で表されるとする。ただし、対象ガスの特徴量は、必ずベクトルとして表現しなければならないわけではない。
【００２４】
ここで、特徴定数θとしては、前述した速度定数βや、速度定数の逆数である時定数τを採用することができる。θとしてβとτを使う場合それぞれについて、式（5）は、以下のように表すことができる。
【００２５】
【数６】
2020026327.tif 000008
【００２６】
なお、図１では、特徴定数として速度定数βが利用されている。そのため、特徴定数の集合Θは、β1 からβm を示している。
【００２７】
＜作用・効果＞
前述したように、センサ１０の検出値に対する分子の寄与は、その分子の種類によって異なると考えられるため、上述した特徴定数の集合Θとそれに対応する寄与値の集合Ξは、対象ガスに含まれる分子の種類やその混合比率に応じて異なるものになると考えられる。よって、特徴定数の集合Θと寄与値の集合Ξとを対応づけた情報は、複数種類の分子が混合されているガスを互いに区別することができる情報、すなわちガスの特徴量として利用することができる。
【００２８】
そこで本実施形態の情報処理装置２０００は、対象ガスをセンサ１０でセンシングすることで得られた時系列データ１４に基づいて、特徴定数の集合Θと、各特徴定数の時系列データ１４に対する寄与を表す寄与値の集合Ξとを算出し、算出した集合ΘとΞとを対応づけた情報を、対象ガスの特徴量として出力する。こうすることで、複数種類の分子が混合しているガスを識別することが可能な特徴量を、そのガスをセンサ１０でセンシングした結果から自動的に生成することができる。
【００２９】
なお、図１を参照した上述の説明は、情報処理装置２０００の理解を容易にするための例示であり、情報処理装置２０００の機能を限定するものではない。以下、本実施形態の情報処理装置２０００についてさらに詳細に説明する。
【００３０】
＜情報処理装置２０００の機能構成の例＞
図３は、実施形態１の情報処理装置２０００の機能構成を例示する図である。情報処理装置２０００は、時系列データ取得部２０２０、算出部２０４０、及び出力部２０６０を有する。時系列データ取得部２０２０は、センサ１０から時系列データ１４を取得する。特徴定数生成部２０３０は、時系列データ１４を用いて、複数の特徴定数θi と、各特徴定数θi の時系列データ１４に対する寄与の大きさを表す寄与値ξi を算出する。出力部２０６０は、特徴定数の集合Θと寄与値の集合Ξとを対応づけた情報を、センサ１０によってセンシングされたガスの特徴量として出力する。
【００３１】
＜情報処理装置２０００のハードウエア構成＞
情報処理装置２０００の各機能構成部は、各機能構成部を実現するハードウエア（例：ハードワイヤードされた電子回路など）で実現されてもよいし、ハードウエアとソフトウエアとの組み合わせ（例：電子回路とそれを制御するプログラムの組み合わせなど）で実現されてもよい。以下、情報処理装置２０００の各機能構成部がハードウエアとソフトウエアとの組み合わせで実現される場合について、さらに説明する。
【００３２】
図４は、情報処理装置２０００を実現するための計算機１０００を例示する図である。計算機１０００は任意の計算機である。例えば計算機１０００は、Personal Computer（PC）やサーバマシンなどの据え置き型の計算機である。その他にも例えば、計算機１０００は、スマートフォンやタブレット端末などの可搬型の計算機である。計算機１０００は、情報処理装置２０００を実現するために設計された専用の計算機であってもよいし、汎用の計算機であってもよい。
【００３３】
計算機１０００は、バス１０２０、プロセッサ１０４０、メモリ１０６０、ストレージデバイス１０８０、入出力インタフェース１１００、及びネットワークインタフェース１１２０を有する。バス１０２０は、プロセッサ１０４０、メモリ１０６０、ストレージデバイス１０８０、入出力インタフェース１１００、及びネットワークインタフェース１１２０が、相互にデータを送受信するためのデータ伝送路である。ただし、プロセッサ１０４０などを互いに接続する方法は、バス接続に限定されない。
【００３４】
プロセッサ１０４０は、CPU（Central Processing Unit）、GPU（Graphics Processing Unit）、FPGA（Field−Programmable Gate Array）などの種々のプロセッサである。メモリ１０６０は、RAM（Random Access Memory）などを用いて実現される主記憶装置である。ストレージデバイス１０８０は、ハードディスク、SSD（Solid State Drive）、メモリカード、又は ROM（Read Only Memory）などを用いて実現される補助記憶装置である。
【００３５】
入出力インタフェース１１００は、計算機１０００と入出力デバイスとを接続するためのインタフェースである。例えば入出力インタフェース１１００には、キーボードなどの入力装置や、ディスプレイ装置などの出力装置が接続される。その他にも例えば、入出力インタフェース１１００には、センサ１０が接続される。ただし、センサ１０は必ずしも計算機１０００と直接接続されている必要はない。例えばセンサ１０は、計算機１０００と共有している記憶装置に時系列データ１４を記憶させてもよい。
【００３６】
ネットワークインタフェース１１２０は、計算機１０００を通信網に接続するためのインタフェースである。この通信網は、例えば LAN（Local Area Network）や WAN（Wide Area Network）である。ネットワークインタフェース１１２０が通信網に接続する方法は、無線接続であってもよいし、有線接続であってもよい。
【００３７】
ストレージデバイス１０８０は、情報処理装置２０００の各機能構成部を実現するプログラムモジュールを記憶している。プロセッサ１０４０は、これら各プログラムモジュールをメモリ１０６０に読み出して実行することで、各プログラムモジュールに対応する機能を実現する。
【００３８】
＜処理の流れ＞
図５は、実施形態１の情報処理装置２０００によって実行される処理の流れを例示するフローチャートである。時系列データ取得部２０２０は時系列データ１４を取得する（Ｓ１０２）。算出部２０４０は、時系列データ１４を用いて、複数の特徴定数θi と、各特徴定数θi に対応する寄与値ξiとを算出する（Ｓ１０４）。出力部２０６０は、特徴定数の集合と寄与値の集合とを対応づけた情報を、対象ガスの特徴量として出力する（Ｓ１０６）。
【００３９】
情報処理装置２０００が図５に示す一連の処理を実行するタイミングは様々である。例えば、情報処理装置２０００は、時系列データ１４を指定する入力操作を受け付け、指定された時系列データ１４について一連の処理を実行する。その他にも例えば、情報処理装置２０００は、時系列データ１４を受信できるように待機しておき、時系列データ１４を受信したこと（すなわち、Ｓ１０２を実行したこと）に応じて、Ｓ１０４以降の処理を実行する。
【００４０】
＜時系列データ１４の取得：Ｓ１０２＞
時系列データ取得部２０２０は時系列データ１４を取得する（Ｓ１０２）。時系列データ取得部２０２０が時系列データ１４を取得する方法は任意である。例えば情報処理装置２０００は、時系列データ１４が記憶されている記憶装置にアクセスすることで、時系列データ１４を取得する。時系列データ１４が記憶されている記憶装置は、センサ１０の内部に設けられていてもよいし、センサ１０の外部に設けられていてもよい。その他にも例えば、時系列データ取得部２０２０は、センサ１０から出力される検出値を順次受信することで、時系列データ１４を得てもよい。
【００４１】
時系列データ１４は、センサ１０が出力した検出値を、センサ１０から出力された時刻が早い順に並べた時系列のデータである。ただし、時系列データ１４は、センサ１０から得られた検出値の時系列データに対して、所定の前処理が加えられたものであってもよい。また、前処理が行われた時系列データ１４を取得する代わりに、時系列データ取得部２０２０が時系列データ１４に対して前処理を行ってもよい。前処理としては、例えば、時系列のデータからノイズ成分を除去するフィルタリングなどを採用することができる。
【００４２】
ここで、時系列データ１４は、センサ１０を対象ガスに曝すことで得られる。ただし、センサを用いてガスに関する測定を行う場合、センサを測定対象のガスに曝す操作と、センサから測定対象のガスを取り除く操作を繰り返すことで、センサから解析対象の時系列データを複数得ることがある。
【００４３】
図６は、センサから得られる複数の時系列データを例示する図である。図６では、立ち上がりの時系列データと立ち下がりの時系列データとを区別しやすいように、立ち上がりの時系列データが実線で表されており、立ち下がりの時系列データが点線で表されている。図６において、期間Ｐ１の時系列データ１４−１と期間Ｐ３の時系列データ１４−３は、センサを測定対象のガスに曝す操作によって得られる。このように、センサを測定対象のガスに曝す操作を行うと、センサの検出値が増加する。そこで、センサを測定対象のガスに曝す操作によって得られる時系列データを、「立ち上がり」の時系列データとも呼ぶ。式（４）における「立ち上がりの場合」は、「時系列データ１４が立ち上がりの時系列データである場合」を意味する。以降の式においても同様である。
【００４４】
一方、期間Ｐ２の時系列データ１４−２と期間Ｐ４の時系列データ１４−４は、センサから測定対象のガスを取り除く操作によって得られる。なお、センサから測定対象のガスを取り除く操作は、例えば、センサをパージガスと呼ばれるガスに曝すことで実現される。このように、センサから測定対象のガスを取り除く操作を行うと、センサの検出値が減少する。そこで、センサから測定対象のガスを取り除く操作によって得られる時系列データを、「立ち下がり」の時系列データとも呼ぶ。式（４）における「立ち下がりの場合」は、「時系列データ１４が立ち下がりの時系列データである場合」を意味する。以降の式においても同様である。
【００４５】
情報処理装置２０００では、センサ１０を対象ガスに曝す操作とセンサ１０から対象ガスを取り除く操作のそれぞれで得られる時系列データ１４とが区別され、それぞれ異なる時系列データ１４として扱われる。例えば図６の例では、期間Ｐ１からＰ４という４つの期間それぞれで得られる時系列データが、それぞれ異なる時系列データ１４として扱われる。そのため、センサ１０を対象ガスに曝す操作とセンサ１０から対象ガスを取り除く操作とを繰り返すことで一連の時系列データを得た場合には、その一連の時系列データを複数の時系列データ１４に分割する必要がある。
【００４６】
センサ１０から得られた一連の時系列データを分割して複数の時系列データ１４を得る方法には、様々な方法を採用することができる。例えば、センサ１０から得られた一連の時系列データを人手で分割することで、複数の時系列データ１４を得る。その他にも例えば、情報処理装置２０００が一連の時系列データを取得し、その時系列データを分割して複数の時系列データ１４を得てもよい。
【００４７】
なお、情報処理装置２０００によって時系列データを分割する方法には、様々な方法を採用できる。例えば、以下のような方法がある。
【００４８】
＜＜（ａ）１階微分を用いる方法＞＞
時系列データ１４において、分割すべき部分では、センサ値の微分が不連続となり、直後に絶対値が最大となる。そこで、１階微分の絶対値が大きくなる点を用いて、時系列データ１４を分割することができる。
【００４９】
＜＜（ｂ）２階微分を用いる方法＞＞
同様に、分割すべき点では微分が不連続となるため、２階微分は無限大に発散する。したがって、２階微分の絶対値が大きくなる点を用いて、時系列データ１４を分割することができる。
【００５０】
＜＜（ｃ）センサから得られるメタデータを利用する方法＞＞
センサの種類によっては、検出値以外のメタデータが提供される。例えば MSS のモジュールでは、測定対象のガス（サンプル）とパージガスの吸引にそれぞれ異なるポンプ（サンプルポンプとパージポンプ）が用意されており、これらを交互にオン／オフすることで、立ち上がりの測定と立ち下がりの測定が行われている。そして、記録される検出値には、ポンプの動作シーケンス（どちらのポンプを利用して得られた検出値であるかを表す情報や、流量のフィードバック制御に用いる流量計測値など）が、時系列情報として付加される。そこで例えば、情報処理装置２０００は、時系列データ１４と共に得られるポンプの動作シーケンスを利用することで、時系列データ１４を分割することができる。
【００５１】
＜＜上記方法の組み合わせ＞＞
（ｃ）の方法では、ポンプが動作してからセンサにガスが届くまでの遅延を加味した修正を加えることが好ましい。そこで例えば、情報処理装置２０００は、（ｃ）の方法で時系列データ１４を仮に複数の区間に分割した後、各区間の中で１階微分の絶対値が最大となる時点を特定し、特定した各時点で時系列データ１４を分割する。
【００５２】
なお、情報処理装置２０００は、センサ１０を対象ガスに曝す操作で得られる時系列データ１４と、センサ１０から対象ガスを取り除く操作で得られる時系列データ１４のうち、いずれか一方のみを利用するように構成されてもよい。
【００５３】
＜特徴定数と寄与値の算出：Ｓ１０４＞
算出部２０４０は、時系列データ１４を用いて、複数の特徴定数θi と、特徴定数θi に対応する寄与値ξi とを算出する（Ｓ１０４）。この処理は、時系列データ１４を、式（５）に示したξi*f(θi) の和に分解することに相当する。
【００５４】
算出部２０４０が時系列データ１４から特徴定数とその寄与値を算出する方法には、様々な方法を利用することができる。以下、まずは特徴定数の算出方法について説明する。なお、特徴定数を算出する過程でその寄与値が算出されるケースもある。このケースについては後述する。
【００５５】
＜＜特徴定数の算出方法１＞＞
例えば算出部２０４０は、時系列データ１４の微分 y'(t) に log をとることで得られる log[y'(t)] の傾きに基づいて、特徴定数の集合を算出する。log[y'(t)] は、以下のようになる。
【数７】
2020026327.tif 000009
【００５６】
なお、log をとるため、y'(t) が常に正であることを前提とする。y'(t) が負になるケースでは、後述する他の方法を利用して、特徴定数を算出する。
【００５７】
図７は、g(t) のグラフを例示する図である。式（８）より、g(t)=log[y'(t)] は、いわゆる log-sum-exp と呼ばれる形式の関数であり、この関数は、積算されている exp のうちの最大値に近い値をとることが知られている。すなわち、式（８）では、g(t)の値が、その時刻 t における ci-βit のうちで最大のものに近い値となる。
【００５８】
ここで、近似的に g(t)=ci-βit とおくと、この関数の傾き g'(t) は速度定数 βi となる。よって、時刻 t における g(t) の傾き g'(t) は、その時刻において ci-βit が最大となる i に対応するβi と略同一になる。また、exp(ci-βit) は単調減少関数であることから、図７に示すように、g(t) の傾きは段階的に切り替わると言える。
【００５９】
そこで算出部２０４０は、各時刻 t における g(t) の傾き g'(t) を算出し、g(t) の定義域（時系列データ１４の測定期間）から、g'(t) が略同一である複数の期間（以下、部分期間）を複数抽出する。そして算出部２０４０は、抽出した各部分期間における g'(t) を、速度定数βi として算出する。例えば算出部２０４０は、各部分期間における g'(t) の統計値（平均値など）を、その部分期間に対応する速度定数として算出する。その他にも例えば、算出部２０４０は、部分期間ごとに、その部分期間に含まれる (t, g(t)) について線形回帰を行い、各部分期間について得られる回帰直線の傾きを、速度定数としてもよい。なお、g'(t) が略同一である期間は、例えば、g'(t) の集合をその値の大きさでクラスタリングすることで特定することができる。
【００６０】
例えば図７において、g(t) の定義域から、３つの部分期間が抽出されている。算出部２０４０は、これら３つの部分期間それぞれに対応する速度定数として、β1、β2、及びβ3 を算出する。このようにして、特徴定数の集合Θ＝{β1, β2, β3} が算出される。
【００６１】
ここで、図７に示す区間２０のように、g(t) の傾きが切り替わる部分に、その前後の傾きのいずれとも異なる傾きが、短い期間現れることがある。このような短い期間だけ現れる傾きについては、特徴定数に含めないようにしてもよい。
【００６２】
例えば、部分期間の長さの最小値を予め定めておく。そして、算出部２０４０は、傾き g'(t) の大きさに基づいて g(t) の定義域を複数の期間に分割し、最小値以上の長さを持つ期間のみを、前述した部分期間として抽出する（すなわち、最小値未満の長さの期間に対応する g'(t) は特徴定数の集合に含めない）。こうすることで、短い期間だけ現れる g'(t) については、特徴定数に含めないようにすることができる。なお、部分期間の長さの最小値は、部分期間に含まれる検出値の数の最小値として表すこともできる。
【００６３】
その他にも例えば、算出部２０４０は、後述する方法で特徴定数の個数 m を決定し、その個数の部分期間を g(t) の定義域から抽出してもよい。この場合、例えば算出部２０４０は、傾き g'(t) の大きさに基づいて g(t) の定義域を複数の期間に分割し、その中から期間の長さについて上位にある m 個の期間を、前述した部分期間として抽出する。
【００６４】
＜＜特徴定数の算出方法２＞＞
前述した算出方法１を採用するためには、y'(t) の log をとる必要があることから、y'(t) の値が常に正である必要がある。これに対し、ここで説明する算出方法２では、前述した g(t) の傾き g'(t) に相当する指標を利用しつつ、y'(t) の値が負になることがあっても利用することができる。
【００６５】
まず、前述した g(t)=log[y'(t)] の傾きである g'(t) は、以下のように表すことができる。
【数８】
2020026327.tif 000010
【００６６】
そこで算出方法２では、y'(t) の log をとらず、y''(t)/y'(t) を g(t) の傾き g'(t) に相当する指標として利用することで、y'(t) が負になることがあるケースも扱えるようにする。
【００６７】
ここで、g'(t) は、(y'(t), y''(t)) というベクトルの向きということができる。さらに言い換えれば、g'(t) は、(y(t), y'(t)) というベクトルの速度ベクトル（時間変化を表すベクトル）の向きであるといえる。
【００６８】
そこで算出部２０４０は、各時刻 t についてベクトル (y(t), y'(t)) を算出し、さらに算出した各ベクトルを利用して、各時刻 t について速度ベクトル (y(t+1)-y(t), y'(t+1)-y'(t)) を算出する。そして算出部２０４０は、算出された複数の速度ベクトルについて、その向きが略同一である複数の部分期間を、y(t) の測定期間から抽出する。なお、y(t+1) は、センサにおいて y(t) の次に得られる検出値である。
【００６９】
そして、算出部２０４０は、抽出した各部分期間における (y(t), y'(t)) の速度ベクトルの向きを、その部分期間に対応する速度定数として算出する。例えば算出部２０４０は、各部分期間における上記速度ベクトルの向きの統計値（平均値など）を、その部分期間に対応する速度定数として算出する。その他にも例えば、算出部２０４０は、部分期間ごとに、その部分期間に含まれる点 (y(t), y'(t)) について線形回帰を行い、各部分期間について得られる回帰直線の傾きを、その部分期間に対応する速度定数としてもよい。なお、y(t) の測定期間の分割は、上記速度ベクトルをその向きでクラスタリングすることで実現することができる。なお、速度ベクトル (y(t+1)-y(t), y'(t+1)-y'(t)) の向きは、{y'(t+1)-y'(t)}/{y(t+1)-y(t)} で表される。
【００７０】
図８は、(y(t), y'(t)) の速度ベクトルを例示する図である。図８では、時系列データ１４から２つの部分期間が得られている。算出部２０４０は、部分期間１における速度ベクトルの向きに基づいて速度定数β1 を算出し、部分期間２における速度ベクトルの向きに基づいて速度定数β2 を算出している。
【００７１】
ここで、図８に示す区間３０のように、速度ベクトルの向きが変化する部分に、その前後の向きのいずれとも異なる向きが、短い期間現れることがある。これは、図７の区間２０と同様である。このような短い期間だけ現れる速度ベクトルの向きについては、特徴定数に含めないようにしてもよい。その具体的な方法には、方法１で説明したように、部分期間の長さの最小値を定めておいたり、特徴定数の個数を決定したりする方法を利用することができる。
【００７２】
＜＜寄与値の算出方法＞＞
ここでは、特徴定数を算出した後に、その特徴定数に対応する寄与値を算出する方法を説明する。算出部２０４０は、算出された特徴定数の集合Θ={θ1,..., θm} に対応する寄与値ξiの集合（すなわち、寄与ベクトル）Ξ={ξ1,..., ξm} をパラメータとして、センサ１０の検出値を予測する予測モデルを生成する。この予測モデル生成する際、観測データである時系列データ１４を利用して寄与ベクトルΞについてパラメータ推定を行うことにより、寄与ベクトルΞを算出することができる。特徴定数として速度定数βを使う場合の予測モデルの一例は、式（６）で表すことができる。また、特徴定数として時定数τを使う場合の予測モデルの一例は、式（７）で表すことができる。
【００７３】
予測モデルのパラメータ推定には、種々の方法を利用することができる。以下、その方法についていくつか例示する。なお、以下の説明では、速度定数βを特徴定数として利用するケースに説明している。時定数τを特徴定数とする場合におけるパラメータ推定の方法は、以下の説明における速度定数βを1/τと読み替えることで実現できる。
【００７４】
＜＜パラメータ推定の方法１＞＞
例えば算出部２０４０は、予測モデルから得られる予測値と、センサ１０から得られた観測値（すなわち、時系列データ１４）とを用いた最尤推定により、パラメータΞを推定する。最尤推定には、例えば最小二乗法を用いることができる。この場合、具体的には、以下の目的関数に従ってパラメータΞを決定する。
【数９】
2020026327.tif 000011
ここで、T は時系列データ１４の長さ（検出値の数）を表す。また、y^(ti) は、時刻 ti の予測値を表す。
【００７５】
上述の目的関数を最小化するベクトルΞは、以下の式（１１）を用いて算出することができる。
【数１０】
2020026327.tif 000012
ここで、ベクトル Y=(y(t0), y(t1),...) である。
【００７６】
そこで、算出部２０４０は、時系列データ Yと特徴定数の集合Θ={β1, β2,...}を上記式（１１）に適用することで、パラメータΞを算出する。
【００７７】
＜＜パラメータ推定の方法２＞＞
上述した最小二乗法において、正則化項を導入して正則化を行ってもよい。例えば以下の式（１２）は、L2 正則化を行う例を示している。 
【数１１】
2020026327.tif 000013
ここで、λは正則化項に与える重みを表すハイパーパラメータである。
【００７８】
この場合、以下の式（１３）に従ってパラメータΞを決定することができる。
【数１２】
2020026327.tif 000014
【００７９】
このような正則化項を導入すると、正則化項を導入しない場合と比較し、行列計算において測定誤差が増幅されてしまうことを抑えることができるため、各寄与値ξi をより正確に算出することができる。また、誤差の増幅を抑えることで、寄与値ξが数値的に安定するため、混合比に対する特徴量のロバスト性が向上する。
【００８０】
なお、前述したように、λはハイパーパラメータであり、予め定めておく必要がある。例えば、テスト測定やシミュレーションを通じて、λの値を決定する。ここで、λの値は、寄与値ξが振動しない程度に小さい値とすることが好ましい。
【００８１】
ここでは、λの値を決定するためのシミュレーションについて説明する。シミュレーションでは、仮想的に「寄与が１の単一の分子」を測定した場合（たとえば、立下りの場合、単一の分子の速度定数をβ0とすれば、y(t)=exp{-β0*t} となる）を考え、この場合の式（１３）による特徴量推定値の結果を観察する。仮想的に、理想的な観測(無限小の測定間隔で無限長時間の測定ができ、観測誤差がゼロ)が可能であるとした場合、仮想的な単一分子のシミュレーションでは、以下のようにβ0だけに尖ったピークがある特徴量が得られ、もとの速度定数β=β0と寄与ξ=1 が完全に再現されるはずである。
【数１３】
2020026327.tif 000015
【００８２】
しかし、実際には理想的な観測は不可能であるため、寄与値のピークが鈍化してしまったり、寄与値が振動してしまったりする。図７は、単一分子について得られる特徴量を例示する図である。この図から、λの値による、ピークの鈍化と振動の増加のトレードオフが分かる。具体的には、λが大きすぎると、振動は減るものの、ピークの幅が大きくなってしまう。ピークの幅が大きくなると、速度定数が近い２つの分子を測定した結果が、１つの大きなピークに見えてしまうため、これらの分子を見分けることが難しくなる。すなわち、感度が下がることになる。一方で、λが小さすぎると、ピークの幅は小さくなるものの、振動が増えてしまう。振動が増えると、後述するように、特徴量のロバスト性が低下してしまう。よって、振動が発生しない（ロバスト性を損なわない）程度にλを小さくするように決定することで、ピークを鋭く（感度を向上）することが好ましいと言える。
【００８３】
シミュレーションの目的は、λを変化させながら、このようなピークの鈍化や振動の発生度合を評価することである。「振動の大きさ」と「ピーク幅」を定量的に測るためには、例えばシミュレーションにより、2つの異なる速度定数β１、β２をそれぞれ持つ２つの仮想的な単一分子の寄与ベクトルΞ１、Ξ２を計算する。そして、これら２つの寄与ベクトルの内積を以下のように計算する。
【数１４】
2020026327.tif 000016
【００８４】
この関数 f(Δv) は、振動しながら減衰する。そこで、この振動のメインローブの幅を「ピーク幅」、サイドローブのレベルを「振動の大きさ」として定量化することができる。そこで、メインローブ幅がなるべく細く、なおかつサイドローブのレベルがなるべく小さくなるようなλの値を選択することにより、λを決定する。
【００８５】
ここで、寄与値の振動を抑えることの利点の１つは、前述したように、特徴量が時定数や速度定数の変化に対してロバストになることである。言い換えれば、特徴量が温度変化に対してロバストになる。以下、その理由を説明する。
【００８６】
温度変化によって時定数や速度定数に変化が生じると、図７や後述する図８に示す寄与値は、Ｘ軸方向に平行移動することになる。寄与値が大きく振動していると、寄与値がＸ軸方向に少し平行移動しただけでも、平行移動前後の寄与ベクトルの距離が大きくなってしまう。すなわち、時定数や速度定数が少し変化しただけで特徴量が大きく変化してしまい、時定数変化や速度定数変化に対する特徴量のロバスト性が低くなる。
【００８７】
これに対し、寄与値の振動が少なければ、平行移動前後の寄与ベクトルの距離は短くなる。これは、時定数や速度定数が少し変化した場合に、特徴量も少しだけ変化することを意味する。すなわち、特徴量のロバスト性が高いことを意味する。よって、特徴量の振動を抑えることで、特徴量のロバスト性が向上すると言える。
【００８８】
なお、最小二乗法における正則化は、前述した L2 正則化には限定されず、L1 正則化などの他の正則化を導入してもよい。
【００８９】
＜＜パラメータ推定の方法３＞＞
この方法では、パラメータΞについて事前分布 P(Ξ) を設定しておく。そして算出部２０４０は、観測値である時系列データ１４を用いた MAP （Maximum a Posteriori）推定により、パラメータΞを決定する。具体的には、以下の目的関数を最大化するパラメータΞを採用する。
【数１５】
2020026327.tif 000017
【００９０】
P(Y｜Ξ)とP(Ξ)は、例えば次のように多変量正規分布で定義される。
【数１６】
2020026327.tif 000018
ここで、N(・|μ,Σ)は、平均μ、共分散Σの多変量正規分布である。また、ベクトル y^=(y^(t1), y^(t2),...)=ΦΞである。σ^2 は観測誤差の分散を表すパラメータである。
【００９１】
ΛはΞの事前分布の共分散行列であり、予め任意の半正定値行列を与えても良いし、後述の方法などにより決定しても良い。
【００９２】
また、P（Y|Ξ）、P（Ξ）は、次のように、ガウス過程(Gaussian Process; GP)によって定めても良い。
【数１７】
2020026327.tif 000019
【００９３】
ここで、GP(ξ(β)|μ(β),Λ(β,β’))は、平均値関数がμ(β)、共分散関数（カーネル関数）がΛ（β，β’）のガウス過程である。また、ガウス過程は連続関数を生成する確率過程であるため、ここでは、ξ(β)は、寄与率をβ(もしくはτ)に関して表した連続関数であり、ベクトルΞは、関数ξ（β）の「β=β1,β2,・・・」における値を配列したベクトルΞ=(ξ(β1),ξ(β2),...)である。この場合、式(１７) は式(１８)の特別な場合とみなすことができ、式(１７) における共分散行列Λの(i,j)成分は、式(１８)の共分散関数Λ(β,β’)の(β,β’)=(β1, β2)における値である。すなわち、式(１７)における行列Λは、所謂ガウス過程におけるグラム行列である。
【００９４】
また、算出部２０４０は、観測値である時系列データ１４を用いたベイズ推定により、パラメータΞを決定しても良い。具体的には、以下の条件付き期待値を計算することにより、パラメータΞを決定する。
【数１８】
2020026327.tif 000020
ここで、E[Ξ|Y］ は、ΞおよびＹが式(１８)の確率分布に従っていると仮定した場合の、条件付き期待値である。
【００９５】
上記目的関数（１４）を最大化する特徴ベクトルΞ、および、上記条件付き期待値（１９）によって求まる特徴ベクトルΞは、いずれも以下の式(２０)によって算出することができる。
【数１９】
2020026327.tif 000021
【００９６】
＜＜＜ハイパーパラメータの定め方＞＞＞
ガウス過程を利用する場合、事前に設定しておくハイパーパラメータとして、a）共分散関数 Λ（β,β’）の形、b）共分散関数のパラメータ、及びc）測定誤差パラメータ σ^2 がある。これらを変えながら、次の手順を実行する。
【００９７】
１．仮想的な速度定数β0の単一分子の測定値をシミュレーションする。
２．シミュレーションした測定値から、寄与値を推定する。
３．推定した寄与値の振動の大きさ・ピーク幅を定量化する。
４．上述の a 〜 c のハイパーパラメータを変更しながら、1〜3を繰り返す。
５．グリッドサーチや最急降下法により、振動が小さく、ピーク幅が狭くなるように a 〜 c のハイパーパラメータを決定する。
【００９８】
なお、特徴量の振動の大きさとピーク幅を定量化する指標には、例えば、上述の関数 f(Δv) のインローブの幅とサイドローブのレベルを用いる。また、そのほかにも、推定したΞを確率分布として見做した際の分散値(二乗分散や絶対値分散)を用いてもよい。これらの分散値は、振動が小さく、ピーク幅が狭くなるほど、小さい値となる。なお、シミュレーションの代わりに、実際の測定（テスト測定）を実施してもよい。
【００９９】
＜特徴定数と共に寄与値が算出されるケース＞
前述したように、特徴定数に対応する寄与値の算出方法の１つとして、最小二乗法を利用する方法がある。算出部２０４０は、特徴定数を算出した後に最小二乗法を利用する代わりに、前述した最小二乗法の目的関数の最小値を特徴定数の集合Θに関して最小化するという組み合わせ最適化問題を解くことにより、特徴定数の集合Θと寄与値の集合Ξの双方を算出してもよい。具体的には、以下の目的関数 h(Θ) を最小化するΘと、最小化された h(Θ) におけるパラメータΞが、算出部２０４０の算出結果となる。
【数２０】
2020026327.tif 000022
【０１００】
ここで、上記の組み合わせ最適化問題を解く具体的な方法には、既存の方法を利用することができる。特に、Ξの要素がすべて正の数であるとき、上記の目的関数は単調減少な優モジュラ関数になるため、たとえば、貪欲アルゴリズムによって精度よく計算することができる。
【０１０１】
＜特徴定数の個数を決定する方法＞
算出部２０４０は、特徴定数の個数を決定してもよい。決定した特徴定数の個数は、例えば前述した、g(t) の定義域や y(t) の測定期間から抽出する部分期間の数を定めるために利用することができる。
【０１０２】
例えば特徴定数の個数は、AIC（Akaike's Information Criterion）や BIC（Bayesian Information Criterion）などの情報量基準を利用して決定することができる。情報量規準の算出には、前述した式（２２）の関数 h(Θ) を利用する。例えば AIC は、以下の様に算出することができる。なお、下記 AIC の導出方法については後述する。
【数２１】
2020026327.tif 000023
ここで、σ^2 は観測誤差の分散である。
【０１０３】
算出部２０４０は、AIC が最小となる整数 K を、特徴定数の個数として決定する。図９は、AIC をグラフで表す図である。図９において、上段は関数 L のグラフであり、下段は AIC のグラフである。AIC が最小になる整数 K は、K* である。そのため、算出部２０４０は、特徴定数の個数を K* とする。
【０１０４】
このように情報量規準を利用して特徴定数の個数を決定することにより、予測モデルのパラメータ（ここでは特徴定数や寄与値の数）の個数を増やすことによって得られる予測精度の向上というメリットと、パラメータの個数を増やすことによってモデルの複雑さが増すというデメリットとのバランスを考慮して、特徴定数の個数を適切に決定することができる。
【０１０５】
式（２４）の AIC を導出する方法について説明する。まず、AIC の定義は以下の通りである。
【数２２】
2020026327.tif 000024
ただし、l は最大尤度であり、p はパラメータ数である。
【０１０６】
パラメータ数 p は、特徴定数と寄与値が K 個ずつであるため、2K である。また、最大尤度は以下のように算出される。
【数２３】
2020026327.tif 000025
【０１０７】
上記 l と p を AIC の定義に代入することで、式（２４）を得ることができる。
【０１０８】
＜特徴量の出力：Ｓ１０６＞
出力部２０６０は、前述した方法で得られた特徴定数の集合Θと寄与値の集合Ξとを対応づけた情報（以下、出力情報）を、ガスの特徴を表す特徴量として出力する（Ｓ１０６）。例えば出力情報は、特徴行列 F を表すテキストデータである。その他にも例えば、出力情報は、特徴定数の集合Θと寄与値の集合Ξとの対応付けを、表やグラフなどでグラフィカルに表現した情報であってもよい。
【０１０９】
図１０は、対象ガスの特徴量をグラフで例示する図である。図１０のグラフは、横軸に時定数τを示し、縦軸に寄与値ξを示している。そして、算出部２０４０によって算出された（τi, ξi）のペアがグラフ上にプロットされている。このように対象ガスの特徴量をグラフィカルな情報で表すことにより、人がガスの特徴を直感的に理解しやすくなる。
【０１１０】
出力情報を出力する具体的な方法は様々である。例えば出力部２０６０は、出力情報を任意の記憶装置に記憶させる。その他にも例えば、出力部２０６０は、出力情報をディスプレイ装置に表示させる。その他に例えば、出力部２０６０は、情報処理装置２０００以外の装置に出力情報を送信してもよい。
【０１１１】
＜特徴定数の集合と寄与値の集合との対応付けを複数算出するケース＞
情報処理装置２０００は、同一の対象ガスについて得られた複数の時系列データ１４それぞれについて、特徴定数の集合Θと寄与値の集合Ξとの対応付けを算出してもよい。そしてこの場合、出力部２０６０は、これら複数の対応付けの集合を、対象ガスの特徴量としてもよい。
【０１１２】
例えば情報処理装置２０００は、立ち上がりの時系列データ１４について、特徴定数の集合Θu と寄与値の集合Ξu とを算出し、これらを対応づけた特徴行列 Fu を生成する。また、情報処理装置２０００は、立ち下がりの時系列データ１４について、特徴定数の集合Θd と寄与値の集合Ξd とを算出し、これらを対応づけた特徴行列Fd を生成する。そして情報処理装置２０００は、生成した特徴行列の組である {Fu, Fd} を対象ガスの特徴量として出力する。
【０１１３】
図１１は、立ち上がりの時系列データ１４と立ち下がりの時系列データ１４それぞれから特徴行列を得るケースを例示する図である。図１１において、立ち上がりの時系列データある時系列データ１４−１から、特徴行列Fu が得られている。また、立ち下がりの時系列データである時系列データ１４−２から、特徴行列Fd が得られている。そこで出力部２０６０は、得られた２つの特徴行列を組み合わせた {Fu, Fd} を、対象ガスの特徴量として出力する。
【０１１４】
なお、出力部２０６０は、立ち上がりの時系列データ１４から得られる特徴行列と、立ち下がりの時系列データ１４から得られる特徴行列とを連結することで得られる１つの行列を、対象ガスの特徴量としてもよい。例えばこの場合、出力部２０６０は、Fu=(ΘuT, ΞuT) とFd=(ΘdT, ΞdT) を連結したFc=(ΘuT, ΞuT, ΘdT, ΞdT) を、対象ガスの特徴量として出力する。なお、連結する特徴行列の行数が互いに異なる場合、行数が少ない方の行列をゼロパディング等の方法で拡張することで、連結する特徴行列の行数を互いに一致させる。
【０１１５】
複数の特徴行列は、立ち上がりの時系列データ１４と立ち下がりの時系列データ１４それぞれから得られるものに限定されない。例えば、特性の異なる複数のセンサ１０それぞれを対象ガスに曝すことで、複数の時系列データ１４を得てもよい。分子をセンサに付着させる場合、センサに対する各分子の付着しやすさは、センサの特性によって異なる。例えば官能膜に分子が付着するタイプのセンサを利用する場合、官能膜の材質によって、その官能膜に対する各分子の付着しやすさが異なる。各分子の離脱しやすさについても同様である。そのため、それぞれ異なる材質の官能膜を持つセンサ１０を用意し、これら複数のセンサ１０それぞれから時系列データ１４を得て解析することで、対象ガスの特徴をより正確に把握することができる。
【０１１６】
そこで情報処理装置２０００は、特性の異なる複数のセンサ１０それぞれから時系列データ１４を取得し、各時系列データ１４について、特徴定数の集合と寄与値の集合とを対応づけた情報を生成する。出力部２０６０は、このようにして得られた複数の情報の組を、対象ガスの特徴量として出力する。
【０１１７】
図１２は、複数のセンサ１０それぞれから時系列データ１４を得ることで、複数の特徴行列を得るケースを例示する図である。この例では、それぞれ特性が異なる３つのセンサ１０−１、センサ１０−２、及びセンサ１０−３が用意されており、それぞれから時系列データ１４−１、時系列データ１４−２、及び時系列データ１４−３が得られている。情報処理装置２０００は、これら複数の時系列データ１４からそれぞれ、特徴行列F1、F2、及び F3 を算出する。そして情報処理装置２０００は、これら３つの特徴行列の組を、対象ガスの特徴量として出力する。なお、前述したように、複数の特徴行列の組を出力する代わりに、これら複数の特徴行列を連結した１つの特徴行列Fc を出力してもよい。
【０１１８】
ここで、特性の異なる複数のセンサ１０は、１つの筐体に収められてもよいし、それぞれ異なる筐体に収められてもよい。前者の場合、例えば、１つのセンサ筐体の中に材質の異なる複数の官能膜を収納し、各官能膜について検出値が得られるように、センサ１０が構成される。
【０１１９】
さらに、図１１で説明した方法と、図１２で説明した方法を組み合わせてもよい。すなわち、情報処理装置２０００は、複数のセンサ１０それぞれから、立ち上がりの時系列データ１４と立ち下がりの時系列データ１４を得て、得られた各時系列データ１４について特徴行列 F を算出し、算出した複数の特徴行列の組や、これらを連結した１つの特徴行列を、対象ガスの特徴量としてもよい。
【０１２０】
＜バイアスを考慮した特徴量の算出＞
センサ１０の検出値には、時間に応じた変化を表さないバイアス項が含まれていることがある。この場合、時系列データ１４は以下のように表される。なお、ここでは特徴定数として、速度定数βを用いている。
【数２４】
2020026327.tif 000026
【０１２１】
バイアスは、例えば、センサ１０のオフセットがずれていることによって生じる。その他にも例えば、バイアスは、対象ガスとパージガスに共通して含まれている成分の寄与（例えば、大気中の窒素や酸素の寄与）により生じる。
【０１２２】
情報処理装置２０００は、時系列データ１４からオフセットを除去する機能を有していてもよい。こうすることで、対象ガスの特徴量をより正確に算出することができる。以下、オフセットを考慮して特徴量を算出する方法について説明する。
【０１２３】
算出部２０４０は、上記式（２７）で表される時系列データ１４の予測モデルを生成することにより、バイアスを考慮して寄与ベクトルΞを算出する。すなわち、算出部２０４０は、式（２７）で表現された予測モデルについて、パラメータΞ及び b の推定を行う。具体的には、算出部２０４０は、目的関数（１０）、（１２）、又は（１６）を、Ξのみならず、b についても最適化することにより、Ξと b を推定する。なお、特徴定数として時定数を用いる場合には、式（２７）においてβk を1/τk に置き換える。
【０１２４】
例えば、目的関数として式(１４)を用いるとする。この場合、算出部２０４０は、Ξ及び b を以下の最適化問題によって計算する。目的関数に（８）や（１０）を用いる場合も同様である。
【数２５】
2020026327.tif 000027
【０１２５】
上記最適化問題の解Ξ及び ｂ は、以下の式によって算出することができる。
【数２６】
2020026327.tif 000028
【０１２６】
このように、バイアス b と寄与ベクトルΞの両方を推定することにより、寄与ベクトルからバイアスの効果が取り除かれ、センサ１０の検定値にバイアスが含まれている場合にも正確に寄与ベクトルを計算することができる。
【０１２７】
なお、出力部２０６０は、特徴行列 F に加え、バイアス b や b0 を出力しても良い。バイアスが上記センサのオフセットのずれによるものである場合、b0の値を用いてセンサのオフセットをキャリブレーションすることができる。
【０１２８】
以上、図面を参照して本発明の実施形態について述べたが、これらは本発明の例示であり、上記以外の様々な構成を採用することもできる。
【０１２９】
上記の実施形態の一部又は全部は、以下の付記のようにも記載されうるが、以下には限られない。
１．  対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得部と、
前記時系列データに対して寄与する複数の特徴定数と、各前記特徴定数の前記時系列データに対する寄与の大きさを表す寄与値とを算出する算出部と、
複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、前記センサによってセンシングされたガスの特徴量として出力する出力部と、を有し、
前記特徴定数は、前記センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である、情報処理装置。
２．  前記算出部は、
前記時系列データの測定期間から複数の部分期間を抽出し、
前記部分期間ごとに、その部分期間における前記検出値の時間変化率の対数に基づいて、前記特徴定数を算出し、
前記部分期間は、その中に含まれる検出値の時間変化率の対数が略同一の期間である、１．に記載の情報処理装置。
３．  前記算出部は、
前記時系列データを用いて、各時刻における検出値及びその時刻における検出値の時間変化率を要素とする時系列のベクトルデータを算出し、
各前記算出した時系列のベクトルデータそれぞれについて速度ベクトルを算出し、
前記速度ベクトルの向きに基づいて、前記時系列データの測定期間から、複数の部分期間を抽出し、
前記部分期間ごとに、その部分期間における前記速度ベクトルの向きに基づいて、前記特徴定数を算出し、
前記部分期間は、その中に含まれる前記速度ベクトルの向きが略同一の期間である、１．に記載の情報処理装置。
４．  前記算出部は、複数の特徴定数それぞれの寄与値をパラメータとする前記センサの検出値の予測モデルについて、前記取得した時系列データを用いたパラメータ推定を行うことで、各寄与値を算出する、１．乃至３．いずれか一つに記載の情報処理装置。
５．  前記算出部は、前記予測モデルから得られる時系列データと前記取得した時系列データとについて、最小二乗法による最尤推定を行うことで、各前記寄与値を算出する、４．に記載の情報処理装置。
６．  前記最小二乗法における最尤推定において、目的関数に正則化項が含まれている、５．に記載の情報処理装置。
７．  前記算出部は、各前記寄与値の事前分布と前記取得した時系列データとを用いた MAP（Maximum a Posteriori）推定又はベイズ推定により、各前記寄与値を算出する、４．に記載の情報処理装置。
８．  前記事前分布は、多変量正規分布又はガウス過程である、７．に記載の情報処理装置。
９．  前記算出部は、前記予測モデルから得られる時系列データと前記取得した時系列データとの二乗誤差を表す目的関数について、その目的関数の最小値を複数の特徴定数に関して最小化することにより、複数の特徴定数及び複数の寄与値を算出する、４．に記載の情報処理装置。
１０．  前記予測モデルには、バイアスを表すパラメータが含まれており、
前記算出部は、前記予測モデルについて寄与値とバイアスそれぞれを表すパラメータを推定する、４．乃至９．いずれか一つに記載の情報処理装置。
１１．  前記時系列データ取得部は、複数の時系列データを取得し、
前記算出部は、複数の時系列データそれぞれについて、特徴定数の集合と寄与値の集合の組みを算出し、
前記出力部は、前記算出された特徴定数の集合と寄与値の集合の組みを複数まとめた情報を、前記対象ガスの特徴量として出力する、１．乃至１０．いずれか一つに記載の情報処理装置。
１２．  前記複数の時系列データは、前記センサを前記対象ガスに曝す際に得られる時系列データと、前記センサから前記対象ガスを取り除く際に得られる時系列データとの双方を含む、１１．に記載の情報処理装置。
１３．  前記複数の時系列データは、特性の異なる複数の前記センサそれぞれから得られる時系列データを含む、１１．に記載の情報処理装置。
【０１３０】
１４．  コンピュータによって実行させる制御方法であって、
対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得ステップと、
前記時系列データに対して寄与する複数の特徴定数と、各前記特徴定数の前記時系列データに対する寄与の大きさを表す寄与値とを算出する算出ステップと、
複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、前記センサによってセンシングされたガスの特徴量として出力する出力ステップと、を有し、
前記特徴定数は、前記センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である、制御方法。
１５．  前記算出ステップにおいて、
前記時系列データの測定期間から複数の部分期間を抽出し、
前記部分期間ごとに、その部分期間における前記検出値の時間変化率の対数に基づいて、前記特徴定数を算出し、
前記部分期間は、その中に含まれる検出値の時間変化率の対数が略同一の期間である、１４．に記載の制御方法。
１６．  前記算出ステップにおいて、
前記時系列データを用いて、各時刻における検出値及びその時刻における検出値の時間変化率を要素とする時系列のベクトルデータを算出し、
各前記算出した時系列のベクトルデータそれぞれについて速度ベクトルを算出し、
前記速度ベクトルの向きに基づいて、前記時系列データの測定期間から、複数の部分期間を抽出し、
前記部分期間ごとに、その部分期間における前記速度ベクトルの向きに基づいて、前記特徴定数を算出し、
前記部分期間は、その中に含まれる前記速度ベクトルの向きが略同一の期間である、１４．に記載の制御方法。
１７．  前記算出ステップにおいて、複数の特徴定数それぞれの寄与値をパラメータとする前記センサの検出値の予測モデルについて、前記取得した時系列データを用いたパラメータ推定を行うことで、各寄与値を算出する、１４．乃至１６．いずれか一つに記載の制御方法。
１８．  前記算出ステップにおいて、前記予測モデルから得られる時系列データと前記取得した時系列データとについて、最小二乗法による最尤推定を行うことで、各前記寄与値を算出する、１７．に記載の制御方法。
１９．  前記最小二乗法における最尤推定において、目的関数に正則化項が含まれている、１８．に記載の制御方法。
２０．  前記算出ステップにおいて、各前記寄与値の事前分布と前記取得した時系列データとを用いた MAP（Maximum a Posteriori）推定又はベイズ推定により、各前記寄与値を算出する、１７．に記載の制御方法。
２１．  前記事前分布は、多変量正規分布又はガウス過程である、２０．に記載の制御方法。
２２．  前記算出ステップにおいて、前記予測モデルから得られる時系列データと前記取得した時系列データとの二乗誤差を表す目的関数について、その目的関数の最小値を複数の特徴定数に関して最小化することにより、複数の特徴定数及び複数の寄与値を算出する、１７．に記載の制御方法。
２３．  前記予測モデルには、バイアスを表すパラメータが含まれており、
前記算出ステップにおいて、前記予測モデルについて寄与値とバイアスそれぞれを表すパラメータを推定する、１７．乃至２２．いずれか一つに記載の制御方法。
２４．  前記時系列データ取得ステップにおいて、複数の時系列データを取得し、
前記算出ステップにおいて、複数の時系列データそれぞれについて、特徴定数の集合と寄与値の集合の組みを算出し、
前記出力ステップにおいて、前記算出された特徴定数の集合と寄与値の集合の組みを複数まとめた情報を、前記対象ガスの特徴量として出力する、１４．乃至２３．いずれか一つに記載の制御方法。
２５．  前記複数の時系列データは、前記センサを前記対象ガスに曝す際に得られる時系列データと、前記センサから前記対象ガスを取り除く際に得られる時系列データとの双方を含む、２４．に記載の制御方法。
２６．  前記複数の時系列データは、特性の異なる複数の前記センサそれぞれから得られる時系列データを含む、２４．に記載の制御方法。
【０１３１】
２７．  １４．乃至２６．いずれか一つに記載の制御方法の各ステップをコンピュータに実行させるプログラム。
【図１】
2020026327.tif 000029
【図２】
2020026327.tif 000030
【図３】
2020026327.tif 000031
【図４】
2020026327.tif 000032
【図５】
2020026327.tif 000033
【図６】
2020026327.tif 000034
【図７】
2020026327.tif 000035
【図８】
2020026327.tif 000036
【図９】
2020026327.tif 000037
【図１０】
2020026327.tif 000038
【図１１】
2020026327.tif 000039
【図１２】
2020026327.tif 000040
【手続補正書】
【提出日】令和3年1月21日(2021.1.21)
【手続補正１】
【補正対象書類名】特許請求の範囲
【補正対象項目名】全文
【補正方法】変更
【補正の内容】
【特許請求の範囲】
【請求項１】
対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得部と、
前記時系列データに対して寄与する複数の特徴定数と、各前記特徴定数の前記時系列データに対する寄与の大きさを表す寄与値とを算出する算出部と、
複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、前記センサによってセンシングされたガスの特徴量として出力する出力部と、を有し、
前記特徴定数は、前記センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である、情報処理装置。
【請求項２】
前記算出部は、
前記時系列データの測定期間から複数の部分期間を抽出し、
前記部分期間ごとに、その部分期間における前記検出値の時間変化率の対数に基づいて、前記特徴定数を算出し、
前記部分期間は、その中に含まれる検出値の時間変化率の対数が略同一の期間である、請求項１に記載の情報処理装置。
【請求項３】
前記算出部は、
前記時系列データを用いて、各時刻における検出値及びその時刻における検出値の時間変化率を要素とする時系列のベクトルデータを算出し、
各前記算出した時系列のベクトルデータそれぞれについて速度ベクトルを算出し、
前記速度ベクトルの向きに基づいて、前記時系列データの測定期間から、複数の部分期間を抽出し、
前記部分期間ごとに、その部分期間における前記速度ベクトルの向きに基づいて、前記特徴定数を算出し、
前記部分期間は、その中に含まれる前記速度ベクトルの向きが略同一の期間である、請求項１に記載の情報処理装置。
【請求項４】
前記算出部は、複数の特徴定数それぞれの寄与値をパラメータとする前記センサの検出値の予測モデルについて、前記取得した時系列データを用いたパラメータ推定を行うことで、各寄与値を算出する、請求項１乃至３いずれか一項に記載の情報処理装置。
【請求項５】
前記算出部は、前記予測モデルから得られる時系列データと前記取得した時系列データとについて、最小二乗法による最尤推定を行うことで、各前記寄与値を算出する、請求項４に記載の情報処理装置。
【請求項６】
前記最小二乗法における最尤推定において、目的関数に正則化項が含まれている、請求項５に記載の情報処理装置。
【請求項７】
前記算出部は、各前記寄与値の事前分布と前記取得した時系列データとを用いた MAP（Maximum a Posteriori）推定又はベイズ推定により、各前記寄与値を算出する、請求項４に記載の情報処理装置。
【請求項８】
前記事前分布は、多変量正規分布又はガウス過程である、請求項７に記載の情報処理装置。
【請求項９】
前記算出部は、前記予測モデルから得られる時系列データと前記取得した時系列データとの二乗誤差を表す目的関数について、その目的関数の最小値を複数の特徴定数に関して最小化することにより、複数の特徴定数及び複数の寄与値を算出する、請求項４に記載の情報処理装置。
【請求項１０】
前記予測モデルには、バイアスを表すパラメータが含まれており、
前記算出部は、前記予測モデルについて寄与値とバイアスそれぞれを表すパラメータを推定する、請求項４乃至９いずれか一項に記載の情報処理装置。
【請求項１１】
前記時系列データ取得部は、複数の時系列データを取得し、
前記算出部は、複数の時系列データそれぞれについて、特徴定数の集合と寄与値の集合の組みを算出し、
前記出力部は、前記算出された特徴定数の集合と寄与値の集合の組みを複数まとめた情報を、前記対象ガスの特徴量として出力する、請求項１乃至１０いずれか一項に記載の情報処理装置。
【請求項１２】
前記複数の時系列データは、前記センサを前記対象ガスに曝す際に得られる時系列データと、前記センサから前記対象ガスを取り除く際に得られる時系列データとの双方を含む、請求項１１に記載の情報処理装置。
【請求項１３】
前記複数の時系列データは、特性の異なる複数の前記センサそれぞれから得られる時系列データを含む、請求項１１に記載の情報処理装置。
【請求項１４】
コンピュータによって実行させる制御方法であって、
対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得ステップと、
前記時系列データに対して寄与する複数の特徴定数と、各前記特徴定数の前記時系列データに対する寄与の大きさを表す寄与値とを算出する算出ステップと、
複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、前記センサによってセンシングされたガスの特徴量として出力する出力ステップと、を有し、
前記特徴定数は、前記センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である、制御方法。
【請求項１５】
請求項１４に記載の制御方法の各ステップをコンピュータに実行させるプログラム。
【国際調査報告】
2020026327.tif 000041
2020026327.tif 000042
2020026327.tif 000043
2020026327.tif 000044