關於matlab作圖,不知如何具體解釋,langrange基函數,求大神教?

稱為一組langrange基,其中

,請用matlab畫出langrange基函數圖


這是一個簡單的Lagrange插值問題,一般就是根據給出的n+1個離散點(x的坐標互異),通過一個n階多項式進行擬合,理論上上一定能找到這樣的一個n階多項式,使得該多項式構造的函數經過所有的樣本點,Lagrange就是這類多項式的一個具體實現。

舉個簡單例子,假設給出這樣一組數據(樣本點){(x_0, y_0), (x_1, y_1), cdots, (x_n, y_n)},其中forall ~i 
eq j,~x_i
eq x_j,那麼一定可以找到這樣一個多項式f(x) = a_0 + a_1x + cdots + a_nx^n使得forall ~i = 0, 1, cdots, n都有f(x_i) = y_i,即該多項式經過所有給出的樣本點。

為什麼這個命題就能成立呢,稍微根據一點線性代數的知識就能推論出來。我們假設這樣一個多項式f(x)存在,那麼其顯然要滿足如下方程:

 left[
 egin{matrix}
   1  x_0  x_0^2  cdots  x_0^n \
   1  x_1  x_1^2  cdots  x_1^n \
   1  vdots  vdots  cdots  vdots\
   1  x_n  x_n^2  cdots  x_n^n \
  end{matrix}
  
ight]
 left(
 egin{matrix}
    a_0 \
    a_1 \
    vdots\
    a_n\
 end{matrix}
 
ight) 
= left(
 egin{matrix}
    y_0 \
    y_1 \
    vdots\
    y_n\
 end{matrix}
 
ight)

我們知道左側矩陣為一個范德蒙德方陣,而x_i彼此互異,該行列式不為零,方程存在唯一解,於是存在性得到證明。

那麼剩下的問題就是如何找出這樣一個多項式(不通過粗暴的解方程實現),而Lagrange插值就是簡單構造出這樣的一個解,其具體構造如下:

L(x) = Sigma_{i=1}^n y_iell_i(x),其中ell_i(x)就是你所說的單位(基)函數,其表達式如下:

ell_j(x) = Pi_{i=0, i 
eq j}^nfrac{x-x_i}{x_j-x_i},顯然有ell_j(x_j)=1ell_j(x_i) = 0, ~forall ~i 
eq j

Lagrange函數就是這些基函數(n階多項式的線性組合),很容易驗證其通過所有樣本點。

最後給出一個具體的算例,假設有四個離散點(-9, 5), (-4, 2), (-1, -2), (7, 9),那麼簡單通過polyfit就能給出其3階多項式的擬合結果(也就是Lagrange插值結果),如下圖中的黑色虛線所示:

而其基函數(乘以y_i)則為其它彩色實線,參考代碼:

%%
sample_x = [-9, -4, -1, 7];
sample_y = [5, 2, -2, 9];
scatter(sample_x, sample_y, "SizeData", 60, "MarkerFaceColor", "b", "MarkerFaceAlpha", .75)
hold on

poly_func = polyfit(sample_x, sample_y, length(sample_x)-1);
predict_x = linspace(min(sample_x)-1, max(sample_x)+1, 100);
predict_y = polyval(poly_func, predict_x);
plot(predict_x, predict_y, "k--", "LineWidth", 2.0)
labels = {"$(x, y)
, "$Sigma ell_iy_i};

syms x
indices = 1:length(sample_x);
for i = indices
tmp_index = find(indices~=i);
lagrange_i = prod(x-sample_x(tmp_index));
lagrange_i = lagrange_i / prod(sample_x(i)-sample_x(tmp_index)) * sample_y(i);
L_i = sym2poly(lagrange_i);
plot(predict_x, polyval(L_i, predict_x), "LineWidth", 1.5)
labels{i+2} = sprintf("$y_%d\ell_%d, i, i");
end

lg = legend(labels, "Location", "northwest");
lg.set("interpreter", "latex", "FontSize", 12)
xlabel("$x, "interp", "latex", "FontSize", 13)
ylabel("$y, "interp", "latex", "FontSize", 13)
title("Demonstration of Lagrange Interpolation", "FontName", "Arial", "FontSize", 14)
grid on

%%
saveas(gcf, "lagrange_interp.png")


yellow說了如何使用符號計算獲得基函數係數,這裡補充一下用數值計算獲得基函數的方法。主要利用poly求出多項式係數,之後用deconv做多項式除法獲得每個基函數的多項式係數:

x = [-9, -4, -1, 7];
y = [ 5, 2, -2, 9];

p = poly(x)";
L = bsxfun(@rdivide,... % 基函數係數
bsxfun(@(~,t)deconv(p,[1 -t]),x",x),... % 基函數分母係數
prod(bsxfun(@minus,x,x")+eye(numel(x)),1)); % 基函數分子
P = L*y"; % 插值函數係數

xi = linspace(x(1)-1,x(end)+1)";
plot(x,y,"ro", .... % 數據點
xi,polyval(P,xi),"--k",... % 插值函數
xi,bsxfun(@power,xi,numel(x)-1:-1:0)*L,... % 基函數
"Linewidth",1);

set([xlabel("$x$") ylabel("$y$")],"Interp","latex","fontsi",15), grid on
legend(["$(x,y)$","$f(x)=sum y_iell_i(x)$",...
sprintfc("$\ell_%d(x)$",1:numel(x))],...
"Interpreter","latex","Location","best","fontsize",12)


我換一個通俗的語言來給你解釋一下

是這樣的,你現在知道一堆離散的點,它們也許像這樣。。

你現在想通過構造一個函數,讓這個函數通過這些點。。

於是這是一種可能

對吧 這是一種可能~

不過是不是只有一種可能呢。。

它還可能像是藍線這樣。。

於是,拉格朗日提出了一種他自己定義的曲線,他的定義就是

通過這n個已知點的n-1次函數成為拉格朗日函數。而且這個n-1次函數是唯一的

還有牛頓插值,插出來的函數是跟拉格朗日一樣的。不同點是,拉格朗日修改數據容易,牛頓添加數據容易。

拉格朗日函數可以通過這個基函數構造,這些基函數的表達式的一個特點就是特別齊整,方便編程,至於你看沒看懂樓上那位同學說的話我就不知道了,要是看不懂可以找點數值分析的書來看,解釋得相對簡單。不過你要是不打算明白,那你只需要知道這只是構造一個通過這些數據點的函數的其中一種方法即可,照套公式就是,反正公式定理就是拿來套的。


wiki上的圖,File:Lagrange polynomial.svgMathematica版

pts={{-9,5},{-4,2},{-1,-2},{7,9}};
{X,Y}=Transpose@pts;
func=Table[Product[If[j!=i,(x-i)/(j-i),1],{i,X}],{j,X}]Y

Plot[{func,Total@func},{x,-10,10},
GridLines-&>Range[-12,12,{2,2}],
Epilog-&>{PointSize@.015,Point@pts},
PlotStyle-&>Append[Automatic{1,1,1,1},{Dashed,Black}],
PlotLegends-&>Table[Style[If[k&<5,Subscript[[ScriptL], k][x],f[x]==[CapitalSigma] Subscript[y, i] Subscript[[ScriptL], i][x]],FontFamily-&>"Times"],{k,5}]]


推薦閱讀:

擬合與插值的區別?
有什麼比較好的方法可以快速準確的開出七次方嗎?

TAG:數學 | 拉格朗日J-LLagrange | MATLAB | 數值分析 |