簡(jiǎn)介

與其它開放源碼語(yǔ)言（比如 Perl 和 Python）相比，PHP 社區(qū)缺少?gòu)?qiáng)有力的工作來(lái)開發(fā)數(shù)學(xué)庫(kù)。

造成這種狀況的一個(gè)原因可能是由于已經(jīng)存在大量成熟的數(shù)學(xué)工具，這可能阻礙了社區(qū)自行開發(fā) PHP 工具的工作。例如，我曾研究過(guò)一個(gè)功能強(qiáng)大的工具 S System，它擁有一組令人印象深刻的統(tǒng)計(jì)庫(kù)，專門被設(shè)計(jì)成用來(lái)分析數(shù)據(jù)集，并且在 1998 年由于其語(yǔ)言設(shè)計(jì)而獲得了 ACM 獎(jiǎng)。如果 S 或者其開放源碼同類 R 僅僅是一個(gè) exec_shell 調(diào)用，那么為何還要麻煩用 PHP 實(shí)現(xiàn)相同的統(tǒng)計(jì)計(jì)算功能呢？有關(guān) S System、它的 ACM 獎(jiǎng)或 R 的更多信息，請(qǐng)參閱相關(guān)參考資料。

難道這不是在浪費(fèi)開發(fā)人員的精力嗎？如果開發(fā) PHP 數(shù)學(xué)庫(kù)的動(dòng)機(jī)是出自節(jié)省開發(fā)人員的精力以及使用最好的工具來(lái)完成工作，那么 PHP 現(xiàn)在的課題是很有意義的。

另一方面，出于教學(xué)動(dòng)機(jī)可能會(huì)鼓勵(lì)對(duì) PHP 數(shù)學(xué)庫(kù)的開發(fā)。對(duì)于大約 10% 的人來(lái)說(shuō)，數(shù)學(xué)是個(gè)值得探索的有趣課題。對(duì)于那些同時(shí)還熟練應(yīng)用 PHP 的人來(lái)說(shuō)，PHP 數(shù)學(xué)庫(kù)的開發(fā)可以增強(qiáng)數(shù)學(xué)學(xué)習(xí)過(guò)程，換句話說(shuō)，不要只 閱讀有關(guān) T 測(cè)試的章節(jié)，還要 實(shí)現(xiàn)一個(gè)能計(jì)算相應(yīng)的中間值并用標(biāo)準(zhǔn)格式顯示它們的類。

通過(guò)指導(dǎo)和訓(xùn)練，我希望證明開發(fā) PHP 數(shù)學(xué)庫(kù)并不是一項(xiàng)很難的任務(wù)，它可能代表一項(xiàng)有趣的技術(shù)和學(xué)習(xí)難題。在本文中，我將提供一個(gè) PHP 數(shù)學(xué)庫(kù)示例，名為 SimpleLinearRegression ，它演示了一個(gè)可以用來(lái)開發(fā) PHP 數(shù)學(xué)庫(kù)的通用方法。讓我們從討論一些通用的原則開始，這些原則指導(dǎo)我開發(fā)這個(gè) SimpleLinearRegression 類。

指導(dǎo)原則

我使用了六個(gè)通用原則來(lái)指導(dǎo) SimpleLinearRegression 類的開發(fā)。

每個(gè)分析模型建立一個(gè)類。 使用逆向鏈接來(lái)開發(fā)類。 預(yù)計(jì)有大量的 getter。 存儲(chǔ)中間結(jié)果。 為詳細(xì)的 API 制定首選項(xiàng)。 盡善盡美并非目標(biāo)。

讓我們更詳細(xì)地逐條研究這些指導(dǎo)方針。　每個(gè)分析模型建立一個(gè)類

每種主要的分析測(cè)試或過(guò)程應(yīng)當(dāng)有一個(gè)名稱與測(cè)試或過(guò)程名相同的 PHP 類，這個(gè)類包含了輸入函數(shù)、計(jì)算中間值和匯總值的函數(shù)和輸出函數(shù)（將中間值和匯總值用文本或圖形格式全部顯示在屏幕上）。

使用逆向鏈接來(lái)開發(fā)類

在數(shù)學(xué)編程中，編碼的目標(biāo)通常是分析過(guò)程（比如 MultipleRegression 、 TimeSeries 或 ChiSquared ）所希望生成的標(biāo)準(zhǔn)輸出值。從解決問(wèn)題的角度出發(fā)，這意味著您可以使用逆向鏈接來(lái)開發(fā)數(shù)學(xué)類的方法。

例如，匯總輸出屏幕顯示了一個(gè)或多個(gè)匯總統(tǒng)計(jì)結(jié)果。這些匯總統(tǒng)計(jì)結(jié)果依賴于中間統(tǒng)計(jì)結(jié)果的計(jì)算，這些中間統(tǒng)計(jì)結(jié)果又可能會(huì)涉及到更深一層的中間統(tǒng)計(jì)結(jié)果，以此類推。這個(gè)基于逆向鏈接的開發(fā)方法導(dǎo)出了下一個(gè)原則。

預(yù)計(jì)有大量的 getter

數(shù)學(xué)類的大部分類開發(fā)工作都涉及到計(jì)算中間值和匯總值。實(shí)際上，這意味著，如果您的類包含許多計(jì)算中間值和匯總值的 getter 方法，您不應(yīng)當(dāng)感到驚訝。

存儲(chǔ)中間結(jié)果

將中間計(jì)算結(jié)果存儲(chǔ)在結(jié)果對(duì)象內(nèi)，這樣您就可以將中間結(jié)果用作后續(xù)計(jì)算的輸入。在 S 語(yǔ)言設(shè)計(jì)中實(shí)施了這一原則。在當(dāng)前環(huán)境下，通過(guò)選擇實(shí)例變量來(lái)表示計(jì)算得到的中間值和匯總結(jié)果，從而實(shí)施了該原則。

為詳細(xì)的 API 制定首選項(xiàng)

當(dāng)為 SimpleLinearRegression 類中的成員函數(shù)和實(shí)例變量制定命名方案時(shí)，我發(fā)現(xiàn)：如果我使用較長(zhǎng)的名稱（類似于 getSumSquaredError 這樣的名稱，而不是 getYY2 ）來(lái)描述成員函數(shù)和實(shí)例變量，那么就更容易了解函數(shù)的操作內(nèi)容和變量所代表的意義。

我沒(méi)有完全放棄簡(jiǎn)寫名稱；但是，當(dāng)我用簡(jiǎn)寫形式的名稱時(shí)，我得設(shè)法提供注釋以完整闡述該名稱的含義。我的看法是：高度簡(jiǎn)寫的命名方案在數(shù)學(xué)編程中很常見，但它們使得理解和證明某個(gè)數(shù)學(xué)例程是否按部就班更為困難，而原本不必造成此種困難。

盡善盡美并非目標(biāo)

這個(gè)編碼練習(xí)的目標(biāo)不是一定要為 PHP 開發(fā)高度優(yōu)化和嚴(yán)格的數(shù)學(xué)引擎。在早期階段，應(yīng)當(dāng)強(qiáng)調(diào)學(xué)習(xí)實(shí)現(xiàn)意義重大的分析測(cè)試，以及解決這方面的難題。

實(shí)例變量

當(dāng)對(duì)統(tǒng)計(jì)測(cè)試或過(guò)程進(jìn)行建模時(shí)，您需要指出聲明哪些實(shí)例變量。

實(shí)例變量的選擇可以通過(guò)說(shuō)明由分析過(guò)程生成的中間值和匯總值來(lái)確定。每個(gè)中間值和匯總值都可以有一個(gè)相應(yīng)的實(shí)例變量，將變量的值作為對(duì)象屬性。

我采用這樣的分析來(lái)確定為清單 1 中的 SimpleLinearRegression 類聲明哪些變量。可以對(duì) MultipleRegression 、 ANOVA 或 TimeSeries 過(guò)程執(zhí)行類似的分析。

清單 1. SimpleLinearRegression 類的實(shí)例變量<?php; // Copyright 2003, Paul Meagher; // Distributed under GPL;class SimpleLinearRegression {; var $n var $X = array(); var $Y = array(); var $ConfInt; var $Alpha; var $XMean; var $YMean; var $SumXX; var $SumXY; var $SumYY; var $Slope; var $YInt; var $PredictedY= array(); var $Error;;;;;= array(); var $SquaredError = array(); var $TotalError; var $SumError; var $SumSquaredError; var $ErrorVariance; var $StdErr; var $SlopeStdErr; var $SlopeVal;// T value of Slope; var $YIntStdErr;;; var $YIntTVal;// T value for Y Intercept var $R; var $RSquared;;; var $DF; // Degrees of Freedom var $SlopeProb// Probability of Slope Estimate var $YIntProb;// Probability of Y Intercept Estimate var $AlphaTVal// T Value for given alpha setting var $ConfIntOfSlope; var $RPath; = '/usr/local/bin/R'// Your path here var $format = '%01.2f'; // Used for formatting output; } ?> 構(gòu)造函數(shù)

SimpleLinearRegression 類的構(gòu)造函數(shù)方法接受一個(gè) X和一個(gè) Y向量，每個(gè)向量都有相同數(shù)量的值。您還可以為您預(yù)計(jì)的 Y值設(shè)置一個(gè)缺省為 95% 的置信區(qū)間（confidence interval）。

構(gòu)造函數(shù)方法從驗(yàn)證數(shù)據(jù)形式是否適合于處理開始。一旦輸入向量通過(guò)了“大小相等”和“值大于 1”測(cè)試，就執(zhí)行算法的核心部分。

執(zhí)行這項(xiàng)任務(wù)涉及到通過(guò)一系列 getter 方法計(jì)算統(tǒng)計(jì)過(guò)程的中間值和匯總值。將每個(gè)方法調(diào)用的返回值賦給該類的一個(gè)實(shí)例變量。用這種方法存儲(chǔ)計(jì)算結(jié)果確保了前后鏈接的計(jì)算中的調(diào)用例程可以使用中間值和匯總值。還可以通過(guò)調(diào)用該類的輸出方法來(lái)顯示這些結(jié)果，如清單 2 所描述的那樣。

清單 2. 調(diào)用類輸出方法

<?php; // Copyright 2003, Paul Meagher; // Distributed under GPL;function SimpleLinearRegression($X, $Y, $ConfidenceInterval='95') { $numX = count($X); $numY = count($Y); if ($numX != $numY) { die('Error: Size of X and Y vectors must be the same.'); }; if ($numX <= 1) {; die('Error: Size of input array must be at least 2.'); } $this->n;;;;= $numX; $this->X;;;;= $X; $this->Y;;;;= $Y; $this->ConfInt = $ConfidenceInterval $this->Alpha= (1 + ($this->ConfInt / 100) ) / 2; $this->XMean= $this->getMean($this->X); $this->YMean= $this->getMean($this->Y); $this->SumXX= $this->getSumXX(); $this->SumYY= $this->getSumYY(); $this->SumXY= $this->getSumXY();;; $this->Slope= $this->getSlope(); $this->YInt;= $this->getYInt(); $this->PredictedY;;;= $this->getPredictedY(); $this->Error= $this->getError(); $this->SquaredError;= $this->getSquaredError(); $this->SumError;;;;;= $this->getSumError(); $this->TotalError;;;= $this->getTotalError();;; $this->SumSquaredError = $this->getSumSquaredError(); $this->ErrorVariance= $this->getErrorVariance(); $this->StdErr; = $this->getStdErr(); $this->SlopeStdErr;;= $this->getSlopeStdErr();;;; $this->YIntStdErr;;;= $this->getYIntStdErr() $this->SlopeTVal;;;;= $this->getSlopeTVal();;; $this->YIntTVal;;;;;= $this->getYIntTVal(); $this->R;;;;= $this->getR();; $this->RSquared;;;;;= $this->getRSquared(); $this->DF;;;= $this->getDF(); $this->SlopeProb;;;;= $this->getStudentProb($this->SlopeTVal, $this->DF); $this->YIntProb;;;;;= $this->getStudentProb($this->YIntTVal, $this->DF); $this->AlphaTVal;;;;= $this->getInverseStudentProb($this->Alpha, $this->DF); $this->ConfIntOfSlope; = $this->getConfIntOfSlope() return true; } ?>

方法名及其序列是通過(guò)結(jié)合逆向鏈接和參考大學(xué)本科學(xué)生使用的統(tǒng)計(jì)學(xué)教科書推導(dǎo)得出的，該教科書一步一步地說(shuō)明了如何計(jì)算中間值。我需要計(jì)算的中間值的名稱帶有“get”前綴，從而推導(dǎo)出方法名。

使模型與數(shù)據(jù)相吻合

SimpleLinearRegression 過(guò)程用于產(chǎn)生與數(shù)據(jù)相吻合的直線，其中直線具有以下標(biāo)準(zhǔn)方程：

y = b + mx

該方程的 PHP 格式看起來(lái)類似于清單 3：

清單 3. 使模型與數(shù)據(jù)相吻合的 PHP 方程

$PredictedY[$i] = $YIntercept + $Slope * $X[$i]

SimpleLinearRegression 類使用最小二乘法準(zhǔn)則推導(dǎo)出 Y 軸截距（Y Intercept）和斜率（Slope）參數(shù)的估計(jì)值。這些估計(jì)的參數(shù)用來(lái)構(gòu)造線性方程（請(qǐng)參閱 清單 3），該方程對(duì) X和 Y值之間的關(guān)系進(jìn)行建模。

使用推導(dǎo)出的線性方程，您就可以得到每個(gè) X值對(duì)應(yīng)的預(yù)測(cè) Y值。如果線性方程與數(shù)據(jù)非常吻合，那么 Y的觀測(cè)值與預(yù)測(cè)值趨近于一致。

如何確定是否非常吻合

SimpleLinearRegression 類生成了相當(dāng)多的匯總值。一個(gè)重要的匯總值是 T統(tǒng)計(jì)值，它可以用來(lái)衡量一個(gè)線性方程與數(shù)據(jù)的 吻合程度。如果非常吻合，那么 T 統(tǒng)計(jì)值往往很大。如果 T 統(tǒng)計(jì)值很小，那么應(yīng)當(dāng)用一個(gè)模型替換該線性方程，該模型假設(shè) Y值的均值是最佳預(yù)測(cè)值（也就是說(shuō)，一組值的均值通常是下一個(gè)觀測(cè)值有用的預(yù)測(cè)值，使之成為缺省模型）。

要測(cè)試 T 統(tǒng)計(jì)值是否大得足以不把 Y值的均值作為最佳預(yù)測(cè)值，您需要計(jì)算獲取 T 統(tǒng)計(jì)值的隨機(jī)概率。如果獲取 T 統(tǒng)計(jì)值的概率很低，那么您可以否定均值是最佳預(yù)測(cè)值這個(gè)無(wú)效假設(shè)，與此相對(duì)應(yīng)，也就確信簡(jiǎn)單線性模型與數(shù)據(jù)非常吻合。

那么，如何計(jì)算 T 統(tǒng)計(jì)值的概率呢？

計(jì)算 T 統(tǒng)計(jì)值概率

由于 PHP 缺少計(jì)算 T 統(tǒng)計(jì)值概率的數(shù)學(xué)例程，因此我決定將此任務(wù)交給統(tǒng)計(jì)計(jì)算包 R（請(qǐng)參閱 參考資料中的 www.r-project.org）來(lái)獲得必要的值。我還想提醒大家注意該包，因?yàn)椋?

R 提供了許多想法，PHP 開發(fā)人員可能會(huì)在 PHP 數(shù)學(xué)庫(kù)中模擬這些想法。 有了 R，可以確定從 PHP 數(shù)學(xué)庫(kù)獲得的值與那些從成熟的免費(fèi)可用的開放源碼統(tǒng)計(jì)包中獲得的值是否一致。 清單 4 中的代碼演示了交給 R 來(lái)處理以獲取一個(gè)值是多么容易。

清單 4. 交給 R 統(tǒng)計(jì)計(jì)算包來(lái)處理以獲取一個(gè)值

<?php; // Copyright 2003, Paul Meagher; // Distributed under GPL;class SimpleLinearRegression {; var $RPath; = '/usr/local/bin/R'// Your path here function getStudentProb($T, $df) {;; $Probability = 0.0;; $cmd = 'echo 'dt($T, $df)' | $this->RPath --slave' $result = shell_exec($cmd);;; list($LineNumber, $Probability) = explode(' ', trim($result)) return $Probability; } function getInverseStudentProb($alpha, $df) { $InverseProbability = 0.0 $cmd = 'echo 'qt($alpha, $df)' | $this->RPath --slave' $result = shell_exec($cmd); list($LineNumber, $InverseProbability) = explode(' ', trim($result)) return $InverseProbability; } } ?>

請(qǐng)注意，這里已經(jīng)設(shè)置了到 R 可執(zhí)行文件的路徑，并在兩個(gè)函數(shù)中使用了該路徑。第一個(gè)函數(shù)根據(jù)學(xué)生的 T 分布返回了與 T 統(tǒng)計(jì)值相關(guān)的概率值，而第二個(gè)反函數(shù)計(jì)算了與給定的 alpha 設(shè)置相對(duì)應(yīng)的 T 統(tǒng)計(jì)值。 getStudentProb 方法用來(lái)評(píng)估線性模型的吻合程度； getInverseStudentProb 方法返回一個(gè)中間值，它用來(lái)計(jì)算每個(gè)預(yù)測(cè)的 Y值的置信區(qū)間。

由于篇幅有限，我不可能逐個(gè)詳細(xì)說(shuō)明這個(gè)類中的所有函數(shù)，因此如果您想搞清楚簡(jiǎn)單線性回歸分析中所涉及的術(shù)語(yǔ)和步驟，我鼓勵(lì)您參考大學(xué)本科學(xué)生使用的統(tǒng)計(jì)學(xué)教科書。

燃耗研究

要演示如何使用該類，我可以使用來(lái)自公共事業(yè)中燃耗（burnout）研究中的數(shù)據(jù)。Michael Leiter 和 Kimberly Ann Meechan 研究了稱為 消耗指數(shù)（Exhaustion Index)的燃耗度量單位和稱之為 集中度（Concentration）的獨(dú)立變量之間的關(guān)系。集中度是指人們的社交接觸中來(lái)自其工作環(huán)境的那部分比例。

要研究他們樣本中個(gè)人的消耗指數(shù)值與集中度值之間的關(guān)系，請(qǐng)將這些值裝入適當(dāng)命名的數(shù)組中，并用這些數(shù)組值對(duì)該類進(jìn)行實(shí)例化。對(duì)類進(jìn)行實(shí)例化后，顯示該類所生成的某些匯總值以評(píng)估線性模型與數(shù)據(jù)的吻合程度。

清單 5 顯示了裝入數(shù)據(jù)和顯示匯總值的腳本：

清單 5. 用于裝入數(shù)據(jù)并顯示匯總值的腳本

<?php; // BurnoutStudy.php // Copyright 2003, Paul Meagher; // Distributed under GPL;include 'SimpleLinearRegression.php'// Load data from burnout study; $Concentration= array(20,60,38,88,79,87,; 68,12,35,70,80,92,; 77,86,83,79,75,81,; 75,77,77,77,17,85,96); $ExhaustionIndex = array(100,525,300,980,310,900,; 410,296,120,501,920,810,; 506,493,892,527,600,855,; 709,791,718,684,141,400,970); $slr = new SimpleLinearRegression($Concentration, $ExhaustionIndex);$YInt;;;= sprintf($slr->format, $slr->YInt); $Slope;;= sprintf($slr->format, $slr->Slope);;;$SlopeTVal = sprintf($slr->format, $slr->SlopeTVal);;;$SlopeProb = sprintf('%01.6f', $slr->SlopeProb);;;?> <table border='1' cellpadding='5'> <tr> <th align='right'>Equation:</th> <td></td> </tr> <tr> <th align='right'>T:</th> <td></td> </tr> <tr> <th align='right'>Prob > T:</th> <td><td> </tr> </table>

通過(guò) Web 瀏覽器運(yùn)行該腳本，產(chǎn)生以下輸出：

Equation: Exhaustion = -29.50 + (8.87 * Concentration); T:; 6.03 Prob > T:; 0.000005

這張表的最后一行指出獲取這樣大 T值的 隨機(jī)概率非常低。可以得出這樣的結(jié)論：與僅僅使用消耗值的均值相比，簡(jiǎn)單線性模型的預(yù)測(cè)能力更好。

知道了某個(gè)人的工作場(chǎng)所聯(lián)系的集中度，就可以用來(lái)預(yù)測(cè)他們可能正在消耗的燃耗程度。這個(gè)方程告訴我們：集中度值每增加 1 個(gè)單位，社會(huì)服務(wù)領(lǐng)域中一個(gè)人的消耗值就會(huì)增加 8 個(gè)單位。這進(jìn)一步證明了：要減少潛在的燃耗，社會(huì)服務(wù)領(lǐng)域中的個(gè)人應(yīng)當(dāng)考慮在其工作場(chǎng)所之外結(jié)交朋友。

這只是粗略地描述了這些結(jié)果可能表示的含義。為全面研究這個(gè)數(shù)據(jù)集的含義，您可能想更詳細(xì)地研究這個(gè)數(shù)據(jù)以確信這是正確的解釋。在下一篇文章中我將討論應(yīng)當(dāng)執(zhí)行其它哪些分析。

您學(xué)到了什么？

其一，要開發(fā)意義重大的基于 PHP 的數(shù)學(xué)包，您不必是一名火箭科學(xué)家。堅(jiān)持標(biāo)準(zhǔn)的面向?qū)ο蠹夹g(shù)，以及明確地采用逆向鏈接問(wèn)題解決方法，就可以相對(duì)方便地使用 PHP 實(shí)現(xiàn)某些較為基本的統(tǒng)計(jì)過(guò)程。

從教學(xué)的觀點(diǎn)出發(fā)，我認(rèn)為：如果只是因?yàn)橐竽谳^高和較低的抽象層次思考統(tǒng)計(jì)測(cè)試或例程，那么這個(gè)練習(xí)是非常有用的。換句話說(shuō)，補(bǔ)充您的統(tǒng)計(jì)測(cè)試或過(guò)程學(xué)習(xí)的一個(gè)好辦法就是將這個(gè)過(guò)程作為算法實(shí)現(xiàn)。

要實(shí)現(xiàn)統(tǒng)計(jì)測(cè)試通常需要超出所給定的信息范圍并創(chuàng)造性地解決和發(fā)現(xiàn)問(wèn)題。對(duì)于發(fā)現(xiàn)對(duì)某個(gè)學(xué)科認(rèn)識(shí)的不足而言，它也是一個(gè)好辦法。

不利的一面，您發(fā)現(xiàn) PHP 對(duì)于取樣分布缺乏內(nèi)在手段，而這是實(shí)現(xiàn)大多數(shù)統(tǒng)計(jì)測(cè)試所必需的。您需要交給 R 來(lái)處理以獲取這些值，但是我擔(dān)心您會(huì)沒(méi)時(shí)間或沒(méi)興趣安裝 R。某些常見概率函數(shù)的本機(jī) PHP 實(shí)現(xiàn)可以解決這個(gè)問(wèn)題。

另一個(gè)問(wèn)題是：該類生成許多中間值和匯總值，但是匯總輸出實(shí)際上沒(méi)有利用這一點(diǎn)。我提供了一些難處理的輸出，但是這既不夠充分也沒(méi)進(jìn)行很好的組織，以致您無(wú)法充分地解釋分析結(jié)果。實(shí)際上，我完全不知道如何可以將輸出方法集成到該類中。這需要得到解決。

最后，要弄明白數(shù)據(jù)，不僅僅是察看匯總值就可以了。您還需要明白各個(gè)數(shù)據(jù)點(diǎn)是如何分布的。最好的辦法之一是將您的數(shù)據(jù)繪制成圖表。再次聲明，我對(duì)這方面不太了解，但是如果要用這個(gè)類來(lái)分析實(shí)際數(shù)據(jù)的話就需要解決這個(gè)問(wèn)題。