下記に示すデータの度数分布表と累積度数表を求めるプログラムを作ってください。
リスト : 身長のデータ our @height = ( 148.7, 149.5, 133.7, 157.9, 154.2, 147.8, 154.6, 159.1, 148.2, 153.1, 138.2, 138.7, 143.5, 153.2, 150.2, 157.3, 145.1, 157.2, 152.3, 148.3, 152.0, 146.0, 151.5, 139.4, 158.8, 147.6, 144.0, 145.8, 155.4, 155.5, 153.6, 138.5, 147.1, 149.6, 160.9, 148.9, 157.5, 155.1, 138.9, 153.0, 153.9, 150.9, 144.4, 160.3, 153.4, 163.0, 150.9, 153.3, 146.6, 153.3, 152.3, 153.3, 142.8, 149.0, 149.4, 156.5, 141.7, 146.2, 151.0, 156.5, 150.8, 141.0, 149.0, 163.2, 144.1, 147.1, 167.9, 155.3, 142.9, 148.7, 164.8, 154.1, 150.4, 154.2, 161.4, 155.0, 146.8, 154.2, 152.7, 149.7, 151.5, 154.5, 156.8, 150.3, 143.2, 149.5, 145.6, 140.4, 136.5, 146.9, 158.9, 144.4, 148.1, 155.5, 152.4, 153.3, 142.3, 155.3, 153.1, 152.3 );
階級 度数 累積度数 ------------------------ 130 - 135 1 1 135 - 140 6 7 140 - 145 12 19 145 - 150 25 44 150 - 155 32 76 155 - 160 17 93 160 - 165 6 99 165 - 170 1 100
階級はデータの範囲を表します。この表では x cm 以上 y cm 未満を x - y で表しています。度数はその階級に出現したデータの個数です。度数を示してある表のことを「度数分布表」といいます。累積度数はその階級までの度数を全部加えたものです。累積度数を示してある表を「累積度数分布表」といいます。
問題11 のデータで、平均値と標準偏差を求めるプログラムを作ってください。データを \(x_1, x_2, \ldots , x_N\) とすると、総計量 \(T\) と平均値 \(M\) は次式で求めることができます。
平均値が同じ場合でも、データの特徴が異なる場合があります。たとえば、A = {4, 4, 5, 5, 5, 6, 6, 6, 7, 7} と B = {1, 2, 3, 4, 5, 6, 7, 8, 9, 10} の平均値は 5.5 になります。A のデータは平均値の近くに集まっていてますが、B のデータはバラバラになっていますね。統計学では、ばらつきの大きさを表すために「分散 (variance)」という値を使います。分散の定義を次に示します。
分散の定義からわかるように、平均値から離れたデータが多いほど、分散の値は大きくなります。逆に、平均値に近いデータが多くなると分散は小さな値になります。そして、分散の平方根が「標準偏差 (SD : standard deviation)」になります。
数値を格納した配列 buff の中から最大値を求める関数 maximum(buff) と最小値を求める関数 minimum(buff) を定義してください。
数値を格納した配列 buff の中から引数 x と等しい要素があるか探索する関数 member(x, buff), 位置を返す関数 position(x, buff), 個数を数える関数 count(x, buff) を定義してください。等しい要素が見つからない場合、member は偽 (0) を、position は -1 を返すものとします。
素因数分解とは、素数でない整数 (合成数) を素数の積の形に書き表すことです。たとえば、12 は 2 * 2 * 3 と素因数分解することができます。素因数分解を行う関数 prime_factor(n) を定義してください。結果は画面 (標準出力) へ出力するものとします。
次の公式を使って平方根の整数部分を求める関数 sqrt_int(n) を定義してください。
式 (1) は、奇数 1 から 2n - 1 の総和は \(n^2\) になることを表しています。式 (2) のように、整数 m の値が \(n^2\) より大きくて \((n + 1)^2\) より小さいのであれば、\(m\) の平方根の整数部分は \(n\) であることがわかります。これは \(m\) から奇数 \(1, 3, 5, \ldots, (2n - 1), (2n + 1)\) を順番に引き算していき、引き算できなくなった時点の \(\frac{2n + 1}{2} = n\) が m の平方根になります。参考文献『平方根計算法』によると、この方法を「めのこ平方」と呼ぶそうです。
sqrtInt(10) =>3 sqrtInt(100) =>10 sqrtInt(1000) =>31 sqrtInt(10000) =>100 sqrtInt(123456789) => 11111
数値を格納した配列 buff を「二分探索 (バイナリサーチ：binary searching)」する関数 binary_search(x, buff) を定義してください。
線形探索の実行時間は要素数 N に比例するので、N が大きくなると時間がかかるようになります。これに対し、二分探索は log2 N に比例する時間でデータを探すことができます。ただし、探索するデータはあらかじめ昇順に並べておく必要があります。これを「ソート (sort)」といいます。二分探索は最初にデータをソートしておかないといけないので、線形探索に比べて準備に時間がかかります。
二分探索の動作を下図に示します。
[11 22 33 44 55 66 77 88 99] key は 66
↑ 66 > 55 後半を探す
11 22 33 44 55 [66 77 88 99] 88 > 66 前半を探す
↑
11 22 33 44 55 [66 77] 88 99 77 > 66 前半を探す
↑
11 22 33 44 55 [66] 77 88 99 66 = 66 発見
↑
図 : 二分探索
二分探索は探索する区間を半分に分割しながら調べていきます。キーが 66 の場合を考えてみましょう。まず区間の中央値 55 とキーを比較します。データが昇順にソートされている場合、66 は中央値 55 より大きいので区間の前半を調べる必要はありません。したがって、後半部分だけを探索すればいいのです。
あとは、これと同じことを後半部分に対して行います。最後には区間の要素が一つしかなくなり、それとキーが一致すれば探索は成功、そうでなければ探索は失敗です。ようするに、探索するデータ数が 1 / 2 ずつ減少していくわけです。上図の場合、線形探索ではデータの比較が 6 回必要になりますが、二分探索であれば 4 回で済みます。また、データ数が 1,000,000 個になったとしても、二分探索を使えば高々 20 回程度の比較で探索を完了することができます。
「ソート (sort)」はある規則に従ってデータを順番に並べ換える操作です。たとえば、データが整数であれば大きい順に並べる、もしくは小さい順に並べます。ソートは昔から研究されている分野で、優秀なアルゴリズムが確立しています。ここでは簡単なソートアルゴリズム (遅いソート) を実際にプログラムしてみましょう。
バブルソート (buble sort) は泡がぶくぶくと浮いてくるように、いちばん小さいデータが後ろから前に浮かび上がってくるアルゴリズムです。隣接する 2 つのデータを比較して、順序が逆であれば入れ換えます。これを順番に後ろから前に行っていけば、いちばん小さなデータは頂上に浮かび上がっているというわけです。先頭が決まったならば、残りのデータに対して同じことを行えば、2 番目には残りのデータの中でいちばん小さいものが浮かび上がります。これをデータ数だけ繰り返せばソートが完了します。
9 5 3 7 6 4 8 交換しない
~~~
9 5 3 7 6 4 8 交換する
~~~
9 5 3 7 4 6 8 交換する
~~~
9 5 3 4 7 6 8 交換しない
~~~
9 5 3 4 7 6 8 交換する
~~~
9 3 5 4 7 6 8 交換する
~~~
3 9 5 4 7 6 8 いちばん小さいデータが決定する
+ 残りのデータに対して同様な操作を行う
図 : バブルソート
配列 buff をバブルソートする関数 buble_sort(buff) を定義してください。
選択ソート (selection sort) は、ソートしていないデータの中から最小値（または最大値）を見つけ、それを先頭のデータと交換する、という手順を繰り返すことでソートを行います。最初は、すべてのデータの中から最小値を探し、それを配列の先頭 buff[0] と交換します。次は、buff[1] 以降のデータの中から最小値を探し、それを buff[1] と交換します。これを繰り返すことでソートすることができます。
[9 5 3 7 6 4 8] 3 と 9 を交換する
+ +
3 [5 9 7 6 4 8] 5 と 4 を交換する
+ +
3 4 [9 7 6 5 8] 9 と 5 を交換する
+ +
3 4 5 [7 6 9 8] 7 と 6 を交換する
+ +
3 4 5 6 [7 9 8] 7 と 7 を交換する
+
3 4 5 6 7 [9 8] 9 と 8 を交換してソート終了
+ +
図 : 選択ソート
配列 buff を選択ソートする関数 select_sort(buff) を定義してください。
単純挿入ソートの考え方はとても簡単です。ソート済みの配列に新しいデータを挿入していくことでソートを行います。次の図を見てください。
[9] 5 3 7 6 4 8 5 を取り出す
[9] * 3 7 6 4 8 5 を[9]の中に挿入する
[5 9] 3 7 6 4 8 9 をひとつずらして先頭に 5 を挿入
[5 9] * 7 6 4 8 3 を取り出して[5 9]の中に挿入する
[3 5 9] 7 6 4 8 先頭に 3 を挿入
[3 5 9] * 6 4 8 7 を取り出して[3 5 9] に挿入
[3 5 7 9] 6 4 8 9 を動かして 7 を挿入
残りの要素も同様に行う
図 : 挿入ソート
最初は先頭のデータひとつがソート済みと考えて、2 番目のデータをそこに挿入することからスタートします。データを挿入するので、そこにあるデータをどかさないといけません。そこで、挿入位置を決めるため後ろから順番に比較するとき、いっしょにデータの移動も行うことにします。
配列 buff を単純挿入ソートする関数 insert_sort(buff) を定義してください。
リスト : 度数分布表と累積度数表
use strict;
use warnings;
our @height = (
・・・省略・・・
);
# 度数分布表
our @freq = (0, 0, 0, 0, 0, 0, 0, 0);
# 累積度数表
our @cum;
# 度数分布表の作成
my $low = 130.0;
my $z = 5.0;
foreach my $x (@height) {
for (my $j = 0; $j < @freq; $j++) {
if ($x < $low + $z * ($j + 1)) {
$freq[$j]++;
last;
}
}
}
# 累積度数表の作成
$cum[0] = $freq[0];
for (my $i = 1; $i < @freq; $i++) {
$cum[$i] = $cum[$i - 1] + $freq[$i];
}
# 表示
for (my $i = 0; $i < @freq; $i++) {
printf("%.1f - %.1f | ", $low + $z * $i, $low + $z * ($i + 1));
printf("%3d %3d\n", $freq[$i], $cum[$i]);
}
配列 @freq が度数分布表、@cum が累積度数表、変数 $low が階級の下限値、$z が階級の幅を表します。@freq は 0 で初期化します。度数分布表の作成は簡単です。最初の foreach で @height の要素を取り出して変数 $x にセットします。次の for ループで階級を求めます。変数 $j が階級を表し、その上限値は $low + $z * ($j + 1) で求めることができます。$x がこの値よりも小さい場合、その要素は階級 $j であることがわかります。$freq[$j] の値をインクリメントして、for ループを脱出します。
累積度数表の作成も簡単です。$cum[0] を $freq[0] で初期化します。あとは、度数分布表の値 $freq[i] と累積度数表の値 $cum[i - 1] を足し算していくだけです。最後に、for ループで 2 つの表を出力します。
実行結果は次のようになります。
$ perl yapp11.pl 130.0 - 135.0 | 1 1 135.0 - 140.0 | 6 7 140.0 - 145.0 | 12 19 145.0 - 150.0 | 25 44 150.0 - 155.0 | 32 76 155.0 - 160.0 | 17 93 160.0 - 165.0 | 6 99 165.0 - 170.0 | 1 100
平均値と標準偏差を求めるプログラムは簡単です。次のリストを見てください。
リスト : 平均値と標準偏差
use strict;
use warnings;
our @height = (
・・・省略・・・
);
# 合計値を求める
sub sum {
my $buff = shift;
my $a = 0;
foreach my $x (@$buff) {
$a += $x;
}
$a;
}
# 平均値と標準偏差
sub meansd {
my $buff = shift;
my $m = sum($buff) / @$buff;
my $s = 0;
foreach my $x (@$buff) {
$s += ($x - $m) * ($x - $m);
}
printf("mean = %.14f, sd = %f\n", $m, sqrt($s / @$buff));
}
sub meansd2 {
my $buff = shift;
my $m = 0;
my $s = 0;
for (my $i = 0; $i < @$buff; $i++) {
my $x = $buff->[$i] - $m;
$m += $x / ($i + 1);
$s += ($i * $x * $x) / ($i + 1);
}
printf("mean = %f, sd = %f\n", $m, sqrt($s / @$buff));
}
meansd(\@height);
meansd2(\@height);
プログラムは簡単なので、説明は不要でしょう。参考文献『Ｃ言語による最新アルゴリズム事典』によると、データを 1 回通読するだけで平均値と標準偏差 (分散) を求めることができるそうです。これを参考にプログラムしたものが関数 meansd2 です。
$ perl yapp12.pl mean = 150.62699999999998, sd = 6.433473 mean = 150.627000, sd = 6.433473
平均値は 150.63 cm で、標準偏差は 6.43 cm になりました。
リスト : 最大値と最小値
use strict;
use warnings;
our @height = (
・・・省略・・・
);
# 最大値を求める
sub maximum {
my $buff = shift;
my $m = $buff->[0];
for (my $i = 1; $i < @$buff; $i++) {
$m = $buff->[$i] if $m < $buff->[$i];
}
$m;
}
# 最小値を求める
sub minimum {
my $buff = shift;
my $m = $buff->[0];
for (my $i = 1; $i < @$buff; $i++) {
$m = $buff->[$i] if $m > $buff->[$i];
}
$m;
}
print maximum(\@height), "\n";
print minimum(\@height), "\n";
$ perl yapp13.pl 167.9 133.7
プログラムは簡単です。@buff の先頭要素を変数 $m に格納します。そして、for ループで 1 番目から順番に要素を取り出して、変数と比較します。$buff->[$i] が $m よりも大きい (または小さい) 場合は変数の値を $buff->[$i] に書き換えます。最後に変数 $m の値を返すだけです。
リスト : 配列の線形探索
use strict;
use warnings;
# $x が @$buff にあるか
sub member {
my ($x, $buff) = @_;
foreach my $y (@$buff) {
return 1 if $x == $y;
}
0;
}
# $x の位置を求める
sub position {
my ($x, $buff) = @_;
for (my $i = 0; $i < @$buff; $i++) {
return $i if $x == $buff->[$i];
}
-1;
}
# $x の個数を求める
sub count {
my ($x, $buff) = @_;
my $c = 0;
foreach my $y (@$buff) {
$c++ if $x == $y;
}
$c;
}
my $a = [1, 2, 3, 1, 2, 3, 4, 1, 2, 3, 4, 5];
if (member(4, $a)) {
print "found\n";
} else {
print "not found\n";
}
if (member(6, $a)) {
print "found\n";
} else {
print "not found\n";
}
print position(1, $a), "\n";
print position(5, $a), "\n";
print position(6, $a), "\n";
print count(1, $a), "\n";
print count(5, $a), "\n";
print count(6, $a), "\n";
$ perl yapp14.pl found not found 0 11 -1 3 1 0
プログラムは簡単です。配列の先頭から順番にデータを比較していくだけです。これを「線形探索 (linear searching)」といいます。member はデータを発見したとき return で真 (1) を返します。for ループを終了したとき、データが見つからなかったので偽 (0) を返します。position はデータが見つかったとき、要素の位置 (添字) を返します。見つからなかった場合は -1 を返します。count はデータを見つけたら変数 $c を +1 します。最後に $c を返すだけです。
リスト : 素因数分解
use strict;
use warnings;
sub prime_factor {
my $n = shift;
# 2 で割り算する
while ($n % 2 == 0) {
print 2, " ";
$n /= 2;
}
# 奇数で割り算する
for (my $i = 3; $i * $i <= $n; $i += 2) {
while ($n % $i == 0) {
print $i, " ";
$n /= $i;
}
}
print $n if $n > 1;
print "\n";
}
prime_factor(1234567890);
prime_factor(1111111111);
最初に 2 で割り算します。それから、奇数で割り算していきます。割り算するときは、その数で割り切れるあいだは割り算を続けることに注意してください。たとえば、27 を素因数分解すると 3 * 3 * 3 になりますが、3 を一回だけしか割り算しないと、結果は 3 * 9 のように素数ではない数が含まれてしまいます。
実行結果は次のようになります。
$ perl yapp15.pl 2 3 3 5 3607 3803 11 41 271 9091
どの数も素数で、2 * 3 * 3 * 5 * 3607 * 3803 = 1234567890 になり、11 * 41 * 271 * 9091 = 1111111111 になります。なお、これはとても単純なアルゴリズムなので、大きな整数の素因数分解には適していません。巨大な合成数の素因数分解はとても難しい問題です。興味のある方は素因数分解について調べてみてください。
リスト : めのこ平方
use strict;
use warnings;
sub sqrt_int {
my $n = shift;
my $m = 1;
for (; $n >= $m; $m += 2) {
$n -= $m
}
int($m / 2);
}
print sqrt_int(10), "\n";
print sqrt_int(100), "\n";
print sqrt_int(1000), "\n";
print sqrt_int(10000), "\n";
print sqrt_int(123456789), "\n";
3 10 31 100 11111
アルゴリズムをそのままプログラムしただけなので、とくに難しいところは無いと思います。この方法はとても簡単ですが、数が大きくなると時間がかかるようになります。そこで、整数を 2 桁ずつ分けて計算していくことにします。次の図を見てください。
整数 6789 を 67 と 89 に分ける 1 + 3 + ... + 15 = 82 < 67 両辺を 100 倍すると 802 < 6700 < 6789 802 = 1 + 3 + ... + 159 (= 2 * 80 - 1) 161 + 163 < (6789 - 6400 = 389) < 161 + 163 + 165
整数 6789 を 67 と 89 に分けます。最初に 67 の平方根を求めます。この場合は 8 になり、82 < 67 を満たします。次に、この式を 100 倍します。すると、802 < 6700 になり、6700 に 89 を足した 6789 も 802 より大きくなります。802 は 1 から 159 までの奇数の総和であることはすぐにわかるので、6789 - 6400 = 389 から奇数 161, 163, ... を順番に引き算していけば 6789 の平方根を求めることができます。
プログラムは次のようになります。
リスト : めのこ平方 (改良版)
use strict;
use warnings;
sub sqrt_int_sub {
my ($n, $m) = @_;
for (; $n >= $m; $m += 2) {
$n -= $m
}
int($m / 2);
}
sub sqrt_int1 {
my $n = shift;
if ($n < 100) {
sqrt_int_sub($n, 1);
} else {
my $m = 10 * sqrt_int1(int($n / 100));
sqrt_int_sub($n - $m * $m, 2 * $m + 1);
}
}
print sqrt_int1(10), "\n";
print sqrt_int1(100), "\n";
print sqrt_int1(1000), "\n";
print sqrt_int1(10000), "\n";
print sqrt_int1(123456789), "\n";
3 10 31 100 11111
sqrt_int1 は $n の平方根の整数部分を求めます。$n が 100 未満の場合は sqrt_int_sub で平方根を求めます。これが再帰呼び出しの停止条件になります。$n が 100 以上の場合は、$n の下位 2 桁を取り除いた値 ($n / 100) の平方根を sqrt_int1 で求め、その値を 10 倍して変数 $m にセットします。そして、sqrt_int_sub で $n - $m * $m から奇数 2 * $m + 1, 2 * $m + 3 ... を順番に引き算していって n の平方根を求めます。
リスト : 二分探索
use strict;
use warnings;
sub binary_search {
my ($x, $buff) = @_;
my $low = 0;
my $high = @$buff - 1;
while ($low <= $high) {
my $mid = $low + int(($high - $low) / 2);
if ($buff->[$mid] == $x) {
return 1;
} elsif ($buff->[$mid] < $x) {
$low = $mid + 1;
} else {
$high = $mid - 1;
}
}
0;
}
my $a = [10, 20, 30, 40, 50, 60, 70, 80, 90];
print binary_search(0, $a), "\n";
print binary_search(10, $a), "\n";
print binary_search(40, $a), "\n";
print binary_search(90, $a), "\n";
print binary_search(100, $a), "\n";
最初に、探索する区間を示す変数 $low と $high を初期化します。$low を 0 に、最後の要素の位置 (@$buff - 1) を $high にセットします。次の while ループで、探索区間を半分ずつに狭めていきます。まず、区間の中央値を求めて変数 $mid にセットします。if 文で $mid の位置にある要素と $x を比較し、等しい場合は探索成功です。return で真 (1) を返します。
$x が大きい場合は区間の後半を調べます。変数 $low に $mid + 1 をセットします。逆に、$x が小さい場合は前半を調べるため、変数 $high に $mid - 1 をセットします。これを二分割できる間繰り返します。$low が $high より大きくなったら分割できないので繰り返しを終了し偽 (0) を返します。
実行結果を示します。
$ perl yapp17.pl 0 1 1 1 0
二分探索はデータを高速に探索することができますが、あらかじめデータをソートしておく必要があります。このため、途中でデータを追加するには、データを挿入する位置を求め、それ以降のデータを後ろへ移動する処理が必要になります。つまり、データの登録には時間がかかるのです。
したがって、二分探索はプログラムの実行中にデータを登録し、同時に探索も行うという使い方には向いていません。途中でデータを追加して探索も行う場合は、他の高速な探索アルゴリズムを検討してみてください。
リスト : バブルソート
use strict;
use warnings;
sub buble_sort {
my $buff = shift;
for (my $i = 0; $i <@$buff; $i++) {
for (my $j = @$buff - 1; $j > $i; $j--) {
if ($buff->[$j] <$buff->[$j - 1]) {
my $temp = $buff->[$j];
$buff->[$j] = $buff->[$j - 1];
$buff->[$j - 1] = $temp;
}
}
}
}
$a = [5,6,4,7,3,8,2,9,1,0];
buble_sort($a);
print "@$a", "\n";
最初のループで @$buff 回だけ繰り返します。2 番目のループで $buff の後ろから前に向かって、確定していないデータを比較していき、もしも順番が逆になっていたら交換します。とても簡単ですね。
実行結果を示します。
$ perl yapp18.pl 0 1 2 3 4 5 6 7 8 9
リスト : 選択ソート
use strict;
use warnings;
sub select_sort {
my $buff = shift;
for (my $i = 0; $i < @$buff - 1; $i++) {
my $min = $buff->[$i];
my $n = $i;
for (my $j = $i + 1; $j < @$buff; $j++) {
if ($buff->[$j] < $min) {
$min = $buff->[$j];
$n = $j;
}
}
$buff->[$n] = $buff->[$i];
$buff->[$i] = $min;
}
}
$a = [5,6,4,7,3,8,2,9,1,0];
select_sort($a);
print "@$a", "\n";
最初のループで変数 $i の値は 0 から @$buff - 1 まで動きます。2 番目のループで $buff-6gt;[$i] から buff->[$@buff - 1] までの中から最小値を探し、それを $buff->[$i] と交換します。最小値とその位置は変数 $min と $n に格納します。2 番目のループを終了したら、$buff->[$i] と $buff->[$n] の値を交換します。
実行結果を示します。
$ perl yapp19.pl 0 1 2 3 4 5 6 7 8 9
選択ソートも実行時間はデータの個数の 2 乗に比例します。
リスト : 単純挿入ソート
use strict;
use warnings;
sub insert_sort {
my $buff = shift;
for (my $i = 1; $i < @$buff; $i++) {
my $j;
my $temp = $buff->[$i];
for ($j = $i - 1; $j >= 0 && $temp < $buff->[$j]; $j--) {
$buff->[$j + 1] = $buff->[$j];
}
$buff->[$j + 1] = $temp;
}
}
$a = [5,6,4,7,3,8,2,9,1,0];
insert_sort($a);
print "@$a", "\n";
最初のループで挿入するデータを選びます。ソート開始時は先頭のデータひとつがソート済みと考えるるので、2 番目のデータ（添字では 1）を取り出して挿入していきます。2 番目のループで挿入する位置を探しています。探索は後ろから前に向かって行います。このとき、挿入位置の検索と同時にデータの移動も行っています。ループが最後まで回れば、そのデータは先頭に挿入されることになります。
実行結果を示します。
$ perl yapp20.pl 0 1 2 3 4 5 6 7 8 9
挿入ソートも実行時間はデータの個数の 2 乗に比例します。
ところで、単純挿入ソートは面白い特徴があり、ソート済みのデータは高速にソートできることが知られています。データがソートされていれば、2 番目のループは繰り返しを行わずに終了するため、最初のループの繰り返し回数でソートが完了することになります。したがって、与えられたデータが大まかにでもソートされていれば、2 番目のループで繰り返す回数が少なくなり、挿入ソートでも高速にソートすることができます。