JP6098257B2 - 信号処理装置、信号処理方法及び信号処理プログラム - Google Patents

信号処理装置、信号処理方法及び信号処理プログラム Download PDF

Info

Publication number
JP6098257B2
JP6098257B2 JP2013054165A JP2013054165A JP6098257B2 JP 6098257 B2 JP6098257 B2 JP 6098257B2 JP 2013054165 A JP2013054165 A JP 2013054165A JP 2013054165 A JP2013054165 A JP 2013054165A JP 6098257 B2 JP6098257 B2 JP 6098257B2
Authority
JP
Japan
Prior art keywords
window
signal
waveform
pulse wave
correlation coefficient
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2013054165A
Other languages
English (en)
Other versions
JP2014176584A (ja
Inventor
一穂 前田
一穂 前田
雅人 阪田
雅人 阪田
大輔 内田
大輔 内田
明大 猪又
明大 猪又
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Fujitsu Ltd
Original Assignee
Fujitsu Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Fujitsu Ltd filed Critical Fujitsu Ltd
Priority to JP2013054165A priority Critical patent/JP6098257B2/ja
Priority to US14/099,457 priority patent/US9962126B2/en
Priority to EP13197941.1A priority patent/EP2777485B1/en
Publication of JP2014176584A publication Critical patent/JP2014176584A/ja
Application granted granted Critical
Publication of JP6098257B2 publication Critical patent/JP6098257B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/72Signal processing specially adapted for physiological signals or for diagnostic purposes
    • A61B5/7235Details of waveform analysis
    • A61B5/7246Details of waveform analysis using correlation, e.g. template matching or determination of similarity
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/0059Measuring for diagnostic purposes; Identification of persons using light, e.g. diagnosis by transillumination, diascopy, fluorescence
    • A61B5/0077Devices for viewing the surface of the body, e.g. camera, magnifying lens
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/0059Measuring for diagnostic purposes; Identification of persons using light, e.g. diagnosis by transillumination, diascopy, fluorescence
    • A61B5/0082Measuring for diagnostic purposes; Identification of persons using light, e.g. diagnosis by transillumination, diascopy, fluorescence adapted for particular medical purposes
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00Measuring for diagnostic purposes; Identification of persons
    • A61B5/02Detecting, measuring or recording pulse, heart rate, blood pressure or blood flow; Combined pulse/heart-rate/blood pressure determination; Evaluating a cardiovascular condition not otherwise provided for, e.g. using combinations of techniques provided for in this group with electrocardiography or electroauscultation; Heart catheters for measuring blood pressure
    • A61B5/024Detecting, measuring or recording pulse rate or heart rate

Landscapes

  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Surgery (AREA)
  • Animal Behavior & Ethology (AREA)
  • Veterinary Medicine (AREA)
  • Public Health (AREA)
  • Physics & Mathematics (AREA)
  • General Health & Medical Sciences (AREA)
  • Biophysics (AREA)
  • Pathology (AREA)
  • Biomedical Technology (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Medical Informatics (AREA)
  • Molecular Biology (AREA)
  • Physiology (AREA)
  • Cardiology (AREA)
  • Artificial Intelligence (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Signal Processing (AREA)
  • Psychiatry (AREA)
  • Measuring Pulse, Heart Rate, Blood Pressure Or Blood Flow (AREA)

Description

本発明は、信号処理装置、信号処理方法及び信号処理プログラムに関する。
信号の波形が類似するか否かを判断する場合に、相関係数が用いられることがある。例えば、相関係数を算出する技術の一例として、次のような心拍計測装置が挙げられる。かかる心拍計測装置は、車両のハンドルに操作者の左右の手に対応して設けた第1の電極及び第2の電極と、リファレンス電位を測定する基準電極とを用いて測定される心電信号から操作者の心拍の間隔を計測するものである。
ここで、上記の心拍計測装置では、心電信号から検出される振幅のピークの時刻を補正する場合に、ハンドルのグリップパターンが相関係数を用いて判別される。すなわち、心拍計測装置は、グリップパターンごとに登録された心電信号の波形と心電信号の波形との間で相関係数を算出することによって類似度を求める。これによって、類似度が最大である心電信号の波形のグリップパターンに対応付けられた補正量にしたがってピークの時刻を補正する。
特開2012−239474号公報
しかしながら、上記の技術では、互いに相似する波形が類似と判定されるので、波形の類似性の判定精度が低下する場合がある。
1つの側面では、波形の類似性の判定精度を向上させることができる信号処理装置、信号処理方法及び信号処理プログラムを提供することを目的とする。
一態様の信号処理装置は、所定の時間幅を持つ第1の窓に含まれる信号の波形と、前記第1の窓内の波形と比較する第2の窓に含まれる信号の波形との間で、前記第1の窓内の振幅値の分散と前記第2の窓内の振幅値の分散との差が大きくなるにしたがって値が小さくなる相関係数を算出する相関係数算出部を有する。
一実施形態によれば、波形の類似性の判定精度を向上させることができる。
図1は、実施例1に係る脈波検出装置の機能的構成を示すブロック図である。 図2は、顔画像の一例を示す図である。 図3は、G信号およびR信号の各信号のスペクトルの一例を示す図である。 図4は、G成分および補正係数kが乗算されたR成分の各信号のスペクトルの一例を示す図である。 図5は、演算後のスペクトルの一例を示す図である。 図6は、波形検出部の機能的構成を示すブロック図である。 図7は、窓幅の設定方法の一例を示す図である。 図8は、脈波波形に含まれる極大点の一例を示す図である。 図9は、ずらし範囲の設定例を示す図である。 図10は、各窓内の波形の一例を示す図である。 図11は、各窓内の波形の一例を示す図である。 図12は、実施例1に係る波形検出処理の手順を示すフローチャートである。 図13は、実施例1に係る相関係数算出処理の手順を示すフローチャートである。 図14は、実施例1及び実施例2に係る信号処理プログラムを実行するコンピュータの一例について説明するための図である。
以下に添付図面を参照して本願に係る信号処理装置、信号処理方法及び信号処理プログラムについて説明する。なお、この実施例は開示の技術を限定するものではない。そして、各実施例は、処理内容を矛盾させない範囲で適宜組み合わせることが可能である。
[脈波検出装置の構成]
まず、本実施例に係る脈波検出装置の機能的構成について説明する。図1は、実施例1に係る脈波検出装置の機能的構成を示すブロック図である。図1に示す脈波検出装置10は、太陽光や室内光などの一般の環境光の下で生体に計測器具を接触させずに、被験者が撮影された画像を用いて被験者の脈波、すなわち心臓の拍動に伴う血液の体積の変動を検出する脈波検出処理を実行するものである。
かかる脈波検出装置10は、上記の脈波検出の一環として、生体が撮影された画像から検出された脈波に対応する周波数帯の信号を脈波検出に用いるか否かを過去の波形と類似性を有するか否か、すなわち波形が周期性を有するか否かによって判断する。かかる類似性を判断するために、脈波検出装置10は、所定の時間幅を持つ基準窓に含まれる信号の波形と、基準窓内の波形と比較するために基準窓よりも時間を遡って設定された比較窓に含まれる信号の波形との間で相関係数を算出する。
このように、本実施例に係る脈波検出装置10は、波形の類似を相関係数を用いて判断する場合に、基準窓内の振幅値の分散と比較窓内の振幅値の分散との差が大きくなるにしたがって値が小さくなる相関係数を算出する信号処理を実行する。このため、本実施例に係る脈波検出装置10では、互いの窓内の波形が相似の関係にある場合でも相関係数が低く算出される。したがって、本実施例に係る脈波検出装置10によれば、波形の類似性の判定精度を向上させることができる。
かかる脈波検出装置10は、一態様として、上記の信号処理がパッケージソフトウェアやオンラインソフトウェアとして提供される信号処理プログラムを所望のコンピュータにインストールさせることによって実装できる。例えば、スマートフォン、携帯電話機やPHS(Personal Handyphone System)などの移動体通信端末のみならず、移動体通信網に接続する能力を持たないデジタルカメラ、タブレット端末やスレート端末を含む携帯端末装置に上記の信号処理プログラムをインストールさせる。これによって、携帯端末装置を脈波検出装置10として機能させることができる。なお、ここでは、脈波検出装置10の実装例として携帯端末装置を例示したが、パーソナルコンピュータを始めとする据置き型の端末装置に信号処理プログラムをインストールさせることもできる。
図1に示すように、脈波検出装置10は、カメラ11と、取得部13と、抽出部14と、波形検出部15と、転換判定部16aと、極値判定部16bと、窓幅設定部17aと、ずらし範囲設定部17bと、相関係数算出部18と、脈波検出部19とを有する。
かかる脈波検出装置10は、図1に示した機能部以外にも既知のコンピュータが有する各種の機能部を有することとしてもかまわない。例えば、脈波検出装置10が据置き端末として実装される場合には、キーボード、マウスやディスプレイなどの入出力デバイスをさらに有することとしてもよい。また、脈波検出装置10がタブレット端末やスレート端末として実装される場合には、タッチパネルをさらに有することとしてもよい。また、脈波検出装置10が移動体通信端末として実装される場合には、アンテナ、移動体通信網に接続する無線通信部、GPS(Global Positioning System)受信機などの機能部をさらに有していてもかまわない。
カメラ11は、CCD(Charge Coupled Device)やCMOS(Complementary Metal Oxide Semiconductor)などの撮像素子を搭載する撮像装置である。例えば、カメラ11には、R(red)、G(green)、B(blue)など3種以上の受光素子を搭載することができる。かかるカメラ11の実装例としては、デジタルカメラやWebカメラを外部端子を介して接続することとしてもよいし、カメラが出荷時から搭載されている場合にはそのカメラを流用できる。なお、ここでは、脈波検出装置10がカメラ11を有する場合を例示したが、ネットワークまたは記憶デバイスを経由して画像を取得できる場合には、必ずしも脈波検出装置10がカメラ11を有さずともよい。
ここで、脈波検出装置10は、上記の脈波検出処理を実行するアプリケーションプログラムがプリインストールまたはインストールされている場合に、カメラ11によって脈波を検出し易い被験者の画像が撮像されるように画像の撮影操作を案内することができる。なお、以下では、上記のアプリケーションプログラムのことを「脈波検出用アプリ」と記載する場合がある。
かかる脈波検出用アプリは、図示しない入力デバイスを介して起動されると、カメラ11を起動する。これを受けて、カメラ11は、カメラ11の撮影範囲に収容された被写体の撮影を開始する。このとき、被験者の顔が映る画像を撮影させる場合には、脈波検出用アプリは、カメラ11が撮影する画像を図示しない表示デバイスに表示しつつ、被験者の鼻を映す目標位置を照準として表示させることもできる。これによって、被験者の眼、耳、鼻や口などの顔パーツの中でも被験者の鼻が撮影範囲の中心部分に収まった画像が撮影できるようにする。そして、脈波検出用アプリは、カメラ11によって被験者の顔が撮影された画像を取得部13へ保存する。なお、以下では、顔が映った画像のことを「顔画像」と記載する場合がある。
取得部13は、画像を取得する処理部である。一態様としては、取得部13は、カメラ11によって撮像された顔画像を取得する。他の一態様としては、取得部13は、顔画像を蓄積するハードディスクや光ディスクなどの補助記憶装置またはメモリカードやUSB(Universal Serial Bus)メモリなどのリムーバブルメディアから画像を取得することもできる。更なる一態様としては、取得部13は、外部装置からネットワークを介して受信することによって顔画像を取得することもできる。なお、取得部13は、CCDやCMOSなどの撮像素子による出力から得られる2次元のビットマップデータやベクタデータなどの画像データを用いて処理を実行する場合を例示したが、1つのディテクタから出力される信号をそのまま取得して後段の処理を実行させることとしてもよい。
抽出部14は、取得部13によって取得された画像から生体領域を抽出する処理部である。一態様としては、顔画像から所定の顔パーツを基準とする生体領域を抽出する。例えば、抽出部14は、顔画像にテンプレートマッチング等の画像処理を実行することによって被験者の眼、耳、鼻や口などの顔パーツのうち特定の顔パーツ、すなわち被験者の鼻を検出する。その上で、抽出部14は、被験者の鼻を中心とし、中心から所定の範囲に含まれる生体領域を抽出する。これによって、被験者の鼻、鼻の周辺に位置する頬の一部の顔中心部分を含んだ生体領域の画像が脈波の検出に使用する画像として抽出される。その後、抽出部14は、原画像から抽出した生体領域の画像を波形検出部15へ出力する。
図2は、顔画像の一例を示す図である。図2には、画像に映る被験者の眼、鼻及び口の一部または全部を含む領域が9つに分割されたブロックが図示されている。図2に示すブロックのうち上段の左及び右のブロックには、被験者の眼が映っている。これらのブロックの画像を検出に用いた場合には、眼の瞬きがノイズとなって心拍数の検出精度の低下を招く場合がある。また、図2に示すブロックのうち下段の3つのブロックには、被験者の口が映っている。これらのブロックの画像を検出に用いた場合には、口の動きがノイズとなって心拍数の検出精度の低下を招く場合がある。一方、図2に示す中段の真ん中のブロック、すなわち斜線の塗りつぶしが図示されたブロックは、眼や口が映るブロックから隔てられており、他のブロックに比べてノイズとなる成分が映っている可能性が低いので、良好な検出結果を期待できる。これらのことから、抽出部14は、原画像から図2に示す中段の真ん中のブロックの画像を生体領域の画像として抽出する。
上記の生体領域の画像を抽出後に、抽出部14は、生体領域に含まれる各画素が持つ画素値に所定の統計処理を実行する。例えば、抽出部14は、生体領域に含まれる各画素が持つ画素値を波長成分ごとに平均する。この他、平均値以外にも、中央値や最頻値を計算することとしてもよく、また、加重平均以外にも任意の平均処理、例えば加重平均や移動平均などを実行することもできる。これによって、生体領域に含まれる各画素が持つ画素値の平均値が当該生体領域を代表する代表値として波長成分ごとに算出される。
波形検出部15は、脈波の検出対象とする領域に含まれる各画素の波長成分別の代表値の信号から、各波長成分の間で脈波が採り得る脈波周波数帯以外の特定周波数帯の成分が互いに相殺された信号の波形を検出する処理部である。
一態様としては、波形検出部15は、画像に含まれる3つの波長成分、すなわちR成分、G成分およびB成分のうち血液の吸光特定が異なるR成分とG成分の2つの波長成分の代表値の時系列データを用いて、波形を検出する。
これを説明すると、顔表面には、毛細血管が流れており、心拍により血管に流れる血流が変化すると、血流で吸収される光量も心拍に応じて変化するため、顔からの反射によって得られる輝度も心拍に伴って変化する。かかる輝度の変化量は小さいが、顔領域全体の平均輝度を求めると、輝度の時系列データには脈波成分が含まれる。ところが、輝度は、脈波以外に体動等によっても変化し、これが、脈波検出のノイズ成分、いわゆる体動アーチファクトとなる。そこで、血液の吸光特性の異なる2種類以上の波長、例えば吸光特性が高いG成分(525nm程度)、吸光特性が低いR成分(700nm程度)で脈波を検出する。心拍は、0.5Hz〜4Hz、1分あたりに換算すれば30bpm〜240bpmの範囲であるので、それ以外の成分はノイズ成分とみなすことができる。ノイズには、波長特性は無い、あるいはあっても極小であると仮定すると、G信号およびR信号の間で0.5Hz〜4Hz以外の成分は等しいはずであるが、カメラの感度差により大きさが異なる。それゆえ、0.5Hz〜4Hz以外の成分の感度差を補正して、G成分からR成分を減算すれば、ノイズ成分は除去されて脈波成分のみを取り出すことができる。
例えば、G成分及びR成分は、下記の式(1)および下記の式(2)によって表すことができる。下記の式(1)における「Gs」は、G信号の脈波成分を指し、「Gn」は、G信号のノイズ成分を指し、また、下記の式(2)における「Rs」は、R信号の脈波成分を指し、「Rn」は、R信号のノイズ成分を指す。また、ノイズ成分は、G成分およびR成分の間で感度差があるので、感度差の補正係数kは、下記の式(3)によって表される。
Ga=Gs+Gn・・・(1)
Ra=Rs+Rn・・・(2)
k=Gn/Rn・・・(3)
感度差を補正してG成分からR成分を減算すると、脈波成分Sは、下記の式(4)となる。これを上記の式(1)及び上記の式(2)を用いて、Gs、Gn、Rs及びRnによって表される式へ変形すると、下記の式(5)となり、さらに、上記の式(3)を用いて、kを消し、式を整理すると下記の式(6)が導出される。
S=Ga−kRa・・・(4)
S=Gs+Gn−k(Rs+Rn)・・・(5)
S=Gs−(Gn/Rn)Rs・・・(6)
ここで、G信号およびR信号は、吸光特性が異なり、Gs>(Gn/Rn)Rsである。したがって、上記の式(6)によってノイズが除去された脈波成分Sを算出することができる。
図3は、G信号およびR信号の各信号のスペクトルの一例を示す図である。図3に示すグラフの縦軸は、信号強度を指し、また、横軸は、周波数(bpm)を指す。図3に示すように、G成分およびR成分は、撮像素子の感度が異なるので、両者の信号強度はそれぞれ異なる。その一方、R成分およびG成分は、いずれにおいても30bpm〜240bpmの範囲外、特に3bpm以上20bpm未満の特定周波数帯でノイズが現れることには変わりはない。このため、図3に示すように、3bpm以上20bpm未満の特定周波数帯に含まれる指定の周波数Fnに対応する信号強度をGn及びRnとして抽出できる。これらGn及びRnによって感度差の補正係数kを導出できる。
図4は、G成分および補正係数kが乗算されたR成分の各信号のスペクトルの一例を示す図である。図4の例では、説明の便宜上、補正係数の絶対値を乗算した結果が図示されている。図4に示すグラフにおいても、縦軸は、信号強度を指し、また、横軸は、周波数(bpm)を指す。図4に示すように、G成分及びR成分の各信号のスペクトルに補正係数kが乗算された場合には、G成分およびR成分の各成分の間で感度が揃う。特に、特定周波数帯におけるスペクトルの信号強度は、大部分においてスペクトルの信号強度が略同一になっている。その一方で、実際に脈波が含まれる周波数の周辺領域400は、G成分およびR成分の各成分の間でスペクトルの信号強度が揃っていない。
図5は、演算後のスペクトルの一例を示す図である。図5では、脈波が現れている周波数帯の視認性を上げる観点から縦軸である信号強度の尺度を大きくして図示している。図5に示すように、G信号のスペクトルから補正係数kの乗算後のR信号のスペクトルが差し引かれた場合には、G成分およびR成分の間での吸光特性の差によって脈波が現れる信号成分の強度が可及的に維持された状態でノイズ成分が低減されていることがわかる。このようにしてノイズ成分だけが除去された脈波波形を検出することができる。
続いて、波形検出部15の機能的構成についてさらに具体的に説明する。図6は、波形検出部15の機能的構成を示すブロック図である。図6に示すように、波形検出部15は、BPF(Band-Pass Filter)152R及び152Gと、抽出部153R及び153Gと、LPF(Low-Pass Filter)154R及び154Gと、算出部155と、BPF156R及び156Gと、乗算部157と、演算部158とを有する。なお、図3〜図5の例では、周波数空間にて脈波を検出する例を説明したが、図6では、周波数成分への変換にかかる時間を削減する観点から、時系列空間にてノイズ成分をキャンセルして脈波を検出する場合の機能的構成を図示している。
例えば、抽出部14から波形検出部15には、生体領域に含まれる各画素が持つR成分の画素値の代表値を信号値とするR信号の時系列データが入力されるとともに、生体領域に含まれる各画素が持つG成分画素値の代表値を信号値とするG信号の時系列データが入力される。このうち、生体領域のR信号は、波形検出部15内のBPF152R及びBPF156Rへ入力されるとともに、生体領域のG信号は、波形検出部15内のBPF152G及びBPF156Gへ入力される。
BPF152R、BPF152G、BPF156R及びBPF156Gは、いずれも所定の周波数帯の信号成分だけを通過させてそれ以外の周波数帯の信号成分を除去するバンドパスフィルタである。これらBPF152R、BPF152G、BPF156R及びBPF156Gは、ハードウェアによって実装されることとしてもよいし、ソフトウェアによって実装されることとしてもよい。
これらBPFが通過させる周波数帯の違いについて説明する。BPF152R及びBPF152Gは、ノイズ成分が他の周波数帯よりも顕著に現れる特定周波数帯の信号成分を通過させる。
かかる特定周波数帯は、脈波が採り得る周波数帯との間で比較することによって定めることができる。脈波が採り得る周波数帯の一例としては、0.5Hz以上4Hz以下である周波数帯、1分あたりに換算すれば30bpm以上240bpm以下である周波数帯が挙げられる。このことから、特定周波数帯の一例としては、脈波として計測され得ない0.5Hz未満及び4Hz超過の周波数帯を採用することができる。また、特定周波数帯は、脈波が採り得る周波数帯との間でその一部が重複することとしてもよい。例えば、脈波として計測されることが想定しづらい0.7Hz〜1Hzの区間で脈波が採り得る周波数帯と重複することを許容し、1Hz未満及び4Hz以上の周波数帯を特定周波数帯として採用することもできる。また、特定周波数帯は、1Hz未満及び4Hz以上の周波数帯を外縁とし、ノイズがより顕著に現れる周波数帯に絞ることもできる。例えば、ノイズは、脈波が採り得る周波数帯よりも高い高周波数帯よりも、脈波が採り得る周波数帯よりも低い低周波数帯でより顕著に現れる。このため、1Hz未満の周波数帯に特定周波数帯を絞ることもできる。また、空間周波数がゼロである直流成分の近傍には、各成分の撮像素子の感度の差が多く含まれるので、3bpm以上60bpm未満の周波数帯に特定周波数帯を絞ることもできる。さらに、人の体の動き、例えば瞬きや体の揺れの他、環境光のチラツキなどのノイズが現れやすい3bpm以上20bpm未満の周波数帯に特定周波数帯を絞ることもできる。
ここでは、一例として、BPF152R及びBPF152Gが特定周波数帯として0.05Hz以上0.3Hz以下の周波数帯の信号成分を通過させる場合を想定して以下の説明を行う。なお、ここでは、特定周波数帯の信号成分を抽出するために、バンドパスフィルタを用いる場合を例示したが、一定の周波数未満の周波数帯の信号成分を抽出する場合などには、ローパスフィルタを用いることもできる。
一方、BPF156R及びBPF156Gは、脈波が採り得る周波数帯、例えば1Hz以上4Hz以下の周波数帯の信号成分を通過させる。なお、以下では、脈波が採り得る周波数帯のことを「脈波周波数帯」と記載する場合がある。
抽出部153Rは、R信号の特定周波数帯の信号成分の絶対強度値を抽出する。例えば、抽出部153Rは、R成分の特定周波数帯の信号成分をべき乗する乗算処理を実行することによって特定周波数帯の信号成分の絶対強度値を抽出する。また、抽出部153Gは、G信号の特定周波数帯の信号成分の絶対強度値を抽出する。例えば、抽出部153Gは、G成分の特定周波数帯の信号成分をべき乗する乗算処理を実行することによって特定周波数帯の信号成分の絶対強度値を抽出する。
LPF154R及びLPF154Gは、特定周波数帯の絶対強度値の時系列データに対し、時間変化に応答させる平滑化処理を実行するローパスフィルタである。これらLPF154R及びLPF154Gは、LPF154Rへ入力される信号がR信号であり、LPF154Gへ入力される信号がG信号である以外に違いはない。かかる平滑化処理によって、特定周波数帯の絶対値強度R´n及びG´nが得られる。
算出部155は、LPF154Gによって出力されたG信号の特定周波数帯の絶対値強度G´nを、LPF154Rによって出力されたR信号の特定周波数帯の絶対値強度R´nで除する除算「G´n/R´n」を実行する。これによって、感度差の補正係数kを算出する。
乗算部157は、BPF156Rによって出力されたR信号の脈波周波数帯の信号成分に算出部155によって算出された補正係数kを乗算する。
演算部158は、乗算部157によって補正係数kが乗算されたR信号の脈波周波数帯の信号成分から、BPF156Gによって出力されたG信号の脈波周波数帯の信号成分を差し引く演算「k*Rs−Gs」を実行する。かかる演算によって得られた信号の時系列データは、脈波の波形に相当する。
このようにして画像の生体領域から検出された脈波波形が転換判定部16aおよび極値判定部16bへ出力される。以下では、脈波波形のi番目の時刻をt(i)、脈波波形の振幅値をv(i)と表す場合がある。さらに、最新の脈波波形のインデックスを「I」と表す場合がある。なお、脈波波形の振幅値が採取される時刻t(i)の間隔、すなわちフレームレートは、δで一定とする。
図1の説明に戻り、転換判定部16aは、波形検出部15によって検出される信号の振幅値の正負が転換するか否かを判定する処理部である。一態様としては、転換判定部16aは、波形検出部15によって脈波波形の振幅値が出力される度に、当該脈波波形の振幅値の正負が転換するか否かを判定する。このとき、転換判定部16aは、下記の式(7)または下記の式(8)のいずれかを満たす場合に、当該振幅値v(i)の正負が転換すると判定する。このように、転換判定部16aは、振幅値の正負が転換する転換点の時刻T(j)=t(I)を図示しない内部メモリ等へ順次記録する。なお、上記の「j」は、何番目の転換点であるかを表し、このうち、最新の転換点のインデックスを「J」と表すこととする。
v(i-1)<0 and v(i)>=0・・・(7)
v(i-1)>=0 and v(i)<0・・・(8)
極値判定部16bは、波形検出部15によって検出される信号の振幅値が極値であるか否かを判定する処理部である。一態様としては、極値判定部16bは、波形検出部15によって脈波波形の振幅値が出力される度に、当該脈波波形の振幅値が出力されるよりも1つ前のインデックスを持つ脈拍波形の振幅値が極大値であるか否かを判定する。このとき、極値判定部16bは、下記の式(9)を満たす場合に、当該振幅値v(i−1)が極大値であると判定する。このように、極値判定部16bは、振幅値が極大値を採る極大点の時刻M(k)=t(I−1)を内部メモリへ順次記録する。なお、上記の「k」は、何番目の極大点であるかを表し、このうち、最新の極大点のインデックスを「K」と表すこととする。
t(i-2)<=t(i-1) and t(i-1)>t(i)・・・(9)
窓幅設定部17aは、基準窓の時間幅を設定する処理部である。かかる基準窓には、2つの窓のうち窓内の波形を脈波の検出に用いるかどうかが判断される方の窓を指し、例えば、最新の時刻t(I)から所定の時間幅を遡った区間が設定される。以下では、窓の時間幅のことを「窓幅」と記載する場合がある。かかる窓幅は、脈波波形の振幅値が取得される度に設定せずともよく、所定の基準窓の決定条件が満たされた場合に絞って窓幅を設定できる。例えば、基準窓の設定条件の一例としては、前回に基準窓の設定から所定のパラメータ、例えば1秒経過したことを契機に、窓幅の設定を開始することができる。他の一例としては、最新のインデックスIが転換点である場合に、窓幅の設定を開始することもできる。
一態様としては、窓幅設定部17aは、所定のパラメータ、例えば3秒間を窓幅Uに設定することができる。なお、窓幅設定部17aは、上記の3秒間以外にも、脈拍数が最も低い人、例えば42bpm以上の時間幅であれば、任意の窓幅Uを設定することができる。
他の一態様としては、窓幅設定部17aは、転換判定部16aによって振幅値の正負が転換すると判定された転換点の数が脈波の周期の整数倍と対応する個数になるまで遡った区間の長さを基準窓の時間幅として設定する。例えば、窓幅設定部17aは、脈波の周期に対応する転換点の個数nに対し、最新の転換点の時刻T(J)と脈波の2周期前の転換点の時刻T(J−2*n)との差を下記の式(10)によって求め、窓幅Uを決定する。
U=T(J)-T(J-2*n)・・・(10)
かかる転換点を用いた窓幅の設定では、脈波波形が脈波周波数帯、例えば42bpm〜240bpmの信号成分を持ち、脈波の測定条件が安定している状況では、脈波と同期して0を中心に正負を繰り返す周期性を利用している。これによって、脈波の測定条件が安定している状況下では、およそ脈波の1周期に対応する転換点の個数nの整数倍、本例では2倍の転換点を基準窓に含めることができる。
図7は、窓幅の設定方法の一例を示す図である。図7に示す縦軸は、振幅であり、横軸は、時間を表す。図7に示すように、脈波の1周期分を遡るには、最新の転換点を除き転換点を2つ遡ることによって脈波の1周期分を遡ることができるので、脈波の2周期分を遡るには4つ(=2*n)遡ればよいことがわかる。このように、上記の式(10)を用いて窓幅Uを設定する場合には、実測の脈波の周期と近似する周期の整数倍を窓幅に設定することができるので、最低の脈拍数、例えば42bpmよりも時間の長さが短い窓幅を設定できる結果、後述の相関係数算出部18における計算量を削減できる。
図1の説明に戻り、ずらし範囲設定部17bは、比較窓をずらす範囲を設定する処理部である。以下では、比較窓をずらす範囲のことを「ずらし範囲」と記載する場合がある。ここで、ずらし範囲設定部17bは、窓幅設定部17aによって基準窓の窓幅が設定されたからといって直ちにずらし範囲を設定する訳ではない。すなわち、ずらし範囲設定部17bは、極値判定部16bによって極大であると判定された極大点のうち窓幅設定部17aによって時間幅が設定された基準窓に含まれる極大点の数が脈波の周期と対応する個数と等しい場合に、ずらし範囲の設定を実行する。
これを具体的に説明すると、まず、ずらし範囲設定部17bは、基準窓に含まれる極大点の数を計数する。例えば、下記の式(11)の条件を満たす整数kの数を計数する。そして、ずらし範囲設定部17bは、基準窓内に含まれる極大点の数M(k)が脈波の1周期に対応する極大点の個数nを超えるか否かを判定する。このとき、基準窓内に脈波の1周期に対応する極大点の個数nよりも多くの極大点が存在する場合には、測定条件が不安定であるとみなすことができる。この場合には、ずらし範囲設定部17bは、ずらし範囲の設定は実行しない。このように、ずらし範囲の設定の開始条件に極大点の数を用いるのは、理想的な脈波波形は1周期あたり一つの極大点を持つという性質を利用し、それよりも多い極大点は測定条件の変化の影響とみなすためである。
T(J)-U<=M(k)<=T(J)・・・(11)
図8は、脈波波形に含まれる極大点の一例を示す図である。図8に示す縦軸は、振幅で指し、横軸は、時間を指す。図8に示すように、脈波の2周期分の窓幅を持つ基準窓に極大点が3つ含まれている。ところが、ノイズを無視できる測定状況下では、脈波の1周期につき極大点は1つしか現れない。よって、図8に示す基準窓内の信号にはピーク検出に支障をきたす程度にノイズが重畳しているとみなし、以降の演算を中止することができる。また、ノイズが重畳した脈波波形によるずらし範囲の設定および相関係数の算出を中止することによって演算量を低減することもできる。
このように、ずらし範囲設定部17bは、基準窓内に含まれる極大点の数M(k)が脈波の1周期に対応する極大点の個数n以下である場合に、ずらし範囲の設定を開始する。
一態様としては、ずらし範囲設定部17bは、一般に知られている人の脈拍数がとりうる範囲、例えば42bpm〜240bpmを基準に設定することができる。例えば、ずらし範囲設定部17bは、ずらし範囲w(秒)を下記の式(12)が表す範囲に設定できる。このとき、下記の式(12)におけるRmin及びRmaxは、下記の式(13)及び下記の式(14)によって表される。このようなずらし範囲を設定するのは、脈波波形が主に被験者の脈波に由来する場合、相関係数算出部において求める相関係数が最大となるずらし範囲は脈拍1回分に相当する時間であるからである。
Rmin≦w≦Rmax・・・(12)
Rmin= 60/240・・・(13)
Rmax= 60/42・・・(14)
他の一態様としては、ずらし範囲設定部17bは、転換点の間隔から脈拍の周期を推定することもできる。所与のパラメータmに対し、最新の転換点の時刻とm*2個前の転換点の時刻との差は、おおよそ脈波m回分であるので、この被験者の現在の脈波1回分の時間は「(T(J)-T(J-2*m))/m(秒)」と推定できる。したがって、所与の0より大きく1より小さいパラメータγに対し、ずらし範囲w(秒)を下記の式(15)〜(17)で表す範囲に設定できる。
Rmin= max(60/240,γ(T(J)-T(J-2*m))/m)・・・(15)
Rmax= min(60/42,(T(J)-T(J-2*m))/m/γ)・・・(16)
Rmin≦w≦Rmax・・・(17)
更なる一態様としては、ずらし範囲設定部17bは、窓幅設定部17aにおいて、n=mに対し、最新の転換点の時刻とn*2個前の転換点の時刻との差を求め、窓幅Uとした場合、ずらし範囲w(秒)を下記の式(18)〜(20)で表す範囲に設定できる。
Rmin= max(60/240,γU/n)・・・(18)
Rmax= min(60/42,U/n/γ)・・・(19)
Rmin≦w≦Rmax・・・(20)
図9は、ずらし範囲の設定例を示す図である。図9には、m=n=2、γ=0.8である場合のずらし範囲が示されている。図9に示す縦軸は、振幅を指し、横軸は、時間を指す。また、図9には、脈波波形を実線で図示するとともに、比較窓の中に基準窓内の脈波波形を破線で重畳して図示している。図9に示すように、人が採り得る最低の脈拍数および最高の脈拍数を収容するように設定された範囲、すなわち60/150〜60/42の範囲よりも、転換点の間隔から設定された範囲、すなわち0.8U/2〜1.25U/2の範囲の方が狭い。このように、転換点の間隔から脈拍の周期を推定することで、ずらし範囲が小さくなる。これによって、被験者の脈拍数から大きく外れたずらし範囲をあらかじめ除外することで不安定な脈波を誤検出することを防ぐとともに、相関係数算出部における計算量も削減できる。なお、ここでは、窓幅Uを用いる場合を例示したが、極大点間隔を用いることもできる。
なお、上記の全ての例において、Rmax>Rminの場合、ずらし範囲が存在しなくなる。この場合には、後述の相関係数算出部18による処理は実行されない。
図1の説明に戻り、相関係数算出部18は、基準窓内の脈波波形と比較窓内の脈波波形との間で相関係数を算出する処理部である。
これを具体的に説明すると、まず、相関係数算出部18は、ずらしインデックスzの候補{z}を求める。かかる候補{z}は、ずらし範囲設定部17bで計算したずらし範囲[Rmin,Rmax]に対し、zδが[Rmin,Rmax]に含まれる全ての整数である。ここでは、{z}の最小値をZmin、{z}の最大値をZmaxと表記することとする。このとき、基準窓に含まれる脈波波形のインデックス数をXとする。上述した例のように基準窓を転換点から設定した場合、かかるXは「T(J)-T(J-2*n)+1」となる。このとき、基準窓に含まれる脈波波形の値を時刻順にx(1),・・・,x(X)は、下記のように記述できる。
x(1)=v(I-X-1)
x(2)=v(I-X-2)

x(X)=v(I)
同様に、基準窓をzずらした比較窓に含まれる脈波波形の値を時刻順にy(1),・・・,y(X)は以下のように記述できる。
y(z,1)=v(I-X-1-z)
y(z,2)=v(I-X-2-z)

y(z,X)=v(I-z)
このように、相関係数算出部18は、ずらしインデックスzの候補{z}を求めた後に、基準窓内の脈波波形と、{z}内の各zについて基準窓をzずらした比較窓内の脈波波形の各々との間で相関係数を算出する。
一態様としては、相関係数算出部18は、ずらしインデックスzをZminからZmaxまで1ずつインクリメントしながら、相関係数S(z)を計算する。かかる相関係数は、基準窓の脈波波形とzずらした比較窓の脈波波形が類似しているほど大きな値を示す指標であり、一例として、下記の式(21)にしたがって相関係数を算出することができる。このようにして全てのずらしインデックスの候補{z}について相関係数を計算したのち、S(z)が最大となるずらしインデックスをZとする。なお、x(1),・・・,x(X)の平均やy(z,1),・・・,y(z,X)の平均は、脈拍波形が0を中心に正負を繰り返す性質があるため、0で代用することとしてもかまわない。
Figure 0006098257
図10は、各窓内の波形の一例を示す図である。図10に示す縦軸は、振幅を指し、横軸は、時間を指す。また、図10には、脈波波形を実線で図示するとともに、比較窓の中に基準窓内の脈波波形を破線で重畳して図示している。図10に示すように、基準窓の脈波波形とzずらした比較窓の脈波波形は、互いの形状がよく類似している。このように、互いの波形の形状が類似する場合には相関係数が高く算出される。
その後、相関係数算出部18は、所定の閾値βに対し、S(Z)とβを比較する。このとき、相関係数算出部18は、S(Z)がβより小さい場合に、基準窓と十分類似した比較窓はみつからなかったとして、脈波検出部19による脈波の検出を実行しない。かかる閾値βには、適切な値を指定することができる。例えば、相関係数であれば、最大が1であることが明らかであるので、0.6〜0.8程度の値を指定することによって波形の周期性を判定する精度を向上できる。
ここで、上記の式(21)のように、共分散を標準偏差で除算することによって相関係数を算出する場合には、基準窓の中で脈波波形の振れ幅が大きく変化すると、互いの窓の振幅値が大きく異なるにもかかわらず、形状が類似する波形、すなわち相似の波形の相関係数が高く算出されてしまう。
図11は、各窓内の波形の一例を示す図である。図11に示す縦軸は、振幅を指し、横軸は、時間を指す。図11に示すように、測定条件の変化が原因となって基準窓の中で脈波波形の振れ幅が大きく変化するケースがある。このように、基準窓の最後に振れ幅が大きくなった場合、時刻をZずらした窓では振れ幅の大きな部分が外れてしまう。このケースでは、上記の式(21)を相関係数の算出に用いた場合、相関係数が非常に大きくなりやすい。これは、相関係数では、両窓の振れ幅を別々に正規化することが原因であり、各窓の右側の振れ幅の大きな部分が、正規化後にほぼ一致してしまっていると考えられる。この場合、各窓の残りの部分(各窓の左側)は相関係数の値にはほとんど影響を与えられない。例えば、図11の例では、相関係数は0.896となってしまう。このため、上記の式(21)を相関係数の算出に用いた場合には、測定条件が不安定なケースでも相関係数が高くなり、不安定な脈波波形を検出に用いているケースが起こりうる。
これらのことから、相関係数算出部18は、下記の式(22)に示すように、共分散を互いの窓の分散のうち大きい方の分散を選択して当該分散で除算することよって補正相関係数を算出する。かかる補正相関係数は、両窓内の脈波波形の分散が等しい場合のみ、相関係数と同じ値となる。両窓内の脈拍波形の分散の差が大きくなるにしたがって相関係数に比べて相対的に小さな値が算出される(値が正の場合)。このため、図11に示すように、両窓内の脈波波形の分散が大きく異なる場合には、補正相関係数は小さくなり、脈波波形の出力を行わなくなる。例えば、図11の例では、相関係数は0.896であるのに対し、補正相関係数は0.197と4分の1となり、閾値βに比べて小さくなる。かかる補正相関係数の算出によって、互いの窓の振幅値が大きく異なるにもかかわらず、形状が類似する波形、すなわち相似の波形の補正相関係数が不当に高く算出される事態を抑制できる。なお、補正相関係数も、最大が1であることが明らかであるので、相関係数と同様にβには0.6〜0.8程度の値を設定すればよいことが多い。
Figure 0006098257
図1の説明に戻り、脈波検出部19は、相関係数算出部18によって出力された波形から脈波を検出する処理部である。一態様としては、脈波検出部19は、被験者の脈拍数や脈波波形など安定した脈波を検出して出力する。すなわち、上記の相関係数算出部18による処理によって時刻t(I)−(U+Zδ)から時刻t(I)までの脈波波形は安定していることが明らかになっているので、その間にn+1拍の脈拍が発生していることがわかる。したがって、脈波検出部19は、下記の式(23)によって脈拍数hを算出することができる。ここで、脈拍数hが人の脈拍数がとりうる範囲、例えば42bpm〜240bpmの範囲から外れている場合、脈拍数および脈拍波形の出力は行わなくてもよい。他の一態様としては、脈波検出部19は、例えば、時刻t(I)−(U+Zδ)から時刻t(I)までの脈波波形を出力することもできる。
h=60/((U+Zδ)/(n+1))[bpm]・・・(23)
このようにして得られる脈拍数や脈波波形は、脈波検出装置10が有する図示しない表示デバイスを始め、任意の出力先へ出力することができる。例えば、脈拍数や脈拍周期のゆらぎから自律神経の働きを診断したり、脈波波形から心疾患等を診断したりする診断プログラムが脈波検出装置10にインストールされている場合には、診断プログラムを出力先とすることができる。また、診断プログラムをWebサービスとして提供するサーバ装置などを出力先とすることもできる。さらに、脈波検出装置10を利用する利用者の関係者、例えば介護士や医者などが使用する端末装置を出力先とすることもできる。これによって、院外、例えば在宅や在席のモニタリングサービスも可能になる。なお、診断プログラムの測定結果や診断結果も、脈波検出装置10を始め、関係者の端末装置に表示させることができるのも言うまでもない。
なお、上記の取得部13、抽出部14、波形検出部15、転換判定部16a、極値判定部16b、窓幅設定部17a、ずらし範囲設定部17b、相関係数算出部18及び脈波検出部19とは、CPU(Central Processing Unit)やMPU(Micro Processing Unit)などに脈波検出プログラムや信号処理プログラムを実行させることによって実現できる。また、上記の各機能部は、ASIC(Application Specific Integrated Circuit)やFPGA(Field Programmable Gate Array)などのハードワイヤードロジックによっても実現できる。
また、上記の内部メモリには、半導体メモリ素子や記憶装置を採用できる。例えば、半導体メモリ素子の一例としては、VRAM(Video Random Access Memory)、RAM(Random Access Memory)、ROM(Read Only Memory)やフラッシュメモリ(flash memory)などが挙げられる。また、メモリではなく、ハードディスク、光ディスクなどの記憶装置を用いることとしてもかまわない。
[処理の流れ]
次に、本実施例に係る脈波検出装置の処理の流れについて説明する。なお、ここでは、脈波検出装置10によって実行される(1)波形検出処理について説明した後に、(2)相関係数算出処理について説明することとする。
(1)波形検出処理
図12は、実施例1に係る波形検出処理の手順を示すフローチャートである。この波形検出処理は、カメラ11から画像が取得する度に処理を起動し、画像が取得されなくなるまで繰り返し実行される処理である。なお、図示しない入力デバイス等を介して中断操作を受け付けた場合には、波形検出処理を中止することもできる。
図12に示すように、被験者が映った画像が取得されると(ステップS101)、抽出部14は、ステップS101で取得された画像から所定の顔パーツ、例えば被験者の鼻を基準とする生体領域の画像を抽出する(ステップS102)。
その上で、抽出部14は、R信号の時系列データをBPF152R及びBPF156Rへ出力するとともに、G信号の時系列データをBPF152G及びBPF156Gへ出力する(ステップS103)。
続いて、BPF152Rは、R信号の特定周波数帯、例えば3bpm以上20bpm未満の周波数帯の信号成分を抽出するとともに、BPF152Gは、G信号の特定周波数帯の信号成分を抽出する(ステップS104A)。
そして、抽出部153Rは、R信号の特定周波数帯の信号成分の絶対強度値を抽出するとともに、抽出部153Gは、G信号の特定周波数帯の信号成分の絶対強度値を抽出する(ステップS105)。
その後、LPF154Rは、R信号の特定周波数帯の絶対強度値の時系列データに対し、時間変化に応答させる平滑化処理を実行するとともに、LPF154Gは、G信号の特定周波数帯の絶対強度値の時系列データに対し、時間変化に応答させる平滑化処理を実行する(ステップS106)。
続いて、算出部155は、LPF154Gによって出力されたG信号の特定周波数帯の絶対値強度G´noiseを、LPF154Rによって出力されたR信号の特定周波数帯の絶対値強度R´noiseで除する除算「G´noise/R´noise」を実行することによって補正係数aを算出する(ステップS107)。
上記のステップS104Aの処理に並行して、BPF156Rは、R信号の脈波周波数帯、例えば42bpm以上240bpm未満の周波数帯の信号成分を抽出するとともに、BPF156Gは、G信号の脈波周波数帯の信号成分を抽出する(ステップS104B)。
その後、乗算部157は、ステップS104Bで抽出されたR信号の脈波周波数帯の信号成分にステップS107で算出された補正係数aを乗算する(ステップS108)。その上で、演算部158は、ステップS108で補正係数aが乗算されたR信号の脈波周波数帯の信号成分から、ステップS104Bで抽出されたG信号の脈波周波数帯の信号成分を差し引く演算「a*Rsignal−Gsignal」を実行する(ステップS109)。
そして、波形検出部15は、演算後の信号の時系列データを脈波波形として転換判定部16a及び極値判定部16bへ出力し(ステップS110)、処理を終了する。
(2)相関係数算出処理
図13は、実施例1に係る相関係数算出処理の手順を示すフローチャートである。この相関係数算出処理は、波形検出部15から脈波波形が入力される限り、繰り返し実行される処理である。
図13に示すように、まず、窓の周期数n、相関係数閾値βとずらし範囲の係数γなどの初期値が各々の機能部によって図示しない内部メモリから取得される(ステップS301)。その後、転換判定部16a及び極値判定部16bは、波形検出部15から脈波波形の時刻および振幅値を取得する(ステップS302)。
そして、転換判定部16aは、波形検出部15によって出力された脈波波形の振幅値の正負が転換するか否かを判定する(ステップS303)。これとともに、極値判定部16bは、波形検出部15によって出力された脈波波形の振幅値が出力されるよりも1つ前のインデックスを持つ脈拍波形の振幅値が極大値であるか否かを判定する(ステップS304)。
続いて、ステップS302で取得された振幅値の正負が転換する場合(ステップS305Yes)には、窓幅設定部17aは、振幅値の正負が転換すると判定された転換点から脈波の周期の整数倍と対応する個数になるまで転換点を遡った区間の長さを基準窓の時間幅として設定する(ステップS306)。
続いて、ずらし範囲設定部17bは、ステップS306で設定された基準窓内に含まれる極大点の数が脈波の1周期に対応する極大点の個数を超えるか否かを判定する(ステップS307)。このとき、基準窓内に含まれる極大点の数が脈波の1周期に対応する極大点の個数nを超える場合(ステップS307Yes)には、測定条件が不安定であることが判明するので、以降の処理を実行せずに、ステップS302に移行する。
一方、基準窓内に含まれる極大点の数が脈波の1周期に対応する極大点の個数を超えない場合(ステップS307No)には、ずらし範囲設定部17bは、ステップS306で設定された基準窓の窓幅U及び窓の周期数nにずらし範囲の係数γを乗算するとともに窓幅U及び窓の周期数nにずらし範囲の係数γの逆数を乗算する計算を並行して実行することによって比較窓にずらし範囲wを設定する(ステップS308)。
このとき、ずらし範囲が存在する場合(ステップS309Yes)には、相関係数算出部18は、ずらしインデックスzをZminからZmaxまで1ずつインクリメントしながら全てのずらしインデックスについて補正相関係数を算出するまで(ステップS311No)、上記の式(22)を用いて、共分散を互いの窓の分散のうち大きい方の分散を選択して当該分散で除算することよって補正相関係数S(z)を算出する(ステップS310)。
その後、全てのずらしインデックスについて補正相関係数S(z)を算出すると(ステップS311Yes)、相関係数算出部18は、補正相関係数S(z)が閾値β以上であるか否かを判定する(ステップS312)。なお、補正相関係数S(z)が閾値βより小さい場合(ステップS312No)には、基準窓と十分類似した比較窓はみつからなかったと判断し、脈波の検出を実行せずにステップS302の処理へ戻る。
ここで、ステップS310で算出された補正相関係数S(z)のうち最大の補正相関係数S(Z)が閾値β以上である場合(ステップS312Yes)には、脈波検出部19は、上記の式(23)を用いて、脈拍数hを算出する(ステップS313)。このとき、脈波検出部19は、脈拍数hが人の脈拍数がとりうる適正範囲、例えば42bpm〜240bpmの範囲から外れている場合(ステップS314No)には、脈拍数および脈拍波形の出力は行わず、ステップS302の処理へ戻る。
一方、脈拍数hが人の脈拍数がとりうる適正範囲である場合(ステップS314Yes)には、脈波検出部19は、脈拍数および脈拍波形を所定の出力先へ出力し(ステップS315)、ステップS302の処理へ戻る。
[実施例1の効果]
上述してきたように、本実施例に係る脈波検出装置10は、波形の類似を相関係数を用いて判断する場合に、基準窓内の振幅値の分散と比較窓内の振幅値の分散との差が大きくなるにしたがって値が小さくなる相関係数を算出する信号処理を実行する。このため、本実施例に係る脈波検出装置10では、互いの窓内の波形が相似の関係にある場合でも相関係数が低く算出される。したがって、本実施例に係る脈波検出装置10によれば、波形の類似性の判定精度を向上させることができる。
さて、これまで開示の装置に関する実施例について説明したが、本発明は上述した実施例以外にも、種々の異なる形態にて実施されてよいものである。そこで、以下では、本発明に含まれる他の実施例を説明する。
[類似性の判定の応用例]
上記の実施例1では、波形の類似性、すなわち脈波の周期性を判定するのに相関係数を用いる場合を例示したが、類似性を判断するための指標は相関係数に限定されない。例えば、脈波検出装置10は、下記の式(24)を用いて、ユークリッド距離を算出したり、下記の式(25)を用いて、正規化されたユークリッド距離を算出することもできる。この場合には、全てのずらし範囲の中で最短距離を持つものが最も周期性が高いずらし範囲Zであると判断できる。このように、ユークリッド距離または正規化されたユークリッド距離を用いる場合にも、最短距離が所定の閾値THよりも大きい場合には、測定条件が不安定であると判断し、当該基準窓内の波形は脈波の検出には使用しない。
Figure 0006098257
Figure 0006098257
[入力信号]
上記の実施例1及び実施例2では、入力信号としてR信号およびG信号の二種類を用いる場合を例示したが、異なる複数の光波長成分を持つ信号であれば任意の種類の信号および任意の数の信号を入力信号とすることができる。例えば、R、G、B、IRおよびNIRなどの光波長成分が異なる信号のうち任意の組合せの信号を2つ用いることもできるし、また、3つ以上用いることもできる。
[他の実装例]
上記の実施例1では、脈波検出装置10が上記の相関係数算出処理をスタンドアローンで実行する場合を例示したが、クライアントサーバシステムとして実装することもできる。例えば、脈波検出装置10は、相関係数算出処理を実行するWebサーバとして実装することとしてもよいし、アウトソーシングによって相関係数算出サービスを始めとするサービスを提供するクラウドとして実装することとしてもかまわない。このように、脈波検出装置10がサーバ装置として動作する場合には、スマートフォンや携帯電話機等の携帯端末装置やパーソナルコンピュータ等の情報処理装置をクライアント端末として収容することができる。これらクライアント端末からネットワークを介して被験者の顔が映った画像が取得された場合に波形検出処理を実行し、その脈波波形を用いて相関係数算出処理を実行した上でその検出結果やその検出結果を用いてなされた診断結果をクライアント端末へ応答することによって脈波検出サービス及び診断サービスを提供できる。
[分散および統合]
また、図示した各装置の各構成要素は、必ずしも物理的に図示の如く構成されていることを要しない。すなわち、各装置の分散・統合の具体的形態は図示のものに限られず、その全部または一部を、各種の負荷や使用状況などに応じて、任意の単位で機能的または物理的に分散・統合して構成することができる。例えば、転換判定部16a、極値判定部16b、窓幅設定部17a、ずらし範囲設定部17bまたは相関係数算出部18を脈波検出装置10の外部装置としてネットワーク経由で接続するようにしてもよい。また、転換判定部16a、極値判定部16b、窓幅設定部17a、ずらし範囲設定部17bまたは相関係数算出部18を別の装置がそれぞれ有し、ネットワーク接続されて協働することで、上記の脈波検出装置10の機能を実現するようにしてもよい。
[信号処理プログラム]
また、上記の実施例で説明した各種の処理は、予め用意されたプログラムをパーソナルコンピュータやワークステーションなどのコンピュータで実行することによって実現することができる。そこで、以下では、図14を用いて、上記の実施例と同様の機能を有する信号処理プログラムを実行するコンピュータの一例について説明する。
図14は、実施例1及び実施例2に係る信号処理プログラムを実行するコンピュータの一例について説明するための図である。図14に示すように、コンピュータ100は、操作部110aと、スピーカ110bと、カメラ110cと、ディスプレイ120と、通信部130とを有する。さらに、このコンピュータ100は、CPU150と、ROM160と、HDD170と、RAM180とを有する。これら110〜180の各部はバス140を介して接続される。
HDD170には、図14に示すように、上記の実施例1で示した転換判定部16a、極値判定部16b、窓幅設定部17a、ずらし範囲設定部17b及び相関係数算出部18と同様の機能を発揮する信号処理プログラム170aが予め記憶される。この信号処理プログラム170aについては、図1に示した各々の転換判定部16a、極値判定部16b、窓幅設定部17a、ずらし範囲設定部17b及び相関係数算出部18の各構成要素と同様、適宜統合又は分離しても良い。すなわち、HDD170に格納される各データは、常に全てのデータがHDD170に格納される必要はなく、処理に必要なデータのみがHDD170に格納されれば良い。
そして、CPU150が、信号処理プログラム170aをHDD170から読み出してRAM180に展開する。これによって、図14に示すように、信号処理プログラム170aは、信号処理プロセス180aとして機能する。この信号処理プロセス180aは、HDD170から読み出した各種データを適宜RAM180上の自身に割り当てられた領域に展開し、この展開した各種データに基づいて各種処理を実行する。なお、信号処理プロセス180aは、図1に示した転換判定部16a、極値判定部16b、窓幅設定部17a、ずらし範囲設定部17b及び相関係数算出部18にて実行される処理、例えば図12〜図13に示す処理を含む。また、CPU150上で仮想的に実現される各処理部は、常に全ての処理部がCPU150上で動作する必要はなく、処理に必要な処理部のみが仮想的に実現されれば良い。
なお、上記の信号処理プログラム170aについては、必ずしも最初からHDD170やROM160に記憶させておく必要はない。例えば、コンピュータ100に挿入されるフレキシブルディスク、いわゆるFD、CD−ROM、DVDディスク、光磁気ディスク、ICカードなどの「可搬用の物理媒体」に各プログラムを記憶させる。そして、コンピュータ100がこれらの可搬用の物理媒体から各プログラムを取得して実行するようにしてもよい。また、公衆回線、インターネット、LAN、WANなどを介してコンピュータ100に接続される他のコンピュータまたはサーバ装置などに各プログラムを記憶させておき、コンピュータ100がこれらから各プログラムを取得して実行するようにしてもよい。
10 脈波検出装置
11 カメラ
13 取得部
14 抽出部
15 波形検出部
16a 転換判定部
16b 極値判定部
17a 窓幅設定部
17b ずらし範囲設定部
18 相関係数算出部
19 脈波検出部

Claims (8)

  1. 心拍または脈波に関する信号を検出する信号検出部と、
    前記信号検出部によって検出された信号を用いて、所定の時間幅を持つ第1の窓に含まれる信号の波形と、前記第1の窓内の波形と比較する第2の窓に含まれる信号の波形との間で、前記第1の窓内の振幅値の分散と前記第2の窓内の振幅値の分散との差が大きくなるにしたがって値が小さくなる相関係数を算出する相関係数算出部と、
    を有することを特徴とする信号処理装置。
  2. 前記信号検出部は、
    生体が撮影された画像から脈波に対応する周波数帯の信号を検出し、
    前記相関係数算出部は、
    前記信号検出部によって検出された信号を用いて、前記第1の窓内の波形と前記第2の窓内の波形との間で前記相関係数を算出することを特徴とする請求項1に記載の信号処理装置。
  3. 前記第1の窓の時間幅を設定する窓幅設定部と、
    前記第2の窓をずらす範囲を設定する範囲設定部とをさらに有し、
    前記相関係数算出部は、
    前記窓幅設定部によって時間幅が設定された第1の窓内の波形と、前記範囲設定部によって設定された範囲内で複数回にわたってずらした第2の窓内の波形の各々との間で前記相関係数を算出することを特徴とする請求項2に記載の信号処理装置。
  4. 前記信号の振幅値の正負が転換するか否かを判定する転換判定部をさらに有し、
    前記窓幅設定部は、
    前記転換判定部によって振幅値の正負が転換すると判定された転換点から前記脈波の周期の整数倍と対応する個数になるまで転換点を遡った区間の長さを前記第1の窓の時間幅として設定することを特徴とする請求項3に記載の信号処理装置。
  5. 前記信号の振幅値が極値であるか否かを判定する極値判定部をさらに有し、
    前記範囲設定部は、
    前記極値判定部によって極大であると判定された極大点のうち前記窓幅設定部によって時間幅が設定された第1の窓に含まれる極大点の数が前記脈波の周期と対応する個数と等しい場合に、前記範囲の設定を実行することを特徴とする請求項3または4に記載の信号処理装置。
  6. 前記範囲設定部は、前記窓幅設定部によって第1の窓に設定された時間幅を用いて、前記範囲を設定することを特徴とする請求項3〜5のいずれか1つに記載の信号処理装置。
  7. コンピュータが、
    心拍または脈波に関する信号を検出し、
    検出された信号を用いて、所定の時間幅を持つ第1の窓に含まれる信号の波形と、前記第1の窓内の波形と比較する第2の窓に含まれる信号の波形との間で、前記第1の窓内の振幅値の分散と前記第2の窓内の振幅値の分散との差が大きくなるにしたがって値が小さくなる相関係数を算出する
    処理を実行することを特徴とする信号処理方法。
  8. コンピュータに、
    心拍または脈波に関する信号を検出し、
    検出された信号を用いて、所定の時間幅を持つ第1の窓に含まれる信号の波形と、前記第1の窓内の波形と比較する第2の窓に含まれる信号の波形との間で、前記第1の窓内の振幅値の分散と前記第2の窓内の振幅値の分散との差が大きくなるにしたがって値が小さくなる相関係数を算出する
    処理を実行させることを特徴とする信号処理プログラム。
JP2013054165A 2013-03-15 2013-03-15 信号処理装置、信号処理方法及び信号処理プログラム Active JP6098257B2 (ja)

Priority Applications (3)

Application Number Priority Date Filing Date Title
JP2013054165A JP6098257B2 (ja) 2013-03-15 2013-03-15 信号処理装置、信号処理方法及び信号処理プログラム
US14/099,457 US9962126B2 (en) 2013-03-15 2013-12-06 Signal processor, signal processing method, and recording medium
EP13197941.1A EP2777485B1 (en) 2013-03-15 2013-12-18 Signal processor, signal processing method, and signal processing program

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2013054165A JP6098257B2 (ja) 2013-03-15 2013-03-15 信号処理装置、信号処理方法及び信号処理プログラム

Publications (2)

Publication Number Publication Date
JP2014176584A JP2014176584A (ja) 2014-09-25
JP6098257B2 true JP6098257B2 (ja) 2017-03-22

Family

ID=49918405

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2013054165A Active JP6098257B2 (ja) 2013-03-15 2013-03-15 信号処理装置、信号処理方法及び信号処理プログラム

Country Status (3)

Country Link
US (1) US9962126B2 (ja)
EP (1) EP2777485B1 (ja)
JP (1) JP6098257B2 (ja)

Families Citing this family (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP6191517B2 (ja) * 2014-03-14 2017-09-06 富士通株式会社 検出装置、検出プログラム、及び検出方法
US10397082B2 (en) * 2014-08-07 2019-08-27 Citrix Systems, Inc. Internet infrastructure measurement method and system adapted to session volume
US9582865B2 (en) * 2014-10-20 2017-02-28 Microsoft Technology Licensing, Llc Visualization for blood flow in skin image data
JP6488722B2 (ja) * 2015-01-26 2019-03-27 富士通株式会社 脈波検出装置、脈波検出方法及び脈波検出プログラム
JP6406171B2 (ja) * 2015-08-25 2018-10-17 トヨタ自動車株式会社 瞬目検知装置
JP6607259B2 (ja) * 2015-11-20 2019-11-20 富士通株式会社 脈波分析装置、脈波分析方法、および脈波分析プログラム
CN105520724A (zh) * 2016-02-26 2016-04-27 严定远 一种测量人体心跳速率和呼吸频率的方法
JP2017158675A (ja) * 2016-03-08 2017-09-14 パナソニックIpマネジメント株式会社 脈拍推定装置、脈拍推定システムおよび脈拍推定方法
KR102588906B1 (ko) * 2017-12-01 2023-10-13 삼성전자주식회사 생체 신호 품질 평가 장치 및 방법
US11457872B2 (en) 2017-12-01 2022-10-04 Samsung Electronics Co., Ltd. Bio-signal quality assessment apparatus and bio-signal quality assessment method
WO2019146025A1 (ja) * 2018-01-24 2019-08-01 富士通株式会社 脈波算出装置、脈波算出方法及び脈波算出プログラム
WO2020070808A1 (ja) * 2018-10-02 2020-04-09 富士通株式会社 脈波算出装置、脈波算出方法及び脈波算出プログラム
CN112836546A (zh) * 2019-11-22 2021-05-25 深圳市理邦精密仪器股份有限公司 检测生理信号质量的方法、装置及电子设备

Family Cites Families (16)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP3029912B2 (ja) * 1992-01-08 2000-04-10 日本コーリン株式会社 脈拍数測定装置
JP3519185B2 (ja) 1995-09-22 2004-04-12 花王株式会社 血流解析装置
FR2794961B1 (fr) * 1999-06-16 2001-09-21 Global Link Finance Procede de determination du decalage temporel entre les instants de passage d'une meme onde de pouls en deux points de mesure distincts d'un reseau arteriel d'un etre vivant et d'estimation de sa pression aortique
US9820658B2 (en) 2006-06-30 2017-11-21 Bao Q. Tran Systems and methods for providing interoperability among healthcare devices
EP1877774A4 (en) 2005-03-25 2011-01-12 Cnoga Holdings Ltd OPTICAL SENSOR DEVICE AND IMAGE PROCESSING UNIT FOR MEASURING CHEMICAL CONCENTRATIONS, CHEMICAL SATURATIONS AND BIOPHYSICAL PARAMETERS
US7725187B1 (en) * 2006-11-06 2010-05-25 Pacesetter, Inc. Motion detection for sensors sensitive to motion induced noise
JP4992043B2 (ja) * 2007-08-13 2012-08-08 株式会社国際電気通信基礎技術研究所 行動識別装置、行動識別システムおよび行動識別方法
US20100152600A1 (en) 2008-04-03 2010-06-17 Kai Sensors, Inc. Non-contact physiologic motion sensors and methods for use
CA2720871A1 (en) * 2008-04-03 2009-10-08 Kai Medical, Inc. Non-contact physiologic motion sensors and methods for use
JP4609539B2 (ja) * 2008-07-04 2011-01-12 トヨタ自動車株式会社 眠気検出装置
JP5582478B2 (ja) * 2009-07-16 2014-09-03 株式会社デルタツーリング 生体状態推定装置及びコンピュータプログラム
JP5234078B2 (ja) * 2010-09-29 2013-07-10 株式会社デンソー 脈波解析装置および血圧推定装置
US10856752B2 (en) * 2010-12-28 2020-12-08 Sotera Wireless, Inc. Body-worn system for continuous, noninvasive measurement of cardiac output, stroke volume, cardiac power, and blood pressure
JP5718126B2 (ja) * 2011-03-31 2015-05-13 沖電気工業株式会社 微細振動特徴量算出装置、微細振動特徴量算出方法及びプログラム
JP5673341B2 (ja) 2011-05-13 2015-02-18 富士通株式会社 心拍計測装置、心拍計測方法及び心拍計測プログラム
US20140073869A1 (en) * 2012-09-11 2014-03-13 Nellcor Puritan Bennett Llc Methods and systems for determining physiological information using status indicators

Also Published As

Publication number Publication date
EP2777485B1 (en) 2017-10-25
EP2777485A1 (en) 2014-09-17
JP2014176584A (ja) 2014-09-25
US20140276114A1 (en) 2014-09-18
US9962126B2 (en) 2018-05-08

Similar Documents

Publication Publication Date Title
JP6098257B2 (ja) 信号処理装置、信号処理方法及び信号処理プログラム
US10292602B2 (en) Blood flow index calculating method, blood flow index calculating apparatus, and recording medium
JP5915757B2 (ja) 脈波検出方法、脈波検出装置及び脈波検出プログラム
JP6256488B2 (ja) 信号処理装置、信号処理方法及び信号処理プログラム
JP6167614B2 (ja) 血流指標算出プログラム、血流指標算出装置および血流指標算出方法
WO2016006027A1 (ja) 脈波検出方法、脈波検出プログラム及び脈波検出装置
US20150366456A1 (en) Pulse wave velocity measurement method
JP6393984B2 (ja) 脈拍計測装置、脈拍計測方法及び脈拍計測プログラム
JP6115263B2 (ja) 脈波検出装置、脈波検出方法及び脈波検出プログラム
JP5920465B2 (ja) バイタルサイン検出方法、バイタルサイン検出装置及びバイタルサイン検出プログラム
Feng et al. Motion artifacts suppression for remote imaging photoplethysmography
US9986923B2 (en) Selecting a region of interest for extracting physiological parameters from a video of a subject
JP6135255B2 (ja) 心拍測定プログラム、心拍測定方法及び心拍測定装置
Lagido et al. Using the smartphone camera to monitor heart rate and rhythm in heart failure patients
Chwyl et al. SAPPHIRE: Stochastically acquired photoplethysmogram for heart rate inference in realistic environments
JPWO2020054122A1 (ja) 情報処理装置、プログラム及び情報処理方法
JP2021146061A (ja) 生体情報取得装置、生体情報取得方法及びプログラム
JP6020015B2 (ja) 脈波検出装置、脈波検出プログラム及び脈波検出方法
JP2016202603A (ja) 生体情報処理システム、プログラム、生体情報処理システムの制御方法
JP6488722B2 (ja) 脈波検出装置、脈波検出方法及び脈波検出プログラム
JP6167615B2 (ja) 血流指標算出プログラム、端末装置および血流指標算出方法
JP7106821B2 (ja) 生体信号解析装置、生体信号解析方法及びプログラム
JP6167849B2 (ja) 脈波検出装置、脈波検出方法及び脈波検出プログラム
Desai et al. Heart rate monitor using mobile camera
Suryasari et al. Illuminance Color Independent in Remote Photoplethysmography for Pulse Rate Variability and Respiration Rate Measurement

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20151007

A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20160921

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20160927

A521 Written amendment

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20161128

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20170124

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20170206

R150 Certificate of patent or registration of utility model

Ref document number: 6098257

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150