JP2009207763A - 生体信号関連性計測装置および生体信号関連性計測方法 - Google Patents

Abstract

【課題】信号主要部分の変化に影響されない生体信号の関連性を検出する装置を提供する。
【解決手段】第１生体情報を計測する手段（１０１）と、第２生体情報を計測する手段（１０５）と、第１生体情報から第１生体信号を抽出する手段（１０２）と、第２生体情報から第２生体信号を抽出する手段（１０６）と、第１生体信号を時間周波数分析し第１分析結果を出力する手段（１０３）と、第２生体信号を時間周波数分析し第２分析結果を出力する手段（１０７）と、第１生体信号に対する第１時間周波数範囲を算出する手段（１０４）と、第２生体信号に対する第２時間周波数範囲を算出する手段（１０８）と、第１時間周波数範囲と第２時間周波数範囲に基づいて第３時間周波数範囲を決定する手段（１０９）と、第１分析結果と第２分析結果と第３時間周波数範囲に基づいて第１生体信号と第２生体信号の関連性を算出する手段（１１０）とを備える。
【選択図】図１

Description

本発明は、生体信号の関連性を計測する装置および方法に関する。

近年、高齢化社会や代謝症候群（メタボリックシンドローム）の問題で、健康の自己管理の重要性が注目を浴びている。健康の自己管理の重要性に関して最も注目されているのは、生活習慣病である。生活習慣病は、生活習慣を変えれば回避や早期回復が可能といわれている。自分の生活習慣に問題があることに気づかせるには、健康状態を短時間で簡単に評価する方法が一番有効だと考えられる。そこで、健康状態の指標として自律神経機能を評価する方法が知られている。

非特許文献１では、血行循環に現れる低周波数ゆらぎに基づいて自律神経機能を評価する方法が提案されている。

体内において、血圧のゆらぎを制御するために、フィードバックループに基づいて自律神経などの働きで心拍数を操作することが常に行われている。体内のシステムにおける遅延のため、血行循環にはマイヤー波と呼ばれる周期約１０秒のゆらぎが現れる。非特許文献１で提案されている方法は、生体状況が急変すると血圧と心拍数のマイヤー波の相互関係が変化するという仮説に基づいている。この方法では、具体的にマイヤー波のゆらぎを心拍数変動と血圧変動から抽出し、相互関係係数を算出し、その相互関係係数の最大値を時系列データとして蓄積する。蓄積された最大相互相関係数の時系列データを用いて自律神経機能の評価を行う。最大相互相関係数が１から離れれば離れるほど、自律神経機能の効率が低下していると評価することができる。

この方法を、非特許文献１の記載に従って、以下、ρｍａｘ法と記す。また、最大相互相関係数を、非特許文献１の記載に従って、以下、ρｍａｘと記す。

安静時にはフィードバックコントロールが乱れなく機能するため、ρｍａｘは高い値になる。図１３は、生体の安静時の生理指標に関する特性の一例を示す図である。図１３（ａ）は、安静時の心拍数のマイヤー波成分と血圧関連情報（脈波拍内積分値）のマイヤー波成分を示すグラフであり、図１３（ｂ）は、遅延を用いた相互相関係数を示すグラフである。相互相関係数の最大値がρｍａｘである。外部刺激などが加えられると、その影響で自律神経機能が崩れ、制御が乱れるためρｍａｘが低下すると考えられている。外部刺激とは、ストレサ（ストレス要因）による刺激である。非特許文献１では、映像酔いを起こす手ぶれ画像を用いてρｍａｘを検証している。図１４は、生体の有刺激時の生理指標に関する特性の一例を示す図である。図１４（ａ）は、外部刺激があった場合の心拍数のマイヤー波成分と血圧関連情報（脈波拍内積分値）のマイヤー波成分を示すグラフであり、図１４（ｂ）は、遅延を用いた相互相関係数を示すグラフである。

以下、図面を参照しながらρｍａｘ法について説明する。
図１５は、ρｍａｘを用いた従来の生体状況評価装置の構成を示すブロック図である。従来の生体情報評価装置２００は、第１の生体計測部２０１と、第１の生体信号抽出部２０２と、第１の周波数成分抽出部２０３と、第２の生体計測部２０４と、第２の生体信号抽出部２０５と、第２の周波数成分抽出部２０６と、相関算出部２０７とを備える。第１の生体計測部２０１は、被験者の第１の生体情報を計測し、第１の生体信号抽出部２０２は、前記第１の生体情報から第１の生体信号を抽出し、第１周波数成分抽出部２０３は、前記第１の生体信号から第１の特定の周波数成分を抽出する。同様に、第２の生体計測部２０４は、被験者の第２の生体情報を計測し、第２の生体信号抽出部２０５は、前記第２の生体情報から第２の生体信号を抽出し、第２周波数成分抽出部２０６は、前記第２の生体信号から第２の特定の周波数成分を抽出する。相関算出部２０７は、前記第１の特定の周波数成分と前記第２の特定の周波数成分の相関値を計算し、出力する。このようにして得られた相関値を用いて、２つの生体信号（前記第１の生体信号および前記第２の生体信号）の関連性を評価する。

図１５に示す従来の生体状況評価装置の動作を、図１６のフローチャートを用いて説明する。例として、第１の生体計測部２０１を脈波測定装置とし、第２の生体計測部２０４を連続血圧計とした場合について説明する。

ステップＳ２０１において、脈波測定装置と連続血圧計でそれぞれ脈波と血圧を測定し、時系列データとして脈波と血圧の波形を記録する。ステップＳ２０２において、ユーザや生体状況評価装置のシステム自体が測定の終了を求めているかどうかを確認する。測定の終了が求められている場合はステップＳ２０３に進み、そうでない場合はステップＳ２０１に戻り、生体計測を継続する。

ステップＳ２０３において、ステップＳ２０２で記録した脈波の波形データから立ち上がり点を検出し、各脈拍の立ち上がり時間と終了時間を記録する。

ステップＳ２０４において、心拍数の変動の時系列データを算出し、記録する。心拍数の算出方法としては、例えば、ステップＳ２０３で記録した各脈拍の立ち上がり時間と終了時間から差分を算出し、逆数を求める方法がある。

ステップＳ２０５において、ステップＳ２０４で算出した心拍数の時系列データからマイヤー波を抽出し、記録する。一般的に、マイヤー波の抽出には、バンドパスフィルタが使用されている。フィルタ処理の方法によっては、ステップＳ２０４で算出した心拍数変動の時系列データを補間する必要がある。例えば、２Ｈｚの周波数でスプライン補間を行うことが望ましい。

ステップＳ２０６において、血圧の変動の時系列データを算出し、記録する。血圧の変動の時系列データの算出方法としては、例えば、ステップＳ２０３で記録した各脈拍の立ち上がり時間と終了時間に基づいて、ステップＳ２０１で記録した血圧の波形を脈拍ごとに分解してピークを求める方法が知られている。

ステップＳ２０７において、ステップＳ２０６で算出した血圧の変動の時系列データからマイヤー波を抽出し、記録する。一般的に、マイヤー波の抽出には、バンドパスフィルタが使用されている。フィルタ処理の方法によっては、ステップＳ２０６で算出した血圧変動の時系列データを補間する必要がある。例えば、２Ｈｚの周波数でスプライン補間を行うことが望ましい。

ステップＳ２０８において、ρｍａｘを算出し、記録する。ρｍａｘは以下の手順で算出する。ステップＳ２０５で抽出した心拍数変動におけるマイヤー波の時系列データと、ステップＳ２０７で抽出した血圧変動におけるマイヤー波の時系列データとから、相互相関係数を算出し、最大相関係数の値を求める。この値をρｍａｘとする。

ステップＳ２０９において、ステップＳ２０８で記録したρｍａｘの時系列データを出力し、動作を終了する。

杉田典大、吉澤誠、田中明、阿部健一、山家智之、仁田新一「血圧−心拍数間の最大相互相関係数を用いた映像刺激の生体影響評価」ヒューマンインタフェース学会論文誌、Ｖｏｌ．４、Ｎｏ．４、ｐｐ．２２７−２３４（２００２年１１月）

しかしながら、上述した非特許文献１に開示されている従来の方法には以下のような問題点がある。この従来の方法では、信号主要部分の周波数が、周波数抽出部が抽出する周波数帯を外れた場合、ノイズの相関を計測する結果となるため、一般に相関値が低下し、生体信号間の関連性（相関）を正しく検出することができなくなる。これを解決するためには、周波数成分抽出部が抽出する周波数帯を拡大することが考えられる。しかし、周波数帯を拡大すると、高周波数にある呼吸の周波数または極低周波数にある血管運動活動を含む可能性が高まり、正しい評価ができなくなるおそれがある。また、広い周波数帯を抽出すると、生体信号関連性計測の感度が低下し、評価が困難になる。

本発明は、このような実情を鑑みてなされたものであり、信号主要部分の変化に影響されない生体信号の関連性を検出する装置および方法を提供する。

本発明の生体信号関連性検出装置は、第１の生体情報を計測する第１の生体情報計測手段と、第２の生体情報を計測する第２の生体情報計測手段と、前記第１の生体情報から第１の生体信号を抽出する第１の生体信号抽出手段と、前記第２の生体情報から第２の生体信号を抽出する第２の生体信号抽出手段と、前記第１の生体信号を時間周波数分析し、第１の分析結果を出力する第１の時間周波数分析手段と、前記第２の生体信号を時間周波数分析し、第２の分析結果を出力する第２の時間周波数分析手段と、前記第１の生体信号に対する第１の時間周波数範囲を算出する第１の時間周波数範囲算出手段と、前記第２の生体信号に対する第２の時間周波数範囲を算出する第２の時間周波数範囲算出手段と、前記第１の時間周波数範囲と前記第２の時間周波数範囲に基づいて第３の時間周波数範囲を決定する時間周波数範囲決定手段と、前記第１の分析結果と前記第２の分析結果と前記第３の時間周波数範囲とに基づいて前記第１の生体信号と前記第２の生体信号の関連性を算出する関連性算出手段とを備えることを特徴とする。

前記第１の生体信号と前記第２の生体信号のうち一方は心拍数であってもよく、他方は血圧または血圧に関連する情報であってもよい。

前記第１の時間周波数範囲算出手段は、各時刻における周波数パワースペクトルにおける予め設定された周波数範囲内で前記第１の生体信号が最大値となる周波数に基づいて前記第１の時間周波数範囲を算出してもよい。

前記第２の時間周波数範囲算出手段は、各時刻における周波数パワースペクトルにおける予め設定された周波数範囲内で前記第２の生体信号が最大値となる周波数に基づいて前記第２の時間周波数範囲を算出してもよい。

前記第１の時間周波数範囲算出手段は、前記第１の生体信号の時間周波数パワースペクトルの勾配に基づいて前記第１の時間周波数範囲を算出してもよい。
前記第２の時間周波数範囲算出手段は、前記第２の生体信号の時間周波数パワースペクトルの勾配に基づいて前記第２の時間周波数範囲を算出してもよい。

前記第１の時間周波数範囲算出手段は、さらに時間周波数領域のしきい値関数にも基づいて前記第１の時間周波数範囲を算出してもよい。
前記第２の時間周波数範囲算出手段は、さらに時間周波数領域のしきい値関数にも基づいて前記第２の時間周波数範囲を算出してもよい。

前記時間周波数領域のしきい値関数は、前記第１の生体信号の各時刻の周波数パワースペクトルにおける最大値の割合関数であってもよい。
前記時間周波数領域のしきい値関数は、前記第２の生体信号の各時刻の周波数パワースペクトルにおける最大値の割合関数であってもよい。

前記時間周波数領域のしきい値関数として分布の逆数関数を用いてもよい。
前記時間周波数領域のしきい値関数は、前記第１の生体信号の時間周波数パワースペクトルにおける特定の主成分に依存する関数であってもよい。

前記時間周波数領域のしきい値関数は、前記第２の生体信号の時間周波数パワースペクトルにおける特定の主成分に依存する関数であってもよい。
前記時間周波数領域のしきい値関数は、過去の前記第１生体信号の時間周波数パワースペクトルにおける周波数の最大値に依存する関数であってもよい。
前記時間周波数領域のしきい値関数は、過去の前記第２生体信号の時間周波数パワースペクトルにおける周波数の最大値に依存する関数であってもよい。

前記時間周波数範囲決定手段は、前記第１の時間周波数範囲と前記第２の時間周波数範囲の論理和に基づいて第３の時間周波数範囲を決定してもよい。
前記時間周波数範囲決定手段は、前記第１の時間周波数範囲と前記第２の時間周波数範囲の論理積に基づいて第３の時間周波数範囲を決定してもよい。
前記時間周波数範囲決定手段は、前記第１の時間周波数範囲と前記第２の時間周波数範囲の排他的論理和に基づいて第３の時間周波数範囲を決定してもよい。

前記関連性算出手段は、時間周波数空間コヒーレンシーを用いて前記関連性を算出してもよい。
前記関連性算出手段は、位相同期を評価する指標を用いて前記関連性を算出してもよい。

前記関連性算出手段は、時間周波数空間コヒーレンシーの各時刻における予め設定された周波数範囲内における最大値を用いて前記関連性を算出してもよい。
前記関連性算出手段は、時間周波数空間コヒーレンシーの各時刻における予め設定された周波数範囲内における平均値を用いて前記関連性を算出してもよい。
前記関連性算出手段は、各時刻における予め設定された周波数範囲内で時間周波数空間コヒーレンシーが最大値となる周波数を用いて前記関連性を算出してもよい。

前記関連性算出手段は、位相同期を評価する指標の各時刻における予め設定された周波数範囲内における最大値を用いて前記関連性を算出してもよい。
前記関連性算出手段は、位相同期を評価する指標の各時刻における予め設定された周波数範囲内における平均値を用いて前記関連性を算出してもよい。
前記関連性算出手段は、各時刻における予め設定された周波数範囲内で位相同期を評価する指標が最大値となる周波数を用いて前記関連性を算出してもよい。

本発明によれば、信号主要部分の変化に影響されずに生体信号の関連性を検出することができる。

図１は、本発明の実施の形態における生体信号関連性計測装置の構成を示すブロック図である。図１において、生体信号関連性検出装置１００は、第１の生体計測部１０１と、第１の生体信号抽出部１０２と、第１の時間周波数分析部１０３と、第１の時間周波数範囲算出部１０４と、第２の生体計測部１０５と、第２の生体信号抽出部１０６と、第２の時間周波数分析部１０７と、第２の時間周波数範囲算出部１０８と、時間周波数範囲決定部１０９と、関連性算出部１１０とを備える。
生体信号関連性検出装置１００は、以下のように２つの生体信号の関連性を算出する。

第１の生体計測部１０１は、被験者の第１の生体情報を計測する。第１の生体情報とは、例えば、脈波や血圧などである。第１の生体信号抽出部１０２は、第１の生体計測部１０１が計測した第１の生体情報から第１の生体信号を抽出する。第１の時間周波数分析部１０３は、第１の生体信号抽出部１０２が抽出した第１の生体信号の各時間における第１の周波数スペクトル分析関数を求める。第１の時間周波数範囲算出部１０４は、第１の生体信号抽出部１０２が抽出した第１の生体信号の瞬時パワースペクトルを求め、第１の生体信号の主要成分の周波数追跡を行い、第１の時間周波数範囲を算出する。

第２の生体計測部１０５は、被験者の第２の生体情報を計測する。第２の生体情報とは、例えば、脈波や血圧などであり、前記第１の生体情報と異なるものである。第２の生体信号抽出部１０６は、第２の生体計測部１０５が計測した第２の生体情報から第２の生体信号を抽出する。第２の時間周波数分析部１０７は、第２の生体信号抽出部１０６が抽出した第２の生体信号の各時間における第２の周波数スペクトル分析関数を求める。第２の時間周波数範囲算出部１０８は、第２の生体信号抽出部１０６が抽出した第２の生体信号の瞬時パワースペクトルを求め、第２の生体信号の主要成分の周波数追跡を行い、第２の時間周波数範囲を算出する。

時間周波数範囲決定部１０９は、第１の時間周波数範囲算出部１０４が算出した第１の時間周波数範囲と、第２の時間周波数範囲算出部１０８が算出した第２の時間周波数範囲とに基づいて、第３の時間周波数範囲を決定する。

関連性算出部１１０は、第１の時間周波数分析部１０３が求めた第１の時間周波数スペクトル分析関数と、第２の時間周波数分析部１０７が求めた第２の周波数スペクトル分析関数と、時間周波数範囲決定部１０９が決定した第３の時間周波数範囲とに基づいて、前記第１の生体信号と前記第２の生体信号の関連性を算出する。

本発明の実施の形態における生体信号関連性計測装置のより具体的な動作について、第１の生体計測部１０１に脈波測定装置、第２の生体計測部１０５に連続血圧計を用いた場合を例としてより詳細に説明する。

第１の生体計測部１０１に用いる脈波測定装置としては、被験者の指先などに装着して使用する光電式指先容量脈波計が一般的である。図２は光電式指先容量脈波計の一般的な構成を示す概略図である。光電式指先容量脈波計１０は、被験者の指先を挿入する凹部の内側に発光部１１と受光部１２が配置されている。発光部１１が発光した光は被験者の指先を経て受光部１２によって受光される。光電式指先容量脈波計は、血流量に応じて光の透過度が変化することを利用し、脈波を測定する。図３は、光電式指先容量脈波計によって得られる脈波の一例を示すグラフである。

第１の生体計測部（脈波測定装置）１０１が測定した脈波の波形を第１の生体信号抽出部１０２に供給する。第１の生体信号抽出部１０２は、脈波の波形から脈拍頻度を測定し、心拍数を抽出する。図４は、脈波の波形から脈拍頻度を測定する方法を説明するグラフである。脈拍間隔とは、脈波の立ち上がり点の間の間隔、すなわち脈拍の時間である。ここでの立ち上がり点とは、脈波速度の最大値となる点の近傍における脈波波形の極小値である。心電図の波形ならば、左心室の興奮により一番大きくスパイク状に現れるＲ波と次のＲ波までの間隔を示すＲ−Ｒ間隔と同等である。

脈拍の時間、すなわち脈拍間隔は、図４に示す脈波波形の終了時間と立ち上がり時間の差をとることによって求まる。なお、立ち上がり時間は、立ち上がり点の時間のことである。また、終了時間は、次の脈拍の立ち上がり時間と同じである。

脈拍間隔（ＲＲ）から下記の式により瞬時心拍数（ＨＲ）が求まる。
HR[bpm]=60/RR[sec]

第１の生体信号抽出部１０２は、このような瞬時心拍数の時系列データを心拍数変動（第１の生体信号）として第１の時間周波数分析部１０３および第１の時間周波数範囲算出部１０４に供給する。

第２の生体計測部１０５に用いる連続血圧計としては、トリメリ法血圧計が一般的である。図５は、トリメリ法血圧計の一般的な構成を示す概略図である。トリメリ法血圧計２０は、制御部２１と、圧力センサアレー２２とを備える。トリメリ法血圧計２０は、橈骨動脈等を圧迫しながら当該動脈上で圧力センサを用い、心拍ごとに圧力波形を取得し、血圧を測定する。図６は、トリメリ法血圧計によって得られる動脈圧力波形の一例を示すグラフである。

第２の生体計測部（連続血圧計）１０５が測定した動脈圧力の波形を第２の生体信号抽出部１０６に供給する。第２の生体信号抽出部１０６は、例えば、図６に示すように動脈圧力の波形からピークを測定し、心臓が血液を送り出すために収縮したときの最高血圧である収縮期血圧として抽出する。第２の生体信号抽出部１０６は、このような収縮期血圧の時系列データを血圧変動（第２の生体信号）として第２の時間周波数分析部１０７および第２の時間周波数範囲算出部１０８に供給する。

第１の時間周波数分析部１０３および第２の時間周波数分析部１０７において、各々、第１の生体信号抽出部１０２が抽出した心拍数変動（第１の生体信号）および第２の生体信号抽出部１０６が抽出した血圧変動（第２の生体信号）の時間周波数分析を行い、各時刻における周波数スペクトル分布関数を求める。周波数スペクトル分布関数を求める場合、例えば、短期間フーリエ変換やウェブレット変換を用いて周波数スペクトルを算出する。高い周波数分解能が得られる、より好ましい方法は、まず、正弦関数や余弦関数を用いて対象の周波数帯を網羅する。次に、窓枠を用いて各時刻のデータに対して前記の正弦関数や余弦関数を用いて重み付き最小二乗回帰を行う。次に、得られた近似関数の重み係数を用いて、前記第１の生体信号および第２の生体信号のパワースペクトルと位相を求める。

第１の時間周波数範囲算出部１０４および第２の時間周波数範囲算出部１０８において、各々、第１の生体信号抽出部１０２が抽出した心拍数変動（第１の生体信号）および第２の生体信号抽出部１０６が抽出した血圧変動（第２の生体信号）に基づいて、第１の時間周波数範囲および第２の時間周波数範囲を算出する。前記第１の生体信号および第２の生体信号の瞬時パワースペクトルを求め、主要成分を追跡するために時間周波数範囲を算出する。

瞬時パワースペクトルを求める方法としては、例えば、窓枠を用いたＹｕｌｅ−Ｗａｌｋｅｒ ＡＲ法やウェブレット変換が知られている。

主要成分の周波数を追跡するための時間周波数範囲の算出方法を以下に具体的に説明する。

まず、周波数範囲を限定することが望ましい。特に、高周波数帯はほとんど呼吸の影響を示しているため、対象範囲に含まない方がよい。安静時の自然な呼吸は０．２５Ｈｚとされている。したがって、周波数範囲の上限を０．２５Ｈｚ以下、例えば０．２０Ｈｚに設定してもよい。

次に、ノイズアーティファクトの影響を防ぐため、各生体信号の瞬時パワースペクトルに時間領域ローパスフィルタ処理を行う。例えば、時間領域移動平均化を適用する。

そして、各生体信号の瞬時パワースペクトルにおける主要成分を強調するために、各時刻においてパワースペクトルを下記のように正規化する。
Sp’(ω,t)=Sp(ω,t)/Max(Sp,t)

ここで、Spは瞬時パワースペクトルであり、Max(Sp,t)はtにおけるSp(ω,t)の最大値を返す関数である。

次に、主要成分を抽出するために各生体信号のピークマスクを生成する。ピークマスクは、例えば、下記の式で求める。
PeakMask(ω,t)=(sign(Sp'(ω,t)-cutoff)+1)/2

ここで、Sp'は正規化された瞬時パワースペクトルであり、cutoffはピークしきい値であり、

であり、

である。

上記の式におけるピークしきい値cutoffを調整することにより、後に説明する関連性算出部１１０が算出する指標の感度を調節することができる。図７は、正規化された瞬間パワースペクトルから生成されるピークマスクのピークしきい値の調整による変化を示す図である。図７の左上に示す正規化瞬時パワースペクトルから、ピークしきい値が０．２５の場合はピークマスクＡが、ピークしきい値が０．５０の場合はピークマスクＢが、ピークしきい値が０．７５の場合はピークマスクＣが生成される。

なおピークしきい値cutoffは、下記のように時間周波数領域の関数を用いてもよい。
PeakMask(ω,t)=(sign(Sp'(ω,t)-cutoff(ω,t))+1)/2

ここで、Sp'は正規化された瞬時パワースペクトルであり、cutoffはピークしきい値であり、

であり、

である。ピークしきい値cutoffに時間周波数領域の関数を用いることによって、様々なピークマスクを生成することができる。例えば、ピークしきい値cutoffの時間周波数領域の関数に分布の逆数関数を用いれば、周波数帯の端の影響を防ぐことができる。

また、呼吸の影響をさらに防ぐために、ピークしきい値cutoffに呼吸の周波数の情報を利用した時間周波数領域の関数を用いることも考えられる。例えば、呼吸の周波数の周辺のcutoffファクタが高くなるような時間周波数領域の関数を用いてもよい。

さらに、ピークマスクにおいて、リッジ（勾配ベクトルに沿った周波数稜線）を抽出したい場合は、ピークの連続性を保つ必要がある。この場合は、例えば、過去のピーク周波数の情報を利用したcutoff時間周波数領域の関数を用いる方法にも適用できる。具体的に、cutoff時間周波数領域の関数としては、過去のピーク周波数を中心にした正規分布の逆数関数が考えられる。

次に、時間周波数範囲決定部１０９は、前記第１の時間周波数範囲と前記第２の時間周波数範囲に基づいて、第３の時間周波数範囲を決定する。第３の時間周波数範囲を決定する方法は、用途に応じて変化する。各生体信号の活動を反映した時間周波数範囲を生成する場合、論理和を用いるとよい。具体的には、第３の時間周波数範囲を以下のように算出する。
PeakMask₁₂(ω,t)=PeakMask₁(ω,t)+PeakMask₂(ω,t)-PeakMask₁(ω,t)×PeakMask₂(ω,t)

ここで、PeakMask₁(ω,t)は第１の時間周波数範囲、PeakMask₂(ω,t)は第２の時間周波数範囲、PeakMask₁₂(ω,t)は第３の時間周波数範囲である。図８は、第３の時間周波数範囲の生成を説明する図である。PeakMask₁(ω,t)とPeakMask₂(ω,t)を重ね合わせることによってPeakMask₁₂(ω,t)が生成される。このようにして第３の時間周波数範囲を決定することにより、生体信号の主成分の追跡を実現し、呼吸や血管運動などの生体ノイズの影響を除去することが可能となる。

関連性算出部１１０は、第１の時間周波数分析部１０３が出力する前記第１の生体信号のパワースペクトルおよび位相と、第２の時間周波数分析部１０７が出力する前記第２の生体信号のパワースペクトルおよび位相と、時間周波数範囲決定部１０９が出力する第３の時間周波数範囲とに基づいて、前記第１の生体信号と前記第２の生体信号の関連性の指標を算出し、出力する。

以下、関連性算出部１１０における関連性指標の算出について具体的に説明する。
まず、第１の時間周波数範囲分析部１０３が出力する前記第１の生体信号のパワースペクトルおよび位相と、第２の時間周波数分析部１０７が出力する前記第２の生体信号のパワースペクトルおよび位相から、前記第１の生体信号と前記第２の生体信号のクロススペクトルを求める。このクロススペクトルは、例えば下記の式によって求められる。

ここで、G₂₁は前記第１の生体信号と前記第２の生体信号のクロススペクトル、Sp₁は前記第１の生体信号のパワースペクトルであり、Sp₂は前記第２の生体信号のパワースペクトルであり、φ₁は前記第１の生体信号の位相であり、φ₂は前記第２の生体信号の位相である。

次に、第１の時間周波数分析部１０３が出力する前記第１の生体信号の位相と、第２の時間周波数分析部１０７が出力する前記第２の生体信号の位相から、前記第１の生体信号と前記第２の生体信号の位相差φ₂₁を求める。具体的には、例えば下記の式のように求める。

次に、前記第１の生体信号のパワースペクトルSp₁と、前記第２の生体信号のパワースペクトルSp₂と、上記のように求めた前記第１の生体信号と前記第２の生体信号のクロススペクトルG₂₁から、前記第１の生体信号と前記第２の生体信号のコヒーレンシーC₂₁を下記の式のように求める。{}_ΔtはΔtの時間帯を対象とした時間領域の平均処理である。

次に、このようにして算出した時間周波数空間におけるコヒーレンシーに、時間周波数範囲決定部１０９が決定した前記第３の時間周波数範囲PeakMask₁₂(ω,t)を適用することにより、前記第１の生体信号と前記第２の生体信号の関連性の正しい評価が可能となる。具体的には、下記の式のように前記第３の時間周波数範囲を適用し、マスクを適用した後の前記第１の生体信号と前記第２の生体信号のコヒーレンシーC₂₁'を得る。
C₂₁'(ω,t)=C₂₁(ω,t)×PeakMask₁₂(ω,t)

被験者に分かりやすい指標を表示するために、時間領域のコヒーレンシー指標を用いるとよい。この時間領域コヒーレンシー指標を得るために、C₂₁'の各時刻における最大値とその周波数を求めるか、平均値を求める。具体的には、例えば下記の式でそれぞれの指標を求める。
Cmax(t)=max_ω(C₂₁'(ω,t))
Ω_C(t)=H(C₂₁'(ω,t))
Cavg(t)={C₂₁'(ω,t)}_Δω

ここで、CmaxはC₂₁'の各時刻における最大値、max_ωは各時刻におけるピークを返す関数、Ω_CはC₂₁'の各時刻における最大値に相当する周波数、Hは各時刻におけるピークに相当する周波数を返す関数、CavgはC₂₁'の各時刻の周波数領域の平均値、{}_ΔωはΔωの周波数帯を対象とした周波数領域の平均処理である。

図９は、相関係数を用いた従来の方法によって得られる指標と、上述した本発明の実施の形態における生体信号関連性計測装置によって得られる指標とを比較するグラフである。下方に示す従来の方法によって得られる指標は非常に不安定である。また、刺激のない期間でも変動していることが分かる。それに対して、上方に示す本発明の生体信号関連性計測装置によって得られる指標は安定している。刺激のあるなし、または、強度に関連した変動を確認することができる。

２つの信号の関連性を示す指標として、コヒーレンシー以外に、位相同期指標を用いてもよい。位相同期指標は下記の式で求められる。

ここで、φ₂₁は前記第１の生体信号と前記第２の生体信号の時間周波数位相差であり、Δφは対象時間帯である。PLV₂₁が１である場合、両信号の位相が完全に同期の状態にあるといえる。PLV₂₁が０である場合、両信号の位相が同期せずに独立して変動していることがいえる。

被験者に分かりやすい指標を表示するために、時間領域の位相同期指標を用いるとよい。この時間領域位相同期指標を得るために、PLV₂₁の各時刻における最大値とその周波数を求めるか、平均値を求める。具体的には、例えば下記の式によってそれぞれの指標を求める。
PLVmax(t)=max_ω(PLV₂₁(ω,t))
Ω_PLV(t)=H(PLV₂₁(ω,t))
PLVavg(t)={PLV₂₁(ω,t)}_Δω

ここで、PLVmaxはPLV₂₁の各時刻における最大値、max_ωは各時刻におけるピークを返す関数、Ω_PLVはPLV₂₁の各時刻における最大値に相当する周波数、Hは各時刻におけるピークに相当する周波数を返す関数、PLVavgはPLV₂₁の各時刻の周波数領域の平均値、{}_ΔωはΔωの周波数帯を対象とした周波数領域の平均処理である。

図１０は、図１に示す生体信号関連性計測装置１００が起動した後に実行されるプログラムの制御動作を示すフローチャートである。このプログラムは、生体情報の計測終了後に生体信号の関連性を算出するような制御構造となっている。

ステップＳ１０１において、第１の生体計測部１０１および第２の生体計測部１０５に被験者の第１および第２の生体情報を測定させる。ステップＳ１０２において、被験者が生体情報の測定の終了を要求しているか否かを判断し、要求している場合はステップＳ１０１に戻って生体情報の測定を継続し、要求している場合はステップＳ１０３に進む。ステップＳ１０３において、第１の生体信号抽出部１０２に、第１の生体計測部１０１が計測した第１の生体情報から第１の生体信号を抽出させ、記録させる。ステップＳ１０４において、第２の生体信号抽出部１０６に、第２の生体計測部１０６が計測した第２の生体情報から第１の生体信号を抽出させ、記録させる。ステップＳ１０５において、第１の時間周波数範囲算出部１０４に、第１の生体信号抽出部１０２が抽出した第１の生体信号に対する第１の時間周波数範囲を算出させる。ステップＳ１０６において、第２の時間周波数範囲算出部１０８に、第２の生体信号抽出部１０６が抽出した第２の生体信号に対する第２の時間周波数範囲を算出させる。ステップＳ１０７において、時間周波数範囲決定部１０９に、第１の時間周波数範囲算出部１０４が算出した第１の時間周波数範囲と第２の時間周波数範囲算出部１０８が算出した第２の時間周波数範囲に基づいて、第３の時間周波数範囲を決定させる。ステップＳ１０８において、第１の時間周波数分析部１０３に、第１の生体信号抽出部１０２が抽出した第１の生体信号の時間周波数分析を行わせ、時間周波数分布関数を算出させ、第２の時間周波数分析部１０７に、第２の生体信号抽出部１０６が抽出した第２の生体信号の時間周波数分析を行わせ、時間周波数分布関数を算出させる。ステップＳ１０９において、関連性算出部１１０に、時間周波数範囲決定部１０９が決定した第３の時間周波数範囲と、第１の時間周波数分析部１０３および第２の時間周波数分析部１０７が算出した時間周波数分布関数とに基づいて、前記第１の生体信号と前記第２の生体信号の関連性を算出させ、内部メモリに記録させる。ステップＳ１１０において、関連性算出部１１０の内部メモリに記録された生体信号関連性のデータを出力させ、動作を終了する。

以上、本発明の実施の形態を説明してきたが、本発明はこのような実施の形態のみに限定されるものではなく、本発明の要旨を逸脱しない範囲における変更や追加があっても本発明に含まれる。

例えば、上記実施の形態では光電容量脈波センサを用いて心拍数変動を抽出したが、心電図に基づいて心拍数変動を抽出してもよい。

また、上記実施の形態では連続血圧計を用いて血圧変動を抽出したが、光電容量脈波センサを用い、脈波信号の特徴に基づいて血圧に関連する情報を抽出してもよい。例えば、図３に示す拍内積分値を算出する方法を用いてもよい。

また、上記実施の形態では２つの生体計測部（第１の生体計測部１０１および第２の生体計測部１０５）を用いたが、１つの生体計測部から複数の生体信号を抽出してもよい。

また、上記実施の形態では、第１の生体信号と第２の生体信号の活動を反映した第３の時間周波数範囲を用いたが、両生体信号の活動が合う時間周波数範囲を用いてもよい。図１１は、このような第３の時間周波数範囲の生成を説明する図である。第１の生体信号に対する第１の時間周波数範囲PeakMask₁と第２の生体信号に対する第２の時間周波数範囲PeakMask₂の共通部分のみが残されたマスクを第３の時間周波数範囲PeakMask₁₂として生成する。具体的には、下記の式を用いればよい。
PeakMask₁₂=PeakMask₁×PeakMask₂

また、上記実施の形態では、第１の生体信号と第２の生体信号の活動を反映した第３の時間周波数範囲を用いたが、各生体信号の独立した活動を反映した第３の時間周波数範囲を用いてもよい。このような第３の時間周波数範囲を生成する場合は、例えば排他的論理和を用いる。図１２は、このような第３の時間周波数範囲の生成を説明する図である。第１の生体信号に対する第１の時間周波数範囲PeakMask₁と第２の生体信号に対する第２の時間周波数範囲PeakMask₂とを重ねて、共通の部分を削除したマスクを第３の時間周波数範囲PeakMask₁₂として生成する。具体的には、下記の式を用いればよい。
PeakMask₁₂=PeakMask₁＋PeakMask₂-2×(PeakMask₁×PeakMask₂)

また、上記実施の形態では、第１の生体信号と第２の生体信号の活動を反映した第３の時間周波数範囲を用いたが、片方の生体信号の活動のみを反映した第３の時間周波数範囲を用いてもよい。具体的には、時間周波数範囲決定部において、第１の時間周波数範囲または第２の時間周波数範囲のいずれか一方を選択し、出力する。

また、上記実施の形態では、第１の時間周波数範囲算出部および第２の時間周波数範囲算出部においてピークマスクを用いて時間周波数範囲を算出したが、リッジを用いる検出アルゴリズムによって時間周波数範囲を算出してもよい。例えば、時間周波数の勾配ベクトルを用いる検出アルゴリズムを使用してもよい。

また、上記実施の形態では、関連性指標としてコヒーレンシーや位相同期指標を用いたが、他の関連性指標や生体評価指標を用いてもよい。例えば、交感神経の活動を示す低周波数ＬＦ（一般的に、０．０５Ｈｚ〜０．１５Ｈｚ）と副交感神経の活動を示す高周波数ＨＦ（一般的に、０．１５Ｈｚ〜０．４Ｈｚ）の割合ＬＦ／ＨＦを用いてもよい。

また、上記実施の形態の生体信号関連性計測装置では、生体情報の計測後に生体信号関連性を算出するような制御構造のプログラムを用いたが、測定される生体情報をリアルタイムに分析するような制御構造のプログラムを用いてもよい。例えば、生体情報から新たな生体信号のデータが抽出されるときに生体信号を同時処理するスレッドを起動させ、関連性指標を更新するような制御構造のプログラムを用いてもよい。

以上説明したように、本発明によれば、従来の相関係数を用いた方法とは異なり、生体信号の活動に変化が生じても生体信号の関連性を正しく算出することができる。また、呼吸に影響されることなく、生体信号の関連性を算出することができる。

本発明の生体信号関連性検出装置は、生体信号の関連性を得ることができること、呼吸の変動の影響を受けないことから、様々な利用が可能であると考えられる。具体的には、自律神経活動の評価ができるので、健康管理モニターとしての利用が考えられる。さらに、呼吸の変動に影響されないため、呼吸が激しくなるような場面でも利用できると考えられる。例えば、生体情報のフィードバックに基づいたスポーツ選手の練習に利用することができる。さらに、小型化することで、いつでもどこでも自分の体調を確認することができるようになる。

本発明は、生体信号関連性検出装置に利用可能である。

本発明の実施の形態における生体信号関連性計測装置の構成を示すブロック図である。 光電式指先容量脈波計の一般的な構成を示す概略図である。 光電式指先容量脈波計によって得られる脈波の一例を示すグラフである。 脈波の波形から脈拍頻度を測定する方法を説明するグラフである。 トリメリ法血圧計の一般的な構成を示す概略図である。 トリメリ法血圧計によって得られる動脈圧力波形の一例を示すグラフである。 正規化された瞬間パワースペクトルから生成されるピークマスクのピークしきい値の調整による変化を示す図である。 第３の時間周波数範囲の生成を説明する図である。 相関係数を用いた従来の方法によって得られる指標と、本発明の実施の形態における生体信号関連性計測装置によって得られる指標とを比較するグラフである。 図１に示す生体信号関連性計測装置１００が起動した後に実行されるプログラムの制御動作を示すフローチャートである。 第３の時間周波数範囲の生成を説明する図である。 第３の時間周波数範囲の生成を説明する図である。 生体の安静時の生理指標に関する特性の一例を示す図であり、（ａ）は安静時の心拍数のマイヤー波成分と血圧関連情報（脈波拍内積分値）のマイヤー波成分を示すグラフであり、（ｂ）は遅延を用いた相互相関係数を示すグラフである。 生体の有刺激時の生理指標に関する特性の一例を示す図であり、（ａ）は外部刺激があった場合の心拍数のマイヤー波成分と血圧関連情報（脈波拍内積分値）のマイヤー波成分を示すグラフであり、（ｂ）は遅延を用いた相互相関係数を示すグラフである。 ρｍａｘを用いた従来の生体状況評価装置の構成を示すブロック図である。 ρｍａｘを用いた従来の生体状況評価装置の処理手順を示すフローチャートである。

符号の説明

１０ 光電式指先容量脈波計
１１ 発光部
１２ 受光部
２０ トノメトリ法血圧計
２１ 制御部
２２ 圧力センサアレー
１００ 生体信号関連性検出装置
１０１ 第１の生体計測部
１０２ 第１の生体信号抽出部
１０３ 第１の時間周波数分析部
１０４ 第１の時間周波数範囲算出部
１０５ 第２の生体計測部
１０６ 第２の生体信号抽出部
１０７ 第２の時間周波数分析部
１０８ 第２の時間周波数範囲算出部
１０９ 時間周波数範囲決定部
１１０ 関連性算出部
２０１ 第１の生体計測部
２０２ 第１の生体信号抽出部
２０３ 第１の時間周波数成分抽出部
２０４ 第２の生体計測部
２０５ 第２の生体信号抽出部
２０６ 第２の時間周波数成分抽出部
２０７ 相関算出部

Claims (28)

1. 第１の生体情報を計測する第１の生体情報計測手段と、
   第２の生体情報を計測する第２の生体情報計測手段と、
   前記第１の生体情報から第１の生体信号を抽出する第１の生体信号抽出手段と、
   前記第２の生体情報から第２の生体信号を抽出する第２の生体信号抽出手段と、
   前記第１の生体信号を時間周波数分析し、第１の分析結果を出力する第１の時間周波数分析手段と、
   前記第２の生体信号を時間周波数分析し、第２の分析結果を出力する第２の時間周波数分析手段と、
   前記第１の生体信号に対する第１の時間周波数範囲を算出する第１の時間周波数範囲算出手段と、
   前記第２の生体信号に対する第２の時間周波数範囲を算出する第２の時間周波数範囲算出手段と、
   前記第１の時間周波数範囲と前記第２の時間周波数範囲に基づいて第３の時間周波数範囲を決定する時間周波数範囲決定手段と、
   前記第１の分析結果と前記第２の分析結果と前記第３の時間周波数範囲とに基づいて前記第１の生体信号と前記第２の生体信号の関連性を算出する関連性算出手段とを備えることを特徴とする生体信号関連性計測装置。
2. 前記第１の生体信号と前記第２の生体信号のうち一方は心拍数であり、他方は血圧または血圧に関連する情報であることを特徴とする請求項１記載の生体信号関連性計測装置。
3. 前記第１の時間周波数範囲算出手段は、各時刻における周波数パワースペクトルにおける予め設定された周波数範囲内で前記第１の生体信号が最大値となる周波数に基づいて前記第１の時間周波数範囲を算出することを特徴とする請求項１または２に記載の生体信号関連性計測装置。
4. 前記第２の時間周波数範囲算出手段は、各時刻における周波数パワースペクトルにおける予め設定された周波数範囲内で前記第２の生体信号が最大値となる周波数に基づいて前記第２の時間周波数範囲を算出することを特徴とする請求項１乃至３のいずれか１項に記載の生体信号関連性計測装置。
5. 前記第１の時間周波数範囲算出手段は、前記第１の生体信号の時間周波数パワースペクトルの勾配に基づいて前記第１の時間周波数範囲を算出することを特徴とする請求項１、２および４のうちいずれか１項に記載の生体信号関連性計測装置。
6. 前記第２の時間周波数範囲算出手段は、前記第２の生体信号の時間周波数パワースペクトルの勾配に基づいて前記第２の時間周波数範囲を算出することを特徴とする請求項１、２、３および５のうちいずれか１項に記載の生体信号関連性計測装置。
7. 前記第１の時間周波数範囲算出手段は、さらに時間周波数領域のしきい値関数にも基づいて前記第１の時間周波数範囲を算出することを特徴とする請求項３又は５に記載の生体信号関連性計測装置。
8. 前記第２の時間周波数範囲算出手段は、さらに時間周波数領域のしきい値関数にも基づいて前記第２の時間周波数範囲を算出することを特徴とする請求項４又は６に記載の生体信号関連性計測装置。
9. 前記時間周波数領域のしきい値関数は、前記第１の生体信号の各時刻の周波数パワースペクトルにおける最大値の割合関数であることを特徴とする請求項７記載の生体信号関連性計測装置。
10. 前記時間周波数領域のしきい値関数は、前記第２の生体信号の各時刻の周波数パワースペクトルにおける最大値の割合関数であることを特徴とする請求項８記載の生体信号関連性計測装置。
11. 前記時間周波数領域のしきい値関数として分布の逆数関数を用いることを特徴とする請求項７に記載の生体信号関連性計測装置。
12. 前記時間周波数領域のしきい値関数として分布の逆数関数を用いることを特徴とする請求項８に記載の生体信号関連性計測装置。
13. 前記時間周波数領域のしきい値関数は、前記第１の生体信号の時間周波数パワースペクトルにおける特定の主成分に依存する関数であることを特徴とする請求項７に記載の生体信号関連性計測装置。
14. 前記時間周波数領域のしきい値関数は、前記第２の生体信号の時間周波数パワースペクトルにおける特定の主成分に依存する関数であることを特徴とする請求項８に記載の生体信号関連性計測装置。
15. 前記時間周波数領域のしきい値関数は、過去の前記第１の生体信号の時間周波数パワースペクトルにおける周波数の最大値に依存する関数であることを特徴とする請求項７に記載の生体信号関連性計測装置。
16. 前記時間周波数領域のしきい値関数は、過去の前記第２の生体信号の時間周波数パワースペクトルにおける周波数の最大値に依存する関数であることを特徴とする請求項８に記載の生体信号関連性計測装置。
17. 前記時間周波数範囲決定手段は、前記第１の時間周波数範囲と前記第２の時間周波数範囲の論理和に基づいて第３の時間周波数範囲を決定することを特徴とする請求項１乃至１６のいずれか１項に記載の生体信号関連性計測装置。
18. 前記時間周波数範囲決定手段は、前記第１の時間周波数範囲と前記第２の時間周波数範囲の論理積に基づいて第３の時間周波数範囲を決定することを特徴とする請求項１乃至１６のいずれか１項に記載の生体信号関連性計測装置。
19. 前記時間周波数範囲決定手段は、前記第１の時間周波数範囲と前記第２の時間周波数範囲の排他的論理和に基づいて第３の時間周波数範囲を決定することを特徴とする請求項１乃至１６のいずれか１項に記載の生体信号関連性計測装置。
20. 前記関連性算出手段は、時間周波数空間コヒーレンシーを用いて前記関連性を算出することを特徴とする請求項１乃至１９のいずれか１項に記載の生体信号関連性計測装置。
21. 前記関連性算出手段は、位相同期を評価する指標を用いて前記関連性を算出することを特徴とする請求項１乃至１９のいずれか１項に記載の生体信号関連性計測装置。
22. 前記関連性算出手段は、時間周波数空間コヒーレンシーの各時刻における予め設定された周波数範囲内における最大値を用いて前記関連性を算出することを特徴とする請求項２０に記載の生体信号関連性計測装置。
23. 前記関連性算出手段は、時間周波数空間コヒーレンシーの各時刻における予め設定された周波数範囲内における平均値を用いて前記関連性を算出することを特徴とする請求項２０に記載の生体信号関連性計測装置。
24. 前記関連性算出手段は、各時刻における予め設定された周波数範囲内で時間周波数空間コヒーレンシーが最大値となる周波数を用いて前記関連性を算出することを特徴とする請求項２０に記載の生体信号関連性計測装置。
25. 前記関連性算出手段は、位相同期を評価する指標の各時刻における予め設定された周波数範囲内における最大値を用いて前記関連性を算出することを特徴とする請求項２１に記載の生体信号関連性計測装置。
26. 前記関連性算出手段は、位相同期を評価する指標の各時刻における予め設定された周波数範囲内における平均値を用いて前記関連性を算出することを特徴とする請求項２１に記載の生体信号関連性計測装置。
27. 前記関連性算出手段は、各時刻における予め設定された周波数範囲内で位相同期を評価する指標が最大値となる周波数を用いて前記関連性を算出することを特徴とする請求項２１に記載の生体信号関連性計測装置。
28. 第１の生体情報を計測するステップと、
    第２の生体情報を計測するステップと、
    前記第１の生体情報から第１の生体信号を抽出するステップと、
    前記第２の生体情報から第２の生体信号を抽出するステップと、
    前記第１の生体信号を時間周波数分析し、第１の分析結果を出力するステップと、
    前記第２の生体信号を時間周波数分析し、第２の分析結果を出力するステップと、
    前記第１の生体信号に対する第１の時間周波数範囲を算出するステップと、
    前記第２の生体信号に対する第２の時間周波数範囲を算出するステップと、
    前記第１の時間周波数範囲と前記第２の時間周波数範囲に基づいて第３の時間周波数範囲を決定するステップと、
    前記第１の分析結果と前記第２の分析結果と前記第３の時間周波数範囲とに基づいて前記第１の生体信号と前記第２の生体信号の関連性を算出するステップとを含むことを特徴とする生体信号関連性計測方法。

Images