最小二乘法原理與MATLAB代碼——線性拟合、多項式拟合、指定函數拟合

最小二乘法

原文連結：原創，轉載請标明出處！

最小二乘法：

已知四個資料點，尋找一條直線（曲線），使得每個點與線的距離的總合最小。

使損失函數Q最小，得到參數β0，β1，對這個兩個參數求偏導數

複現了下面文章中最小二乘法的不同類型的拟合方法

文章連結： https://wenku.baidu.com/view/1071a6db6f1aff00bed51e66.html

1、線性拟合

已知一組資料點，滿足函數y=ax+b，利用最小二乘法的原理求解a和b。

圖中數字10是有10個資料點

matlab代碼

%% 線性拟合 Linear fitting
x=[0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1];
y=[2.3201,2.6470,2.9070,3.2885,3.6008,3.9090,4.2147,4.5191,4.8232,5.1275];
figure(1);
plot(x,y,'o');
hold on;
[a,b]=Linear_fitting(x,y);

t=0.1:0.1:1;
z=a*t+b;
plot(t,z);

function [a,b]=Linear_fitting(x,y)
n=size(x,2);
sum_xy=sum(x.*y);
sum_x=sum(x);
sum_y=sum(y);
sum_xx=sum(x.*x);

a=(n*sum_xy-sum_x*sum_y)/(n*sum_xx-sum_x^2);
b=(sum_xy-sum_xx*sum_y/sum_x)/(sum_x-n*sum_xx/sum_x);

% a=(n*sum(x.*y)-sum(x)*sum(y))/(n*sum(x.*x)-sum(x)*sum(x));
% b=(sum(x.*y)-sum(x.*x)*sum(y)/sum(x))/(sum(x)-n*sum(x.*x)/sum(x));
end
           

2、多項式拟合

已知一組資料點，滿足多項式函數

利用最小二乘法的原理求解參數a0,a1,…,ak。

得到損失函數

為了目标函數最優，對每個參數求偏導數，并令k個偏導數為0.

利用矩陣乘法形式XA=Y求解A

講解連結：最小二乘法—多項式拟合非線性函數 - 簡書 (jianshu.com)

MATLAB最小二乘法 - 凱魯嘎吉 - 部落格園 (cnblogs.com)

matlab代碼

x=[0.1732;0.1775;0.1819;0.1862;0.1905;0.1949;0.1992;0.2035;0.2079;0.2122;0.2165;0.2208;0.2252;0.2295;0.2338;0.2384];
y=[-3.41709;-4.90887;-6.09424;-6.95362;-7.63729;-8.12466;-8.37153;-8.55049;-8.61958;-8.65326;-8.60021;-8.52824;-8.43502;-8.32234;-8.20419;-8.04472];
m=3;
figure(2);
plot(x,y,'o');
hold on;
p=Polynomial_fitting(x,y,m);
Polynomial_fitting_plot(p,x,m);%繪制圖像

function P=Polynomial_fitting(x,y,m)  %x,y為序列長度相等的資料向量，m為拟合多項式次數
a=zeros(2*m+1,1);
for i=0:2*m
    a(i+1,1)=sum(x.^i);
end
b=zeros(m+1,1);
for i=0:m
    b(i+1,1)=sum((x.^i).*y);
end
A=zeros(m+1,m+1);
for i=1:m+1
    for j=1:m+1
        A(i,j)=a(i+j-1,1);
    end
end
p=A\b;
P=fliplr(p'); %左右翻轉矩陣
%https://www.cnblogs.com/kailugaji/p/6932482.html
end

function Polynomial_fitting_plot(p,t,m)
z=zeros(1,size(t,2));
for i=0:m
    z=z+p(m-i+1)*t.^i;
end
plot(t,z);
end
           

3、指定非線性函數拟合

已知一組資料點，滿足函數

y = a ( 1 ) x 2 + a ( 2 ) s i n ( x ) + a ( 3 ) x 3 y=a(1)x^2+a(2)sin(x)+a(3)x^3 y=a(1)x2+a(2)sin(x)+a(3)x3

利用最小二乘法求解參數a(1),a(2),a(3)

matlab 中提供一個lsqcurvefit函數

x = [3.6 7.7 9.3 4.1 8.6 2.8 1.3 7.9 10.0 5.4];
y = [16.5 150.6 263.1 24.7 208.5 9.9 2.7 163.9 325.0 54.3];
a=[1 1 1]; %初始估計值
[p,resnorm]=lsqcurvefit(@fun,a,x,y); %指定拟合曲線    resnorm殘差的平方範數，傳回為非負實。
figure(3);
plot(x,y,'o');
hold on;
t=1:1:10;
z=fun(p,t);
plot(t,z);

function F=fun(a,x)
F=a(1)*x.^2+a(2)*sin(x)+a(3)*x.^3;

%https://wenku.baidu.com/view/1071a6db6f1aff00bed51e66.html
end
           

原文連結：原創，轉載請标明出處！