JP2004081513A - 心拍変動解析装置、心拍変動解析方法、および心拍変動解析用プログラム - Google Patents
心拍変動解析装置、心拍変動解析方法、および心拍変動解析用プログラム Download PDFInfo
- Publication number
- JP2004081513A JP2004081513A JP2002246635A JP2002246635A JP2004081513A JP 2004081513 A JP2004081513 A JP 2004081513A JP 2002246635 A JP2002246635 A JP 2002246635A JP 2002246635 A JP2002246635 A JP 2002246635A JP 2004081513 A JP2004081513 A JP 2004081513A
- Authority
- JP
- Japan
- Prior art keywords
- heartbeat
- signal
- subject
- information
- heart rate
- 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.)
- Granted
Links
Images
Classifications
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/24—Detecting, measuring or recording bioelectric or biomagnetic signals of the body or parts thereof
- A61B5/316—Modalities, i.e. specific diagnostic methods
- A61B5/318—Heart-related electrical modalities, e.g. electrocardiography [ECG]
- A61B5/346—Analysis of electrocardiograms
- A61B5/349—Detecting specific parameters of the electrocardiograph cycle
- A61B5/352—Detecting R peaks, e.g. for synchronising diagnostic apparatus; Estimating R-R interval
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
- A61B5/00—Measuring for diagnostic purposes; Identification of persons
- A61B5/02—Detecting, 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/024—Detecting, measuring or recording pulse rate or heart rate
- A61B5/02405—Determining heart rate variability
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Cardiology (AREA)
- Medical Informatics (AREA)
- Surgery (AREA)
- Biophysics (AREA)
- Pathology (AREA)
- Engineering & Computer Science (AREA)
- Biomedical Technology (AREA)
- Heart & Thoracic Surgery (AREA)
- Veterinary Medicine (AREA)
- Molecular Biology (AREA)
- Physics & Mathematics (AREA)
- Animal Behavior & Ethology (AREA)
- General Health & Medical Sciences (AREA)
- Public Health (AREA)
- Physiology (AREA)
- Measurement And Recording Of Electrical Phenomena And Electrical Characteristics Of The Living Body (AREA)
- Measuring Pulse, Heart Rate, Blood Pressure Or Blood Flow (AREA)
- Measuring And Recording Apparatus For Diagnosis (AREA)
Abstract
【課題】ローパスフィルタの設計が容易で、ノイズの混入が少なく高精度に心拍ゆらぎを計測する。
【解決手段】心電図情報から被検体の心拍を表す心拍信号を心拍信号算出手段１２で算出し、この心拍信号から被検体の呼吸状態を反映させた呼吸信号を呼吸信号演算手段１３により演算する。これら心拍信号と呼吸信号とに基づいて被検体の呼吸状態を反映させた被検体の心拍のゆらぎ情報を心拍ゆらぎ情報提供手段１４から提供する。
【選択図】 図５
【解決手段】心電図情報から被検体の心拍を表す心拍信号を心拍信号算出手段１２で算出し、この心拍信号から被検体の呼吸状態を反映させた呼吸信号を呼吸信号演算手段１３により演算する。これら心拍信号と呼吸信号とに基づいて被検体の呼吸状態を反映させた被検体の心拍のゆらぎ情報を心拍ゆらぎ情報提供手段１４から提供する。
【選択図】 図５
Description
【０００１】
【発明の属する技術分野】
本発明は、心電図情報から心拍のゆらぎ成分を求めるための心拍変動解析装置、心拍変動解析方法、および心拍変動解析用プログラムに関する。
【０００２】
【従来の技術】
従来、この種の心拍変動解析方法には、ＣＤＭ（Ｃｏｍｐｌｅｘ Ｄｅｍｏｄｕｌａｔｉｏｎ）法がある。このＣＤＭ法は、図１に示すように心電図計測装置１により計測した被検体の心電図情報のピーク間隔をピーク間隔検出器２により検出してＲ−Ｒ間隔データを算出し、このＲ−Ｒ間隔データに０．３Ｈｚ正弦波発生器３から発生した０．３Ｈｚの正弦波を変調器４で乗算し、直交した中間周波数成分Ｉ（ｉｎ Ｐｈａｓｅ）、Ｑ（ｑｕａｄｒａｔｕｒｅ）信号を生成する。
【０００３】
これらＩ，Ｑ信号から復調器５において帯域幅０．１５Ｈｚのローパスフィルタにより復調し、心拍ゆらぎ出力装置６から心拍のゆらぎ情報を出力している。
【０００４】
ここで、心拍のゆらぎは心拍変動であり、この心拍変動とは心周期のゆらぎのうち、洞結節に対する自律神経入力のゆらぎを原因とする心周期の変動であり、心拍変動を解析するには、心電図の１拍ごとのピーク間隔（以下、Ｒ−Ｒ間隔ともいう。）を計測するのが一般的である。
【０００５】
次に、ＣＤＭ法により心拍ゆらぎの呼吸成分の大きさを求めるための数学的手法について説明する。
【０００６】
まず、呼吸周波数をｆｒ、位相をφ、中間周波数をｆｗとする。そして、心拍ゆらぎを含むＲ−Ｒ間隔の信号を
【０００７】
【数１】
ｙ＝Ａ×ｓｉｎ（２π×ｆｒ×ｔ＋φ）
とする。ここで、Ａは心拍ゆらぎの呼吸性成分の大きさ、ｆｒは呼吸の周波数、ｔは時間である。
【０００８】
そして、０．３Ｈｚ正弦波発生器３から発生したローカル信号は、
【０００９】
【数２】
ｓｉｎ（２π×ｆｗ×ｔ）
【００１０】
【数３】
ｃｏｓ（２π×ｆｗ×ｔ）
また、中間周波数信号Ｉ，Ｑは、それぞれＲ−Ｒ間隔の信号と上記ローカル信号とを変調器４において乗算すると、
【００１１】
【数４】
Ｉ＝ｙ×ｓｉｎ（２π×ｆｗ×ｔ）
＝ｓｉｎ（２π×ｆｒ×ｔ＋φ）×Ａ×ｓｉｎ（２π×ｆｒ×ｔ）
【００１２】
【数５】
Ｑ＝ｙ×ｃｏｓ（２π×ｆｗ×ｔ）
＝ｓｉｎ（２π×ｆｒ×ｔ＋φ）×Ａ×ｃｏｓ（２π×ｆｗ×ｔ）
これらの式は三角関数の積和公式より
【００１３】
【数６】
ｓｉｎ（ｘ）×ｓｉｎ（ｙ）
＝１／２×（ｃｏｓ（ｘ−ｙ）−ｃｏｓ（ｘ＋ｙ））
【００１４】
【数７】
ｃｏｓ（ｘ）×ｓｉｎ（ｙ）
＝１／２×（ｓｉｎ（ｘ−ｙ）−ｓｉｎ（ｘ＋ｙ））
となるので、
【００１５】
【数８】
Ｉ＝１／２×Ａ×（ｃｏｓ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
−ｃｏｓ（２π×ｆｗ×ｔ＋φ＋ｆｒ×２π×ｔ）
【００１６】
【数９】
Ｑ＝１／２×Ａ×（ｓｉｎ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
−ｓｉｎ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
となる。
【００１７】
ここで、呼吸の周波数ｆｒは０．１５〜０．４５Ｈｚであるので、０．３Ｈｚ正弦波発生器３から発生したローカル信号ｆｗを０．３Ｈｚとすると、ｆｗ−ｆｒは−０．１５Ｈｚ〜＋０．１５Ｈｚ、ｆｗ＋ｆｒは０．４５〜０．７５Ｈｚである。すなわち、復調器５において通過帯域０．１５Ｈｚ、遮断帯域０．４５〜０．７５Ｈｚのローパスフィルタを用いればよいことになる。
【００１８】
中間周波数信号Ｉ，Ｑに帯域幅０．１５Ｈｚのローパスフィルタを通すことにより、式の第２項が消去でき、それぞれＩＸ，ＱＸとすると、
【００１９】
【数１０】
ＩＸ＝１／２×Ａ×（ｃｏｓ（２π×（ｆｗ−ｆｒ）×ｔ＋φ）
【００２０】
【数１１】
ＱＸ＝１／２×Ａ×（ｓｉｎ（２π×（ｆｗ−ｆｒ）×ｔ＋φ）
となる。
【００２１】
三角関数の公式より
【００２２】
【数１２】
ｓｉｎ（ｘ）×ｓｉｎ（ｘ）＋ｃｏｓ（ｘ）×ｃｏｓ（ｘ）＝１
であるので、
【００２３】
【数１３】
（ＩＸ）２＋（ＱＸ）２＝（１／２×Ａ）２
故に
【００２４】
【数１４】
Ａ＝２×（（ＩＸ）２＋（ＱＸ）２）０．５
である。
【００２５】
すなわち、図２（Ａ）〜図２（Ｅ）は従来の信号を示す波形図である。入力信号として振幅１、周波数０．２Ｈｚの正弦波を心拍変動とする。図２（Ａ）は入力信号、図２（Ｂ）は上述した０．３Ｈｚのローカル信号、図２（Ｃ）は上記入力信号と上記ローカル信号を掛け合わせた中間周波数信号（Ｉ，Ｑ）である。
【００２６】
また、図２（Ｄ）は中間周波数（Ｉ，Ｑ）にローパスフィルタを通した信号（ＩＸ，ＱＸ）、図２（Ｅ）はローパスフィルタを通した信号（ＩＸ，ＱＸ）より計算した振幅成分である。
【００２７】
図３は従来の信号と折り返しノイズとの関係を示す説明図である。図３に示すように、呼吸の周波数ｆｒは０．１５〜０．４５Ｈｚであるので、０．３Ｈｚ正弦波発生器３から発生したローカル信号ｆｗが０．３Ｈｚであることから、信号は０〜０．１５Ｈｚの範囲であり、折り返しノイズは０．４５〜０．７５Ｈｚである。この折り返しノイズを除去するローパスフィルタは、図３の破線に示す特性となる。すなわち、従来のローパスフィルタは、通過帯域が０．１５Ｈｚであって、遮断帯域が０．４５〜０．７５Ｈｚである。
【００２８】
図４（Ａ），（Ｂ）は実際に用いるローパスフィルタの特性を示す波形図である。図４（Ａ）に示すように、減衰量は遮断帯域０．４５Ｈｚで−２２ｄＢ程度である。
【００２９】
【発明が解決しようとする課題】
ところで、上述した従来のＣＤＭ法を用いた心拍変動解析方法では、図３に示すようにローパスフィルタの通過帯域が０．１５Ｈｚ、遮断帯域が０．４５〜０．７５Ｈｚである。すなわち、ローパスフィルタの通過帯域幅が０．１５Ｈｚと大きいため、図２（Ｅ）に示すように振幅信号にノイズが重畳し、呼吸以外の周波数成分が混入するという課題が一例として挙げられる。
【００３０】
また、呼吸が乱れた場合、従来の心拍変動解析方法では、非定常性のノイズが混入し、呼吸周波数以外に混入するノイズが大きいという課題もあり、その場合、呼吸の乱れを検出する手段がなかった。そして、算出した心拍変動成分の信頼性が不明であるという課題も一例として挙げられる。
【００３１】
さらに、従来の心拍変動解析方法では、図４（Ａ），（Ｂ）に示すようにローパスフィルタが遮断周波数においても十分な減衰量を得ることができないという課題が一例として挙げられる。
【００３２】
本発明が解決しようとする課題としては、ローパスフィルタの設計が容易で、ノイズの混入が少なく高精度に心拍ゆらぎを計測可能な心拍変動解析装置、心拍変動解析方法、および心拍変動解析用プログラムを提供することが一例として挙げられる。
【００３３】
【課題を解決するための手段】
上記課題を解決するため、請求項１記載の発明は、被検体の心電図情報を検出する心電図情報検出手段と、前記心電図情報から前記被検体の心拍を表す心拍信号を算出する心拍信号算出手段と、前記心拍信号から前記被検体の呼吸状態を反映させた呼吸信号を演算する呼吸信号演算手段と、前記心拍信号と前記呼吸信号とに基づいて前記被検体の呼吸状態を反映させた前記被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供手段と、を備えたことを特徴とする。
【００３４】
上記課題を解決するため、請求項７記載の発明は、被検体の心電図情報を検出する心電図情報検出工程と、前記心電図情報から前記被検体の心拍を表す心拍信号を算出する心拍信号算出工程と、前記心拍信号から前記被検体の呼吸状態を反映させた呼吸信号を演算する呼吸信号演算工程と、前記心拍信号と前記呼吸信号とに基づいて前記被検体の呼吸状態を反映させた前記被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供工程と、を備えたことを特徴とする。
【００３５】
上記課題を解決するため、請求項１３記載の発明は、心拍変動解析装置に含まれるコンピュータを、被検体の心電図情報を検出する心電図情報検出手段、前記心電図情報から前記被検体の心拍を表す心拍信号を算出する心拍信号算出手段、前記心拍信号から前記被検体の呼吸状態を反映させた呼吸信号を演算する呼吸信号演算手段、前記心拍信号と前記呼吸信号とに基づいて前記被検体の呼吸状態を反映させた前記被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供手段、として機能させることを特徴とする。
【００３６】
【発明の実施の形態】
以下、本願の実施形態を図面に基づいて説明する。
【００３７】
図５は本願に係る心拍変動解析装置の第１実施形態の構成を示すブロック図である。なお、本実施形態の心拍変動解析装置は、例えば車両の通信ナビゲーションシステムなどに搭載され、被検体である運転者の心拍ゆらぎを解析し、その解析結果から精神状態や体調に合わせて提供する音楽や走行ルートを選択するために用いられる。この場合、運転者の心電図情報は、ハンドルセンサを通して心電図計測装置により検出される。
【００３８】
図５に示すように、本実施形態の心拍変動解析装置は、被検体の心電図情報を検出する心電図情報検出手段の一例としての心電図計測装置１１と、この心電図計測装置１１の心電図情報から被検体の心拍を表す心拍信号を算出する心拍信号算出手段の一例としてのピーク間隔検出器１２と、このピーク間隔検出器１２により算出された心拍信号から被検体の呼吸状態を反映させた呼吸信号を演算する呼吸信号演算手段の一例としての自己相関演算器１３と、ピーク間隔検出器１２により算出された心拍信号と自己相関演算器１３により演算された呼吸信号とに基づいて被検体の呼吸状態を反映させた被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供手段の一例としての心拍ゆらぎ情報提供装置１４とから大略構成されている。
【００３９】
ピーク間隔検出器１２は、心電図情報を表す心電図波形のピーク間隔データを算出する手段であって、このピーク間隔データは心電図情報のＲ波相互間の間隔（Ｒ−Ｒ間隔）データである。
【００４０】
自己相関演算器１３は、心電図波形のピーク間隔データ、すなわちＲ−Ｒ間隔データの自己相関値を演算するものである。
【００４１】
心拍ゆらぎ情報提供装置１４は、心電図波形のピーク間隔データの自己相関値から被検体の呼吸数を算出する呼吸数算出手段の一例としての呼吸数算出器１５と、ピーク間隔データを、呼吸数に基づいて周波数を変えた変調信号で直交変更する変調手段の一例としての変調器１６と、この変調器１６により変調された変調信号を復調する復調手段の一例としての復調器１７と、この復調器１７による復調結果を心拍ゆらぎ情報として出力する心拍ゆらぎ情報出力手段の一例としての心拍ゆらぎ情報出力装置１８と、を備えて構成されている。
【００４２】
自己相関演算器１３により演算されたＲ−Ｒ間隔データの自己相関値は、信頼度情報出力手段の一例としての信頼度情報出力装置１９に入力され、この信頼度情報出力装置１９は被検体の呼吸状態の信頼度を表す情報として出力し、自己相関値が閾値０．４より大きい場合には、呼吸が安定し信頼度が高いと判断される。
【００４３】
上記ピーク間隔検出器１２、自己相関演算器１３、呼吸数算出器１５、変調器１６、および復調器１７はＣＰＵ２０から構成され、このＣＰＵ２０はメモリ２１との間でデータの授受が実行され、ＣＰＵ２０内の各機器はメモリ２１に予め記憶されているプログラムに従って動作する。
【００４４】
次に、本実施形態の心拍変動解析処理を図６のフローチャートに従って説明する。
【００４５】
まず、心電図計測装置１１により計測した心電図情報を読込み、その心電図情報をメモリ２１に記憶する（ステップＳ１）。次いで、ステップＳ２で心電図波形のピーク間隔データであるＲ−Ｒ間隔データの算出タイミングであるか否かを判断し、Ｒ−Ｒ間隔データの算出タイミングである場合（ステップＳ２；ＹＥＳ）にはステップＳ３に進む一方、Ｒ−Ｒ間隔データの算出タイミングでない場合（ステップＳ２；ＮＯ）にはステップＳ１に戻り、Ｒ−Ｒ間隔データの算出タイミングになるまで待機状態となる。
【００４６】
ステップＳ３では、ピーク間隔検出器１２によりＲ−Ｒ間隔データを算出するとともに、そのＲ−Ｒ間隔データをメモリ２１に記憶する。そして、ピーク間隔検出器１２により算出したＲ−Ｒ間隔データの自己相関値を自己相関演算器１３により求めてメモリ２１に記憶（ステップＳ４）した後、呼吸数算出器１５により呼吸周波数を算出する（ステップＳ５）。ここで、ピーク間隔のデータ（不等間隔データ）をスプライン補間して１０Ｈｚでリサンプリングを行う。
【００４７】
さらに、呼吸周波数に応じて変調信号の周波数を変えた後、その変えた周波数を用いて変調器１６により直交変調処理を実行する（ステップＳ６，Ｓ７）。すなわち、ステップＳ７では、Ｒ−Ｒ間隔データに呼吸周波数の位相の９０度異なる正弦波を掛け合わせ、周波数変調を行う。次のステップＳ８の復調処理では、周波数変調した２つの変調波を掛け合わせ、直流帯域通過、遮断帯域０．３〜０．９Ｈｚのローパスフィルタを通して復調する。そして、心拍ゆらぎ成分を取り出して心拍ゆらぎ出力装置１８から出力する（ステップＳ９）。
【００４８】
次いで、自己相関値を閾値０．４と比較し、その比較結果に応じた出力態様で表示する（ステップＳ１０，Ｓ１１）。すなわち、ステップＳ１０，Ｓ１１では、自己相関値が閾値０．４より大きい場合に呼吸が安定していると判断し、心拍ゆらぎの信頼度が高いと判断し、その旨を信頼度情報出力装置１９から出力表示する。一方、自己相関値が閾値０．４より小さい場合には、呼吸が不安定であると判断し、心拍ゆらぎの信頼度が低いと判断し、その旨を信頼度情報出力装置１９から出力表示する。つまり、信頼度情報出力装置１９では、心拍ゆらぎの信頼度が高い場合と低い場合とで異なるように画像や音声にて出力表示する。
【００４９】
これら一連の処理が終了したら、全体の処理を終了するかを判断し、終了の場合（ステップＳ１２；ＹＥＳ）には、全体の処理を終了する一方、終了しない場合（ステップＳ２；ＮＯ）にはステップＳ１に戻り、上記一連の処理を再度実行する。
【００５０】
次に、本実施形態により心拍ゆらぎの呼吸成分の大きさを求めるための数学的手法について説明する。
【００５１】
まず、呼吸周波数をｆｒ、位相をφ、中間周波数をｆｗとする。そして、心拍ゆらぎを含むＲ−Ｒ間隔の信号を
【００５２】
【数１５】
ｙ＝Ａ×ｓｉｎ（２π×ｆｒ×ｔ＋φ）
とする。ここで、Ａは心拍ゆらぎの呼吸性成分の大きさ、ｆｒは呼吸の周波数、ｔは時間である。
【００５３】
また、中間周波数信号Ｉ，Ｑは、それぞれＲ−Ｒ間隔の信号とローカル信号とを変調器１６において乗算すると、
【００５４】
【数１６】
Ｉ＝１／２×Ａ×（ｃｏｓ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
−ｃｏｓ（２π×ｆｗ×ｔ＋φ＋ｆｒ×２π×ｔ）
【００５５】
【数１７】
Ｑ＝１／２×Ａ×（ｓｉｎ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
−ｓｉｎ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
呼吸の周波数ｆｒはローカル信号ｆｗと等しいので、
【００５６】
【数１８】
Ｉ＝１／２×Ａ×（φ）−ｃｏｓ（２×２π×ｆｗ×ｔ＋φ）
【００５７】
【数１９】
Ｑ＝１／２×Ａ×（φ）−ｓｉｎ（２×２π×ｆｗ×ｔ＋φ）
２×ｆｗは、０．３〜０．９Ｈｚである。したがって、本実施形態では、直流成分が通過し、遮断帯域０．３〜０．９Ｈｚのローパスフィルタを用いればよいことになる。
【００５８】
中間周波数信号Ｉ，Ｑに直流帯域通過のローパスフィルタを通すことにより式の第２項が消去でき、それぞれＩＸ，ＱＸとすると、
【００５９】
【数２０】
ＩＸ＝１／２×Ａ×（ｃｏｓ（２π×（ｆｗ−ｆｒ）×ｔ＋φ）
【００６０】
【数２１】
ＱＸ＝１／２×Ａ×（ｓｉｎ（２π×（ｆｗ−ｆｒ）×ｔ＋φ）
三角関数の公式より
【００６１】
【数２２】
ｓｉｎ（ｘ）×ｓｉｎ（ｘ）＋ｃｏｓ（ｘ）×ｃｏｓ（ｘ）＝１
であるので、
【００６２】
【数２３】
（ＩＸ）２＋（ＱＸ）２＝（１／２×Ａ）２
故に
【００６３】
【数２４】
Ａ＝２×（（ＩＸ）２＋（ＱＸ）２）０．５
である。
【００６４】
すなわち、図７（Ａ）〜図７（Ｅ）は本実施形態の信号を示す波形図である。入力信号として振幅１、周波数０．２Ｈｚの正弦波を心拍変動とする。図７（Ａ）は入力信号、図７（Ｂ）は０．２Ｈｚのローカル信号、図７（Ｃ）は上記入力信号と上記ローカル信号を掛け合わせた中間周波数信号（Ｉ，Ｑ）である。
【００６５】
また、図７（Ｄ）は中間周波数信号（Ｉ，Ｑ）をローパスフィルタに通した信号（ＩＸ，ＱＸ）、図７（Ｅ）はローパスフィルタに通した信号（ＩＸ，ＱＸ）より計算した振幅成分である。したがって、本実施形態では、図７（Ｅ）に示すように計算した振幅成分にノイズが重畳していないことが判る。
【００６６】
図８は本実施形態の信号と折り返しノイズとの関係を示す説明図である。図８に示すように、直流信号成分は０〜０．０５Ｈｚの範囲であり、折り返しノイズは０．３〜０．９Ｈｚの範囲である。この折り返しノイズを除去する復調器１７のローパスフィルタは、図８の破線に示す特性となる。これにより、本実施形態では、通過帯域幅が０．０５Ｈｚと狭く、設計の容易なローパスフィルタを用いることができる。
【００６７】
図９（Ａ），（Ｂ）は本実施形態で用いるローパスフィルタの特性を示す波形図である。本実施形態では、図９（Ａ）に示すように遮断周波数においても十分な減衰効果が得られる。
【００６８】
次に、自己相関による呼吸数と相関値の算出方法を説明する。図１０は心拍変動の信号を示す波形図である。図１０は、１０Ｈｚサンプリングで６７ポイントのデータである。この入力信号の自己相関係数を調べると図１１になる。図１１に示すように遅延時間３１ポイントで自己相関係数が０．４５である。すなわち、呼吸の周期３．１秒、周波数は１／３．１＝０．３２Ｈｚであり、自己相関係数が閾値０．４より大きいので、呼吸が安定していることになる。
【００６９】
このように本実施形態によれば、被検体の心電図情報を検出する心電図計測装置１１と、この心電図計測装置１１の心電図情報から被検体の心拍を表す心拍信号を算出するピーク間隔検出器１２と、このピーク間隔検出器１２により算出された心拍信号と呼吸信号とに基づいて被検体の呼吸状態を反映させた被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供装置１４とを備えたことにより、ローパスフィルタの設計が容易で、ノイズの混入が少なく、呼吸の計測が不要で、高精度に心拍ゆらぎを計測することができる。
【００７０】
また、本実施形態によれば、ピーク間隔検出器１２は、心電図情報を表す心電図波形のピーク間隔データを算出する手段であって、このピーク間隔データは心電図情報のＲ波相互間の間隔（Ｒ−Ｒ間隔）データであることにより、高精度に心拍ゆらぎを計測することができる。
【００７１】
さらに、本実施形態によれば、心拍ゆらぎ情報提供装置１４は、心電図波形のピーク間隔データの自己相関値から被検体の呼吸数を算出する呼吸数算出器１５と、ピーク間隔データを、呼吸数に基づいて周波数を変えた変調信号で直交変更する変調器１６と、この変調器１６により変調された変調信号を復調する復調器１７と、この復調器１７による復調結果を心拍ゆらぎ情報として出力する心拍ゆらぎ情報出力装置１８とを備えたことにより、復調器１７に用いるローパスフィルタの設計が容易で、ノイズの混入が少なく、呼吸の計測が不要で、高精度に心拍ゆらぎを計測することができる。
【００７２】
さらにまた、本実施形態によれば、自己相関演算器１３は、心電図波形のピーク間隔データ、すなわちＲ−Ｒ間隔データの自己相関値を演算し、この自己相関値は、信頼度情報出力装置１９に入力され、この信頼度情報出力装置１９は被検体の呼吸状態の信頼度を表す情報として出力することにより、計測した心拍ゆらぎの信頼度を計測することができる。
【００７３】
図１２は本願に係る心拍変動解析装置の第２実施形態の構成を示すブロック図である。なお、前記第１実施形態と同一の部分には、同一の符号を付して説明を省略する。
【００７４】
前記第１実施形態は、自己相関関数を用いて呼吸周波数および相関値を算出したが、本実施形態は、他の周波数解析方法であるＦＦＴ（高速フーリェ変換）２２にてスペクトルを調べて、ピーク周波数を算出して中間周波数としている。
【００７５】
このように、ＦＦＴ２２にてスペクトルを分析し、ピーク周波数を算出して中間周波数としても前記第１実施形態と同様の作用および効果が得られる。なお、ＦＦＴ２２に代えてＭＥＭ（最大エントロピー法），ＡＲ（自己回帰）モデルを用いてスペクトルを分析すれば、ＦＦＴよりも短いデータから分解能の高いスムーズなスペクトルが得られるなどの効果を奏する。
【００７６】
なお、上記各実施形態は、本発明を車両の通信ナビゲーションシステムなどに搭載して被検体である運転者の心拍ゆらぎを解析するために用いたが、これに限らず、例えばスポーツジムにおける各種運動器具に搭載して心拍ゆらぎを解析するようにしてもよい。
【００７７】
また、図６に示したフローチャートに対応するプログラムを、フレキシブルディスクまたはハードディスクなどの情報記録媒体に記憶させておき、これを汎用のマイクロコンピュータなどに読み出して実行させることで、当該マイクロコンピュータを上記実施形態に係るＣＰＵ２０として機能させることも可能である。
【図面の簡単な説明】
【図１】従来の心拍変動解析装置の構成を示すブロック図である。
【図２】（Ａ）〜（Ｅ）は従来の信号を示す波形図である。
【図３】従来の信号と折り返しノイズとの関係を示す説明図である。
【図４】（Ａ），（Ｂ）は従来例で用いるローパスフィルタの特性を示す波形図である。
【図５】本願に係る心拍変動解析装置の第１実施形態の構成を示すブロック図である。
【図６】第１実施形態の心拍変動解析処理を示すフローチャートである。
【図７】（Ａ）〜（Ｅ）は第１実施形態の信号を示す波形図である。
【図８】第１実施形態の信号と折り返しノイズとの関係を示す説明図である。
【図９】（Ａ），（Ｂ）は第１実施形態で用いるローパスフィルタの特性を示す波形図である。
【図１０】第１実施形態の心拍変動の信号を示す波形図である。
【図１１】第１実施形態の入力信号の自己相関係数を示す波形図である。
【図１２】本願に係る心拍変動解析装置の第２実施形態の構成を示すブロック図である。
【符号の説明】
１１ 心電図計測装置
１２ ピーク間隔検出器
１３ 自己相関演算器
１４ 心拍ゆらぎ情報提供装置
１５ 呼吸数算出器
１６ 変調器
１７ 復調器
１８ 心拍ゆらぎ情報出力装置
１９ 信頼度情報出力装置
２０ ＣＰＵ
２１ メモリ
２２ ＦＦＴ
【発明の属する技術分野】
本発明は、心電図情報から心拍のゆらぎ成分を求めるための心拍変動解析装置、心拍変動解析方法、および心拍変動解析用プログラムに関する。
【０００２】
【従来の技術】
従来、この種の心拍変動解析方法には、ＣＤＭ（Ｃｏｍｐｌｅｘ Ｄｅｍｏｄｕｌａｔｉｏｎ）法がある。このＣＤＭ法は、図１に示すように心電図計測装置１により計測した被検体の心電図情報のピーク間隔をピーク間隔検出器２により検出してＲ−Ｒ間隔データを算出し、このＲ−Ｒ間隔データに０．３Ｈｚ正弦波発生器３から発生した０．３Ｈｚの正弦波を変調器４で乗算し、直交した中間周波数成分Ｉ（ｉｎ Ｐｈａｓｅ）、Ｑ（ｑｕａｄｒａｔｕｒｅ）信号を生成する。
【０００３】
これらＩ，Ｑ信号から復調器５において帯域幅０．１５Ｈｚのローパスフィルタにより復調し、心拍ゆらぎ出力装置６から心拍のゆらぎ情報を出力している。
【０００４】
ここで、心拍のゆらぎは心拍変動であり、この心拍変動とは心周期のゆらぎのうち、洞結節に対する自律神経入力のゆらぎを原因とする心周期の変動であり、心拍変動を解析するには、心電図の１拍ごとのピーク間隔（以下、Ｒ−Ｒ間隔ともいう。）を計測するのが一般的である。
【０００５】
次に、ＣＤＭ法により心拍ゆらぎの呼吸成分の大きさを求めるための数学的手法について説明する。
【０００６】
まず、呼吸周波数をｆｒ、位相をφ、中間周波数をｆｗとする。そして、心拍ゆらぎを含むＲ−Ｒ間隔の信号を
【０００７】
【数１】
ｙ＝Ａ×ｓｉｎ（２π×ｆｒ×ｔ＋φ）
とする。ここで、Ａは心拍ゆらぎの呼吸性成分の大きさ、ｆｒは呼吸の周波数、ｔは時間である。
【０００８】
そして、０．３Ｈｚ正弦波発生器３から発生したローカル信号は、
【０００９】
【数２】
ｓｉｎ（２π×ｆｗ×ｔ）
【００１０】
【数３】
ｃｏｓ（２π×ｆｗ×ｔ）
また、中間周波数信号Ｉ，Ｑは、それぞれＲ−Ｒ間隔の信号と上記ローカル信号とを変調器４において乗算すると、
【００１１】
【数４】
Ｉ＝ｙ×ｓｉｎ（２π×ｆｗ×ｔ）
＝ｓｉｎ（２π×ｆｒ×ｔ＋φ）×Ａ×ｓｉｎ（２π×ｆｒ×ｔ）
【００１２】
【数５】
Ｑ＝ｙ×ｃｏｓ（２π×ｆｗ×ｔ）
＝ｓｉｎ（２π×ｆｒ×ｔ＋φ）×Ａ×ｃｏｓ（２π×ｆｗ×ｔ）
これらの式は三角関数の積和公式より
【００１３】
【数６】
ｓｉｎ（ｘ）×ｓｉｎ（ｙ）
＝１／２×（ｃｏｓ（ｘ−ｙ）−ｃｏｓ（ｘ＋ｙ））
【００１４】
【数７】
ｃｏｓ（ｘ）×ｓｉｎ（ｙ）
＝１／２×（ｓｉｎ（ｘ−ｙ）−ｓｉｎ（ｘ＋ｙ））
となるので、
【００１５】
【数８】
Ｉ＝１／２×Ａ×（ｃｏｓ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
−ｃｏｓ（２π×ｆｗ×ｔ＋φ＋ｆｒ×２π×ｔ）
【００１６】
【数９】
Ｑ＝１／２×Ａ×（ｓｉｎ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
−ｓｉｎ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
となる。
【００１７】
ここで、呼吸の周波数ｆｒは０．１５〜０．４５Ｈｚであるので、０．３Ｈｚ正弦波発生器３から発生したローカル信号ｆｗを０．３Ｈｚとすると、ｆｗ−ｆｒは−０．１５Ｈｚ〜＋０．１５Ｈｚ、ｆｗ＋ｆｒは０．４５〜０．７５Ｈｚである。すなわち、復調器５において通過帯域０．１５Ｈｚ、遮断帯域０．４５〜０．７５Ｈｚのローパスフィルタを用いればよいことになる。
【００１８】
中間周波数信号Ｉ，Ｑに帯域幅０．１５Ｈｚのローパスフィルタを通すことにより、式の第２項が消去でき、それぞれＩＸ，ＱＸとすると、
【００１９】
【数１０】
ＩＸ＝１／２×Ａ×（ｃｏｓ（２π×（ｆｗ−ｆｒ）×ｔ＋φ）
【００２０】
【数１１】
ＱＸ＝１／２×Ａ×（ｓｉｎ（２π×（ｆｗ−ｆｒ）×ｔ＋φ）
となる。
【００２１】
三角関数の公式より
【００２２】
【数１２】
ｓｉｎ（ｘ）×ｓｉｎ（ｘ）＋ｃｏｓ（ｘ）×ｃｏｓ（ｘ）＝１
であるので、
【００２３】
【数１３】
（ＩＸ）２＋（ＱＸ）２＝（１／２×Ａ）２
故に
【００２４】
【数１４】
Ａ＝２×（（ＩＸ）２＋（ＱＸ）２）０．５
である。
【００２５】
すなわち、図２（Ａ）〜図２（Ｅ）は従来の信号を示す波形図である。入力信号として振幅１、周波数０．２Ｈｚの正弦波を心拍変動とする。図２（Ａ）は入力信号、図２（Ｂ）は上述した０．３Ｈｚのローカル信号、図２（Ｃ）は上記入力信号と上記ローカル信号を掛け合わせた中間周波数信号（Ｉ，Ｑ）である。
【００２６】
また、図２（Ｄ）は中間周波数（Ｉ，Ｑ）にローパスフィルタを通した信号（ＩＸ，ＱＸ）、図２（Ｅ）はローパスフィルタを通した信号（ＩＸ，ＱＸ）より計算した振幅成分である。
【００２７】
図３は従来の信号と折り返しノイズとの関係を示す説明図である。図３に示すように、呼吸の周波数ｆｒは０．１５〜０．４５Ｈｚであるので、０．３Ｈｚ正弦波発生器３から発生したローカル信号ｆｗが０．３Ｈｚであることから、信号は０〜０．１５Ｈｚの範囲であり、折り返しノイズは０．４５〜０．７５Ｈｚである。この折り返しノイズを除去するローパスフィルタは、図３の破線に示す特性となる。すなわち、従来のローパスフィルタは、通過帯域が０．１５Ｈｚであって、遮断帯域が０．４５〜０．７５Ｈｚである。
【００２８】
図４（Ａ），（Ｂ）は実際に用いるローパスフィルタの特性を示す波形図である。図４（Ａ）に示すように、減衰量は遮断帯域０．４５Ｈｚで−２２ｄＢ程度である。
【００２９】
【発明が解決しようとする課題】
ところで、上述した従来のＣＤＭ法を用いた心拍変動解析方法では、図３に示すようにローパスフィルタの通過帯域が０．１５Ｈｚ、遮断帯域が０．４５〜０．７５Ｈｚである。すなわち、ローパスフィルタの通過帯域幅が０．１５Ｈｚと大きいため、図２（Ｅ）に示すように振幅信号にノイズが重畳し、呼吸以外の周波数成分が混入するという課題が一例として挙げられる。
【００３０】
また、呼吸が乱れた場合、従来の心拍変動解析方法では、非定常性のノイズが混入し、呼吸周波数以外に混入するノイズが大きいという課題もあり、その場合、呼吸の乱れを検出する手段がなかった。そして、算出した心拍変動成分の信頼性が不明であるという課題も一例として挙げられる。
【００３１】
さらに、従来の心拍変動解析方法では、図４（Ａ），（Ｂ）に示すようにローパスフィルタが遮断周波数においても十分な減衰量を得ることができないという課題が一例として挙げられる。
【００３２】
本発明が解決しようとする課題としては、ローパスフィルタの設計が容易で、ノイズの混入が少なく高精度に心拍ゆらぎを計測可能な心拍変動解析装置、心拍変動解析方法、および心拍変動解析用プログラムを提供することが一例として挙げられる。
【００３３】
【課題を解決するための手段】
上記課題を解決するため、請求項１記載の発明は、被検体の心電図情報を検出する心電図情報検出手段と、前記心電図情報から前記被検体の心拍を表す心拍信号を算出する心拍信号算出手段と、前記心拍信号から前記被検体の呼吸状態を反映させた呼吸信号を演算する呼吸信号演算手段と、前記心拍信号と前記呼吸信号とに基づいて前記被検体の呼吸状態を反映させた前記被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供手段と、を備えたことを特徴とする。
【００３４】
上記課題を解決するため、請求項７記載の発明は、被検体の心電図情報を検出する心電図情報検出工程と、前記心電図情報から前記被検体の心拍を表す心拍信号を算出する心拍信号算出工程と、前記心拍信号から前記被検体の呼吸状態を反映させた呼吸信号を演算する呼吸信号演算工程と、前記心拍信号と前記呼吸信号とに基づいて前記被検体の呼吸状態を反映させた前記被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供工程と、を備えたことを特徴とする。
【００３５】
上記課題を解決するため、請求項１３記載の発明は、心拍変動解析装置に含まれるコンピュータを、被検体の心電図情報を検出する心電図情報検出手段、前記心電図情報から前記被検体の心拍を表す心拍信号を算出する心拍信号算出手段、前記心拍信号から前記被検体の呼吸状態を反映させた呼吸信号を演算する呼吸信号演算手段、前記心拍信号と前記呼吸信号とに基づいて前記被検体の呼吸状態を反映させた前記被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供手段、として機能させることを特徴とする。
【００３６】
【発明の実施の形態】
以下、本願の実施形態を図面に基づいて説明する。
【００３７】
図５は本願に係る心拍変動解析装置の第１実施形態の構成を示すブロック図である。なお、本実施形態の心拍変動解析装置は、例えば車両の通信ナビゲーションシステムなどに搭載され、被検体である運転者の心拍ゆらぎを解析し、その解析結果から精神状態や体調に合わせて提供する音楽や走行ルートを選択するために用いられる。この場合、運転者の心電図情報は、ハンドルセンサを通して心電図計測装置により検出される。
【００３８】
図５に示すように、本実施形態の心拍変動解析装置は、被検体の心電図情報を検出する心電図情報検出手段の一例としての心電図計測装置１１と、この心電図計測装置１１の心電図情報から被検体の心拍を表す心拍信号を算出する心拍信号算出手段の一例としてのピーク間隔検出器１２と、このピーク間隔検出器１２により算出された心拍信号から被検体の呼吸状態を反映させた呼吸信号を演算する呼吸信号演算手段の一例としての自己相関演算器１３と、ピーク間隔検出器１２により算出された心拍信号と自己相関演算器１３により演算された呼吸信号とに基づいて被検体の呼吸状態を反映させた被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供手段の一例としての心拍ゆらぎ情報提供装置１４とから大略構成されている。
【００３９】
ピーク間隔検出器１２は、心電図情報を表す心電図波形のピーク間隔データを算出する手段であって、このピーク間隔データは心電図情報のＲ波相互間の間隔（Ｒ−Ｒ間隔）データである。
【００４０】
自己相関演算器１３は、心電図波形のピーク間隔データ、すなわちＲ−Ｒ間隔データの自己相関値を演算するものである。
【００４１】
心拍ゆらぎ情報提供装置１４は、心電図波形のピーク間隔データの自己相関値から被検体の呼吸数を算出する呼吸数算出手段の一例としての呼吸数算出器１５と、ピーク間隔データを、呼吸数に基づいて周波数を変えた変調信号で直交変更する変調手段の一例としての変調器１６と、この変調器１６により変調された変調信号を復調する復調手段の一例としての復調器１７と、この復調器１７による復調結果を心拍ゆらぎ情報として出力する心拍ゆらぎ情報出力手段の一例としての心拍ゆらぎ情報出力装置１８と、を備えて構成されている。
【００４２】
自己相関演算器１３により演算されたＲ−Ｒ間隔データの自己相関値は、信頼度情報出力手段の一例としての信頼度情報出力装置１９に入力され、この信頼度情報出力装置１９は被検体の呼吸状態の信頼度を表す情報として出力し、自己相関値が閾値０．４より大きい場合には、呼吸が安定し信頼度が高いと判断される。
【００４３】
上記ピーク間隔検出器１２、自己相関演算器１３、呼吸数算出器１５、変調器１６、および復調器１７はＣＰＵ２０から構成され、このＣＰＵ２０はメモリ２１との間でデータの授受が実行され、ＣＰＵ２０内の各機器はメモリ２１に予め記憶されているプログラムに従って動作する。
【００４４】
次に、本実施形態の心拍変動解析処理を図６のフローチャートに従って説明する。
【００４５】
まず、心電図計測装置１１により計測した心電図情報を読込み、その心電図情報をメモリ２１に記憶する（ステップＳ１）。次いで、ステップＳ２で心電図波形のピーク間隔データであるＲ−Ｒ間隔データの算出タイミングであるか否かを判断し、Ｒ−Ｒ間隔データの算出タイミングである場合（ステップＳ２；ＹＥＳ）にはステップＳ３に進む一方、Ｒ−Ｒ間隔データの算出タイミングでない場合（ステップＳ２；ＮＯ）にはステップＳ１に戻り、Ｒ−Ｒ間隔データの算出タイミングになるまで待機状態となる。
【００４６】
ステップＳ３では、ピーク間隔検出器１２によりＲ−Ｒ間隔データを算出するとともに、そのＲ−Ｒ間隔データをメモリ２１に記憶する。そして、ピーク間隔検出器１２により算出したＲ−Ｒ間隔データの自己相関値を自己相関演算器１３により求めてメモリ２１に記憶（ステップＳ４）した後、呼吸数算出器１５により呼吸周波数を算出する（ステップＳ５）。ここで、ピーク間隔のデータ（不等間隔データ）をスプライン補間して１０Ｈｚでリサンプリングを行う。
【００４７】
さらに、呼吸周波数に応じて変調信号の周波数を変えた後、その変えた周波数を用いて変調器１６により直交変調処理を実行する（ステップＳ６，Ｓ７）。すなわち、ステップＳ７では、Ｒ−Ｒ間隔データに呼吸周波数の位相の９０度異なる正弦波を掛け合わせ、周波数変調を行う。次のステップＳ８の復調処理では、周波数変調した２つの変調波を掛け合わせ、直流帯域通過、遮断帯域０．３〜０．９Ｈｚのローパスフィルタを通して復調する。そして、心拍ゆらぎ成分を取り出して心拍ゆらぎ出力装置１８から出力する（ステップＳ９）。
【００４８】
次いで、自己相関値を閾値０．４と比較し、その比較結果に応じた出力態様で表示する（ステップＳ１０，Ｓ１１）。すなわち、ステップＳ１０，Ｓ１１では、自己相関値が閾値０．４より大きい場合に呼吸が安定していると判断し、心拍ゆらぎの信頼度が高いと判断し、その旨を信頼度情報出力装置１９から出力表示する。一方、自己相関値が閾値０．４より小さい場合には、呼吸が不安定であると判断し、心拍ゆらぎの信頼度が低いと判断し、その旨を信頼度情報出力装置１９から出力表示する。つまり、信頼度情報出力装置１９では、心拍ゆらぎの信頼度が高い場合と低い場合とで異なるように画像や音声にて出力表示する。
【００４９】
これら一連の処理が終了したら、全体の処理を終了するかを判断し、終了の場合（ステップＳ１２；ＹＥＳ）には、全体の処理を終了する一方、終了しない場合（ステップＳ２；ＮＯ）にはステップＳ１に戻り、上記一連の処理を再度実行する。
【００５０】
次に、本実施形態により心拍ゆらぎの呼吸成分の大きさを求めるための数学的手法について説明する。
【００５１】
まず、呼吸周波数をｆｒ、位相をφ、中間周波数をｆｗとする。そして、心拍ゆらぎを含むＲ−Ｒ間隔の信号を
【００５２】
【数１５】
ｙ＝Ａ×ｓｉｎ（２π×ｆｒ×ｔ＋φ）
とする。ここで、Ａは心拍ゆらぎの呼吸性成分の大きさ、ｆｒは呼吸の周波数、ｔは時間である。
【００５３】
また、中間周波数信号Ｉ，Ｑは、それぞれＲ−Ｒ間隔の信号とローカル信号とを変調器１６において乗算すると、
【００５４】
【数１６】
Ｉ＝１／２×Ａ×（ｃｏｓ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
−ｃｏｓ（２π×ｆｗ×ｔ＋φ＋ｆｒ×２π×ｔ）
【００５５】
【数１７】
Ｑ＝１／２×Ａ×（ｓｉｎ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
−ｓｉｎ（２π×ｆｗ×ｔ＋φ−ｆｒ×２π×ｔ）
呼吸の周波数ｆｒはローカル信号ｆｗと等しいので、
【００５６】
【数１８】
Ｉ＝１／２×Ａ×（φ）−ｃｏｓ（２×２π×ｆｗ×ｔ＋φ）
【００５７】
【数１９】
Ｑ＝１／２×Ａ×（φ）−ｓｉｎ（２×２π×ｆｗ×ｔ＋φ）
２×ｆｗは、０．３〜０．９Ｈｚである。したがって、本実施形態では、直流成分が通過し、遮断帯域０．３〜０．９Ｈｚのローパスフィルタを用いればよいことになる。
【００５８】
中間周波数信号Ｉ，Ｑに直流帯域通過のローパスフィルタを通すことにより式の第２項が消去でき、それぞれＩＸ，ＱＸとすると、
【００５９】
【数２０】
ＩＸ＝１／２×Ａ×（ｃｏｓ（２π×（ｆｗ−ｆｒ）×ｔ＋φ）
【００６０】
【数２１】
ＱＸ＝１／２×Ａ×（ｓｉｎ（２π×（ｆｗ−ｆｒ）×ｔ＋φ）
三角関数の公式より
【００６１】
【数２２】
ｓｉｎ（ｘ）×ｓｉｎ（ｘ）＋ｃｏｓ（ｘ）×ｃｏｓ（ｘ）＝１
であるので、
【００６２】
【数２３】
（ＩＸ）２＋（ＱＸ）２＝（１／２×Ａ）２
故に
【００６３】
【数２４】
Ａ＝２×（（ＩＸ）２＋（ＱＸ）２）０．５
である。
【００６４】
すなわち、図７（Ａ）〜図７（Ｅ）は本実施形態の信号を示す波形図である。入力信号として振幅１、周波数０．２Ｈｚの正弦波を心拍変動とする。図７（Ａ）は入力信号、図７（Ｂ）は０．２Ｈｚのローカル信号、図７（Ｃ）は上記入力信号と上記ローカル信号を掛け合わせた中間周波数信号（Ｉ，Ｑ）である。
【００６５】
また、図７（Ｄ）は中間周波数信号（Ｉ，Ｑ）をローパスフィルタに通した信号（ＩＸ，ＱＸ）、図７（Ｅ）はローパスフィルタに通した信号（ＩＸ，ＱＸ）より計算した振幅成分である。したがって、本実施形態では、図７（Ｅ）に示すように計算した振幅成分にノイズが重畳していないことが判る。
【００６６】
図８は本実施形態の信号と折り返しノイズとの関係を示す説明図である。図８に示すように、直流信号成分は０〜０．０５Ｈｚの範囲であり、折り返しノイズは０．３〜０．９Ｈｚの範囲である。この折り返しノイズを除去する復調器１７のローパスフィルタは、図８の破線に示す特性となる。これにより、本実施形態では、通過帯域幅が０．０５Ｈｚと狭く、設計の容易なローパスフィルタを用いることができる。
【００６７】
図９（Ａ），（Ｂ）は本実施形態で用いるローパスフィルタの特性を示す波形図である。本実施形態では、図９（Ａ）に示すように遮断周波数においても十分な減衰効果が得られる。
【００６８】
次に、自己相関による呼吸数と相関値の算出方法を説明する。図１０は心拍変動の信号を示す波形図である。図１０は、１０Ｈｚサンプリングで６７ポイントのデータである。この入力信号の自己相関係数を調べると図１１になる。図１１に示すように遅延時間３１ポイントで自己相関係数が０．４５である。すなわち、呼吸の周期３．１秒、周波数は１／３．１＝０．３２Ｈｚであり、自己相関係数が閾値０．４より大きいので、呼吸が安定していることになる。
【００６９】
このように本実施形態によれば、被検体の心電図情報を検出する心電図計測装置１１と、この心電図計測装置１１の心電図情報から被検体の心拍を表す心拍信号を算出するピーク間隔検出器１２と、このピーク間隔検出器１２により算出された心拍信号と呼吸信号とに基づいて被検体の呼吸状態を反映させた被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供装置１４とを備えたことにより、ローパスフィルタの設計が容易で、ノイズの混入が少なく、呼吸の計測が不要で、高精度に心拍ゆらぎを計測することができる。
【００７０】
また、本実施形態によれば、ピーク間隔検出器１２は、心電図情報を表す心電図波形のピーク間隔データを算出する手段であって、このピーク間隔データは心電図情報のＲ波相互間の間隔（Ｒ−Ｒ間隔）データであることにより、高精度に心拍ゆらぎを計測することができる。
【００７１】
さらに、本実施形態によれば、心拍ゆらぎ情報提供装置１４は、心電図波形のピーク間隔データの自己相関値から被検体の呼吸数を算出する呼吸数算出器１５と、ピーク間隔データを、呼吸数に基づいて周波数を変えた変調信号で直交変更する変調器１６と、この変調器１６により変調された変調信号を復調する復調器１７と、この復調器１７による復調結果を心拍ゆらぎ情報として出力する心拍ゆらぎ情報出力装置１８とを備えたことにより、復調器１７に用いるローパスフィルタの設計が容易で、ノイズの混入が少なく、呼吸の計測が不要で、高精度に心拍ゆらぎを計測することができる。
【００７２】
さらにまた、本実施形態によれば、自己相関演算器１３は、心電図波形のピーク間隔データ、すなわちＲ−Ｒ間隔データの自己相関値を演算し、この自己相関値は、信頼度情報出力装置１９に入力され、この信頼度情報出力装置１９は被検体の呼吸状態の信頼度を表す情報として出力することにより、計測した心拍ゆらぎの信頼度を計測することができる。
【００７３】
図１２は本願に係る心拍変動解析装置の第２実施形態の構成を示すブロック図である。なお、前記第１実施形態と同一の部分には、同一の符号を付して説明を省略する。
【００７４】
前記第１実施形態は、自己相関関数を用いて呼吸周波数および相関値を算出したが、本実施形態は、他の周波数解析方法であるＦＦＴ（高速フーリェ変換）２２にてスペクトルを調べて、ピーク周波数を算出して中間周波数としている。
【００７５】
このように、ＦＦＴ２２にてスペクトルを分析し、ピーク周波数を算出して中間周波数としても前記第１実施形態と同様の作用および効果が得られる。なお、ＦＦＴ２２に代えてＭＥＭ（最大エントロピー法），ＡＲ（自己回帰）モデルを用いてスペクトルを分析すれば、ＦＦＴよりも短いデータから分解能の高いスムーズなスペクトルが得られるなどの効果を奏する。
【００７６】
なお、上記各実施形態は、本発明を車両の通信ナビゲーションシステムなどに搭載して被検体である運転者の心拍ゆらぎを解析するために用いたが、これに限らず、例えばスポーツジムにおける各種運動器具に搭載して心拍ゆらぎを解析するようにしてもよい。
【００７７】
また、図６に示したフローチャートに対応するプログラムを、フレキシブルディスクまたはハードディスクなどの情報記録媒体に記憶させておき、これを汎用のマイクロコンピュータなどに読み出して実行させることで、当該マイクロコンピュータを上記実施形態に係るＣＰＵ２０として機能させることも可能である。
【図面の簡単な説明】
【図１】従来の心拍変動解析装置の構成を示すブロック図である。
【図２】（Ａ）〜（Ｅ）は従来の信号を示す波形図である。
【図３】従来の信号と折り返しノイズとの関係を示す説明図である。
【図４】（Ａ），（Ｂ）は従来例で用いるローパスフィルタの特性を示す波形図である。
【図５】本願に係る心拍変動解析装置の第１実施形態の構成を示すブロック図である。
【図６】第１実施形態の心拍変動解析処理を示すフローチャートである。
【図７】（Ａ）〜（Ｅ）は第１実施形態の信号を示す波形図である。
【図８】第１実施形態の信号と折り返しノイズとの関係を示す説明図である。
【図９】（Ａ），（Ｂ）は第１実施形態で用いるローパスフィルタの特性を示す波形図である。
【図１０】第１実施形態の心拍変動の信号を示す波形図である。
【図１１】第１実施形態の入力信号の自己相関係数を示す波形図である。
【図１２】本願に係る心拍変動解析装置の第２実施形態の構成を示すブロック図である。
【符号の説明】
１１ 心電図計測装置
１２ ピーク間隔検出器
１３ 自己相関演算器
１４ 心拍ゆらぎ情報提供装置
１５ 呼吸数算出器
１６ 変調器
１７ 復調器
１８ 心拍ゆらぎ情報出力装置
１９ 信頼度情報出力装置
２０ ＣＰＵ
２１ メモリ
２２ ＦＦＴ
Claims (18)
- 被検体の心電図情報を検出する心電図情報検出手段と、
前記心電図情報から前記被検体の心拍を表す心拍信号を算出する心拍信号算出手段と、
前記心拍信号から前記被検体の呼吸状態を反映させた呼吸信号を演算する呼吸信号演算手段と、
前記心拍信号と前記呼吸信号とに基づいて前記被検体の呼吸状態を反映させた前記被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供手段と、
を備えたことを特徴とする心拍変動解析装置。 - 請求項１記載の心拍変動解析装置において、
前記心拍信号算出手段は、前記心電図情報を表す心電図波形のピーク間隔データを算出する手段であることを特徴とする心拍変動解析装置。 - 請求項２記載の心拍変動解析装置において、
前記ピーク間隔データは、前記心電図情報のＲ波相互間の間隔データであることを特徴とする心拍変動解析装置。 - 請求項１記載の心拍変動解析装置において、
前記心拍ゆらぎ情報提供手段は、前記心電図波形のピーク間隔データの自己相関値から前記被検体の呼吸数を算出する呼吸数算出手段と、
前記ピーク間隔データを、前記呼吸数に基づいて周波数を変えた変調信号で直交変更する変調手段と、
前記変調手段により変調された変調信号を復調する復調手段と、
前記復調手段による復調結果を心拍ゆらぎ情報として出力する心拍ゆらぎ情報出力手段と、
を備えたことを特徴とする心拍変動解析装置。 - 請求項１記載の心拍変動解析装置において、
前記呼吸信号演算手段は、前記心電図波形のピーク間隔データの自己相関値を演算する自己相関値演算手段であることを特徴とする心拍変動解析装置。 - 請求項５記載の心拍変動解析装置において、
前記自己相関値を、前記被検体の呼吸状態の信頼度を表す情報として出力する信頼度情報出力手段を備えることを特徴とする心拍変動解析装置。 - 被検体の心電図情報を検出する心電図情報検出工程と、
前記心電図情報から前記被検体の心拍を表す心拍信号を算出する心拍信号算出工程と、
前記心拍信号から前記被検体の呼吸状態を反映させた呼吸信号を演算する呼吸信号演算工程と、
前記心拍信号と前記呼吸信号とに基づいて前記被検体の呼吸状態を反映させた前記被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供工程と、
を備えたことを特徴とする心拍変動解析方法。 - 請求項７記載の心拍変動解析方法において、
前記心拍信号算出方法は、前記心電図情報を表す心電図波形のピーク間隔データを算出する工程であることを特徴とする心拍変動解析方法。 - 請求項８記載の心拍変動解析方法において、
前記ピーク間隔データは、前記心電図情報のＲ波相互間の間隔データであることを特徴とする心拍変動解析方法。 - 請求項７記載の心拍変動解析方法において、
前記心拍ゆらぎ情報提供工程は、前記心電図波形のピーク間隔データの自己相関値から前記被検体の呼吸数を算出する呼吸数算出工程と、
前記ピーク間隔データを、前記呼吸数に基づいて周波数を変えた変調信号で直交変更する変調工程と、
前記変調手段により変調された変調信号を復調する復調工程と、
前記復調手段による復調結果を心拍ゆらぎ情報として出力する心拍ゆらぎ情報出力工程と、
を備えたことを特徴とする心拍変動解析方法。 - 請求項７記載の心拍変動解析方法において、
前記呼吸信号演算工程は、前記心電図波形のピーク間隔データの自己相関値を演算する自己相関値演算工程であることを特徴とする心拍変動解析方法。 - 請求項１１記載の心拍変動解析方法において、
前記自己相関値を、前記被検体の呼吸状態の信頼度を表す情報として出力する信頼度情報出力工程を備えることを特徴とする心拍変動解析方法。 - 心拍変動解析装置に含まれるコンピュータを、
被検体の心電図情報を検出する心電図情報検出手段、
前記心電図情報から前記被検体の心拍を表す心拍信号を算出する心拍信号算出手段、
前記心拍信号から前記被検体の呼吸状態を反映させた呼吸信号を演算する呼吸信号演算手段、
前記心拍信号と前記呼吸信号とに基づいて前記被検体の呼吸状態を反映させた前記被検体の心拍のゆらぎ情報を提供する心拍ゆらぎ情報提供手段、
として機能させることを特徴とする心拍変動解析用プログラム。 - 請求項１３記載の心拍変動解析用プログラムにおいて、
前記心拍信号算出手段は、前記心電図情報を表す心電図波形のピーク間隔データを算出する手段として機能させることを特徴とする心拍変動解析用プログラム。 - 請求項１４記載の心拍変動解析用プログラムにおいて、
前記ピーク間隔データは、前記心電図情報のＲ波相互間の間隔データとして機能させることを特徴とする心拍変動解析用プログラム。 - 請求項１３記載の心拍変動解析用プログラムにおいて、
前記心拍ゆらぎ情報提供手段は、前記心電図波形のピーク間隔データの自己相関値から前記被検体の呼吸数を算出する呼吸数算出手段、
前記ピーク間隔データを、前記呼吸数に基づいて周波数を変えた変調信号で直交変更する変調手段、
前記変調手段により変調された変調信号を復調する復調手段、
前記復調手段による復調結果を心拍ゆらぎ情報として出力する心拍ゆらぎ情報出力手段、
として機能させることを特徴とする心拍変動解析用プログラム。 - 請求項１３記載の心拍変動解析用プログラムにおいて、
前記呼吸信号演算手段は、前記心電図波形のピーク間隔データの自己相関値を演算する自己相関値演算手段として機能させることを特徴とする心拍変動解析用プログラム。 - 請求項１７記載の心拍変動解析用プログラムにおいて、
前記自己相関値を、前記被検体の呼吸状態の信頼度を表す情報として出力する信頼度情報出力手段として機能させることを特徴とする心拍変動解析用プログラム。
Priority Applications (4)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2002246635A JP3946108B2 (ja) | 2002-08-27 | 2002-08-27 | 心拍変動解析装置、心拍変動解析方法、および心拍変動解析用プログラム |
| EP03254954A EP1394687B1 (en) | 2002-08-27 | 2003-08-08 | Apparatus and method for analyzing heart-rate variability based on electrocardiogram information |
| DE60306856T DE60306856T2 (de) | 2002-08-27 | 2003-08-08 | Vorrichtung und Verfahren zur Analyse von Herzfrequenzschwankungen auf der Grundlage von Elektrokardiogramminformationen |
| US10/647,589 US7104959B2 (en) | 2002-08-27 | 2003-08-26 | Apparatus and method for analyzing heart-rate variability based on electrocardiogram information |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2002246635A JP3946108B2 (ja) | 2002-08-27 | 2002-08-27 | 心拍変動解析装置、心拍変動解析方法、および心拍変動解析用プログラム |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| JP2004081513A true JP2004081513A (ja) | 2004-03-18 |
| JP3946108B2 JP3946108B2 (ja) | 2007-07-18 |
Family
ID=31492536
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP2002246635A Expired - Fee Related JP3946108B2 (ja) | 2002-08-27 | 2002-08-27 | 心拍変動解析装置、心拍変動解析方法、および心拍変動解析用プログラム |
Country Status (4)
| Country | Link |
|---|---|
| US (1) | US7104959B2 (ja) |
| EP (1) | EP1394687B1 (ja) |
| JP (1) | JP3946108B2 (ja) |
| DE (1) | DE60306856T2 (ja) |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2010085100A (ja) * | 2008-09-29 | 2010-04-15 | Toto Ltd | 人体検知装置及びそれを備えた小便器 |
| JP2010214016A (ja) * | 2009-03-18 | 2010-09-30 | Iwate Medical Univ | 体表面心電図を解析し、ｔ波交互脈または心房細動波に関する２次元機能図を生成する心電図解析装置 |
| JP2015042267A (ja) * | 2008-04-16 | 2015-03-05 | コーニンクレッカ フィリップス エヌ ヴェ | 睡眠／覚醒状態評価方法及びシステム |
| WO2022071342A1 (ja) | 2020-09-30 | 2022-04-07 | ダイキン工業株式会社 | 温度推定装置、空調制御装置、空調制御システム |
Families Citing this family (17)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| DE102006053728A1 (de) * | 2006-08-28 | 2008-03-06 | Biotronik Crm Patent Ag | Bestimmung von Schwellenwerten für RR-Schlagfilter |
| US20080161651A1 (en) * | 2006-12-27 | 2008-07-03 | Cardiac Pacemakers, Inc. | Surrogate measure of patient compliance |
| TW201121501A (en) * | 2009-12-29 | 2011-07-01 | V Tac Technology Co Ltd | Examination method of exercise heart rate. |
| JP5709017B2 (ja) * | 2010-02-15 | 2015-04-30 | 国立大学法人九州大学 | 被験体状態解析用信号のピーク周波数測定システム |
| FI126631B (en) * | 2014-07-28 | 2017-03-15 | Murata Manufacturing Co | Method and apparatus for monitoring stress |
| AT516204B1 (de) | 2014-08-29 | 2016-07-15 | Peter R Mmag Msc Hauschild | Verfahren und Anordnung zur Analyse der Interaktion von elektromagnetischen Hochfrequenz-Immissionen mit vegetativen Regulationsmechanismen eines Testsubjekts |
| DE112015005804T5 (de) * | 2014-12-24 | 2017-10-19 | Asahi Kasei Kabushiki Kaisha | Atemzustandsschätzvorrichtung, tragbare Vorrichtung, an Körper tragbare Vorrichtung, Programm, Medium, Atemzustandsschätzverfahren und Atemzustandsschätzer |
| JP6687645B2 (ja) * | 2016-02-04 | 2020-04-28 | 日本電信電話株式会社 | 生体信号処理方法および生体信号処理装置 |
| CN105982659A (zh) * | 2016-03-17 | 2016-10-05 | 厦门鑫奥力电器有限公司 | 一种用于测量心率的耳夹及其控制方法 |
| AT16253U1 (de) | 2017-09-04 | 2019-05-15 | Dr Peter R Hauschild Mmag | Messgeräte-/Ladestations-Anordnung mit Akkulade- und Datenübertragungsschnittstelle |
| AT520401B1 (de) | 2017-09-04 | 2021-04-15 | R Hauschild Mmag Dr Peter | Portables Messgerät zum nichtinvasiven Erfassen bioelektrischer Signale |
| CN107928660B (zh) * | 2017-10-25 | 2022-07-26 | 厦门传芯健康科技有限公司 | 一种基于心电图的健康预警系统 |
| AT523903B1 (de) | 2020-06-12 | 2023-02-15 | Matissek Msc Michael | Umwelt-/Gesundheitsmonitoringsystem zur Erfassung von elektromagnetischen Feldern und der Beeinflussung vegetativer Regulationsmechanismen des menschlichen Organismus |
| US11723568B2 (en) * | 2020-09-10 | 2023-08-15 | Frictionless Systems, LLC | Mental state monitoring system |
| CN113747391B (zh) * | 2021-11-05 | 2022-03-18 | 比亚迪股份有限公司 | 监测轨道车辆驾驶员健康的方法和系统、管理中心和设备 |
| CN114403892B (zh) * | 2022-01-24 | 2024-04-02 | 中山大学 | 基于心率变异性交互性调节呼吸反馈的方法、装置及设备 |
| CN117100237B (zh) * | 2023-10-18 | 2024-02-06 | 众保健康科技服务(济南)有限公司 | 一种居家养老用智能监测系统 |
Family Cites Families (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US5042497A (en) * | 1990-01-30 | 1991-08-27 | Cardiac Pacemakers, Inc. | Arrhythmia prediction and prevention for implanted devices |
| US5902250A (en) * | 1997-03-31 | 1999-05-11 | President And Fellows Of Harvard College | Home-based system and method for monitoring sleep state and assessing cardiorespiratory risk |
| US6480733B1 (en) * | 1999-11-10 | 2002-11-12 | Pacesetter, Inc. | Method for monitoring heart failure |
| WO2003073930A1 (en) * | 2002-03-01 | 2003-09-12 | Christine Ross | Novel utilization of heart rate variability in animals |
| TWI245622B (en) * | 2003-12-26 | 2005-12-21 | We Gene Technologies Inc | Auto analyzing method and device thereof for autonomic nervous |
-
2002
- 2002-08-27 JP JP2002246635A patent/JP3946108B2/ja not_active Expired - Fee Related
-
2003
- 2003-08-08 DE DE60306856T patent/DE60306856T2/de not_active Expired - Fee Related
- 2003-08-08 EP EP03254954A patent/EP1394687B1/en not_active Expired - Fee Related
- 2003-08-26 US US10/647,589 patent/US7104959B2/en not_active Expired - Fee Related
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2015042267A (ja) * | 2008-04-16 | 2015-03-05 | コーニンクレッカ フィリップス エヌ ヴェ | 睡眠／覚醒状態評価方法及びシステム |
| JP2010085100A (ja) * | 2008-09-29 | 2010-04-15 | Toto Ltd | 人体検知装置及びそれを備えた小便器 |
| JP2010214016A (ja) * | 2009-03-18 | 2010-09-30 | Iwate Medical Univ | 体表面心電図を解析し、ｔ波交互脈または心房細動波に関する２次元機能図を生成する心電図解析装置 |
| WO2022071342A1 (ja) | 2020-09-30 | 2022-04-07 | ダイキン工業株式会社 | 温度推定装置、空調制御装置、空調制御システム |
Also Published As
| Publication number | Publication date |
|---|---|
| DE60306856D1 (de) | 2006-08-31 |
| EP1394687B1 (en) | 2006-07-19 |
| US7104959B2 (en) | 2006-09-12 |
| JP3946108B2 (ja) | 2007-07-18 |
| EP1394687A2 (en) | 2004-03-03 |
| DE60306856T2 (de) | 2006-11-23 |
| EP1394687A3 (en) | 2005-05-25 |
| US20040092835A1 (en) | 2004-05-13 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| JP3946108B2 (ja) | 心拍変動解析装置、心拍変動解析方法、および心拍変動解析用プログラム | |
| US11622693B2 (en) | Method and apparatus for non-contact fast vital sign acquisition based on radar signal | |
| JP6552013B2 (ja) | マルチレートｅｃｇ処理を使用するｒ−ｒ間隔の測定 | |
| US20160106328A1 (en) | Determining arterial pulse transit time from time-series signals obtained at proximal and distal arterial sites | |
| US8655436B2 (en) | Heart rate meter and heart beat detecting method | |
| JP6056389B2 (ja) | 心拍推定装置、心拍推定方法及びプログラム | |
| CN108113706A (zh) | 一种基于音频信号的心率监测方法、装置及系统 | |
| JP2004121625A (ja) | 脈波検出装置及びフーリエ変換処理装置 | |
| EP3678553B1 (en) | Diagnosis of pathologies using infrasonic signatures | |
| US20190038164A1 (en) | Biological signal processing method and biological signal processing apparatus | |
| JPH0389174A (ja) | 電子測定装置及び周波数推定方法 | |
| Wu et al. | A new approach for analysis of heart rate variability and QT variability in long-term ECG recording | |
| Rong et al. | Smart homes: See multiple heartbeats through wall using wireless signals | |
| Wang et al. | IQ signal demodulation for noncontact vital sign monitoring using a CW Doppler radar: A review | |
| JPH0548128B2 (ja) | ||
| Kong et al. | Heart rate tracking using a wearable photoplethysmographic sensor during treadmill exercise | |
| JPH11128185A (ja) | 心拍変動解析方法および装置 | |
| JP4765100B2 (ja) | 位相測定装置、方法、プログラムおよび記録媒体 | |
| JP2016182165A (ja) | 生体信号処理装置、生体信号処理プログラム、生体信号処理プログラムを記録したコンピュータ読み取り可能な記録媒体、生体信号処理方法 | |
| JP2723337B2 (ja) | 呼吸計測装置 | |
| JP2000296118A (ja) | 生体信号解析方法および装置 | |
| EP2938247B1 (en) | Method and apparatus for reducing motion artifacts in ecg signals | |
| US11911136B2 (en) | System and method for calculating cardiac pulse transit or arrival time information | |
| JP6152806B2 (ja) | 生体情報測定方法 | |
| KR20210110724A (ko) | 힐버트 변환을 이용한 심전도 신호의 크기 측정 장치, 방법 및 컴퓨터로 독출 가능한 기록 매체 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20050725 |
|
| A977 | Report on retrieval |
Free format text: JAPANESE INTERMEDIATE CODE: A971007 Effective date: 20070130 |
|
| A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20070206 |
|
| A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20070320 |
|
| 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: 20070410 |
|
| A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20070410 |
|
| R150 | Certificate of patent or registration of utility model |
Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20110420 Year of fee payment: 4 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20120420 Year of fee payment: 5 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20130420 Year of fee payment: 6 |
|
| LAPS | Cancellation because of no payment of annual fees |