(57)【要約】【課題】脈波信号から高精度に脈波を解析する。【解決手段】生体の脈波を表す脈波信号のうち、当該脈波の基本周波数を含む解析帯域内の成分を観測信号として抽出する成分抽出部と、前記基本周波数で変動する第１周期信号を、前記観測信号に乗算することで、同相信号を生成する第１処理部と、前記第１周期信号に直交し、前記基本周波数で変動する第２周期信号を、前記観測信号に乗算することで、直交信号を生成する第２処理部と、前記同相信号および前記直交信号で規定され、時間軸上の第１時点における信号値が示す第１偏角と、前記同相信号および前記直交信号で規定され、前記第１時点とは異なる第２時点における信号値が示す第２偏角との差分を、前記基本周波数に対する周波数の変動を表す周波数変動指標として算定する差分算定部と、前記生体の脈周期に関する周期指標を前記周波数変動指標から算定する指標算定部とを具備する脈波解析装置。【選択図】図３
2019118669.tif 000002
【特許請求の範囲】
【請求項１】
生体の脈波を表す脈波信号のうち、当該脈波の基本周波数を含む解析帯域内の成分を観測信号として抽出する成分抽出部と、
前記基本周波数で変動する第１周期信号を、前記観測信号に乗算することで、同相信号を生成する第１処理部と、
前記第１周期信号に直交し、前記基本周波数で変動する第２周期信号を、前記観測信号に乗算することで、直交信号を生成する第２処理部と、
前記同相信号および前記直交信号で規定され、時間軸上の第１時点における信号値が示す第１偏角と、前記同相信号および前記直交信号で規定され、前記第１時点とは異なる第２時点における信号値が示す第２偏角との差分を、前記基本周波数に対する周波数の変動を表す周波数変動指標として算定する差分算定部と、
前記生体の脈周期に関する周期指標を前記周波数変動指標から算定する指標算定部と
を具備する脈波解析装置。
【請求項２】
前記脈波信号の解析により前記基本周波数を特定する周波数解析部を具備する
請求項１の脈波解析装置。
【請求項３】
前記差分算定部は、前記同相信号の信号値と前記直交信号の信号値との関係を表す複素関数について、前記第１時点と前記第２時点との間における自己相関の偏角を、前記周波数変動指標として算定する
請求項１または請求項２の脈波解析装置。
【請求項４】
前記周期指標に応じて前記解析帯域の帯域幅を制御する帯域制御部を具備する
請求項１から請求項３の何れかの脈波解析装置。
【請求項５】
生体の脈波を表す脈波信号のうち、当該脈波の基本周波数を含む解析帯域内の成分を観測信号として抽出し、
前記基本周波数で変動する第１周期信号を、前記観測信号に乗算することで、同相信号を生成し、
前記第１周期信号に直交し、前記基本周波数で変動する第２周期信号を、前記観測信号に乗算することで、直交信号を生成し、
前記同相信号および前記直交信号で規定され、時間軸上の第１時点における信号値が示す第１偏角と、前記同相信号および前記直交信号で規定され、前記第１時点とは異なる第２時点における信号値が示す第２偏角との差分を、前記基本周波数に対する周波数の変動を表す周波数変動指標として算定し、
前記生体の脈周期に関する周期指標を前記周波数変動指標から算定する
脈波解析方法。
【請求項６】
生体の脈波を表す脈波信号のうち、当該脈波の基本周波数を含む解析帯域内の成分を観測信号として抽出する成分抽出部、
前記基本周波数で変動する第１周期信号を、前記観測信号に乗算することで、同相信号を生成する第１処理部、
前記第１周期信号に直交し、前記基本周波数で変動する第２周期信号を、前記観測信号に乗算することで、直交信号を生成する第２処理部、
前記同相信号および前記直交信号で規定され、時間軸上の第１時点における信号値が示す第１偏角と、前記同相信号および前記直交信号で規定され、前記第１時点とは異なる第２時点における信号値が示す第２偏角との差分を、前記基本周波数に対する周波数の変動を表す周波数変動指標として算定する差分算定部、および、
前記生体の脈周期に関する周期指標を前記周波数変動指標から算定する指標算定部
としてコンピューターを機能させるプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、生体の脈波を解析するための技術に関する。
【背景技術】
【０００２】
生体の脈波に関する指標を算定する技術が従来から提案されている。例えば特許文献１には、生体の脈波の大きさに応じた信号（容積脈波）から、心電図のＲ−Ｒ間隔の変動（すなわち心拍変動）に相当するａ−ａ間隔の変動を特定する構成が開示されている。
【先行技術文献】
【特許文献】
【０００３】
【特許文献１】特開２０１０−０５１８２２号公報
【発明の概要】
【発明が解決しようとする課題】
【０００４】
しかし、特許文献１の技術のように、容積脈波から生体の脈波を解析する場合、生体の体動に起因したノイズや装置の回路ノイズの影響を解析結果が受けやすいという問題がある。以上の事情を考慮して、本発明の好適な態様は、脈波信号から高精度に脈波を解析することを目的とする。
【課題を解決するための手段】
【０００５】
以上の課題を解決するために、本発明の好適な態様に係る脈波解析装置は、生体の脈波を表す脈波信号のうち、当該脈波の基本周波数を含む解析帯域内の成分を観測信号として抽出する成分抽出部と、前記基本周波数で変動する第１周期信号を、前記観測信号に乗算することで、同相信号を生成する第１処理部と、前記第１周期信号に直交し、前記基本周波数で変動する第２周期信号を、前記観測信号に乗算することで、直交信号を生成する第２処理部と、前記同相信号および前記直交信号で規定され、時間軸上の第１時点における信号値が示す第１偏角と、前記同相信号および前記直交信号で規定され、前記第１時点とは異なる第２時点における信号値が示す第２偏角との差分を、前記基本周波数に対する周波数の変動を表す周波数変動指標として算定する差分算定部と、前記生体の脈周期に関する周期指標を前記周波数変動指標から算定する指標算定部とを具備する。以上の態様では、脈波信号のうち脈波の基本周波数を含む解析帯域内の成分が観測信号として抽出され、当該観測信号を利用して脈周期に関する周期指標が算定されるから、例えば脈波信号（つまり解析帯域外の成分も含む信号）から周期指標を算定する構成と比較して、生体の体動に起因したノイズや装置の回路ノイズの影響を低減した周期指標を算定することができる。したがって、脈波信号から周波数変動指標を高精度に算定することができる。
【０００６】
本発明の好適な態様において、前記脈波信号の解析により前記基本周波数を特定する周波数解析部を具備する。以上の態様では、周波数変動指標の算定に利用される脈波信号の解析により基本周波数が特定されるから、周波数変動指標の算定に利用する脈波信号を基本周波数の特定に流用できるという利点がある。
【０００７】
本発明の好適な態様において、前記差分算定部は、前記同相信号の信号値と前記直交信号の信号値との関係を表す複素関数について、前記第１時点と前記第２時点との間における自己相関の偏角を、前記周波数変動指標として算定する。以上の態様では、同相信号の信号値と直交信号の信号値との関係を表す複素関数について、第１時点と第２時点との間における自己相関の偏角が、周波数変動指標として算定されるから、第１偏角と第２偏角とをそれぞれ算定することが不要になる。
【０００８】
本発明の好適な態様において、前記周期指標に応じて前記解析帯域の帯域幅を制御する帯域制御部を具備する。以上の態様では、周期指標に応じて解析帯域の帯域幅が制御されるから、周期指標の変動する範囲に応じた適切な解析帯域が設定される。
【０００９】
本発明の好適な態様に係る脈波解析方法は、生体の脈波を表す脈波信号のうち、当該脈波の基本周波数を含む解析帯域内の成分を観測信号として抽出し、前記基本周波数で変動する第１周期信号を、前記観測信号に乗算することで、同相信号を生成し、前記第１周期信号に直交し、前記基本周波数で変動する第２周期信号を、前記観測信号に乗算することで、直交信号を生成し、前記同相信号および前記直交信号で規定され、時間軸上の第１時点における信号値が示す第１偏角と、前記同相信号および前記直交信号で規定され、前記第１時点とは異なる第２時点における信号値が示す第２偏角との差分を、前記基本周波数に対する周波数の変動を表す周波数変動指標として算定し、前記生体の脈周期に関する周期指標を前記周波数変動指標から算定する。
【００１０】
本発明の好適な態様に係るプログラムは、生体の脈波を表す脈波信号のうち、当該脈波の基本周波数を含む解析帯域内の成分を観測信号として抽出する成分抽出部、前記基本周波数で変動する第１周期信号を、前記観測信号に乗算することで、同相信号を生成する第１処理部、前記第１周期信号に直交し、前記基本周波数で変動する第２周期信号を、前記観測信号に乗算することで、直交信号を生成する第２処理部、前記同相信号および前記直交信号で規定され、時間軸上の第１時点における信号値が示す第１偏角と、前記同相信号および前記直交信号で規定され、前記第１時点とは異なる第２時点における信号値が示す第２偏角との差分を、前記基本周波数に対する周波数の変動を表す周波数変動指標として算定する差分算定部、および、前記生体の脈周期に関する周期指標を前記周波数変動指標から算定する指標算定部としてコンピューターを機能させる。
【図面の簡単な説明】
【００１１】
【図１】本発明の第１実施形態に係る脈波解析装置の側面図である。
【図２】脈波解析装置の機能に着目した構成図である。
【図３】制御装置の機能に着目した構成図である。
【図４】強度スペクトルを示すグラフである。
【図５】観測周波数の時間変化を表すグラフである。
【図６】制御装置が実行する処理のフローチャートである。
【図７】第２実施形態に係る脈波解析装置の機能に着目した構成図である。
【図８】実製品に付与される検出信号を示すグラフである。
【図９】第３実施形態に係る脈波解析装置の使用例を示す模式図である。
【図１０】第３実施形態に係る脈波解析装置の他の使用例を示す模式図である。
【図１１】変形例における脈波解析装置の構成図である。
【図１２】変形例における脈波解析装置の構成図である。
【図１３】変形例における脈波解析装置の構成図である。
【発明を実施するための形態】
【００１２】
＜第１実施形態＞
図１は、本発明の第１実施形態に係る脈波解析装置１００の側面図である。第１実施形態の脈波解析装置１００は、被験者の心拍変動（Heart Rate Variability）に関する指標（以下「心拍変動指標」という）を算定する生体測定機器である。心拍変動は、生体の脈波の周期（以下「脈周期」という）の時間的な変動である。第１実施形態の脈波解析装置１００は、被験者の身体のうち特定の部位（以下「測定部位」という）Ｍの脈周期から心拍変動指標を測定する。以下の説明では、被験者の手首または上腕を測定部位Ｍとして例示する。
【００１３】
脈波解析装置１００は、測定部位Ｍに装着される。第１実施形態の脈波解析装置１００は、図１に例示される通り、筐体部１２とベルト１４とを具備する腕時計型の携帯機器である。脈波解析装置１００は、測定部位Ｍにベルト１４を巻回することで被験者の身体に装着される。
【００１４】
図２は、脈波解析装置１００の機能に着目した構成図である。図２に例示される通り、第１実施形態の脈波解析装置１００は、制御装置２０と記憶装置２２と表示装置２４と検出装置３０とを具備する。制御装置２０および記憶装置２２は、筐体部１２の内部に設置される。図１に例示される通り、表示装置２４（例えば液晶表示パネル）は、例えば筐体部１２における測定部位Ｍとは反対側の表面に設置され、測定結果を含む各種の画像を制御装置２０による制御のもとで表示する。
【００１５】
図２の検出装置３０は、測定部位Ｍの状態に応じた検出信号（脈波信号の例示）Ｓa(ｔ)を生成する光学センサーモジュールである。図２に例示される通り、第１実施形態の検出装置３０は、発光部３１と受光部３２と駆動回路３３と出力回路３４とを具備する。発光部３１および受光部３２は、例えば筐体部１２において測定部位Ｍに対向する位置（典型的には測定部位Ｍに接触する表面）に設置される。なお、駆動回路３３および出力回路３４の一方または双方を検出装置３０とは別体の外部回路として設置することも可能である。
【００１６】
発光部３１は、測定部位Ｍに光を照射する光源である。第１実施形態の発光部３１は、例えばインコヒーレントな光を出射するＬＥＤ（Light Emitting Diode）等の発光素子が発光部３１として好適に利用される。駆動回路３３は、制御装置２０の制御のもとで発光部３１を発光させる。なお、相異なる波長の光を出射する複数の発光素子で発光部３１を構成してもよい。発光部３１から測定部位Ｍに入射した光は、測定部位Ｍの内部を通過しながら拡散反射を繰返したうえで筐体部１２側に出射する。具体的には、測定部位Ｍの内部に存在する動脈（例えば、上腕動脈、橈骨動脈または尺骨動脈）等の血管と血管内の血液とを通過した光が測定部位Ｍから筐体部１２側に出射する。
【００１７】
受光部３２は、測定部位Ｍから到来する光を受光する。第１実施形態の受光部３２は、測定部位Ｍから到達する光の強度の時間変化を表す検出信号（つまり生体を通過した光の強度の時間変化を表す信号）を生成する。例えば、受光部３２は、受光強度に応じた電荷を発生するフォトダイオード（ＰＤ：Photo Diode）等の受光素子を含む。測定部位Ｍの内部の血管は、拍動と同等の周期で反復的に拡張および収縮するから、測定部位Ｍからの受光強度に応じて受光部３２が生成する検出信号は、測定部位Ｍの血管の血流量の変動に対応した周期的な変動成分を含む脈波信号（脈波を表す信号）である。すなわち、検出信号には、脈周期で変動する成分が含まれている。
【００１８】
第１実施形態の検出装置３０は、発光部３１と受光部３２とが測定部位Ｍに対して片側に位置する反射型の光学センサーである。ただし、発光部３１と受光部３２とが測定部位Ｍを挟んで反対側に位置する透過型の光学センサーを検出装置３０として利用してもよい。
【００１９】
出力回路３４は、受光部３２が受光した光の強度に応じた検出信号Ｓa(ｔ)を生成する。具体的には、出力回路３４は、受光部３２が生成した検出信号を増幅する増幅回路（図示略）と、増幅後の検出信号を所定のサンプリング周期毎にアナログからデジタルに変換することで検出信号Ｓa(ｔ)を生成するＡ/Ｄ変換器（図示略）とを具備する。出力回路３４が生成した検出信号Ｓa(ｔ)は、制御装置２０に供給される。符号ｔは、Ａ/Ｄ変換における時間軸上の任意の１個のサンプリング点を意味する。
【００２０】
図２の制御装置２０は、ＣＰＵ（Central Processing Unit）またはＦＰＧＡ（Field-Programmable Gate Array）等の演算処理装置であり、脈波解析装置１００の全体を制御する。第１実施形態の制御装置２０は、記憶装置２２に記憶されたプログラム（アプリケーションプログラム）を実行することで、心拍変動指標を算定するための機能（周波数解析部６２，成分抽出部６４，演算処理部６６，変動算定部６８）を実現する。記憶装置２２は、例えば不揮発性の半導体メモリーで構成され、制御装置２０が実行するプログラムと制御装置２０が使用する各種のデータとを記憶する。なお、制御装置２０の機能を複数の集積回路に分散した構成、または、制御装置２０の一部または全部の機能を専用の電子回路で実現した構成も採用され得る。また、図２では制御装置２０と記憶装置２２とを別個の要素として図示したが、記憶装置２２を内包する制御装置２０を例えばＡＳＩＣ（Application Specific Integrated Circuit）等により実現することも可能である。
【００２１】
図３は、制御装置２０の機能に着目した構成図である。周波数解析部６２は、検出信号Ｓa(ｔ)の解析により、検出信号Ｓa(ｔ)が表す脈波の基本周波数ｆ0を特定する。基本周波数ｆ0は、被験者の脈拍数に対応する。具体的には、周波数解析部６２は、検出装置３０が生成した検出信号Ｓa(ｔ)から周波数に関する強度スペクトル（例えばパワースペクトル）Ｘを算定し、当該強度スペクトルＸから基本周波数ｆ0を特定する。強度スペクトルＸは、各周波数ｆにおける検出信号Ｓa(ｔ)の信号成分の信号強度（パワーまたは振幅）Ｇ(f)の分布である。強度スペクトルＸの算定には、高速フーリエ変換（ＦＦＴ：Fast Fourier Transform）等の公知の周波数解析が任意に採用され得る。図４は、周波数解析部６２が算定した強度スペクトルＸを示すグラフである。図４に例示される通り、強度スペクトルＸには、脈波を表す複数のピークが含まれる。強度スペクトルＸに含まれる複数のピークのうち、強度Ｇ(f)が最大となるピークの周波数が脈波の基本周波数ｆ0として特定される。図４では、１Ｈｚが基本周波数ｆ0として特定される場合が例示されている。周波数解析部６２による基本周波数ｆ0の特定は、時間軸上の単位期間（フレーム）毎に反復的に実行される。単位期間は、例えば拍動の１拍分に相当する期間より充分に長い時間長（例えば１６秒）である。
【００２２】
図３の成分抽出部６４は、検出信号Ｓa(ｔ)のうち、脈波の基本周波数ｆ0を含む帯域（以下「解析帯域」という）Ｂ内の成分を観測信号Ｓb(ｔ)として抽出する。図４に例示される通り、解析帯域Ｂは、例えば基本周波数ｆ0を中心周波数とした所定の帯域幅（例えば１Ｈｚ）の通過帯域である。解析帯域Ｂ内の成分を選択的に通過させて、他の成分を除去するような周波数応答を有するバンドパスフィルターが成分抽出部６４として利用される。つまり、解析帯域Ｂ外に存在するノイズ（例えば体動に起因したノイズや脈波解析装置１００の回路ノイズ）の影響を低減した観測信号Ｓb(ｔ)が成分抽出部６４により生成される。
【００２３】
成分抽出部６４で生成される観測信号Ｓb(ｔ)には、脈周期で変動する成分が含まれる。脈周期は、観測信号Ｓb(ｔ)において相前後する２個のピーク（心臓から測定部位Ｍに到達する進行波に相当するピーク）間の間隔である。被験者の脈周期は、時間軸上で変動し得る。図５は、脈周期の逆数である周波数（以下「観測周波数」）ｆの時間変化を表すグラフである。図５から把握される通り、観測周波数ｆは、脈周期の変動に応じて経時的に変化する。図５には、基本周波数ｆ0に対する観測周波数ｆの変動量Δｆが図示されている。基本周波数ｆ0と観測周波数ｆとの差分が変動量Δｆである。すなわち、周波数解析部６２が特定する基本周波数ｆ0は、観測周波数ｆの変動量Δｆにおける基準値である。時間軸上の任意の時点ｔにおける観測周波数ｆは、基本周波数ｆ0と基本周波数ｆ0に対する変動量Δｆとの和（ｆ0＋Δｆ）で表現されるから、観測信号Ｓb(ｔ)は、以下の数式(1)で表現される。数式(1)から把握される通り、観測信号Ｓb(ｔ)は、観測周波数ｆ(＝ｆ0＋Δｆ)を係数とした正弦関数で表現される。
【数１】
2019118669.tif 000003
【００２４】
図３の演算処理部６６は、成分抽出部６４が抽出した観測信号Ｓb(ｔ)から、生体の脈周期Ｐを算定する。観測信号Ｓb(ｔ)に対する直交検波により脈周期Ｐを算定する。具体的には、演算処理部６６は、信号処理部９１と第１処理部９３と第２処理部９５と差分算定部９７と指標算定部９８とを具備する。
【００２５】
信号処理部９１は、基本周波数ｆ0で変動する第１周期信号Ｚi(ｔ)と、基本周波数ｆ0で変動する第２周期信号Ｚq(ｔ)とを生成する。第１周期信号Ｚi(ｔ)と第２周期信号Ｚq(ｔ)とは相互に直交する。具体的には、第１周期信号Ｚi(ｔ)は、以下の数式(2)の余弦関数で表現され、第２周期信号Ｚq(ｔ)は、以下の数式(3)の正弦関数で表現される。具体的には、信号処理部９１は、第１周期信号Ｚi(ｔ)を生成する発振器と、第１周期信号Ｚi(ｔ)を９０°だけ移相することで第２周期信号Ｚq(ｔ)を生成する移相器とを含んで構成される。
【数２】
2019118669.tif 000004
【数３】
2019118669.tif 000005
【００２６】
第１処理部９３は、信号処理部９１が生成した第１周期信号Ｚi(ｔ)を観測信号Ｓb(ｔ)に乗算することで、同相信号Ｓi(ｔ)を生成する。第１実施形態の第１処理部９３は、乗算部９３１と低域通過部９３３とを具備する。乗算部９３１は、数式(1)で表現される観測信号Ｓb(ｔ)と、数式(2)で表現される第１周期信号Ｚi(ｔ)とを乗算する。乗算部９３１による乗算結果（数式(1)と数式(2)との積で表現される信号）には、基本周波数ｆ0と変動量Δｆとの和成分(ｆ0＋Δｆ＋ｆ0)および差成分(ｆ0＋Δｆ−ｆ0)が含まれる。低域通過部９３３は、乗算部９３１の乗算結果から、以下の数式(4)で表現される同相信号Ｓi(ｔ)を生成する。数式(4)から理解される通り、同相信号Ｓi(ｔ)は、変動量Δｆを係数とした余弦関数で表現される。乗算部９３１の乗算結果のうち、同相信号Ｓi(ｔ)を表す成分を通過させ、同相信号Ｓi(ｔ)を表す成分よりも高域側の成分を除去するローパスフィルターが低域通過部９３３として利用される。具体的には、乗算部９３１の乗算結果から和成分(ｆ0＋Δｆ＋ｆ0)が低域通過部９３３により除去される。
【数４】
2019118669.tif 000006
【００２７】
第２処理部９５は、信号処理部９１が生成した第２周期信号Ｚq(ｔ)を観測信号Ｓb(ｔ)に乗算することで、直交信号Ｓq(ｔ)を生成する。第１実施形態の第２処理部９５は、乗算部９５１と低域通過部９５３とを具備する。乗算部９５１は、数式(1)で表現される観測信号Ｓb(ｔ)と、数式(3)で表現される第２周期信号Ｚq(ｔ)とを乗算する。乗算部９５１による乗算結果（数式(1)と数式(3)との積で表現される信号）には、基本周波数ｆ0と変動量Δｆとの和成分(ｆ0＋Δｆ＋ｆ0)および差成分(ｆ0＋Δｆ−ｆ0)が含まれる。低域通過部９５３は、乗算部９５１の乗算結果から、以下の数式(5)で表現される直交信号Ｓq(ｔ)を生成する。数式(5)から理解される通り、直交信号Ｓq(ｔ)は、変動量Δｆを係数とした正弦関数で表現される。乗算部９５１の乗算結果のうち、直交信号Ｓq(ｔ)を表す成分を通過させ、直交信号Ｓq(ｔ)を表す成分よりも高域側の成分を除去するローパスフィルターが低域通過部９５３として利用される。具体的には、乗算部９５１の乗算結果から和成分(ｆ0＋Δｆ＋ｆ0)が低域通過部９５３により除去される。
【数５】
2019118669.tif 000007
【００２８】
差分算定部９７は、同相信号Ｓi(ｔ)と直交信号Ｓq(ｔ)とから、基本周波数ｆ0に対する観測周波数ｆの変動を表す指標（以下「周波数変動指標」という）を算定する。具体的には、差分算定部９７は、同相信号Ｓi(ｔ)および直交信号Ｓq(ｔ)で規定され、時間軸上の第１時点ｔ1における信号値が示す第１偏角θ1と、同相信号Ｓi(ｔ)および直交信号Ｓq(ｔ)で規定され、第１時点ｔ1とは異なる第２時点ｔ2における信号値が示す第２偏角θ2との差分を周波数変動指標として算定する。第１時点ｔ1は、Ａ／Ｄ変換器のサンプリング周期Ｔのｎ個分に相当する時間だけ時点ｔ0よりも後方の時点（ｔ0＋ｎＴ）である。第２時点ｔ2は、サンプリング周期Ｔの（ｎ＋１）個分に相当する時間だけ時点ｔ0よりも後方の時点（ｔ0＋(ｎ＋１)Ｔ）である。すなわち、第１時点ｔ1と第２時点ｔ2とは、時間軸上で相前後するサンプリング点である。時点ｔ0は時間軸上の特定の時点である。
【００２９】
ここで、数式(4)の同相信号Ｓi(ｔ)の信号値と数式(5)の直交信号Ｓq(ｔ)の信号値とは、以下の数式(6)の複素関数で表現することが可能である（オイラーの公式）。数式(6)は、時点ｔにおける偏角が２πΔｆｔで表現される複素関数である。
【数６】
2019118669.tif 000008
第１実施形態の差分算定部９７は、数式(6)で表現される複素関数について、第１時点ｔ1における第１偏角θ1と、第２時点ｔ2における第２偏角θ2との差分を周波数変動指標として算定する。具体的には、数式(6)で表現される複素関数について、第１時点ｔ1と第２時点ｔ2との間における自己相関Ｃ(ｔ)の偏角ＡrgＣ(ｔ)が周波数変動指標として算定される。自己相関Ｃ(ｔ)は、以下の数式(7)で表現される。数式(7)の２πΔｆｔ1が第１偏角θ1であり、数式(7)の２πΔｆｔ2が第２偏角θ2である。数式(7)における自己相関Ｃ(ｔ)の偏角ＡrgＣ(ｔ)である２πΔｆＴが、第１偏角θ1と第２偏角θ2との差分を表す。
【数７】
2019118669.tif 000009
【００３０】
図３の指標算定部９８は、差分算定部９７が算定した周波数変動指標（ＡrgＣ(ｔ)）から脈周期Ｐを算定する。脈周期Ｐは、数式(7)の関係を考慮すると、以下の数式(8)の演算により算定される。
【数８】
2019118669.tif 000010
数式(8)から理解さる通り、脈周期Ｐは、基本周波数ｆ0に対する変動量Δｆの割合で表現される。指標算定部９８は、検出信号Ｓa(ｔ)のサンプリング周期Ｔ毎に脈周期Ｐを算定する。つまり、指標算定部９８は、時間軸上の相異なる時点に対応する複数の脈周期Ｐの時系列（すなわち心拍変動）を特定する。
【００３１】
変動算定部６８は、指標算定部９８が特定した心拍変動から心拍変動指標を算定する。例えば心拍変動の標準偏差を表すＳＤＮＮ(standard deviation of NN intervals)や、交感神経のバランスを表したＬＦ(Low Frequency)／ＨＦ（High Frequency）比が心拍変動指標として例示される。ＳＤＮＮおよびＬＦ／ＨＦ比の算定には、公知の技術が任意に採用される。表示装置２４は、変動算定部６８が算定した心拍変動指標を表示する。
【００３２】
図６は、制御装置２０が実行する処理（以下「脈波解析処理」という）のフローチャートである。例えば、被験者からの指示を契機として脈波解析処理が実行される。図６の脈波解析処理は、周波数解析部６２による基本周波数ｆ0の特定後に、時間軸上のサンプリング周期Ｔ毎に実行される。脈波解析処理を開始すると、成分抽出部６４は、検出信号Ｓa(ｔ)のうち脈波の基本周波数ｆ0を含む解析帯域Ｂ内の成分を観測信号Ｓb(ｔ)として抽出する（Ｓ1）。すなわち、検出信号Ｓa(ｔ)のうち解析帯域Ｂ外に存在するノイズの影響を低減した観測信号Ｓb(ｔ)が生成される。信号処理部９１は、数式(1)で表現される第１周期信号Ｚi(ｔ)と、数式(2)で表現される第２周期信号Ｚq(ｔ)とを生成する（Ｓ2）。第１処理部９３は、第１周期信号Ｚi(ｔ)を観測信号Ｓb(ｔ)に乗算することで、前述の数式(4)で表現される同相信号Ｓi(ｔ)を生成する（Ｓ3）。第２処理部９５は、第２周期信号Ｚq(ｔ)を観測信号Ｓb(ｔ)に乗算することで、前述の数式(5)で表現される直交信号Ｓq(ｔ)を生成する（Ｓ4）。なお、ステップＳ3とステップＳ4との順番の先後は不問である。差分算定部９７は、同相信号Ｓi(ｔ)と直交信号Ｓq(ｔ)とから周波数変動指標を算定する（Ｓ5）。数式(4)の同相信号Ｓi(ｔ)の信号値と数式(5)の直交信号Ｓq(ｔ)の信号値との関係を表す数式(6)の複素関数について、第１時点ｔ1と第２時点ｔ2との間における自己相関Ｃ(ｔ)の偏角ＡrgＣ(ｔ)が周波数変動指標として算定される。自己相関Ｃ(ｔ)は、前述の数式(7)で表現される。指標算定部９８は、周波数変動指標から脈周期Ｐを算定する。脈周期Ｐの算定には、前述の数式(8)が利用される。
【００３３】
以上の説明から理解される通り、第１実施形態では、検出信号Ｓa(ｔ)のうち脈波の基本周波数ｆ0を含む解析帯域Ｂ内の成分が観測信号Ｓb(ｔ)として抽出され、当該観測信号Ｓb(ｔ)を利用して脈周期Ｐが算定される。以上の構成によれば、例えば検出信号Ｓa(ｔ)から（つまり解析帯域Ｂ以外の成分も含む信号から）脈周期Ｐを算定する構成と比較して、生体の体動に起因したノイズや装置の回路ノイズの影響を低減した脈周期Ｐを算定することができる。したがって、検出信号Ｓa(ｔ)から脈周期Ｐを高精度に算定することができる。ひいては、脈周期Ｐから高精度に周波数変動指標を算定することが可能である。
【００３４】
また、上述の通り、周波数解析部６２による基本周波数ｆ0の特定と、成分抽出部６４による観測信号Ｓb(ｔ)の抽出には共通の検出信号Ｓa(ｔ)が利用される。すなわち、周波数変動指標の算定に利用される検出信号Ｓa(ｔ)の解析により基本周波数ｆ0が特定される。以上の構成によれば、周波数変動指標の算定に利用する検出信号Ｓa(ｔ)を基本周波数ｆ0の特定に流用できるという利点がある。
【００３５】
＜第２実施形態＞
本発明の第２実施形態を説明する。なお、以下に例示する各形態において作用または機能が第１実施形態と同様である要素については、第１実施形態の説明で使用した符号を流用して各々の詳細な説明を適宜に省略する。
【００３６】
第１実施形態では、成分抽出部６４における解析帯域Ｂの帯域幅は固定である。それに対して、第２実施形態では、成分抽出部６４における解析帯域Ｂの帯域幅は可変である。図７は、第２実施形態に係る制御装置２０の機能に着目した構成図である。第２実施形態の制御装置２０は、第１実施形態の制御装置２０が実現した機能（周波数解析部６２，成分抽出部６４，演算処理部６６，変動算定部６８）に加えて、帯域制御部６９の機能を実現する。
【００３７】
帯域制御部６９は、指標算定部９８が算定した脈周期Ｐに応じて、成分抽出部６４における解析帯域Ｂの帯域幅を制御する。具体的には、脈周期Ｐの時系列（つまり心拍変動）の散布度（脈周期Ｐのバラツキを表す指標）に応じて、解析帯域Ｂの帯域幅が制御される。心拍変動の散布度が大きい場合には帯域幅が広くなるように制御され、心拍変動の散布度が小さい場合には帯域幅が狭くなるように制御される。成分抽出部６４は、帯域制御部６９により設定された帯域幅の解析帯域Ｂを利用して、検出信号Ｓa(ｔ)から観測信号Ｓb(ｔ)を抽出する。
【００３８】
第２実施形態においても第１実施形態と同様の効果が実現される。第２実施形態では、脈周期Ｐに応じて、成分抽出部６４における解析帯域Ｂの帯域幅が制御されるから、解析帯域Ｂの帯域幅が固定である第１実施形態の構成と比較して、脈周期Ｐの変動する範囲に応じた適切な解析帯域Ｂが設定される。
【００３９】
第１実施形態および第２実施形態では、周波数解析部６２が１個の単位期間について特定した基本周波数ｆ0を、成分抽出部６４における解析帯域Ｂの設定と、同相信号Ｓi(ｔ)および直交信号Ｓq(ｔ)の周波数の設定とに適用する構成を便宜的に例示した。実際には、時間軸上の単位期間毎に周波数解析部６２が基本周波数ｆ0を反復的に特定し、基本周波数ｆ0の特定毎に、成分抽出部６４における解析帯域Ｂの設定と、第１周期信号Ｚi(ｔ)および第２周期信号Ｚq(ｔ)の生成とに当該基本周波数ｆ0が適用される。すなわち、基本周波数ｆ0は順次に更新される。
【００４０】
以上に説明した通り、本発明の好適な態様においては、検出信号Ｓa(ｔ)から順次に特定される基本周波数ｆ0に応じて成分抽出部６４の解析帯域Ｂを設定する構成（以下「構成Ａ」という）を採用する。構成Ａを採用することで実際の脈波解析装置（以下「実製品」という）から観測できる挙動を以下に説明する。
【００４１】
図８に例示される通り、信号値が周期的に変動する検出信号Ｓa(ｔ)を実製品に付与した場合を想定する。検出信号Ｓa(ｔ)の周波数は、時点ｔa1まで１Ｈｚに維持され、時点ｔa1において１Ｈｚから２Ｈｚに不連続に変化し、時点ｔa1の後の時点ｔb1において２Ｈｚから１Ｈｚに不連続に変化する。
【００４２】
実製品では、所定長の単位期間毎に基本周波数ｆ0が更新される。したがって、検出信号Ｓa(ｔ)の周波数が１Ｈｚから２Ｈｚに変動する時点ｔa1では基本周波数ｆ0は未だ更新されず、時点ｔa1から時間が経過した時点ｔa2（すなわち単位期間の端点）が到来することで基本周波数ｆ0が更新される。したがって、時点ｔa1から時点ｔa2までの期間内においては、検出信号Ｓa(ｔ)の周波数は２Ｈｚであるけれども、成分抽出部６４の解析帯域Ｂの中心周波数は依然として１Ｈｚに維持されるという状況になる。すなわち、観測信号Ｓb(ｔ)は、本来の周期的な変動である２Ｈｚの成分は低減されてノイズが優勢な信号となる。したがって、実製品Ａにおいては、指標算定部９８が算定する脈周期Ｐは、ノイズに起因した異常値となる。同様の理由により、図８から理解される通り、時点ｔb1から時点ｔb2の間においても脈周期Ｐは異常値となる。なお、検出信号Ｓa(ｔ)は、実際には生体の脈波を表す信号であるから、周波数が不連続に変動することはない。したがって、実際の測定の場面では、図８に例示したようなノイズは観測されず、脈周期Ｐを高精度に測定できるという所期の効果は実現される。
【００４３】
他方、構成Ａを採用することなく脈周期Ｐを測定する構成（例えば成分抽出部６４を省略した構成）では、検出信号Ｓa(ｔ)の周波数の変化には迅速に追従して脈周期Ｐが測定される。したがって、時点ｔa1や時点ｔb1の直後の期間において脈周期Ｐは異常値とならない。以上の説明から理解される通り、検出信号Ｓa(ｔ)の周波数が不連続に変化した直後に脈周期Ｐが異常値となる場合、実製品は構成Ａを採用している可能性が高いと言える。なお、構成Ａを採用しない構成において、検出信号Ｓa(ｔ)の周波数が不連続に変化する時点の直後において異常値の発生を回避できるとは言っても、構成Ａを採用していない以上、脈周期Ｐの測定値にはノイズに起因した誤差が発生する。したがって、脈周期Ｐの高精度な測定という観点からすると、構成Ａを採用した構成が、構成Ａを採用しない構成と比較して有利であるという事情は変わらない。なお、実製品に付与する検出信号Ｓa(ｔ)は、図８の例示に限定されない。例えば周波数が経時的に変化するチャープ信号を検出信号Ｓa(ｔ)として利用してもよい。
【００４４】
なお、以上の説明では、脈周期Ｐに着目したが、実製品が構成Ａを採用しているか否かを判断するための指標は、脈周期Ｐに限定されない。脈周期Ｐが異常値を取る場合は、脈周期Ｐの時系列（心拍変動）から算定される心拍変動指標も異常値をとり得る。したがって、心拍変動指標を利用して構成Ａを採用しているか否かを判断してもよい。具体的には、時点ｔa1から時点ｔa2までの期間や時点ｔb1から時点ｔb2までの期間において心拍変動指標が異常値となる場合には構成Ａを採用している可能性が高い。
【００４５】
＜第３実施形態＞
図９は、第３実施形態における脈波解析装置１００の使用例を示す模式図である。図９に例示される通り、脈波解析装置１００は、相互に別体で構成された検出ユニット７１と表示ユニット７２とを具備する。検出ユニット７１は、前述の各形態で例示した検出装置３０を具備する。図１０には、被験者の上腕に装着される形態の検出ユニット７１が例示されている。図１０に例示される通り、被験者の手首に装着される形態の検出ユニット７１も好適である。
【００４６】
表示ユニット７２は、前述の各形態で例示した表示装置２４を具備する。例えば携帯電話機またはスマートフォン等の情報端末が表示ユニット７２の好適例である。ただし、表示ユニット７２の具体的な形態は任意である。例えば、被験者が携帯可能な腕時計型の情報端末、または、脈波解析装置１００の専用の情報端末を表示ユニット７２として利用してもよい。
【００４７】
検出信号Ｓa(ｔ)から心拍変動指標を算定するための要素（以下「演算部」という）は、例えば表示ユニット７２に搭載される。演算部は、図３に例示された要素（周波数解析部６２，成分抽出部６４，演算処理部６６，変動算定部６８）を包含する。演算部は、例えば表示ユニット７２に搭載される。検出ユニット７１の検出装置３０が生成した検出信号Ｓa(ｔ)が有線または無線で表示ユニット７２に送信される。表示ユニット７２の演算部は、検出信号Ｓa(ｔ)から心拍変動指標を算定して表示装置２４に表示する。
【００４８】
なお、演算部を検出ユニット７１に搭載してもよい。演算部は、検出装置３０が生成した検出信号Ｓa(ｔ)から心拍変動指標を算定し、当該心拍変動指標を表示するためのデータを表示ユニット７２に有線または無線で送信する。表示ユニット７２の表示装置２４は、検出ユニット７１から受信したデータが示す心拍変動指標を表示する。
【００４９】
＜変形例＞
以上に例示した各形態は多様に変形され得る。具体的な変形の態様を以下に例示する。以下の例示から任意に選択された２以上の態様を適宜に併合することも可能である。
【００５０】
（１）前述の各形態の脈波解析装置１００は変動算定部６８を具備する構成を例示したが、変動算定部６８を脈波解析装置１００から省略してもよい。例えば指標算定部９８が算定した脈周期Ｐが表示装置２４に表示される。また、脈波解析装置１００が算定した脈周期Ｐを脈波解析装置１００とは別体の情報端末（例えば携帯電話機またはスマートフォン）に送信して、当該情報端末が脈周期Ｐから心拍変動指標を算定してもよい。
【００５１】
（２）前述の各形態では、指標算定部９８は、生体の脈周期Ｐを算定したが、指標算定部９８が算定する指標は脈周期Ｐに限定されない。例えば脈周期Ｐの段階を複数の候補（例えば、異常／高目／通常、など）から特定して被験者に報知することも可能である。または、基準となる脈周期Ｐ（例えば基本周波数ｆ0の逆数に相当する周期）からの変動量を算定して被験者に報知してもよい。以上の説明から理解される通り、指標算定部９８が算定する指標は、生体の脈周期Ｐに関する周期指標として包括的に表現される。
【００５２】
（３）前述の各形態では、数式(6)の複素関数について、第１時点ｔ1と第２時点ｔ2との間における自己相関Ｃ(ｔ)の偏角ＡrgＣ(ｔ)を周波数変動指標として算定したが、周波数変動指標の算定方法は以上の例示に限定されない。数式(6)の複素関数について、第１偏角θ1と第２偏角θ2とをそれぞれ算定し、当該第１偏角θ1と第２偏角θ2との差分を周波数変動指標として算定してもよい。ただし、自己相関Ｃ(ｔ)の偏角ＡrgＣ(ｔ)を周波数変動指標として算定する前述の各形態によれば、第１偏角θ1と第２偏角θ2とをそれぞれ算定することが不要になる。
【００５３】
（４）前述の各形態では、周波数解析部６２による基本周波数ｆ0の特定と、成分抽出部６４による観測信号Ｓb(ｔ)の抽出とに共通の検出信号Ｓa(ｔ)を利用したが、周波数解析部６２による基本周波数ｆ0の特定と、成分抽出部６４による観測信号Ｓb(ｔ)の抽出とに別個の検出信号Ｓa(ｔ)を利用してもよい。検出装置３０を利用しない方法で基本周波数ｆ0を取得してもよい。例えば、脈波解析装置１００とは別体の情報端末から基本周波数ｆ0を取得して周期指標の特定に利用してもよい。
【００５４】
（５）前述の各形態では、単体の機器として構成された脈波解析装置１００を例示したが、以下の例示の通り、脈波解析装置１００の複数の要素は相互に別体の装置として実現され得る。なお、以下の説明では、検出信号Ｓa(ｔ)から心拍変動指標を算定する要素を「演算部２７」と表記する。演算部２７は、例えば、図３に例示された要素（周波数解析部６２，成分抽出部６４，演算処理部６６，変動算定部６８）を包含する。なお、第２実施形態の構成では、帯域制御部６９も演算部２７に含まれる。
【００５５】
前述の各形態では、検出装置３０を具備する脈波解析装置１００を例示したが、図１１に例示される通り、検出装置３０を脈波解析装置１００とは別体とした構成も想定される。検出装置３０は、例えば被験者の手首や上腕等の測定部位Ｍに装着される可搬型の光学センサーモジュールである。脈波解析装置１００は、例えば携帯電話機またはスマートフォン等の情報端末で実現される。腕時計型の情報端末で脈波解析装置１００を実現してもよい。検出装置３０が生成した検出信号Ｓa(ｔ)が有線または無線で脈波解析装置１００に送信される。脈波解析装置１００の演算部２７は、検出信号Ｓa(ｔ)から心拍変動指標を算定して表示装置２４に表示する。以上の説明から理解される通り、検出装置３０は脈波解析装置１００から省略され得る。
【００５６】
前述の各形態では、表示装置２４を具備する脈波解析装置１００を例示したが、図１２に例示される通り、表示装置２４を脈波解析装置１００とは別体とした構成も想定される。脈波解析装置１００の演算部２７は、検出信号Ｓa(ｔ)から心拍変動指標を算定し、当該心拍変動指標を表示するためのデータを表示装置２４に送信する。表示装置２４は、専用の表示機器であってもよいが、例えば、携帯電話機もしくはスマートフォン等の情報端末、または、被験者が携帯可能な腕時計型の情報端末に搭載されてもよい。脈波解析装置１００の演算部２７が算定した心拍変動指標は、有線または無線により表示装置２４に送信される。表示装置２４は、脈波解析装置１００から受信した心拍変動指標を表示する。以上の説明から理解される通り、表示装置２４は脈波解析装置１００から省略され得る。
【００５７】
図１３に例示される通り、検出装置３０および表示装置２４を脈波解析装置１００（演算部２７）とは別体とした構成も想定される。例えば、脈波解析装置１００（演算部２７）が、携帯電話機やスマートフォン等の情報端末に搭載される。
【００５８】
（６）前述の各形態では、筐体部１２とベルト１４とを具備する腕時計型の脈波解析装置１００を例示したが、脈波解析装置１００の具体的な形態は任意である。例えば、被験者の身体に貼付可能なパッチ型、被験者の耳部に装着可能な耳装着型、被験者の指先に装着可能な指装着型（例えば着爪型）、または、被験者の頭部に装着可能な頭部装着型など、任意の形態の脈波解析装置１００が採用され得る。
【００５９】
（７）前述の各形態では、被験者の心拍変動指標を表示装置２４に表示したが、心拍変動指標を被験者に報知するための構成は以上の例示に限定されない。例えば、心拍変動指標を音声で被験者に報知することも可能である。被験者の耳部に装着可能な耳装着型の脈波解析装置１００においては、心拍変動指標を音声で報知する構成が特に好適である。また、心拍変動指標を被験者に報知することは必須ではない。例えば、脈波解析装置１００が算定した心拍変動指標を通信網から他の通信装置に送信してもよい。また、脈波解析装置１００の記憶装置２２や脈波解析装置１００に着脱可能な可搬型の記録媒体に心拍変動指標を格納してもよい。
【００６０】
（８）前述の各形態に係る脈波解析装置１００は、前述の例示の通り、制御装置２０とプログラムとの協働により実現される。本発明の好適な態様に係るプログラムは、コンピューターが読取可能な記録媒体に格納された形態で提供されてコンピューターにインストールされ得る。また、配信サーバーが具備する記録媒体に格納されたプログラムを、通信網を介した配信の形態でコンピューターに提供することも可能である。記録媒体は、例えば非一過性（non-transitory）の記録媒体であり、ＣＤ-ＲＯＭ等の光学式記録媒体（光ディスク）が好例であるが、半導体記録媒体または磁気記録媒体等の公知の任意の形式の記録媒体を包含し得る。なお、非一過性の記録媒体とは、一過性の伝搬信号（transitory, propagating signal）を除く任意の記録媒体を含み、揮発性の記録媒体を除外するものではない。
【符号の説明】
【００６１】
１００…脈波解析装置、１２…筐体部、１４…ベルト、２０…制御装置、２２…記憶装置、２４…表示装置、３０…検出装置、３１…発光部、３２…受光部、３３…駆動回路、３４…出力回路、６２…周波数解析部、６４…成分抽出部、６６…演算処理部、６８…変動算定部、６９…帯域制御部、９１…信号処理部、９３…第１処理部、９５…第２処理部、９７…差分算定部、９８…指標算定部、９３１，９５１…乗算部、９３３，９５３…低域通過部。
【図１】
2019118669.tif 000011
【図２】
2019118669.tif 000012
【図３】
2019118669.tif 000013
【図４】
2019118669.tif 000014
【図５】
2019118669.tif 000015
【図６】
2019118669.tif 000016
【図７】
2019118669.tif 000017
【図８】
2019118669.tif 000018
【図９】
2019118669.tif 000019
【図１０】
2019118669.tif 000020
【図１１】
2019118669.tif 000021
【図１２】
2019118669.tif 000022
【図１３】
2019118669.tif 000023