【解決手段】生体信号解析装置1は、生体信号を取得する生体信号取得部211と、生体信号取得部211で取得された複数の生体信号のゆらぎを表す1次元の特徴パラメータを算出する演算部213と、を備える。また、演算部213は、複数の生体信号を混合した尺度混合モデルを用いて、特徴パラメータを算出する。
【発明の概要】
【発明が解決しようとする課題】
【0004】
特許文献1の診断支援装置では、複数の部位で測定された脳波を個々に分析して、判別情報を生成する。したがって、個々の脳波に含まれるノイズの影響を受けやすく、判定精度は低いものとなる。
【0005】
また、特許文献1の診断支援装置では脳波を周波数解析して被験者の状態を判断する。特許文献1のように、生体信号の周波数又は振幅の大きさに着目した解析方法では、被験者の身体的特徴による個体差や生体信号に含まれるノイズの影響を受けやすいため解析精度が低い。
【0006】
本発明は、上述の事情に鑑みてなされたものであり、複数の生体信号を混合して被験者の身体の状態を分析することにより、判定精度の高い生体信号解析装置及び生体信号解析方法を提供することを目的とする。
【課題を解決するための手段】
【0007】
上記目的を達成するために、この発明の第1の観点に係る生体信号解析装置は、
生体信号を取得する生体信号取得部と、
前記生体信号取得部で取得された複数の生体信号のゆらぎを表す1次元の特徴パラメータを算出する演算部と、を備え、
前記演算部は、
前記複数の生体信号を混合した尺度混合モデルを用いて、前記特徴パラメータを算出する。
【0008】
また、前記尺度混合モデルは、多次元分散分布モデルであり、
前記特徴パラメータは、前記多次元分散分布モデルの自由度を含む、
こととしてもよい。
【0009】
また、前記多次元分散分布モデルを表す分散共分散行列の分布として、逆ウィシャート分布を用いる、
こととしてもよい。
【0010】
また、前記生体信号は脳波である、
こととしてもよい。
【0011】
また、前記生体信号に基づいて、異常脳波を伴う疾患を検出する、
こととしてもよい。
【0012】
また、前記生体信号を、複数の周波数フィルタを用いて、所定の周波数帯域の信号に分解するフィルタ部を備え、
前記演算部は、
各周波数帯域に分解された前記生体信号について特徴パラメータを算出する、
こととしてもよい。
【0013】
また、前記演算部は、
所定の時間幅を有する窓を設定し、窓の開始時刻を所定のずらし幅でずらしながら、設定された窓の時間幅の前記生体信号について特徴パラメータを算出する、
こととしてもよい。
【0014】
また、本発明の第2の観点に係る生体信号解析方法は、
複数の生体信号を取得するステップと、
取得された前記複数の生体信号のゆらぎを表す1次元の特徴パラメータを算出するステップと、を含み、
前記特徴パラメータは、前記複数の生体信号を混合した尺度混合モデルを用いて算出される。
【発明の効果】
【0015】
本発明の生体信号解析装置及び生体信号解析方法によれば、複数の生体信号を混合して、生体信号のゆらぎを表す1次元パラメータを算出するので、容易にかつ精度よく被験者の状態変化を分析することが可能である。
【発明を実施するための形態】
【0017】
以下、図を参照しつつ、本発明の実施の形態に係る生体信号解析装置1について説明する。
図1のブロック図に示すように、生体信号解析装置1は、センサ11、制御ユニット20を備える。
【0018】
センサ11は、取得する生体信号の種類に合わせて選択される検出器である。取得する生体信号は、例えば脳波、脈拍等であり、特に限定されない。
【0019】
制御ユニット20は、例えばコンピュータ装置であり、
図1に示すように、制御部21、記憶部22、表示部23、入力部24を備える。
【0020】
制御部21は、CPU(Central Processing Unit)、ROM(Read Only Memory)、RAM(Random Access Memory)、水晶発振器等から構成されており、生体信号解析装置1の動作を制御するとともに、センサ11で取得された生体信号に基づいて特徴パラメータを算出し、被験者の身体状態を判定する。制御部21は、制御部21のROM、記憶部22等に記憶されている各種動作プログラム及びデータをRAMに読み込んでCPUを動作させることにより、
図1に示される制御部21の各機能を実現させる。これにより、制御部21は、生体信号取得部211、フィルタ部212、演算部213として動作する。
【0021】
生体信号取得部211は、センサ11を制御して、被験者の生体信号を取得する。また、取得した生体信号をフィルタ部212へ送信する。本実施の形態では、複数のセンサ11を用いて、同時に複数の生体信号が取得される。
【0022】
フィルタ部212は、特定の周波数帯域の信号を通過させる周波数フィルタを備え、生体信号取得部211で取得した生体信号を特定の周波数帯域の信号へ変換する。
【0023】
演算部213は、フィルタ部212から受信した複数の生体信号に基づいて、被験者の身体の状態を示す1次元の特徴パラメータを算出する。特徴パラメータは、複数の生体信号のゆらぎを示すものであり、ゆらぎの大きさ、すなわち時間的な変化の大きさによって被験者の身体状態の変化を表す。特徴パラメータの詳細な算出方法については、後述する。
【0024】
記憶部22は、ハードディスク、フラッシュメモリ等の不揮発性メモリであり、複数の生体信号から特徴パラメータを算出する演算アルゴリズム、算出された特徴パラメータ等を記憶する。
【0025】
表示部23は、コンピュータ装置である制御ユニット20に備えられた表示用デバイスであり、例えば液晶パネルである。表示部23は、センサ11で取得された生体信号、演算部213で算出された特徴パラメータ等を表示する。
【0026】
入力部24は、生体信号取得の開始、終了指示、特徴パラメータ算出のための条件変更等を入力するための入力デバイスである。入力部24は、制御ユニット20に備えられたキーボード、マウス等である。
【0027】
続いて、生体信号解析装置1を用いた生体信号解析方法について、
図2のフローチャートを参照しつつ、具体的に説明する。本実施の形態では、複数チャンネルの脳波を生体信号として取得し、被験者が癲癇の発作を生じているか否かを解析する場合を例として説明する。
【0028】
医師等のユーザは、
図3に示すように、被験者の頭部に複数のセンサ11を取り付ける。そして、入力部24を操作して、脳波の取得を開始する。生体信号取得部211は、脳波の取得開始指示により、センサ11を制御して、被験者の脳波を取得する(ステップS11)。生体信号取得部211は、取得した脳波データをフィルタ部212へ送信する(ステップS12)。
【0029】
フィルタ部212は、受信した各チャンネルの脳波データを、周波数フィルタによって、所定の周波数帯域の信号に変換する(ステップS13)。本実施の形態に係るフィルタ部212は、取得された脳波データを、α帯域(8〜12Hz)、β帯域(13〜24Hz)、γ帯域(25〜100Hz)、δ帯域(1〜3Hz)、θ帯域(4〜7Hz)の5つの帯域の信号に分解する。そして、演算部213は、フィルタ部212から受信した各周波数帯域の信号について、脳波のゆらぎを表す1次元パラメータである特徴パラメータを算出する。これにより、癲癇発作の種類等について、より詳細な分析が可能となる。
【0030】
以下、特徴パラメータの算出方法について説明する。演算部213は、分解された各帯域の脳波について、同一の時刻から所定の時間幅を持つ窓を設定する(ステップS14)。
【0031】
演算部213は、分解された帯域ごとに、信号のゆらぎを表す特徴パラメータを算出する(ステップS15)。より具体的には、演算部213は、まず、ステップS14で設定した窓内の、各チャンネルの脳波を混合して解析するための尺度混合モデルを生成する。本実施の形態に係る尺度混合モデルでは、センサ11であるL極(Lは自然数)の電極で計測された脳波xをx∈R
L、平均ゼロの確率変数と考え、分散共分散行列Σ∈R
L×Lの多変量正規分布を用いて表現する。ここで、Σは確率変数であり、Σの分布は自由度ν∈R
1と尺度行列Ψ∈R
L×Lによって特徴付けられる。
【0032】
上記の分散共分散行列Σが与えられた下での脳波xの条件付き分布P(x|Σ)は、以下の式(1)に示す多変量正規分布で表される。
【数1】
【0033】
分散共分散行列Σの分布として、以下の式(2)に示す、多変量正規分布の共役事前分布である逆ウィシャート分布P(Σ;ν,Ψ)を仮定する。
【数2】
【0034】
ここで、ν、Ψは、逆ウィシャート分布を決定するパラメータであり、νは自由度、Ψは尺度行列と呼ばれる。
【0035】
xの周辺分布P(x)は、以下の式(3)、(4)で表される。
【数3】
【0036】
ただし、Δはマハラノビス距離の平方で、以下の式(5)で表される。
【数4】
【0037】
上記の式(1)、(3)より、P(x)は、分散共分散行列の異なる多変量正規分布を無限個足し合わせた分布である。本実施の形態では、この多次元分散分布モデルP(x)を尺度混合モデルとして用いる。
【0038】
続いて、計測された脳波xから、分散共分散行列Σのパラメータとして定義される特徴パラメータを推定する方法について説明する。
【0039】
本実施の形態では、観測されたNサンプルの脳波データX={x
n∈R
L;n=1,・・・、N}に基づいて分散共分散行列Σのパラメータν、Ψを推定する。観測されたサンプル系列が尤も生起され易いパラメータは、次式で表す周辺分布の尤度を最大化することにより推定できる。
【数5】
【0040】
本実施の形態では、EM(Expectation-Maximization)アルゴリズムを用いて上記の推定を行う。
【0041】
簡単のため、式(4)において、パラメータをν’=ν−L+1,Ψ’=Ψ/(ν−L+1)と置き換えると、周辺尤度P(x
n)は以下の式(7)のように書き換えることができる。
【数6】
【0042】
ただし、Δ’は以下の式(8)で表される。
【数7】
ここで、式(7)はStudentのt分布と等価である。
【0043】
式(7)、(8)を用いて、以下の手順(a)〜(d)で周辺尤度P(x
n)が最大となるパラメータν’、Ψ’を推定する。
(a)各パラメータに任意の初期値を与える。
(b)現在のパラメータ値を使って、以下の式(9)に示す重みパラメータw
nを計算する(Eステップ)。
【数8】
【0044】
(c)手順(b)で求めたw
nを用いて、以下の式(10)により、Ψ’を更新する(Mステップ)。
【数9】
【0045】
ν’については、対数周辺尤度を最大化するν’を線形探索によって求める。
【数10】
【0046】
(d)更新後のパラメータを用いて対数周辺尤度を算出し、収束基準を満たすまで、上記の手順(b)及び(c)を反復する。最後に、推定されたt分布のパラメータであるν’、Ψ’を、元の分散共分散行列Σのパラメータν、Ψに変換する。
【0047】
上記の手順(a)〜(d)により、計測した脳波xから分散共分散行列Σのパラメータを推定する。
【0048】
推定されたパラメータであるνは、t分布の自由度に対応し、分布のガウス性を特徴付ける。本実施の形態の尺度混合モデルにおいては、分散に相当するパラメータの確率的ゆらぎによって、分布のガウス性が変化すると考えられる。すなわち、観測された脳波xから自由度νを推定することにより、脳波xの分散共分散行列Σのゆらぎを評価することができる。以上のようにして、脳波のゆらぎを表す1次元パラメータである特徴パラメータとして自由度νを算出する。
【0049】
特徴パラメータとしての自由度νを算出した後、解析を継続する場合(ステップS16のYES)、窓の開始時刻を所定のずらし幅で遅らせて(ステップS17)、ステップS15に戻り、特徴パラメータである自由度νの算出を繰り返す。
【0050】
解析を継続しない場合(ステップS16のNO)、制御部21は、算出された自由度νを時系列データとして記憶部22へ記憶させる。また、制御部21は、算出された自由度νの時系列データを、被験者の身体状態の変化を表す特徴パラメータとして、表示部23に表示させる(ステップS18)。そして、生体信号解析装置1は、解析処理を終了する。
【0051】
図4(A)、(B)は、癲癇患者の脳波データの例であり、
図4(A)は左半球の脳波データ、
図4(B)は右半球の脳波データを表している。
図5(A)、(B)は、
図4(A)、(B)の各データを上述の解析方法によって解析した解析結果を示す図である。
図5(A)、(B)では、自由度νの逆数である1/νを特徴パラメータとして、周波数帯域ごとに時系列のコンター図として示している。計測チャンネル数は、
図3に示す16チャンネルであり、計測時間は50±2秒、サンプリング周波数は500Hz、窓幅は10秒、ずらし幅は1秒である。
【0052】
図4(A)、(B)に示す脳波は、被験者が焦点発作を発生している場合のデータであり、
図5(A)、(B)に示すように、発作の発生を明確に判別できることがわかる。
【0053】
図6(A)、(B)及び
図7(A)、(B)は別の癲癇患者の脳波の例である。
図6(A)、(B)は、欠神発作を発生している場合のデータであり、
図7(A)、(B)に示すように、発作の種類が異なる場合でも、本発明に係る生体信号解析により、発作の発生を明確に判別できることがわかる。
【0054】
図8(A)、(B)及び
図9(A)、(B)は、さらに別の癲癇患者の脳波の例である。
図8(A)、(B)は、ミオクロニー発作を発生している場合のデータであり、
図9(A)、(B)に示すように、発作が上記2つの発作の場合と異なる周波数帯域、異なる大きさの信号として表れる場合であっても、複数の周波数フィルタを適用することで、明確に判別できることがわかる。
【0055】
以上、説明したように、本発明に係る生体信号解析装置及び生体信号解析方法では、複数の生体信号を混合して、生体信号のゆらぎを表す1次元パラメータを算出するので、容易に、かつ精度よく被験者の身体状態の変化を分析することが可能である。
【0056】
すなわち、本発明の解析方法によって1次元パラメータとして解析できるので、被験者の頭部に取り付けられた複数のセンサ11で取得される各チャンネルの脳波を医師が観察し、発作が生じているか否かを経験的に判定する方法に比べて、癲癇発作の発生等、被験者の身体状態の変化をより容易に分析することができる。
【0057】
また、複数の生体信号のゆらぎを解析できるので、個々の生体信号の振幅変化、周波数変化を解析する従来の方法に比べて、より精度よく癲癇発作等被験者の身体状態の変化を分析することができる。
【0058】
例として、γ帯域の信号に着目して、本発明の方法によって発作を識別した場合の識別率と、二乗平均平方根(RMS:Root Mean Square)により振幅変化を特徴量として発作を識別した場合の識別率とを比較する。
図10に示すように、本発明の方法による識別率は、二重平均平方根を用いた場合の識別率よりも高いことがわかる。
【0059】
また、本実施の形態では、センサ11で取得された生体信号をフィルタ部212の周波数フィルタによって、異なる周波数帯域を持つ複数の信号に分割して解析している。これにより、身体状態の変化の種類によって、その変化の表れやすい周波数帯域が異なる場合でも、適切に状態変化の発生を分析することができる。
【0060】
本実施の形態では、式(7)〜(11)に示すEMアルゴリズムによって特徴パラメータを推定することとしたが、これに限られない。特徴パラメータは、式(6)に示すモデルの周辺尤度を最大化可能な任意のアルゴリズムで推定できる。また、周辺尤度最大化以外のハイパーパラメータ推定法を適用することも可能である。
【0061】
本実施の形態では、被験者の脳波を生体信号として解析したが、これに限られない。例えば、脈拍、呼吸、発汗、心電図等の生体信号であってもよい。また、脳波と脈拍、呼吸と心電図等、複数の種類の生体信号を混合して解析することとしてもよい。これにより、身体の様々な部位の状態変化や複合的な状態変化を分析することができる。
【0062】
また、本実施の形態では、生体信号である被験者の脳波に基づいて、癲癇の発作を検出することとしたが、これに限られない。例えば、脳波を含む生体信号に基づいて、意識障害を伴う脳炎、糖尿病性昏睡、脳血管障害、脳腫瘍、アルツハイマー型痴呆等の異常脳波を伴う疾患を、検出することとしてもよい。また、被験者の覚醒状態、睡眠状態の別や睡眠深度などの脳状態を分析することとしてもよい。
【0063】
また、本実施の形態では、センサ11で取得した生体信号である脳波を、制御ユニット20で逐次解析することとしたが、これに限られない。例えば、過去に取得された生体信号を制御ユニット20の生体信号取得部211に読み込んで、解析を実行し、現在の生体信号の解析結果と比較することとしてもよい。