㈠ MATLAB代碼對一組數據用最小二乘法處理。急急急~
7.1.1 分段線性插值
所謂分段線性插值就是通過插值點用折線段連接起來逼近原曲線,這也是計算機繪制圖形的基本原理。實現分段線性插值不需編制函數程序,MATLAB自身提供了內部函數interp1其主要用法如下:
interp1(x,y,xi) 一維插值
◆ yi=interp1(x,y,xi)
對一組點(x,y) 進行插值,計算插值點xi的函數值。x為節點向量值,y為對應的節點函數值。如果y 為矩陣,則插值對y 的每一列進行,若y 的維數超出x 或 xi 的維數,則返回NaN。
◆ yi=interp1(y,xi)
此格式默認x=1:n ,n為向量y的元素個數值,或等於矩陣y的size(y,1)。
◆ yi=interp1(x,y,xi,』method』)
method用來指定插值的演算法。默認為線性演算法。其值常用的可以是如下的字元串。
● nearest 線性最近項插值。
● linear 線性插值。
● spline 三次樣條插值。
● cubic 三次插值。
所有的插值方法要求x是單調的。x 也可能並非連續等距的。
正弦曲線的插值示例:
>> x=0:0.1:10;
>> y=sin(x);
>> xi=0:0.25:10;
>> yi=interp1(x,y,xi);
>> plot(x,y,』0』,xi,yi)
則可以得到相應的插值曲線(讀者可自己上機實驗)。
Matlab也能夠完成二維插值的運算,相應的函數為interp2,使用方法與interpl基本相同,只是輸入和輸出的參數為矩陣,對應於二維平面上的數據點,詳細的用法見Matlab聯機幫助。
7.1.2 最小二乘法擬合
在科學實驗的統計方法研究中,往往要從一組實驗數據中尋找出自變數x 和因變數y之間的函數關系y=f(x) 。由於觀測數據往往不夠准確,因此並不要求y=f(x)經過所有的點 ,而只要求在給定點上誤差按照某種標准達到最小,通常採用歐氏范數作為誤差量度的標准。這就是所謂的最小二乘法。在MATLAB中實現最小二乘法擬合通常採用polyfit函數進行。
函數polyfit是指用一個多項式函數來對已知數據進行擬合,我們以下列數據為例介紹這個函數的用法:
>> x=0:0.1:1;
>> y=[ -0.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2 ]
為了使用polyfit,首先必須指定我們希望以多少階多項式對以上數據進行擬合,如果我們指定一階多項式,結果為線性近似,通常稱為線性回歸。我們選擇二階多項式進行擬合。
>> P= polyfit (x, y, 2)
P=
-9.8108 20.1293 -0.0317
函數返回的是一個多項式系數的行向量,寫成多項式形式為:
為了比較擬合結果,我們繪制兩者的圖形:
>> xi=linspace (0, 1, 100); %繪圖的X-軸數據。
>> Z=polyval (p, xi); %得到多項式在數據點處的值。
當然,我們也可以選擇更高冪次的多項式進行擬合,如10階:
>> p=polyfit (x, y, 10);
>> xi=linspace (0, 1,100);
>> z=ployval (p, xi);
讀者可以上機繪圖進行比較,曲線在數據點附近更加接近數據點的測量值了,但從整體上來說,曲線波動比較大,並不一定適合實際使用的需要,所以在進行高階曲線擬合時,「越高越好」的觀點不一定對的。
7.2 符號工具箱及其應用
在數學應用中,常常需要做極限、微分、求導數等運算,MATLAB稱這些運算為符號運算。MATLAB的符號運算功能是通過調用符號運算工具箱(Symbolic Math Toolbox)內的工具實現,其內核是借用Maple數學軟體的。MATLAB的符號運算工具箱包含了微積分運算、化簡和代換、解方程等幾個方面的工具,其詳細內容可通過MATLAB系統的聯機幫助查閱,本節僅對它的常用功能做簡單介紹。
7.2.1 符號變數與符號表達式
MATLAB符號運算工具箱處理的對象主要是符號變數與符號表達式。要實現其符號運算,首先需要將處理對象定義為符號變數或符號表達式,其定義格式如下:
格式1: sym (『變數名』) 或 sym (『表達式』)
功能: 定義一個符號變數或符號表達式。
例如:
>> sym (『x』) % 定義變數x為符號變數
>> sym(『x+1』) % 定義表達式x+1為符號表達式
格式2: syms 變數名1 變數名2 …… 變數名n
功能: 定義變數名1、變數2 ……、變數名 n為符號變數。
例如:
>> syms a b x t % 定義a,b, x,t 均為符號變數
7.2.2 微積分運算
1、極限
格式:limit (f, t, a, 『left』 or 『right』)
功能:求符號變數t 趨近a 時,函數f 的(左或右)極限。『left』 表示求左極限,『right』 表示求右極限,省略時表示求一般極限;a省略時變數t 趨近0; t省略時默認變數為x ,若無x則尋找(字母表上)最接近字母x 的變數。
例如:求極限的命令及結果為:
>> syms x t
>> limit ((1+2*t/x)^(3*x) , x, inf )
ans=
exp(6*t)
再如求函數x / |x| ,當時的左極限和右極限,命令及結果為:
>> syms x
>> limit(x/abs(x), x, 0, 』left』) ans = -1
>> limit(x/abs(x),x, 0, 』right』) ans = 1
2、導數
格式: diff (f,t,n)
功能: 求函數f 對變數 t的n 階導數。當n省略時,默認 n=1;當t省略時,默認變數x, 若無x時則查找字母表上最接近字母x 的字母。
例如:求函數f=a*x^2+b*x+c對變數 x的一階導數, 命令及結果為
>> syms a b c x
>> f=a*x^2+b*x+c;
>> diff(f)
ans=
2*a*x+b
求函數f 對變數b的一階導數(可看作求偏導), 命令及結果為
>> diff(f,b) ans=x
求函數f 對變數x的二階導數, 命令及結果為
>> diff(f,2) ans=2*a
3、積分
格式: int(f,t,a,b)
功能: 求函數f 對變數 t從a 到b的定積分. 當a和b省略時求不定積分;當t省略時, 默認變數為(字母表上)最接近字母x的變數。
例如:求函數f=a*x^2+b*x+c對變數x不定積分, 命令及結果為
>> syms a b c x
>> f=a*x^2+b*x+c;
>> int(f)
ans=
1/3*a*x^3+1/2*b*x^2+c*x
求函數f 對變數b不定積分, 命令及結果為
>> int(f,b)
ans=
a*x^2*b+1/2*b^2*x+c*b
求函數f 對變數x 從 1到5的定積分, 命令及結果為
>> int(f,1,5)
ans=
124/3*a+12*b+4*c
4、級數求和
格式: symsum (s,t,a,b)
功能:求表達式s中的符號變數t從第a項到第b項的級數和。
例如: 求級數的前三項的和, 命令及結果為
>> symsum(1/x,1,3) ans=11/6
7.2.3 化簡和代換
MATLAB符號運算工具箱中,包括了較多的代數式化簡和代換功能,下面僅舉出部分常見運算。
simplify 利用各種恆等式化簡代數式
expand 將乘積展開為和式
factor 把多項式轉換為乘積形式
collect 合並同類項
horner 把多項式轉換為嵌套表示形式
例如:進行合並同類項執行
>> syms x
>> collect(3*x^3-0.5*x^3+3*x^2)
ans=
5/2*x^3+3*x^2)
進行因式分解執行
>> factor(3*x^3-0.5*x^3+3*x^2)
ans=
1/2*x^2*(5*x+6)
7.2.4 解方程
1、代數方程
格式:solve (f,t)
功能:對變數t 解方程f=0,t 預設時默認為x 或最接近字母x 的符號變數。
例如:求解一元二次方程f=a*x^2+b*x+c的實根,
>> syms a b c x
>> f=a*x^2+b*x+c;
>> solve (f,x)
ans=
[1/2/a*(-b+(b^2-4*a*c)^ (1/2))]
[1/2/a*(-b-(b^2-4*a*c)^ (1/2))]
2、微分方程
格式:dsolve(『s』, 』s1』, 』s2』,…, 』x』)
其中s為方程;s1,s2,……為初始條件,預設時給出含任意常數c1,c2,……的通解;x為自變數,預設時默認為t 。
例如:求微分方程的通解
>> dsolve(『Dy=1+y^2』)
ans=
tan(t+c1)
7.3 優化工具箱及其應用
在工程設計、經濟管理和科學研究等諸多領域中,人們常常會遇到這樣的問題:如何從一切可能的方案中選擇最好、最優的方案,在數學上把這類問題稱為最優化問題。這類問題很多,例如當設計一個機械零件時如何在保證強度的前提下使重量最輕
或用量最省(當然偷工減料除外);如何確定參數,使其承載能力最高;在安排生產時,如何在現有的人力、設備的條件下,合理安排生產,使其產品的總產值最高;在確定庫存時如何在保證銷售量的前提下,使庫存成本最小;在物資調配時,如何組織運輸使運輸費用最少。這些都屬於最優化問題所研究的對象。
MATLAB的優化工具箱被放在toolbox目錄下的optim子目錄中,其中包括有若干個常用的求解函數最優化問題的程序。MATLAB的優化工具箱也在不斷地完善。不同版本的MATLAB,其工具箱不完全相同。在MATLAB5.3版本中,對優化工具箱作了全面的改進。每個原有的常用程序都重新編制了一個新的程序。除fzero和fsolve外都重新起了名字。這些新程序使用一套新的控制演算法的選項。與原有的程序相比,新程序的功能增強了。在MATLAB5.3和6.0版本中,原有的優化程序(除fzero和fsolve外)仍然保留並且可以使用,但是它們遲早會被撤消的。鑒於上述情況,本書將只介紹那些新的常用的幾個優化程序。
㈡ 關於最小二乘法的c語言程序
#include <stdio.h>
void main ()
{
int num,i;
float x,y,l,m,n,p,a,b;
i=1;
l=0.0;
m=0.0;
n=0.0;
p=0.0;
printf ("請輸入你想計算的x,y的個數:");
scanf("%d",&num);
if (num>=1)
{
while (i<=num);
{
printf("請輸入x的值");
scanf ("%lf",&x);
printf("請輸入y的值");
scanf ("%lf",&y);
l+=x;
m+=y;
n+=x*y;
p+=x*x;
i++;
}
a=(num*n-l*m)/(num*p-l*l);
b=(p*m-n*l)/(num*p-l*l);
printf("最小二乘法所算得的斜率和截距分別為%f和%f\n",a,b);
}
else printf("mun"輸入有誤!);
}
㈢ 求非線性最小二乘法的演算法及代碼.
main ()
{
int n,i,flag2;
char flag1='y';
float ar[50],br[50],x,y,xe,ye,xye,xxe,sx,sy,sxy,sxx,a,b;
printf ("\n歡迎使用最小二乘法數據處理程序\n");
printf ("\n說明:本程序運行結果保留小數點後三位\n");
for (;flag1=='y'||flag1=='Y';)
{printf ("\n請輸入您要處理的數據的組數:");
printf ("\n*****提示:本程序定義一對x,y值為一組數據:");
scanf ("%d",&n);
if (n>50) {printf ("\n對不起,本程序暫時無法處理50組以上的數據");
continue;
}
printf ("\n請選擇您的數據的處理方式:");
printf ("\n\t1.y與x為一次線性關系");
printf ("\n\t2.y與x的2次為線性關系\n");
scanf ("%d",&flag2);
if (flag2>2||flag2<1) {printf ("\n對不起,您的輸入不正確\n");continue;}
for (i=0;i<n;i++)
{printf ("\n請輸入第%2d個x的值\tx%2d=",i+1,i+1);
scanf (" %f",&ar[i]);
printf ("\n請輸入對應的y的值:\ty%2d=",i+1);
scanf (" %f",&br[i]);
}
if (flag2!=1)
{for (i=0;i<n;i++)
br[i]=br[i]/(ar[i]*ar[i]);
}
sx=sy=sxx=sxy=0;
for (i=0;i<n;i++)
{sx=sx+ar[i];
sy=sy+br[i];
sxx=sxx+ar[i]*ar[i];
sxy=sxy+ar[i]*br[i];
}
xe=sx/n;
ye=sy/n;
xye=sxy/n;
xxe=sxx/n;
b=(xye-xe*ye)/(xxe-xe*xe);
a=ye-b*xe;
printf ("\n對您輸入的數據的處理已經完成,結果如下:");
printf ("\n\ta=%8.3f\n\tb=%8.3f\n",a,b);
printf ("\nb即為擬合直線的斜率,a為截距\n");
printf ("\n*********************************************************\n");
printf ("\n是否繼續使用本程序處理數據?(y/n)?");
scanf (" %c",&flag1);
if (flag1=='y'||flag1=='Y') continue;
else if (flag1=='n'||flag1=='N') break;
else {printf ("\n***操作非法,本程序將關閉***\n");
exit (0);
}
}
exit (0);
}