JP6488722B2 - 脈波検出装置、脈波検出方法及び脈波検出プログラム - Google Patents

Description

本発明は、脈波検出装置、脈波検出方法及び脈波検出プログラムに関する。

被験者が撮影された画像から血液の体積の変動、いわゆる脈波を検出する技術が知られている。例えば、脈波の一態様である心拍数を検出する技術の一例としては、心電計の電極を生体に装着することによって計測された心電波形のピーク、例えばＰ波やＲ波などを用いて心拍数を検出する心電図法が挙げられる。他の一例としては、指や耳たぶなどの末梢血管に赤外線を照射し、その反射光が血流および吸光特性によって周期的に変動する光学的な変化から心拍とほぼ等価な脈拍を検出する光電脈波法が挙げられる。

これら心電図法や光電脈波法を用いる場合には、電極を生体に装着したり、あるいは生体に感光面を密着させたりするので、計測器具が生体に接触しないと検出が困難である上、計測器具を装着した状態で日常を生活するのは煩わしい場合がある。

このことから、生体に計測器具が接触しない状態で脈波を測定するために、被験者が撮影された画像を用いた心拍測定方法が提案されている。かかる心拍測定方法では、カメラによって被験者の顔が撮影された画像の信号成分に対し、例えば、ＦＦＴ（Fast Fourier Transform）や最大エントロピー法によってパワースペクトルを求め、かかるパワースペクトルの周波数のピーク値から心拍数を得る。

特開２００２−３３０９３５号公報 特開２０１３−１３８７０５号公報

しかしながら、上記の技術では、脈波の周期変動に追従できない場合がある。

すなわち、上記の心拍測定方法は、ＦＦＴや最大エントロピー法を用いて、一定区間における周波数の分布を計算し、どの周波数成分が最大であるのかを求めるものである。それ故、実際の脈拍数が一拍ごとに変動する場合、その変動に追従することができない。

１つの側面では、本発明は、脈波の周期変動に追従できる脈波検出装置、脈波検出方法及び脈波検出プログラムを提供することを目的とする。

一態様の脈波検出装置は、画像を取得する取得部と、前記画像から脈波波形を検出する波形検出部と、前記脈波波形との間で互いの一致度が大きくなるように、所定の周期的なモデル波形であって周期変動を定義するパラメータを含むモデル波形に設定されるパラメータを更新する更新部と、前記パラメータの更新後のモデル波形と前記脈波波形との一致度が所定の閾値以上である場合に、前記脈波波形を出力する出力部とを有する。

脈波の周期変動に追従できる。

図１は、実施例１に係る脈波検出装置の機能的構成を示すブロック図である。 図２は、顔画像の一例を示す図である。 図３は、波形検出部の機能的構成を示すブロック図である。 図４は、初期のモデル波形の一例を示す図である。 図５は、パラメータの更新結果の一例を示す図である。 図６は、実施例１に係る波形検出処理の手順を示すフローチャートである。 図７は、実施例１に係るパラメータの更新処理の手順を示すフローチャートである。 図８は、実施例１及び実施例２に係る脈波検出プログラムを実行するコンピュータのハードウェア構成例を示す図である。

以下に添付図面を参照して本願に係る脈波検出装置、脈波検出方法及び脈波検出プログラムについて説明する。なお、この実施例は開示の技術を限定するものではない。そして、各実施例は、処理内容を矛盾させない範囲で適宜組み合わせることが可能である。

［脈波検出装置の構成］
まず、本実施例に係る脈波検出装置の機能的構成について説明する。図１は、実施例１に係る脈波検出装置の機能的構成を示すブロック図である。図１に示す脈波検出装置１０は、太陽光や室内光などの一般の環境光の下で生体に計測器具を接触させずに、被験者が撮影された画像を用いて被験者の脈波、すなわち心臓の拍動に伴う血液の体積の変動を検出する脈波検出処理を実行するものである。

かかる脈波検出装置１０は、上記の脈波検出の一環として、生体が撮影された画像から検出された脈波波形との間で一致度が最大となるように、モデル波形の周期変動に関するパラメータを更新し、その周期変動から脈波波形や脈拍数を求める。これによって、脈波の周期変動に追従することを目指す。

かかる脈波検出装置１０は、一態様として、上記の脈波検出処理がパッケージソフトウェアやオンラインソフトウェアとして提供される信号処理プログラムを所望のコンピュータにインストールさせることによって実装できる。例えば、スマートフォン、携帯電話機やＰＨＳ（Personal Handyphone System）などの移動体通信端末のみならず、移動体通信網に接続する能力を持たないデジタルカメラ、タブレット端末やスレート端末を含む携帯端末装置に上記の脈波検出プログラムをインストールさせる。これによって、携帯端末装置を脈波検出装置１０として機能させることができる。なお、ここでは、脈波検出装置１０の実装例として携帯端末装置を例示したが、パーソナルコンピュータを始めとする据置き型の端末装置に脈波検出プログラムをインストールさせることもできる。

図１に示すように、脈波検出装置１０は、カメラ１１と、取得部１３と、抽出部１４と、波形検出部１５と、転換判定部１６と、初期設定部１７と、更新部１８と、脈波検出部１９とを有する。

かかる脈波検出装置１０は、図１に示した機能部以外にも既知のコンピュータが有する各種の機能部を有することとしてもかまわない。例えば、脈波検出装置１０が据置き端末として実装される場合には、キーボード、マウスやディスプレイなどの入出力デバイスをさらに有することとしてもよい。また、脈波検出装置１０がタブレット端末やスレート端末として実装される場合には、タッチパネルをさらに有することとしてもよい。また、脈波検出装置１０が移動体通信端末として実装される場合には、アンテナ、移動体通信網に接続する無線通信部、ＧＰＳ（Global Positioning System）受信機などの機能部をさらに有していてもかまわない。

カメラ１１は、ＣＣＤ（Charge Coupled Device）やＣＭＯＳ（Complementary Metal Oxide Semiconductor）などの撮像素子を搭載する撮像装置である。

例えば、カメラ１１には、Ｒ（red）、Ｇ（green）、Ｂ（blue）など３種以上の受光素子を搭載することができる。かかるカメラ１１の実装例としては、デジタルカメラやＷｅｂカメラを外部端子を介して接続することとしてもよいし、カメラが出荷時から搭載されている場合にはそのカメラを流用できる。なお、ここでは、脈波検出装置１０がカメラ１１を有する場合を例示したが、ネットワークまたは記憶デバイスを経由して画像を取得できる場合には、必ずしも脈波検出装置１０がカメラ１１を有さずともよい。

ここで、脈波検出装置１０は、上記の脈波検出処理を実行するアプリケーションプログラムがプリインストールまたはインストールされている場合に、カメラ１１によって脈波を検出し易い被験者の画像が撮像されるように画像の撮影操作を案内することができる。なお、以下では、上記のアプリケーションプログラムのことを「脈波検出用アプリ」と記載する場合がある。

かかる脈波検出用アプリは、図示しない入力デバイスを介して起動されると、カメラ１１を起動する。これを受けて、カメラ１１は、カメラ１１の撮影範囲に収容された被写体の撮影を開始する。このとき、被験者の顔が映る画像を撮影させる場合には、脈波検出用アプリは、カメラ１１が撮影する画像を図示しない表示デバイスに表示しつつ、被験者の鼻を映す目標位置を照準として表示させることもできる。これによって、被験者の眼、耳、鼻や口などの顔パーツの中でも被験者の鼻が撮影範囲の中心部分に収まった画像が撮影できるようにする。そして、脈波検出用アプリは、カメラ１１によって被験者の顔が撮影された画像を取得部１３へ保存する。なお、以下では、顔が映った画像のことを「顔画像」と記載する場合がある。

取得部１３は、画像を取得する処理部である。

一実施形態として、取得部１３は、カメラ１１によって撮像された顔画像を取得する。他の実施形態として、取得部１３は、顔画像を蓄積するハードディスクや光ディスクなどの補助記憶装置またはメモリカードやＵＳＢ（Universal Serial Bus）メモリなどのリムーバブルメディアから画像を取得することもできる。更なる実施形態として、取得部１３は、外部装置からネットワークを介して受信することによって顔画像を取得することもできる。なお、取得部１３は、ＣＣＤやＣＭＯＳなどの撮像素子による出力から得られる２次元のビットマップデータやベクタデータなどの画像データを用いて処理を実行する場合を例示したが、１つのディテクタから出力される信号をそのまま取得して後段の処理を実行させることとしてもよい。

抽出部１４は、取得部１３によって取得された画像から生体領域を抽出する処理部である。

一実施形態として、抽出部１４は、顔画像から所定の顔パーツを基準とする生体領域を抽出する。例えば、抽出部１４は、顔画像にテンプレートマッチング等の画像処理を実行することによって被験者の眼、耳、鼻や口などの顔パーツのうち特定の顔パーツ、すなわち被験者の鼻を検出する。その上で、抽出部１４は、被験者の鼻を中心とし、中心から所定の範囲に含まれる生体領域を抽出する。これによって、被験者の鼻、鼻の周辺に位置する頬の一部の顔中心部分を含んだ生体領域の画像が脈波の検出に使用する画像として抽出される。その後、抽出部１４は、原画像から抽出した生体領域の画像を波形検出部１５へ出力する。

図２は、顔画像の一例を示す図である。図２には、画像に映る被験者の眼、鼻及び口の一部または全部を含む領域が９つに分割されたブロックが図示されている。図２に示すブロックのうち上段の左及び右のブロックには、被験者の眼が映っている。これらのブロックの画像を検出に用いた場合には、眼の瞬きがノイズとなって心拍数の検出精度の低下を招く場合がある。また、図２に示すブロックのうち下段の３つのブロックには、被験者の口が映っている。これらのブロックの画像を検出に用いた場合には、口の動きがノイズとなって心拍数の検出精度の低下を招く場合がある。一方、図２に示す中段の真ん中のブロック、すなわち斜線の塗りつぶしが図示されたブロックは、眼や口が映るブロックから隔てられており、他のブロックに比べてノイズとなる成分が映っている可能性が低いので、良好な検出結果を期待できる。これらのことから、抽出部１４は、原画像から図２に示す中段の真ん中のブロックの画像を生体領域の画像として抽出する。

上記の生体領域の画像を抽出後に、抽出部１４は、生体領域に含まれる各画素が持つ画素値に所定の統計処理を実行する。例えば、抽出部１４は、生体領域に含まれる各画素が持つ画素値を波長成分ごとに平均する。この他、平均値以外にも、中央値や最頻値を計算することとしてもよく、また、相加平均以外にも任意の平均処理、例えば加重平均や移動平均などを実行することもできる。これによって、生体領域に含まれる各画素が持つ画素値の平均値が当該生体領域を代表する代表値として波長成分ごとに算出される。

波形検出部１５は、脈波の検出対象とする領域に含まれる各画素の波長成分別の代表値の信号から、各波長成分の間で脈波が採り得る脈波周波数帯以外の特定周波数帯の成分が互いに相殺された信号の波形を検出する処理部である。

一実施形態として、波形検出部１５は、画像に含まれる３つの波長成分、すなわちＲ成分、Ｇ成分およびＢ成分のうち血液の吸光特定が異なるＲ成分とＧ成分の２つの波長成分の代表値の時系列データを用いて、波形を検出する。

これを説明すると、顔表面には、毛細血管が流れており、心拍により血管に流れる血流が変化すると、血流で吸収される光量も心拍に応じて変化するため、顔からの反射によって得られる輝度も心拍に伴って変化する。かかる輝度の変化量は小さいが、顔領域全体の平均輝度を求めると、輝度の時系列データには脈波成分が含まれる。ところが、輝度は、脈波以外に体動等によっても変化し、これが、脈波検出のノイズ成分、いわゆる体動アーチファクトとなる。そこで、血液の吸光特性の異なる２種類以上の波長、例えば吸光特性が高いＧ成分(５２５ｎｍ程度)、吸光特性が低いＲ成分(７００ｎｍ程度)で脈波を検出する。心拍は、０．５Ｈｚ〜４Ｈｚ、１分あたりに換算すれば３０ｂｐｍ〜２４０ｂｐｍの範囲であるので、それ以外の成分はノイズ成分とみなすことができる。ノイズには、波長特性は無い、あるいはあっても極小であると仮定すると、Ｇ信号およびＲ信号の間で０．５Ｈｚ〜４Ｈｚ以外の成分は等しいはずであるが、カメラの感度差により大きさが異なる。それゆえ、０．５Ｈｚ〜４Ｈｚ以外の成分の感度差を補正して、Ｇ成分からＲ成分を減算すれば、ノイズ成分は除去されて脈波成分のみを取り出すことができる。

例えば、Ｇ成分及びＲ成分は、下記の式（１）および下記の式（２）によって表すことができる。下記の式（１）における「Ｇｓ」は、Ｇ信号の脈波成分を指し、「Ｇｎ」は、Ｇ信号のノイズ成分を指し、また、下記の式（２）における「Ｒｓ」は、Ｒ信号の脈波成分を指し、「Ｒｎ」は、Ｒ信号のノイズ成分を指す。また、ノイズ成分は、Ｇ成分およびＲ成分の間で感度差があるので、感度差の補正係数ｋは、下記の式（３）によって表される。

Ｇａ＝Ｇｓ＋Ｇｎ・・・（１）
Ｒａ＝Ｒｓ＋Ｒｎ・・・（２）
ｋ＝Ｇｎ／Ｒｎ・・・（３）

感度差を補正してＧ成分からＲ成分を減算すると、脈波成分Ｓは、下記の式（４）となる。これを上記の式（１）及び上記の式（２）を用いて、Ｇｓ、Ｇｎ、Ｒｓ及びＲｎによって表される式へ変形すると、下記の式（５）となり、さらに、上記の式（３）を用いて、ｋを消し、式を整理すると下記の式（６）が導出される。

Ｓ＝Ｇａ−ｋＲａ・・・（４）
Ｓ＝Ｇｓ＋Ｇｎ−ｋ（Ｒｓ＋Ｒｎ）・・・（５）
Ｓ＝Ｇｓ−（Ｇｎ／Ｒｎ）Ｒｓ・・・（６）

ここで、Ｇ信号およびＲ信号は、吸光特性が異なり、Ｇｓ＞（Ｇｎ／Ｒｎ）Ｒｓである。したがって、上記の式（６）によってノイズが除去された脈波成分Ｓを算出することができる。

図３は、波形検出部１５の機能的構成を示すブロック図である。図３に示すように、波形検出部１５は、ＢＰＦ（Band-Pass Filter）１５２Ｒ及び１５２Ｇと、抽出部１５３Ｒ及び１５３Ｇと、ＬＰＦ（Low-Pass Filter）１５４Ｒ及び１５４Ｇと、算出部１５５と、ＢＰＦ１５６Ｒ及び１５６Ｇと、乗算部１５７と、演算部１５８とを有する。なお、図３には、時系列空間にてノイズ成分をキャンセルして脈波を検出する場合の機能的構成が例示されている。

例えば、抽出部１４から波形検出部１５には、生体領域に含まれる各画素が持つＲ成分の画素値の代表値を信号値とするＲ信号の時系列データが入力されるとともに、生体領域に含まれる各画素が持つＧ成分の画素値の代表値を信号値とするＧ信号の時系列データが入力される。このうち、生体領域のＲ信号は、波形検出部１５内のＢＰＦ１５２Ｒ及びＢＰＦ１５６Ｒへ入力されるとともに、生体領域のＧ信号は、波形検出部１５内のＢＰＦ１５２Ｇ及びＢＰＦ１５６Ｇへ入力される。

ＢＰＦ１５２Ｒ、ＢＰＦ１５２Ｇ、ＢＰＦ１５６Ｒ及びＢＰＦ１５６Ｇは、いずれも所定の周波数帯の信号成分だけを通過させてそれ以外の周波数帯の信号成分を除去するバンドパスフィルタである。これらＢＰＦ１５２Ｒ、ＢＰＦ１５２Ｇ、ＢＰＦ１５６Ｒ及びＢＰＦ１５６Ｇは、ハードウェアによって実装されることとしてもよいし、ソフトウェアによって実装されることとしてもよい。

これらＢＰＦが通過させる周波数帯の違いについて説明する。ＢＰＦ１５２Ｒ及びＢＰＦ１５２Ｇは、ノイズ成分が他の周波数帯よりも顕著に現れる特定周波数帯の信号成分を通過させる。

かかる特定周波数帯は、脈波が採り得る周波数帯との間で比較することによって定めることができる。脈波が採り得る周波数帯の一例としては、０．５Ｈｚ以上４Ｈｚ以下である周波数帯、１分あたりに換算すれば３０ｂｐｍ以上２４０ｂｐｍ以下である周波数帯が挙げられる。このことから、特定周波数帯の一例としては、脈波として計測され得ない０．５Ｈｚ未満及び４Ｈｚ超過の周波数帯を採用することができる。また、特定周波数帯は、脈波が採り得る周波数帯との間でその一部が重複することとしてもよい。例えば、脈波として計測されることが想定しづらい０．７Ｈｚ〜１Ｈｚの区間で脈波が採り得る周波数帯と重複することを許容し、１Ｈｚ未満及び４Ｈｚ以上の周波数帯を特定周波数帯として採用することもできる。また、特定周波数帯は、１Ｈｚ未満及び４Ｈｚ以上の周波数帯を外縁とし、ノイズがより顕著に現れる周波数帯に絞ることもできる。例えば、ノイズは、脈波が採り得る周波数帯よりも高い高周波数帯よりも、脈波が採り得る周波数帯よりも低い低周波数帯でより顕著に現れる。このため、１Ｈｚ未満の周波数帯に特定周波数帯を絞ることもできる。また、空間周波数がゼロである直流成分の近傍には、各成分の撮像素子の感度の差が多く含まれるので、３ｂｐｍ以上６０ｂｐｍ未満の周波数帯に特定周波数帯を絞ることもできる。さらに、人の体の動き、例えば瞬きや体の揺れの他、環境光のチラツキなどのノイズが現れやすい３ｂｐｍ以上２０ｂｐｍ未満の周波数帯に特定周波数帯を絞ることもできる。

ここでは、一例として、ＢＰＦ１５２Ｒ及びＢＰＦ１５２Ｇが特定周波数帯として３ｂｐｍ以上２０ｂｐｍ未満の周波数帯の信号成分を通過させる場合を想定して以下の説明を行う。なお、ここでは、特定周波数帯の信号成分を抽出するために、バンドパスフィルタを用いる場合を例示したが、一定の周波数未満の周波数帯の信号成分を抽出する場合などには、ローパスフィルタを用いることもできる。

一方、ＢＰＦ１５６Ｒ及びＢＰＦ１５６Ｇは、脈波が採り得る周波数帯、例えば４２ｂｐｍ以上２４０ｂｐｍ未満の周波数帯の信号成分を通過させる。なお、以下では、脈波が採り得る周波数帯のことを「脈波周波数帯」と記載する場合がある。

抽出部１５３Ｒは、Ｒ信号の特定周波数帯の信号成分の絶対強度値を抽出する。例えば、抽出部１５３Ｒは、Ｒ成分の特定周波数帯の信号成分をべき乗する乗算処理を実行することによって特定周波数帯の信号成分の絶対強度値を抽出する。また、抽出部１５３Ｇは、Ｇ信号の特定周波数帯の信号成分の絶対強度値を抽出する。例えば、抽出部１５３Ｇは、Ｇ成分の特定周波数帯の信号成分をべき乗する乗算処理を実行することによって特定周波数帯の信号成分の絶対強度値を抽出する。

ＬＰＦ１５４Ｒ及びＬＰＦ１５４Ｇは、特定周波数帯の絶対強度値の時系列データに対し、時間変化に応答させる平滑化処理を実行するローパスフィルタである。これらＬＰＦ１５４Ｒ及びＬＰＦ１５４Ｇは、ＬＰＦ１５４Ｒへ入力される信号がＲ信号であり、ＬＰＦ１５４Ｇへ入力される信号がＧ信号である以外に違いはない。かかる平滑化処理によって、特定周波数帯の絶対値強度Ｒ´ｎ及びＧ´ｎが得られる。

算出部１５５は、ＬＰＦ１５４Ｇによって出力されたＧ信号の特定周波数帯の絶対値強度Ｇ´ｎを、ＬＰＦ１５４Ｒによって出力されたＲ信号の特定周波数帯の絶対値強度Ｒ´ｎで除する除算「Ｇ´ｎ／Ｒ´ｎ」を実行する。これによって、感度差の補正係数ｋを算出する。

乗算部１５７は、ＢＰＦ１５６Ｒによって出力されたＲ信号の脈波周波数帯の信号成分に算出部１５５によって算出された補正係数ｋを乗算する。

演算部１５８は、乗算部１５７によって補正係数ｋが乗算されたＲ信号の脈波周波数帯の信号成分から、ＢＰＦ１５６Ｇによって出力されたＧ信号の脈波周波数帯の信号成分を差し引く演算「ｋ＊Ｒｓ−Ｇｓ」を実行する。かかる演算によって得られた信号の時系列データは、脈波の波形に相当する。

このようにして画像の生体領域から検出された脈波波形が転換判定部１６へ出力される。以下では、脈波波形のｉ番目の時刻をｔ（ｉ）、脈波波形の振幅値をｖ（ｉ）と表す場合がある。さらに、最新の脈波波形のインデックスを「Ｉ」と表す場合がある。なお、脈波波形の振幅値が採取される時刻ｔ（ｉ）の間隔、すなわちフレームレートは、δで一定とする。

図１の説明に戻り、転換判定部１６は、波形検出部１５によって検出される信号の振幅値の正負が転換するか否かを判定する処理部である。

一実施形態として、転換判定部１６は、波形検出部１５によって脈波波形の振幅値が出力される度に、当該脈波波形の振幅値の符号、すなわち正負が転換するか否かを判定する。このとき、転換判定部１６は、下記の式（７）または下記の式（８）のいずれかを満たす場合に、当該振幅値ｖ（ｉ）の正負が転換すると判定する。このように、転換判定部１６は、振幅値の正負が転換する転換点の時刻Ｔ（ｊ）＝ｔ（Ｉ）を図示しない内部メモリ等へ順次記録する。なお、上記の「ｊ」は、何番目の転換点であるかを表し、このうち、最新の転換点のインデックスを「Ｊ」と表すこととする。

v(i-1)<0 and v(i)>=0・・・（７）
v(i-1)>=0 and v(i)<0・・・（８）

初期設定部１７は、モデル波形を定義するパラメータの初期設定を行う処理部である。

ここで言う「モデル波形」とは、波形検出部１５により出力される脈波波形にフィッティングさせる模式的な波形を指す。かかるモデル波形は、人の脈波由来成分に類似した波形で、周期的な波形であることが好ましい。以下では、モデル波形の一例として、周期変動ｓｉｎ波を用いる場合を例示するが、周期が変動する三角波など他の周期変動波形であってもよく、また、複数の周期変動波形の和で構成される合成波であっても構わない。

例えば、周期変動ｓｉｎ波として、時刻に対し、周期が単調に短くなるか、あるいは長くなるモデルを例に挙げる。これは、想定する窓幅、例えば数秒間程度の比較的短い時間の場合、人の脈拍数はおおよそ単調に変動することが多いからである。さらに、ここでは、周期が線形に変動する、下記の式（９）の周期変動ｓｉｎ波をモデル波形として例示する。

ｓｉｎ（α_ｎｔ_ｉ ^２＋β_ｎｔ_ｉ＋γ_ｎ）・・・（９）

上記の式（９）における「α」は、モデル波形の周期変動を指し、「β」は、モデル波形の窓開始時点における周期を指し、また、「γ」は、位相を指す。なお、「ｎ」は、モデル波形を定義する「α」、「γ」及び「β」を含むパラメータが繰り返し演算により更新される回数を指す。

ここで、上記の式（９）において、周期は、「α_ｎｔ_ｉ＋β_ｎ」であり、α_ｎが０以外である場合には、モデル波形の周期が線形に変動することがわかる。例えば、αが正の値にである場合には、時間経過にしたがって周期が単調に増加するモデル波形を定義でき、αが負の値である場合には、時間経過にしたがって周期が単調に減少するモデル波形が定義できることがわかる。

例えば、初期設定部１７は、モデル波形の周期数ｍを設定すると共に、上記のパラメータα、β及びγの初期値を設定する。図４は、初期のモデル波形の一例を示す図である。図４には、波形検出部１５が出力する脈波波形を破線で示すと共に、初期のモデル波形が実線で示されている。図４に示すグラフの縦軸は、振幅を指し、グラフの横軸は、時間を指す。図４に示すように、波形検出部１５が出力する脈波波形から転換点Ｊが転換判定部１６により判定された場合、初期設定部１７は、所定の周期数ｍ、図示の例ではｍ＝３を設定し、当該転換点Ｊが検出された時刻からｍ＊２個前、すなわち６個前の転換点まで遡った時刻をモデル波形とフィッティングさせる窓幅に設定する。すなわち、初期設定部１７は、最新の転換点Ｔ（Ｊ）を窓の終了時刻に設定すると共に、Ｔ（Ｊ−ｍ＊２）が窓の開始時刻に設定する。その上で、初期設定部１７は、図４に示す通り、上記の窓幅にｍ周期分のｓｉｎ波がちょうど収まるように、β_０及びγ_０を設定する。なお、α_０には、初期値の一例としてゼロを設定することとし、初期のモデル波形には周期変動がない状態でパラメータの更新が実行されることとする。

更新部１８は、モデル波形のパラメータを更新する処理部である。

一実施形態として、更新部１８は、波形検出部１５が出力する脈波波形のうち初期設定部１７により設定された窓内の部分波形とモデル波形との一致度が大きくなるように、パラメータα、β及びγを繰り返し演算により更新する。

例えば、窓内の脈波波形のデータ点を（ｔ１，ｘ１），（ｔ２，ｘ２），・・・，（ｔＮ，ｘＮ）とするとき、ｎ反復目のモデル波形は３つのパラメータα_ｎ、β_ｎ及びγ_ｎを持つｓｉｎ（α_ｎｔ_ｉ ^２＋β_ｎｔ_ｉ＋γ_ｎ）で表される。このとき、評価値（一致度）の一例は、窓内の脈波波形とｓｉｎ波の内積（平均を０と仮定した相関係数）により、下記の式（１０）で表すことができる。なお、上記の各データ点における振幅値ｘ１，ｘ２，・・・，ｘＮには、波形検出部１５が出力する振幅値をそのまま用いることもできるが、一例として、−１から＋１までの値に正規化することにより、振幅値が大きく互いの波形が相似でない場合に一致度が高く算出される事態を抑制できる。

そして、上記の式（１０）でｆ_ｎが最大となるα_ｎ、β_ｎ及びγ_ｎは、下記の式（１１）〜式（１３）を満たす。

よって、上記の式（１１）〜式（１３）を満たすα_ｎ、β_ｎ及びγ_ｎを計算することができればよい。このようなα_ｎ、β_ｎ及びγ_ｎを求める場合、一例として、ニュートン法を適用し、下記の式（１４）、式（１５）及び式（１７）を交互に繰り返して計算することにより求めることができる。ただし、下記の式（１５）における「ｇ（β_ｎ）」は、下記の式（１６）を満たし、下記の式（１７）における「ｈ（α_ｎ）」は、下記の式（１８）を満たすこととする。

このように、更新部１８は、パラメータα_ｎ、β_ｎ及びγ_ｎが収束するように、式（１４）、式（１５）及び式（１７）の計算を繰り返すことにより、パラメータα_ｎ、β_ｎ及びγ_ｎをパラメータα_ｎ＋１、β_ｎ＋１及びγ_ｎ＋１へ更新する計算を繰り返す。

ここで、更新部１８は、反復回数が所定の回数になることを条件に、上記のパラメータの更新を終了することができる。例えば、脈波波形に含まれるノイズが小さい場合、脈波波形が周期変動ｓｉｎ波に元々近いはずなので、すぐに正しいα_ｎ、β_ｎ及びγ_ｎに収束する可能性が高い。この場合、反復回数は小さい数に固定しておけばよい。なお、パラメータの更新の終了条件は、これに限定されず、一致度ｆ_ｎが所定の閾値以上になるまでパラメータの更新を繰り返すこともできる。

図５は、パラメータの更新結果の一例を示す図である。図５の上部には、初期のモデル波形、すなわち反復回数ｎ＝０のモデル波形が示される一方で、図５の下部には、パラメータ更新後のモデル波形、すなわち反復回数ｎ＝上限Ｔｈのモデル波形が示されている。図５の上部に示す初期のパラメータα_０、β_０及びγ_０を式（１４）、式（１５）及び式（１７）にしたがって更新を繰り返すことにより、図５の下部に示すパラメータα_ｎ、β_ｎ及びγ_ｎを求める。すなわち、図５の例では、波形検出部１５が出力する脈波波形にフィッティングするように、モデル波形の周期を単調減少させるパラメータα_ｎ、β_ｎ及びγ_ｎが求められた更新結果が示されている。

脈波検出部１９は、更新部１８によりパラメータが更新されたモデル波形から脈波を検出する処理部である。

一実施形態として、脈波検出部１９は、更新部１８によりパラメータの更新が所定の反復回数にわたって行われた場合、当該パラメータα_ｎ、β_ｎ及びγ_ｎを上記の式（１０）へ代入することにより、脈波波形とモデル波形の一致度を算出する。その上で、脈波検出部１９は、一致度が所定の閾値以上であるか否かを判定する。このとき、一致度が閾値未満である場合、脈波検出部１９は、脈波、すなわち脈波波形や脈拍数を出力しない。一方、一致度が閾値以上である場合、脈波検出部１９は、波形検出部１５が出力する脈波波形または更新後のモデル波形を出力したり、更新後のモデル波形から脈拍数を求めて出力する。例えば、脈拍数は、更新後のモデル波形が含む極大点ごとに当該極大点の時刻ｔｉにおいてその時点での周期であるα_ｎｔ_ｉ＋β_ｎを用いて下記の式（１９）から計算できる。

（α_ｎｔ_ｉ＋β_ｎ）／２π×６０・・・（１９）
ただし、α_ｎｔ_ｉ ^２＋β_ｎｔ_ｉ＋γ_ｎが２πで除算された余りがπ／２であるｔ_ｉ

このようにして得られる脈拍数や脈波波形は、脈波検出装置１０が有する図示しない表示デバイスを始め、任意の出力先へ出力することができる。例えば、脈拍数や脈拍周期のゆらぎから自律神経の働きをＬＦ（Low Frequency）／ＨＦ（Hi Frequency）法を用いて診断したり、脈波波形から心疾患等を診断したりする診断プログラムが脈波検出装置１０にインストールされている場合には、診断プログラムを出力先とすることができる。また、診断プログラムをＷｅｂサービスとして提供するサーバ装置などを出力先とすることもできる。さらに、脈波検出装置１０を利用する利用者の関係者、例えば介護士や医者などが使用する端末装置を出力先とすることもできる。これによって、院外、例えば在宅や在席のモニタリングサービスも可能になる。なお、診断プログラムの測定結果や診断結果も、脈波検出装置１０を始め、関係者の端末装置に表示させることができるのも言うまでもない。

なお、上記の取得部１３、抽出部１４、波形検出部１５、転換判定部１６、初期設定部１７、更新部１８及び脈波検出部１９は、ＣＰＵ（Central Processing Unit）やＭＰＵ（Micro Processing Unit）などに脈波検出プログラムや信号処理プログラムを実行させることによって実現できる。また、上記の各機能部は、ＡＳＩＣ（Application Specific Integrated Circuit）やＦＰＧＡ（Field Programmable Gate Array）などのハードワイヤードロジックによっても実現できる。

また、上記の内部メモリには、半導体メモリ素子や記憶装置を採用できる。例えば、半導体メモリ素子の一例としては、ＶＲＡＭ（Video Random Access Memory）、ＲＡＭ（Random Access Memory)やフラッシュメモリ（flash memory）などが挙げられる。また、メモリではなく、ハードディスク、光ディスクなどの記憶装置を用いることとしてもかまわない。

［処理の流れ］
次に、本実施例に係る脈波検出装置の処理の流れについて説明する。なお、ここでは、脈波検出装置１０によって実行される（１）波形検出処理について説明した後に、（２）パラメータの更新処理について説明することとする。

（１）波形検出処理
図６は、実施例１に係る波形検出処理の手順を示すフローチャートである。この波形検出処理は、カメラ１１から画像が取得する度に処理を起動し、画像が取得されなくなるまで繰り返し実行される処理である。なお、図示しない入力デバイス等を介して中断操作を受け付けた場合には、波形検出処理を中止することもできる。

図６に示すように、被験者が映った画像が取得されると（ステップＳ１０１）、抽出部１４は、ステップＳ１０１で取得された画像から所定の顔パーツ、例えば被験者の鼻を基準とする生体領域の画像を抽出する（ステップＳ１０２）。

その上で、抽出部１４は、Ｒ信号の時系列データをＢＰＦ１５２Ｒ及びＢＰＦ１５６Ｒへ出力するとともに、Ｇ信号の時系列データをＢＰＦ１５２Ｇ及びＢＰＦ１５６Ｇへ出力する（ステップＳ１０３）。

続いて、ＢＰＦ１５２Ｒは、Ｒ信号の特定周波数帯、例えば３ｂｐｍ以上２０ｂｐｍ未満の周波数帯の信号成分を抽出するとともに、ＢＰＦ１５２Ｇは、Ｇ信号の特定周波数帯の信号成分を抽出する（ステップＳ１０４Ａ）。

そして、抽出部１５３Ｒは、Ｒ信号の特定周波数帯の信号成分の絶対強度値を抽出するとともに、抽出部１５３Ｇは、Ｇ信号の特定周波数帯の信号成分の絶対強度値を抽出する（ステップＳ１０５）。

その後、ＬＰＦ１５４Ｒは、Ｒ信号の特定周波数帯の絶対強度値の時系列データに対し、時間変化に応答させる平滑化処理を実行するとともに、ＬＰＦ１５４Ｇは、Ｇ信号の特定周波数帯の絶対強度値の時系列データに対し、時間変化に応答させる平滑化処理を実行する（ステップＳ１０６）。

続いて、算出部１５５は、ＬＰＦ１５４Ｇによって出力されたＧ信号の特定周波数帯の絶対値強度Ｇ´_noiseを、ＬＰＦ１５４Ｒによって出力されたＲ信号の特定周波数帯の絶対値強度Ｒ´_noiseで除する除算「Ｇ´_noise／Ｒ´_noise」を実行することによって補正係数ａを算出する（ステップＳ１０７）。

上記のステップＳ１０４Ａの処理に並行して、ＢＰＦ１５６Ｒは、Ｒ信号の脈波周波数帯、例えば４２ｂｐｍ以上２４０ｂｐｍ未満の周波数帯の信号成分を抽出するとともに、ＢＰＦ１５６Ｇは、Ｇ信号の脈波周波数帯の信号成分を抽出する（ステップＳ１０４Ｂ）。

その後、乗算部１５７は、ステップＳ１０４Ｂで抽出されたＲ信号の脈波周波数帯の信号成分にステップＳ１０７で算出された補正係数ａを乗算する（ステップＳ１０８）。その上で、演算部１５８は、ステップＳ１０８で補正係数ａが乗算されたＲ信号の脈波周波数帯の信号成分から、ステップＳ１０４Ｂで抽出されたＧ信号の脈波周波数帯の信号成分を差し引く演算「ａ＊Ｒ_signal−Ｇ_signal」を実行する（ステップＳ１０９）。

そして、波形検出部１５は、演算後の信号の時系列データを脈波波形として転換判定部１６へ出力し（ステップＳ１１０）、処理を終了する。

（２）パラメータの更新処理
図７は、実施例１に係るパラメータの更新処理の手順を示すフローチャートである。この処理は、顔画像が撮像されるフレーム周波数に対応するサンプリング間隔で波形検出部１５から脈波波形を形成する最新の振幅値が入力される度に実行される。

図７に示すように、転換判定部１６は、波形検出部１５によって出力される脈波波形の最新のデータ点に関する振幅値および時刻を取得する（ステップＳ３０１）。そして、転換判定部１６は、ステップＳ３０１で取得された最新のデータ点が前後で振幅値の正負が転換する転換点であるか否かを判定する（ステップＳ３０２）。

このとき、転換点である場合（ステップＳ３０２Ｙｅｓ）、初期設定部１７は、モデル波形の周期数ｍを設定すると共に、上記のパラメータα、β及びγの初期値を設定することにより、初期のモデル波形を定義する各パラメータを初期設定する（ステップＳ３０３）。

その上で、更新部１８は、一例として、ニュートン法により導出される、上記の式（１４）、式（１５）及び式（１７）にしたがってパラメータα、β及びγを更新する（ステップＳ３０４）。

そして、パラメータα、β及びγの更新が所定の反復回数にわたって行われるまで（ステップＳ３０５Ｎｏ）、更新部１８は、ステップＳ３０４に戻り、パラメータα、β及びγの更新を繰り返し実行する。

一方、パラメータα、β及びγの更新が所定の反復回数にわたって行われた場合（ステップＳ３０５Ｙｅｓ）、脈波検出部１９は、所定の反復回数にわたって更新されたパラメータα_ｎ、β_ｎ及びγ_ｎを上記の式（１０）へ代入することにより、脈波波形とモデル波形の一致度を算出する（ステップＳ３０６）。

ここで、脈波検出部１９は、一致度が所定の閾値以上であるか否かを判定する（ステップＳ３０７）。このとき、一致度が閾値未満である場合（ステップＳ３０７Ｎｏ）、脈波、すなわち脈波波形や脈拍数を出力しない。一方、一致度が閾値以上である場合（ステップＳ３０７Ｙｅｓ）、脈波検出部１９は、波形検出部１５が出力する脈波波形または更新後のモデル波形を出力したり、あるいは更新後のモデル波形から脈拍数を求めて出力し（ステップＳ３０８）、処理を終了する。

なお、図７に示すフローチャートでは、ステップＳ３０４でモデル波形のパラメータが更新される度に、窓にモデル波形がｍ周期分入るよう、窓幅を更新することもできる。これによって、例えば、初期のモデル波形に比べて更新後のモデル波形の周期が長い場合に、窓内に含まれる周期数が小さくなる結果、脈拍数を誤検出する可能性があるが、かかる誤検出を低減できる。

また、図７に示すフローチャートでは、最終的なモデル波形のパラメータが所与の範囲に含まれている場合のみ脈拍数を計算することとしてもかまわない。例えば、上記の式（１９）により得られる脈拍数が一般的な脈拍数範囲である４２ｂｐｍ〜２００ｂｐｍに含まれる場合に絞って脈拍数を計算したり、位相パラメータγ_ｎの絶対値が所与の閾値（例えば、π／４）より小さい場合に絞って脈拍数を計算する。これによって、ノイズを含む脈拍波形から誤った脈拍数を計算する可能性を低減させることもできる。

［効果の一側面］
上述してきたように、本実施例に係る脈波検出装置１０は、生体が撮影された画像から検出された脈波波形との間で一致度が最大となるように、モデル波形の周期変動に関するパラメータを更新し、その周期変動から脈波波形や脈拍数を求める。一側面として、一拍ごとに脈拍の周期が変動する場合でも、脈波波形にモデル波形をフィッティングさせることができる。したがって、本実施例に係る脈波検出装置１０によれば、脈波の周期変動に追従できる。他の側面として、モデル波形の周期を一定に固定する場合よりもグラフ形状の自由度を向上させることができる。それ故、脈拍数の算出精度も高めることができる。

さて、これまで開示の装置に関する実施例について説明したが、本発明は上述した実施例以外にも、種々の異なる形態にて実施されてよいものである。そこで、以下では、本発明に含まれる他の実施例を説明する。

［入力信号］
例えば、上記の実施例１では、入力信号としてＲ信号およびＧ信号の二種類を用いる場合を例示したが、異なる複数の光波長成分を持つ信号であれば任意の種類の信号および任意の数の信号を入力信号とすることができる。例えば、Ｒ、Ｇ、Ｂ、ＩＲおよびＮＩＲなどの光波長成分が異なる信号のうち任意の組合せの信号を２つ用いることもできるし、また、３つ以上用いることもできる。

［他の実装例］
上記の実施例１では、脈波検出装置１０が上記一連の処理をスタンドアローンで実行する場合を例示したが、クライアントサーバシステムとして実装することもできる。例えば、脈波検出装置１０は、図６及び図７に示す処理を実行するＷｅｂサーバとして実装することとしてもよいし、アウトソーシングによって脈波検出サービスを始めとするサービスを提供するクラウドとして実装することとしてもかまわない。このように、脈波検出装置１０がサーバ装置として動作する場合には、スマートフォンや携帯電話機等の携帯端末装置やパーソナルコンピュータ等の情報処理装置をクライアント端末として収容することができる。これらクライアント端末からネットワークを介して被験者の顔が映った画像が取得された場合に波形検出処理を実行し、その脈波波形を用いてパラメータの更新処理を実行した上でその検出結果やその検出結果を用いてなされた診断結果をクライアント端末へ応答することによって脈波検出サービス及び診断サービスを提供できる。

［分散および統合］
また、図示した各装置の各構成要素は、必ずしも物理的に図示の如く構成されておらずともよい。すなわち、各装置の分散・統合の具体的形態は図示のものに限られず、その全部または一部を、各種の負荷や使用状況などに応じて、任意の単位で機能的または物理的に分散・統合して構成することができる。例えば、取得部１３、抽出部１４、波形検出部１５、転換判定部１６、初期設定部１７、更新部１８または脈波検出部１９を脈波検出装置１０の外部装置としてネットワーク経由で接続するようにしてもよい。また、取得部１３、抽出部１４、波形検出部１５、転換判定部１６、初期設定部１７、更新部１８または脈波検出部１９を別の装置がそれぞれ有し、ネットワーク接続されて協働することで、上記の脈波検出装置１０の機能を実現するようにしてもよい。

［脈波検出プログラム］
また、上記の実施例で説明した各種の処理は、予め用意されたプログラムをパーソナルコンピュータやワークステーションなどのコンピュータで実行することによって実現することができる。そこで、以下では、図８を用いて、上記の実施例と同様の機能を有する脈波検出プログラムを実行するコンピュータの一例について説明する。

図８は、実施例１及び実施例２に係る脈波検出プログラムを実行するコンピュータのハードウェア構成例を示す図である。図８に示すように、コンピュータ１０００は、操作部１１００ａと、スピーカ１１００ｂと、カメラ１１００ｃと、ディスプレイ１２００と、通信部１３００とを有する。さらに、このコンピュータ１０００は、ＣＰＵ１５００と、ＲＯＭ１６００と、ＨＤＤ１７００と、ＲＡＭ１８００とを有する。これら１１００〜１８００の各部はバス１４００を介して接続される。

ＨＤＤ１７００には、図８に示すように、上記の実施例１で示した取得部１３、抽出部１４、波形検出部１５、転換判定部１６、初期設定部１７、更新部１８及び脈波検出部１９と同様の機能を発揮する脈波検出プログラム１７００ａが記憶される。この脈波検出プログラム１７００ａは、図１に示した取得部１３、抽出部１４、波形検出部１５、転換判定部１６、初期設定部１７、更新部１８及び脈波検出部１９の各構成要素と同様、統合又は分離してもかまわない。すなわち、ＨＤＤ１７００には、必ずしも上記の実施例１で示した全てのデータが格納されずともよく、処理に用いるデータがＨＤＤ１７００に格納されればよい。

このような環境の下、ＣＰＵ１５００は、ＨＤＤ１７００から脈波検出プログラム１７００ａを読み出した上でＲＡＭ１８００へ展開する。この結果、脈波検出プログラム１７００ａは、図８に示すように、脈波検出プロセス１８００ａとして機能する。この脈波検出プロセス１８００ａは、ＲＡＭ１８００が有する記憶領域のうち脈波検出プロセス１８００ａに割り当てられた領域にＨＤＤ１７００から読み出した各種データを展開し、この展開した各種データを用いて各種の処理を実行する。例えば、脈波検出プロセス１８００ａが実行する処理の一例として、図６や図７に示す処理などが含まれる。なお、ＣＰＵ１５００では、必ずしも上記の実施例１で示した全ての処理部が動作せずともよく、実行対象とする処理に対応する処理部が仮想的に実現されればよい。

なお、上記の脈波検出プログラム１７００ａは、必ずしも最初からＨＤＤ１７００やＲＯＭ１６００に記憶されておらずともかまわない。例えば、コンピュータ１０００に挿入されるフレキシブルディスク、いわゆるＦＤ、ＣＤ−ＲＯＭ、ＤＶＤディスク、光磁気ディスク、ＩＣカードなどの「可搬用の物理媒体」に各プログラムを記憶させる。そして、コンピュータ１０００がこれらの可搬用の物理媒体から各プログラムを取得して実行するようにしてもよい。また、公衆回線、インターネット、ＬＡＮ、ＷＡＮなどを介してコンピュータ１０００に接続される他のコンピュータまたはサーバ装置などに各プログラムを記憶させておき、コンピュータ１０００がこれらから各プログラムを取得して実行するようにしてもよい。

１０ 脈波検出装置
１１ カメラ
１３ 取得部
１４ 抽出部
１５ 波形検出部
１６ 転換判定部
１７ 初期設定部
１８ 更新部
１９ 脈波検出部

Claims (7)

1. 画像を取得する取得部と、
   前記画像から脈波波形を検出する波形検出部と、
   前記脈波波形と、所定の周期的なモデル波形であって周期変動をパラメータに含むモデル波形との一致度に基づいて、前記パラメータの更新を前記パラメータが収束するまで反復する更新部と、
   反復後のパラメータを有するモデル波形と前記脈波波形との一致度が所定の閾値以上である場合に、前記脈波波形を出力する出力部と
   を有することを特徴とする脈波検出装置。
2. 前記反復後のモデル波形のパラメータを用いて脈拍数を算出する算出部をさらに有することを特徴とする請求項１に記載の脈波検出装置。
3. 前記モデル波形として、周期が時間に対して単調に変動する周期変動正弦波を用いることを特徴とする請求項１または２に記載の脈波検出装置。
4. 前記モデル波形として、周期が時間に対して線形に変動する周期変動正弦波を用いることを特徴とする請求項３に記載の脈波検出装置。
5. 前記モデル波形は、前記周期変動と共に周期および位相をパラメータに含むことを特徴とする請求項１〜４のいずれか１つに記載の脈波検出装置。
6. コンピュータが、
   画像を取得し、
   前記画像から脈波波形を検出し、
   前記脈波波形と、所定の周期的なモデル波形であって周期変動をパラメータに含むモデル波形との一致度に基づいて、前記パラメータの更新を前記パラメータが収束するまで反復し、
   反復後のパラメータを有するモデル波形と前記脈波波形との一致度が所定の閾値以上である場合に、前記脈波波形を出力する、
   処理を実行することを特徴とする脈波検出方法。
7. コンピュータに、
   画像を取得し、
   前記画像から脈波波形を検出し、
   前記脈波波形と、所定の周期的なモデル波形であって周期変動をパラメータに含むモデル波形との一致度に基づいて、前記パラメータの更新を前記パラメータが収束するまで反復し、
   反復後のパラメータを有するモデル波形と前記脈波波形との一致度が所定の閾値以上である場合に、前記脈波波形を出力する、
   処理を実行させることを特徴とする脈波検出プログラム。

Images