1. 程式人生 > >概率論中指數分佈介紹及C++11中std::exponential_distribution的使用

概率論中指數分佈介紹及C++11中std::exponential_distribution的使用

指數分佈:在深度學習中,我們經常會需要一個在x=0點處取得邊界點(sharp point)的分佈。為了實現這一目的,我們可以使用指數分佈(exponential distribution):

p(x;λ)= λlx≥0exp(-λx)

指數分佈使用指示函式(indicator function) lx≥0來使得當x取負值時的概率為零。

指數分佈:在概率論和統計學中,指數分佈(Exponential distribution)是一種連續概率分佈。指數分佈可以用來表示獨立隨機事件發生的時間間隔,比如旅客進入機場的時間間隔、打進客服中心電話的時間間隔、中文維基百科新條目出現的時間間隔等等。

一個指數分佈的概率密度函式是:

其中λ>0是分佈的一個引數,常被稱為率引數(rate parameter)。即每單位時間發生該事件的次數。指數分佈的區間是[0,∞)。如果一個隨機變數X呈指數分佈,則可以寫作:X∽Exponential(λ)。

指數分佈的累積分佈函式可以寫成:

隨機變數X(X的率引數是λ)的期望值是:E[X]=1/λ,比方說:如果你平均每個小時接到2次電話,那麼你預期等待每一次電話的時間是半個小時。X的方差是:D[X]=1/λ2,X的偏離係數是:V[X]=1.

以下簡單介紹Laplace分佈、Dirac分佈、經驗分佈、混合分佈:

Laplace分佈:一個聯絡緊密的概率分佈是Laplace分佈(Laplace distribution),它允許我們在任意一點μ處設定概率質量的峰值:

拉普拉斯分佈:在概率論與統計學中,拉普拉斯分佈是以皮埃爾-西蒙·拉普拉斯的名字命名的一種連續概率分佈。由於它可以看作是兩個不同位置的指數分佈背靠背拼接在一起,所以它也叫做雙指數分佈。兩個相互獨立同概率分佈指數隨機變數之間的差別是按照指數分佈的隨機時間布朗運動,所以它遵循拉普拉斯分佈。

如果隨機變數的概率密度函式分佈為:

那麼它就是拉普拉斯分佈。其中,μ是位置引數,b>0是尺度引數。如果μ=0,那麼,正半部分恰好是尺度為1/2的指數分佈。

拉普拉斯分佈的概率密度函式讓我們聯想到正態分佈,但是,正態分佈是用相對於μ平均值的差的平方來表示,而拉普拉斯概率密度函式用相對於平均值的差的絕對值來表示。因此,拉普拉斯分佈的尾部比正態分佈更加平坦。

Boost 1.64.0庫中實現了對Laplace分佈的實現,在math/distributions/laplace.hpp檔案中,對應類為laplace_distribution。

Dirac分佈:在一些情況下,我們希望概率分佈中的所有質量都集中在一個點上。這可以通過Dirac delta函式(Dirac delta function)δ(x)定義概率密度函式來實現:p(x)= δ(x-μ)

Dirac delta函式被定義成在除了0以外的所有點的值都為0,但是積分為1。Dirac delta函式不像普通函式一樣對x的每一個值都有一個實數值的輸出,它是一種不同型別的數學物件,被稱為廣義函式(generalized function),廣義函式是依據積分性質定義的數學物件。我們可以把Dirac delta函式想成一系列函式的極限點,這一系列函式把除0以外的所有點的概率密度越變越小。

通過把p(x)定義成δ函式左移-μ個單位,我們得到了一個在x=μ處具有無限窄也無限高的峰值的概率質量。

狄拉克δ函式:在科學和數學中,狄拉克δ函式或簡稱δ函式是在實數線上定義的一個廣義函式或分佈。它在除零以外的點上都等於零,且其在整個定義域上的積分等於1。δ函式有時可看作是在原點處無限高、無限細,但是總面積為1 的一個尖峰,在物理上代表了理想化的質點或點電荷的密度。

從純數學的觀點來看,狄拉克δ函式並非嚴格意義上的函式,因為任何在擴充套件實數線上定義的函式,如果在一個點以外的地方都等於零,其總積分必須為零。δ函式只有在出現在積分以內的時候才有實質的意義。根據這一點,δ函式一般可以當做普通函式一樣使用。它形式上所遵守的規則屬於運算微積分的一部分,是物理學和工程學的標準工具。

δ函式的圖形通常可以視為整條x軸和正y軸。雖然稱為函式,但δ函式並非真正的函式,至少它的值域不在實數以內。

Dirac分佈經常作為經驗分佈(empirical distribution)的一個組成部分出現:

經驗分佈將概率密度1/m賦給m個點x(1),…,x(m)中的每一個,這些點是給定的資料集或者取樣的集合。只有在定義連續型隨機變數的經驗分佈時,Dirac delta函式才是必要的。對於離散型隨機變數,情況更加簡單:經驗分佈可以被定義成一個Multinoulli分佈,對於每一個可能的輸入,其概率可以簡單地設為在訓練集上那個輸入值的經驗頻率(empirical frequency)。

當我們在訓練集上訓練模型時,我們可以認為從這個訓練集上得到的經驗分佈指明瞭我們取樣來源的分佈。關於經驗分佈另外一個重要的觀點是,它是訓練資料的似然最大的那個概率密度函式。

Empirical distribution function: In statistics, an empirical distribution function is the distribution function associated with the empirical measure of a sample. This cumulative distribution function is a step function that jumps up by 1/n at each of the n data points. Its value at any specified value of the measured variable is the fraction of observations of the measured variable that are less than or equal to the specified value.

The empirical distribution function is an estimate of the cumulative distribution function that generated the points in the sample. It converges with probability 1 to that underlying distribution, according to the Glivenko–Cantelli theorem.A number of results exist to quantify the rate of convergence of the empirical distribution function to the underlying cumulative distribution function.

分佈的混合:通過組合一些簡單的概率分佈來定義新的概率分佈也是很常見的。一種通用的組合方法是構造混合分佈(mixture distribution)。混合分佈由一些元件(component)分佈構成。每次實驗,樣本是由哪個元件分佈產生的取決於從一個Multinoulli分佈中取樣的結果:這裡P(c)是對各元件的一個Multinoulli分佈。

一個非常強大且常見的混合模型是高斯混合模型(Gaussian Mixture Model),它的元件p(x|c=i)是高斯分佈。每個元件都有各自的引數,均值μ(i)和協方差矩陣∑(i)。有一些混合可以有更多的限制。和單個高斯分佈一樣,高斯混合模型有時會限制每個元件的協方差矩陣為對角的或者各向同性的(標量乘以單位矩陣)。

除了均值和協方差以外,高斯混合模型的引數指明瞭給每個元件i的先驗概率(prior probability)αi=P(c=i)。”先驗”一詞表明瞭在觀測到x之前傳遞給模型關於c的信念(belief)。作為對比,P(c|x)是後驗概率(posterior probability),因為它是在觀測到x之後進行計算的。高斯混合模型是概率密度的萬能近似器(universal approximator),在這種意義上,任何平滑的概率密度都可以用具有足夠多元件的高斯混合模型以任意精度來逼近。

Mixture distribution: In probability and statistics, a mixture distribution is the probability distribution of a random variable that is derived from a collection of other random variables as follows: first, a random variable is selected by chance from the collection according to given probabilities of selection, and then the value of the selected random variable is realized. The underlying random variables may be random real numbers, or they may be random vectors (each having the same dimension), in which case the mixture distribution is a multivariate distribution.

以下是對C++11中的指數分佈模板類std::exponential_distribution的介紹:

C++11中引入了指數分佈(Exponential distribution)模板類std::exponential_distribution.指數分佈是一個連續概率分佈,產生隨機非負浮點數,可以用來表示獨立隨機事件發生的時間間隔。
std::exponential_distribution:Random number distribution that produces floating-point values according to an exponential distribution, which is described by the following probability density function: p(x|λ)= λe-λx, x>0
This distribution produces random numbers where each value represents the interval between two random events that are independent but statistically defined by a constant average rate of occurrence (its lambda, λ).The distribution parameter, lambda, is set on construction. To produce a random value following this distribution, call its member function operator().Its analogous discrete distribution is the geometric_distribution. 

 下面是從其它文章中copy的std::exponential_distribution測試程式碼,詳細內容介紹可以參考對應的reference:

#include "exponential_distribution.hpp"
#include <iostream>
#include <random>
#include <string>
#include <chrono>
#include <thread>
#include <iomanip>
#include <map>

///////////////////////////////////////////////////////////
// reference: http://www.cplusplus.com/reference/random/exponential_distribution/
int test_exponential_distribution_1()
{
{
	const int nrolls = 10000;  // number of experiments
	const int nstars = 100;    // maximum number of stars to distribute
	const int nintervals = 10; // number of intervals

	std::default_random_engine generator;
	std::exponential_distribution<double> distribution(3.5);

	int p[nintervals] = {};

	for (int i = 0; i<nrolls; ++i) {
		double number = distribution(generator);
		if (number<1.0) ++p[int(nintervals*number)];
	}

	std::cout << "exponential_distribution (3.5):" << std::endl;
	std::cout << std::fixed; std::cout.precision(1);

	for (int i = 0; i<nintervals; ++i) {
		std::cout << float(i) / nintervals << "-" << float(i + 1) / nintervals << ": ";
		std::cout << std::string(p[i] * nstars / nrolls, '*') << std::endl;
	}
}

{ // (1)、exponential_distribution::exponential_distribution: Construct exponential distribution,
//   Constructs an exponential_distribution object, adopting the distribution parameters specified either by lambda or by object parm.
//   (2)、exponential_distribution::lambda: Returns the parameter lambda (λ) associated with the exponential_distribution.
//   This parameter represents the number of times the random events are observed by interval, on average.
//   This parameter is set on construction.
//   (3)、exponential_distribution::max: Maximum value
//   Returns the least upper bound of the range of values potentially returned by member operator().
//   (4)、exponential_distribution::min: Minimum value
//   Returns the greatest lower bound of the range of values potentially returned by member operator(),
//   which for exponential_distribution is always zero.
//   (5)、exponential_distribution::operator(): Generate random number
//   Returns a new random number that follows the distribution's parameters associated to the object (version 1)
//   or those specified by parm (version 2).

	// construct a trivial random generator engine from a time-based seed:
	int seed = std::chrono::system_clock::now().time_since_epoch().count();
	std::default_random_engine generator(seed);

	std::exponential_distribution<double> distribution(1.0);

	std::cout << "ten beeps, spread by 1 second, on average: " << std::endl;
	for (int i = 0; i<10; ++i) {
		double number = distribution(generator);
		std::chrono::duration<double> period(number);
		std::this_thread::sleep_for(std::chrono::seconds(1)/*period*/);
		std::cout << "beep!" << std::endl;
	}

	std::cout << "lambda: " << distribution.lambda() << std::endl;
	std::cout << "max: " << distribution.max() << std::endl;
	std::cout << "min: " << distribution.min() << std::endl;
}

{ // exponential_distribution::param: Distribution parameters
//   The first version(1) returns an object with the parameters currently associated with the distribution object.
//   The second version(2) associates the parameters in object parm to the distribution object.
	std::default_random_engine generator;
	std::exponential_distribution<double> d1(0.8);
	std::exponential_distribution<double> d2(d1.param());

	// print two independent values:
	std::cout << d1(generator) << std::endl;
	std::cout << d2(generator) << std::endl;
}

{ // exponential_distribution::reset: Resets the distribution, so that subsequent uses of the object do not depend on values already produced by it.
//   This function may have no effect if the library implementation for this distribution class produces independent values.
	std::default_random_engine generator;
	std::exponential_distribution<double> distribution(1.0);

	// print two independent values:
	std::cout << distribution(generator) << std::endl;
	distribution.reset();
	std::cout << distribution(generator) << std::endl;
}

	return 0;
}

////////////////////////////////////////////////////////////////////
// reference: http://en.cppreference.com/w/cpp/numeric/random/exponential_distribution
int test_exponential_distribution_2()
{
	std::random_device rd;
	std::mt19937 gen(rd());

	// if particles decay once per second on average,
	// how much time, in seconds, until the next one?
	std::exponential_distribution<> d(1);

	std::map<int, int> hist;
	for (int n = 0; n<10000; ++n) {
		++hist[2 * d(gen)];
	}
	for (auto p : hist) {
		std::cout << std::fixed << std::setprecision(1)
			<< p.first / 2.0 << '-' << (p.first + 1) / 2.0 <<
			' ' << std::string(p.second / 200, '*') << '\n';
	}

	return 0;
}

////////////////////////////////////////////////////////////////
// reference: https://msdn.microsoft.com/en-us/library/bb982914.aspx
static void test(const double l, const int s)
{
	// uncomment to use a non-deterministic generator
	//    std::random_device gen;
	std::mt19937 gen(1701);

	std::exponential_distribution<> distr(l);

	std::cout << std::endl;
	std::cout << "min() == " << distr.min() << std::endl;
	std::cout << "max() == " << distr.max() << std::endl;
	std::cout << "lambda() == " << std::fixed << std::setw(11) << std::setprecision(10) << distr.lambda() << std::endl;

	// generate the distribution as a histogram  
	std::map<double, int> histogram;
	for (int i = 0; i < s; ++i) {
		++histogram[distr(gen)];
	}

	// print results  
	std::cout << "Distribution for " << s << " samples:" << std::endl;
	int counter = 0;
	for (const auto& elem : histogram) {
		std::cout << std::fixed << std::setw(11) << ++counter << ": "
			<< std::setw(14) << std::setprecision(10) << elem.first << std::endl;
	}
	std::cout << std::endl;
}
int test_exponential_distribution_3()
{
	double l_dist = 0.5;
	int samples = 10;

	std::cout << "Use CTRL-Z to bypass data entry and run using default values." << std::endl;
	std::cout << "Enter a floating point value for the 'lambda' distribution parameter (must be greater than zero): ";
	//std::cin >> l_dist;
	std::cout << "Enter an integer value for the sample count: ";
	//std::cin >> samples;

	test(l_dist, samples);

	return 0;
}
GitHubhttps://github.com/fengbingchun/Messy_Test