1. 程式人生 > >C語言下泊松分佈以及指數分佈隨機數生成器實現

C語言下泊松分佈以及指數分佈隨機數生成器實現

最近實驗室的專案需要實現模擬檔案訪問序列,要求單位時間內的資料請求次數符合泊松分佈,而兩次請求見的時間間隔符合指數分佈。沒辦法只好重新撿起已經丟掉多時的概率知識。於是也就有了這篇關於在C語言下符合泊松分佈和指數分佈的隨機數生成器的實現。

泊松分佈

在實際的事例中,當某一事件,比如進站乘客數量,電話交換機接收到的通話請求以固定的瞬時速率λ獨立且隨機地出現時,就可以認為該事件在單位時間內發生的次數符合泊松分佈。

首先必須由二項分佈引出:

如果做一件事情成功的概率是 p 的話,那麼獨立嘗試做這件事情 n 次,成功次數的分佈就符合二項分佈。展開來說,在做的 n 次中,成功次數有可能是 0 次、1 次 …… n次。成功 i 次的概率是:

( n 中選出 i 項的組合數) * p ^ i * (1-p)^ (n-i)

以上公式很容易推導,用一點概率學最基本的知識就夠了。因為每一特定事件成功的概率是 p ,不成功的概率是 1-p 。i 次成功的事件可以任意分佈在總共的 n 次嘗試中。把它們乘起來就是恰好成功 i 次的概率。

當我們把二項分佈推而廣之後,就可以得到波鬆分佈。

可以這樣考慮,在一個特定時間內,某件事情會在任意時刻隨機發生(前提是,每次發生都是獨立的,且跟時間無關)。當我們把這個時間段分成非常小的時間片構成時,可以認為,每個時間片內,該事件可能發生,也可能不發生。幾乎可以不考慮發生多於一次的情況(因為時間片可被分的足夠小)。

當時間片分的越小,該時間片內發生這個事件的概率 p 就會成正比的減少。即:特定時間段被分成的時間片數量 n 與每個時間片內事件發生的概率 p 的乘積 n*p 為一個常數。這個常數表示了該事件在指定時間段發生的頻度。

回過頭來再來看這段時間內,指定事件恰好發生 i 次的概率是多少?代入上面推匯出來的公式得到:

n * (n-1)... (n-i+1) / i! * p^i * (1-p) ^ (n-i) => np(np-p)...(np-ip+p) / i! * ((1-p) ^ (-1/p))^(-np) / (1-p) ^i

當 n 趨向無窮大時,p 趨向 0 。而此時 (1-p)^(-1/p) 趨向 e 。

上面這個公式可以劃簡為

其中λ 事件發生的平均速率。

有了這些知識,那麼寫出符合泊松分佈的隨機數生成器程式也就不難了:

int possion(int Lambda)
{
	int  k = 0;
	long double p = 1.0;
	long double l = exp(-Lambda);
        srand((unsigned)time(NULL));
        //printf("%1.5Lf\n", l);
	while(p>=l)
	{
		double u = U_Random();
		p *= u;
		k++;
	}
	return k-1;
}

double U_Random()
{
	double f;
	f = (float)(rand() %100);
	return f/100;
}

指數分佈

指數分佈是一種連續的概率分佈,通常可以用來描述連續的獨立隨機事件發生的時間間隔。其概率密度函式為:

公式

其中λ是事件發生的平均速率。

對於較為簡單的概率分佈,其隨機數生成器可以通過其對應的累積分佈函式的逆函式構造。具體的證明可以參考這篇文章:點選開啟連結。對於指數分佈來說,其累積分佈函式為:

那麼其逆函式則為:


所以其隨機數生成器程式為:

double randomExponential(double lambda)
{
	double pv = 0.0;
	pv = (double)(rand()%100)/100;
	while(pv == 0)
	{
		pv = (double)(rand() % 100)/100;
	}
	pv = (-1  / lambda)*log(1-pv);
	return pv;
}