数値解析法及び演習第十二回

[1]モンテカルロ法

乱数を用いると多重積分を容易に行うことができる．例としてN次元単位球の体積を求めてみよう．

乱数をN個用意する．各乱数 x_i (1≤i≤N) は0から1の値をとる．この乱数の組をW回生成する．生成した乱数のうち，

$\begin{displaymath} x_1^2+x_2^2+x_3^2+....+x_N^2\le 1 \end{displaymath}$

(1)

を満たす組の数をwとすると，N次元球の体積は近似的に

$\begin{displaymath} V_{N}=2^N\times {w\over W} \end{displaymath}$

(2)

と求めることができる．一方で解析解は

$\begin{displaymath} V_{2m+1}={2(2\pi)^m\over (2m+1)!!} \quad m=0,1,2,... \end{displaymath}$

(3)

$\begin{displaymath} V_{2m}={\pi^m\over m!} \quad m=1,2,3,... \end{displaymath}$

(4)

となる．

練習問題1

N=2の場合について，W=10, 10², 10³, 10⁴と変化させて誤差の変化を観察せよ．
W=10⁴とし，N=1～10までについて体積を算出せよ．体積はどう変化したか？

[3]最小自乗法II

それでは逆行列のアルゴリズムを使う練習を行う．

またまた最小自乗法をとりあげる．いままでに行ったのは，フィッティングする関数を

$\begin{displaymath} y=a+bx \end{displaymath}$

(7)

と線形に仮定してきた．今度は多項式

$\begin{displaymath} y=\sum_{i=1}^{n}a_ix^{i-1}=a_1+a_2x+...+a_nx^{n-1} \end{displaymath}$

(8)

を仮定する．ここで実験データがmこあるとすると，多項式の項数はm以下でなければならない(式の数が足りなくなるから)．

ここで係数aを決定するのが仕事となる．以前と同じように，関数のフィッティングは

$\begin{displaymath} f=\sum_{j=1}^m(y_j-\sum_{i=1}^{n}a_ix_j^{i-1})^2 \end{displaymath}$

(9)

を最小にすることで行う．そのためにはfを各aで偏微分してゼロにすれば良い．

$\begin{displaymath} {\partial f\over \partial a_k}=0 \quad 1\le k \le n \end{displaymath}$

(10)

これを計算すると

$\begin{displaymath} \sum_{i=1}^{n}a_i\sum_{j=1}^m x_j^{i+k-2}=\sum_{j=1}^mx_j^{k-1}y_j \quad 1\le k \le n \end{displaymath}$

(11)

という連立方程式を解くことになる．行列の形で書き替えると

$\begin{displaymath} \left( \begin{array}{cccc} m & \sum_{j=1}^m x_j & \ldots & \... ...c}{\dotfill}\\ \sum_{j=1}^m x_j^{n-1} y_j \end{array}\right) \end{displaymath}$

(12)

となる．

練習問題5

逆行列をもとめるアルゴリズムを用いて

x  y

1.D0 3.D0
3.D0 6.D0
5.D0 5.5D0
7.D0 5.D0
9.D0 7.D0

を二次式，三次式でフィッティングせよ．

練習問題6

このデータはある曲面から作成したものである．データを保存してgnuplotを実行し，

splot 'kyokumen.dat'

で曲面の形状を確認せよ．曲面の形状を

$\begin{displaymath} z=ax^2+bxy+cy^2 \end{displaymath}$

(13)

と仮定し，データをフィッティングするa,b,cを決定せよ．a+cは平均曲率，4ac-bはガウス曲率と呼ばれる．先ほどと同様に，

$\begin{displaymath} g=\sum_{j=1}^m(z_j-(ax_j^2+bx_jy_j+cy_j^2))^2 \end{displaymath}$

(14)

を最小にするために

$\begin{displaymath} {\partial g\over \partial a}={\partial g\over \partial b}={\partial g\over \partial c}=0 \end{displaymath}$

(15)

から構成される連立方程式を解けばよい．

[4] 多変量の相関

いくつかの特性に対するデータが得られたときそれをどう味わうかを考える．ここでは多変量の相関を考える．例として以下のデータをみてみよう．

$\begin{displaymath} \begin{array}{cccc} & {\rm math} & {\rm phys}& {\rm chem} \... ... C} & 45. & 68. & 74.\\ {\rm D} & 56. & 75. & 76. \end{array}\end{displaymath}$

(16)

このように，サンプル(Aさん，Bさん．．)を行にとり，変量(数学，理科．．)を列にとって作った測定値の行列をデータ行列と呼ぶ．多変量解析の目的は, どの変量とどの変量が相互に関係しているのかを統計的に推測することである.

一般に異なる変量は異なる次元を持つ．そのためまず各変量をおなじように扱えるようにするために標準化という作業をおこなう．これは，平均値がゼロ，分散が1となるような変換をデータに施すことである．ある変量の平均値は

$\begin{displaymath} \bar{x}_j={1\over n}\sum_{i=1}^nx_{ij} \end{displaymath}$

(17)

であたえられる．分散（標準偏差の二乗）は

$\begin{displaymath} \sigma_j^2={1\over n}\sum_{i=1}^n(x_{ij}-\bar{x}_j)^2 \end{displaymath}$

(18)

となる．データの標準化は

$\begin{displaymath} z_{ij}={x_{ij}-\bar{x}_j\over \sigma_j} \end{displaymath}$

(19)

の変換によっておこなわれる．この変換をおこなった行列を標準化行列と呼ぶ．標準化の結果, 平均値が0, 分散が1となる.

$\begin{displaymath} Z= \left( \begin{array}{cccc} z_{11} & z_{12} & \ldots & z_{... ...fill}\\ z_{n1} & z_{n2} & \ldots & z_{nm} \end{array}\right) \end{displaymath}$

(20)

ここでこの行列の転置行列をかけると

$\displaystyle R={1\over n}Z^{T}Z$	$\textstyle =$	$\displaystyle {1\over n} \left( \begin{array}{cccc} z_{11} & z_{21} & \ldots & ... ...icolumn{4}{c}{\dotfill}\\ z_{n1} & z_{n2} & \ldots & z_{nm} \end{array}\right)$
	$\textstyle =$	$\displaystyle {1\over n}\left( \begin{array}{cccc} \sum_{i=1}^{n} z_{i1}^2 & \s... ...um_{i=1}^{n} z_{im}z_{i2} & \ldots & \sum_{i=1}^{n} z_{im}^2 \end{array}\right)$	(21)

となり各要素は二つの変量の間の相関係数を与える．変量jと変量kの相関係数は

$\begin{displaymath} r_{jk}={\sigma_{jk}\over \sigma_j\sigma_k} ={\sum_{i=1}^n(x_... ..._{i=1}^n(x_{ij}-\bar{x}_j)^2\sum_{i=1}^n(x_{ik}-\bar{x}_k)^2}} \end{displaymath}$

(22)

であたえられる．二つの変量の相関がよいほど相関係数の値の絶対値は1に近くなる．行列Rの要素はそれぞれ二つの変量の相関係数のため，Rは相関行列と呼ばれる．それではでて来た相関係数はどう味わえばよいのか？それには次の表を用いる. Nはデータの数, rは相関係数である. それぞれの数値は, 無相関のN個のデータがあった時にr以上の相関係数が出現する確率（パーセント）である. この数値が小さい程, 二つのデータには相関があると思ってよい. 5以下なら有意, 1以下なら高度に有意であるという.

					r
N	0.0	0.1	0.2	0.3	0.4	0.5	0.6	0.7	0.8	0.9	1.0
3	100	94	87	81	74	67	59	51	41	29	0
6	100	85	70	56	43	31	21	12	6	1	0
10	100	78	58	40	25	14	7	2	0.5		0
20	100	67	40	20	8	2	0.5	0.1			0
50	100	49	16	3	0.4						0

無相関な二つの量がN組あったとき, 相関係数がr以上になる確率. 空白は0.05%以下.

今日の宿題 (12/23日まで)

以下の問題に回答せよ.

「健康診断データ」，「大手私鉄データ」，「株式，為替レートデータ」，「体力測定データ」を解析し，どの項目とどの項目に高度に有意な相関があるか指摘せよ．

日程表へ戻る << 　

$\displaystyle 5.x+2.y+6.z$	$\textstyle =$	$\displaystyle 7.$
$\displaystyle 2.x+2.8y+4.4z$	$\textstyle =$	$\displaystyle 7.8$	(6)
$\displaystyle 4.x+3.y+8.2z$	$\textstyle =$	$\displaystyle 6.$

数値解析法及び演習第十二回