%最小二乘法插值
%用法:
%　　1.將該文件放在work文件夾下
%    2.在matlab輸入命令zuixiao(n); (n為逼近精度，n>0整數)

function zuixiao(N) 
format long;

X=[3 4 5 6 7 8 9];                       %插值數據x坐標,x內有k個元素,N的最大值為k
F=[2.01 2.98 3.50 5.02 5.47 6.02 7.05];  %插值數據y坐標

syms x;
% Aa=B;
A=zeros(N);
B=zeros(N,1);
a=zeros(N,1);
for i=1:1:N %構成逼近基函數1  x  x^2...
    R{i}=x^(i-1); 
end
for k=1:1:N %構成逼近基函數1  x  x^2...方程系數矩陣A
    for j=1:1:N
        for i=1:1:length(X)
            x=X(i);  
            A(k,j)=A(k,j)+eval(R{k}*R{j});
        end
    end
end
for k=1:1:N %構成系數矩陣B
        for i=1:1:length(X)
            x=X(i);           
            B(k)=B(k)+F(i)*eval(R{k});
        end
end
a=inv(A)*B;%求解逼近基函數1  x  x^2...的系數a;
syms x;
fun=0;
for k=1:1:N%獲得最小二乘法插值函數fun
    fun=fun+a(k)*R{k};
end
fun                       %輸出獲得最小二乘法插值函數fun

plot(X,F,'o');%畫出插值節點
hold on

%畫出最小二乘法插值函數fun
i=linspace(X(1),X(length(X)),100);
for k=1:1:length(i)
    x=i(k);
    y_fun(k)=eval(fun);
end
plot(i,y_fun);
hold on