国产色视频在线,99久久中文字幕组,久久日本亚洲

上節(jié)課主要介紹了特征值與特征向量的概念，低階矩陣的特征值可以通過(guò)列出特征方程求解，高階矩陣則可以通過(guò)冪法與反冪法迭代求解出最大特征值與最小特征值（模），要求出矩陣的全部特征值則需要借助矩陣的 QR分解來(lái)將矩陣相似化為一個(gè)上三角矩陣，相似化過(guò)程不改變矩陣的特征值，因此轉(zhuǎn)化后的上三角矩陣的對(duì)角線元素即為原矩陣的特征值。本節(jié)課主要介紹曲線擬合（機(jī)器學(xué)習(xí)中的線性回歸算法）與插值算法。

1. 擬合與插值（Fitting, Interpolation）

曲線擬合是給定一個(gè)數(shù)據(jù)集，構(gòu)建一個(gè)數(shù)學(xué)方程能夠最好的匹配這些數(shù)據(jù)點(diǎn)
插值則是構(gòu)建一個(gè)多項(xiàng)式方程來(lái)完美地契合給定的數(shù)據(jù)集，當(dāng)數(shù)據(jù)點(diǎn)比較少的時(shí)候，一個(gè)簡(jiǎn)單的多項(xiàng)式方程就可以實(shí)現(xiàn)插值，而數(shù)據(jù)點(diǎn)比較多的時(shí)候，更多地采用的是樣條插值，在各個(gè)區(qū)間構(gòu)建不同的多項(xiàng)式方程

Fitting and Interpolation.png

2. 線性曲線擬合

線性曲線擬合是使用如下線性方程（一元一階多項(xiàng)式，對(duì)應(yīng)機(jī)器學(xué)習(xí)線性回歸算法中只有一個(gè)特征的情況）來(lái)擬合給定的數(shù)據(jù)集： $f(x)=a_1x+a_0$ 若數(shù)據(jù)集中只有兩個(gè)點(diǎn)，則該線性方程可以準(zhǔn)確求解出來(lái)，通過(guò)這兩個(gè)點(diǎn)；當(dāng)數(shù)據(jù)集超過(guò)兩個(gè)點(diǎn)時(shí)，該線性方程無(wú)法通過(guò)所有的數(shù)據(jù)點(diǎn)，此時(shí)需要找出擬合數(shù)據(jù)點(diǎn)最好的那個(gè)方程。

image.png

如何衡量擬合的準(zhǔn)確度？對(duì)數(shù)據(jù)點(diǎn)，定義和的差值為殘差：

image.png

一個(gè)可選的衡量測(cè)度是計(jì)算殘差和：然而這種測(cè)度有時(shí)候會(huì)不好用，因?yàn)楫?dāng)有些點(diǎn)殘差為很大的正數(shù)，而有些點(diǎn)殘差為很小的負(fù)數(shù)時(shí)，將會(huì)給出一個(gè)接近于的殘差和，

image.png

另外一個(gè)可選的測(cè)度是殘差的絕對(duì)值和：該測(cè)度總是正的，可以用來(lái)衡量擬合的精度，但無(wú)法用來(lái)確定哪個(gè)是最好的擬合曲線，因?yàn)椴煌臄M合曲線有可能給出相同的殘差絕對(duì)值和。

image.png

因此，一個(gè)誤差函數(shù)應(yīng)該既能夠衡量擬合曲線的好壞，又能夠用來(lái)確定出唯一地?cái)M合曲線，目前廣泛采用的是殘差的平方和：采用線性最小二乘回歸算法就可以找出使得上述誤差最小的那條擬合曲線。給定數(shù)據(jù)集，是關(guān)于的非線性方程，結(jié)合微積分相關(guān)概念，可知取最小值的點(diǎn)滿(mǎn)足：化簡(jiǎn)有：其解為：如果定義則有
MATLAB實(shí)現(xiàn)：

function [a1,a0] = LinearRegression(x, y)
% LinearRegration calculates the coefficients a1 and a0 of the linear
% equation y = a1*x + a0 that best fit n data points.
% Input variables:
% x    A vector with the coordinates x of the data points.
% y    A vector with the coordinates y of the data points.
% Output variable:
% a1   The coefficient a1.
% a0   The coefficient a0.

nx = length(x);
ny = length(y);
if nx ~= ny
    disp('ERROR: The number of elements in x must be the same as in y.')
    a1 = 'Error';
    a0 = 'Error';
else
Sx = sum(x);
Sy = sum(y);
Sxy = sum(x.*y);
Sxx = sum(x.^2);
a1 = (nx*Sxy - Sx*Sy)/(nx*Sxx - Sx^2);
a0 = (Sxx*Sy - Sxy*Sx)/(nx*Sxx - Sx^2);
end

3. 非線性曲線擬合

很多實(shí)際應(yīng)用中，使用線性方程的擬合效果不夠好，采用非線性擬合則能夠達(dá)到比較好的效果。

image.png

3.1 常見(jiàn)的幾種非線性方程

（這里同線性擬合類(lèi)似，只考慮單變量的情況，對(duì)應(yīng)線性回歸算法中的一個(gè)特征）： $y=bx^m$ $y=be^{mx}$ $y=\frac{1}{mx+b}$ 對(duì)非線性擬合，很容易將其改寫(xiě)為線性形式： $ln(y)=mln(x)+ln(b)$ $ln(y) = mx+ln(b)$ $\frac{1}{y}=mx+b$ 將數(shù)據(jù)集進(jìn)行相應(yīng)的特征變換后，采用線性最小二乘法求解相應(yīng)的參數(shù)。然而，對(duì)于給定的數(shù)據(jù)集，如何判斷使用哪一種非線性方程進(jìn)行擬合呢？一個(gè)可行的方案是首先通過(guò)繪制數(shù)據(jù)集中的點(diǎn)（當(dāng)特征很多時(shí)，很難可視化），確定是選擇線性還是非線性進(jìn)行擬合，如果選擇非線性擬合，具體選擇哪一種非線性方程則要取決于其背后的物理現(xiàn)象和數(shù)學(xué)原理。同樣的，可以通過(guò)以特定的方式繪制數(shù)據(jù)點(diǎn)并檢查這些點(diǎn)是否符合直線來(lái)預(yù)測(cè)所提出的非線性函數(shù)是否具有很好的擬合能力。
MATLAB實(shí)現(xiàn)：

texp = 2:2:30;
vexp = [9.7 8.1 6.6 5.1 4.4 3.7 2.8 2.4 2.0 1.6 1.4 1.1 0.85 0.69 0.6];
vexpLOG = log(vexp);
R = 5E6;
[a1,a0] = LinearRegression(texp, vexpLOG)
b = exp(a0)
C = -1/(R*a1)
t = 0:0.5:30;
v = b*exp(a1*t);
plot(t,v,texp,vexp,'or')

image.png

3.2 二次及高階多項(xiàng)式擬合

多項(xiàng)式是指具有如下形式的函數(shù)（單變量）： $f(x)=a_nx^n+a_{n-1}x^{n-1}+...+a_1x+a_0$ $n$ 稱(chēng)作多項(xiàng)式的階數(shù)。如果給定數(shù)據(jù)集中有 $n$ 個(gè)點(diǎn)，則可以確定 $n-1$ 個(gè)從1階到 $n-1$ 階的擬合多項(xiàng)式，其中 $n-1$ 階多項(xiàng)式可以實(shí)現(xiàn)通過(guò)每一個(gè)數(shù)據(jù)點(diǎn)，但并不是多項(xiàng)式的階數(shù)越高越好，這取決于實(shí)際問(wèn)題，當(dāng)數(shù)據(jù)集點(diǎn)本身就具有噪聲污染時(shí)，提高擬合的階數(shù)并沒(méi)有意義。（同機(jī)器學(xué)習(xí)中的過(guò)擬合，overfitting類(lèi)似）

image.png

對(duì)多項(xiàng)式擬合，采用多項(xiàng)式回歸算法來(lái)確定各個(gè)參數(shù)。給定個(gè)數(shù)據(jù)點(diǎn)，想要使用階多項(xiàng)式對(duì)其進(jìn)行擬合：以為例，首先定義總誤差：要讓總誤差最小，取對(duì)應(yīng)的導(dǎo)數(shù)為即可：將其更改為線性系統(tǒng)形式：
該線性系統(tǒng)的解即為多項(xiàng)式的各個(gè)參數(shù)。
MATLAB實(shí)現(xiàn)：

clear, clc
x = 0:0.4:6;
%x=1:4
y = [0 3 4.5 5.8 5.9 5.8 6.2 7.4 9.6 15.6 20.7 26.7 31.1 35.6 39.3 41.5];
n = length(x);
m = 4;
for i = 1:2*m
    xsum(i) = sum(x.^(i));
end
% Beginning of Step 3
a(1,1) = n;
b(1,1)= sum(y);
for j= 2:m+1
    a(1,j)=xsum(j-1);
end
for i = 2:m+1
    for j=1:m+1
        a(i,j)= xsum(j+i-2);
    end
    b(i,1) = sum(x.^(i-1).*y);
end
% Step 4
p=(a\b)'
for i = 1:m+1
    Pcoef(i) = p(m+2-i);
end
epsilon=0:0.1:6;
stressfit=polyval(Pcoef,epsilon);
plot(x,y,'ro',epsilon,stressfit,'k','linewidth',2)
xlabel('Strain','fontsize',20)
ylabel('Stress (MPa)','fontsize',20)
%%%%%%%%%%%%%%%%%%%%
>> help polyval
polyval - 多項(xiàng)式計(jì)算

  此 MATLAB 函數(shù) 返回在 x 處計(jì)算的 n 次多項(xiàng)式的值。輸入?yún)?shù) p 是長(zhǎng)為 n+1 的向量，
其元素是按要計(jì)算的多項(xiàng)式降冪排序的系數(shù)。

    y = polyval(p,x)
    [y,delta] = polyval(p,x,S)
    y = polyval(p,x,[],mu)
    [y,delta] = polyval(p,x,S,mu)

    See also polyder, polyfit, polyint, polyvalm

    Reference page for polyval
    Other functions named polyval

image.png

上述多項(xiàng)式擬合可以擴(kuò)展到更一般的情況： $f(x)=c_1f_1(x)+c_2f_2(x)+...+c_mf_m(x)$ 定義總誤差為 $E=\sum^n_{i=1} [y_i-\sum^m_{j=1}c_jf_j(x_i)]^2$ 令偏導(dǎo)為 $0$ $\frac{\partial E}{\partial c_k} = 2\sum^n_{i=1}(y_i-\sum^m_{j=1}C_jf_j(x_i))(-f_k(x_i))=0,k=1,...,m$ 改為線性系統(tǒng)形式，有 $\sum^n_{i=1}\sum^m_{j=1}c_{j}f_j(x_i)f_k(x_i)=\sum^n_{i=1}y_if_k(x_i),k=1,...,m$ 矩陣形式為 $\begin{bmatrix}\sum^n_{i=1}f^2_1(x_i) & \sum^n_{i=1}f_1(x_i)f_2(x_i) & \cdots & \sum^n_{i=1}f_m(x_i)f_1(x_i) \\ \sum^n_{i=1}f_1(x_i)f_2(x_i) & \sum^n_{i=1}f^2_2(x_i) & \cdots & \sum^n_{i=1}f_m(x_i)f_2(x_i)\\ \vdots & \vdots & \ddots & \vdots\\ \sum^n_{i=1}f_1(x_i)f_m(x_i) & \sum^n_{i=1}f_2(x_i)f_m(x_i) & \cdots & \sum^n_{i=1}f^2_m(x_i) \end{bmatrix}\begin{bmatrix}c_1\\c_2\\ \vdots \\ c_m \end{bmatrix}=\begin{bmatrix}\sum^n_{i=1}y_if_1(x_i)\\\sum^n_{i=1}y_if_2(x_i) \\ \vdots \\ \sum^n_{i=1}y_if_m(x_i) \end{bmatrix}$
MATLAB實(shí)現(xiàn)：
$f(x)=\frac{A}{x}+\frac{Be^{-2x^2}}{x}$

clear all
x = [0.6 0.8 0.85 0.95 1.0 1.1 1.2 1.3 1.45 1.6 1.8];
y = [0.08 0.06 0.07 0.07 0.07 0.06 0.06 0.06 0.05 0.05 0.04];
a(1,1) = sum(1./x.^2);
a(1,2) = sum(exp(-2*x.^2)./x.^2);
a(2,1) = a(1,2);
a(2,2) = sum(exp(-4*x.^2)./x.^2);
b(1,1) = sum(y./x);
b(2,1) = sum((y.*exp(-2*x.^2))./x);
AB = a\b
xfit = 0.6:0.02:1.8;
yfit = AB(1)./xfit + AB(2)*exp(-2*xfit.^2)./xfit;
plot(x,y,'o',xfit,yfit)

image.png

4. 多項(xiàng)式插值

如前所述，插值是讓擬合的曲線通過(guò)所有給定的數(shù)據(jù)點(diǎn)（總誤差為0，機(jī)器學(xué)習(xí)中的 $E_{in}$ 為 $0$ ），對(duì) $n$ 個(gè)點(diǎn)，總是能找到一個(gè) $n-1$ 階的多項(xiàng)式函數(shù)通過(guò)所有的點(diǎn)： $P(x) = a_0+a_1x+...+a_mx^m=y$ 矩陣形式 $\begin{bmatrix}x_1^m & x_1^{m-1} & \cdots & x_1 & 1\\ x_2^m & x_2^{m-1} & \cdots & x_2 & 1\\ \vdots & \vdots & \ddots & \vdots & \vdots\\ x_n^m & x_n^{m-1} & \cdots & x_n & 1 \end{bmatrix}\begin{bmatrix}a_m \\ a_{m-1}\\ \vdots \\ a_0\end{bmatrix}=\begin{bmatrix}y_1 \\ y_2\\ \vdots \\ y_n\end{bmatrix}$

4.1 多項(xiàng)式插值的拉格朗日形式：

image.png

考慮兩個(gè)點(diǎn) $(x_1,y_1)$ , $(x_2,y_2)$ ，其一階拉格朗日多項(xiàng)式具有如下形式： $y=a_1(x-x_2)+a_2(x-x_1)$ 令該函數(shù)通過(guò)這兩個(gè)點(diǎn)，得 $a_1=\frac{y_1}{x_1-x_2}, a_2=\frac{y_2}{x_2-x_1}$ 因此多項(xiàng)式為： $y=\frac{y_1(x-x_2)}{x_1-x_2}+\frac{y_2(x-x_1)}{x_2-x_1}$ 類(lèi)似地，對(duì)三個(gè)點(diǎn) $(x_1,y_1),(x_2,y_2),(x_3,y_3)$ ，可得其二階拉格朗日具有形式 $y=\frac{(x-x_2)(x-x_3)}{(x_2-x_1)(x_3-x_1)}y_1+\frac{(x-x_1)(x-x_3)}{(x_1-x_2)(x_3-x_2)}y_2+\frac{(x-x_1)(x-x_2)}{(x_1-x_3)(x_1-x_3)}y_3$ 因此，拉格朗日多項(xiàng)式的通解形式為： $f(x)=\sum^n_{i=1}y_i\prod_{j\neq i}\frac{(x-x_j)}{(x_i-x_j)}$ MATLAB實(shí)現(xiàn)：

function Yint = LagrangeINT(x,y,Xint)
% LagrangeINT fits a Lagrange polynomial to a set of given points and
% uses the polynomial to determines the interpolated value of a point.
% Input variables:
% x  A vector with the x coordinates of the given points.
% y  A vector with the y coordinates of the given points.
% Xint  The x coordinate of the point to be interpolated.
% Output variable:
% Yint  The interpolated value of Xint.

n = length(x);
for i = 1:n
    L(i) = 1;
    for j = 1:n
        if j ~= i
            L(i)= L(i)*(Xint-x(j))/(x(i)-x(j));
        end
    end
end
y.*L
Yint = sum(y.*L);
%%%%%%%%
>> x = [1 2 4 5 7];
>> y = [52 5 -5 -40 10];
>> Yinterpolated = LagrangeINT(x,y,3)

ans =

   -5.7778    2.6667   -4.4444   13.3333    0.2222

Yinterpolated =

    6.0000

>> i=1;
>> for x1=0:0.1:10
y1(i)=LagrangeINT(x,y,x1);
i = i+1;
end
>> x1 = 0:0.1:10;
>> plot(x,y,'*',x1,y1)

image.png

4.2 牛頓多項(xiàng)式插值

上述拉格朗日多項(xiàng)式插值有一個(gè)缺陷，當(dāng)數(shù)據(jù)集增加一個(gè)插值點(diǎn)時(shí)，每個(gè)系數(shù)都需要重新計(jì)算，而牛頓多項(xiàng)式插值則不需要重新計(jì)算： $f(x)=a_1+a_2(x-x_1)+a_3(x-x_1)(x-x_2)+...+a_n(x-x_1)(x-x_2)...(x-x_{n-1})$

image.png

以?xún)蓚€(gè)數(shù)據(jù)點(diǎn)為例 , 容易計(jì)算出當(dāng)增加一個(gè)數(shù)據(jù)點(diǎn)時(shí)，可以看到 , ，這兩個(gè)參數(shù)都沒(méi)有發(fā)生變化，因此只需計(jì)算，通過(guò)得，一直進(jìn)行下去，就可以得到牛頓多項(xiàng)式插值的通解形式：
第一層：其中
第二層：其中
第三層：其中

image.png

MATLAB實(shí)現(xiàn)：

function Yint = NewtonsINT(x,y,Xint)
% NewtonsINT fits a Newtons polynomial to a set of given points and
% uses the polynomial to determines the interpolated value of a point.
% Input variables:
% x  A vector with the x coordinates of the given points.
% y  A vector with the y coordinates of the given points.
% Xint  The x coordinate of the point to be interpolated.
% Output variable:
% Yint  The interpolated value of Xint.

n = length(x);
a(1) = y(1);
for i = 1:n-1
    divDIF(i,1) = (y(i+1)-y(i))/(x(i+1)-x(i));
end
for j = 2:n-1
    for i = 1:n-j
        divDIF(i,j) = (divDIF(i+1,j-1) - divDIF(i,j-1))/(x(j+i) - x(i));
    end
end
for j=2:n
    a(j)=divDIF(1,j-1);
end
Yint=a(1);
xn=1;
for k=2:n
    xn=xn*(Xint-x(k-1));
    Yint=Yint+a(k)*xn;
end

5. 樣條插值

5.1 線性樣條插值

當(dāng)數(shù)據(jù)點(diǎn)比較少時(shí)，使用低階多項(xiàng)式就可以實(shí)現(xiàn)插值，然而當(dāng)數(shù)據(jù)點(diǎn)比較多時(shí)，插值的多項(xiàng)式需要很高的階數(shù)，這會(huì)造成插值點(diǎn)不準(zhǔn)確（同回歸中預(yù)測(cè)不準(zhǔn)確）。因此，這個(gè)時(shí)候可以采用區(qū)間插值的方法，即每?jī)蓚€(gè)或三個(gè)相鄰點(diǎn)之間，構(gòu)建具有不同參數(shù)的插值多項(xiàng)式（階數(shù)一致）。

線性插值就是每?jī)蓚€(gè)相鄰的點(diǎn)之間構(gòu)建各自的線性函數(shù)進(jìn)行插值。此時(shí)，第一個(gè)點(diǎn) $(x_1,y_1)$ 與第二個(gè)點(diǎn) $(x_2,y_2)$ 之間的插值函數(shù)(拉格朗日插值公式)為： $f_1(x)=\frac{x-x_2}{x_1-x_2}y_1+\frac{x-x_1}{x_2-x_1}y_2$ ，同樣地，很容易得到第 $i$ 個(gè)點(diǎn) $(x_i,y_i)$ 與第 $i+1$ 個(gè)點(diǎn)之間的插值函數(shù)為： $f_i(x)=\frac{x-x_{i+1}}{x_i-x_{i+1}}y_i+\frac{x-x_i}{x_{i+1}-x_i}y_{i+1},i=1,...,n-1$

image.png

可以看到，線性插值是一種連續(xù)插值，因?yàn)橄噜彽膬蓚€(gè)多項(xiàng)式在同一節(jié)點(diǎn)的值相同，但在該節(jié)點(diǎn)處的斜率（導(dǎo)數(shù)）是不連續(xù)的。
MATLAB實(shí)現(xiàn)：

function Yint = LinearSpline(x, y, Xint)
% LinearSpline calculates interpolation using linear splines.
% Input variables:
% x    A vector with the coordinates x of the data points.
% y    A vector with the coordinates y of the data points.
% Xint The x coordinate of the interpolated point. 
% Output variable:
% Yint The y value of the interpolated point.

n = length(x);
for i = 2:n
    if Xint < x(i)
        break
    end
end
Yint = (Xint - x(i))*y(i-1)/(x(i-1)-x(i)) + (Xint - x(i-1))*y(i)/(x(i)-x(i-1));

5.2 二次樣條插值

給定 $n$ 個(gè)數(shù)據(jù)點(diǎn)，二次樣條插值是在 $n-1$ 個(gè)區(qū)間內(nèi)各自構(gòu)建一個(gè)二次曲線插值數(shù)據(jù)： $f_i(x)=a_ix^2+b_ix+c_i,i=1,2,...,n-1$ 可以看到每個(gè)多項(xiàng)式有 $3$ 個(gè)待定系數(shù)，共有 $3n-3$ 個(gè)待定系數(shù)，因此求解規(guī)則包括：

每個(gè)二次曲線在區(qū)間兩端點(diǎn) $x_i$ 與 $x_{i+1}$ 的值為 $y_i,y_{i+1}$ ： $f_i(x_i)=a_ix^2_i+b_ix_i+c_i=y_i,i=1,2,...,n-1$ $f_{i}(x_{i+1})=a_{i}x^2_{i+1}+b_ix_{i+1}+c_{i}=y_{i+1},i=1,2,...,n-1$ 共構(gòu)建出 $2n-2$ 個(gè)方程
每?jī)蓚€(gè)相鄰區(qū)間的節(jié)點(diǎn)處，對(duì)應(yīng)多項(xiàng)式的斜率相同（保證導(dǎo)數(shù)連續(xù)）： $f'_i(x_{i+1})=2a_ix_i+b_i=f'_{i+1}(x_{i+1})=2a_{i+1}x_{i+1}+b_{i+1},i=1,2,...,n-2$ 共構(gòu)建 $n-2$ 個(gè)方程
上述條件共構(gòu)建了 $3n-4$ 個(gè)方程，想要求解 $3n-3$ 個(gè)參數(shù)，還少一個(gè)條件，因此，通常假設(shè)第一個(gè)點(diǎn)或者最后一個(gè)點(diǎn)的二階導(dǎo)數(shù)為零： $f''_1(x_1)=2a_1=0$ 即第一個(gè)或最后一個(gè)多項(xiàng)式為線性多項(xiàng)式

image.png

5.3 三次樣條插值

對(duì) $n$ 個(gè)數(shù)據(jù)點(diǎn)，三次樣條插值是在 $n-1$ 個(gè)區(qū)間內(nèi)各自構(gòu)建一個(gè)三次多項(xiàng)式進(jìn)行插值。三次樣條插值包含標(biāo)準(zhǔn)形式與拉格朗日形式?？紤]如下標(biāo)準(zhǔn)形式： $f_i(x)=a_ix^3+b_ix^2+c_ix+d_i,i=1,2,...,n-1$ 共有 $4n-4$ 個(gè)待定參數(shù)。求解規(guī)則：

每個(gè)三次曲線在區(qū)間兩端點(diǎn) $x_i$ 與 $x_{i+1}$ 的值為 $y_i,y_{i+1}$ ： $f_i(x_i)=a_ix^3_i+b_ix^2_i+c_ix_i+d_i=y_i,i=1,2,...,n-1$ $f_{i}(x_{i+1})=a_{i}x^3_{i+1}+b_ix^2_{i+1}+c_{i}x_{i+1}+d_i=y_{i+1},i=1,2,...,n-1$ 共構(gòu)建出 $2n-2$ 個(gè)方程
每?jī)蓚€(gè)相鄰區(qū)間的節(jié)點(diǎn)處，對(duì)應(yīng)多項(xiàng)式的斜率相同（保證一階導(dǎo)數(shù)連續(xù)）： $f'_i(x_{i+1})=3a_ix^2_i+2b_ix_i+c_i=f'_{i+1}(x_{i+1})=3a_{i+1}x^2_{i+1}+2b_{i+1}x_{i+1}+c_{i+1},i=1,2,...,n-2$ 共構(gòu)建 $n-2$ 個(gè)方程
每?jī)蓚€(gè)相鄰區(qū)間的節(jié)點(diǎn)處，對(duì)應(yīng)多項(xiàng)式的二階導(dǎo)數(shù)相同（當(dāng)從一條曲線轉(zhuǎn)換為下一條曲線時(shí)，斜率的變化是連續(xù)的）： $f''_i(x_{i+1})=6a_ix_i+2b_i=f''_{i+1}(x_{i+1})=6a_{i+1}x_{i+1}+2b_{i+1},i=1,2,...,n-2$ 共構(gòu)建 $n-2$ 個(gè)方程
上述條件共構(gòu)建了 $4n-6$ 個(gè)方程，還差兩個(gè)條件，通常要求第一個(gè)點(diǎn)和最后一個(gè)點(diǎn)的二階導(dǎo)數(shù)為零(自然三次樣條插值，natural cubic splines)： $6a_1x_1+2b_1=0$ $6a_{n-1}x_n+2b_{n-1}=0$

image.png

三次樣條插值的標(biāo)準(zhǔn)形式易于理解，使用簡(jiǎn)單，但需要求解個(gè)線性方程，計(jì)算復(fù)雜?？紤]三次樣條的拉格朗日形式：

從三次多項(xiàng)式的二階導(dǎo)數(shù)出發(fā)，三次多項(xiàng)式的二階導(dǎo)數(shù) $f''_i(x)$ 是一個(gè)線性函數(shù)，采用拉格朗日多項(xiàng)式插值形式（點(diǎn) $(x_i,f''_i(x_i)),(x_{i+1},f''_i(x_{i+1}))$ ）： $f''_i(x)=\frac{x-x_{i+1}}{x_i-x_{i+1}}f''_i(x_i)+\frac{x-x_{i}}{x_{i+1}-x_i}f''_i(x_{i+1})$ 積分可得： $f'_i(x)=\frac{(x-x_{i+1})^2}{2(x_i-x_{i+1})}f''_i(x_i)+\frac{(x-x_i)^2}{2(x_{i+1}-x_{i})}f''_i(x_{i+1})+C_1$ 再次積分，有 $f_i(x)=\frac{(x-x_{i+1})^3}{6(x_i-x_{i+1})}f''_i(x_i)+\frac{(x-x_{i})^3}{6(x_{i+1}-x_{i})}f''_i(x_{i+1})+C_1x+C_2$ 令 $f_i(x_i)=y_i,f_i(x_{i+1})=y_{i+1}$ 解得
$C_1=\frac{y_{i+1}-y_i}{x_{i+1}-x_i}+\frac{x_{i+1}-x_{i}}{6}f''_i(x_i)-\frac{x_{i+1}-x_{i}}{6}f''_i(x_{i+1})$
$C_2=\frac{y_ix_{i+1}-y_{i+1}x_i}{x_{i+1}-x_i}-\frac{(x_{i+1}-x_i)}{6}f''_i(x_i)x_{i+1}+\frac{(x_{i+1}-x_i)}{6}f''_i(x_{i+1})x_i$ 因此 $\begin{aligned}f_i(x)&=\frac{(x_{i+1}-x)^3}{6(x_{i+1}-x_i)}f''_i(x_i)+\frac{(x-x_{i})^3}{6(x_{i+1}-x_{i})}f''_i(x_{i+1})+C_1x+C_2\\ &=\frac{(x_{i+1}-x)^3}{6(x_{i+1}-x_i)}f''_i(x_i)+\frac{(x-x_{i})^3}{6(x_{i+1}-x_{i})}f''_i(x_{i+1})\\ &+[\frac{y_i}{x_{i+1}-x_i}-\frac{(x_{i+1}-x_i)f''_i(x_i)}{6}](x_{i+1}-x)\\ & +[\frac{y_{i+1}}{x_{i+1}-x_i}-\frac{(x_{i+1}-x_i)f''_i(x_{i+1})}{6}](x-x_{i})\end{aligned},i=1,2,...,n-1$ 定義 $h_i=x_{i+1}-x_i,a_i=f''(x_i)$ ，則 $f_i(x)=\frac{a_i}{6h_i}(x_{i+1}-x)^3+\frac{a_{i+1}}{6h_i}(x-x_i)^3+[\frac{y_i}{h_i}-\frac{h_ia_i}{6}](x_{i+1}-x)+[\frac{y_{i+1}}{h_i}-\frac{h_ia_{i+1}}{6}](x-x_{i})$ 進(jìn)一步利用 $f'_i(x_{i+1})=f'_{i+1}(x_{i+1}),i=1,...,n-2$ ，可得： $6[\frac{y_{i+2}-y_{i+1}}{h_{i+1}}-\frac{y_{i+1}-y_i}{h_i}]=h_ia_i+2(h_i+h_{i+1})a_{i+1}+h_{i+1}a_{i+2}$ 注意 $a_1,a_2,...,a_{n-1},a_n$ 是待定數(shù)，然而這里只有 $n-2$ 個(gè)方程，還需要兩個(gè)條件，利用標(biāo)準(zhǔn)形式中的最后一個(gè)條件，即第一個(gè)點(diǎn)和最后一個(gè)點(diǎn)的二階導(dǎo)數(shù)為零，可得 $a_1=0,a_n=0$ 。

MATLAB實(shí)現(xiàn)：

function Yint = CubicSplines(x,y,xint)
n = length(x);
a1=0;
an=0;
A = zeros(n-2,n-2);
b = zeros(n-2,1);
for i=1:n-2
    h(i) = x(i+1)-x(i);
    h(i+1) = x(i+2)-x(i+1);
    if i == 1
       A(i,i) = 2*(h(i)+h(i+1));
       A(i,i+1) =  h(i+1);
       b(i,1) = 6*((y(i+2)-y(i+1))/h(i+1)-(y(i+1)-y(i))/h(i));
    elseif i == n-2
        A(i,i-1) = h(i);
        A(i,i) = 2*(h(i)+h(i+1));
        b(i,1) = 6*((y(i+2)-y(i+1))/h(i+1)-(y(i+1)-y(i))/h(i));
    else
        A(i,i-1) = h(i); 
        A(i,i) = 2*(h(i)+h(i+1));
        A(i,i+1) = h(i+1);
        b(i,1) = 6*((y(i+2)-y(i+1))/h(i+1)-(y(i+1)-y(i))/h(i));
    end
end
 p = A\b;
 a(1) = 0;
 a(n) = 0;
 for i=1:n-2
    a(1+i) = p(i);
 end
 
 for i = 2:n
    if xint <= x(i)
        break
    end
 end
  h(i) = x(i)-x(i-1);
  Yint = a(i-1)/(6*h(i))*(x(i)-xint)^3+a(i)/(6*h(i))*(xint-x(i-1))^3+(y(i-1)/h(i)-a(i-1)*h(i)/6)*(x(i)-xint)+(y(i)/h(i)-a(i)*h(i)/6)*(xint-x(i-1));
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
x = [5 7.5 10 12.5 15 17.5 20 22.5 25 27.5 30];
h = [3.3 7.5 41.8 51.8 61 101.1 132.9 145.5 171.4 225.8 260.9];


xint = 5:0.5:30;
l = length(xint);
yint = zeros(l,1);
for i = 1:l
    yint(i) = CubicSplines(x,h,xint(i));
end
figure
subplot(1,2,1)
plot(x,h,'k*',xint,yint,'r-')
title('cubic')

p = interp1(x,h,xint,'spline');
subplot(1,2,2)
plot(x,h,'k*',xint,p,'r-')
title('spline')

image.png

6. 總結(jié)

本節(jié)課主要介紹了曲線的擬合與插值，與機(jī)器學(xué)習(xí)中的回歸算法相似。給定一組數(shù)據(jù)集，擬合是找到最好的曲線來(lái)匹配數(shù)據(jù)集中的點(diǎn)（線性回歸，非線性回歸），插值則是找到通過(guò)所有數(shù)據(jù)點(diǎn)的曲線。線性回歸是最簡(jiǎn)單的曲線擬合算法，非線性回歸可以通過(guò)特征變換轉(zhuǎn)化為線性回歸，而其中的高階多項(xiàng)式擬合還可以進(jìn)一步擴(kuò)展到更一般的情形，求解方法均是基于最優(yōu)點(diǎn)處的導(dǎo)數(shù)為零得到。對(duì)于插值，任意 $n$ 個(gè)點(diǎn)都可以找到一個(gè) $n-1$ 階的多項(xiàng)式函數(shù)通過(guò)所有數(shù)據(jù)點(diǎn)，求解方法有拉格朗日法和牛頓法，拉格朗日法計(jì)算簡(jiǎn)單，但增加數(shù)據(jù)點(diǎn)時(shí)，所有參數(shù)均需重新計(jì)算，而牛頓法則不需要重新計(jì)算，但計(jì)算較為復(fù)雜。另一方面，當(dāng)數(shù)據(jù)點(diǎn)比較少時(shí)，多項(xiàng)式插值較為簡(jiǎn)單，但對(duì)于數(shù)據(jù)點(diǎn)較多的情況，多項(xiàng)式的階數(shù)需要很高，計(jì)算復(fù)雜，并且插值時(shí)并不可靠（過(guò)擬合）。樣條插值是在每個(gè)區(qū)間內(nèi)構(gòu)建各自的多項(xiàng)式插值函數(shù)，一般有線性樣條插值、二次樣條插值以及三次樣條插值。其中為求解簡(jiǎn)便，三次樣條插值除了標(biāo)準(zhǔn)形式外，還有一個(gè)便于求解的拉格朗日形式。

色偷偷精品伊人,欧洲久久精品,欧美综合婷婷骚逼,国产AV主播,国产最新探花在线,九色在线视频一区,伊人大交九欧美,1769亚洲,黄色成人av

數(shù)值計(jì)算day6-曲線擬合與插值

數(shù)值計(jì)算day6-曲線擬合與插值

1. 擬合與插值（Fitting, Interpolation）

2. 線性曲線擬合

3. 非線性曲線擬合

3.1 常見(jiàn)的幾種非線性方程

3.2 二次及高階多項(xiàng)式擬合

4. 多項(xiàng)式插值

4.1 多項(xiàng)式插值的拉格朗日形式：

4.2 牛頓多項(xiàng)式插值

5. 樣條插值

5.1 線性樣條插值

5.2 二次樣條插值

5.3 三次樣條插值

6. 總結(jié)

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容

色偷偷精品伊人,欧洲久久精品,欧美综合婷婷骚逼,国产AV主播,国产最新探花在线,九色在线视频一区,伊人大交九 欧美,1769亚洲,黄色成人av

數(shù)值計(jì)算day6-曲線擬合與插值

1. 擬合與插值（Fitting, Interpolation）

2. 線性曲線擬合

3. 非線性曲線擬合

3.1 常見(jiàn)的幾種非線性方程

3.2 二次及高階多項(xiàng)式擬合

4. 多項(xiàng)式插值

4.1 多項(xiàng)式插值的拉格朗日形式：

4.2 牛頓多項(xiàng)式插值

5. 樣條插值

5.1 線性樣條插值

5.2 二次樣條插值

5.3 三次樣條插值

6. 總結(jié)

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容

色偷偷精品伊人,欧洲久久精品,欧美综合婷婷骚逼,国产AV主播,国产最新探花在线,九色在线视频一区,伊人大交九欧美,1769亚洲,黄色成人av