用 PHP 实现的简单线性回归

本文的目标是向您展示如何使用第 1 部分中讨论的 SimpleLinearRegression 类来构建一个重要的数据研究工具。

简要回顾：概念

简单线性回归建模背后的基本目标是从成对的 X值和 Y值（即 X和 Y测量值）组成的二维平面中找到最吻合的直线。一旦用 最小方差法找到这条直线，就可以执行各种统计测试，以确定这条直线与观测到的 Y值的偏离量吻合程度。

线性方程（ y = mx + b）有两个参数必须根据所提供的 X和 Y数据估算出来，它们是斜率（ m）和 y 轴截距（ b）。一旦估算出这两个参数，就可以将观测值输入线性方程，并观察方程所生成的 Y预测值。

要使用最小方差法估算出 m和 b参数，就要找到 m 和 b 的估计值，使它们对于所有的 X值得到的 Y值的观测值和预测值最校观测值和预测值之差称为误差（ y i- (mx i+ b) ），并且，如果对每个误差值都求平方，然后求这些残差的和，其结果是一个被称为 预测平方差的数。使用最小方差法来确定最吻合的直线涉及寻找使预测方差最小的 m和 b的估计值。

可以用两种基本方法来找到满足最小方差法的估计值 m和 b。第一种方法，可以使用数值搜索过程设定不同的 m和 b值并对它们求值，最终决定产生最小方差的估计值。第二种方法是使用微积分找到用于估算 m和 b 的方程。我不打算深入讨论推导出这些方程所涉及的微积分，但我确实在 SimpleLinearRegression 类中使用了这些分析方程，以找到 m和 b 的最小平方估计值（请参阅 SimpleLinearRegression 类中的 getSlope() 和 getYIntercept 方法）。

即使拥有了可以用来找到 m和 b的最小平方估计值的方程，也并不意味着只要将这些参数代入线性方程，其结果就是一条与数据良好吻合的直线。这个简单线性回归过程中的下一步是确定其余的预测方差是否可以接受。

可以使用统计决策过程来否决“直线与数据吻合”这个备择假设。这个过程基于对 T 统计值的计算，使用概率函数求得随机大的观测值的概率。正如第 1 部分所提到的， SimpleLinearRegression 类生成了为数众多的汇总值，其中一个重要的汇总值是 T 统计值，它可以用来衡量线性方程与数据的吻合程度。如果吻合良好，则 T 统计值往往是一个较大的值；如果 T 值很小，就应该用一个缺省模型代替您的线性方程，该模型假定 Y值的平均值是最佳预测值（因为一组值的平均值通常可以是下一个观测值的有用的预测值）。

要测试 T 统计值是否大到可以不用 Y值的平均值作为最佳预测值，需要计算随机获得 T 统计值的概率。如果概率很低，那就可以不采用平均值是最佳预测值这一无效假设，并且相应地可以确信简单线性模型是与数据良好吻合的。（有关计算 T 统计值概率的更多信息，请参阅第 1 部分。）

回过头讨论统计决策过程。它告诉您何时不采用无效假设，却没有告诉您是否接受备择假设。在研究环境中，需要通过理论参数和统计参数来建立线性模型备择假设。

您将构建的数据研究工具实现了用于线性模型（T 测试）的统计决策过程，并提供了可以用来构造理论和统计参数的汇总数据，这些参数是建立线性模型所需要的。数据研究工具可以归类为决策支持工具，供知识工作者在中小规模的数据集中研究模式。

从学习的角度来看，简单线性回归建模值得研究，因为它是理解更高级形式的统计建模的必由之路。例如，简单线性回归中的许多核心概念为理解多次回归（Multiple Regression）、要素分析（Factor Analysis）和时间序列（Time Series）等建立了良好的基矗

简单线性回归还是一种多用途的建模技术。通过转换原始数据（通常用对数或幂转换），可以用它来为曲线数据建模。这些转换可以使数据线性化，这样就可以使用简单线性回归来为数据建模。所生成的线性模型将被表示为与被转换值相关的线性公式。

回页首

概率函数

在前一篇文章中，我通过交由 R 来求得概率值，从而避开了用 PHP 实现概率函数的问题。我对这个解决方案并非完全满意，因此我开始研究这个问题：开发基于 PHP 的概率函数需要些什么。

我开始上网查找信息和代码。一个两者兼有的来源是书籍 Numerical Recipes in C 中的概率函数。我用 PHP 重新实现了一些概率函数代码（ gammln.c 和 betai.c 函数），但我对结果还是不满意。与其它一些实现相比，其代码似乎多了些。此外，我还需要反概率函数。

幸运的是，我偶然发现了 John Pezzullo 的 Interactive Statistical Calculation。John 关于 概率分布函数的网站上有我需要的所有函数，为便于学习，这些函数已用 JavaScript 实现。

我将 Student T 和 Fisher F 函数移植到了 PHP。我对 API 作了一点改动，以便符合 Java 命名风格，并将所有函数嵌入到名为 Distribution 的类中。该实现的一个很棒的功能是 doCommonMath 方法，这个库中的所有函数都重用了它。我没有花费力气去实现的其它测试（正态测试和卡方测试）也都使用 doCommonMath 方法。

这次移植的另一个方面也值得注意。通过使用 JavaScript，用户可以将动态确定的值赋给实例变量，譬如：

var PiD2 = pi() / 2

在 PHP 中不能这样做。只能把简单的常量值赋给实例变量。希望在 PHP5 中会解决这个缺陷。

请注意 清单 1中的代码并未定义实例变量 — 这是因为在 JavaScript 版本中，它们是动态赋予的值。

清单 1. 实现概率函数

<?php // Distribution.php // Copyright John Pezullo // Released under same terms as PHP. // PHP Port and OO'fying by Paul Meagher class Distribution { function doCommonMath($q, $i, $j, $b) { $zz = 1; $z = $zz; $k = $i; while($k <= $j) { $zz = $zz * $q * $k / ($k - $b); $z = $z + $zz; $k = $k + 2; } return $z; } function getStudentT($t, $df) { $t = abs($t); $w = $t / sqrt($df); $th = atan($w); if ($df == 1) { return 1 - $th / (pi() / 2); } $sth = sin($th); $cth = cos($th); if( ($df % 2) ==1 ) { return 1 - ($th + $sth * $cth * $this->doCommonMath($cth * $cth, 2, $df - 3, -1)) / (pi()/2); } else { return 1 - $sth * $this->doCommonMath($cth * $cth, 1, $df - 3, -1); } } function getInverseStudentT($p, $df) { $v = 0.5; $dv = 0.5; $t = 0; while($dv > 1e-6) { $t = (1 / $v) - 1; $dv = $dv / 2; if ( $this->getStudentT($t, $df) > $p) { $v = $v - $dv; } else { $v = $v + $dv; } } return $t; } function getFisherF($f, $n1, $n2) { // implemented but not shown } function getInverseFisherF($p, $n1, $n2) { // implemented but not shown } } ?>

回页首

图形输出

迄今为止，您已经实现的输出方法都以 HTML 格式显示汇总值。它也适合于用 GIF、JPEG 或 PNG 格式显示这些数据的分布图（scatter plot）或线图（line plot）。

与其亲自编写生成线图和分布图的代码，我认为最好使用名为 JpGraph的基于 PHP 的图形库。JpGraph 正由 Johan Persson 积极开发，其 项目网站这样描述它：

无论是对于只有最少代码的“以快捷但不恰当方式获得的”图形，还是对于需要非常细粒度控制的复杂专业图形，JpGraph 都可以使它们的绘制变得简单。JpGraph 同样适用于科学和商业类型的图形。

JpGraph 分发版中包含大量可以根据特定需求进行定制的示例脚本。将 JpGraph 用于数据研究工具非常简单，只需找到功能与我的需求类似的示例脚本，然后对该脚本进行改写以满足我的特定需求即可。

清单 3中的脚本是从样本数据研究工具（ explore.php）中抽取的，它演示了如何调用该库以及如何将来自于 SimpleLinearRegression 分析的数据填入 Line 和 Scatter 类。这段代码中的注释是 Johan Persson 编写的（JPGraph 代码库的文档化工作做得很好）。

清单 3. 来自于样本数据研究工具 explore.php 的函数的详细内容

<?php // Snippet extracted from explore.php script include ("jpgraph/jpgraph.php"); include ("jpgraph/jpgraph_scatter.php"); include ("jpgraph/jpgraph_line.php"); // Create the graph $graph = new Graph(300,200,'auto'); $graph->SetScale("linlin"); // Setup title $graph->title->Set("$title"); $graph->img->SetMargin(50,20,20,40); $graph->xaxis->SetTitle("$x_name","center"); $graph->yaxis->SetTitleMargin(30); $graph->yaxis->title->Set("$y_name"); $graph->title->SetFont(FF_FONT1,FS_BOLD); // make sure that the X-axis is always at the // bottom at the plot and not just at Y=0 which is // the default position $graph->xaxis->SetPos('min'); // Create the scatter plot with some nice colors $sp1 = new ScatterPlot($slr->Y, $slr->X); $sp1->mark->SetType(MARK_FILLEDCIRCLE); $sp1->mark->SetFillColor("red"); $sp1->SetColor("blue"); $sp1->SetWeight(3); $sp1->mark->SetWidth(4); // Create the regression line $lplot = new LinePlot($slr->PredictedY, $slr->X); $lplot->SetWeight(2); $lplot->SetColor('navy'); // Add the pltos to the line $graph->Add($sp1); $graph->Add($lplot); // ... and stroke $graph_name = "temp/test.png"; $graph->Stroke($graph_name); ?> <img src='<?php echo $graph_name ?>' vspace='15'> ?>

回页首

数据研究脚本

该数据研究工具由单个脚本（ explore.php）构成，该脚本调用 SimpleLinearRegressionHTML 类和 JpGraph 库的方法。

该脚本使用了简单的处理逻辑。该脚本的第一部分对所提交的表单数据执行基本验证。如果这些表单数据通过验证，则执行该脚本的第二部分。

该脚本的第二部分所包含的代码用于分析数据，并以 HTML 和图形格式显示汇总结果。 清单 4中显示了 explore.php脚本的基本结构：

清单 4. explore.php 的结构

<?php // explore.php if (!empty($x_values)) { $X = explode(",", $x_values); $numX = count($X); } if (!empty($y_values)) { $Y = explode(",", $y_values); $numY = count($Y); } // display entry data entry form if variables not set if ( (empty($title)) OR (empty($x_name)) OR (empty($x_values)) OR (empty($y_name)) OR (empty($conf_int)) OR (empty($y_values)) OR ($numX != $numY) ) { // Omitted code for displaying entry form } else { include_once "slr/SimpleLinearRegressionHTML.php"; $slr = new SimpleLinearRegressionHTML($X, $Y, $conf_int); echo "<h2>$title</h2>"; $slr->showTableSummary($x_name, $y_name); echo "<br><br>"; $slr->showAnalysisOfVariance(); echo "<br><br>"; $slr->showParameterEstimates($x_name, $y_name); echo "<br>"; $slr->showFormula($x_name, $y_name); echo "<br><br>"; $slr->showRValues($x_name, $y_name); echo "<br>"; include ("jpgraph/jpgraph.php"); include ("jpgraph/jpgraph_scatter.php"); include ("jpgraph/jpgraph_line.php"); // The code for displaying the graphics is inline in the // explore.php script. The code for these two line plots // finishes off the script: // Omitted code for displaying scatter plus line plot // Omitted code for displaying residuals plot } ?>