JP2022071812A - 生体信号分析装置、コンピュータプログラム及び記録媒体 - Google Patents

生体信号分析装置、コンピュータプログラム及び記録媒体 Download PDF

Info

Publication number
JP2022071812A
JP2022071812A JP2021107109A JP2021107109A JP2022071812A JP 2022071812 A JP2022071812 A JP 2022071812A JP 2021107109 A JP2021107109 A JP 2021107109A JP 2021107109 A JP2021107109 A JP 2021107109A JP 2022071812 A JP2022071812 A JP 2022071812A
Authority
JP
Japan
Prior art keywords
waveform
time
measurement
frequency
data
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.)
Pending
Application number
JP2021107109A
Other languages
English (en)
Inventor
悦則 藤田
Yoshinori Fujita
由美 小倉
Yumi Ogura
慎一郎 前田
Shinichiro Maeda
重行 小島
Shigeyuki Kojima
正博 堀川
Masahiro Horikawa
良香 延廣
Yoshika Nobuhiro
成彦 金子
Shigehiko Kaneko
正生 吉栖
Masao Yoshizumi
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.)
Delta Tooling Co Ltd
Original Assignee
Delta Tooling Co 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 Delta Tooling Co Ltd filed Critical Delta Tooling Co Ltd
Priority to US18/250,827 priority Critical patent/US20230389814A1/en
Priority to PCT/JP2021/039924 priority patent/WO2022092243A1/ja
Priority to EP21886360.3A priority patent/EP4238491A1/en
Publication of JP2022071812A publication Critical patent/JP2022071812A/ja
Pending legal-status Critical Current

Links

Images

Abstract

【課題】 生体状態、特に、心臓の動きに起因する状態を把握できるようにする。【解決手段】 本発明によれば、心尖拍動と心音との境界周波数を求めることができる。これにより、生体音や体内振動の集合体である生体信号から、心音はもとより、心尖拍動を知ることが可能となった。この境界周波数が心拍数の2次関数で表し、それに基づく相関データを利用して、心拍数又は境界周波数をリアルタイムで知ることができるようになった。その結果、例えば、心拍数がわかれば、境界周波数の大小を知ることができ、心臓に負荷がかかっているかどうか等、従来と比較してより広範な観点からの健康状態の推定を行いやすくすることができる。【選択図】 図12

Description

本発明は、生体状態の推定に関する技術に関し、特に、心尖拍動に関連する指標を用いた生体信号分析装置、コンピュータープログラム及び記録媒体に関する。
心尖拍動は、左室心尖部の運動を反映し、心臓、主として左室機能の評価に用いられる。心尖拍動は、心機図(Mechanocardiogram:MCG)法により記録される心尖拍動図(Apex cardiogram: ACG)により客観化されるのが一般であるが、心尖拍動図は、心臓と胸壁の間の種々の組織を介して胸壁の動きを捉えているため、記録巧拙の影響を受けやすい。また、心機図記録装置は比較的大型のものが多く、測定にあたっても防音室で行わなければならないなどの制約がある。
この点に鑑み特許文献1では、体表面から心尖部近傍に向けて超音波信号を入射させる超音波発信部と、反射超音波を受信する受信センサを備えた心尖拍動図用の生体情報測定装置を開示しており、防音室等の特殊な環境でなくても容易に心尖拍動図を作成できるとしている。
一方、本発明者らは、特許文献2~5等において、人の背部の体表面に生じる振動を非拘束で捉え、その振動を解析して人の状態を推定する技術を提案している。人の背部の体表面に生じる振動は、心臓と大動脈等の生体内の振動が伝播したものであり、心房及び心室の収縮期及び拡張期の情報や、循環の補助ポンプとなる血管壁の弾力情報及び反射波の情報を含んでいる。
特許文献2では、体表面を介して伝播する振動(生体信号)から抽出した1Hz近傍の背部体表脈波の時系列波形に所定の時間幅を適用してスライド計算を行って周波数傾きの時系列波形を求め、その変化の傾向から、例えば、振幅が増幅傾向にあるか、減衰傾向にあるかなどによって生体状態の推定を行っている。また、生体信号を周波数解析し、予め定めたULF帯域(極低周波帯域)からVLF帯域(超低周波帯域)に属する機能調整信号、疲労受容信号及び活動調整信号に相当する各周波数のパワースペクトルを求め、各パワースペクトルの時系列変化から人の状態を判定することも開示している。
特許文献3~4では、恒常性維持機能レベルを判定する手段を開示している。また、特許文献5では、生体信号の音・振動情報に対応した固有振動数を含む固有振動子を備えた共鳴層を具備する音・振動情報収集機構を開示している。
また、非特許文献1には、心音が不規則振動で構成されることが示されている。非特許文献2には、上記の心尖拍動図の波形を分類して診断に利用することが開示されている。非特許文献3には、心周期のゆらぎに関する説明がある。
WO2012/165660号公報 特開2011-167362号公報 特開2014-117425号公報 特開2014-223271号公報 特開2016-26516号公報
Nobuhiro, Y. et al. Development of Physical Condition Fluctuation Prediction Model Using Trunk Biosignals, Proceeding of Asia-Pacific Vibration Conference 2019 (2019) Yoshikawa,J.et al. Physical Examination in Cardiology, Bunkodo Shoten Co.,Ltd.,79-95 (2005) Kobayashi M., Musha T.: 1/f fluctuation of heart beat period. IEEE Transactions on biomedical engineering. BME, 29, 45-457 (1982)
特許文献1に開示のものは、192個の超音波振動子を用いた装置であり、12個を一組としてそれらを順次ずらしてスキャンして心尖拍動図を作成する。また、より高い精度の心尖拍動図を作成するために、192個の超音波振動子を一列とし、それを複数列用いた装置も開示している。特許文献1では、超音波振動子をこのように多数使用しなければ所望の心尖拍動図を得ることは難しい。よって、製品としては非常に高価なものにならざるを得ず、実用的ではない。また、鮮明な心尖拍動図を作成するためには、肋骨の位置の検出も考慮する必要があり、発信周波数も限定され、また、受信素子の位置も限られ、それらの条件を満たす装置としなければならず、その点も、構造の複雑化、製造コストの増加につながる。
また、特許文献1では、従来の心機図記録装置が複数人数で取り扱わなければならなかったのに対し、一人でも測定できるということを利点として述べている。しかし、超音波振動子をスキャンするために、医師や技師の操作が必要であり、測定に手間がかかると共に、スキャンしている間のデータしか得られない。すなわち、測定対象者の心尖拍動を長時間に亘って捉えることには適していない。
一方、特許文献2~5に開示した体表面に生じる振動を非拘束で捉えるセンサは、三次元立体編物、三次元立体編物を取り囲むフィルム、マイクロフォン等からなり、人の体に接触させておくだけで体表面を介した生体信号を取得できる。すなわち、人の体に取り付ければ、医師等が何らの操作をしなくても生体信号データが得られる。しかし、特許文献2~5は、主に自律神経機能や心拍変動の解析を行っており、左室心尖部の運動である心尖拍動に関する解析は行われていない。
また、特許文献2~5で用いたセンサは、人の背部の体表面から伝播される20~100Hzの音響振動を受動的に捉えるのに適しており、10Hz以下の体内振動である心尖拍動を抽出するのは適していない。
また、非特許文献3の心尖拍動図を利用して診断を行うには、心尖拍動図の波形を読み取る専門的知識が必要となると共に、判定者によって判定結果に差が生じる場合もあり、判定精度や迅速性の点で改善の余地がある。従って、心尖拍動図を利用した診断はあくまで補助的に利用されているに過ぎない。さらに、心尖拍動に関連する情報として、このような心尖拍動図は、補助的ながらも診断に利用されているが、心尖拍動図以外の心尖拍動に関連する指標を用いて人の様々な健康状態(心疾患等の疾病の有無、疾病の特定、体調の良否等)を推定することは行われていない。
本発明は上記に鑑みなされたものであり、簡易で安価に製造可能であると共に、長時間に亘る測定も可能で、生体信号検出センサにより収集される生体信号を分析して、心尖拍動に関する情報を取り出すことを可能とし、さらには、人の様々な健康状態を推定することもできる生体信号分析装置、コンピュータプログラム及び記録媒体を提供することを課題とする。
上記課題を解決するため、本発明者らは、三次元立体編物とマイクロフォン等を用いたもので、0.5~80Hzの生体信号データを取得できる生体信号検出センサを開発した。一方、生体信号検出センサは、体表面から伝播される音響振動を受動的に捉えるだけであり、捕捉される生体信号データには、様々な生体音や体内振動が含まれている。収集可能な周波数帯域を低周波帯域まで広げた結果、心尖拍動に関するデータが含まれるようになったものの、心音のデータに隠れ、心尖拍動のデータのみを抽出することが困難であるという新たな課題に直面した。そこで、本発明者らはかかる課題を解決すべく鋭意研究し、本発明を完成するに至った。
すなわち、本発明の生体信号分析装置は、
体表面を介して生体信号検出センサにより得られる生体信号データを周波数解析する周波数解析手段と、
前記周波数解析手段から得られる前記生体信号データの周波数解析結果から、前記生体信号データ中、心尖拍動により生じる振動と心音により生じる振動との境界周波数を求める境界周波数特定手段と
を有することを特徴とする。
前記境界周波数特定手段は、
前記周波数解析結果中、調和振動と不規則振動との境界となるパワースペクトルの急変部を求め、この急変部を基準に前記境界周波数を特定する手段を含むことが好ましい。
前記境界周波数特定手段は、
前記周波数解析結果に、同時に測定した心音データの周波数解析結果を加味して、前記パワースペクトルの急変部を求める手段を含むことが好ましい。
前記境界周波数特定手段は、
前記生体信号データ及び前記心音データの各周波数解析結果をそれぞれ加算平均処理した波形を対数差分法を用いて両対数軸表示する両対数軸表示手段と、両対数軸表示した波形からゆらぎの変化点を求め、このゆらぎの変化点を前記パワースペクトルの急変部として特定する急変部を特定する手段と
を有することが好ましい。
前記周波数解析手段として、短時間フーリエ変換が適用され、
前記境界周波数特定手段は、前記短時間フーリエ変換の解析結果から前記パワースペクトルの急変部を求める手段を含むことが好ましい。
前記周波数解析手段は、前記短時間フーリエ変換の解析結果を、時間、周波数及びパワースペクトルの変動の程度を示す画像データに出力する手段を含み、
前記境界周波数特定手段は、前記画像データから、前記パワースペクトルの急変部を求める手段を含むことが好ましい。
前記境界周波数に関連する相関データが記憶された相関データ記憶部と、
前記相関データに照合して、測定対象者の測定時の健康状態を推定する測定時状態推定手段と
をさらに有することが好ましい。
前記相関データが、前記境界周波数と心拍数との相関データであり、
前記測定時状態推定手段は、測定対象者の測定時の心拍数を、前記境界周波数と心拍数との相関データに照合して、前記測定対象者の測定時の前記境界周波数を推定する手段を含むことが好ましい。
前記相関データが、前記境界周波数と心拍変動のゆらぎ特性との相関データであり、
前記測定時状態推定手段は、前記測定対象者の測定時の前記境界周波数を、前記境界周波数と心拍変動のゆらぎ特性との相関データに照合して、前記測定対象者の測定時の心拍変動のゆらぎ特性を推定する手段を含むことが好ましい。
さらに、前記生体信号データを、前記境界周波数特定手段により特定される前記境界周波数を上限値に設定してフィルタリングし、前記心尖拍動により生じる振動の波形を求める心尖拍動波形抽出手段と、
前記心尖拍動により生じる振動の波形を分類する波形分類手段と
を有すると共に、
前記測定時状態推定手段は、前記心尖拍動波形抽出手段により求められる、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形を、前記波形分類手段により分類し、当該測定時の波形分類結果のデータをもとに、前記測定対象者の測定時の健康状態を推定する手段を含むことが好ましい。
前記相関データ記憶部には、前記波形分類結果と健康状態との相関データが予め記憶されており、
前記測定時状態推定手段は、前記心尖拍動波形抽出手段により求められる、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形を、前記波形分類手段により分類し、当該測定時の波形分類結果を、予め記憶された前記波形分類結果と健康状態との相関データに照合して、前記測定対象者の測定時の健康状態を推定する手段を含むことが好ましい。
前記心尖拍動波形抽出手段により出力され、前記波形分類手段における分類対象となる波形情報と前記健康状態とを教師データとして用い、前記波形情報から前記健康状態を推定する推定モデルを機械学習により生成するモデル作成手段を有し、
前記測定時状態推定手段は、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形情報を入力として、前記モデル作成手段により作成された推定モデルを用い、取得した測定時の波形情報から、測定時の前記健康状態の値を出力する手段を含むことが好ましい。
前記波形分類手段は、次の(A)及び(B)の手段:
(A)フーリエ級数展開を用いた数学的アプローチにより前記波形を分類する手段、
(B)次の(1)~(4)の工学的アプローチにより得られるデータを1又は2以上組み合わせて前記波形を分類する手段、
(1)前記心尖拍動により生じる振動の波形の所定時間範囲における時間波形の形状を示すデータ、
(2)前記時間波形の前記周波数解析手段により得られる周波数とパワースペクトルを横軸と縦軸にとったグラフのデータ、及び、
(3)前記時間波形の前記短時間フーリエ変換の解析結果である時間、周波数及びパワースペクトルの変動の程度を示す所定時間範囲における画像データ、
(4)コレログラムにより抽出される波形の種類と周期に関する所定時間範囲におけるデータ、
を用い、
健常者の前記心尖拍動により生じる振動の波形を5つに分類し、
前記測定対象者の測定時の前記健康状態を、前記5つに分類された波形分類結果のいずれかに対応するか否かにより、健康状態か否かを推定することが好ましい。
前記生体信号検出センサは、
三次元立体編物と、
前記三次元立体編物の周囲を密閉的に被覆する収容フィルムと、
前記収容フィルムの外側に配設されるマイクロフォンと
前記マイクロフォンをカバーするケースと、前記ケース内で前記マイクロフォンへの外乱の混入抑制機能を果たす外乱混入抑制部材と
を有して構成されることが好ましい。
前記外乱混入抑制部材がゲルであることが好ましい。
また、本発明は、体表面を介して生体信号検出センサにより得られる生体信号データを処理し、コンピュータを、生体信号分析装置として機能させるコンピュータプログラムであって、
前記生体信号データを周波数解析する手順と、
その周波数解析結果から、前記生体信号データ中、心尖拍動により生じる振動と心音により生じる振動との境界周波数を特定する手順と
を前記コンピュータに実行させるコンピュータプログラムを提供する。
前記境界周波数を特定する手順では、
前記周波数解析結果中、調和振動と不規則振動との境界となるパワースペクトルの急変部を求め、この急変部を基準に前記境界周波数を特定することが好ましい。
前記境界周波数を特定する手順では、
前記周波数解析結果に、同時に測定した心音データの周波数解析結果を加味して、前記パワースペクトルの急変部を求めることが好ましい。
前記境界周波数を特定する手順では、
前記生体信号データ及び前記心音データの各周波数解析結果をそれぞれ加算平均処理した波形を対数差分法を用いて両対数軸表示し、両対数軸表示した波形からゆらぎの変化点を求め、このゆらぎの変化点を前記パワースペクトルの急変部として特定することが好ましい。
前記周波数解析する手順では、短時間フーリエ変換が適用され、
前記境界周波数を特定する手順では、前記短時間フーリエ変換の解析結果から前記パワースペクトルの急変部を求めることが好ましい。
前記周波数解析する手順では、前記短時間フーリエ変換の解析結果を、時間、周波数及びパワースペクトルの変動の程度を示す画像データに出力し、
前記境界周波数を特定する手順では、前記画像データから、前記パワースペクトルの急変部を求めることが好ましい。
相関データ記憶部に記憶された前記境界周波数に関連する相関データに照合して、測定対象者の測定時の健康状態を推定する手順を実行することが好ましい。
前記相関データが、前記境界周波数と心拍数との相関データであり、
前記測定時の状態を推定する手順では、測定対象者の測定時の心拍数を、前記境界周波数と心拍数との相関データに照合して、前記測定対象者の測定時の前記境界周波数を推定することが好ましい。
前記相関データが、前記境界周波数と心拍変動のゆらぎ特性との相関データであり、
前記測定時の状態を推定する手順では、前記測定対象者の測定時の前記境界周波数を、前記境界周波数と心拍変動のゆらぎ特性との相関データに照合して、前記測定対象者の測定時の心拍変動のゆらぎ特性を推定することが好ましい。
さらに、前記生体信号データを、前記境界周波数特定手段により特定される前記境界周波数を上限値に設定してフィルタリングし、前記心尖拍動により生じる振動の波形を求める手順と、
前記心尖拍動により生じる振動の波形を分類する手順と
をコンピュータに実行させ、
前記測定時の状態を推定する手順では、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形を、前記波形を分類する手順の実行により分類し、その波形分類結果をもとに、前記測定対象者の測定時の健康状態を推定することが好ましい。
前記相関データ記憶部には、前記波形分類結果と健康状態との相関データが記憶されており、
前記測定時の状態を推定する手順では、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形を、前記心尖拍動により生じる振動の波形を分類する手順の実行により分類し、その測定時の波形分類結果を、前記波形分類結果と健康状態との相関データに照合して、前記測定対象者の測定時の健康状態を推定することが好ましい。
前記心尖拍動波形を抽出する手順の実行により出力され、前記波形を分類する手順における分類対象となる波形情報と前記健康状態とを教師データとして用い、前記波形情報から前記健康状態を推定する推定モデルを機械学習により生成し、
前記測定時の状態を推定する手順が、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形情報を入力として、前記推定モデルを用い、取得した測定時の波形情報から、測定時の前記健康状態の値を出力することが好ましい。
前記波形を分類する手順は、次の(A)及び(B)の手順:
(A)フーリエ級数展開を用いた数学的アプローチにより前記波形を分類する手順、
(B)次の(1)~(4)の工学的アプローチにより得られるデータを1又は2以上組み合わせて前記波形を分類する手順、
(1)前記心尖拍動により生じる振動の波形の所定時間範囲における時間波形の形状を示すデータ、
(2)前記時間波形の前記周波数解析手段により得られる周波数とパワースペクトルを横軸と縦軸にとったグラフのデータ、及び、
(3)前記時間波形の前記短時間フーリエ変換の解析結果である時間、周波数及びパワースペクトルの変動の程度を示す所定時間範囲における画像データ、
(4)コレログラムにより抽出される波形の種類と周期に関する所定時間範囲におけるデータ、
を用い、
健常者の前記心尖拍動により生じる振動の波形を5つに分類し、
前記測定対象者の測定時の前記健康状態を、前記5つに分類された波形分類結果のいずれかに対応するか否かにより、健康状態か否かを推定することが好ましい。
また、本発明は、上記のコンピュータプログラムが記録された記録媒体を提供する。
本発明によれば、心尖拍動と心音との境界周波数を求めることができる。この心尖拍動に関連する新たな指標である境界周波数は、心拍変動のゆらぎ特性と相関があり、境界周波数を用いた健康状態(心疾患等の疾病の有無、疾病の特定、体調の良否等)の推定を行うことができる。また、本発明によれば、この境界周波数が心拍数との関数で表すことができるため、それに基づく相関データを利用して、境界周波数をリアルタイムで知ることができる。
また、本発明では、上記の境界周波数を用いることで求められる心尖拍動により生じる振動の波形を分類し、各波形と人の健康状態(心疾患等の疾病の有無、疾病の特定、体調の良否等)を関連づけることで、健康状態の推定を行うことができる。それにより、医師等の専門家が行っていた心尖拍動図の波形の評価を自動的に行わせることができる。特に、心尖拍動の波形情報と健康状態とを教師データとして用い、波形情報から健康状態を推定する推定モデルを機械学習により生成する構成とすることで、健康状態の推定精度を向上させることができる。
図1(a)は、本発明の第1の実施形態にかかる生体信号検出センサ(4SR)の示した外観斜視図であり、図1(b)は、エアパックとゲルパックとを分離して示した外観斜視図であり、図1(c)は、断面図である。 図2は、本発明の第1の実施形態にかかる生体信号検出センサ(4SR)の構造を、従来の生体信号検出センサ(3SR)との比較で示した図である。 図3(a)は、スピーカーから出力され、実施形態の生体信号検出センサ(4SR)と従来の生体信号検出センサ(3SR)に入力される心音図(Phonocardiogram:PCG)の機械的複製音(PCG複製入力)と、三次元立体編物が伝達する音響振動データを示す図であり、図3(b)は、PCG複製入力の共鳴波形と変調されたPCG複製入力を構成する周波数成分を示す図であり、図3(c)は、実施形態の生体信号検出センサ(4SR)と従来の生体信号検出センサ(3SR)のばね特性及び減衰特性を示したリサージュ図形であり、図3(d)はボード線図を示し、図3(e)は実施形態の生体信号検出センサ(4SR)における三次元立体編物が収容フィルム内に密封収容された構造(sealed air pack)のばね定数を示した図である。 図4は、PCGの波形を入力したときの増幅効果を示した図であり、(a)は、実施形態の生体信号検出センサ(4SR)と従来の生体信号検出センサ(3SR)の出力波形を、(b)は、それぞれの周波数解析結果を、(c)は、4SR/3SRのゲインを示す。 図5は、APW(Acoustic Pulse Wave)の増幅効果を示した図であり、(a)は、実施形態の生体信号検出センサ(4SR)と従来の生体信号検出センサ(3SR)の出力波形を、(b)は、それぞれの周波数解析結果を、(c)は、4SR/3SRのゲインを示す。 図6は、第1の実施形態の生体信号検出センサ(4SR)と、ゲルパックを除いたエアパックのみの構造体との比較実験を説明するための図である。 図7(a)~(c)は、4SRの実験結果を示し、図7(d)~(f)は、エアパック1Aのみの構造体(4SR without gel pack)の実験結果を示す。図7(g)は、(d)の拡大図、図7(h)は、(e)の拡大図である。図7(i)は、ボード線図を示し、図7(j)は、4SRとエアパック1Aのみの構造体の荷重-たわみ特性を示し、図7(k)は、両者の振幅比(4SR/4SR without gel pack)を示した図である。 図8(a)~(d)は、平均心拍数58/minの被験者のF-APW、R-APW、L-APWとPCGの3秒間の時系列波形計測結果を示した図である。 図9(a)~(d)は、平均心拍数71/minの被験者のF-APW、R-APW、L-APWとPCGの3秒間の時系列波形計測結果を示した図である。 図10(a)~(d)は、平均心拍数79/minの被験者のF-APW、R-APW、L-APWとPCGの3秒間の時系列波形計測結果を示した図である。 図11(a)~(d)は、平均心拍数93/minの被験者のF-APW、R-APW、L-APWとPCGの3秒間の時系列波形計測結果を示した図である。 図12は、本発明の第1の実施形態に係る生体信号分析装置の構造を模式的に示した図である。 図13は、境界周波数を抽出する工程を説明するためのフローチャートである。 図14(a)~(d)は、平均心拍数58/min、71/min、79/min、93/minの各被験者のCAB(Cardiac Apex Beat)及びCAS(Cardiac Acoustic Sound)の周波数解析結果を示した図である。 図15は、平均心拍数58/minの被験者から得られたFront CAB(0.5-BF)、Front CAB(5-BF)、Rear CAB(0.5-BF)、Lumbar CAB(0.5-BF)、PCG、Front CAS(BF-50)、Rear CAS(BF-50)、及びLumbar CAS(BF-50)の各処理波形を示した図である。 図16は、平均心拍数71/minの被験者から得られたFront CAB(0.5-BF)、Front CAB(5-BF)、Rear CAB(0.5-BF)、Lumbar CAB(0.5-BF)、PCG、Front CAS(BF-50)、Rear CAS(BF-50)、及びLumbar CAS(BF-50)の各処理波形を示した図である。 図17は、平均心拍数79/minの被験者から得られたFront CAB(0.5-BF)、Front CAB(5-BF)、Rear CAB(0.5-BF)、Lumbar CAB(0.5-BF)、PCG、Front CAS(BF-50)、Rear CAS(BF-50)、及びLumbar CAS(BF-50)の各処理波形を示した図である。 図18は、平均心拍数93/minの被験者から得られたFront CAB(0.5-BF)、Front CAB(5-BF)、Rear CAB(0.5-BF)、Lumbar CAB(0.5-BF)、PCG、Front CAS(BF-50)、Rear CAS(BF-50)、及びLumbar CAS(BF-50)の各処理波形を示した図である。 図19(a)~(d)は、平均心拍数58/min、71/min、79/min、93/minの4名の被験者について、ECG、F-APW、PCGをまとめて示した図である。 図20(a)~(d)は、平均心拍数58/min、71/min、79/min、93/minの4名の被験者について、ECG、Front CAB、Front CASを示した図である。 図21(a)~(d)は、平均心拍数58/min、71/min、79/min、93/minの4名の被験者について、Front CAB及びFront CASの各パワースペクトル密度(PSD)を示した図である。 図22は、計測時間10秒間の平均心拍数93/minの被験者のSTFT解析結果を示した図であり、(a)はF-APWの解析結果を、(b)がPCGの解析結果を示す。 図23(a)は、計測時間10秒間の平均心拍数93/minの被験者のF-APW×PCG-1のSTFT解析結果を示した図であり、図23(b)は、図23(a)の拡大図である。 図24(a)は、計測時間10秒間の平均心拍数93/minの被験者のF-APW×PCGのSTFT解析結果を示した図であり、図24(b)は、F-APW×PCG-1、F-APW、PCG、F-APW×PCGのパワースペクトルを示した図である。 図25(a)は、平均心拍数58/minの被験者のF-APW×PCG-1、F-APW、PCG、F-APW×PCGのパワースペクトルを示した図であり、図25(b)は、F-APW×PCG-1のSTFT解析結果を示した図である。 図26(a)は、平均心拍数71/minの被験者のF-APW×PCG-1、F-APW、PCG、F-APW×PCGのパワースペクトルを示した図であり、図26(b)は、F-APW×PCG-1のSTFT解析結果を示した図である。 図27(a)は、平均心拍数79/minの被験者のF-APW×PCG-1、F-APW、PCG、F-APW×PCGのパワースペクトルを示した図であり、図27(b)は、F-APW×PCG-1のSTFT解析結果を示した図である。 図28(a)は50名の被験者のBF-HR相関図であり、図28(b)は、高調波消失次数とHRから相関図である。 図29(a)は、心拍変動のゆらぎ特性と心拍数の相関を示した図であり、図29(b)は、心拍変動のゆらぎ特性とBFとの相関を示した図であり、図29(c)は、ある1名の被験者の心拍数とBFとの相関を示した図であり、図29(d)は、他の1名の被験者の心拍数とBFとの相関を示した図であり、図29(e)は、さらに他の1名の被験者の心拍数とBFとの相関を示した図である。 図30は、第2の実施形態に係る生体信号分析装置の構成を示した模式図である。 図31は、第2の実施形態に係る生体信号分析装置の他の態様を説明するめの模式図である。 図32(a)~(e)は、Front CAB(0.5-BF)の代表的な波形分析結果であって、疑似正弦波の事例を説明するための図である。 図33(a)~(e)は、Front CAB(0.5-BF)の代表的な波形分析結果であって、三角波の事例を説明するための図である。 図34(a)~(e)は、Front CAB(0.5-BF)の代表的な波形分析結果であって、余弦波の事例を説明するための図である。 図35(a)~(e)は、Front CAB(0.5-BF)の代表的な波形分析結果であって、ガウス波の事例を説明するための図である。 図36(a)~(e)は、Front CAB(0.5-BF)の代表的な波形分析結果であって、疑似全波整流波の事例を説明するための図である。 図37は、2種のSTFT、パワースペクトル、Correlogramを用いて波形検討を行った結果を示した図である。 図38(a)~(e)は、5人の被験者の計測部位別の周波数解析結果を示した図である。 図39は、計測部位別の振幅比を示した図である。 図40は、安静時呼吸抑制時下で心拍変動が生じた際のBF-HR相関図である。 図41(a)は仰臥位におけるBF-HR相関図であり、図41(b)は左側臥位におけるBF-HR相関図であり、図41(c)は右側臥位におけるBF-HR相関図である。
以下、図面に示した本発明の実施形態に基づき、本発明をさらに詳細に説明する。
(第1の実施形態)
・生体信号検出センサ
まず、図1(a)~(c)に基づき、本実施形態で用いた生体信号検出センサ1の構成を説明する。本実施形態の生体信号検出センサ1は、エアパック1Aとゲルパック1Bとの積層構造からなる。エアパック1Aは、三次元立体編物(3Dネット)10及び該三次元立体編物(3Dネット)10を密閉的に収容する収容フィルム20とを有して構成される。ゲルパック1Bは、ケース40内にマイクロフォン30が固定配置され、マイクロフォン30の周囲にゲル50が充填されている。
三次元立体編物10は、互いに離間して配置された一対のグランド編地同士を連結糸で結合することにより形成されている。各グランド編地は、例えば、繊維を撚った糸から、ウェール方向及びコース方向のいずれの方向にも連続したフラットな編地組織(細目)に形成したり、ハニカム状(六角形)のメッシュを有する編地組織に形成したりすることができる。連結糸は、一方のグランド編地と他方のグランド編地とが所定の間隔を保持するように、三次元立体編物に所定の剛性を付与している。従って、面方向に張力が付与されることにより、三次元立体編物を構成する対向するグランド編地の糸、あるいは、対向するグランド編地間を連結する連結糸を弦振動させることが可能となる。それにより、生体信号である心臓・血管系の音・振動によって弦振動が生じ、三次元立体編物の面方向に伝播される。
三次元立体編物のグランド編地を形成する糸又は連結糸の素材としては、種々のものを用いることができるが、例えば、ポリプロピレン、ポリエステル、ポリアミド、ポリアクリロニトリル、レーヨン等の合成繊維や再生繊維、ウール、絹、綿等の天然繊維が挙げられる。上記素材は単独で用いてもよいし、これらを任意に併用することもできる。好ましくは、ポリエチレンテレフタレート(PET)、ポリブチレンテレフタレート(PBT)などに代表されるポリエステル系繊維、ナイロン6、ナイロン66などに代表されるポリアミド系繊維、ポリエチレン、ポリプロピレンなどに代表されるポリオレフィン系繊維、あるいはこれらの繊維を2種類以上組み合わせたものである。また、グランド糸又は連結糸の糸形状も限定されるものではなく、丸断面糸、異形断面糸、中空糸等のいずれでもよい。さらに、カーボン糸、金属糸等を使用することもできる。
使用可能な三次元立体編物としては、例えば、以下のようなものを用いることができる。
(a) 製品番号:49013D(住江織物(株)製)、厚さ10mm
材質:
表側のグランド編地・・・450デシテックス/108fのポリエチレンテレフタレート繊維仮撚加工糸の2本の撚り糸
裏側のグランド編地・・・450デシテックス/108fのポリエチレンテレフタレート繊維仮撚加工糸の2本の撚り糸
連結糸・・・・・・・・・350デシテックス/1fのポリトリメチレンテレフタレートモノフィラメント
(b)製品番号:AKE70042(旭化成(株)製)、厚さ7mm
(c)製品番号:T28019C8G(旭化成(株)製)、厚さ7mm
三次元立体編物10は、収容フィルム20により被覆されている。収容フィルム20は、本実施形態では、合成樹脂製の2枚のフィルム21,22を用いてなり、これを三次元立体編物10の表面及び裏面を被覆するように配置し、両者の周縁部を溶着等により固着している。これにより、三次元立体編物10は、収容フィルム20内に密閉的に収容される。なお、フィルム21,22の周縁部を固着する際には、フィルム21,22によって三次元立体編物10を厚み方向に若干押圧されるように固着することが好ましい。三次元立体編物10の張力が高まり、該三次元立体編物10を構成する糸の弦振動がより生じやすくなる。
収容フィルム20の外側には、ケース40が取り付けられ、そのケース40内にマイクロフォン30が配設されている。ケース40内であって、マイクロフォン30の周囲は、外乱混入抑制部材としてのゲル50が充填されている。ケース40は、合成樹脂製で、マイクロフォン30に伝播される音響振動の外部への拡散を防ぐ機能を有し、ゲル50により、外部振動がマイクロフォン30によって捕捉されることを抑制する。なお、マイクロフォン30には検出した音響振動データを搬送するコード30aが接続されている。
生体信号検出センサ1は、人の各種部位、例えば、背部、胸部、腰部などに当接して使用される。体表面の振動が収容フィルム20及び三次元立体編物10に伝播されてマイクロフォン30によって捕捉されるものであるが、皮膚表面に直接貼着する場合に限らず、衣服の表面に取り付けて用いることができる。
ここで、本実施形態で採用した新規の生体信号検出センサ1の検出性能を、本発明者らが特許文献2等で開示している従来の生体信号検出センサと比較して説明する。なお、以下において、本実施形態の生体信号検出センサ1を4SR(Sound Sensing System using Stochastic Resonance)と称し、従来の生体信号検出センサ1000を3SR(Sound Sensing System using Resonance)と称する。
(心音の周波数帯域での4SRセンシング性能評価)
図2に4SR(生体信号検出センサ1)と3SR(従来の生体信号検出センサ1000)のセンシング性能を比較するための実験装置の概略と部品構成を示す。中央の実験装置の左右の各3枚の写真は、加速度センサ(a-1)とそれぞれ3SRと4SRの構成部品(a-2~a-6)を示し、実験装置下部には3SR(a-7)と4SR(a-8)の概観写真を示す。また、3SR、4SRに用いられている厚さ10mmの3次元立体編物 (以後場合により、3Dネットと呼ぶ)の断面を写真(a-4)に示す。写真(Section:a-a, Section:b-b)は、3SRと4SR用の3Dネットのカット断面である。
まず、4SRの振幅増幅の基本原理を3Dネットの断面の写真(a-4)を用いて説明する。本実験で用いた3Dネットの上下面はマルチフィラメント糸で作られた基布層で、真ん中のモノフィラメント糸で作られたパイルは上下の基布層に編み込んであり、マルチフィラメント糸と摩擦結合している。被験者の体重が3Dネット上面に付加されると、エアパックに内蔵された3DネットのX型に編み込まれたパイル(以後、X型パイルと呼ぶ)がたわみ、体動によりパイル同士が擦れ、音と振動が生じる。この音と振動が確率共鳴現象を生じさせる。そして、X型パイルはたわむことで張力が発生し、固有振動子が構成される。張力作用下のX型パイルの固有振動数は20Hz近傍にある。そこで、密閉エアパックである4SRは確率共鳴と弦の振動という二つの共鳴機構により、0.5~80Hzの振動を増大させるものと考えられる。
次に、心音に対する増幅能について検討する。図3(a)に示したように、検討に用いた心音は、PCG上I音の分裂間隔を0.03秒にしたI音を用いた。この心音は心拍数56/minの被験者から計測されたもので、フィルター処理を行った後、スピーカーによって機械的に複製した。スピーカーはFOSTEX社製P1000Kをエンクロージャー無しで用いた。P1000Kは径が10cmのコーン型フルレンジユニットである。P1000Kの振動板に600gのウェイトを載せ振動板の重さを1005gにし、最低共振周波数を82Hzから52Hzに低く設定した。これにより52Hz~16kHzが再生周波数帯域となった。最低共振周波数f0は、
Figure 2022071812000002
で求められる。s0はばねの強さを示し、m0は振動板の重さを示す。
スピーカーから出力された機械的複製音(以後、Reproduction PCGと呼ぶ)は、加速度センサの上面に置かれた3Dネットを経由した後、3SRと4SRの各マイクロフォンで計測される。3SRのマイクロフォン1010は、3Dネット1020とともにエラストマフィルム1030とビーズ1040で包まれる。減衰は、スリット内を流れる空気の流動によって、内蔵された3Dネットがもつばね力によりマイクロフォン周りの空気が流動することによって発生する。その大きさは写真(a-7)に示すようにマイクロフォン付属のコード1010a周辺部のスリット幅で決められる。3SRは、エラストマフィルム1030内にマイクロフォン1010を配置しているため、コード1010aをエラストマフィルム1030外へ引き出さざるを得ず、その部位において必ず隙間が生じる。一方、4SRは、三次元立体編物(3Dネット)10が収容フィルム20内に収容されたエアパック1Aを、マイクロフォン30が配置されているゲルパック1Bとは分離した構成であるため、マイクロフォン30のコード30aはゲルパック1Bから外部に引き出されるため、エアパック1Aの密閉性がそれにより阻害されることがない。すなわち、4SRは、三次元立体編物(3Dネット)10の全周がエラストマフィルム(収容フィルム20)で密閉されており、4SRに負荷がかかると均一の空気圧がエラストマフィルムにかかることになるため、エラストマフィルムに張力が生じ、音響振動が伝播しやすい環境が作られる(図3(c)参照)。
なお、図3(a)は、上記のようにPCG複製入力と3Dネットが伝える音響振動データを示すが、3Dネットが伝える音響振動データは、3Dネットの固有振動子の共鳴効果による波形に続いて、変調されたPCG複製入力からなる。図3(c)のリサージュ図形の面積は、振動の1サイクル中に消費するエネルギーを示し、減衰容量から求めた3SRの減衰比は0.53で、4SRは0.20であった。
(心尖拍動の周波数帯域での4SRセンシング性能評価)
健常者のACGでは、収縮早期のみに認められる持続時間の短い拍動(収縮期波は持続の短い尖鋭な陽性波とそれに続く下行脚(dicrotic limb)からなり、左心房の収縮に由来するA波は痕跡のように認められる)でtappingと形容され、ACGで良好な記録をとることが難しい。そこで触診でtappingと形容される心尖拍動を触知できる座位姿勢の被験者(心拍数:62/min)を選定し、3SRと4SRを用いてF-APW(胸部前面から計測した生体情報:Front Acoustic Pulse Wave)を計測し、データロガーに記録し、フィルター処理を行って以下の解析を行った。
データロガーに記録した入力波形を周波数解析し、周波数に対する4SRと3SRのパワースペクトルの比(PSD4SR/PSD3SR)をゲインとし、センシング性能をゲインで評価する。4SRは左第5肋間、胸部正中線から左6cmの胸部前面に配置し、3SRは左10cmの胸部前面(心電図のV4誘導付近に相当)に配置する。PCGならびに4SRセンシングシステム用の各マイクロフォンの設置位置は、肺胞呼吸音が混入し易い部位であるため、呼吸を止めた状態で、F-APWを10秒間計測する。なお、計測実験の手順は安静時の呼吸を5分間継続して2秒間の吸気の後、呼吸を止めた。被験者には40歳代の心循環系に基礎疾患のない健康な男性を選定した。
(APW計測)
次に、APWの計測を行う。計測は、全て座位姿勢で行う。心尖拍動を含むと考えられるF-APWを捉えるための4SRは、左第5肋間、左鎖骨中線上(V4付近)で、胸部正中線から左10cmの胸部前面に配置し、R-APW(胸部後部から計測した生体情報:Rear Acoustic Pulse Wave)のための4SRは、胸部前面のF-APWのための4SRと同じ高さで背部正中線から左6cmの胸部後面に配置する。腰部のL-APW(腰部から計測した生体情報:Lumber Acoustic Pulse Wave)を捉える4SRは臍部の真後ろである第3~4腰椎部正中に配置し、PCG用マイクロフォンは心尖部に置く。また、ECG(Electrocardiogram:心電図)は、II誘導を取得する。
なお、胸部前部の心尖部マイクロフォンセンサ配置は、I音、III音、IV音の至適記録部位であるが、II音については良好な記録をとることが難しい。II音は、収縮終期に心室へ逆流しようとする動脈内血液が半月弁の過伸展、ついで動脈壁に反動を生じさせる(cardiohemic system)ことで発生すると言われており、心尖部では減弱する。
PCGならびに4SRセンシングシステム用の各マイクロフォンの設置位置は、肺胞呼吸音が混入し易い部位であるため、呼吸を止めた状態での10秒間の計測結果を用いて計算する。なお、安静時のAPW計測実験は、a)自然呼吸を5分間継続して行った状態、b)2秒間の吸気の後30秒間呼吸を止めた状態、c)30秒間の自然呼吸の状態、d)2秒間の吸気の後30秒間呼吸を止めた状態の4つの条件で行う。
被験者は、比較的胸壁の薄い集団と比較的胸壁の厚い集団の異なる特徴を持つ健常な男女70人とした。
比較的胸壁が薄い健常人で構成された50人は、年齢は20~60歳代に分布し、男性が78%だった。胸壁が薄いと心尖拍動の特徴である正弦波ないし三角波を観測するのが容易となる。この対象者を用いて、後述する5つの波形タイプの定義と、心尖拍動に起因する振動と心音に起因する振動との境界周波数(Boundary Frequency:BF)を決定した。
比較的胸壁の厚い20歳代で筋肉質な被験者20人を用いて、先行して決定した5つの波形タイプをもとに、解析方法の適合性評価と拡張性の検討を行った。
比較的胸壁が薄い健常人で構成された50人の臨床プロフィールは、次のとおりである。平均年齢 38.2(±10.3)歳、男性39名、女性11名、BMI値平均 23(±3.4)、そのうち高血圧症4名、喫煙者10名である。
一方、比較的胸壁の厚い20人の臨床プロフィールは、平均年齢 23.4(±3.6)歳、男性20名、BMI値平均 20.8(±2.4)、そのうち高血圧症1名、喫煙者5名である。
解析対象データは、呼吸を止めた最初の30秒間のうちRR変動が15%以内の部分から10秒間のデータを用い、心拍数はその10秒間のデータでの平均値とした。
(4SRセンシングシステムの性能の評価結果とAPW計測結果)
図3(a)中で黒の細線で示した波形がReproduction PCGで、グレーの太線で示した波形が加速度センサで計測したAcceleration PCGである。Reproduction PCGは、I音(A波とB波)と30Hz近傍の波形(C波)、およびII音(D波、E波、F波)で構成され、Acceleration PCGは、Reproduction PCGの共鳴波形(a波とb波、およびd波とe波) とそれに続く振幅と周波数が変化した波形(c波とf波)で構成される。図3(b)の周波数解析結果で示す様に60Hz以上の周波数帯域にあるA波とB波およびD波とE波は、振幅が大きくなって、a波とb波、およびd波とe波になっており、30Hz近傍のC波とF波は、ヘテロダインにより振幅と周波数が変調されたc波とf波になった。その結果、30Hz近傍ではAcceleration PCG のパワースペクトルは、Reproduction PCGのものよりも小さくなった。
次に、センシングシステムを構成する機械振動系に含まれる減衰性能について検討する。図3(c)は、3SRと4SRの各センサがもつ減衰特性を示すリサージュ図形である。リサージュ図形は加圧板を110×110mmにして振幅±1.0mm、加振周波数1.34Hzの入力で、サーボパルサーを用いて描いた。予備圧縮力は652Nである。リサージュ図形で囲まれる面積は、振動の1サイクル中に消費するエネルギーで、減衰容量である。これをWと表すと、Wは、
Figure 2022071812000003
となる。4SR と3SRの動的特性を比較すると、3SRは4SRとほぼ同等のばね特性に加えて、エアの流入流出による減衰の影響が顕著であることが分かる。
また、振動変位を評価する場合に用いられる振幅倍率Z/Yは、相対変位をZとし、心尖拍動の変位をYとすれば、
Figure 2022071812000004
で与えられる。ここに4SRと3SRの各センサの機械振動系の固有角振動数は、
Figure 2022071812000005
であり、減衰比は、
Figure 2022071812000006
である。
図3(d)は、式(3)を用いて計算した4SRと3SRの振幅倍率曲線を示す。4SRの各センサのばね定数を比較したものである。
Figure 2022071812000007
の値を代入した。
なお、図3(e)は、φ15mmの加圧板で図2(a-8)で示したセンサ(4SR)をデジタルフォースゲイジ(RZ-20)で計測して得た荷重-たわみ特性である。心尖拍動の計測を想定してばね定数を計算した。胸部前部から計測する場合は、センサの押し付け力を小さくしているために、たわみ量が1mmの場合のばね定数を用いた。胸部後部からの計測は体重がセンサにかかるため、たわみ量が3SRの設計値である4~5mmの場合のばね定数を用いた。
4SRはセンシングシステムの固有振動数を低く設定してあるため、Z/Y≒1の範囲で測定すると4SRは0.45Hz以上の周波数帯域が測定でき、3SRは2.10Hz以上が測定できることが分かる。なお、減衰比をζ=0.7とした場合、4SRは0.99Hz以上が測定でき、3SRは3.15Hz以上が測定できる。このように、減衰比を変化させるとメカニカルフィルターの特性を調整することが可能である。
心臓聴診は、50~100Hz近傍の高い周波数範囲にある心音の特徴をつかみ、解釈し、それを確認することにあり、高い周波数の音響振動を記録することが重要である。一方、心尖拍動は心音よりも低い周波数波範囲内にある収縮早期のみに認められる持続の短い拍動であり、高い周波数の計測に力点を置くと良好な心尖拍動図を記録するのが難しくなる。このように、心音と心尖拍動ではセンシングシステムの要求仕様が異なるため、センシングシステムの定量評価を個別に評価した。図4(a)は、Acceleration PCGを入力にして、3SRと4SRの各センシングシステムが捉えた音響振動情報(以後、3SR data、4SR dataと呼ぶ)を比較したものである。また、図4(b)は、周波数解析結果を示す。3SR data及び4SR dataは、確率共鳴と固有振動子による共鳴の影響を受けた音響振動情報で、両者を比較すると同じ周波数帯域でパワースペクトルに増減が見られる。図4(c)は、4SR/3SRのゲイン(PSD4SR/PSD3SR)を示す。10Hz近傍の確率共鳴の効果と20Hz近傍の固有振動子の共鳴効果及び20~80Hzのセンサのもつ機械的振動特性の差がゲインの変化(最大11.6dBと最小8.0dB)となって表れ、図4(a)で示す波形の最大振幅差は4.3倍になった。なお、80Hz以上のゲインの急低下は4SRセンシングシステムに取り付けられたゲルの減衰効果によるものと考えられる。
図5(a)に、スピーカに変えて人を音源としたときの3SR dataと4SR dataの時系列波形を示し、図5(b)に周波数解析結果を示す。図5(c)に、4SR/3SRのゲイン(PSD4SR/PSD3SR)を示す。3SRセンシングシステムに内在する減衰は1.34Hz以上で機能し、1Hz近傍のゲインが13.4dBであるのに、1.34~10Hz近傍のゲインは4.9dBまで低下している。波形で見ると減衰の影響は顕著で、図5(a)に示されるように収縮期波の最大振幅は6.8倍となって表れた。
なお、3SRのエアダンピングが、心音に関係する不規則振動に与える影響は小さいものの、図5(b)に示すように13Hz以下の心尖拍動波のパワースペクトルには影響を与え、その値は小さくなった。また、図3(c)に示すダンピング効果は、主に基本調波と2次高調波で顕著に認められた。これは、図3(d)のシミュレーション結果とも一致し、3SRでは心尖拍動波は描記できないことが分かった。
振動伝達機構について考察すると、エラストマフィルム(4SRでは符号20の収容フィルム、3SRでは符号1030のエラストマフィルム)自体にも張力が生じるため、張力変動を介して、4SRでは1Hz近傍の低周波振動は内部の3Dネットに伝達される。その結果、13Hz以上の振動は、20Hz近傍で生じる固有振動子の効果で、図5(b)に示される13~30Hz間のほぼ同等のパワースペクトルとして3SRと4SRで計測された。その振動成分が3Dネットのパイルを伝わるものと考えられる。また、図5(b)で心音の最小周波数と心尖拍動波の高調波成分が心音に混じる周波数帯が10Hz近傍にあることが分かった。そこで、4SRの方が、考察対象周波数の全帯域にわたり、3SRよりも感度が高いことが分かった。
従来、3SRの解析対象周波数は10~30Hzに設定していたが、以上の実験結果より4SRが計測対象とするF-APWの振動数は、マイクロフォンの性能保証範囲が0.1Hz以上であることを考慮し、0.5~80Hzまでとし、心音については、大部分の音が低い可聴域に集まっていることで、100Hz以上ではなく、25~45Hzないし40~80Hzを解析対象周波数とすることとした。また、サンプリング周波数は、3SRでは200Hzとしていたが、4SRでは1000Hzに設定した。
ここで、本実施形態の生体信号検出センサ1(4SR)と、4SRからゲルパック1Bを取り除き、エアパック1Aのみの構造体(4SR without gel pack)で音響振動を測定した場合との比較実験を行った。図6は、その実験の様子を示した図であり、ウレタンベース上、図の右側に生体信号検出センサ1(4SR)を、左側にエアパック1Aのみの構造体(4SR without gel pack)を配置して実験を行った。その他の実験条件は、上記の3SRと4SRの比較実験と同様である。図7にその結果を示す。
図7(a)~(c)は、4SRの実験結果を示す。図7(d)~(f)は、エアパック1Aのみの構造体(4SR without gel pack)の実験結果を示す。(a),(d)は、人を音源としたときに捉えられた音響振動情報の時系列波形を示し、(b),(e)は、(a),(d)の波形の振幅の中心をゼロ点として示した時系列波形であり、(c),(f)は、Acceleration PCGを入力して捉えられた音響振動情報の時系列波形を示す。図7(g)は、(d)の拡大図、図7(h)は、(e)の拡大図である。
これらの図から明らかなように、4SRは、エアパック1Aのみの構造体(4SR without gel pack)よりも、音響振動情報の検出感度が高い。特に、人を音源とした場合、4SRでは、±0.2[V]程度の振幅変化が生じ、音響振動情報を明確に捉えることができていたのに対し、エアパック1Aのみの構造体では振幅変化を捉えることは非常に困難であった。また、図7(j)に示したように、4SRは、エアパック1Aのみの構造体と比較して、荷重-たわみ特性における比較でヒステリシスロスが小さかった。
また、マイクロフォン30は、合成樹脂製のケース40の底面に、薄板の円板状の支持部材(図6の符号30bの部材)上に配置されているが、ケース40は、厚さ0.1~0.2mm程度で、円板状の支持部材30bよりも剛性が低い。そのため、ケース40は、音響振動に反応して振動しやすく、音響振動を共鳴させ、いわばスピーカーにおける振動板とエンクロージャーのような役割を果たし、マイクロフォン30による検出感度を高め、振幅を大きくする機能を有している。図7(k)は、両者の振幅比(4SR/4SR without gel pack)を示したものである。この図から、4SRは、心尖拍動のような周波数の極めて低い音響振動に対して検出感度が高いだけでなく、心音の領域(25~80Hz)においても検出感度が高いことがわかる。
また、胸部に直接取り付けたマイクロフォンにより計測される、体幹から出力された微小音響・振動(micro acoustic vibration(MAV))と、4SRのから計測されるF-APW(4SRを構成する上記部材の機械振動系の応答、並びに、固有振動子と確率共鳴による共鳴の影響を受けた音響振動情報)の各時間波形を周波数解析し、振幅倍率を計算したものがF-APW/MAVのゲインとなるが、F-APW/MAVのゲインは、心尖拍動の基本調波の周波数帯域で42dBとなり、高調波の周波数帯域で30dB強となり、固有振動子の周波数帯域である20~30Hz近傍の心音の最小周波数帯域で20dBとなった。機械振動系の応答と確率共鳴の効果は、最大で40dB超で最小10dBとなって表れ、振幅差は200倍弱となった。F-APWは1~1.5Hzの帯域の波形が観測され、心音と心尖拍動の情報が含まれた。MAVとF-APWのTime Lagは0.12秒であった。
図8~図11に、平均心拍数58、71、79、93/minの被験者4名についてのF-APW、R-APW、L-APWとPCGの3秒間の時系列波形計測結果を示す。これら3種類のAPWには心音と心尖拍動の情報が含まれる。F-APWは1~1.5Hzの帯域の波形、R-APWとL-APWでは5~7Hzの帯域の波形が観測された。
・生体信号分析装置
次に、本実施形態の生体信号検出センサ1から得られるデータを処理するコンピュータプログラムが設定されたコンピュータ機能を有する生体信号分析装置100について図12に基づき説明する。
生体信号分析装置100は、生体信号検出センサ1によって取得される生体信号の時系列データを処理して心尖拍動の波形を得、さらに生体状態、すなわち、人の様々な健康状態(心疾患等の疾病の有無、疾病の特定、体調の良否等)を推定する。生体信号分析装置100は、コンピュータ(パーソナルコンピュータ、機器に組み込まれるマイクロコンピュータ等も含む)からなり、生体信号検出センサ1のマイクロフォン30から送信される生体信号の時系列データを受信する。そして、受信した時系列データを用いて所定の処理を行う周波数解析手段110、境界周波数特定手段120を有している。
また、記憶部には、予め測定して構築された、心尖拍動に起因する振動と心音に起因する振動との境界周波数と、心拍数との相関データが記憶された相関データ記憶部150を有している。
より詳細には、生体信号分析装置100は、周波数解析手段110、境界周波数特定手段120として機能する手順を実行させるコンピュータプログラムが記憶部(当該コンピュータ(生体信号分析装置100)としての内蔵のハードディスク等の記録媒体のほか、リムーバブルの各種記録媒体、通信手段で接続された他のコンピュータの記録媒体等も含む)に記憶されている。なお、生体信号分析装置100は、周波数解析手段110、境界周波数特定手段120を実現するコンピュータプログラムが組み込まれた1以上の記憶回路を有する電子回路を用いて実現することもできる。
また、コンピュータプログラムは、記録媒体に記憶させて提供することができる。コンピュータプログラムを記憶した記録媒体は、非一過性の記録媒体であっても良い。非一過性の記録媒体は特に限定されないが、例えば フレキシブルディスク、ハードディスク、CD-ROM、MO(光磁気ディスク)、DVD-ROM、メモリカードなどの記録媒体が挙げられる。また、通信回線を通じてコンピュータプログラムをコンピュータに伝送してインストールすることも可能である。
ここで、生体信号検出センサ1は、上記のように心音を捉えることができることから、得られた生体信号データを受信したならば、周波数解析手段110は、体表面を介して生体信号検出センサ1により得られる生体信号データを周波数解析する。
境界周波数特定手段120は、周波数解析手段110から得られる生体信号データの周波数解析結果から、生体信号データ中、心尖拍動に起因する振動と心音に起因する振動との境界周波数(Boundary Frequency:BF)を求める。境界周波数特定手段120は、周波数解析手段110により得られる周波数解析結果中、調和振動と不規則振動との境界となるパワースペクトルの急変部を求め、この急変部を基準に境界周波数を特定する手段を含む。
本実施形態によれば、境界周波数が求められることで、生体信号から心尖拍動の波形を心音の波形とは区別して得ることが可能となる。
また、本実施形態の生体信号分析装置100は、境界周波数そのものを判定指標として用いている。すなわち、後述の実験から明らかなように、境界周波数は心拍数との2次関数で表され、心拍数を基に両者の相関データから境界周波数を知ることができる。境界周波数と心拍数との相関を見ると、効率のよい状態の心臓の心拍数は約70~80/minである。これは、心臓の運動負荷が少ない状態であり、その場合の境界周波数は最小となっている。よって、心拍数を測定して相関データを用いて得られる測定時の境界周波数により、その時点の心臓がどのような状態であるかが推定できる。
本実施形態の生体信号分析装置100は、図12において想像線で示したように、境界周波数と心拍数との相関データが記憶された相関データ記憶部150を有し、相関データ記憶部150に記憶された相関データに照らし、測定対象者の心拍数をもとに境界周波数を求め、測定対象者の測定時の境界周波数を推定する測定時状態推定手段130を有する構成とすることが好ましい。
測定時状態推定手段130によって状態を推定する際には,心拍数の情報が必要となるが、心拍数の情報は、図12に想像線で示した心拍数演算手段140により、心拍数を求めることができる。この場合には、生体信号検出センサ1と、その生体信号データを処理して心拍数を求める心拍数演算手段140をあわせて心拍数測定部として機能する。なお、心拍数の測定は、心電図法、心音図法、光電脈波法、血圧計法等があり、それらの測定装置を心拍数測定部として用いることも可能である。但し、本実施形態の生体信号検出センサ1は、人の体表面に当接するだけで、非拘束で測定でき、測定のしやすさの点で優れている。
相関データ記憶部150に記憶される相関データは、予め測定された心尖拍動に起因する振動と心音に起因する振動との境界周波数と、心拍数との相関データからなる。この相関データは、境界周波数と、境界周波数を求めた際の心拍数との相関を示した2次関数で構築されている。よって、測定した心拍数をこの相関データ記憶部150の記憶部に記憶された相関データに照合することで境界周波数を知ることができる。
なお、相関データ記憶部150に記憶される境界周波数は、予め、多数の生体信号データを解析して求められる。境界周波数は、境界周波数特定手段120により特定される上記の急変部として求められる。
急変部の求め方としては、まず、体表面を介して生体信号検出センサ1により得られる生体信号データを周波数解析し、その周波数解析結果に、同時に測定した心音データの周波数解析結果を加味し、生体信号データ中、心尖拍動に起因する振動と前記心音に起因する振動との境界周波数を抽出して求めることができる。
具体的には、周波数解析手段110により求められる生体信号データの周波数解析結果を加算平均処理し、その波形を周波数とパワースペクトルを用いて両対数軸表示すると共に、心音データの周波数解析結果を加算平均処理してその波形を周波数とパワースペクトルを用いて両対数軸表示し、両対数軸表示した2つの波形の対数差分の波形からゆらぎの変化点を求め、このゆらぎの変化点を心尖拍動の高調波成分が極めて小さくなる実質的な消失点とし、消失点の周波数を上記の急変部に相当するものとして境界周波数を求めることができる。
また、周波数解析手段110として、短時間フーリエ変換(STFT : Short Time Fourier Transform)の手法を適用することができる。この場合、境界周波数特定手段120は、短時間フーリエ変換の解析結果から上記のパワースペクトルの急変部を求める。短時間フーリエ変換は、好ましくは、その解析結果を、時間、周波数及びパワースペクトルの変動の程度を示す画像データに出力する。これにより、境界周波数特定手段120は、短時間フーリエ変換の画像データから、パワースペクトルの急変部を求めることができる。具体的な急変部の特定、すなわち、境界周波数を求める手段については後述する。
ここで、生体信号検出センサ1は、人の体の背部、胸部、腰部等に取り付けられ、体内の音、振動を捉えるものであり、得られる生体信号は、種々の生体音、体内振動の集合体である。一方、心尖拍動は、心音の波形に隠れた振動であり、それらを分離することは難しい。しかしながら、本実施形態では、心尖拍動と波形と心音の波形とを区別する上記の境界周波数(Boundary Frequency:BF)を見出し、心拍数に基づいた心尖拍動に関連する生体状態の推定を可能としている。
以下、周波数解析手段110により周波数解析し、境界周波数特定手段120により境界周波数を特定する方法について詳述する。本実施形態においては、上記のように、両対数軸表示した波形を用い、その波形において高調波調和振動が実質的に消失する周波数に基づきゆらぎの変化点を求め、このゆらぎの変化点をパワースペクトルの急変部とする手法と、短時間フーリエ変換の解析結果からパワースペクトルの急変部を求める手法の2つがある。
[境界周波数(BF)の抽出法]
(APWの高調波調和振動消失周波数に基づくBF抽出法)
心音は、非特許文献1に示されているように不規則振動であるが、心尖拍動波は、後述の実験結果より、主に基本調波と高調波からなる調和振動で構成されることが明らかになった。そこで、胸部前部に取り付けた4SRセンシングシステムから得られるAPWの波形を周波数解析し、加算平均処理を行い、高調波調和振動と不規則振動のパワースペクトルの変化点、すなわちパワースペクトルの急変部を見つけることに着目した。心尖拍動波のCAB(Cardiac Apex Beat)の高調波成分は周波数が高くなるにつれて、そのパワースペクトルが小さくなる。一方、心音のCAS(Cardiac Acoustic Sound)の不規則振動系のパワースペクトルは、周波数に依存せず、パワースペクトルの大きさが変化する。ここでは、このCABの高調波成分のパワースペクトルが小さくなり、そしてCASの変動挙動が変化する周波数を境界周波数:Boundary Frequencyと呼び、以下では、これをBFと呼ぶ。なお、心拍変動や血圧変動が生じるとAPWの波形が瞬時に変動するが、APW波形の瞬時変動の要因の一つに、PCGの卓越周波数の変動があると考えられる。
被験者の個人差により、体幹の音響振動の伝達特性に関してもインピーダンスが異なり、BFを見つけにくい場合もある。そこで、対数差分法を適用して心尖拍動波の高調波成分が消失する周波数(高調波成分のパワースペクトルが非常に小さくなり、実質的に無視できるレベルの周波数)を見つける。APWの加算平均処理は時間窓8.2秒で、90%オーバーラップでフーリエ変換を行う。加算平均処理された波形を両対数軸表示で表す。両対数軸表示のパワースペクトル波形は、加法定理が使えるためである。なお、心音にはPCGの波形を用いる。対数差分法を適用することで、PCGの0.5~20Hzの周波数帯域のパワースペクトルの絶対値が大きくなり、かつパワースペクトルの波形を逆位相にできる。これにより、心尖拍動の高調波成分と心音の不規則振動のパワースペクトルを小さくし、BFを見つけ易くする。
(短時間フーリエ変換を用いたBF抽出法)
F-APW、PCG、APW×PCG、APW×PCG-1にパワースペクトルの短時間フーリエ変換(STFT : Short Time Fourier Transform)を適用し、ゆらぎの変化点を可視化する。本解析で用いる短時間フーリエ変換の式は、
Figure 2022071812000008
となる。
x(t)は計測信号データであり、W(t-τ)は窓関数である。周波数分解能を0.5Hz以下にするために、窓関数の点数は212=4096とする。4096点はサンプリング周波数を1000Hzにとると4.096秒となり、4.096秒は0.244Hzとなる。これにより、周波数分解能を0.5Hz以下の0.244Hzに設定できる。そして、窓関数の移動幅は50点に設定した。tは計測時間を示し、τは窓関数の移動時間による時系列データを示す。計測信号データに窓関数を掛け合わせることで、時間と共に変化する周波数情報を作ることができる。なお、周波数情報は50/1000=0.05秒毎に生じる。
次に、上記二つのBF抽出法に基づいて決定したBFを用いて、F-APWから心尖拍動波形(Front CAB)と心音(Front CAS)を描記し、同様にR-APW、L-APWからRear CAB、Rear CASおよびLumbar CAB、Lumbar CASを描記する。また予備実験と予備的検討で、BFに関与しているのがHRであることが判明したので、両者の相関関係を調べるためにBFとHRの相関図を50名分のデータを用いて作成する。
なお、高調波成分が消失する周波数に基づくBF抽出法は周波数からBFを決定する手法で、一方、短時間フーリエ変換を用いたBF抽出法は、パワースペクトルを閾値としてBFを決定する手法である。具体的には、前者は周波数を変化させてパワースペクトル高さが下降から上昇に転ずる周波数をBFと定める方法で、後者はBF付近にあるパワースペクトル群の平均値の70%を閾値に設定し、BFを定める方法である。
ここで図19~図21に、HR(心拍数):58/min、71/min、79/min、93/minの4名の被験者について、ECG、F-APW、PCG、Front CAB、Front CASの比較をまとめて示す。
HR:58/minの被験者事例では、Front CABによって計測された収縮期波がACGのLate systolic bulgeを呈しているが、ECGとPCG(Front CAS)からは拡張型心筋症は考えにくく、重力の影響で心尖拍動がSustainedなパターンに近づいたと考えられる。
HR:71/minの被験者事例では、ECGではP波がはっきりとせず、小さくて見えない。Front CABにはACGで確認される左房収縮に由来するA波(心房性隆起波)が確認でき、Front CASにもPCGのI音より早く低周波の波があり、Ectopic rhythm(異所性調律)ではない。
HR:79/minの被験者事例では、Front CAB によって計測された波形はSustainedなパターンを呈すが、A波とRapid Filling Wave (RFW)の波形とFront CASから、左室肥大と収縮期に左室からのリークがある病例(僧房弁逆流等)ではないことが示唆される。
HR:93/minの被験者事例では、ECGで低電位を示し安静状態で高い心拍数を示すが、低周波成分・高周波成分をもつFront CABとI音が亢進も減弱していないFront CASから、基礎心疾患の疑いは消える。
(心拍変動のゆらぎ特性と心拍数およびBFとの相関)
両対数軸表示のパワースペクトルの傾きをここではゆらぎ係数と呼び、急峻なゆらぎ係数の変化点を、以後、ゆらぎの変化点と呼ぶ。ゆらぎの変化点では、複雑性の喪失やBreakpointが生じる。ゆらぎ係数はパワースペクトルP(f)と周波数(f)で表され、両対数軸表示でlog P(f)とlog(f)となり、log P(f) とlog(f)の間に直線関係がある場合は、
Figure 2022071812000009
で表される。
ECGから求めた心拍変動のゆらぎ特性を評価するためAPW計測の実験プロトコルにより計測した実験データを用いる。なお、ゆらぎ特性は、式(5)を用いて計算する。得られるnの値はマイナスの数値(およそ-0.5~-1.5)となるが、非特許文献3によれば、心周期の1/fゆらぎと定義されているため、絶対値で表記する。
心拍変動のゆらぎを評価するため、無呼吸、自然呼吸、無呼吸の各30秒間、計90秒間のデータを用いる。ここでの心拍数は、90秒間の平均値を用いる。
[各抽出法によるBFの抽出結果]
(F-APWの高調波調和振動消失周波数からのBF抽出の結果)
図13は、F-APW×PCG-1を用いてのBFを抽出する手順を示したフローチャートであり、次の(1)から(10)は図13中の丸囲み数字に該当する。以下、計算の手順と着眼点を説明する。
(1) 4SR dataからフィルタ処理を行わず、F-APWを求め、周波数解析する。以下ではこれをF-APWと表す。
(2) F-APWに対して加算平均処理を行う(時間窓8.2秒で、90%オーバーラップ)。
(3) PCGから波形を抽出し、周波数解析する。以下ではこれをPCGと表す。
(4) (2)と同様の加算平均処理(時間窓8.2秒, 90%オーバーラップ)を行う。
(5) (2)と(4)で求められたF-APW とPCGのパワースペクトルF-APW、PCGを用いてF-APW×PCG -1を生成し、Breakpointを見つける。
(6) F-APW×PCG-1 波形上にあるBreakpointから、ゆらぎをもつ心尖拍動波の高調波消失周波数とホワイトノイズとなる心音の不規則振動成分出現周波数を同定する。
(7) F-APW×PCG -1上の高調波消失点からBFラインを引き、BFラインとF-APWとの交点をBFとする。
(8) 0.5Hz~BFがCAB帯域(以後、CAB (0.5-BF))となり、BF~50HzがCAS帯域(以後、CAS(BF-50))となる。
(9) F-APWに0.5Hz~BFのバンドパスフィルタを適用し、CABの高調波成分がBF近傍で減少していることを確認する。
(10) F-APWにBF~50Hzのバンドパスフィルタを適用し、BF近傍では周波数に依らずパワースペクトルが変化していることを確認する。
図14は、両軸線形表示で心拍数(HR)別(58~93/min)にまとめた CABとCASの周波数解析結果である。この図から、BFでのCABの高調波成分の消失、BF近傍からCASの不規則振動の出現、およびBFの HR依存性が理解でき、本発明の手法が、CABとCASを分離してその特徴を可視化できることがわかる。
次に、図15~図18に、APWの高調波調和振動消失周波数に基づくBF抽出法によって描記されたCAB、CAS波形を示す。平均心拍数58、71、79、93/minの被験者から得られたF-APW、R-APW、L-APWから抽出されたFront CAB(0.5-BF)、Front CAS(5-BF)、Rear CAB(0.5-BF)、Lumbar CAB(0.5-BF)および PCGならびにFront CAS(BF-50)、Rear CAS(BF-50)、Lumbar CAS(BF-50)の各処理波形を示す。描記時間は1.5秒間である。括弧内は適応したフィルターの周波数帯域を示す。各CAB、各CAS波形の縦軸、横軸のレンジは計測部位毎で統一されている。4人の被験者のFront CAB(0.5-BF)は、基本調波が主となる心尖拍動波が描記されていることがわかる。
得られた結果の考察から、Front CAB(5-BF)は、複数の5~10Hzの波形の組み合わせで形成されていることが分かった。一方、Rear CAB(0.5-BF) 、Lumbar CAB(0.5-BF)については、Front CAB(0.5-BF)に示す波形より周波数の高い波形が主体の調和振動波形となり、Front CAB(0.5-BF)とは異なる波形になった。また、心音は、複数の振動成分と複数の振動波形が合成された高周波波形からなる不規則振動であった。
(短時間フーリエ変換を用いたBF抽出の結果)
図22~図27は、計測時間10秒間の平均心拍数が93/min、58/min、71/min、79/minの被験者のF-APW、PCGの周波数解析によるパワースペクトルのSTFT解析結果を示す。
STFTの図は、縦軸は周波数を示し、横軸はSTFTの出力時刻を表す。色の濃淡はパワースペクトルの変動の程度を示し、色の濃淡のレンジを変えることで式(5)で表されるゆらぎ係数の大小も表現できる。パワースペクトルの変動の大小はRed, Green, Blue(RGB)の濃淡で示される。パワースペクトルの変動が少ない箇所やゆらぎ係数が急変しない部位や小振幅の不規則振動では、面状のRed、Green、Blueが出現し、逆にパワースペクトルの変動が大きい所や、ゆらぎ係数が急変する部位はRed、Green、Blueが線状に出現する。
上記の対数差分(F-APW×PCG-1)は、図23(a),(b)に示される。図24(b)は周波数解析結果を示し、Breakpointに対応するBFは15Hzとなる。
図24(a)では、0.5~30Hzの周波数帯域にある高調波成分と不規則振動を強調するために、APWとPCGの対数和分(APW×PCG)を計算した。BFを境として心尖拍動波と不規則振動である心音の両者が強調されたSTFTとなった。このため、14~17HzのBlueとYellowが主体の帯域の中にBFが出現した。
APW×PCG-1とAPW×PCG、これら二つの線状描記ラインから、心尖拍動の高調波と心音の不規則振動の変動挙動が検出され、心尖拍動の高調波消失の周波数と心音出現周波数、およびBFを同定することができた。
次に、図22~図27の各図の考察のための着眼点を説明する。
図22~図24は計測時間10秒間の平均心拍数が93/minの被験者のSTFT解析結果を示す。
図22(a)は、APWが1.6Hzを基本調波とした9次までの高調波を成分にもつことを示す。本図のBFはパワースペクトルの急変部をbreakpointとするため、BFはBlueで示される14~17Hzの間に存在する。
図22(b)は、10~50Hzの間でパワースペクトルの変化が大きい振動成分をもつのがPCGで、図22(a)と併せて見ると5次の高調波から不規則振動までを捉えていることが分かる。
図23(a)はAPW×PCG-1で、BFはパワースペクトルの変動が少ない帯域での急変部をbreakpointとするため、BFはRed、Blue、Redで示される14~17Hzの間に存在する。図23(b)はAPW×PCG-1の拡大図で15Hz近傍にbreakpointを示す線状のRed lineが出現した。
図24(a)のAPW×PCGは、1.49Hzの基本調波がBlue Lineとなり、2次高調波がGreen Lineで、3次高調波から15Hzまでの高調波がRed lineとなり、基本調波より高調波の方がパワースペクトルが大きくなっていることが分かる。そして15Hz近傍のRed-Blue-Redで示されたBlue範囲にbreakpointがあり、Redが始まる15~20Hzの間から不規則振動があることを示した。
図24(b)は、APW×PCG-1、APW、PCG、APW×PCGのパワースペクトルを示す。APW×PCG-1にbreakpointがあり、breakpointは15Hzとなり、BFは15Hzに生じた。
図25~図28は、上記で説明したAPW×PCG-1とSTFTを用いたBF抽出法から求めたBFである。被験者の心拍数は58、71、79/minである。
図25(a)は、図22(a)の事例とは逆でbreakpointが山型のパワースペクトルの中に現れた事例で、図25(b)においてBFはRed-Green-Red-Blueで挟まれた帯域の中で線状のRedで現れた。
図26(a)は、図22(a)に示されるbreakpointが山型のパワースペクトルの稜線部に現れた事例で、図26(b)においてBFはRed-Blue-Red-Blueで挟まれた線状のRedで現れた。
図27(a),(b)は、パワースペクトルのゆらぎが変化した事例で、図27(a)の周波数解析結果からは判定が難しいが、図27(b)のSTFT法であればBFを判定できることがわかる。図27(b)では、BFはRed-Blueの境界線に線状のRedで現れた。
図28(a),(b)は、高調波消失周波数からのBF抽出法により求めると共に、STFTを用いたBF抽出法によっても確認した50名の被験者のBFと心拍数との相関図(BF-HR相関図)である。図28(a)は、BFの周波数を縦軸に、HRを横軸にとって示したもので、両者の関係を2次関数で近似すると、その関係は、y = 0.0173x2 - 2.5847x + 107.6111で表され、HRが75/minのときに10Hz近傍で最小値をとる。このときR2 = 0.8242となった。本図を用いると、PCGを計測する必要が無く、BFを知ることができる。
また、BFは年齢、SBP、DBP、BMIとの相関はなかった。さらに、20名の若年対照群におけるBFの検討において、STFTによって求められたBFと BF-HR曲線によって決定されるBFの間に強い相関がみられた(R2=0.85)。すなわち、本実施形態の生体信号検出センサ1(4SR)を用いて、体幹から出力される微小音響・振動(micro acoustic vibration (MAV))を採取しているため、生体信号分析装置100により分析することで、従来特定が困難であった「胸壁の厚い集団」の被験者からもBF、心尖拍動波の特定や、波形分類(正弦波ないし三角波)を行うことができる。
図28(b)は、高調波消失次数とHRから相関図を求めたものである。次数を縦軸に、HRを横軸にとり、両者の関係を2次関数で近似すると、その関係は、y = 0.0203x2 -3.2258x + 135.4587で表され、HRが79/minのときに10Hz近傍で最小値をとる。このときR2 = 0.9479となった。なお、図28(a),(b)は、安静状態でRR変動が15%以内の計測時間10秒間のデータで構成されている。
よって、図28(a),(b)の相関図及びその2次関数を相関データ記憶部150に予め記憶させておくことにより、測定時状態推定手段130が、心拍数から速やかに測定時の境界周波数を求めることができ、その境界周波数に応じて、測定時の人の心臓の状態や健康状態を速やかに判定し、出力することができる。
(心拍変動のゆらぎ特性と心拍数およびBFとの関係)
APW計測実験のECGから求めた心拍変動のゆらぎ特性と心拍数の相関を50人の被験者において検討した。図29(a)に示すように両者の関係は2次関数で近似され、HRが72/minのときに1/fゆらぎ(|n|=1)となった。一方、HRが50または90/min付近のときは、1/f2ゆらぎ(|n|=2)となった。
図29(b)は、心拍変動のゆらぎ特性とBFの相関を検討したものである。両者の関係は1次関数で近似された。BFが8~12Hz付近において1/fゆらぎ(|n|=1)となり、BFがおよそ20Hz以上で1/f2ゆらぎ(|n|=2)となった。
図29(c)~(e)は、同一被験者で90秒間の実験データを9等分して、安静時呼吸抑制時下で心拍変動が生じたときのBFの変動を求め、BF-HR相関図の回帰曲線上に描記したものである。3名の被験者は回帰曲線上を心拍変動に伴いBFが移動した。
安静状態で心拍数が70~80/minの場合に心拍変動のゆらぎが1/fとなっており、これをBFに照らすと、心拍変動のゆらぎが1/fとなるのは上記のBFが8~12Hz付近のときであり、BFがこの付近の値を示すときに心臓が効率の良い運動を行っていると推定できる。よって、図29(b)の相関図を相関データ記憶部150に記憶させておくことにより、測定時状態推定手段130が、境界周波数特定手段120により得られたBFから速やかに測定時の心拍変動のゆらぎ特性を推定できる。このとき、ゆらぎ特性が1/fと推定できれば、心臓の動きがよいこと、すなわち心臓の動きに関して健康であることが推定でき、BFの値が1/fゆらぎを示す8~12Hz付近から外れるほど、加齢の影響、疾病の影響等が大きいと推定できる。
(第2の実施形態)
次に、本発明の第2の実施形態について図30に基づいて説明する。本実施形態では、生体信号分析装置100が、上記の第1の実施形態の構成に加え、心尖拍動波形抽出手段160、波形分類手段170を有している。
心尖拍動波形抽出手段160は、生体信号データを、境界周波数特定手段120により特定される境界周波数(BF)を上限値に設定してフィルタリングし、心尖拍動により生じる振動の波形を求める。本実施形態では、上記のように0.5Hz~BFのバンドパスフィルターを用いてフィルタリングしている。なお、第1の実施形態の図15~図18において、心尖拍動波形(CAB)及び心音波形(CAS)を描記している。その際の心尖拍動波形は、本実施形態の心尖拍動波形抽出手段160と同様の手法により、0.5Hz~BFのバンドパスフィルターにより求めたものである。また、心音波形は、BF~50Hzのバンドパスフィルターを適用して求めている。
波形分類手段170は、心尖拍動波形抽出手段160により出力される心尖拍動により生じる振動の波形(心尖拍動波形)を分類する。本実施形態では、(A)フーリエ級数展開を用いた数学的アプローチにより前記波形を分類する手段と、(B)工学的アプローチにより得られるデータから波形を分類する手段とを併用して分類する。詳細は後述するが、2つのアプローチを併用して、本実施形態では、健常者の心尖拍動波形を5つに分類している。
また、工学的アプローチにより得られるデータから波形を分類する手段では、次の4つのうちの少なくとも1つの判定要素を用いて波形を分類する。分類精度の向上のためには、好ましくは、2つ以上の判定要素を組み合わせて分類する。
(1)前記心尖拍動により生じる振動の波形の所定時間範囲における時間波形の形状を示すデータ(図32~図36の各(b)の図)
(2)前記時間波形の前記周波数解析手段により得られる周波数とパワースペクトルを横軸と縦軸にとったグラフのデータ(図32~図36の各(c)の図)
(3)前記時間波形の前記短時間フーリエ変換の解析結果である時間、周波数及びパワースペクトルの変動の程度を示す所定時間範囲における画像データ(図32~図36の各(d)の図)
(4)コレログラムにより抽出される波形の種類と周期に関する所定時間範囲におけるデータ(図32~図36の各(e)の図)
上記のうち、まず、健常者の時間波形が、フーリエ級数展開を用いた数学的アプローチにより5つの波形に分類されることを求め、その実証として、工学的アプローチにより、振幅、周波数、時相を分析して、5つの波形に対応するデータを、Front CAB(0.5-BF)の時間波形の形状データ(図32~図36の各(b)の図)、パワースペクトルを横軸と縦軸にとったグラフのデータ(図32~図36の各(c)の図)、STFTの解析結果の図形データ(図32~図36の各(d)の図)、コレログラムにより抽出されるデータ(図32~図36の各(e)の図)から求めている。なお、この波形の評価法についてはさらに後述する。
波形分類手段170により分類された波形分類結果の情報は、相関データ記憶部150に記憶される。ここで、波形分類結果の情報は、健康状態(心疾患等の疾病の有無、疾病の特定、体調の良否等)に関する情報に紐付けられて記憶されている。これは、できるだけ多くの人の心尖拍動波形を測定し、その時の健康状態に対応させてデータ化されたものである。
測定時状態推定手段130は、心尖拍動波形抽出手段160により求められる、測定対象者の測定時の心尖拍動により生じる振動の波形について、波形分類手段170を用いて分類させ、当該測定時の波形分類結果を、相関データ記憶部150にアクセスして、予め記憶された波形分類結果の情報に照合する。それにより、測定時の波形分類結果に対応する健康状態が出力される。
なお、図31に示したように、心尖拍動波形抽出手段160により出力され、波形分類手段170において分類対象となる心尖拍動波形の波形情報と健康状態とを教師データとして用い、波形情報から健康状態を推定する推定モデルを機械学習により生成するモデル作成手段180を有する構成とすることもできる。また、波形分類手段170による分類結果の情報(分類情報)を加味して機械学習を行わせることにより、推定モデルの精度が向上する。
この場合、測定時状態推定手段130は、測定対象者の測定時の心尖拍動により生じる振動の波形情報を入力として、モデル作成手段180により作成された推定モデルを用い、取得した測定時の波形情報から、測定時の健康状態の値を出力する。このような機械学習を組み込んだ構成とすることにより、波形情報からの健康状態の推定精度が向上する。なお、学習対象の波形情報の特徴量は、図32~図36に示した各波形の構築している振幅、周波数、時相等となるが、図32~図36の各(c)に示したSTFTの解析結果としての画像データは情報量が多く、学習対象の特徴量を抽出しやすい。よって、波形情報の把握には、STFTの画像データを学習対象に含めることが好ましい。また、図32~図36の各(c)の上図及び下図は、0.5Hz~BFの心尖拍動波形を、1.95Hzと7.81Hzの2つに切り分けたものであり、これにより学習対象の特徴量がより増加し、推定精度の向上に寄与できる。
なお、第2の実施形態のその他の構成は、生体信号検出センサ1を含め、第1の実施形態と同様である。
以下、波形の分類、評価する方法についてさらに具体的に説明する。
[Front CABの時間波形の評価法]
(Front CABの時間波形分析)
次に、Front CAB波形の評価法について説明する。
心血管系の動態について非侵襲的に評価し、病態及び治療の指標となる理学的診断ツールの一つに上記の心尖拍動図(ACG)が知られているが、診察の体位は仰臥位が絶対的なものではなく、座位や左側臥位がとられる場合がある。そこで、仰臥位、左側臥位、座位を網羅した、計測体位別心尖拍動正常パターンに関する情報があると、PCG、ECGの波形分析や時相分析との関係について検討することができる。
そこで、胸壁上に置かれたマイクロフォンで測定されたF-APW波形を周波数解析し、0.5HzからBFの周波数帯域を取り出したもの、すなわちFront CAB(0.5-BF)を生成し、臨床的に知られる心尖拍動波のものと比較する。波形分析の着眼点としては、視診でも触診でも収縮早期のみに認められる持続の短い拍動であるE波のtappingと拍動の持続が長く力強いsustainedに注目する。ここではFront CAB(0.5-BF)の波形の形状とACGのピークが出現する時相を用いて波形分析を行う。
Front CAB(0.5-BF)での各時間波形の振幅の評価は、収縮期の時相で行う。左心房の収縮に由来するA波並びにA波の痕跡をFront CAB(0.5-BF)で確認する。A波の後のボトム値をACGのC pointとし、ACGのE波に相当するピーク値までの全振幅を計測し、収縮期波の振幅とする。このACGのC pointとE波に至る全振幅を、Front CAB(5-BF)、Rear CAB(0.5-BF)、Lumbar CAB(0.5-BF)の時間波形 についても計算する。計測部位別の時間波形の振幅比較を行う場合は、Front CAB(0.5-BF)を基準値にして行う。各波形毎の比をもって振幅比とする。時間波形の振幅のボトム値は心電図のP波とR波の間に存在する。時間波形の振幅のピーク値はR波以降、PCGのII音の開始点までに存在する。次に、被験者50名分のFront CAB(0.5-BF)の時間波形の分類を行う。
(短時間フーリエ変換を用いた高調波成分分析法)
次に、短時間フーリエ変換を用いた波形勾配分析法と高調波成分分析法について説明する。
Front CAB(0.5-BF)の時間波形を構成する高調波を対象にして、式(4)を用いて時間波形のどの部位に高調波の成分が重畳されているかを検討した。周波数は、1~4Hzと4~16Hzに分けて分析する。1~4Hzの波形分析のために窓点数は、29=512とした。512点はサンプリング周波数1000Hzにより0.512秒となり、計算上では、1.95Hzのパワースペクトルの成分が出現する。ハニング窓関数の移動幅は51点に設定する。その時には、周波数情報は51/1000=0.051秒毎に作り出される。
高調波によって構成されるE波を捉えるための窓点数は、27=128とした。128点はサンプリング周波数1000Hzにより0.128秒となり、7.8Hzでパワースペクトルの成分が出現する。分解能を0.12秒間(7.8Hz)とした理由は、心拍数90/minの場合の収縮期の持続時間を0.26秒間と想定し、収縮期の持続時間のピークを捉える検出周波数の必要分解能を0.13秒間としたことによる。この場合には、ハニング窓関数の移動幅を12点に設定した。その時には、周波数情報は12/1000=0.012秒毎に作り出される。
(自己相関関数を用いた振動の峻別)
自己相関関数(Autocorrelation function:ACF)を用いてCorrelogramを描き、波形の峻別を行う。Correlogramの縦軸はACFで、横軸は時間である。それぞれの時間間隔TについてのACFを計算し、Correlogramから波形の種類と周期を捉える。なお、ACF、すなわちR(τ)は、
Figure 2022071812000010
で表される。τは時間間隔、ラグ (Lag)を示し、Tは全区間を示す。
正規化自己相関関数は、τ=0で最大値1をとり、その値は2乗平均値となる。ACF値が0.2以下の場合、自己相関を持たないことに対する95%以上の信頼区間にある。振動の種類分析は、10秒間のデータで判定する。
[Front CABの時間波形の評価結果]
座位における正常パターンを50人の被験者から計測されたFront CAB(0.5-BF)を用いて、分類と検定を行った。図32~図36がその結果を示す。Front CAB(0.5-BF)の波形は、図32~図36の各(b)の図に示され、図32~図36の各(a)の図に示したECGの1秒間の時間波形に対応する。Front CAB(0.5-BF)の波形は、5つの波形に分類さる。心尖拍動図のA波(心房収縮波:矢頭)及びE波(*印)に始まる心室収縮期波に類した波形が認められる。
5つの波形パターンの分類手法はフーリエ級数展開に基づくものである。ただし、数式表現は直観的な理解には不向きである。そこで、Front CABの波形の考察にあたり、振幅に注目する時にはPSD(パワースペクトル密度)、振幅と位相に注目する時にはSTFT、周期性を持つかどうかを調べる時にはCorrelogramを用いることで、Front CABの波形の特徴を分類した。
本明細書では、分類した5つの波形を疑似正弦波(S波)、疑似三角波群(T波(三角波)、C波(余弦波)、G波(ガウス波))、疑似全波整流波(R波)と名付けた。このうちS波とG波は、収縮期波の持続時間から捉えると視診・触診で言われているsustained、T波、C波、R波はtappingが対応する。なお、これらの波形は、完全な正弦波、三角波、全波整流波ではないものの、その波形形状は類似していることから、「疑似」を接頭語に付している。
波形分析のために必要なものは、基本調波の波形とbottomとpeakの各近傍にある各波形の振幅と周波数、および位相である。
図32~図36の各(c)では、周波数解析によりパワースペクトル密度(PSD)を求め、振幅と周波数の関係を分類し、図32~図36の各(e)は、Correlogram により振動の種類を分類した。これらのパラメータにより、Front CAB(0.5-BF)の時間波形は5つの代表的な波形に分類された。
より詳細には、PSDによると以下のような波形パターンの分類が可能となる。
まず、Front CABのS波が示す収縮期波E波は、高周波成分の寄与が小さく正弦波に近くなる。すなわち、基本周波数成分のみが卓越する(図に示した「X」のパターン)。
次に、Front CABのT波・C波・G波が示す収縮期波E波は、三角波に収束していく過程の中にあるもので、三角波として分類できる。すなわち、奇数次の周波数成分が卓越し、かつ、その振幅は、次数が高くなるにつれて、次数の二乗の逆数に比例して小さくなる(図に示した「Y」のパターン)。
Front CABのR波は、基本周波数成分の振幅が小さく、その整数倍の周波数成分の振幅との差が少ないため、三角波の部分が正弦波の全波整流波形に置き換わる。対応するPSDは、偶数次の周波数成分が卓越し、かつ、高次になるにつれて振幅は小さくなる(図に示した「Z」のパターン)。
STFTは、位相情報の影響を調べる際に適しており、PSDとSTFTを重ね合わせることにより、波形分類の精度が向上する。これにCorrelogram を用いて、A波とE波の波形間の周期性を検討することで、正常な波形のバリアントとして、任意の波形を振幅、位相、周期性の観点から分類することができる。
図32~図36の各(d)の図のSTFT図は、図32~図36の各(c)に現れる基本調波から5Hzまでのパワースペクトル情報を、中心周波数1.95Hzと7.81Hz の2つの周波数で切り分けて、パワースペクトルに関する情報量を2倍にして、パターンマッチングを行いやすくしたものである((d)の上図が1.95Hz、下図が7.8Hz)。図32~図36の各(c)の基本調波並びに1~5Hz間の高調波の情報は、STFTの1~5Hz間に描記される。なお、時間波形の情報はRed、Green、Blueによる等高線で表した。
波形I(S波)は基本調波の出現する時間帯が3箇所あり、そして5~10Hzの成分が出現する時間帯が2箇所あり、心房性隆起波と収縮期波の出現間隔が短く、矩形波であることを示す。波形II(T波)は5Hz以下の成分の出現を示すSTFTが2.3~2.5秒間に集中し、5~10Hzの成分は2.5~2.7秒間に集中し、おおよそ収縮期の時相にこれらの成分が収まり、三角波であることを示す。同様に、波形III(C波)も収縮期の時相に収まり、心房性隆起波と収縮期波までの間隔が短く、三角波であることを示す。波形IV(G波)は波5~10Hzのピーク値間隔が波形I(S波)に近くなり、S波とT波を組み合わせたような形状で、全体としては三角波の特徴を示す。波形V(R波)は5Hz以下のSTFTのピーク値が右下がりの等高線を示し基本調波から高調波のパワースペクトルが同じ大きさにあることが示され、さらに5Hz以上10Hz以下のピーク値が3点以上あり、そして出現の時相が5Hz以下の成分と同一となり、全波整流波の特徴を示した。
図32~図36の各(e)のCorrelogram により振動の種類(A波とE波の波形間の周期性)を自己相関係数が0.5以上、0.2以下を基準値とし、正弦波のノッチの出方で、分類した。そこで、Front CAB(0.5-BF)の Correlogramの波形は、正弦波であるi、正弦波にノッチが1つ入る波形であるii、2つ以上入る波形であるvが抽出される。残りの波形は、包絡線が片減衰のiiiと両減衰のivに分けられた。
[波形分類の理論的根拠と客観的分類基準]
ここで、フーリエ級数展開を用いた心尖拍動波への数学的アプローチについて詳述する。これらの数学的アプローチにより健常者の心尖拍動波が5つに分類されることが示され、上記のPSD、STFT及びCorrelogramはこの5つに分類を工学的アプローチによる裏付けとなっている。
まず、任意の周期関数をf(t)とすると、f(t)は、次のフーリエ級数に展開できる。
ただし、周期をTとして、ω=2π/Tである。
Figure 2022071812000011
S波は、その波形形状が正弦波に近いため、式(7)で、1次成分が卓越する状況にあたる。ただし、若干の高調波成分も含まれているため、STFT解析結果からは、基本周波数が主体となり、これに高次周波数成分が伴ったものと解釈できる。
これは、収縮期の持続時間が長い、Sustainedの所見に合致する。
T波は、三角波に近いため、以下のようにフーリエ級数に展開できると考えられる。
Figure 2022071812000012
式(8)よりnの次数が高くなるにつれて、急激に展開係数が小さくなってゆくという特徴があることが分かる。
また、スペクトルとして発現する次数は、奇数次であることが分かる。
これによってPSDの形が説明できる。
Correlogramで示される自己相関係数がt=0でピークを持ち、それ以外の時間では自己相関係数び絶対値が1より小さくなることは、正弦波とは別の高調波の波形が影響を及ぼしていることから説明可能である。
これは、正弦波に、収縮早期の持続の短い尖鋭な心室収縮に由来する波形の重畳を意味し、収縮早期の持続時間が短いTappingの所見に合致する。
C波の心室収縮に由来する波形についてもT波と同様な考察が可能である。T波との違いは、式(8)のフーリエ級数展開において、C波は、T波よりも1次成分の卓越度が高いことである。
STFT解析の結果を見ると、基本周波数と第2高調波が主体となり、第7高調波までの周波数成分で心室収縮に由来する波形が形作られ、Correlogramは、T波と同様に、収縮早期の持続時間が短いTappingの所見に合致する。
G波は、STFT解析から、S波にT波を組み合わせたものとして考察でき、1次成分が卓越する。T波よりもS波に近いため、収縮期の持続時間が長く、Sustainedの所見に合致する。
R波の心室収縮に由来する波形は、波形自体は、単峰ではなく多峰であるが、全波整流波形と類似しているとの仮説を立てると、次のようなフーリエ級数に展開できる。
Figure 2022071812000013
式(9)は、2次成分から始まり、これに偶数次の高調波が重畳したものと考えられる。nが大きくなるとともに展開係数は小さくなる。なお、R 波には1次成分も含まれているため、PSDの図に示すような形をとるものと考えられる。
Correlogramの分析から、T波と同様に1次成分以外が卓越している理由が分かる。これは収縮早期の持続時間が短いTappingの所見に合致する。
以上のとおり、Correlogram からは正弦波の基本周期を持つS波とG波に分類され、高調波が顕在化するT波、C波、R波に分類される。また、ラグに対する自己相関関数の差に注目すると、Correlogramの+側でS波とG波が区別でき、Correlogramの-側でT波とC波が区別でき、多峰性を示すノッチによりR波が区別できる。
すなわち、Correlogramによる分類は、S波とG波は正弦波的な周期性に着目し、他方T波、C波、R波はひずみ波的な周期性に着目したものである。収縮期波の周期性はCorrelogramの+側の自己相関関数で分類され、心房性隆起波と他の波形の周期性はCorrelogramの-側の自己相関関数で分類される。S波とG波は、ラグに対する自己相関の差からさらに仕分け、T波、C波、R波は-側の自己相関係数の値(上記の例では、0,-0.2,-0.5)で仕分けることができた。そして心房性隆起波と他の波形の凸凹は描かれた自己相関関数の形(ノッチ)を利用した。R波は多峰であるため、ノッチが2つ以上となった。
図32~図36の各(b)~(d)の図を用いて、50事例をパターンマッチングして波形を検定した。その結果、疑似正弦波(S波)、三角波(T波)、余弦波(C波)、ガウス波(G波)、疑似全波整流波(R波)の発生比率は、それぞれ4/50=8%、8/50=16%、15/50=30%、11/50=22%、12/50=24%となった。
図37は、図32~図36の各(a)の図で決められた波形検定を用いずに、図32~図36の各(b)~(e)の図に示されているパワースペクトル、2種のSTFT、Correlogramを用いて波形検討を行った結果を示す。図中の太字は、波形検定ができた数値を示す。STFTとCorrelogramで、50件中、三角波の1件以外の43件は波形検定ができた。パワースペクトルは、50件中38件の波形検定は可能となった。なお、P値はSTFT検定でP=2.2×10-16、パワースペクトル検定でP=7.9×10-5、Correlogram検定でP=2.2×10-16となった。Wiener-Khintchineの定理に従うと、STFTとCorrelogram (ACF) の両者は逆フーリエ変換の関係にあることから、両者のP値が一致したことは、解析手法の妥当性を支持する結果となった。したがって、この3つの検定を合わせて行うことで、Front CAB(0.5-BF)の波形検定の精度向上が可能となり、従来のACGの波形の形状により行っていた評価を、高い精度で分類できる可能性が見出され、4SRセンシングシステムが、正常心と肥大・拡大心の鑑別等の異常所見を見つけるために必要な波形解像度を有していることがわかった。
よって、本実施形態の測定時状態推定手段130は、波形分類手段170により得られる測定対象者の心尖拍動波の分類結果が、上記の健常者の波形分類である5つのパターンに対応する場合には、健常と判断してその結果を出力し、それらのいずれにもあてはまらない場合には、心疾患に関する異常所見を有していると判断し、その結果を出力する。
図38(a)~(e)には、図32~図36の各(a)の図で定義された5つの波形の代表事例として5人の被験者のFront CAB(0.5-BF)、Front CAB(5-BF)、Rear CAB(0.5-BF)、Lumbar CAB(0.5-BF)の周波数解析結果を示す。図39に計測部位別のCABの振幅比の平均値と標準偏差と周波数解析結果を示す。振幅比はFront CAB(0.5-BF)に対する計測部位別CABの振幅の大きさの違いから計測の容易性を表す指標で、上記の(Front CABの時間波形分析)の項で示した収縮期波の左室内圧上昇開始点から心室収縮波の最大振幅値を用いて計算する。Front CAB(5-BF)の振幅比はFront CAB(5-BF)とFront CAB(0.5-BF)の比から求められる。同様にRear CAB(0.5-BF)の振幅比はRear CAB(0.5-BF)とFront CAB(0.5-BF) の比、Lumbar CAB(0.5-BF)の振幅比はLumbar CAB(0.5-BF)とFront CAB(0.5-BF)の比から求められる。ここにFront CAB(0.5-BF)に対してFront CAB(5-BF)は0.5Hz~5Hzの周波数成分を含まないため振幅比は0.2となった。Rear CAB(0.5-BF)は振幅比が1.1となり、Lumbar CAB(0.5-BF)は振幅比は2.6となった。平均値と標準偏差は50人の被験者の振幅比を用いて計算した。
次に、図38(a)~(e)のパワースペクトルと図39の振幅について検討する。Rear CAB(0.5-BF)とLumbar CAB(0.5-BF)では、5Hz近傍のパワースペクトルが増大した。具体的には、Rear CAB(0.5-BF)は5Hz近傍が1.8~4.8倍となり、Lumbar CAB(0.5-BF)は5Hz近傍が3.7~7.9倍となった。5Hz近傍の高調波のパワースペクトルの差は、大動脈の長さに起因するもので、大動脈と圧力脈動による液柱共振現象が生じたことが示唆された。
以上より、本発明によれば、具体的に次のような効果が得られた。
従来の3SRは車両積載を目的として開発したセンシングシステムで、外部振動をエアパックの持つエアダンピング効果により減衰させ、20Hz以上の周波数成分を含む心音のI音、II音と心拍変動の周期成分を取り出していた。一方、本発明の4SRは室内で使用されることを目的として開発したもので、1~2Hz近傍にある心尖拍動の周期と振幅を捉えるためのセンサとして適していることがわかった。
F-APWは、心尖拍動波と近似した波形となった。そして、F-APWから心尖拍動波と心音の各波形を生成できた。R-APWとL-APWは、圧力脈動による共振現象が起因して異なる波形になったと考えられる。但し、F-APWと異なる波形となったR-APWとL-APWも、振動解析を適用することで心尖拍動を含む心臓の運動に関する振動情報(CAB)と心音(CAS)を抽出でき、臨床的意義のあるものとなる可能性が示された。
F-APWの高調波消失周波数からのBF抽出法と STFT法を用いたBF抽出法を組み合わせることにより、BF-心拍数の相関図が作成され、実際の計測結果からBFは心拍数の2次関数で表すことができ、PCGを用いることなく心拍数からリアルタイムでF-APWのデータのみでBFを決定することが可能となった。そしてSTFT法による可視化法を用いると心臓の振動や大動脈の圧力脈動の振動分析のみならず、変化の度合いを活用した診断ができる可能性が示唆された。
心拍変動のゆらぎの実験結果から安静状態で心拍数が70~80/minの場合に1/fとなり、効率の良い運動を行っていることが示唆された。心拍変動のゆらぎの実験結果からBFが10Hz付近のときに心拍変動のゆらぎが1/fとなり効率の良い運動を行っていることが示唆された。この特性を活用すれば、薬物が循環器機能に及ぼす変化を定量化することができる可能性がある。
上記のように心拍数(この場合は10秒間)からBFが決定される。そのBFから得られる心拍変動のゆらぎ(90秒間)の値が、心拍数(90秒間)から直接得られる心拍変動のゆらぎ(90秒間)の値と傾向が一致した。これらの関係式を代入して求めた結果は、BFの妥当性を支持するものとなった。
さらに、APWに対して0.5HzからBFの帯域でバンドパスフィルター処理をすることで、ACG近似の波形が描記でき、周波数解析、STFT法を用いた時系列解析法、Correlogramを用いた波形分析法から心尖拍動を表すFront CAB(0.5-BF)を大きく5つの波形に分類することができた。
安静時呼吸抑制時下で心拍変動が生じると、BF-HR相関図の回帰曲線上をBFが移動する(図40参照)ため、上記で説明した振動解析法とBF-HR相関図から心臓の運動解析と圧力脈動、および心臓の運動の様子や心室内血流を評価できる可能性が示唆された。
また、仰臥位および左右側臥位でも座位と同じBF-HR相関図が得られた(図41参照)。したがって、BF-HR相関図は姿勢にとらわれることなく、一つの線図で推定できる可能性が示唆され、4SRを用いた計測法は容易にデータ取得が可能であるため実用性が高く、視診・触診所見の客観化への応用が可能となり、APWは心血管系の動態について非侵襲的に評価し、病態及び治療の指標となる理学的診断ツールとなる可能性が出てきた。
本発明によれば、二つの周波数帯にまたがる信号をBFを特定して分離している。
すなわち、APWから心尖拍動と心音に対応する信号が分離できたことにより、1つの生理的信号を観測することで、心臓活動の異なる種類の動作の結果を定量的に認識することが可能となり、心臓状態に関する有用な相補情報を得ることができる。また、信号分離によって波形分解能が高くなることで、例えば、二つの周波数帯の境界となる周波数とその近傍のFront CABの高調波成分とFront CASの卓越周波数成分が、心循環系の新たなパラメータとなる可能性がある。
また、Front CABは心尖拍動図の収縮期波の圧力情報を有し、Front CASは心周期に関する時相情報を有し、圧力と時相情報を組み合わせることにより、相補情報から診断目的に使えると考えられる。
例えば、Front CASから心音I音が特定でき、収縮期波の上行脚のI音に対応するノッチは、図32~図36の各(b),(c)に示したSTFTで広い範囲でも狭い範囲でも観察できる。0~20Hz帯域のウインド幅を変更することにより観察できる範囲が変えられ、コントラストを調整することでノッチが判別し易くなる。収縮期の圧力波のノッチの有無の判定やこのノッチが心臓活動の異常原因の特定に使え、診断に活用できる。心臓の振動情報であるFront CABと心臓の音響情報であるFront CASは、主に心房・心室・弁・血液・大動脈の力学的特性を示す。そしてその相補性の結果は、血圧や心拍数、および呼吸数のゆらぎとなって表れる。また、薬物治療で、患者に投薬する間に起きる変化を定量化することもできる。例えば左心室不全の診断において、警戒すべき左心室の力学的特性を同じ年齢層の健常人のものと比較することでも左心室不全の存在を診断する際に役立つ。さらに、導通経路や速度の特徴、および電位などのECG情報と関連付けられることにより、心筋の細胞の健康状態の推定ができる可能性もある。Front CABとFront CASの情報から、心臓状態の監視や機能評価に関連するデータが提供できる可能性もある。その結果、心室の大きさ、心臓の壁の厚さ並びに心臓の力学的特性のモデル化ができ、運動が数式で記述できる可能性もある。
そして、大動脈に沿ってデータを計測することで、大動脈の血行動態の中で心尖拍動ないし心音に由来するものに分離でき、大動脈の血行動態に関する新たな知見を得ることができる可能性もある。
上記のことから明らかなように、0.5~80Hzの音響振動情報(Acoustic Pulse Wave :APW)を抽出し、所定のフィルタ処理を適用することで、ACG相当波形(Cardiac Apex Beat: CAB)と心音相当波形(Cardiac Acoustic Sound: CAS)を生成できることがわかる。
また、CABとCASは、発生源と伝達経路が異なり、これら二つの振動が重畳したものがAPWである。収縮能が心筋固有の特性で、前負荷や後負荷から独立していると考えると、収縮期のCABの波形は、収縮の際の化学的刺激やホルモンの影響に関する情報を有する。収縮期に心室から駆出される血液の運動情報を有するものとしてCASを考えると、APWは前負荷と後負荷、および収縮能の情報を有するものとなる。
よって、APW、CAB、CASには多くの音・振動・圧に関する情報があり、ECGの持つ分極現象と合わせて観察することで、心臓の故障診断ができる。
そして、APW・CAB・CASのリアルタイム計測により、臨床的意義のある音振動・圧力情報が加わり、機械の故障診断に用いられてきた工学的アプローチが健康診断へ適用可能であることもわかった。また、APWを計測する本発明の4SRシステムは長期にわたり繰返し非侵襲的に生体情報を取得することができるPhysical diagnostic tool となり、上記の機械学習を含むAIを使用した分析を含め広範な応用が期待される。
1 生体信号検出センサ(4SR)
10 三次元立体編物
20 収容フィルム
30 マイクロフォン
40 ケース
50 ゲル
100 生体信号分析装置
110 周波数解析手段
120 境界周波数特定手段
130 測定時状態推定手段
140 心拍数演算手段
150 相関データ記憶部
160 心尖拍動波形抽出手段
170 波形分類手段
180 モデル作成手段
1000 生体信号検出センサ(3SR)

Claims (31)

  1. 体表面を介して生体信号検出センサにより得られる生体信号データを周波数解析する周波数解析手段と、
    前記周波数解析手段から得られる前記生体信号データの周波数解析結果から、前記生体信号データ中、心尖拍動により生じる振動と心音により生じる振動との境界周波数を求める境界周波数特定手段と
    を有することを特徴とする生体信号分析装置。
  2. 前記境界周波数特定手段は、
    前記周波数解析結果中、調和振動と不規則振動との境界となるパワースペクトルの急変部を求め、この急変部を基準に前記境界周波数を特定する手段を含む請求項1記載の生体信号分析装置。
  3. 前記境界周波数特定手段は、
    前記周波数解析結果に、同時に測定した心音データの周波数解析結果を加味して、前記パワースペクトルの急変部を求める手段を含む請求項2記載の生体信号分析装置。
  4. 前記境界周波数特定手段は、
    前記生体信号データ及び前記心音データの各周波数解析結果をそれぞれ加算平均処理した波形を対数差分法を用いて両対数軸表示する両対数軸表示手段と、両対数軸表示した波形からゆらぎの変化点を求め、このゆらぎの変化点を前記パワースペクトルの急変部として特定する急変部を特定する手段と
    を有する請求項3記載の生体信号分析装置。
  5. 前記周波数解析手段として、短時間フーリエ変換が適用され、
    前記境界周波数特定手段は、前記短時間フーリエ変換の解析結果から前記パワースペクトルの急変部を求める手段を含む請求項2又は3記載の生体信号分析装置。
  6. 前記周波数解析手段は、前記短時間フーリエ変換の解析結果を、時間、周波数及びパワースペクトルの変動の程度を示す画像データに出力する手段を含み、
    前記境界周波数特定手段は、前記画像データから、前記パワースペクトルの急変部を求める手段を含む請求項5記載の生体信号分析装置。
  7. 前記境界周波数に関連する相関データが記憶された相関データ記憶部と、
    前記相関データに照合して、測定対象者の測定時の健康状態を推定する測定時状態推定手段と
    をさらに有する請求項1~6のいずれか1に記載の生体信号分析装置。
  8. 前記相関データが、前記境界周波数と心拍数との相関データであり、
    前記測定時状態推定手段は、測定対象者の測定時の心拍数を、前記境界周波数と心拍数との相関データに照合して、前記測定対象者の測定時の前記境界周波数を推定する手段を含む請求項7記載の生体信号分析装置。
  9. 前記相関データが、前記境界周波数と心拍変動のゆらぎ特性との相関データであり、
    前記測定時状態推定手段は、前記測定対象者の測定時の前記境界周波数を、前記境界周波数と心拍変動のゆらぎ特性との相関データに照合して、前記測定対象者の測定時の心拍変動のゆらぎ特性を推定する手段を含む請求項7記載の生体信号分析装置。
  10. さらに、前記生体信号データを、前記境界周波数特定手段により特定される前記境界周波数を上限値に設定してフィルタリングし、前記心尖拍動により生じる振動の波形を求める心尖拍動波形抽出手段と、
    前記心尖拍動により生じる振動の波形を分類する波形分類手段と
    を有すると共に、
    前記測定時状態推定手段は、前記心尖拍動波形抽出手段により求められる、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形を、前記波形分類手段により分類し、当該測定時の波形分類結果のデータをもとに、前記測定対象者の測定時の健康状態を推定する手段を含む請求項7記載の生体信号分析装置。
  11. 前記相関データ記憶部には、前記波形分類結果と健康状態との相関データが予め記憶されており、
    前記測定時状態推定手段は、前記心尖拍動波形抽出手段により求められる、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形を、前記波形分類手段により分類し、当該測定時の波形分類結果を、予め記憶された前記波形分類結果と健康状態との相関データに照合して、前記測定対象者の測定時の健康状態を推定する手段を含む請求項10記載の生体信号分析装置。
  12. 前記心尖拍動波形抽出手段により出力され、前記波形分類手段における分類対象となる波形情報と前記健康状態とを教師データとして用い、前記波形情報から前記健康状態を推定する推定モデルを機械学習により生成するモデル作成手段を有し、
    前記測定時状態推定手段は、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形情報を入力として、前記モデル作成手段により作成された推定モデルを用い、取得した測定時の波形情報から、測定時の前記健康状態の値を出力する手段を含む請求項10記載の生体信号分析装置。
  13. 前記波形分類手段は、次の(A)及び(B)の手段:
    (A)フーリエ級数展開を用いた数学的アプローチにより前記波形を分類する手段、
    (B)次の(1)~(4)の工学的アプローチにより得られるデータを1又は2以上組み合わせて前記波形を分類する手段、
    (1)前記心尖拍動により生じる振動の波形の所定時間範囲における時間波形の形状を示すデータ、
    (2)前記時間波形の前記周波数解析手段により得られる周波数とパワースペクトルを横軸と縦軸にとったグラフのデータ、及び、
    (3)前記時間波形の前記短時間フーリエ変換の解析結果である時間、周波数及びパワースペクトルの変動の程度を示す所定時間範囲における画像データ、
    (4)コレログラムにより抽出される波形の種類と周期に関する所定時間範囲におけるデータ、
    を用い、
    健常者の前記心尖拍動により生じる振動の波形を5つに分類し、
    前記測定対象者の測定時の前記健康状態を、前記5つに分類された波形分類結果のいずれかに対応するか否かにより、健康状態か否かを推定する請求項10~12のいずれか1に記載の生体信号分析装置。
  14. 前記画像データが、複数の周波数を基準として切り分けられ、所定時間範囲において複数生成されている請求項13記載の生体信号分析装置。
  15. 前記生体信号検出センサは、
    三次元立体編物と、
    前記三次元立体編物の周囲を密閉的に被覆する収容フィルムと、
    前記収容フィルムの外側に配設されるマイクロフォンと
    前記マイクロフォンをカバーするケースと、前記ケース内で前記マイクロフォンへの外乱の混入抑制機能を果たす外乱混入抑制部材と
    を有して構成される請求項1~14のいずれか1に記載の生体信号分析装置。
  16. 前記外乱混入抑制部材がゲルである請求項15記載の生体信号分析装置。
  17. 体表面を介して生体信号検出センサにより得られる生体信号データを処理し、コンピュータを、生体信号分析装置として機能させるコンピュータプログラムであって、
    前記生体信号データを周波数解析する手順と、
    その周波数解析結果から、前記生体信号データ中、心尖拍動により生じる振動と心音により生じる振動との境界周波数を特定する手順と
    を前記コンピュータに実行させるコンピュータプログラム。
  18. 前記境界周波数を特定する手順では、
    前記周波数解析結果中、調和振動と不規則振動との境界となるパワースペクトルの急変部を求め、この急変部を基準に前記境界周波数を特定する請求項17記載のコンピュータプログラム。
  19. 前記境界周波数を特定する手順では、
    前記周波数解析結果に、同時に測定した心音データの周波数解析結果を加味して、前記パワースペクトルの急変部を求める請求項18記載のコンピュータプログラム。
  20. 前記境界周波数を特定する手順では、
    前記生体信号データ及び前記心音データの各周波数解析結果をそれぞれ加算平均処理した波形を対数差分法を用いて両対数軸表示し、両対数軸表示した波形からゆらぎの変化点を求め、このゆらぎの変化点を前記パワースペクトルの急変部として特定する請求項19記載のコンピュータプログラム。
  21. 前記周波数解析する手順では、短時間フーリエ変換が適用され、
    前記境界周波数を特定する手順では、前記短時間フーリエ変換の解析結果から前記パワースペクトルの急変部を求める請求項18又は19記載のコンピュータプログラム。
  22. 前記周波数解析する手順では、前記短時間フーリエ変換の解析結果を、時間、周波数及びパワースペクトルの変動の程度を示す画像データに出力し、
    前記境界周波数を特定する手順では、前記画像データから、前記パワースペクトルの急変部を求める請求項21記載のコンピュータプログラム。
  23. 相関データ記憶部に記憶された前記境界周波数に関連する相関データに照合して、測定対象者の測定時の健康状態を推定する手順を実行する請求項17~22のいずれか1に記載のコンピュータプログラム。
  24. 前記相関データが、前記境界周波数と心拍数との相関データであり、
    前記測定時の状態を推定する手順では、測定対象者の測定時の心拍数を、前記境界周波数と心拍数との相関データに照合して、前記測定対象者の測定時の前記境界周波数を推定する請求項23記載のコンピュータプログラム。
  25. 前記相関データが、前記境界周波数と心拍変動のゆらぎ特性との相関データであり、
    前記測定時の状態を推定する手順では、前記測定対象者の測定時の前記境界周波数を、前記境界周波数と心拍変動のゆらぎ特性との相関データに照合して、前記測定対象者の測定時の心拍変動のゆらぎ特性を推定する請求項23記載のコンピュータプログラム。
  26. さらに、前記生体信号データを、前記境界周波数特定手段により特定される前記境界周波数を上限値に設定してフィルタリングし、前記心尖拍動により生じる振動の波形を求める手順と、
    前記心尖拍動により生じる振動の波形を分類する手順と
    をコンピュータに実行させ、
    前記測定時の状態を推定する手順では、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形を、前記波形を分類する手順の実行により分類し、その波形分類結果をもとに、前記測定対象者の測定時の健康状態を推定する請求項23記載のコンピュータプログラム。
  27. 前記相関データ記憶部には、前記波形分類結果と健康状態との相関データが記憶されており、
    前記測定時の状態を推定する手順では、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形を、前記心尖拍動により生じる振動の波形を分類する手順の実行により分類し、その測定時の波形分類結果を、前記波形分類結果と健康状態との相関データに照合して、前記測定対象者の測定時の健康状態を推定する請求項26記載のコンピュータプログラム。
  28. 前記心尖拍動波形を抽出する手順の実行により出力され、前記波形を分類する手順における分類対象となる波形情報と前記健康状態とを教師データとして用い、前記波形情報から前記健康状態を推定する推定モデルを機械学習により生成し、
    前記測定時の状態を推定する手順が、前記測定対象者の測定時の前記心尖拍動により生じる振動の波形情報を入力として、前記推定モデルを用い、取得した測定時の波形情報から、測定時の前記健康状態の値を出力する請求項26記載のコンピュータプログラム。
  29. 前記波形を分類する手順は、次の(A)及び(B)の手順:
    (A)フーリエ級数展開を用いた数学的アプローチにより前記波形を分類する手順、
    (B)次の(1)~(4)の工学的アプローチにより得られるデータを1又は2以上組み合わせて前記波形を分類する手順、
    (1)前記心尖拍動により生じる振動の波形の所定時間範囲における時間波形の形状を示すデータ、
    (2)前記時間波形の前記周波数解析手段により得られる周波数とパワースペクトルを横軸と縦軸にとったグラフのデータ、及び、
    (3)前記時間波形の前記短時間フーリエ変換の解析結果である時間、周波数及びパワースペクトルの変動の程度を示す所定時間範囲における画像データ、
    (4)コレログラムにより抽出される波形の種類と周期に関する所定時間範囲におけるデータ、
    を用い、
    健常者の前記心尖拍動により生じる振動の波形を5つに分類し、
    前記測定対象者の測定時の前記健康状態を、前記5つに分類された波形分類結果のいずれかに対応するか否かにより、健康状態か否かを推定する請求項26~28のいずれか1に記載のコンピュータプログラム。
  30. 前記画像データが、複数の周波数を基準として切り分けられ、所定時間範囲において複数生成されている請求項29記載のコンピュータプログラム。
  31. 前記請求項17~30のいずれか1に記載のコンピュータプログラムが記録された記録媒体。
JP2021107109A 2020-10-28 2021-06-28 生体信号分析装置、コンピュータプログラム及び記録媒体 Pending JP2022071812A (ja)

Priority Applications (3)

Application Number Priority Date Filing Date Title
US18/250,827 US20230389814A1 (en) 2020-10-28 2021-10-28 Biological signal analysis device, computer program, and recording medium
PCT/JP2021/039924 WO2022092243A1 (ja) 2020-10-28 2021-10-28 生体信号分析装置、コンピュータプログラム及び記録媒体
EP21886360.3A EP4238491A1 (en) 2020-10-28 2021-10-28 Biological signal analysis device, computer program, and recording medium

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2020180963 2020-10-28
JP2020180963 2020-10-28

Publications (1)

Publication Number Publication Date
JP2022071812A true JP2022071812A (ja) 2022-05-16

Family

ID=81593859

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2021107109A Pending JP2022071812A (ja) 2020-10-28 2021-06-28 生体信号分析装置、コンピュータプログラム及び記録媒体

Country Status (1)

Country Link
JP (1) JP2022071812A (ja)

Similar Documents

Publication Publication Date Title
EP3078948B1 (en) Acoustic and vibration information accumulation mechanism, acoustic and vibration sensing system, and computer program
JP6073799B2 (ja) 冠状動脈疾患を診断するための周波数電力を検出するシステム
Jiang et al. A cardiac sound characteristic waveform method for in-home heart disorder monitoring with electric stethoscope
JP6131404B2 (ja) 心臓の機能不全と異常を表す情報を測定するための方法及び器械
US20220211323A1 (en) Method and apparatus for non-invasive detection of physiological and patho-physiological sleep conditions
WO2020262563A1 (ja) 健康監視装置、コンピュータプログラム、記録媒体及び生体信号測定装置
JP6097495B2 (ja) 生体状態分析装置及びコンピュータプログラム
EP1968434A2 (en) A method for detecting cardiovascular problems using micro or nano vibrations
Wang et al. Heart sound measurement and analysis system with digital stethoscope
US7517319B2 (en) Method and system for analyzing cardiovascular sounds
Feng et al. Non-contact home health monitoring based on low-cost high-performance accelerometers
WO2022092243A1 (ja) 生体信号分析装置、コンピュータプログラム及び記録媒体
WO2020166260A1 (ja) 体調判定装置、コンピュータプログラム及び記録媒体
JP2022071812A (ja) 生体信号分析装置、コンピュータプログラム及び記録媒体
WO2019139155A1 (ja) 血圧推定装置、血圧推定方法、コンピュータプログラム及び記録媒体
JP7278566B2 (ja) 生体状態推定装置及び生体状態推定システム
Park et al. Ballistocardiography
JP2022071786A (ja) 心尖拍動検出装置、コンピュータプログラム及び記録媒体
Azad et al. Spatial distribution of seismocardiographic signals
Fujita et al. Extraction of apex beat waveform from acoustic pulse wave by sound sensing system using stochastic resonance
JP2024004423A (ja) 血圧推定装置、コンピュータプログラム及び記録媒体
Pouladian Identification of brachial artery characteristic using a cuff for noninvasive detection of coronary artery disease

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20240507