数値計算法

数値計算法

2008 4/23

林田 清

（大阪大学大学院理学研究科）

実験データの統計処理 その１

„

誤差について

„

母集団と標本

„

平均値と標準偏差

„

誤差伝播

„

最尤法

„

平均値につく誤差

誤差

(Error)：真の値からのずれ

„

測定誤差

„ 物差しが曲がっていた „ 測定する対象が室温が低いため縮んでいた „ 1gの単位までしかデジタル表示されない計りで1g以下 „ 計りの目盛りを読み取る角度によって値が異なる „

統計誤差

„ 放射線源を検出器で測定したときの計数率 „ テレビの視聴率 „

偶然誤差

(Random Error)と系統誤差(Systematic

Error)

測定値ｘの分布 (ヒストグラム） 0 1 2 3 4 5 6 0 1 2 3 4 5 6 7 8 9 10 ｘ 頻度

測定値の分布

„

n個の測定値 x

₁

, x

₂

, …, x

_n

の分布

„ 例えば „ 1本の棒の長さをn人の人が同じものさしを 使って測定する „ （同じ設計で製作した）n本の棒の長さ „ 1個の放射線源について１分間あたりの放射 線の検出個数をn回測定する „ ある振り子の振動周期をn回測定する „ 分布の広がりが誤差を表す

母集団と標本

„

母集団

„

同じ条件で無限回の測定を繰り返したときの測

定値の分布

_{(極限頻度分布）}

„ 実際には無限回の測定は不可能 „ 極限頻度分布は存在すると仮定する „

測定は母集団から標本を採取する操作

„

採集された標本から母集団の分布を推定す

るのが統計的解析

„

真の値は不可知

平均値、標準偏差

„

その他、中央値、最頻値

1 2 1 2 2 1 1 2 2 1

n

, ,....,

1

1

(

)

1

1

lim

1

lim

(

)

n n i i n i i n i n i n i n i

x x

x

x

x

n

s

x

x

n

x

n

x

n

σ

μ

σ

μ

= = →∞ ₌ →∞ ₌

≡

≡

−

−

=

=

−

∑

∑

∑

∑

２

回の（独立な）測定

各々の誤差は

標本の平均値

標本の分散（=標準偏差 ）　

母集団の平均　

母集団の分散

(

)

(

)

(

)

1 2 2 2 2 2 2 2 2 2 2 1 1 ( ) 1 1 2 2 1 1 ( ) ( ) 2 2 1 ( ) ( ) 2( )( ) 2 1 ( 1) 1 1 ( ) ( 1) 2 n i i i j i j ij i j i j i j i j i j n n ij n ij i j i j n i i x x n x x x x x x x x x x x x x x x x x x x x s n n x x n n = = = ≠ = ≡ ⎛ ⎛ + ⎞⎞ ⎛ ⎛ + ⎞⎞ Δ ≡ _⎜ −_{⎜ ⎟⎟} +_⎜ −_{⎜ ⎟⎟} ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ = − = − − − = − + − − − − Δ = Δ − = − −

∑

∑ ∑

∑

平均　　　　　　　 二項間の分散の和　 の平均　　　　　

(

)

( )

(

)

2 2 1 2 1 1 1 2 1 ( ) 2( )( ) 1 ( ) ( 1) ( 1) 1 ( ) ( 1) n j i j j n n n i i j i i j n i i x x x x x x n x x x x x x n n n n x x n = = = = = + − − − − = − − − − − − = − −

∑

∑

∑

∑

∑

標本の分散（標準偏差

2

_{） （なぜ}

_{n-1で割るのか？）}

„

（不偏）分散

s

_n2 „

標準偏差

s

_n

誤差伝播１

2 2 1 2 2 1 2 2 2 2 ( , ,...) 1 lim ( ) ( ) ( ) 1 lim ( ) ( ) 1 lim ( ) ( ) 2( )( n x i n i i i i n x i i n i i i i i n x f u v x x n x x x x u u v v u v x x u u v v n u v x x u u v v u u v v n u v σ σ →∞ ₌ →∞ ₌ →∞ = ⎡ ⎤ = _⎢ − _⎥ ⎣ ⎦ ∂ ∂ ⎛ ⎞ ⎛ ⎞ − ≈ − _{⎜ ⎟}+ − _{⎜ ⎟} + ∂ ∂ ⎝ ⎠ ⎝ ⎠ ⎡ _⎛ ∂ _{⎞ ⎛} ∂ _⎞ ⎤ ≈ _⎢ − _{⎜ ⎟}+ − _{⎜ ⎟} + _⎥ ∂ ∂ ⎝ ⎠ ⎝ ⎠ ⎣ ⎦ ∂ ∂ ⎛ ⎞ ⎛ ⎞ = − _{⎜ ⎟} + − _{⎜ ⎟} + − − ∂ ∂ ⎝ ⎠ ⎝ ⎠

∑

∑

" "

[

]

1 2 2 2 2 1 1 1 2 2 2 2 2 ) 1 1 lim ( ) , lim ( ) 1 lim ( )( ) (covariance) n i n n u i v i n n i i n uv i i n i x u v uv x x u v u u v v n n u u v v n x x u v σ σ σ σ σ σ σ = →∞ ₌ →∞ ₌ →∞ ₌ ⎡ _{⎛ ∂ ⎞⎛ ∂ ⎞ +} ⎤ ⎢ _{⎜ ∂ ⎟⎜ ∂ ⎟} ⎥ ⎝ ⎠⎝ ⎠ ⎢ ⎥ ⎣ ⎦ ⎡ ⎤ ⎡ ⎤ = _⎢ − _⎥ = _⎢ − _⎥ ⎣ ⎦ ⎣ ⎦ ⎡ ⎤ ≡ _⎢ − − _⎥ ⎣ ⎦ ∂ ∂ ⎛ ⎞ ⎛ ⎞ ≈ _{⎜ ⎟} + _{⎜ ⎟} + + ∂ ∂ ⎝ ⎠ ⎝ ⎠

∑

∑

∑

∑

" " 共分散 x x u v ∂ ∂ ⎛ ⎞⎛ _{⎞ +} ⎜ _∂ ⎟⎜ _∂ ⎟ ⎝ ⎠⎝ ⎠ " „ 測定値u,vの関数としてｘが定 義されているとき、ｘの誤差は u,vの測定誤差からどう計算（伝 播）されるか

誤差伝播２

[

]

2 1 2 2 2 2 2 1 lim n ( )( ) uv i i n i x u v uv u u v v n x x x x u v u v σ σ σ σ σ →∞ ₌ ⎡ ⎤ ≡ _⎢ − − _⎥ ⎣ ⎦ ∂ ∂ ∂ ∂ ⎛ ⎞ ⎛ ⎞ ⎛ ⎞⎛ ⎞ ≈ _{⎜ ⎟} + _{⎜ ⎟} + + _{⎜ ⎟⎜ ⎟}+ ∂ ∂ ∂ ∂ ⎝ ⎠ ⎝ ⎠ ⎝ ⎠⎝ ⎠

∑

" " „

uとｖが独立のとき（相関がないとき）、共分

散

σ

_uv

はゼロ

2 2 2 2 2 x u v x x u v σ ≈ σ ⎛_⎜ ∂ ⎞_⎟ +σ ⎛_⎜ ∂ ⎞_⎟ + ∂ ∂ ⎝ ⎠ ⎝ ⎠ "

誤差伝播３

„

足し算、引き算 。。。誤差は同じ

„

バックグランドの引き算で誤差が大きくなる

„

かけ算

„

相対誤差の大きい成分が全体の誤差を決める

2 2 2 x u v x u v x u v σ σ σ = + = − = + あるいは　　 　 2 2 2 2 2 2 2 2 2 2 2 u v x u v x uv v u u v u v σ σ σ σ σ = ⎛ ⎞ = + = _⎜ + _⎟ ⎝ ⎠ 　

平均値の誤差

_{(Error)、不確かさ(Uncertainty)}

„

測定を

N回繰り返して平均を取ることで、（偶

然）誤差を１

_{/√nに小さくできる}

1 2 1 2 2 2 2 1

n

, ,....,

1

1

1

n n i i n x x i

x x

x

x

x

n

n

n

σ

σ

σ

σ

σ

= =

≡

=

=

∑

∑

誤差伝播則を使

回の（独立な）測定

各々の誤差は

標本平均値

標本平均値の誤差

は

うと

最尤法

(Maximum Likelihood

Method)

1 2 2

n

, ,....,

μ

Gauss)

(

)

1

1

exp

2

2

'

'

n i i i i i i

x x

x

x

x

dx

dQ

Pdx

x

P

σ

μ

σ

σ

π

μ

μ

μ

+

=

⎡

_⎛

₋

_⎞

⎤

≡

_⎢

−

_⎜

_⎟

_⎥

⎝

⎠

⎢

⎥

⎣

⎦

回の（独立な）測定

で、

母集団が平均値 　　標準偏差 の正規（

分布の場合

１回の測定で　 から

　の間の値を観測する確率は

　は不可知、推定値を とする。

尤度が最大になる はどう決められるか

最尤法２

最尤法（正規分布の場合の例） 2 1 2 1 2 1 ' ' ' 1 1 ( ') exp 2 2 , ,..., ( ') ( ') ' 1 1 exp 2 2 ( ') ' i i i n n i i n n i i x x P n x x x P P x P μ σ σ μ μ σ σ π μ μ μ σ σ π μ μ μ = = = ⎡ _{⎛ − ⎞} ⎤ = _⎢− _{⎜ ⎟ ⎥} ⎝ ⎠ ⎢ ⎥ ⎣ ⎦ = ⎡ ₋ ⎤ ⎛ ⎞ ⎛ ⎞ = _{⎜ ⎟ ⎢}− _{⎜ ⎟ ⎥} ⎝ ⎠ ⎝ ⎠ ⎢_⎣ ⎥_⎦

∏

∑

平均値 、標準偏差 の正規分布を仮定すると を観測する確率は 回の測定で を観測する確率（尤度）は を最大にする が最も確からしい の推定値 考え方： 最も確率の高い標本分布（測定 値の組）が実現されているはず

最尤法３

„

最も確からしい母集団平均

(mean)の推定値は加

算平均

_(average)

2 1 2 1 1 ( ') ' 1 2 ' 0 ' 1 ' n i i n i i n i i P X x X x dX d x x n

μ

μ

σ

μ

μ

σ

μ

= = = − ⎛ ⎞ = _{⎜ ⎟} ⎝ ⎠ − ⎛ ⎞ = − _{⎜ ⎟} = ⎝ ⎠ = =

∑

∑

∑

を最大にすることは次の を最小にするのと同じ

誤差が異なるデータの場合

（重みつき平均）

2 1 1 2 ₂ 2 2 1 1 2 ' ₂ ' 1 1 ( ') exp 2 2 ' ( / ) ' ' 1 0 ' ' 2 _{(1/ )} 1 ' (1/ ) i i n n i i i _i i n n i i i i i _i i _{i i} i x x P x x x d d μ

σ

μ

μ

σ

σ

π

μ

σ

μ

μ

_μ

μ

σ

σ

σ

μ

σ

σ

= = = = ⎡ ⎤ ⎛ ⎞ _{⎛ − ⎞} ⎢ ⎥ = _{⎜⎜ ⎟⎟} − _{⎜ ⎟} ⎢ _{⎝ ⎠} ⎥ ⎝ ⎠ _{⎣ ⎦} ⎛ ⎞ ⎛ − ⎞ − = − _{⎜ ⎟} = = ⎜ ⎟ ⎝ ⎠ ⎝ ⎠ =

∑

∏

∑

∑

∑

_∑

∑

各測定値 につく誤差が異なる の場合 の最尤推定値は より また推定値 に関する誤差は

どうやって誤差を評価するか？

„

例えば使用説明書に書いてある測定器の精度を使

用するのは一般には不十分

„ Conservativeな測定値の範囲を示すには有効 „

一般に系統誤差を評価するのは困難

„ 全く独立な実験を行い結果を比較する „

測定を

同じ条件で

複数回繰り替えすことができる場

合は測定値の（標本）標準偏差が（偶然）誤差の推

定値を与える

„ 最尤法を使い誤差を推定することもできる „

統計誤差の場合、理論的に推定できることがある

„ 例）放射線源を決まった時間だけ計測する際の計数ｘは ポアソン分布に従う。 この場合統計誤差は√x。

レポート課題１（締め切りは

5/14）

„ 平均値と標準偏差を求めるプログラム „ 入力：データの数、データ „ データは以下の１０個（例えばある月の最高気温（℃）１０日分） „ 35.5,25.0,34.2,34.6,22.8,27.7,30.6,26.8,23.0,39.3 „ 出力：（標本）平均値、標準偏差 „ ソースプログラムと出力結果をメイルの本文にして khclass@ess.sci.osaka-u.ac.jp までメイルせよ。 実行形 式は添付しないこと。 メイルのタイトルは report1_学籍番 号とすること。 „ 他の人のソースプログラムと結果をそのままコピーするのは ダメ。

データ入力

„

Cの場合

int n; double a; scanf("%d %lf",&n,&a); „

Fortranの場合

c23456 integer n real*8 a read(*,*) n,a 1. 可能であれば、データをファイルに書き込み、ファイルから読 み出すようにプログラムしてみよう。使用するのは_fopen, fscanf, fclose。 リダイレクト利用するのも一案。 2. ファイルに含まれるデータの数が任意の場合にも対応できる ようにできるか？

繰り返しループ

„ Cの場合 int i,n; double a,sum; scanf("%d",&n); sum=0;

for(i=0; i<n; i++) { scanf("%lf",&a); sum = sum+a； } „ Fortranの場合 c23456 integer i,n real*8 a,sum read(*,*) n sum=0 do 10 i=1,n read(*,*) a sum = sum+a 10 continue

配列の利用

„ Cの場合

int i,n;

double a,amat[100]; scanf("%d",&n);

for(i=0; i<n; i++ ) { scanf("%lf",&a); amat[i] = a； } „ Fortranの場合 c23456 integer i,n real*8 a,amat(100) read(*,*) n sum=0 do 10 i=1,n read(*,*) a amat(i)=a 10 continue