最小二乘法曲線擬合

除了梯度下降法去擬合曲線,最小二乘法也是另一個方法娄琉。和梯度下降不斷迭代求出極值不同的是曲横,最小二乘法是直接求導計算出參數(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;
}

最后編輯于
?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者
  • 序言:七十年代末,一起剝皮案震驚了整個濱河市现喳,隨后出現(xiàn)的幾起案子凯傲,更是在濱河造成了極大的恐慌,老刑警劉巖嗦篱,帶你破解...
    沈念sama閱讀 211,194評論 6 490
  • 序言:濱河連續(xù)發(fā)生了三起死亡事件冰单,死亡現(xiàn)場離奇詭異,居然都是意外死亡灸促,警方通過查閱死者的電腦和手機诫欠,發(fā)現(xiàn)死者居然都...
    沈念sama閱讀 90,058評論 2 385
  • 文/潘曉璐 我一進店門,熙熙樓的掌柜王于貴愁眉苦臉地迎上來浴栽,“玉大人荒叼,你說我怎么就攤上這事〉浼Γ” “怎么了被廓?”我有些...
    開封第一講書人閱讀 156,780評論 0 346
  • 文/不壞的土叔 我叫張陵,是天一觀的道長椿每。 經(jīng)常有香客問我伊者,道長,這世上最難降的妖魔是什么间护? 我笑而不...
    開封第一講書人閱讀 56,388評論 1 283
  • 正文 為了忘掉前任亦渗,我火速辦了婚禮,結(jié)果婚禮上汁尺,老公的妹妹穿的比我還像新娘法精。我一直安慰自己,他們只是感情好痴突,可當我...
    茶點故事閱讀 65,430評論 5 384
  • 文/花漫 我一把揭開白布搂蜓。 她就那樣靜靜地躺著,像睡著了一般辽装。 火紅的嫁衣襯著肌膚如雪帮碰。 梳的紋絲不亂的頭發(fā)上,一...
    開封第一講書人閱讀 49,764評論 1 290
  • 那天拾积,我揣著相機與錄音殉挽,去河邊找鬼丰涉。 笑死,一個胖子當著我的面吹牛斯碌,可吹牛的內(nèi)容都是我干的一死。 我是一名探鬼主播,決...
    沈念sama閱讀 38,907評論 3 406
  • 文/蒼蘭香墨 我猛地睜開眼傻唾,長吁一口氣:“原來是場噩夢啊……” “哼投慈!你這毒婦竟也來了?” 一聲冷哼從身側(cè)響起冠骄,我...
    開封第一講書人閱讀 37,679評論 0 266
  • 序言:老撾萬榮一對情侶失蹤伪煤,失蹤者是張志新(化名)和其女友劉穎,沒想到半個月后猴抹,有當?shù)厝嗽跇淞掷锇l(fā)現(xiàn)了一具尸體带族,經(jīng)...
    沈念sama閱讀 44,122評論 1 303
  • 正文 獨居荒郊野嶺守林人離奇死亡,尸身上長有42處帶血的膿包…… 初始之章·張勛 以下內(nèi)容為張勛視角 年9月15日...
    茶點故事閱讀 36,459評論 2 325
  • 正文 我和宋清朗相戀三年蟀给,在試婚紗的時候發(fā)現(xiàn)自己被綠了蝙砌。 大學時的朋友給我發(fā)了我未婚夫和他白月光在一起吃飯的照片。...
    茶點故事閱讀 38,605評論 1 340
  • 序言:一個原本活蹦亂跳的男人離奇死亡跋理,死狀恐怖择克,靈堂內(nèi)的尸體忽然破棺而出,到底是詐尸還是另有隱情前普,我是刑警寧澤肚邢,帶...
    沈念sama閱讀 34,270評論 4 329
  • 正文 年R本政府宣布,位于F島的核電站拭卿,受9級特大地震影響骡湖,放射性物質(zhì)發(fā)生泄漏。R本人自食惡果不足惜峻厚,卻給世界環(huán)境...
    茶點故事閱讀 39,867評論 3 312
  • 文/蒙蒙 一响蕴、第九天 我趴在偏房一處隱蔽的房頂上張望。 院中可真熱鬧惠桃,春花似錦浦夷、人聲如沸。這莊子的主人今日做“春日...
    開封第一講書人閱讀 30,734評論 0 21
  • 文/蒼蘭香墨 我抬頭看了看天上的太陽。三九已至呐馆,卻和暖如春肥缔,著一層夾襖步出監(jiān)牢的瞬間,已是汗流浹背汹来。 一陣腳步聲響...
    開封第一講書人閱讀 31,961評論 1 265
  • 我被黑心中介騙來泰國打工续膳, 沒想到剛下飛機就差點兒被人妖公主榨干…… 1. 我叫王不留怒见,地道東北人。 一個月前我還...
    沈念sama閱讀 46,297評論 2 360
  • 正文 我出身青樓姑宽,卻偏偏與公主長得像,于是被迫代替她去往敵國和親闺阱。 傳聞我的和親對象是個殘疾皇子炮车,可洞房花燭夜當晚...
    茶點故事閱讀 43,472評論 2 348

推薦閱讀更多精彩內(nèi)容