幾種隨機算法的實現原理

轉自：http://www.androidstar.cn/%E5%87%A0%E7%A7%8D%E9%9A%8F%E6%9C%BA%E7%AE%97%E6%B3%95%E7%9A%84%E5%AE%9E%E7%8E%B0%E5%8E%9F%E7%90%86/

　　在日常工作中，經常需要使用隨機算法。比如面對大量的數據， 需要從其中隨機選取一些數據來做分析。 又如在得到某個分數后， 為了增加隨機性， 需要在該分數的基礎上， 添加一個擾動， 並使該擾動服從特定的概率分布（偽隨機）。本文主要從這兩個方面出發， 介紹一些算法， 供大家參考。

　　首先假設我們有一個使用的隨機函數float frand()， 返回值在(0, 1)上均勻分布。大多數的程序語言庫提供這樣的函數。 在其他的語言如C/C++中， 可以通過間接方法得到。如 frand()= ((float)rand() ) / RAND_MAX;

1, 隨機選取數據

　　假設我們有一個集合A(a_1,…,a_n), 對於數m，0≤m≤n, 如何從集合A中等概率地選取m個元素呢？
　　通過計算古典概率公式可以得到， 每個元素被選取的概率為m/n。 如果集合A里面的元素本來就具有隨機性， 每個元素在各個位置上出現的概率相等， 並且只在A上選取一次數據，那么直接返回A的前面m個元素就可以了， 或者可以采取每隔k個元素取一個等類似的方法。這樣的算法局限很大， 對集合A的要求很高， 因此下面介紹兩種其他的算法。

1.1 假設集合A中的元素在各個位置上不具有隨機性， 比如已經按某種方式排序了，那么我們可以遍歷集合A中的每一個元素a_i, 0<=n 根據一定的概率選取ai。如何選擇這個概率呢？

設m’為還需要從A中選取的元素個數， n’為元素a_i及其右邊的元素個數， 也即n’=(n-i+1)。那么選取元素a_i的概率為 m’/n’。
由於該算法的證明比較繁瑣， 這里就不再證明。
我們簡單計算一下前面兩個元素(2<=m<=n)各被選中的概率。
1) 設p(a_i=1)表示a_i被選中的概率。顯而易見， p(a_1=1)=m/n, p(a_1=0)為(n-m)/n;
2)第二個元素被選中的概率為
p(a_2=1)= p(a_2=1,a_1=1)+p(a_2=1,a_1=0)
= p(a_1=1)*p(a_2=1│a_1=1)+ p(a_1=0)* p(a_2=1│a_1=0)
= m/n * (m-1)/(n-1) + (n-m)/n*m/(n-1)
= m/n

我們用c++語言， 實現了上述算法

template<class T>
bool getRand(const vector vecData, int m, vector& vecRand)
{
        int32_t nSize = vecData.size();
    if(nSize  < m || m < 0)
        return false;
    vecRand.clear();
        vecRand.reserve(m);
    for(int32_t i = 0, isize = nSize; i < isize ; i++){
            float fRand = frand();
            if(fRand <=(float)(m)/nSize){
                vecRand.push_back(vecData[i]);
                m--;
            }
        nSize --;
    }
    return true;
}

　　利用上述算法， 在m=4， n=10， 選取100w次的情況下， 統計了每個位置的數被選取的概率

位置 概率
1    0.399912
2    0.400493
3    0.401032
4    0.399447
5    0.399596
6    0.39975
7    0.4
8    0.399221
9    0.400353
10    0.400196

　　還有很多其他算法可以實現這個功能。比如對第i個數， 隨機的從a_i, …, a_n中， 取一個數和a_i交換。這樣就不單獨介紹了。

1.2 在有些情況下，我們不能直接得到A的元素個數。比如我們需要從一個很大的數據文件中隨機選取幾條數據出來。在內存不充足的情況下，為了知道我們文件中數據的個數， 我們需要先遍歷整個文件，然后再遍歷一次文件利用上述的算法隨機的選取m個元素。

　　又或者在類似hadoop的reduce方法中， 我們只能得到數據的迭代器。我們不能多次遍歷集合， 只能將元素存放在內存中。 在這些情況下， 如果數據文件很大， 那么算法的速度會受到很大的影響， 而且對reduce機器的配置也有依賴。

這個時候，我們可以嘗試一種只遍歷一次集合的算法。

1. 取前m個元素放在集合A’中。
2. 對於第i個元素(i>m), 使i在 m/i的概率下， 等概率隨機替換A’中的任意一個元素。直到遍歷完集合。
3. 返回A’

下面證明在該算法中，每一個元素被選擇的概率為m/n.

1. 當遍歷到到m+1個元素時， 該元素被保存在A’中的概率為 m/(m+1), 前面m個元素被保存在A’中的概率為 1- (m/m+1 * 1/m) = m/m+1
2. 當遍歷到第i個元素時，設前面i-1個元素被保存在A’中的概率為 m/(i-1)。根據算法， 第i個元素被保存在A’中的概率為m/i , 前面i-1各個元素留在A’中的概率為 m/(i-1) * (1-(m/i* 1/m) = m/i;
3. 通過歸納，即可得到每個元素留在A’中的概率為 m/n;

我們在類似 hadoop的reduce函數中， 用java實現該算法。

public void reduce(TextPair key, Iterator value, OutputCollector collector, int m)
{
    Text[] vecData = new Text[m];
    int nCurrentIndex = 0;
    while(value.hasNext()){
        Text tValue = value.next();
        if(nCurrentIndex < m){
           vecData[nCurrentIndex] = tValue;
        }
        else if(frand() < (float)m / (nCurrentIndex+1)) {
           int nReplaceIndex = (int)(frand() * m);
           vecData[nReplaceIndex] = tValue;
        }
        nCurrentIndex ++;
    }
 //collect data
…….
}

　　利用上述算法，在m=4， n=10， 經過100w次選取之后， 計算了每個位置被選擇的選擇的概率

位置 概率
1    0.400387
2    0.400161
3    0.399605
4    0.399716
5    0.400012
6    0.39985
7    0.399821
8    0.400871
9    0.400169
10    0.399408

2. 隨機數的計算

在搜索排序中，有些時候我們需要給每個搜索文檔的得分添加一個隨機擾動， 並且讓該擾動符合某種概率分布（這還算是隨機擾動嗎O(∩_∩)O）。假設我們有一個概率密度函數f(x), min<=x<=max, 並且有

那么可以利用f(x)和frand設計一個隨機計算器r(frand()), 使得r(frand())返回的數據分布， 符合概率密度函數f(x)。
令

那么函數

符合密度函數為f(x)的分布。
下面對這個以上的公式進行簡單的證明：

由於g(x)是單調函數， 並且x在[0,1]上均勻分布，那

由於上述公式太復雜， 計算運算量大， 在線上實時計算的時候通常采用線性差值的方法。
算法為：
1）在offline計算的時候， 設有數組double A[N+1]；對於所有的i, 0<=i<=N, 令

2）在線上實時計算的時候，令

f = frand()，
lindex = (int) (f* N);
rindex = lindex +1;

那么線性插值的結果為 A[lindex]*(A[rindex]-f) + A[rindex] * (f – A[lindex])
我們做了一組實驗，令f(x)服從標准正太分布N(0,1)， N=10000， 並利用該算法取得了200*N個數。對這些數做了個簡單的統計， 得到x軸上每個小區間的概率分布圖。

3后記
　　在日常工作中， 還有其他一些有趣的算法。比如對於top 100w的query， 每個query出現的頻率不一樣， 需要從這100w個query， 按照頻率越高， 概率越高的方式隨機選擇query。限於篇幅， 就不一一介紹了。

幾種隨機算法比較