(57)【要約】      （修正有）【課題】短い観測時間で心拍変動の周波数解析を精度良く行う方法を提供する。【解決手段】心拍変動の周波数解析方法であって、所定数のＲ波発生時刻からＲＲ間隔を算出するステップＳ２と、各前記ＲＲ間隔の逆数を瞬時心拍数とする時間関数を定義するステップＳ３と、前記時間関数に対するフーリエ変換を適用して連続スペクトルを生成するステップＳ４と、を含む周波数解析方法。最新の所定数のＲ波発生時刻について、上記ステップを出来る限り短い周期で繰り返し実行することにより、心拍変動の経時的変化を高い時間分解能で解析する。【選択図】図１
2019092697.tif 000002
【特許請求の範囲】
【請求項１】
心拍変動の周波数解析方法であって、
所定数のＲ波発生時刻からＲＲ間隔を算出するステップと、
各前記ＲＲ間隔の逆数を瞬時心拍数とする時間関数を定義するステップと、
前記時間関数に対するフーリエ変換を適用して連続スペクトルを生成するステップと、
を含む周波数解析方法。
【請求項２】
心拍変動の経時的変化を解析する方法であって、下記ステップ（ａ）乃至（ｃ）を繰り返し実行することを特徴とする方法。
（ａ）最新の所定数のＲ波発生時刻からＲＲ間隔を算出するステップ。
（ｂ）各前記ＲＲ間隔の逆数を瞬時心拍数とする時間関数を定義するステップ。
（ｃ）前記時間関数に対するフーリエ変換を適用して連続スペクトルを生成するステップ。
【請求項３】
コンピュータに、請求項１または２に記載の方法の各ステップを実行させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、心拍変動の周波数解析方法に関し、より詳細には、短い観測時間で心拍変動の周波数解析を精度良く行う方法に関する。
【背景技術】
【０００２】
従来、心拍変動(ＨＲＶ：Heart Rate Variability)の周波数スペクトル特性が自律神経系の指標として用いられており、その周波数解析においては、不等時間間隔のＲＲ間隔を等時間間隔に補正した離散的な時系列データに対して離散フーリエ変換（ＤＦＴ／ＦＦＴ）を実行する方法が一般的である。
【０００３】
一方、心拍変動の周波数解析を利用して対象者の状態（心理状態、精神状態、身体状態など）を判定することが種々検討されている（例えば、特許文献１、２）。この場合、対象者の状態の変化を高い時間分解能で判定しようとすれば、心拍データの観測時間をできるだけ短くする必要があるが、一定の時間データ長を基底とする離散フーリエ変換では、観測時間の逆数が周波数分解能となるため、観測時間を短くするほど周波数分解能が低下して解析精度が落ちるというトレードオフの問題が生じる。
【先行技術文献】
【特許文献】
【０００４】
【特許文献１】特開２０１６−２２０８２号公報
【特許文献１】特開２０１７−５１４９６号公報
【発明の概要】
【発明が解決しようとする課題】
【０００５】
本発明は、上記従来技術における課題に鑑みてなされたものであり、短い観測時間で心拍変動の周波数解析を精度良く行う方法を提供することを目的とする。
【課題を解決するための手段】
【０００６】
本発明者は、短い観測時間で心拍変動の周波数解析を精度良く行う方法につき鋭意検討した結果、以下の構成に想到し、本発明に至ったのである。
【０００７】
すなわち、本発明によれば、心拍変動の周波数解析方法であって、所定数のＲ波発生時刻からＲＲ間隔を算出するステップと、各前記ＲＲ間隔の逆数を瞬時心拍数とする時間関数を定義するステップと、前記時間関数に対するフーリエ変換を適用して連続スペクトルを生成するステップと、を含む周波数解析方法が提供される。
【発明の効果】
【０００８】
上述したように、本発明によれば、短い観測時間で心拍変動の周波数解析を精度良く行う方法が提供される。
【図面の簡単な説明】
【０００９】
【図１】本実施形態の心拍変動の解析方法のフローチャート。
【図２】瞬時心拍数の時間関数を示す図。
【図３】心拍信号を模したサンプルデータの周波数解析結果を示す図。
【発明を実施するための形態】
【００１０】
以下、本発明を図面に示した実施の形態をもって説明するが、本発明は、図面に示した実施の形態に限定されるものではない。
【００１１】
本発明は、心拍数を基底とする心拍変動の周波数解析方法を開示する。以下、本発明の実施形態である心拍変動の周波数解析方法の手順を図１に示すフローチャートに基づいて説明する。なお、以下では、便宜的に、心拍数「６」を基底とする場合を例にとって説明を行う。
【００１２】
まず、ステップＳ１では、対象者の心拍信号から、時間的に連続する所定数のＲ波のピーク値の発生時刻（以下、Ｒ波発生時刻という）を取得する。
【００１３】
本例では、心拍数「６」を基底とするので、図２（ａ）に示すように、６個のＲ波発生時刻ｔ０〜ｔ５を取得する。
【００１４】
続くステップＳ２では、検出した６個のＲ波発生時刻から５個のＲＲ間隔を算出する。
【００１５】
本例では、図２（ｂ）に示すように、時刻ｔ１と時刻ｔ０の差分としてＲＲ間隔１（以下、ＲＲ１という）を算出し、時刻ｔ２と時刻ｔ１の差分としてＲＲ間隔２（以下、ＲＲ２という）を算出し、時刻ｔ３と時刻ｔ２の差分としてＲＲ間隔３（以下、ＲＲ３という）を算出し、時刻ｔ４と時刻ｔ３の差分としてＲＲ間隔４（以下、ＲＲ４という）を算出し、時刻ｔ５と時刻ｔ４の差分としてＲＲ間隔５（以下、ＲＲ５という）を算出する。
【００１６】
続くステップＳ３では、算出したＲＲ間隔の逆数を瞬時心拍数とする時間関数を定義する。具体的には、隣り合う２つのＲ波発生時刻の間隔（ＲＲ間隔）の逆数を当該２つのＲ波発生時刻の間の各時刻の瞬時心拍数と推定して、瞬時心拍数の時間関数を定義する。
【００１７】
この例では、ＲＲ１の逆数（ＲＲ１）−１を時刻ｔ０と時刻ｔ１の間の各時刻の瞬時心拍数Ａ１と推定し、ＲＲ２の逆数（ＲＲ２）−１を時刻ｔ１と時刻ｔ２の間の各時刻の瞬時心拍数Ａ２と推定し、ＲＲ３の逆数（ＲＲ３）−１を時刻ｔ２と時刻ｔ３の間の各時刻の瞬時心拍数Ａ３と推定し、ＲＲ４の逆数（ＲＲ４）−１を時刻ｔ３と時刻ｔ４の間の各時刻の瞬時心拍数Ａ４と推定し、ＲＲ５の逆数（ＲＲ５）−１を時刻ｔ４と時刻ｔ５の間の各時刻の瞬時心拍数Ａ５と推定することにより、瞬時心拍数の時間関数として、図２（ｃ）に示すような階段関数が定義される。
【００１８】
続くステップＳ４では、定義した瞬時心拍数の時間関数に対するフーリエ変換を適用して連続スペクトルを生成する。なお、ここでいう「フーリエ変換」とは、連続的な関数に対するフーリエ変換(FT：Fourier Transform)を意味し、高速フーリエ変換（FFT：Fast Fourier Transform）を含む離散フーリエ変換（DFT：Discrete Fourier Transform）を意味しないことに留意されたい。
【００１９】
ここで、本実施形態では、図２（ｃ）に示す瞬時心拍数の時間関数のフーリエ変換を以下の手順で行う。
【００２０】
瞬時心拍数の時間関数のt0〜t1区間をフーリエ変換した関数F1(ω)は、下記式（１）で表される。
【００２１】
【数１】
2019092697.tif 000003
【００２２】
同様に、瞬時心拍数の時間関数のt1〜t2区間、t2〜t3区間、t3〜t4区間、t4〜t5間の各々をフーリエ変換した関数F2(ω)、F3(ω)、F4(ω)、F5(ω)は、下記式（２）で表される。
【００２３】
【数２】
2019092697.tif 000004
【００２４】
したがって、フーリエ変換の線形性より、t0〜t5区間をフーリエ変換した関数F(ω)は、下記式（３）で表される。
【００２５】
【数３】
2019092697.tif 000005
【００２６】
そして、F(ω)の振幅スペクトルA(ω)は、下記式（４）で表される。
【００２７】
【数４】
2019092697.tif 000006
【００２８】
つまり、本実施形態では、予め、瞬時心拍数の時間関数（階段関数）のフーリエ変換を実施して一般式を導出しておき、当該一般式にA1〜A5、t0〜t5などのデータを入力して演算を行うことにより連続関数を生成する。ここで、上記式（４）において、a1〜a5、R1〜R5、I1〜I5のそれぞれの値は、定数であるA1〜A5、t0〜t5、およびωで構成され、A(ω)は、ω≠0の区間で連続する値を取り、連続スペクトルとなる。
【００２９】
以上、説明したように、本実施形態によれば、心拍変動の周波数解析結果として連続スペクトルが得られるので、短い観測時間でも必要な周波数分解能が確保される。さらに、最新の所定数のＲ波発生時刻について、図１に示すステップＳ１〜Ｓ４を出来る限り短い周期で繰り返し実行すれば、心拍変動の経時的変化を高い時間分解能で解析することが可能になる。
【００３０】
なお、図１に示した各ステップをコンピュータに実行させるためのプログラムは、任意のプログラム言語で記述することができ、任意の記録媒体に格納して頒布することができ、任意のネットワークを介して伝送することができる。
【００３１】
以上、本発明について実施形態をもって説明してきたが、本発明は上述した実施形態に限定されるものではなく、当業者が推考しうるその他の実施態様の範囲内において、本発明の作用・効果を奏する限り、本発明の範囲に含まれるものである。
【実施例】
【００３２】
以下、本発明の心拍変動の解析方法について、実施例を用いてより具体的に説明を行なうが、本発明は、後述する実施例に限定されるものではない。
【００３３】
下記表１に示す内容の心拍信号を模したサンプルデータを作成し、実施例として、本発明の解析方法による周波数解析を行った。
【００３４】
【表１】
2019092697.tif 000007
【００３５】
併せて、比較例として、同じサンプルデータについて、MATLABのFFT関数を使用して周波数解析を行った（サンプリング周波数10Hz）。
【００３６】
図３（ａ）は比較例の解析結果を示す。図３（ａ）に示すように、比較例では、３種類の周波数成分（0.090Hz、0.095Hz、0.10Hz）の振幅ピーク値について最大で２割程度の差が生じる結果となった。これは、比較例の周波数分解能が約0.026Hz（観測時間≒0.39secの逆数）であることから、LF成分に対応するプロット周波数が0.078Hz、0.104Hzとなるところ、サンプルデータのLF成分0.090Hzがその中間にあたるために乖離が大きくなったことが原因と考えられる。
【００３７】
一方、図３（ｂ）は実施例の解析結果を示す。図３（ｂ）に示すように、実施例では、３種類の周波数成分（0.090Hz、0.095Hz、0.10Hz）のそれぞれの振幅ピーク値が概ね等しかった。この結果から、本発明の解析方法によれば、約39秒という短い観測時間で精度良い結果が得られることが示された。
【図１】
2019092697.tif 000008
【図２】
2019092697.tif 000009
【図３】
2019092697.tif 000010