伽爾頓板實驗小球分布的研究

晉宏營 劉美云

（榆林學院能源工程學院，陜西 榆林 719000）

1 引言

伽爾頓板實驗可以形象地說明大數目隨機事件中的統計規律，以及統計規律中伴隨的漲落現象.伽爾頓板裝置是在一塊豎直木板的上部規則地釘上許多釘子，木板的下部用豎直隔板隔成許多等寬的狹槽，從板頂漏斗形的入口處可以投入小球，板前覆蓋玻璃，以使小球留在狹槽內［1］.實驗表明：當從入口處投入一個小球時，小球最后落入哪個狹槽是偶然的；當投入大量小球時，可看到最后落入各狹槽的小球數目不相同，在中央的槽內小球數目最多，離中央越遠的槽內小球越少；當小球數目較多時，重復該實驗，每次得到的小球分布彼此近似地重合［1，2］.伽爾頓板實驗中大量小球的分布服從一定的統計規律，近似于正態分布，但由于伽爾頓板左右側面的阻擋限制，該分布的范圍與正態分布有差別，不是從-∞到＋∞.伽爾頓板實驗中小球分布的函數解析式是什么，教材和其他文獻中沒有給出［1～5］.我們使用最大熵原理研究了伽爾頓板實驗中小球的分布，得到了小球分布的概率密度函數解析表達式；我們還在計算機上對該實驗進行了蒙特卡羅模擬，并把模擬得到的小球分布與導出的小球理論分布進行了比較.

2 伽爾頓板實驗中小球理論分布的推導

最大熵原理是統計物理中的一個基本原理，它指出：一個宏觀系統的信息熵（廣義熵）在一組約束條件下趨于約束極大值.按照此原理，對于一個宏觀系統，如果我們選擇合適的約束條件，利用拉格朗日乘子法等方法計算其信息熵的約束極大值，原則上可以求出該系統的分布［6］.作為自然界的一個基本規律，最大熵原理已在很多領域得到廣泛應用［6～8］，下面我們使用最大熵原理對伽爾頓板實驗中小球分布的具體表達式進行推導.

圖1 伽爾頓板實驗裝置

伽爾頓板實驗裝置如圖1所示，圖中黑點代表釘子，下面是狹槽.設入口處相對于狹槽的高度為h，以板底中心為坐標原點，沿板底為x軸建立坐標系，見圖1，原點到板底兩端的距離均為L.設小球落在坐標x處的概率密度為f（x），即落在區間x—x＋dx之間的概率為f（x）dx，由概率歸一化條件，可得

根據最大熵原理，信息熵S定義為

從入口處投入小球，則小球在下落過程中先后與許多釘子碰撞，最后落入某一狹槽.設各個小球落在板底的位置到入口處的距離平方的平均值為C，則有

按照最大熵原理，小球在板底的分布應使得信息熵S在約束條件式（1）和式（3）下取得極大值，這類約束極值問題可使用拉格朗日乘子法解決.根據拉格朗日乘子法，引入函數：

式中，α是由約束條件式（1）引入的拉格朗日乘子；β是由約束條件式（3）引入的拉格朗日乘子.

由δF［f（x）］＝0，可計算得信息熵S在約束條件（1）和式（3）下取極大值的概率密度函數f（x）為

把式（5）代入式（1），得

由式（6）得

與不同x值對應的誤差函數erf（x）的值可從一般積分表所附的誤差函數表中直接查出［1］.

把式（7）代回式（5），即可得伽爾頓板實驗中，小球分布的概率密度函數為

這樣我們便使用最大熵原理推導出了伽爾頓板實驗中小球理論分布的具體函數表達式.

3 計算機模擬伽爾頓板實驗

計算機模擬實驗是科學研究的重要手段，它可以克服真實實驗中遇到的許多困難，彌補實驗儀器不足的缺陷［9］.使用計算機模擬伽爾頓板實驗可以方便地改變實驗參數，便于反復進行多次實驗，并快速得到實驗結果.

我們使用Matlab語言編寫了計算機模擬程序，對伽爾頓板實驗進行了蒙特卡羅模擬，模擬的伽爾頓板實驗裝置形狀如圖1所示，釘子的總行數和每行的釘子個數均可調整.設共有m行釘子，奇數行的釘子數相同為2n-1個，偶數行的釘子數相同為2n個.從上向下統計行數，最上邊的釘子為第一行，且第一行中間的那個釘子正對入口處，往下每行釘子交錯排開.以第一行中間的釘子為坐標原點，沿著第一行釘子為x軸，向右為x軸正方向，豎直向下為y軸的正方向，建立坐標系.規定同一行中相鄰兩個釘子的距離為1，相鄰的兩行距離也是1，奇數行兩端的釘子與板邊的距離為1，偶數行兩端的釘子與板邊的距離為0.5；規定第一行中間釘子的坐標為（0，1），從入口處落下的小球第一次都和坐標為（0，1）的釘子相碰，即小球落到第一行時的坐標都為（0，1），在隨后的下落過程中每個小球依次與下面的每行釘子中的一個釘子相碰，具體和哪個釘子相碰，由randn函數生成的一個正態分布隨機數決定.模擬中小球總數取為N，定義一個N行×2列的矩陣，用來存放這N個小球的位置；矩陣中的每一行代表一個小球，矩陣的第一列用來存放小球位置的x坐標，第二列用來存放小球位置的y坐標.

現以小球從第一行下落到第二行為例說明一下模擬過程.當生成的隨機數0≤randn＜n時，小球下落到第二行時位于它在第一行位置的右邊，具體下落到哪個位置是這樣規定的：當0≤randn＜1時，小球的橫坐標加0.5，縱坐標加1；當1≤randn＜2時，小球的橫坐標加1.5，縱坐標加1；……；當n-1≤randn＜n時，小球的橫坐標加n-0.5，縱坐標加1.當生成的隨機數-n≤randn＜0時，小球下落到第二行時位于它在第一行位置的左邊，具體規定如下：當-1≤randn＜0時，小球的橫坐標加-0.5，縱坐標加1；當-2≤randn＜-1時，小球的橫坐標加-1.5，縱坐標加1；……；當-n≤randn＜-n＋1時，小球的橫坐標加-n＋0.5，縱坐標加1.當生成的隨機數n≤randn≤3n時，小球下落到第二行的位置規定為：當0≤randn-n＜1時，小球的橫坐標為n-0.5，縱坐標加1；當1≤randn-n＜2時，小球的橫坐標為n-1.5，縱坐標加1；……；當2n-1≤randn-n≤2n時，小球的橫坐標為-n＋0.5，縱坐標加1.當生成的隨機數-3n≤randn＜-n時，小球下落到第二行的位置規定為：當-1≤randn＋n＜0時，小球的橫坐標為-n＋0.5，縱坐標加1；當-2≤randn＋n＜-1時，小球的橫坐標為-n＋1.5，縱坐標加1；……；當-2n≤randn＋n＜-2n＋1時，小球的橫坐標為n-0.5，縱坐標加1.當生成的隨機數randn＞3n或randn＜-3n時，拋棄該隨機數，令計算機再重新生成一個.小球從第二行下落到第三行等繼續下落的過程依次類推.

在下面的模擬中，我們采用了如下設置：共有24行釘子（m＝24），奇數行的釘子數為21個（n＝11），偶數行的釘子數為22個.我們取了1 000 000個小球（N＝1 000 000），對它們在伽爾頓板中的下落過程進行了模擬，得到了小球頻數按照落點位置分布的統計直方圖，見圖2.從圖中可看出，在正對小球入口處（中央位置）的小球數目最多，離中央位置越遠處的小球數目越少，分布形狀近似于正態分布，但與正態分布有差別之處，正態分布的范圍是從-∞到＋∞，而伽爾頓板實驗中小球的分布由于受到板左右側面的阻擋限制，分布范圍是從板的左邊緣到右邊緣.

圖2 小球落點頻數統計直方圖

4 理論推導結果與計算機模擬結果的比較

為了便于比較，我們把上述模擬得到的小球頻數除以小球總數，轉換為頻率，從而得到了小球頻率按照小球落點位置分布的數據，利用這些數據作出了小球頻率按照小球落點位置分布的條形圖，見圖3.

圖3 小球概率分布理論曲線與計算機模擬的小球頻率條形圖的比較

把頻率取自然對數，可得到小球頻率的自然對數與小球落點位置間關系的數據.然后使用最小二乘法，對小球頻率的自然對數與小球落點位置間關系的數據進行二次曲線擬合，擬合得到的二次曲線方程為

把我們導出的小球理論分布的概率密度函數式（9）兩邊取自然對數，得

式（10）與式（11）比較，得到

由上式可計算得

函數（14）為我們導出的小球理論分布概率密度函數的解析表達式，將此函數的關系曲線與小球頻率按照小球落點位置分布的條形圖作在同一張圖上，見圖3.由圖3可見，理論導出的小球分布關系曲線（細線）與模擬得到的小球分布頻率的條形圖符合得較好.這說明我們使用最大熵原理導出的伽爾頓板實驗中小球理論分布函數解析式較好地符合了實際情況.

5 結論

本文使用最大熵原理研究了伽爾頓板實驗小球的分布，推導出了小球落點分布的概率密度函數解析表達式，見式（9）.接著使用蒙特卡羅方法對伽爾頓板實驗進行了計算機模擬.通過把導出的小球理論分布概率密度函數解析式作成曲線，并把此函數曲線與小球頻率按落點位置分布的條形圖作在一起進行比較，發現二者符合得較好，這說明導出的小球分布函數解析式較為成功.

［1］李椿.熱學 第二版［M］.北京：高等教育出版社，2008

［2］郝志峰，謝國瑞，汪國強.概率論與數理統計 第二版［M］.北京：高等教育出版社，2009

［3］彭芳麟.伽爾頓板實驗的計算機模擬［J］.大學物理，2005，24（1）：45～49

［4］廖旭，任學藻.用二項式分布研究伽爾頓板實驗的分布曲線［J］.實驗科學與技術，2006，（1）：79～81

［5］聶燕.高爾頓釘板試驗的算法實現及分析［J］.中國民航飛行學院學報，2008，19（3）：62～64

［6］Banavar J R，Maritan A，Volkov I.Applications of the principle of maximum entropy：from physics to ecology［J］.JournalofPhysics：Condensed Matter，2010，22：063101

［7］Plastino A，Curado E M F.Equivalence between maximum entropy principle and enforcing dU＝TdS［J］.PhysicalReviewE，2005，72：047103

［8］Jin H Y，Luo L F，Zhang L R.Using estimative reaction free energy to predict splice sites and their flanking competitors［J］.Gene，2008，424（1-2）：115～120

［9］彭芳麟.計算物理基礎［M］.北京：高等教育出版社，2010