zoukankan      html  css  js  c++  java
  • 最小二乘法

    原文:http://blog.csdn.net/jairuschan/article/details/7517773/

    概念

    最小二乘法多项式曲线拟合,根据给定的m个点,并不要求这条曲线精确地经过这些点,而是曲线y=f(x)的近似曲线y= φ(x)。

    原理

    [原理部分由个人根据互联网上的资料进行总结,希望对大家能有用]

         给定数据点pi(xi,yi),其中i=1,2,…,m。求近似曲线y= φ(x)。并且使得近似曲线与y=f(x)的偏差最小。近似曲线在点pi处的偏差δi= φ(xi)-y,i=1,2,...,m。 

    常见的曲线拟合方法:

         1.使偏差绝对值之和最小

         

         2.使偏差绝对值最大的最小

         

         3.使偏差平方和最小

         

         按偏差平方和最小的原则选取拟合曲线,并且采取二项式方程为拟合曲线的方法,称为最小二乘法。

    推导过程:

         1. 设拟合多项式为:

              

         2. 各点到这条曲线的距离之和,即偏差平方和如下:

              

         3. 为了求得符合条件的a值,对等式右边求ai偏导数,因而我们得到了: 

              

              

                             .......

              

         4. 将等式左边进行一下化简,然后应该可以得到下面的等式:

              

              

                         .......

              

         5. 把这些等式表示成矩阵的形式,就可以得到下面的矩阵:

              

         6. 将这个范德蒙得矩阵化简后可得到:

              

         7. 也就是说X*A=Y,那么A = (X'*X)-1*X'*Y,便得到了系数矩阵A,同时,我们也就得到了拟合曲线。

    实现

    运行前提:

    1. Python运行环境与编辑环境;
    2. Matplotlib.pyplot图形库,可用于快速绘制2D图表,与matlab中的plot命令类似,而且用法也基本相同。

    代码:

     1 # coding=utf-8
     2 
     3 '''
     4 作者:Jairus Chan
     5 程序:多项式曲线拟合算法
     6 '''
     7 import matplotlib.pyplot as plt
     8 import math
     9 import numpy
    10 import random
    11 
    12 fig = plt.figure()
    13 ax = fig.add_subplot(111)
    14 
    15 #阶数为9阶
    16 order=9
    17 
    18 #生成曲线上的各个点
    19 x = numpy.arange(-1,1,0.02)
    20 y = [((a*a-1)*(a*a-1)*(a*a-1)+0.5)*numpy.sin(a*2) for a in x]
    21 #ax.plot(x,y,color='r',linestyle='-',marker='')
    22 #,label="(a*a-1)*(a*a-1)*(a*a-1)+0.5"
    23 
    24 #生成的曲线上的各个点偏移一下,并放入到xa,ya中去
    25 i=0
    26 xa=[]
    27 ya=[]
    28 for xx in x:
    29     yy=y[i]
    30     d=float(random.randint(60,140))/100
    31     #ax.plot([xx*d],[yy*d],color='m',linestyle='',marker='.')
    32     i+=1
    33     xa.append(xx*d)
    34     ya.append(yy*d)
    35 
    36 '''for i in range(0,5):
    37     xx=float(random.randint(-100,100))/100
    38     yy=float(random.randint(-60,60))/100
    39     xa.append(xx)
    40     ya.append(yy)'''
    41 
    42 ax.plot(xa,ya,color='m',linestyle='',marker='.')
    43 
    44 
    45 #进行曲线拟合
    46 matA=[]
    47 for i in range(0,order+1):
    48     matA1=[]
    49     for j in range(0,order+1):
    50         tx=0.0
    51         for k in range(0,len(xa)):
    52             dx=1.0
    53             for l in range(0,j+i):
    54                 dx=dx*xa[k]
    55             tx+=dx
    56         matA1.append(tx)
    57     matA.append(matA1)
    58 
    59 #print(len(xa))
    60 #print(matA[0][0])
    61 matA=numpy.array(matA)
    62 
    63 matB=[]
    64 for i in range(0,order+1):
    65     ty=0.0
    66     for k in range(0,len(xa)):
    67         dy=1.0
    68         for l in range(0,i):
    69             dy=dy*xa[k]
    70         ty+=ya[k]*dy
    71     matB.append(ty)
    72  
    73 matB=numpy.array(matB)
    74 
    75 matAA=numpy.linalg.solve(matA,matB)
    76 
    77 #画出拟合后的曲线
    78 #print(matAA)
    79 xxa= numpy.arange(-1,1.06,0.01)
    80 yya=[]
    81 for i in range(0,len(xxa)):
    82     yy=0.0
    83     for j in range(0,order+1):
    84         dy=1.0
    85         for k in range(0,j):
    86             dy*=xxa[i]
    87         dy*=matAA[j]
    88         yy+=dy
    89     yya.append(yy)
    90 ax.plot(xxa,yya,color='g',linestyle='-',marker='')
    91 
    92 ax.legend()
    93 plt.show()
  • 相关阅读:
    LeetCode(Weekly Contest 187)题解
    2021 暑期实习招聘思之未雨绸缪
    LeetCode(Weekly Contest 185)题解
    构建私有的 docker registry
    LeetCode(Weekly Contest 182)题解
    Win10图标变白
    远程连接linux上的mysql
    IDEA复制粘贴html文件打不进war包问题
    Git命令
    数据源收集
  • 原文地址:https://www.cnblogs.com/hpuwangjunling/p/4374454.html
Copyright © 2011-2022 走看看