文档章节

曲线拟合最小二乘法算法设计及Scala和C++代码实现

安艳青
 安艳青
发布于 2017/01/09 20:20
字数 1165
阅读 1457
收藏 11

一、曲线拟合最小二乘原理

    设实验数据集为,构造的拟合曲线为,则在插值节点处产生的误差为,若想使能够描述原曲线,则应使最小。常使用下列三种标准来度量误差的大小:

                        

    曲线拟合最小二乘法选择的是使用向量2-范数最为总体误差的度量标准。

    具体地,对给定一组实验数据拟构造的曲线为,其中是基函数是当为极小值时的解。

    设多元函数,令 ,具体如下:

               

                

    设,则上式可以写为:把它改写成矩阵形式为:

                        

    这是关于系数的线性方程组,也称为正则方程。由于线性无关,故上述线性方程组的系数矩阵行列式不为零,因此方程组有唯一解

二、算法设计

    曲线拟合最小二乘法算法设计主要分成两部分,首先应得到线性方程组系数的增广矩阵,然后再求解线性方程组,方程组的求解可用高斯消去法得到。

    多项式型经验公式的计算步骤:

    (1)输入实验数据:

    (2)构造正则方程:

             右端项:

             系数矩阵:

    (3)用高斯消去法求解正则方程;

    (4)输出经验公式。

    具体算法描述如下图:

                            

    其中,n为插值节点下标上限,m为拟合曲线最高次方,X,Y为插值节点向量。

三、Scala实现

    在IDEA里实现这个算法,具体代码如下:

package com.hzznv.ls

import breeze.numerics.pow

/**
  * Created by YanqingAn on 2017/1/6.
  */
object LeastSquare_2 {

  def main(args: Array[String]): Unit = {

    val x = new Array[Double](20)
    val y = new Array[Double](20)
    val a = Array.ofDim[Double](20, 20)
    var t: Double = 0

    var m: Int = 0
    var n: Int = 0

    n = readLine("请输入n:").toInt
    val xy:String = readLine("请输入xy:")
    for(i <- 0 to 2*n-1){
      if(i % 2 == 0){
        x(i/2) = xy.split(" ")(i).toDouble
      }else if(i % 2 == 1) {
        y((i-1)/2) = xy.split(" ")(i).toDouble
      }
    }
    m = readLine("请输入m:").toInt
    println("输入的插值节点是:")
    for(i <- 0 to n-1){
      print(x(i) + " ")
    }
    println()
    for(i <- 0 to n-1){
      print(y(i) + " ")
    }
    println()
    for(i <- 0 to m){
      for(j <- 0 to m){
        for(k <- 0 to n-1){
          a(i)(j) += pow(x(k), i+j)
        }
      }
      for(k <- 0 to n-1){
        a(i)(m+1) += y(k) * pow(x(k), i)
      }
    }
    println("增广矩阵是:")
    for(i <- 0 to m){
      for(j <- 0 to m+1){
        print(" " + a(i)(j))
      }
      println()
    }

    for(k <- 0 to m-1){
      for(i <- k+1 to m){
        t = -a(i)(k) / a(k)(k)
        for(j <- k+1 to m+1){
          a(i)(j) += a(k)(j) * t
        }
      }
    }

    for(i <- (0 to m).reverse){
      for(j <- i+1 to m){
        a(i)(m+1) -= a(i)(j) * a(j)(m+1)
      }
      a(i)(m+1) /= a(i)(i)
      a(i)(m+1) = f"${a(i)(m+1)}%1.5f".toDouble
    }

    print("拟合曲线是: y=" + a(0)(m+1))
    for(i <- 1 to m){
      if(a(i)(m+1) >= 0){
        print("+")
      }
      print(a(i)(m+1))
      for(j <- 0 until i){
        print("*x")
      }
    }
  }
}


    执行过程和结果如下:

四、附C++代码

3 C++实现 在 Vi s ua l
 c++6. 0里实现这个算法 ,具体代码 如下:
 
#i ncl ude<i ost ream. h>
#i ncl ude” mat h. h”
 
int  main (){
 double x[ 20], y[ 20], a[ 20] [ 20]= { 0}, t ;
 int  i , j, m, n, k;
 cout <<”请输入 n”<<e nd l :
 cin>>n;
 
 cout <<”请输入 X Y”<<endl :
 for( i -1; i <=n; i ++)
   cin>>x[ i ]>>y[ i ];
 cout <<”请输入 m”<<endl :
  cin>>m;
 cout <<endl <<”输入的插值节点是:”<<endl ;
for( i =1; i <=n. i ++)
  cout <<x[ i ]<<” ”;
cout <<endl ;
for( i =1; i <=n; i ++)
cout <<y[ i ]<<” ”;
cout <<endl ;
for( i =0; i <=m; i ++)
 {
 for( j =0; j <=m ; j ++)
 for( k=1; k<=n; k++)
 a[ i ] [ j ]+=pow ( x[ k ], i +j );
 
for( k=1; k<=n; k++)
 a[ i ] [ j ]+=y[ k] * pow (x[ k ], i );
}
 cout <<endl <<”增广矩阵是:”<<endl ;
 
for( i =0; i <=m; i ++)
 
{ for( j =0; j <=m+1; j ++)
 co ut <<” ”<<a[ i ] [ j ];
 cout <<endl ; }
for( k=0; k<=m-1 ; k++)
 for( i =k+l ; i <=m; i ++)
 { t =-a[ i ] [ k]/ a[ k ] [ k ] k ;
 for( j =k +1 ; j <= m+1 ; j ++ )
 a[ i ] [ j ]+=a[ k ] [ j ] * t ;
 }
 for(i=m;i>=0;i--){

    for(j=i+1;j<=m;j++)

    a[i][m+1] -= a[i][j] * a[j][m+1];

    a[i][m+1] /= a[i][i];}

cout<<endl<<" 拟合曲线是:";

cout<<" y="<<a[0][m+1];

for(i=1;i<=m;i++){

    if(a[i][m+1]>0) cout<<" +";

    cout<<a[i][m+1];

    for(j=0;j<i;j++)

    cout<<" *x";

}

五、结束语

    曲线拟合最小二乘算法在插值问题的数值算法中是比较经典的算法,它应用的领域也比较宽泛,比如在轨道曲线拟合、人口预测、煤矿瓦斯预警等问题上都有突出贡献。

 

© 著作权归作者所有

下一篇: 机器学习
安艳青
粉丝 2
博文 6
码字总数 9107
作品 0
杭州
程序员
私信 提问
加载中

评论(4)

库特
库特
scala写得不那么 "scala",所以用scala实现是伪命题,换汤不换药
库特
库特
scala写得不那么 "scala",所以用scala实现是伪命题,换汤不换药
安艳青
安艳青 博主

引用来自“IdleMan”的评论

研究数据的杰出女青年:bowtie:
😊刚开始入门儿旅程的数据小菜鸟~
IdleMan
IdleMan
研究数据的杰出女青年:bowtie:
SP++3.0 发布,欢迎大家使用

消息来自 Jerry 的博客: SP++ (Signal Processing in C++) 是一个关于信号处理与数值计算的开源C++程序库,该库提供了信号处理与数值计算中常用算法的C++实现。SP++中所有算法都以C++类模板...

红薯
2011/02/12
5.3K
4
SP++3.0已发布,欢迎大家使用(同心协力,共创开源)

SP++ (Signal Processing in C++) 是一个关于信号处理与数值计算的开源C++程序库,该库提供了信号处理与数值计算中常用算法的C++实现。SP++中所有算法都以C++类模板方法实现,以头文件形式组...

张明
2011/02/12
11.5K
56
C语言/C++编程学习:算法之排序:折半插入法

C语言是面向过程的,而C++是面向对象的 C和C++的区别: C是一个结构化语言,它的重点在于算法和数据结构。C程序的设计首要考虑的是如何通过一个过程,对输入(或环境条件)进行运算处理得到...

小辰带你看世界
2018/05/13
0
0
泛型编程与设计新思维

作者: 徐景周 转帖: http://www.vckbase.com/document/viewdoc/?id=955 前言 永远记住,编写代码的宗旨在于简单明了,不要使用语言中的冷僻特性,耍小聪明,重要的是编写你理解的代码,理解...

ValueError
2011/01/12
357
1
Scala 不是改良的 Java,你会考虑使用 Scale 吗?

Scala 本身就基于 Java 平台,却要来跟 Java 比较,这有点不靠谱。 Scala编程语言拥有所有Java的语言特征,而且还支持所有的新兴的有趣的概念,例如闭包,higher-kinded类型,内联XML。如果你...

绿悠悠
2010/07/28
3.6K
10

没有更多内容

加载失败,请刷新页面

加载更多

OSChina 周六乱弹 —— 早上儿子问我他是怎么来的

Osc乱弹歌单(2019)请戳(这里) 【今日歌曲】 @凉小生 :#今日歌曲推荐# 少点戾气,愿你和这个世界温柔以待。中岛美嘉的单曲《僕が死のうと思ったのは (曾经我也想过一了百了)》 《僕が死の...

小小编辑
今天
2.4K
15
Excption与Error包结构,OOM 你遇到过哪些情况,SOF 你遇到过哪些情况

Throwable 是 Java 中所有错误与异常的超类,Throwable 包含两个子类,Error 与 Exception 。用于指示发生了异常情况。 Java 抛出的 Throwable 可以分成三种类型。 被检查异常(checked Exc...

Garphy
今天
41
0
计算机实现原理专题--二进制减法器(二)

在计算机实现原理专题--二进制减法器(一)中说明了基本原理,现准备说明如何来实现。 首先第一步255-b运算相当于对b进行按位取反,因此可将8个非门组成如下图的形式: 由于每次做减法时,我...

FAT_mt
昨天
40
0
好程序员大数据学习路线分享函数+map映射+元祖

好程序员大数据学习路线分享函数+map映射+元祖,大数据各个平台上的语言实现 hadoop 由java实现,2003年至今,三大块:数据处理,数据存储,数据计算 存储: hbase --> 数据成表 处理: hive --> 数...

好程序员官方
昨天
61
0
tabel 中含有复选框的列 数据理解

1、el-ui中实现某一列为复选框 实现多选非常简单: 手动添加一个el-table-column,设type属性为selction即可; 2、@selection-change事件:选项发生勾选状态变化时触发该事件 <el-table @sel...

everthing
昨天
21
0

没有更多内容

加载失败,请刷新页面

加载更多

返回顶部
顶部