Octave のプログラミング要素 スクリプトファイルの作成

5
S 科 環境情報学 ’05
Octave のプログラミング要素
繰り返し Octave では次の２つの構文により，繰り返し
処理（ループ）を行うことができます．
• while 条件式 実行文 end[while]
• for 要素指定 実行文 end[while]
while ループは，C 言語とほとんど変わりませんが．for
ループはベクトルを用いて要素を柔軟に指定することがで
きます．典型的な例は，次のようなもので，k の値を変化
させて，その 2 乗を表示します．
for k = 1:0.5:10
k**2
endfor
与える要素は等差数列でなくてもよいので，
for k = [0 10 2 3 7 pi]
k**2
endfor
などと指定することも可能です．また，for ループを入れ子
にすることも他のプログラミング言語と同様に可能です．
例えば，2 次元行列の要素を 0 で初期化する場合には，素
直に
$ emacs foo.oct
‡‡ Memo ‡‡
行列を初期化する方法として，for ループを用いる
☞
方法および最後の要素を初期化する方法を用いた場合の実
行時間を比較する次のスクリプトを emacs で作成し，実
行してみなさい．
N = 100;
M = 200;
tic
for n = 1:N
for m = 1:M
A(n,m) = 0;
endfor
endfor
toc
tic; B(N,M) = 0; toc
7 tic と toc は，その間に記述された文の実行時間を表示します．
Box-Müller 変換により生成した正規乱数と，randn()
☞
により生成した正規乱数を，ヒストグラムを描いて比較す
る次のスクリプトを動かしてみなさい．
N = 100;
M = 200;
for n = 1:N
for m = 1:M
A(n,m) = 0;
endfor
endfor
と記述できます．しかし，このような for ループは，要素が
増える度に行列の再定義を実行するようなもので，octave
では大変手間がかかります．単に
N = 100;
M = 200;
A(N,M) = 0;
とすべきです．
スクリプトファイルの作成
長いプログラムを対話的に実行すると，間違えた場合
にやり直しが大変です．そのような場合には，スクリプト
ファイル（名前を foo.oct とします) を，エディターで作
成し，シェルのコマンドラインから octave に直接
$ octave -q foo.oct
と実行させます．
Unix 系 OS では，Emacs（イーマックス）が普及してい
ますので，それを使ってください．foo.oct を作成するに
は，シェルのコマンドラインから次のように起動します．
目次へ
### boxmuller-ts.oct ###
axis([-5,5,0,0.25])
X = randn(50000);
[X1,X2] = boxmuller(50000);
c = [-5:0.5:5];
hold on
hist(X,c,1);
hist(X1,c,1);
hist(X2,c,1);
pause
7 関数ファイル（boxmuller.m）を先に作成しておきなさい．
### boxmuller.m ###
function [r1,r2] = boxmuller(N)
R1 = rand(N,1);
R2 = rand(N,1);
r1 = sqrt(-2*log(R1)).*cos(2*pi*R2);
r2 = sqrt(-2*log(R1)).*sin(2*pi*R2);
endfunction
6
S 科 環境情報学 ’05
3
2 つの確率変数 X, Y の平均値と分散については次の性
質があります．ここに, a, b, c は定数です．
確率分布
ある変数 X が確率変数であるとは，それがとりうる値
x と x となる確率 p が同時に定められた場合をいいます．
確率変数 X が値 x をとる確率を改めて記号 P (X = x) で
表わし，これを確率分布と呼びます．また，確率変数 X が
ある値 x までをとる確率 F (x) = P (X ≤ x) を分布関数と
呼びます．
サイコロの目のように確率変数が離散的な値である場合
についてはこの講義では省略して，身長や体重のように確
率変数が連続的な場合について統計学的な特性をまとめて
おきましょう．
連続的な場合には，任意の a, b (a < b) に対して X が区
間 [a, b] の間にある確率は，確率密度関数f (x) の積分とし
て次のように与えられます．
Z
P (a ≤ x ≤ b) =
b
f (x) dx
(9)
a
E(aX + bY + c) = aE(X) + bE(Y ) + c
(17)
V (aX + b) = a2 V (X)
(18)
また X, Y が独立の場合には
E(X · Y ) = E(X) · E(Y )
(19)
V (X + Y ) = V (X) + V (Y )
(20)
が成り立ちます．
3.1.1
一様分布
区間 (0, 1) の一様乱数を与える確率密度関数 f (x) は，
一様という定義より，一定値 f0 であるから，
Z
もちろん，確率密度関数 f (x) は
Z
f (x) dx = 1, f (x) ≥ 0
P (0 < x < 1) =
Z
1
dx = f0
0
(10)
Ω
1
dx = f0 × 1 = 1
0
すなわち f (x) = 1 (0 < x < 1) です．したがって，平均
値と分散は
という性質があります．したがって連続な場合の分布関数
F (x) = P (X ≤ x) は，
Z x
F (x) = P (X ≤ x) =
(11)
f (x)dx
¸1
1
x2
= =µ
2
2
0
0
µ ¶2
Z 1
1
1 1
1
V (X) =
x2 dx −
= − =
2
3 4
12
0
Z
E(X) =
xmin
·
1
x dx =
と表現され，以下の基本的な性質を満たします．
‡‡ Memo ‡‡
( i ) F (xmin ) = 0, F (xmax ) = 1
( ii ) F (x) は単調非減少関数である
(iii) 右連続である
☞
区間 (0, 1) の一様分布の平均値と分散を求めて，理
論的な値と比較しなさい．
さらに，次の関係があることも自明です．
P (a < x ≤ b) = F (b) − F (a)
(12)
また F (x) が微分可能であれば，次の関係が成立します．
f (x) =
3.1
dF (x)
dx
(13)
連続的な確率分布の平均値と分散
確率変数 X の定義域を Ω，その確率密度関数を f (x) と
するとき，平均値 E(X) および分散 V (X) は次のように定
義されます．
Z
(14)
µ = E(X) =
xf (x) dx
Ω
Z
(15)
σ 2 = V (X) = (x − µ)2 f (x) dx
ZΩ
=
x2 f (x) dx − µ2
(16)
Ω
目次へ