JP6696639B2 - 周波数測定器、周波数測定方法及び周波数測定プログラム - Google Patents

周波数測定器、周波数測定方法及び周波数測定プログラム Download PDF

Info

Publication number
JP6696639B2
JP6696639B2 JP2016035780A JP2016035780A JP6696639B2 JP 6696639 B2 JP6696639 B2 JP 6696639B2 JP 2016035780 A JP2016035780 A JP 2016035780A JP 2016035780 A JP2016035780 A JP 2016035780A JP 6696639 B2 JP6696639 B2 JP 6696639B2
Authority
JP
Japan
Prior art keywords
data
frequency
distance data
mechanical vibrations
independent component
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2016035780A
Other languages
English (en)
Other versions
JP2017148401A (ja
Inventor
義彦 桑原
義彦 桑原
大樹 松本
大樹 松本
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Shizuoka University NUC
Original Assignee
Shizuoka University NUC
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 Shizuoka University NUC filed Critical Shizuoka University NUC
Priority to JP2016035780A priority Critical patent/JP6696639B2/ja
Publication of JP2017148401A publication Critical patent/JP2017148401A/ja
Application granted granted Critical
Publication of JP6696639B2 publication Critical patent/JP6696639B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Measuring Pulse, Heart Rate, Blood Pressure Or Blood Flow (AREA)
  • Measuring And Recording Apparatus For Diagnosis (AREA)

Description

本発明は周波数測定器、周波数測定方法及び周波数測定プログラムに係り、特に振幅の異なる複数の機械的振動のリアルタイム測定に応用することに好適な周波数測定器、周波数測定方法及び周波数測定プログラムに関する。
高齢化社会を迎え、寝たきり老人の健康管理のために、心拍・呼吸状態を安価かつ簡便な装置で監視する必要性が増加している。このため、定在波レーダモジュールで人体の呼吸による胸の上下動を測定して人体の異常検出する人体異常検知装置が提案されている(特許文献1参照。)。この特許文献1の別の実施例には心拍数の測定も可能であるかの記載も見いだせる。
しかし、特許文献1には、寝たきり老人の呼吸と心拍の両方を同時に測定できる旨の明示の記載はない。仮に、特許文献1に記載された発明で、寝たきり老人の呼吸と心拍の両方を同時にリアルタイムで計測するとすれば、呼吸数測定用と心拍数測定用の2台の定在波レーダモジュールが必要になる。よって、特許文献1に記載された発明では、装置が大がかりになるので、寝たきり老人の健康管理のための安価かつ簡便な装置としては適当ではない。
そもそも、特許文献1に心拍数の測定も可能であるとの希望的な記載は見いだせるものの、特許文献1の図25及び図26に呼吸周期を定在波レーダで測定した結果が示されているのみであって、心拍数の測定結果等の現実のデータの開示は特許文献1のどこを探してもない。現実には、胸の上下動に比べて心拍や脈拍の変位は1/100程度であるので、心拍の測定は胸の上下動にマスクされてしまい、振幅の異なる複数の機械的振動のリアルタイム測定は、非常に困難である。
特開2014−217453号公報
本発明は、被測定対象に同時に発生し、振幅の異なる複数の機械的振動であっても、単一距離データとして1台の距離データ測定装置で測定し、単一距離データから複数の機械的振動を分離してそれぞれの周波数をリアルタイムで測定できる周波数測定器、この周波数測定器を用いた周波数測定方法及びこの周波数測定方法をコンピュータに実行させる周波数測定プログラムを提供することを目的とする。
上記目的を達成するために、本発明の第1の態様は、(a) 同一位置で振動している周波数の異なる複数の機械的振動を、単一距離データとして測定する距離データ測定装置と、(b)単一距離データを、複数の機械的振動と同数の複数の疑似信号流データに分離する複数信号流生成手段と、(c)複数の疑似信号流データをそれぞれフーリエ変換するフーリエ変換装置と、(d)フーリエ変換によって得られた各周波数成分について独立成分分析をして、複数の機械的振動のそれぞれの周波数を取得する独立成分分析手段と、を備える周波数測定器であることを要旨とする。
本発明の第2の態様は、(a)距離データ測定装置によって、同一位置で振動している周波数の異なる複数の機械的振動を、単一距離データとして測定するステップと、(b)複数信号流生成手段によって、単一距離データを、複数の機械的振動と同数の複数の疑似信号流データに分離するステップと、(c)フーリエ変換装置によって、複数の疑似信号流データをそれぞれフーリエ変換するステップと、(d)独立成分分析手段によって、フーリエ変換によって得られた各周波数成分について独立成分分析をして、複数の機械的振動のそれぞれの周波数を取得するステップと、を含む周波数測定方法であることを要旨とする。
本発明の第3の態様は、(a)複数信号流生成手段に、距離データ測定装置が測定した同一位置で振動している周波数の異なる複数の機械的振動を単一距離データとして入力させ、この単一距離データを、複数の機械的振動と同数の複数の疑似信号流データに分離させる命令と、(b)フーリエ変換装置によって、複数の疑似信号流データをそれぞれフーリエ変換する命令と、(c)独立成分分析手段によって、フーリエ変換によって得られた各周波数成分について独立成分分析をして、複数の機械的振動のそれぞれの周波数を取得する命令と、を含む一連の命令をコンピュータシステムを構成している信号分離演算処理装置に実行させる周波数測定プログラムであることを要旨とする。
本発明によれば、被測定対象に同時に発生し、振幅の異なる複数の機械的振動であっても、単一距離データとして1台の距離データ測定装置で測定し、単一距離データから複数の機械的振動を分離してそれぞれの周波数をリアルタイムで測定できる周波数測定器、この周波数測定器を用いた周波数測定方法及びこの周波数測定方法をコンピュータに実行させる周波数測定プログラムを提供することができる。
本発明第1の実施の形態に係る周波数測定器の一態様としての構成の概略を説明する模式的なブロック図である。 第1の実施の形態に係る周波数測定器を、寝たきり老人の健康管理として、心拍・呼吸状態を非接触バイタルモニタリング法で管理する状況を例示的に説明する模式図である。 第1の実施の形態に係る周波数測定器の定在波レーダを用い、横隔膜の振動を介して心拍・呼吸状態を非接触でバイタルモニタリングする場合の原理を説明する模式図である。 第1の実施の形態に係る周波数測定器に用いられる独立成分分析手段の論理的なハードウェア資源の概略構造を説明するブロック図である。 図4に示した独立成分分析手段を論理的なハードウェア資源として用いた場合の高速独立成分分析(FastICA)のアルゴリズムの一例の概略を説明するフローチャトである。 第1の実施の形態に係る周波数測定器の定在波レーダによる測定における合成波の電力の距離依存性を説明する模式図である。 第1の実施の形態に係る周波数測定器の定在波レーダによる測定における周波数掃引幅と発信波の中心周波数との関係、並びに合成波の電力の実数部と虚数部等を説明する模式図である。 第1の実施の形態に係る周波数測定器が測定した、呼吸成分と心拍成分が重畳した単一の距離データの変位の時間変化を示す図である。 図8の単一の距離データをフーリエ変換した後のスペクトルを示す図である。 図8の単一の距離データから呼吸と心拍のスペクトルを分離したことを示す図である。 図8の単一の距離データからデータの開始点を3サンプル点分ずらして実線と破線で示した2本の疑似信号流データに分離したことを示す模式図である。 図11に示した2本の疑似信号流データのうちの1本として実線の曲線の方をフーリエ変換した後のスペクトルを示す図である。 図12のフーリエ変換後のスペクトルに対し独立成分分析を行った結果得られたスペクトルに、高調波成分があることを示す図である。 高調波成分をカットする補正をして、図12のフーリエ変換後のスペクトルに対し独立成分分析を行った結果のスペクトルを示す図である。 図15(a)は振動シミュレータで生成した呼吸に対応する変位スペクトル、図15(b)は振動シミュレータで生成した心拍に対応する変位スペクトル、図15(c)は振動シミュレータで生成した呼吸と心拍の合成変位スペクトルを示す図である。 図15(c)の波形をフーリエ変換した後のスペクトルを示す図である。 図16のスペクトルに対し、高調波成分をカットして独立成分分析をした結果を示す図である。 単一距離データを複数の疑似信号流データに分離する際の遅延量の定義を説明するために、観測信号の時間変化を模式的に示した波形カーブである。 図18に示した遅延時間を規定するサンプル点の数が1で、分離に用いる遅延時間が最適値より短い場合の、独立成分分析による周波数測定の結果を示す図である。 遅延時間を規定するサンプル点の数が2で、分離に用いる遅延時間が最適値の範囲に含まれる場合の、独立成分分析による周波数測定の結果を示す図である。 遅延時間に含まれるサンプル点の数が3で、分離に用いる遅延時間が最適値の範囲に含まれる場合の、独立成分分析による周波数測定の結果を示す図である。 遅延時間に含まれるサンプル点の数が4で、分離に用いる遅延時間が最適値の範囲に含まれる場合の、独立成分分析による周波数測定の結果を示す図である。 遅延時間に含まれるサンプル点の数が5で、分離に用いる遅延時間が最適値より長い場合の、独立成分分析による周波数測定の結果を示す図である。 遅延時間に含まれるサンプル点の数が8で、分離に用いる遅延時間が最適値より長い場合の、独立成分分析による周波数測定の結果を示す図である。 遅延時間に含まれるサンプル点の数が10で、分離に用いる遅延時間が最適値より長い場合の、独立成分分析による周波数測定の結果を示す図である。 遅延時間に含まれるサンプル点の数が12で、分離に用いる遅延時間が最適値より長い場合の、独立成分分析による周波数測定の結果を示す図である。
次に、図面を参照して、本発明の第1の実施の形態を説明する。以下の図面の記載において、同一又は類似の部分には同一又は類似の符号を付している。ただし、図面は模式的なものであり、比率等は現実のものとは異なることに留意すべきである。したがって、具体的な技術的内容は以下の説明を参酌して判断すべきものである。又、図面相互間においても互いの寸法の関係や比率が異なる部分が含まれていることは勿論である。
又、以下に示す第1の実施の形態は、本発明の技術的思想を具体化するための装置や方法を例示するものであって、本発明の技術的思想は、構成部品の構造、配置等を下記のものに特定するものでない。
例えば、以下に示す第1の実施の形態では独立した機械的振動が2つの場合について例示的に説明するが、独立成分分析が対象とする複数の機械的振動は、3つ以上の複数であっても構わないことは独立成分分析の原理上、勿論である。よって、本発明の技術的思想は、特許請求の範囲に記載された請求項が規定する技術的範囲内において、種々の変更を加えることができる。
(第1の実施の形態)
本発明の第1の実施の形態に係る周波数測定器は、図1に示すように、同一位置で振動している周波数の異なる複数の機械的振動(S1,S2)を単一距離データとして測定する距離データ測定装置2と、単一距離データを複数の機械的振動(S1,S2)と同数の複数の疑似信号流データに分離する複数信号流生成手段31と、複数の疑似信号流データをそれぞれフーリエ変換するフーリエ変換装置32と、フーリエ変換によって得られた各周波数成分について独立成分分析をして複数の機械的振動(S1,S2)のそれぞれの周波数を取得する独立成分分析手段33を備える。
第1の実施の形態に係る周波数測定器の複数信号流生成手段31は、単一距離データを構成するデータ点の内から時間の異なるデータ点を選び、それぞれを複数の疑似信号流データの開始データとしている。このため、複数信号流生成手段31が分離する複数の疑似信号流データは、互いに一定時間遅延した時間的な関係にある。又、複数信号流生成手段31が分離する複数の疑似信号流データの流れの本数は、複数の機械的振動(S1,S2)と同数である。
距離データ測定装置2は、単一距離データを電磁波の信号として単一の受信アンテナ29で受信する定在波レーダである。図3では、複数の機械的振動(S1,S2)が被測定対象1である人間の横隔膜が発生する振動に重畳されているモデルを例示している。図3は複数の機械的振動(S1,S2)が、被測定対象1である人間の横隔膜を介して検出される呼吸成分の振動S1と、心臓が発生する心拍成分の振動S2とを含むモデルであり、第1の実施の形態に係る周波数測定器は、呼吸成分の振動S1と心拍成分の振動S2とを独立成分分析手段33によって分離する心拍・呼吸モニタ装置として機能している。
第1の実施の形態に係る周波数測定器は、図2に示したような、被測定対象1としての寝たきり老人の健康管理のための心拍・呼吸状態の監視をする心拍・呼吸モニタ装置として応用可能である。図2の模式図ではベッドの上に寝た被測定対象1としての寝たきり老人の腹部に定在波レーダ(距離データ測定装置)2を用いた非接触バイタルモニタリングを実現する心拍・呼吸モニタ装置が提供できる。
第1の実施の形態に係る周波数測定器によれば、被測定対象1からの反射波VR(S1,S2)とレーダからの入射波VTにより生じる定在波を観測することで、呼吸成分の振動S1と心拍成分の振動S2の検出を行う。呼吸成分の振動S1に比べ、体表変位が小さいため心拍成分の振動S2を検出するのは困難であったが、図1に示す独立成分分析手段33によって独立成分分析を使用することにより、呼吸成分と心拍成分の混ざった単一距離データの信号から呼吸成分の振動S1と心拍成分の振動S2を分離することが可能になる。
図1に示す独立成分分析手段33は、複数の信号源が統計的に独立であることを前提に、単一の観測データである距離データのみを用いて、呼吸成分の振動S1と心拍成分の振動S2との元の信号源を推定する独立成分分析法を実行する。即ち、第1の実施の形態に係る周波数測定器では、便宜上、図2及び図3に示したような、元の信号源が2つの最も簡単な場合について例示的に説明するが、元の信号源が3つ以上、N個であっても構わないことは勿論である。
第1の実施の形態に係る周波数測定器では信号源の数N=2の場合であるが、統計的に独立なN個(≧2)の信号源が、N個観測される場合となる一般的な混合モデルでは、N個の観測信号(混合行列)の行列x(t) = [x1(t), x2(t), …, xN(t)]Tは、混合行列Aを用いると、独立なN個の信号源の行列s(t) = [s1(t), s2(t), …, sN(t)]Tに対して(1)式で表される関係となる:

x(t) = As(t) ……(1)

ここで、混合行列Aは、信号源s(t)が観測点に到達するまでの伝達特性をもつ未知の定数の行列である。
混合モデルに対する分離モデルでは、分離信号の行列u(t) = [u1(t), u2(t), …, uN(t)]T、分離行列Wを用いて(2)式で表される:

u(t) = Wx(t) ……(2)

混合信号 x(t)だけをデータとして使用し、分離信号が統計的に独立となるように分離行列Wを逐次更新しながら、分離信号 u(t)を生成するのが独立成分分析法(ICA法)のアルゴリズムである。分離モデルの式を周波数領域に変換すると(3)式を得る:

u(ω,ts )=W(ω)x(ω,ts) ……(3)
第1の実施の形態に係る周波数測定器では、フーリエ変換装置32がフーリエ変換して得られるスペクトログラムx(ω,ts)の各周波数帯に、独立成分分析手段33が周波数領域における独立成分分析法(周波数領域ICA)を適用して分離行列を求める。第1の実施の形態に係る周波数測定器では信号源の数N=2の場合について、独立成分分析手段33が高速独立成分分析(FastICA)のアルゴリズムを用いて分離行列を求める。
図1に示すように、定在波レーダ(距離データ出力装置)2は、被測定対象1に対しマイクロ波、ミリ波又はサブミリ波の電磁波を照射し、被測定対象1から反射された電磁波を検出するRF送受信装置25と、このRF送受信装置25を制御し、RF送受信装置25が検出した電磁波の信号を処理し、距離データとして出力するワンチップ・プロセッサ(集積回路)21を備えている。
このため、RF送受信装置25は、マイクロ波等の電磁波を発振する発振器251と、発振器251が発振した電磁波を被測定対象1に対し照射する送信アンテナ28と、被測定対象1から反射された電磁波を受信する受信アンテナ29と、受信アンテナ29からの電気信号を位相の異なる電気信号に分離するミキサ252と、ミキサ252で分離された位相の異なる電気信号をそれぞれ検波する第1検波ダイオード253及び第2検波ダイオード254とを備える。発振器251は発信波の中心周波数 f0 =3GHz〜3THzの電磁波を発信するが、距離データの分解能(変位分解能)を考慮すると、発信波の中心周波数 f0 =20GHz以上の高周波であることが好ましい。
なお、図1に示した定在波レーダ(距離データ出力装置)2には、送信アンテナ28と受信アンテナ29とが、それぞれ別の空中線である場合が模式的に表示されているが、例示に過ぎない。送信アンテナ28と受信アンテナ29とを同一の空中線として構成しても構わないことは定在波レーダの技術としては明らかな事項であり、送信アンテナ28と受信アンテナ29とを同一アンテナとした方が、第1の実施の形態に係る周波数測定器の構造がより簡略化できるので、工業的な意味では好ましい。又、送信アンテナ28と受信アンテナ29は素子長がλ/2のパッチアンテナを30素子程度配列して構成しても、他の構成でも構わない。
図3において、送信波VT の振幅と周波数をA、f、また光速をcとすると送信波VT は(4)式で表せる:

Figure 0006696639
ただし、図3において送信波VT (f、x)の照射方向に定義されるX軸上の任意の一点をx=0として図6のx1とx2を定義している。
また、周波数掃引幅fwを用いると、発振器251の発信波の中心周波数 f0に対し、図7に示すように、周波数をf = f0 - fw/2からf = f0 + fw/2までステップ状に変化させる。n番目の被測定対象1の距離をdn、x軸上の任意の点における送信波VT に対する反射波VRの大きさの比をγn、位相差をφnとすれば、被測定対象1からの反射波VR は(5)式のようになる:

Figure 0006696639
第1検波ダイオード253及び第2検波ダイオード254からの信号は、図1に示すように、それぞれ第1検波増幅器27及び第2検波増幅器26を介して、ワンチップ・プロセッサ(集積回路)21のA/D変換装置24に入力される。
ワンチップ・プロセッサ(集積回路)21は、A/D変換装置24が出力したディジタル信号を処理する制御・演算処理回路22を備える。制御・演算処理回路22から距離データの信号が信号分離演算処理装置3に出力されると同時に、D/A変換装置23で変換されたアナログ信号が変調信号増幅器25を介して、RF送受信装置25の発振器251に出力される。
図1に示すように、制御・演算処理回路22は、D/A変換装置23にFM変調された周波数制御電圧を変調信号として出力する変調信号生成手段221と、A/D変換装置24が出力したディジタル信号を処理する信号処理手段222と、信号処理手段222の出力を用いて距離データを演算する距離データ演算手段223と、距離データ演算手段223の演算結果を図8に示したような時間変化のデータとして信号分離演算処理装置3に出力する距離データ出力手段224を備える。
制御・演算処理回路22の信号処理手段222は、図3及び図6に示したx=x1での位置における合成波の電力p(f、x1 )を実数部、π/2シフトしたx=x2での位置における合成波の電力p(f、x2 )を虚数部とおき、γn≪1, f = f0 - fd, (-fw/2 < fd < +fw/2)として、(6)式を用いて、解析信号pa(fd)を周波数のシフト量fdの関数として求める:

Figure 0006696639
解析信号pa(fd )は周期関数であり、解析信号pa(fd )の周期は被測定対象1との距離に逆比例の関係がある。このため、制御・演算処理回路22の距離データ演算手段223が、解析信号pa(fd )をフーリエ変換することにより、距離スペクトルP(x)が得られる。そして、この距離スペクトルP(x)のピークの位置から、被測定対象1までの距離が(7)式を用いて求めることができる:

Figure 0006696639
信号分離演算処理装置3は、通常のコンピュータシステムの中央演算処理装置(CPU)と同様な機能を有し、通常のコンピュータシステムと同様に、信号分離演算処理装置3の演算結果や演算の途中のデータを記憶するデータ記憶装置36が信号分離演算処理装置3に接続されている。信号分離演算処理装置3には通常のコンピュータシステムと同様に、入力装置37及び出力装置35が接続されている。
図1に示すように信号分離演算処理装置3は、1本の距離データから独立成分分析が対象とする複数本の疑似信号流データを生成する複数信号流生成手段31と、複数信号流生成手段31が生成した複数の疑似信号流データのそれぞれをフーリエ変換するフーリエ変換装置32と、フーリエ変換されたそれぞれの周波数成分に対して、独立成分分析をする独立成分分析手段33と、独立成分分析手段33が分析し分離した複数の機械的振動を出力装置35に出力させる独立成分表示手段34とを備えるマイクロプロセッサである。ここで、複数信号流生成手段31、フーリエ変換装置32、独立成分分析手段33、独立成分表示手段34は、それぞれ信号分離演算処理装置3を構成する論理的なハードウェア資源であり、必ずしも専用の回路が現実に存在することを前提とするものではない。
独立成分分析手段33は、(3)式で定義される「周波数領域分離モデル」の式に対し、高速独立成分分析のアルゴリズムを実行する。このため、独立成分分析手段33は、図4に示すように、前処理としての白色化処理を実行する白色化手段331と、非線形関数を導出する非線形関数導出手段332と、解空間をユニタリ行列に限定する解空間限定手段333と、ユニタリ行列を初期値に設定する初期値設定手段334と、分離信号を計算する分離信号計算手段335と、ニュートン法による更新と降下を実行する更新手段336と、グラム・シュミットの直交化を実行する直交化手段337と、ノルムの正規化を実行する正規化手段338と、各分離信号の成分iについて収束したか否かを判定する最小化判定手段339を、それぞれ論理的なハードウェア資源として備える。
信号分離演算処理装置3には、マイクロチップとして実装されたマイクロプロセッサ(MPU)等を使用してコンピュータシステムを構成することが可能である。又、コンピュータシステムを構成する信号分離演算処理装置3として、算術演算機能を強化し信号処理に特化したデジタルシグナルプロセッサ(DSP)や、メモリや周辺回路を搭載し組込み機器制御を目的としたマイクロコントローラ(マイコン)等を用いてもよい。或いは、現在の汎用コンピュータのメインCPUを信号分離演算処理装置3に用いてもよい。更に、信号分離演算処理装置3の一部の構成又はすべての構成をフィールド・プログラマブル・ゲート・アレイ(FPGA)のようなプログラマブル・ロジック・デバイス(PLD)で構成してもよい。
PLDによって、信号分離演算処理装置3の一部又はすべてを構成した場合は、データ記憶装置36は、PLDを構成する論理ブロックの一部に含まれるメモリブロック等のメモリ要素として構成することができる。更に、信号分離演算処理装置3は、CPUコア風のアレイとPLD風のプログラム可能なコアを同じチップに搭載した構造でもよい。このCPUコア風のアレイは、あらかじめPLD内部に搭載されたハードマクロCPUと、PLDの論理ブロックを用いて構成したソフトマクロCPUを含む。つまりPLDの内部においてソフトウェア処理とハードウェア処理を混在させた構成でもよい。
又、各ハードウェア資源を構成するマイクロプロセッサがコア化され、そのコアの複数が一つのマイクロチップに実装されてプロセッサ・パッケージ内に収められたマルチコア型CPUでも構わない。信号分離演算処理装置3は第1の実施の形態に係る周波数測定器が構成するコンピュータシステム内での演算を行なう中心的機能をなす。
第1の実施の形態に係る周波数測定器の信号分離演算処理装置3は、図示を省略したプログラム記憶装置に格納されたプログラムに規定された命令列を順に読み込んで解釈・実行することで独立成分分析等に必要な情報の加工を行う。図1では詳細な図示を省略しているが、通常のコンピュータシステムと同様に、信号分離演算処理装置3はバスと呼ばれる信号線を介してデータ記憶装置36等の主記憶装置や入出力回路に接続され、何段階かの入出力回路を介して入力装置37及び出力装置35と情報のやりとりを行う。図1ではデータ記憶装置36を模式的に表示しているが、実際には、図示を省略した補助記憶装置や表示装置、通信装置などの周辺機器に接続されてデータやプログラムなど情報のやりとりを行うことができる。
図1の信号分離演算処理装置3の独立成分分析手段33における高速独立成分分析(FastICA)のアルゴリズムは、例えば、図5に示したフローチャートに示したような手順で実行できる。即ち、先ず、図4に示した独立成分分析手段33の白色化手段331が、図5のフローチャートのステップS11において、独立成分分析の前処理として、主成分分析等を用いて白色化処理を実行する。即ち、相関行列が単位行列となるように変換して、観測信号を無相関とする。
そして、ステップS12において、独立成分分析手段33の非線形関数導出手段332が、分離信号yiの非線形関数G(yi)を導出する。非線形関数G(yi)は、対象とする事象の情報量となる。ステップS13において、独立成分分析手段33の解空間限定手段333が分離行列の解空間をユニタリ行列Uに限定する。解空間をユニタリ行列Uに限定することにより高速な独立成分分析が可能になる。更に、ステップS14において、独立成分分析手段33の初期値設定手段334がユニタリ行列Uを初期値に設定する。
その後、ステップS15において、独立成分分析手段33の分離信号計算手段335が、分離信号を計算する。そして、ステップS16において、独立成分分析手段33の更新手段336がニュートン法による更新と降下を実行する。即ち、非線形関数G(yi)の1階微分及び2階微分を求めて、行列の空間での自然な最急降下方向を決める。そして、最急降下方向に勾配法を適用して更新と降下を実行し、非線形関数G(yi)の最適化を実行する。
更に、ステップS17において、独立成分分析手段33の直交化手段337がグラム・シュミットの直交化を実行する。ステップS18において、独立成分分析手段33の正規化手段338がノルムを1に正規化する。ステップS19において、独立成分分析手段33の最小化判定手段339が各iについて、分離信号yiの非線形関数G(yi)の期待値が最小化されたかを判定する。即ち、個々の分離信号yiのエントロピーが最小化され、個々の分離信号yiが収束したか否かを判定する。
ステップS19において分離信号yiが収束していないと判定された場合は、ステップS15に戻り、分離信号yiの非線形関数G(yi)の期待値が最小化されるまで、ステップS15〜ステップS19までの処理を繰り返す。ステップS19において分離信号yiが収束したと判定された場合は、高速独立成分分析のアルゴリズムが終了となる。
図1のワンチップ・プロセッサ(集積回路)21の距離データ出力手段224は、離散データとして30〜100msのデータ間隔で、距離データ演算手段223が演算した距離データを、信号分離演算処理装置3の複数信号流生成手段31に対して図8に示すような周期的な波形の信号流を形成するように、逐次出力する。距離データ演算手段223が演算した結果である距離データは、一旦データ記憶装置36に格納された後、複数信号流生成手段31がデータ記憶装置36から距離データを読み出すようにしても良い。
図9は、図8の単一の距離データをフーリエ変換した後のスペクトルであるが、図9のスペクトルからは呼吸成分の振動S1は読めても、振動振幅の小さい心拍成分の振動S2を読むことができない。第1の実施の形態に係る周波数測定器においては、信号分離演算処理装置3の複数信号流生成手段31は、図8に示したような単一の距離データを、図11に実線と破線で示したような、変位の時間変化を示す2本の疑似信号流データの流れに分離する。このとき、複数信号流生成手段31は、図11に示した単一距離データを構成するデータ点の内から時間の異なる2点のデータ点を選び、それぞれを実線と破線で示した2本の疑似信号流データの開始データとしている。図11では実線と破線を区別して見やすくするため、破線のカーブが実線のカーブより3サンプル点分進んでいるかのように模式的に示してはいるが、実線と破線は本来、単一の距離データを示す同一の1本の曲線である。
図11に破線で示した疑似信号流データのカーブと実線で示した疑似信号流データのカーブとを分離する遅延量(遅延時間)には最適値がある。 図18に示す観測信号の時間変化を示す波形カーブは、図11の破線のカーブと実線のカーブに対応する2本の疑似信号流データを分離する際に用いられる遅延量(遅延時間)の定義を示している。
シミュレータを用いて観測信号である単一距離データを2つの疑似信号流データに分離する際の遅延量(遅延時間)の長さを変えて、それぞれ独立成分分析を行い、図18で定義した「分離に必要な遅延量(遅延時間)」の最適値を決定した。シミュレーションでは、観測信号(単一の距離データ)を141サンプル点分の長さで測定し、図18に示すx1=128サンプル点長、x2=128サンプル点長の同一の信号長の信号流を2つ取り出すこととした。そして、図19〜図26に示すように、遅延時間内に含まれるサンプル点の数を1,2,3,4,5,8,10,12と順に変えて、様々な遅延時間の長さで独立成分分析を行い、最適値の範囲に属するサンプル点の個数を求めた。
図19〜図26に示す独立成分分析の実験では、呼吸成分の振動S1の周波数を0.2Hz±10%の範囲で,心拍成分の振動S2の周波数を1.3Hz±3%の範囲で一周期ごとに変化させランダムで変化を行うようにした。また、後述する「補正有り」の条件での測定を行った。
図19から図26に示す観測信号の測定は連続で行い、図19〜図26では独立成分分析の結果分離した2つの周波数スペクトルのパワーが最も高い周波数を示している。 図19〜図26に破線で示す真値は揺らぎの中心周波数を示し、実線で示す測定値と比較している。
分離に用いられる遅延時間に含まれるサンプル点の数が、図19では1の場合、図20では2の場合、図21では3の場合、図22では4の場合、図23では5の場合、図24では8の場合、図25では10の場合、及び図26では12の場合である。図19〜図26に示すように、遅延時間を規定するサンプル点の数を変えて、単一距離データ(観測信号)をそれぞれ2つの疑似信号流データに分離して独立成分分析を行った結果、分離に必要な遅延量(遅延時間)の最適値は、図20〜図22に示した2〜4サンプル点分の長さであることが理解できる。
一方、図19が示すとおり、2つの疑似信号流データに分離する際の遅延時間が図20〜図22に示した最適値よりも短い場合の独立成分分析では、心拍成分の振動S2の周波数検出の精度が下がっていることが分かる。又、図23〜図26の心拍成分の振動S2の周波数のバラツキから分かるように、遅延時間が図20〜図22に示した最適値よりも長い場合の独立成分分析でも、周波数検出の精度が下がっていることが分かる。
図8及び図11に示す距離データは、発振器251が中心周波数 f0 =24GHzの信号を発振し、周波数掃引幅fw = 56.25MHzとして定在波レーダ(距離データ出力装置)2が取得したデータである。図8及び図11は、振幅が2mmで0.33Hzの呼吸成分と、振幅が0.02mmで1.3Hzの心拍成分が重畳した変位を示す距離データである。
定在波レーダ(距離データ出力装置)2の送信アンテナ28からは、例えば半値幅が水平方向12°、垂直方向25°となるビーム幅の電磁波が出射される。送信アンテナ28から出射される電磁波のサイドローブは、例えば水平方向15dB、垂直方向15dBに抑制されていれば良い。図2に示したような測定環境であれば、送信アンテナ28から出射される電磁波の出力は、例えば、4mW程度あれば十分である。電磁波の出力が4mW程度あると、100m程度の測定レンジが実現できる。
例えば、距離データ出力手段224がデータ間隔70msで、単一の距離データとなる信号流を複数信号流生成手段31に時々刻々送信すると、複数信号流生成手段31は3サンプル点分の長さ、即ち図11において実線で示した疑似信号流データに対して70×3=210msだけ遅延した破線の疑似信号流データをなす一連の離散データを生成する。この結果、複数信号流生成手段31は、N=2個の信号源の行列s(t) = [s1(t), s2(t)]Tに対応させてN=2個の観測信号の行列x(t) = [x1(t), x2(t)]Tを、単一の距離データを構成している信号流から生成する。
図11に示した計測変位の波形は信号分離演算処理装置3における独立成分分析の対象となるN=2個の信号源の信号であるが、独立成分分析手段33に入力される前に、それぞれ信号分離演算処理装置3のフーリエ変換装置32によって、フーリエ変換される。図11において実線で示した疑似信号流データをフーリエ変換装置32でフーリエ変換したスペクトルを図12に示す。図11において破線で示した疑似信号流データをフーリエ変換装置32でフーリエ変換した結果も、単一の距離データを構成している信号流から生成されているので、図12に示したスペクトルと同一である。
図12のフーリエ変換後のスペクトルを元に、独立成分分析手段33が独立成分分析を行った結果が図13のスペクトルである。実験前に心拍計により心拍を求めた結果1.19Hzであったが、図13のスペクトルに示すように、1.19Hzの付近には、心拍以外に呼吸の第2次高調波のスペクトルも現れていることが分かる。
独立成分分析手段33が分離した呼吸の周波数スペクトルにおいて、もっとも大きい電力を持つスペクトルの周波数が呼吸の基本周波数で,これから0.8Hz 以上まで現れるスペクトルは高調波成分と考えられることから、この間の周波数成分をカットして補正を行う必要がある。
独立成分分析手段33が高調波成分をカットする補正をして独立成分分析を行った結果を図14に示す。0.33Hzに呼吸のスペクトルが、1.17Hzに心拍のスペクトルがはっきり現れ、統計的に独立なN個(≧2)の信号源である、呼吸成分の振動S1と心拍成分の振動S2の周波数を検出することができた。
既に述べたとおり、 図8及び図11に示した距離データには呼吸成分の振動S1と心拍成分の振動S2が重畳されている。 図8及び図11では、呼吸成分の振動S1は振幅が2mmであるのに対し、心拍成分の振動S2は振幅が0.02mmであり、心拍成分の振動S2は呼吸成分の振動S1の振幅の1/100である。このため、従来技術では呼吸成分の振動S1に心拍成分の振動S2がマスクされてしまい、心拍成分の振動S2をリアルタイムで測定するのは、非常に困難であった。
しかしながら、第1の実施の形態に係る周波数測定器によれば、被測定対象1の人間に同時に発生し、振幅が100倍異なる2つの機械的振動である呼吸成分の振動S1と心拍成分の振動S2を1台の距離データ測定装置2で測定しても、図10及び図14に示すように、距離データ測定装置2が出力する単一距離データから呼吸成分の振動S1と心拍成分の振動S2を分離して、呼吸成分の振動S1と心拍成分の振動S2のそれぞれの周波数をリアルタイムで測定できる。
呼吸と心拍の振動を模擬する振動シミュレータの波形を図15(a)〜(c)に示す。振動シミュレータは、厚さ2mm、8cm×12cmのアルミ板が、振幅が100倍異なる心拍と呼吸の2つの振動を同時に模擬するような機構を有している。図15(a)は振動シミュレータで生成した呼吸に対応する0.22Hzの変位スペクトルで、図15(b)は振動シミュレータで生成した心拍に対応する1.2Hzの変位スペクトルである。図15(c)は振動シミュレータで生成した呼吸と心拍の合成変位スペクトルである。
図16は、図15(c)の計測波形をフーリエ変換装置32でフーリエ変換したスペクトルを示す。被測定対象1が人間である場合の測定と同様に、被測定対象1をアルミ板として0.8Hz 以上まで現れるスペクトルをカットした信号を使って独立成分分析手段33で独立成分分析をした結果を図17に示す。図17に示すように、振動シミュレータを用いて、被測定対象1をアルミ板として機械的に発生した信号を合成して観測信号を生成した場合でも、振幅が100倍異なる2つの独立した信号源の周波数を、他方の振動にマスクされずにリアルタイムで検出することができたことが分かる。
(その他の実施の形態)
上記のように、本発明は1つの実施の形態によって記載したが、この開示の一部をなす論述及び図面は本発明を限定するものであると理解すべきではない。この開示から当業者には様々な代替実施の形態、実施例及び運用技術が明らかとなろう。
例えば、既に述べた第1の実施の形態の説明においては、被測定対象1を人間とし、人間の横隔膜を介して検出される呼吸成分の振動S1と心拍成分の振動S2とを独立成分分析の対象としたモデルについて説明したが、被測定対象1は人間に限定されるものではない。又、独立成分分析で分離される振動は、呼吸成分の振動S1と心拍成分の振動S2に限定されるものでもなく、独立成分分析が対象とする独立成分は2に限定されるものでもない。
本発明によれば、人間以外の種々の被測定対象1に対しても、同様に、振幅が100倍程度等大きく異なっている測定環境や測定条件において、複数の機械的振動を1台の距離データ測定装置2で測定し、距離データ測定装置2が出力する単一距離データから複数の機械的振動を分離して、複数の機械的振動のそれぞれの周波数をリアルタイムで測定できる。
更に、上記の第1の実施の形態の説明においては振幅が100倍も異なっている場合について例示的に説明したが、発明の技術的思想を説明する上の便宜に過ぎない。本発明によれば、振幅が大きく異なっていない測定環境や測定条件であっても、同様に、種々の被測定対象1の複数の機械的振動を1台の距離データ測定装置2で測定し、複数の機械的振動のそれぞれの周波数をリアルタイムで測定できることは、独立成分分析の原理から容易に理解できるであろう。
このように、本発明はここでは記載していない様々な実施の形態等を含むことは勿論である。したがって、本発明の技術的範囲は上記の説明から妥当な特許請求の範囲に係る発明特定事項によってのみ定められるものである。
本発明は高齢化社会における寝たきり老人の心拍・呼吸状態を非接触バイタルモニタリング法によって監視する医療機器の技術分野や、複雑な振動のリアルタイム測定の技術分野に利用可能である。
1…被測定対象
2…距離データ測定装置
21…ワンチップ・プロセッサ(集積回路)
22…演算処理回路
221…変調信号生成手段
222…信号処理手段
223…距離データ演算手段
224…距離データ出力手段
251…発振器
252…ミキサ
253…検波ダイオード
254…検波ダイオード
23…A変換装置
24…D変換装置
25…変調信号増幅器
25…RF送受信装置
26…検波増幅器
27…検波増幅器
28…送信アンテナ
29…受信アンテナ
3…信号分離演算処理装置
31…複数信号流生成手段
32…フーリエ変換装置
33…独立成分分析手段
34…独立成分表示手段
35…出力装置
36…データ記憶装置
37…入力装置
331…白色化手段
332…非線形関数導出手段
333…解空間限定手段
334…初期値設定手段
335…分離信号計算手段
336…更新手段
337…直交化手段
338…正規化手段
339…最小化判定手段

Claims (7)

  1. 同一位置で振動している周波数の異なる複数の機械的振動を、変位の時間変化を示す単一距離データとして測定する距離データ測定装置と、
    前記単一距離データを構成するデータ点内から時間の異なるデータ点を選び、それぞれを複数の疑似信号流データの開始データとして分離して、前記単一距離データを、前記複数の機械的振動と同数で互いに一定時間遅延した時間的な関係にある、前記複数の疑似信号流データに分離する複数信号流生成手段と、
    前記複数の疑似信号流データをそれぞれフーリエ変換するフーリエ変換装置と、
    前記フーリエ変換によって得られた各周波数成分について独立成分分析をして、前記複数の機械的振動のそれぞれの周波数を取得する独立成分分析手段と、
    を備えることを特徴とする周波数測定器。
  2. 前記単一距離データを、電磁波の信号として単一のアンテナで受信することを特徴とする請求項に記載の周波数測定器。
  3. 前記距離データ測定装置が、定在波レーダであることを特徴とする請求項に記載の周波数測定器。
  4. 前記複数の機械的振動が、人間の横隔膜が発生する振動に重畳されていることを特徴とする請求項1〜のいずれか1項に記載の周波数測定器。
  5. 前記複数の機械的振動が、人間の横隔膜を介して検出される呼吸の振動と、心臓が発生する振動とを含むことを特徴とする請求項に記載の周波数測定器。
  6. 距離データ測定装置によって、同一位置で振動している周波数の異なる複数の機械的振動を、変位の時間変化を示す単一距離データとして測定するステップと、
    前記単一距離データを構成するデータ点内から時間の異なるデータ点を選び、それぞれを複数の疑似信号流データの開始データとして分離して、複数信号流生成手段によって、前記単一距離データを、前記複数の機械的振動と同数で互いに一定時間遅延した時間的な関係にある、前記複数の疑似信号流データに分離するステップと、
    フーリエ変換装置によって、前記複数の疑似信号流データをそれぞれフーリエ変換するステップと、
    独立成分分析手段によって、前記フーリエ変換によって得られた各周波数成分について独立成分分析をして、前記複数の機械的振動のそれぞれの周波数を取得するステップと、
    を含むことを特徴とする周波数測定方法。
  7. 複数信号流生成手段に、距離データ測定装置が測定した同一位置で振動している周波数の異なる複数の機械的振動を変位の時間変化を示す単一距離データとして入力させ、前記単一距離データを構成するデータ点内から時間の異なるデータ点を選び、それぞれを複数の疑似信号流データの開始データとして分離して、前記単一距離データを、前記複数の機械的振動と同数で互いに一定時間遅延した時間的な関係にある、前記複数の疑似信号流データに分離させる命令と、
    フーリエ変換装置によって、前記複数の疑似信号流データをそれぞれフーリエ変換する命令と、
    独立成分分析手段によって、前記フーリエ変換によって得られた各周波数成分について独立成分分析をして、前記複数の機械的振動のそれぞれの周波数を取得する命令と、
    を含む一連の命令をコンピュータシステムを構成している信号分離演算処理装置に実行させることを特徴とする周波数測定プログラム。
JP2016035780A 2016-02-26 2016-02-26 周波数測定器、周波数測定方法及び周波数測定プログラム Active JP6696639B2 (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP2016035780A JP6696639B2 (ja) 2016-02-26 2016-02-26 周波数測定器、周波数測定方法及び周波数測定プログラム

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2016035780A JP6696639B2 (ja) 2016-02-26 2016-02-26 周波数測定器、周波数測定方法及び周波数測定プログラム

Publications (2)

Publication Number Publication Date
JP2017148401A JP2017148401A (ja) 2017-08-31
JP6696639B2 true JP6696639B2 (ja) 2020-05-20

Family

ID=59739397

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2016035780A Active JP6696639B2 (ja) 2016-02-26 2016-02-26 周波数測定器、周波数測定方法及び周波数測定プログラム

Country Status (1)

Country Link
JP (1) JP6696639B2 (ja)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2019120539A (ja) * 2017-12-28 2019-07-22 有限会社起福 トイレブース使用状況報知システム

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2002102187A (ja) * 2000-09-29 2002-04-09 Matsushita Electric Ind Co Ltd 生体検出装置
JP3873793B2 (ja) * 2002-03-29 2007-01-24 日本電気株式会社 顔メタデータ生成方法および顔メタデータ生成装置
JP4405343B2 (ja) * 2004-08-23 2010-01-27 株式会社デンソー 心拍計測装置
WO2014002276A1 (ja) * 2012-06-29 2014-01-03 富士通株式会社 バイタルサイン検出方法、バイタルサイン検出装置及びバイタルサイン検出プログラム
JP5383890B1 (ja) * 2012-10-22 2014-01-08 沖電気工業株式会社 物体数検知装置及びプログラム
JP6123450B2 (ja) * 2013-04-16 2017-05-10 富士通株式会社 生体情報取得装置、方法及びプログラム
JP2014217453A (ja) * 2013-05-02 2014-11-20 斎藤 光正 定在波レーダーによる人体異常検知装置及びその利用方法

Also Published As

Publication number Publication date
JP2017148401A (ja) 2017-08-31

Similar Documents

Publication Publication Date Title
US11464419B2 (en) Methods for training a model for use in radio wave based health monitoring
US10746863B2 (en) Target extraction system, target extraction method, information processing apparatus, and control method and control program of information processing apparatus
JP2014033914A (ja) 超音波画像診断装置及び超音波画像診断方法
CN104545859B (zh) 生理量测的感测系统与方法
JP6696639B2 (ja) 周波数測定器、周波数測定方法及び周波数測定プログラム
JP2015109960A (ja) 超音波診断装置、超音波診断装置の制御器及び超音波診断装置の制御方法
JP6484290B2 (ja) 近傍界測定装置及び近傍界測定方法
JP2013160508A (ja) レーダ装置評価システムおよびレーダ装置評価方法
JP2014039702A5 (ja)
CN107595260A (zh) 非接触式体征检测方法、装置、存储介质及其计算机设备
US10709423B2 (en) Ultrasound processing apparatus and method
Gunasekara Contactless estimation of breathing rate using UWB radar
EP2944265A1 (en) Ultrasonic diagnostic device
JP2015045631A (ja) 周波数変動信号解析装置
JP6747600B2 (ja) 心拍測定装置
JP2014153206A (ja) 信号処理装置及び信号処理方法
US20170258437A1 (en) Ultrasonic diagnostic device
JP5260897B2 (ja) 超音波診断装置
JP7327868B1 (ja) 生体状態測定装置、生体状態測定方法、プログラム及び生体状態測定システム
Burov et al. Acoustic tomography of the nonlinear parameter by a small number of transducers
WO2023037613A1 (ja) 変位検知装置及び方法
WO2024038629A1 (ja) 生体状態測定装置、生体状態測定方法、プログラム及び生体状態測定システム
JP2010216940A (ja) 信号処理方法及び信号処理装置
JP6161028B2 (ja) 映像装置及び映像形成方法
JP2018108142A (ja) 超音波診断装置及び制御方法

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20190109

A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20191031

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20191112

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20191130

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: 20200324

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20200409

R150 Certificate of patent or registration of utility model

Ref document number: 6696639

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150

R250 Receipt of annual fees

Free format text: JAPANESE INTERMEDIATE CODE: R250

R250 Receipt of annual fees

Free format text: JAPANESE INTERMEDIATE CODE: R250