除了梯度下降法去擬合曲線,最小二乘法也是另一個方法娄琉。和梯度下降不斷迭代求出極值不同的是曲横,最小二乘法是直接求導計算出參數(shù)。數(shù)學的問題還是參考數(shù)學相關的資料吧踩蔚,同時可以@refer http://blog.csdn.net/jairuschan/article/details/7517773
下面直接上代碼。
/*
* 利用最小二乘法求擬合參數(shù)
*
* @auther menmei
* @date 2017/03/27
*/
/*
* Least Square Procedure
* @fomula w = (x'x)-1x'y
*/
/*
* 獲取矩陣的轉(zhuǎn)秩
*
* @param Array
* @return Array
*
*/
function getMatrixT($origin){
$lines = count($origin);
if($lines <= 0) return False;
$rowsArr = $origin[0];
if(!is_array($rowsArr)) return False;
$rows = count($origin[0]);
$return = [];
for($i = 0; $i < $lines; $i ++){
for($j = 0; $j < $rows; $j ++){
$return[$j][$i] = $origin[$i][$j];
}
}
return $return;
}
/*
* 矩陣相乘
*
* @param Array $a
* @param Array $b
*
* @condition rows(A) = lines(B)
* @return Array $return
*/
function matrixMulti($a, $b){
$linesA = count($a);
$rowsA = count($a[0]);
$linesB = count($b);
$rowsB = count($b[0]);
if($rowsA != $linesB) return False;
//echo $rows;
//if($linesA == 0 || $rowsB == 0) return False;
for($i = 0; $i < $linesA; $i++){
for($j = 0; $j < $rowsB; $j ++){
$num = 0;
for($z = 0; $z < $rowsA; $z ++){
$num += $a[$i][$z] * $b[$z][$j];
}
$return[$i][$j] = $num;
}
}
return $return;
}
/*
*求矩陣的模
* @refer http://blog.csdn.net/shanshanpt/article/details/16820325
*/
function calculate_A( $src, $n )
{
$result = 0.0;
if( $n == 1 )
{
return $src[0][0];
}
for( $i = 0; $i < $n; ++$i )
{
for( $j = 0; $j < $n - 1; ++$j )
{
for( $k = 0; $k < $n - 1; ++$k )
{
$x = $j + 1;
$y = $k >= $i ? $k + 1 : $k;
$tmp[$j][$k] = $src[$x][$y];
}
}
$t = calculate_A( $tmp, $n - 1 );
if( $i % 2 == 0 )
{
$result += $src[0][$i] * $t;
}
else
{
$result -= $src[0][$i] * $t;
}
}
return $result;
}
/*
* 伴隨矩陣
* @refer http://blog.csdn.net/shanshanpt/article/details/16820325
*/
function getAdjoint( $origin )
{
$n = count($origin);
if( $n == 1 )
{
$dst[][] = 1;
return;
}
for( $i = 0; $i < $n; ++$i )
{
for( $j = 0; $j < $n; ++$j )
{
for( $k = 0; $k < $n - 1; ++$k )
{
for( $t = 0; $t < $n - 1; ++$t )
{
$x = $k >= $i ? $k + 1 : $k ;
$y = $t >= $j ? $t + 1 : $t;
$tmp[$k][$t] = $origin[$x][$y];
}
}
$ret[$j][$i] = calculate_A( $tmp, $n - 1 );
if( ( $i + $j ) % 2 == 1 )
{
$ret[$j][$i] = -1*$ret[$j][$i];
}
$ret[$j][$i] = 1/260*$ret[$j][$i];
}
}
return $ret;
}
/***************** EXAMPLE ******************/
$a = [[1, 2], [1, 6], [1, 9], [1, 13]];
$y = [[4],[8],[12],[21]];
$aT = getMatrixT($a);
$aTa = matrixMulti($aT, $a);
$t = getAdjoint($aTa);
$ret = matrixMulti($t, $aT);
$ret = matrixMulti($ret, $y);
var_dump($ret);
//$dataset = [[1,4],[2,5],[5,1],[4,2]];
//$dataret = [19,26, 19, 20];
//$expect = [10, 10];
//$step = 0.001;
//$times = 1000000;
求伴隨矩陣和矩陣的模的兩個函數(shù) 自己寫的不太好枚粘,運行起來比代碼里用的函數(shù)慢多了馅闽。。就不放出來了,這里用的是參考一篇文章的代碼福也。原代碼是C語言寫的局骤。謝謝[shanshanpt] http://blog.csdn.net/shanshanpt/article/details/16820325。
還有一些相關的函數(shù)暴凑。
/*
* 獲得逆序數(shù)
* @auther menmei
* @date 2017/03/27
*/
function getInversionNumber($arr){
$Total = count($arr);
$inversionNum = 0;
$exist = [];
for($i = 0; $i < $Total; $i ++){
if($arr[$i] >= 1){
$tmp = range(0, $arr[$i] - 1);
//在tmp中有 但exist沒有的峦甩。
$diff = array_diff($tmp, $exist);
$inversionNum += count($diff);
}
$exist[] = $arr[$i];
}
return $inversionNum;
}