2019/6/17

ベイズ更新の例題

梁のヤング率推定
信頼性工学特論
レポート課題例

2019/6/17

概要• 梁の固有振動数の計測データからヤング率を求める• 固有振動数の計測データはモンテカルロ法で生成する• 実測データではなく、シミュレートする• ベイズ更新によりヤング率の確率分布を推定する
𝐸: ヤング率、𝜌: 密度ℎ𝐿𝑏

2019/6/17

ベイズ更新• 𝐷：データ、𝐻：パラメータに対して
尤度

事前分布

𝑝 𝐷𝐻 𝑝 𝐻
𝑝 𝐻𝐷 =𝑝 𝐷

事後分布• 𝐷：計測した固有振動数𝑓• 𝐻：対数正規分布するヤング率𝐸の平均𝜇 や分散𝜎
• 今回はヤング率𝐸の平均𝜇 だけを考える2019/6/17

数値実験によるデータ生成
MaterialMember

𝐸

𝜌𝑏ℎ𝐿

[N/m2][kg/m3]

[m]

[m]

[m]

Average

2.05×1011

7850

0.75010.00201

1.829

S.D.

2.00×109

50.00150.00001

0.001𝐼𝐴

[m4]

[m2]

5.08×10-100.001508

Average2019/6/17

乱数発生
Excel関数を用いて材料・部材パラメータを生成=LOGNORM.INV(RAND(),平均,分散)
Sample12345678
910Erho2.07.E+11 7849.735
2.06.E+11 7853.153
2.05.E+11 7848.674
2.07.E+11 7852.375
2.03.E+11 7852.570
2.06.E+11 7855.607
2.01.E+11 7845.315
2.04.E+11 7853.192
2.05.E+11 7845.399
2.06.E+11 7845.064

IAbhL
7.498.E-01 2.005.E-03 1.829471
7.462.E-01 2.023.E-03 1.8306097.474.E-01 2.010.E-03 1.829343
7.527.E-01 2.015.E-03 1.829732
7.482.E-01 2.022.E-03 1.828159
7.506.E-01 2.008.E-03 1.829068
7.501.E-01 2.003.E-03 1.829525
7.507.E-01 2.027.E-03 1.828176
7.492.E-01 2.001.E-03 1.829364
7.515.E-01 1.998.E-03 1.830082

freq

Noisefreq

2019/6/17

乱数発生𝑏ℎ𝐼=12Sample
12345678910

𝐴 = 𝑏ℎ𝜋𝑓=𝐿ErhoIAbh
L
2.07.E+11 7849.735 5.035E-10 1.503.E-03 7.498.E-01 2.005.E-03 1.829471
2.06.E+11 7853.153 5.148E-10 1.510.E-03 7.462.E-01 2.023.E-03 1.830609
2.05.E+11 7848.674 5.055E-10 1.502.E-03 7.474.E-01 2.010.E-03 1.829343
2.07.E+11 7852.375 5.134E-10 1.517.E-03 7.527.E-01 2.015.E-03 1.829732
2.03.E+11 7852.570 5.152E-10 1.513.E-03 7.482.E-01 2.022.E-03 1.828159
2.06.E+11 7855.607 5.063E-10 1.507.E-03 7.506.E-01 2.008.E-03 1.829068
2.01.E+11 7845.315 5.025E-10 1.503.E-03 7.501.E-01 2.003.E-03 1.829525
2.04.E+11 7853.192 5.212E-10 1.522.E-03 7.507.E-01 2.027.E-03 1.828176
2.05.E+11 7845.399 5.005E-10 1.499.E-03 7.492.E-01 2.001.E-03 1.829364
2.06.E+11 7845.064 4.993E-10 1.501.E-03 7.515.E-01 1.998.E-03 1.830082

freq
1.40
1.40
1.39
1.40
1.39
1.39
1.37
1.40
1.39
1.39

𝐸𝐼
𝜌𝐴
Noise

freq

2019/6/17

乱数発生
𝐷：データ

観測データには更にノイズを加える

Sample12345678910ErhoIAbhL
2.07.E+11 7849.735 5.035E-10 1.503.E-03 7.498.E-01 2.005.E-03 1.829471
2.06.E+11 7853.153 5.148E-10 1.510.E-03 7.462.E-01 2.023.E-03 1.830609
2.05.E+11 7848.674 5.055E-10 1.502.E-03 7.474.E-01 2.010.E-03 1.829343
2.07.E+11 7852.375 5.134E-10 1.517.E-03 7.527.E-01 2.015.E-03 1.829732
2.03.E+11 7852.570 5.152E-10 1.513.E-03 7.482.E-01 2.022.E-03 1.828159
2.06.E+11 7855.607 5.063E-10 1.507.E-03 7.506.E-01 2.008.E-03 1.829068
2.01.E+11 7845.315 5.025E-10 1.503.E-03 7.501.E-01 2.003.E-03 1.829525
2.04.E+11 7853.192 5.212E-10 1.522.E-03 7.507.E-01 2.027.E-03 1.828176
2.05.E+11 7845.399 5.005E-10 1.499.E-03 7.492.E-01 2.001.E-03 1.829364
2.06.E+11 7845.064 4.993E-10 1.501.E-03 7.515.E-01 1.998.E-03 1.830082

freq
1.40
1.40
1.39
1.40
1.39
1.39
1.37
1.40
1.39
1.39

Noise0.000455-0.001333
0.000392
0.000811
0.000858-0.001364-0.000106
0.001370-0.002126-0.000589

freq
1.396
1.401
1.392
1.401
1.396
1.391
1.373
1.403
1.385
1.385

2019/6/17

尤度• 𝑝 𝐷 𝐻 ：𝐻のとき𝐷が得られる確率• 正規分布確率を仮定（分散：𝑠 ）

𝑝 𝑓𝜇=

12𝜋𝑠

𝑓−
exp −𝜋𝐿𝜇 𝐼
𝜌𝐴2𝑠尤度は求めたいパラメータ（ここでは𝐸の確率分布の平均𝜇 ）の関数になる

2019/6/17

尤度とパラメータ
尤度は得られたデータ𝑓が、ある𝜇 の値に対してどれくらい尤もらしいかを表す
𝜋𝐿

𝜇 =3が3のとき

𝜇 𝐼
𝜌𝐴

𝑓は3になりやすい𝜇 =2が2のとき

31

𝑓は2になりやすい
2が1のとき𝜇 =1𝑓は1になりやすい
112312𝑓

3𝑓23𝑓

3.00.E+112.90.E+112.80.E+112.70.E+112.60.E+112.50.E+112.40.E+112.30.E+112.20.E+112.10.E+112.00.E+111.90.E+111.80.E+111.70.E+111.60.E+112019/6/17

事前分布の仮定• 𝜇 に関する事前情報は無いので一様分布と定める
1.200

1.000

0.800

0.600

0.400

0.200

0.000

𝜇2019/6/17

事後分布の計算• 尤度：𝑝 𝑓 𝜇=は𝜇 の関数

exp −• 𝑓は1つ目のデータ＝1.396等を代入する• 尤度関数の標準偏差： 𝑠は適当に設定する（たとえば0.2とした）
2.50.E+00
2.00.E+00
1.50.E+00
1.00.E+00
5.00.E-01

3.00.E+112.90.E+112.80.E+112.70.E+112.60.E+112.50.E+112.40.E+112.30.E+112.20.E+112.10.E+112.00.E+111.90.E+111.80.E+111.70.E+110.00.E+00

1.60.E+11• 尤度関数

𝜇2.50.E+00

1.000
2.00.E+00

0.800
1.50.E+00

0.600

0.400
1.00.E+00

0.200
5.00.E-01

0.000
0.00.E+00

𝑝 𝜇=1𝑝 𝜇 𝑓 =12𝜋𝑠
1.396 −

exp −𝜋𝐿2𝑠
𝜇 𝐼
𝜌𝐴

1.60.E+111.70.E+111.80.E+111.90.E+112.00.E+112.10.E+112.20.E+112.30.E+112.40.E+112.50.E+112.60.E+112.70.E+112.80.E+112.90.E+113.00.E+111.2001.60.E+111.70.E+111.80.E+111.90.E+112.00.E+112.10.E+112.20.E+112.30.E+112.40.E+112.50.E+112.60.E+112.70.E+112.80.E+112.90.E+113.00.E+111.60.E+111.70.E+111.80.E+111.90.E+112.00.E+112.10.E+112.20.E+112.30.E+112.40.E+112.50.E+112.60.E+112.70.E+112.80.E+112.90.E+113.00.E+112019/6/17ベイズ更新• ヤング率𝐸の平均𝜇 の確率分布
𝑝 𝜇 𝑓 =
12𝜋𝑠

exp −1.396 −
𝜋𝐿2𝑠

2.00.E+00

1.50.E+00

1.00.E+00

0.00.E+00

5.00.E-01

𝜇 𝐼
𝜌𝐴×12𝜋𝑠exp −1.401 −

4.50.E+00

4.00.E+00

3.50.E+00

3.00.E+00

2.50.E+00・・・・・・

𝜋𝐿2𝑠

𝜇 𝐼
𝜌𝐴

3.00.E+112.90.E+112.80.E+112.70.E+112.60.E+112.50.E+112.40.E+112.30.E+112.20.E+112.10.E+112.00.E+111.90.E+111.80.E+111.70.E+111.60.E+11

2019/6/17

結果• 10個の計測データについてベイズ更新を繰り返す
1.20.E+03

1.00.E+03

8.00.E+02

6.00.E+02

4.00.E+02

2.00.E+02

0.00.E+00• 数値実験における𝜇 の正解値：2.05×1011……精度良好