確率の宿題について

[osusi]

初めて質問させていただきます。
学校の宿題でC言語を使ってプログラミングをしてシミュレーションをするのが出たのですが、
どうも詰まってしまって、今回メールさせていただきました。
問題 最初に箱が一つ(箱１)あって、その中に玉が１個入っている。
次の時間に玉が１個増えるが、それは箱１に入るか、
新しい箱(箱２)が出現して、その中に玉が１個入るかのどちらかである。
仮に、後者の場合になったとすると、 箱が２つ、それぞれの箱の中に玉が1つ入っている。
次の時間に、玉が1個増えるが、それは箱１もしくは箱２に入るか、新たな箱(箱３)が出現してそこに玉が入るかである。
と、これを１０００回ほど繰り返して、 箱が何個あって、その中にそれぞれ玉が何個入ってるか、シミュレーションして表示せよと言う宿題です。
玉が既存の箱の中のどこにはいるか、または、新しい箱が出現してそこに入るかの確率は、各箱の中の玉の数に依存しています。
その確率はこんな感じです。 箱がｋ個(A1～Aｋ)あって、それぞれの箱の中に玉が何個か入っている(合計の玉はｎ個)。
次の時間に球が一個増えるが、それがA1～Aｋの箱に入るか、それとも新しい箱A(ｋ＋1)に入るか、
どの箱が選ばれるか、その確率レートは
箱１(玉n1個)→n1＋１個の場合　確率n1/(n+θ)
箱２(玉n2個)→n2＋１個の場合　確率n2/(n+θ) 　　　・ 　　　・ 　　　・
箱i (玉ni個)→ni＋１個の場合　確率ni/(n+θ)
　・ 　　　・ 　　　・
箱ｋ(玉nk個)→nk＋１個の場合　確率nk/(n+θ)
新しい箱が出現する確率　θ/(n+θ)
この確率を全部足すと１になります。 (θは定数です)
今プログラムを考えているのですが、 箱は配列を１００個用意して？？？。。。 １０００回の繰
り返しはfor文使って。。。 どこの箱にするかは、箱１から箱ｋまでの確率と新たな箱の出現の確率を０から１までの長さの直線にして、 それか
ら、０から１の間をとる乱数を発生させて、それがどこの箱の確率の幅に当てはまるかをプログラム(if文で？)して。。。 だけど、for文で書くか
ら毎回箱の数が増える。どうやって確率の幅のプログラムをかいたらいいのか？？？？ とC言語の基礎の本を読みながら格闘しています。 にっち
もさっちも行かないので、お知恵を貸してください。 お願いします。

[管理人]

写真１
画像クリックで拡大します。
(写真２は次の投稿)

[管理人]

osusi様改めまして、こんにちは＾＾
質問の方、掲示板にうつさせていただきました。
サンプルプログラムが書けたので紹介しますね。
まず、サンプルプログラムはこちらにございます。

http://dixq.net/board/program/num.cpp

写真１はサンプルプログラムの実行結果で、
代入過程が示されています。
写真２は最終結果の出力の様子です。
1000回目まで表示した後に、出力するようになっている様子を
お見せするため、上ちょっと切れていますけど、
実際に見ていただきたいのは
box[0]からの箱の中身です。
最終出力はbox番号とboxの中の玉の個数が示されています。
ではサンプルプログラムのアルゴリズムを紹介する前にちょっと
色々整理してみましょう。

確率的に処理をする方法はいろいろありますが、
私が考えた確率の出し方は次の通りです。
例えば4/9の確率的に処理をする方法は
x = rand()%9;
if(0>=x && 4<x)
,
この条件に当てはまらない(箱を作る時など)はelseを用いずに
if(4>=x && 9<x)
,
と書きました。
では具体的に箱とたまについて

今
box[0]と
box[1]の２つの箱がある状態だとします。

次の瞬間、
box[0]に玉が入る確率はbox[0]/(トータル個数+θ)
box[1]に玉が入る確率はbox[1]/(トータル個数+θ)
新しく箱が作られる確率はθ/(トータル個数+θ)
ですね。

変数は次のような物を使いました。
int i,j,ram,num=0,boxnum=1,box[1000];

iは1000回までループさせるforで使うカウンタで、
この値は同時に玉のトータル個数を表します。

jは0～箱があるまで処理をする時のforで使うカウンタです。
上記の箱の例では0～1を表すわけですね。
もし箱の数が４つなら0～3を表します。

ramは乱数を格納する変数です。
ram=rand()%(SITA+i);
と書くことで、iはトータル個数ですから確率の分母が作成できますね。
サンプルプログラムではこの値にちょっと都合で+1していますが
ifでnumやiと比較するためで、深い意味はありません。
(for文のiが0ではなく1から始まっているせいです)

numは、確率的処理を行うための分子の役割をします。
先ほど確率的処理の説明をしましたね。
このように
x = rand()%9;
if(0>=x && 4<x)
,
この条件に当てはまらない(箱を作る時など)はelseを用いずに
if(4>=x && 9<x)
,
この0～4の値を変数にしたのがnumです。
if(0>=x && 4<x)
は、１つめの箱に４個玉が入っているときの条件ですから
num=0;
box[0]=4;
if(num>=x && num+box[0]<x)
こういう意味になりますね。

boxnumという変数は今何個箱があるかという変数です。
boxnumの値が３なら３つ箱がある事を示します。

box[/url]は箱の配列で要素は玉の個数を代入します。
1000個もいりませんけど、一応1000個用意しました。


次に処理についての説明をします。

#define SITA 5
はθの値を５で指定してあります。
θの値を変更したい時はこの5の数字を任意に変更してください。
#define文の意味がお分かりでないばあいは聞いてください。

srand((unsigned) time(NULL));
は、乱数を取得する時の初期化を行っています。
意味はわからなくて結構です。
今の時間を初期値に、乱数がいつも違う値になるようにしているのです。

ではサンプルプログラムを上から見ていきましょう。

実際にはこちらのプログラムをご使用ください。

http://dixq.net/board/program/num.cpp

では解説します。

#include <stdio.h>
//乱数を扱うためのライブラリ
#include <stdlib.h>
//乱数を初期化するためのライブラリ
#include <time.h>
//θを5に指定。変更可能
#define SITA 5

void main(){
//乱数の初期値に今の秒数を指定
　　srand((unsigned) time(NULL));
　　int i,j,ram,num=0,boxnum=1,box[1000];
//箱が作られた時には玉が１つはいっているはずだから予め箱全てに1をいれておく。
　　for(i=0;i<1000;i++)
　　　　box=1;
//ここから1000回試行を行います
　　for(i=1;i<=1000;i++){
//乱数代入。+1の意味に深い意味はありません。
　　　　ram=rand()%(SITA+i)+1;
//今の箱の数だけ試行します
　　　　for(j=0;j<boxnum;j++){
//確率的処理です。意味は上記参照。
//num～num+boxとする事で、実際には分子がboxである事になります。
　　　　　　if(ram>num && ram<=num+box[j]){
//もし条件にあえばあった箱の中身の個数を1増やします。
　　　　　　　　box[j]++;
//printf文は省略します。
　　　　　　　　printf("..(省略)
//条件にあえばブレイクし、次の試行に移ります。
　　　　　　　　break;
　　　　　　}
//確率の分子の役割をするnumに箱の中身の個数分プラス
　　　　　　num+=box[j];
//もし新しい箱を作成する条件にあえば
　　　　　　if(ram>i && ram<=SITA+i){
//箱の個数を増やす。
　　　　　　　　boxnum++;
　　　　　　　　printf("..(省略)
　　　　　　　　break;
　　　　　　}
　　　　}
//確率の分子の役割をするnumを０に戻す
　　　　num=0;
　　}
//1000回終わったら最終結果を出力。
　　for(i=0;i<boxnum;i++)
　　　　printf("box[%d] = %d\n",i,box);
}

これで全ての処理が終わりです。
何か出力や計算で間違いがあったら教えてください。

こういうアルゴリズムを考える問題って特に人が書いたプログラムを読むのはだるいですから
アルゴリズムの説明だけきいて、自分でプログラムを書くほうがホントはいいんですが
サンプルプログラムがないと結果がつかめないかとおもい、
一応、かいて見ました。
人の書いたプログラムってみにくいですからね。
読んでも読んでもわからない場合は無理によもうとせず、
アルゴリズムを考えて自分でかいて見てください。

また、わからない部分があったら気軽に聞いてくださいね☆
がんばってください＾＾

[osusi]

ありがとうございます。
凄い早いお返事、非常にありがたいです。

このプログラムと解説をじっくり読んで自分で理解してみます。
そこで、わいてきた疑問・質問等があったらまた伺わせていただきますので、そのときはよろしくお願いします。

[管理人]

落ち着いてみたら特に難しいプログラムではないですよ。
乱数の初期化とprint文をはぶいてスマートに書くとこんな風にかけます。


void main(){
　　int i,j,ram,num=0,boxnum=1,box[1000];
　　for(i=0;i<1000;i++) box=1;
　　for(i=1;i<=1000;i++){
　　　　ram=rand()%(SITA+i)+1;
　　　　for(j=0;j<boxnum;j++){
　　　　　　if(ram>num && ram<=num+box[j]){
　　　　　　　　box[j]++; break;}
　　　　　　if(ram>i && ram<=SITA+i){
　　　　　　　　boxnum++; break;}
　　　　　　num+=box[j];
　　　　}
　　　　num=0;
　　}
　　for(i=0;i<boxnum;i++) printf("box[%d] = %d\n",i,box);
}

[osusi]

管理人さんが作成したプログラム、だいたい理解できました
(*´д｀*)

で、新たに質問があります。

確率レートを、
箱１(玉n1個) →n1+1個の場合 (n1－r)/(n+θ)
箱２(玉n2個)→n2+1個の場合 (n2－r)/(n+θ)
箱i(玉ni個)→ni+1個の場合 (ni－r)/(n+θ)
箱k(玉nk個)→nk+1個の場合 (nk－r)/(n+θ)

新たな箱が出現する確率レート (θ＋kr)/(n+θ)

合計の確率は１。

とし、新たにｒを導入した場合

プログラムは、
サンプルプログラムを使うとどこをかえればいいんでしょうか？

#define でｒを定義して、
それをif文での条件に導入したらいいと思うんですが、
どうもわからないので・・・。

教えてくださいm(__)m

[管理人]

やはり１から１０まで私が書いたのがよくなかったですかね・・。
他人の書いたプログラムを見るのは苦労するものです。

確率を計算したい時は確率をa/bとした時

x=rand()%b;
if(xが0からaの時)

としてやればいいんです。
例えば確率12/33で処理を行いたい場合、
x=rand()%33;
if(x>=0 && x<12)
としてやればその確率で処理が行われます。
このやり方が一番簡単ですね。

もしも分子が小数点になった場合は整数に直しましょう。
例えばrが0.5の時、n1が１なら小数になりますよね。
つまり
0.5/(n+θ)となってしまいます。この時は２倍してやりましょう。
1/2*(n+θ)として

if(rand()%(2*(n+θ))==0)

とする事でn1=1,r=0.5の時の箱１のP1が求まります。
確率の分母を乱数をわったあまり「%」をする事で処理が出来ます。
以前のサンプルは引用しない方がいいと思います。
上記のアルゴリズムで自分でかいて見るのが一番よくわかります。
頑張ってください。

[osusi]

なるほど！
そうやればｒが少数でも対応できますね。

ありがとうございましたm(__)m

[管理人]

こんにちは、osusi様＾＾

解決されましたか？
昨日思ったのですが、もし小数点第２位までｒを下げると、使えなくなります。
といいますのも
乱数が生成できる数は
0～32767
までなんです。
小数点第２位を整数にするには100倍しますね。そのうえ１～1000までカウントアップしますから、
乱数は10万まで必要になります。
そういう場合はrand2()等という名前で他に関数を作っておき、

void main(){
int x;
x=rand2();
}

void rand2(){
return (rand()+rand()+rand()+rand()%(100000-32767*3));
}

このようにすることでrand2は0～10万までの乱数が生成できます。
32767を3回たしても若干10万には届かず、98301になるので4回目のrandは丁度10万になるように調整したものです。

return (rand()+rand()+rand());はrand()*3としてはいけません。
なぜなら例えば2という乱数は生成されなくなります。
rand()が0なら0、rand()が1なら3になってしまいます。
もし10万以上の乱数を作る場合はfor文を使った方がいいでしょう。

あらゆる状態に常に対応できるようなプログラムを書くべきですが、
私の頭の範囲では、こういう改善策しか思いつきません。。
参考にして頂けたら幸いです。

[osusi]

何とか解決できそうです。
この乱数の扱い方も参考にしてみます。

ありがとうございます(　＾ω＾)

[管理人]

はい、頑張ってくださいね☆

[osusi]

ども。お世話になってますm(__)m

いま上の宿題を少し自分で拡張させているんですが、
そこで質問があるのでカキコしました。

この箱モデルでは、１０００回確率のさいころをふって、
箱[k]に何個ボールが入っているかを出力させてますよね。

拡張問題として
「箱[k]に何個ボールが入っているか」の他に、「１０００回の確率のさいころをふった結果、箱自体が何箱あるかを出す」。
そして、この試行を１００回行わせて、別ファイルに保存して平均・標準偏差をだそうかなぁと考えました。

で、管理人さんが考えてくれたプログラムに少し足して出来たプログラムがこれです。↓↓↓

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#define SITA 2
void main(){

FILE *fw,*fw1;

srand((unsigned) time(NULL));
int i,j,z,ram,num=0,boxnum=1,box[1000];

fw=fopen("箱の番号とその中の個数.txt","w+");
fw1=fopen("箱の数.txt","w+");

for(z=0;z<100;z++){
for(i=0;i<1000;i++)
box=1;
for(i=1;i<=1000;i++){
ram=rand()%(SITA+i)+1;
for(j=0;j<boxnum;j++){
if(ram>num && ram<=num+box[j]){
box[j]++;

break;
}
num+=box[j];
if(ram>i && ram<=SITA+i){
boxnum++;

break;
}
}
num=0;
}



for(i=0;i<boxnum;i++){


printf("box[%d] = %d\n",i,box);
fprintf(fw,"%d\t %d\n",i,box);

}

fprintf(fw1,"%d\n",boxnum);

}

fclose(fw);
fclose(fw1);
}


試行を１００回繰り返すのは、新たにfor文を付け加えました(zのところ)。
また、箱の数を数えるプログラムの書き方が分からないので、計算後のboxnumの値が箱の数と等しい(あってるかな？)と思い、それを付け加えました。

で、テストとして、今１００回繰り返しではなく、３・４回の繰り返しでやっているんですが、どうも出た結果がおかしいのです。

繰り返しの数が多いほど、箱の数が増えてしまう結果になってしまいます。

繰り返しが１回だと、箱の数が大体１５箱ぐらいでおさまるのですが、回数を４回とかにすると、４回目の試行の箱の数が４０とか５０とかそれぐらいになってしまいます。

このプログラムで何かおかしいところってありますか？

[管理人]

boxnumが箱の数ですが、
１回の試行が全て完了した時にこの箱の数を０で初期化していないからだと思います。
１００回試行するのでしたら１回ごとに試行が終わるたびに箱の数を０に戻してやらないと
２回目に試行が始まった時、１回目の箱の個数から始まってしまいます。
ですから箱の数が段々増えるのでしょう。
まだプログラムを10秒位みただけなんでよくわかりませんが＾＾；

ちゃんと見て見ますね。

[管理人]

うん、やっぱそうだと思いますよ。
でも０で初期化じゃないですね、
最初箱が１つあるはずなんで、１で初期化しないといけませんでした。
最後、箱の数を書き出した後かまたはforで一番最初に
boxnum=1;
を書けば解決します。

[osusi]

出来ました～(　＾ω＾)

ありがとうございますm(__)m

一つでもプログラムが違うと変な値がはじき出されるから、
C言語は大変です( ´д｀).∵