計測自動制御学会東北支部 第 307 回研究集会(2017.2.27)
資料番号 307-4

ガイド波による骨密度測定で用いる
線形多孔質モデルが適用可能な空隙率の検討
Investigation of porosity suitable for the linear porous model used
for bone density estimation with guided waves○小野聡 齊藤玄敏
Satoshi Ono, Hiroyuki Saitou
弘前大学大学院
Hirosaki Univ
キーワード：

ガイド波(Guided wave)， アレイ信号処理(Array signal processing)，速度分散(Dispersion)，多孔質体(Porous solid)，センブランス解析(Semblance analysis)，最尤推定法(Maximum likelihood estimation)

連絡先： 〒036-8561

青森県弘前市文京町 3 番地 弘前大学大学院 理工学研究科

齊藤玄敏，E-mail: saitoh@hirosaki-u.ac.jp

＿＿＿＿＿
1． 緒言

成された導波路に沿って，長手方向に減衰するこ

近年，平均寿命が延びたことや生活習慣の変化となく伝播する特性を持っている．また境界からにより，骨粗鬆症患者は増加傾向にある．この予1 波長程度の深さの力学的性質を反映することか防にはまめに検査を受信することで早期に症状のら，皮膚近くに存在する骨であれば，その特性を進行度を把握することが必要となる．

評価することが出来ると思われる．

超音波骨密度測定の利点は，計測装置が軽量・

ただし，ガイド波には速度分散性があるため，安価，操作に専門的技術を必要としない，大人数速度を正確に測定するためにはアレイ計測を行うに対する頻繁な検査に適している，などがある．

必要がある．このため，測定部にはセンサを展開

一般的な超音波骨密度検査法は踵骨にパルスをするのに十分な大きさが求められる．本研究では透過させる方式で，骨内伝搬速度を計測し，それ測定部として下腿脛骨部を想定している．皮膚やを指標に骨密度をデータベースから割り出す．問

脂肪などの生体組織が他の場所より薄いためで測

題点としては，計測波形が 1 波形のみで，波形の定部位として適していると考える．

立ち上がりを正確に同定することが難しいこと，

斉藤ら 2)は空隙率の異なる同一素材の多孔質材

また，実体波は高密度領域を選択的に伝播するた料を用いて，空隙率の違いによるガイド波の速度め，骨密度の低下を正しくとらえることが出来な分散性を調べ，その特性から多孔質丸棒の実体波

いこと，などが挙げられる．

速度を推定し，試料の密度推定を行った．この研

この問題を解決するために，我々の研究室では究から空隙率が 50%を超えると密度推定に大きなガイド波を用いた新たな骨密度検査法の開発に取誤差が生じることがわかった．この原因として，り組んでいる 1)2)．ガイド波は物理境界によって形線形多孔質モデルの適用限界とせん断波速度の推1定誤差が考えられる．

X N−1 (ω, τ, p)は解析窓で切り出した波形のフーリ
エスペクトル，𝑒 𝑗𝜔𝑝𝑧𝑁−1 は位相補正項，＋はエルミ

本研究では，多孔質材料の圧縮波とせん断波の計測と速度推定を行い，すべて実測した圧縮波，ート行列を表す．

せん断波速度，密度から試料の弾性係数を算出す

次に最尤推定法について説明をする．計測デーる．また，線形多孔質モデルにおける縦弾性係数

タに対して，あるパラメータ値が得られる確率をと剛性率間の比例定数を算出し，線形多孔質モデ尤度とみなし，尤度が最大になるようなパラメー

ルの適用限界について確認する．

タを探索する方法である 3)．R(ω, τ, p)を式(2)として周波数センブランス法と同様の解法で行う．

2． 解析方法
本研究では速度解析法として，周波数センブラM
^




R (, , p )  
R jj ( , , p ) 

 j1




ンス法(FSM)と最尤推定法(MLM)を用いる．
センブランス法は，加算平均波形のエネルギー^
R ij ( , , p ) ^ 1/ 2
^ 1/ 2 R ii ( , , p ) R jj (, , p ) と各波形エネルギーの平均値と比が最大になる位相速度を検索する方法である．センブランス値は次式で定義されている．
 N 1 x k ( t , Z k ) t    pZ  k  0  pZ k  TS(, p)  1/ M(2)3．計測実験

2

実際に使用する試料に圧縮波・せん断波を入射しk  pZ k  T N 1N

  x (t, Z )k2

波形を計測する．計測に使用した試料にはソリッ(1)k

ドリジットポリウレタンフォームを使用した．こt    pZ k k  0

ここに， p は位相速度の逆数，t は時刻，  は計の素材は 96.0~99.0%に細かな独立気泡を一様に算開始時刻，N は波形データ数，k は波形の番号，含んでいる．空隙率の異なる 3 種類

zk は送受信機関距離，T は解析窓の幅，xk (t, zk )は(56.5%,42.0%,27.5%)の試料を用意した．寸法が縦波形データを表す．一方，センブランスの周波数4cm×横 4cm で長さが 2cm，3cm，4cm，5cm，表現は次式で与えられる．

6cm，7cm，8 ㎝の長方体である．これに周期 2µsの両振れの矩形波を入射信号として，試料に伝播

(, , p)  e  (, p)R 1 (, , p)e(, p)

させ，計測に使用する試料の長さを変えて，各空ここに，隙率につき 7 本の波形を取得した．実験の様子をFig.1 に示す．

R (, , p)  X(, , p)X (, , p) X 0 (, , p) X (, , p)   X (, , p)
N
1 e jpZ
e(, , p)   
e jpZ


0N 1




Fig.1 超音波トランスデューサ及び実験試料(左：圧縮波計測時，右：せん断波計測時)

2圧縮波速度の推定には試料長さ 2~8cm の波形データを，せん断波速度の推定には試料長さ 2~4cmの波形データを用い，周波数センブランス法と最尤推定法により速度を推定した．Fig5~8 に解析結果を示す．縦軸速度，横軸周波数として表示をしている．白い領域はセンブランス値が高く，黒い領域はセンブランス値が低い．これらの図から切り出した波形群の相関ピークの速度を抽出した．
圧縮波の速度推定結果を表 1，せん断波の速度推定結果を表 2 に示している．せん断波の解析結果
で Fig8 のように同じ時間位置に 2 つの速度のピークが存在したので速度の高い算出結果を( )で表 2に示した．
Fig.2 計測波形(圧縮波)(空隙率：42.0%，パルス幅：2µs)

Fig.4 ガイド波波形図(空隙率 42.0% )Table1 周波数センブランス法による速度推定結果

Fig.3 計測波形(せん断波)(空隙率：42.0%，パルス幅：2[µs] )計測した圧縮波を Fig.2 に，せん断波を Fig．

空隙率27.5%

42.0%

56.5%

圧縮波速度[m/s]1763

16141432

せん断波速度[m/s]

850784738

3 に示す．
ガイド波についても計測実験を行い，波形デ

Table2 最尤推定法による速度推定結果ータを取得した．ソリッドリジットポリウレタン
フォームの中実丸棒にパルス幅 10[µs]の方形パルスを入力信号として，発振点から 100~150[mm]離れた位置でセンサの位置を 10[mm] ずつずらして， 6 本の波形を取得した．(Fig.4)
3

空隙率27.5%

42.0%

56.5%

圧縮波速度[m/s]1763

16141432

せん断波速度[m/s]

850784738(959)

(874)

( － )

Velocity[m/s]

1900

10.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0

1200

500
0.045

0.06250.08

Time[ms]
Fig.5 周波数センブランス解析結果(圧縮波)
Fig.8 最尤推定法解析結果(せん断波)

G 0  E 0
G 0E0
3(3G 0  E 0 )K0 
E  E01

Fig.6 周波数センブランス解析結果(せん断波)
1900

Velocity[m/s]11 P(11K 0  5E 0 )( 9 K 0  E 0 ) P

(3)
(4)

4 K 0 ( 27 K 0  5E 0 )

(5)

1 P
G  G06( K 0  2G 0 ) P1
9 K 0  8G 0

(6)

1 P
K  K0
3K 0 P
1
4G 0

(7)

0.9

4G3


(8)

G


(9)

K

0.8

Cp 

0.7
0.6
0.5

1600

Cs 

0.4
0.30.2

ρ:多孔質材料の密度，P:空隙率，α: 比例定数

0.1
0E0 :中実材料の縦弾性係数，

1300
0.05

0.06

0.07

E:多孔質材料の縦弾性係数，

Time[ms]

G0 :中実材料の剛性率，Fig.7 最尤推定法解析結果(圧縮波)

G:多孔質材料の剛性率，
K0 :中実材料の体積弾性率

4． 弾性係数の算出
密度推定の評価には多孔質体の力学的性質を表

K:多孔質材料の体積弾性率す式を用いる（式(3)~(9)）．中実材料の弾性係数とCp :P 波速度，Cs :S 波速度

求める多孔質材料の空隙率を式に代入することで多孔質材料の弾性係数について求める 4)．

式(3)に示したように剛性率 G と縦弾性係数 E4は比例関係にある．密度推定に利用するときには周波数センブランス法で推定した速度で算出しこの比例定数は，空隙率の高い試料で計測を行った理論分散曲線を①，最尤推定法により推定したたとしても一定であると仮定している．各空隙率

表 2 の圧縮波速度と表 2( )のせん断波速度を用いで圧縮波・せん断波速度を求めたのでそれぞれのて算出した理論分散曲線を②とする．

空隙率で速度から弾性係数を算出し，その比例定比較を行うために計測したガイド波について周数が一定になるか確認した．周波数センブランス

波数センブランス解析と最尤推定法を用いて解析法による速度推定結果より算出した弾性系数を表をし，周波数ごとに抽出したピークを Fig9a と3，最尤推定法による速度推定結果より算出した弾Fig9b に示す．

性系数を表 4 に示す．
表 2 のせん断波速度を用いて比例定数を算出したために表 3，表 4 の結果に差が出ている．2 つのF2速度が出ているためどちらが正しい速度か確かめ

L1

るために圧縮波・せん断波の速度から理論分散曲線を導出し，ガイド波速度推定の結果と比較をすることにした．せん断波速度が正しければ抽出したピークと理論分散曲線は重なり合うと考える．

F1

5．多孔質丸棒の理論分散曲線

P 波:1614m/s
S 波:784m/s

理論分散曲線の算出には，中実丸棒の固体導波路上における円柱座標系の波動方程式からもとめた L モード，F モードの特性方程式を使用する．

Fig.9a ガイド波解析結果とHTable3

周波数センブランス法による速度推定

理論分散曲線①との比較値より算出した弾性係数
空隙率

27.5%

42.0%

56.5%

剛性率[GPa]

0.971

0.894

0.864

縦弾性係数[GPa]

2.729

2.564

2.455

比例定数

0.35600.34880.35192000
2000

F218001800

H
HH HHH
HHH

1600
1600
1400
1400

L112001200

Table4 最尤推定法法による速度推定値より

1000
1000

算出した弾性係数

800
800
600
600

F1

空隙率

27.5%

42.0%

56.5%

400
400

剛性率[GPa]

0.857

0.883

0.910

200
200

(1.261)(1.146)

( － )

00

2.462

2.542

2.561

(3.347)(3.097)( － )Fig.9b ガイド波解析結果と0.3482

0.34720.3554

理論分散曲線②との比較(0.3766)(0.3699)( － )縦弾性係数[GPa]比例定数

S 波:874m/s

005P 波:1622m/s

101020203030404050506060707080809090100100

は理論分散曲線と抽出したピークがほぼ重なり合っている一方．②では理論分散曲線からずれた位置でピークが抽出されている．②を算出するときのせん断波速度は正しいせん断波ではないために，ピーク位置がずれる．に基づく多孔質丸棒の密度推定”，計測自動制御学会東北支部(2013)[3] Kai Hsu,Arthur B.Baggeroer”Application of the maximum-likelihood
method(MLM) for sonic velocity logging”[4] 近 藤 連 一 “ 多 孔 材 料 ‐ 性 質 と 利 用 ,pp165-183,1973”

6．結言
多孔質材料の圧縮波とせん断波計測行い，最尤推定法で正しい圧縮波とせん断波速度について求め，縦弾性係数と剛性率の比例定数を算出したところ，空隙率 56.5%でわずかであるが，0.01 ほど他の空隙率の場合での比例定数より増加した．周波数センブランス法による速度推定結果もあるが，空隙率 27.5%で速度が最尤推定法の結果よりも 50m/s 速い速度で算出している．
他の空隙率では比例定数が近い値で出ているのにも関わらず，値が大きく異なるのでせん断波の速度推定を正しく行えていない可能性がある。
そのため比例定数の変化を確かめるのには使用できないことがわかっている．
空隙率が高い試料を計測した時，比例定数が変動すると，多孔質体の力学的性質を表す式の比例定数は一定であるという仮定が成り立たなくなる．空隙率が高くなる試料の密度推定結果に誤差が生じるのは多孔質体の力学的性質を表す式によるものだと確認された．

謝辞
本研究は科学研究費補助金(課題番号：24620002)の支援を受けて行われたものである．
この場を借りて関係者各位に謝意を表する．

参考文献[1] 鎌田正洋，齊藤玄敏，”フィルタバンク法による骨を伝播するガイド波の速度分散特性の評
価”，計測自動制御学会東北支部(2011)[2] 斉藤匠，齊藤玄敏，”ガイド波の速度分散性6