数据分析MATLAB相关算法

数据分析-MATLAB相关算法

一.非线性方程求解

1.二分法

m文件：

function[c/err/yc]=bisect(f,a,b,delta)

%f是所要求解的函数

%a,b是区间的左右限

%delta是误差界

%c是近似解

%yc是函数f在c上的值，err为误差

ifnargin<4

delta=le-5;

end

ya=feval('f',a);

yb=feval('f'/b);

ifyb==O

c=b,return/end

ifya*yb>0,

disp('(a,b)不是有根区间')；

return,end

maxl=l+round((log(b-a)-log(delta))/log(2));

fork=l:maxl

c=(a+b)/2;

yc=feval('f',c);

ifyc==O

a=c;

b=c;

break,

elseifyb*yc>0

b=c;

yb=yc;

else

a=c;

ya=c;

end

if(b-a)<delta

break,

end

,,

k/c=(a+b)/2,err=abs(b-a)/yc=feval(f/c);

end

2.不动点迭代

m文件

function[pO/k,err,p]=fixpt(g,pO,tol.maxl)

%g是给定的迭代函数

%pO是初值

%p是不动点的近似值

P(l)=pO；

fork=2:maxl

P(k)=feval('g;P(k-l));

k/err=abs(P(k)-P(k-l))

p=P(k);

if(err<tol)

break;end

ifk==maxl

disp('超过了迭代最大次数；

end

end

P

选定的迭代函数：

functiony=g(x)

y=sin(x)/x;

二.方阵的特征值与特征向量

1幕法

function[mxbiaozhi]=mifa(A,jingdu,cishu)

%累法求矩阵最大特征值，其中

%m为绝对值最大的特征值，x为对应最大特征值的特征向量

%biaozhi表明迭代是否成功

ifnargin<3

cishu=100;

end

ifnargin<2

jingdu=le-5;

end

n=length(A);

x=ones(n,l);

biaozhi='迭代失败！'；

k=0;ml=0;

whilek<=cishu

v=A*x;

[vmax,k]=max(abs(v));

m=v(k);

x=v/m;

ifabs(m-ml)<jingdu

biaozhi='迭代成功！';break;

end

ml=m;

k=k+l;

end

三.常微分方程求解

L改进欧拉法解微分方程

functionyout=gaijinoula(f,xO,yO,xn,n)

%定义输入输出

x=zeros(l,n+l);y=zeros(l,n+l);

x(l)=xO;y(l)=yO;

h=(xn-xO)/n;

fori=l:n

x(i+l)=x(i)+h;

zO=y(i)+h*feval(f,x(i)/y(i));

y(i+l)=y(i)+(feval(f,x(i),y(i))+feval(f,x(i+l),z0))*h/2;

end

'

shuchu=[x'/y]

fprintfCx(i)y(i)1)

举例：

所要求解的函数:

functionDy=f(x,y)

Dy=x+y;

2.四阶龙格库塔法解微分方程

functionyout=xin(bianliang)

%定义输入输出

clearall

xO=O;xn=l;yO=l;h=O.l;

%设置初始值、区间和步长

[y,x]=lgkt4j(xO,xn,yO,h);

%四阶龙格库塔法

n=length(x);

fprintfCix(i)y(i)\n');

%输出格式

fori=l:n

fprintf('%2d%3.3f%4.4f\n',i,x(i),y(i));

end

function[y,x]=lgkt4j(x0,xn,y0,h)

x=xO:h:xn;

%设置区间

n=length(x);

yl=x;

yl⑴=yO;

fori=l:n

Kl=f(x(i),yl(i));

K2=f(x(i)+h/2,yl(i)+h/2*Kl);

K3=f(x(i)+h/2,yl(i)+h/2*K2);

K4=f(x(i)+h,yl(i)+h*K3);

yl(i+l)=yl(i)+h/6*(Kl+2*K2+2*K3+K4);

end

y=yi；

举例：

functionDy=f(x,y)

Dy=y-2*x/y;

执行结果为：

ix(i)y(i)

10.0001.0000

20.1001.0954

30.2001.1832

40.3001.2649

50.4001.3416

60.5001.4142

70.6001.4832

80.7001.5492

90.8001.6125

100.9001.6733

111.0001.7321

四.插值法

1拉格朗日插值

function[c,l]=lglr(x,y)

%乂为门个节点的横坐标组成的向量，y为纵坐标组成的向量

%c为插值函数的系数组成的向量

%输出为差值多项式的系数

w=length(x);

n=w-l;

l=zeros(w,w);

fork=l:n+l

v=l;

forj=l:n+l

ifk-=j

v=conv(v,poly(x(j)))/(x(k)-x(j));

end

end

l(k,:)=v;

end

c=y*l;

五.数据拟合

L最小二乘拟合

functionc=zxrc(x,y,m)

%x是数据点横坐标，y数据点纵坐标

%m要构造的多项式的系数，c是多项式由高到低的系数所组成的向量

n=length(x);

b=zeros(l:m+l);

f=zeros(n,m+l);

fork=l:m+l

A

f(:fk)=x'.(k-l);

end

a=f'*f;

b=f'*y';

c=a\b;

c=flipud(c);

六.矩阵相关程序

1.求矩阵的行列式

functiond=hanglieshi(a)

%求任意输入矩阵的行列式

clearall;

a=input('输入矩阵a=');

d=l;

n=size(a);%方阵的行(或者列)数

fork=l:n-l

设矩阵的主元

e=a(kfk);%

fori=k:n%求出矩阵的全主元

forj=k:n

ifabs(a(i,j))>e

e=a(i,j);

P=i；

q=j；

elsec=0;

end

end

end

forj=k:n%行交换

t=a(k,j);

a(kj)=a(p,j);

a(p,j)=t;

end

ifp~=k%判断行列式是否换号

d=d*(-l);

elsed=d;

end

fori=k:n%列交换

t=a(i,k);

a(i,k)=a(i,q);

a(i,q)=t;

end

ifq~=k%判断行列式是否换号

d=d*(-l);

elsed=d;

end

ifa(k,k)-=O

fori=k+l:n%消元

r=a(i,k)/a(k,k);

forj=k+l:n

a(ij)=a(ij)-r*a(kfj);

end

end

elsed=d;

end

end

fori=l:n%求行列式

d=d*a(i,i);

end

disp('矩阵a的行列式为：’)

2.矩阵的换行

functionc=huanhang(a)

%实现矩阵换行

clearall;

a=input('输入矩阵a=');

[m,n]=size(a);

forj=l:n

t=a(l,j);

a(l,j)=a(2fj);

a(2fj)=t;

end

c=a;

disp('换行后矩阵a变为：’)

c

3.列主元消元法解方程

functiond=jiefang(a)

%列主元消元法解方程

clearall;

a二input('输入矩阵a=');

[rowfcolumn]=size(a);

fori=l:column%每一列的列标

m(i)=i;

s(i)=O;

x(i)=O;

end

fork=l:row・l%最后一行不用比较

e=a(kfk);

p=k;

q=k;

fori=k:row

forj=k:column-l

ifabs(a(ifj))>abs(e)

e=a(i,j);

p=i；

q=j；

elsec=0;

end

end

end

t=m(k);%换列标记

m(k)=m(q);

m(q)=t;

fori=l:row%列交换

t=a(i,k);

a(i,k)=a(i,q);

a(Lq)=t；

end

forj=k:column%行变换

t=a(k,j);

a(k,j)=a(p,j);

a(PJ)=t；

end

ifa(k,k)==0%消元

dispf非唯一解')

elsefori=k+l:row

r=a(i,k)/a(k,k);

forj=k:column

a(ij)=a(i,j)-r*a(kj);

end

end

end

end

ifa(row,row)==0

disp('非唯一解')

else

s(row)=a(row,column)/a(row,row);

s(row)

q=m(row);

x(q)=s(row);

fori=row-l:l

forj=i+l:row

s(i)=s(i)+a(i,j)*x(i)；

end

s(i)=[a(i,column)-s(i)]/a(i,i);

q=m(i);

x(q)=s(i);

end

end

fori=l:row

x(i)

end

end

4.两矩阵相乘

functiond=chengfa(A,B)

%实现两个矩阵相乘

clearall;

A=input('输入矩阵A=');

B二input('输入矩阵B=')

[mn]=size(A);

[nbp]=size(B);

C=zeros(m,p);

ifn~=nb

disp('不满足矩阵相乘条件')

elsefori=l:m

forj=l:p

d=0;

fork=l:n

d=d+A(i,k)*B(k,j);

end

C(ifj)=d;

end

end

di$p('矩阵AB结果为一)

C

End

5.矩阵元素最大值及下标

functiond=xunzhuyuan(a)

%求一个矩阵的最大元素及其下标

clearall