|
|
EDA365欢迎您登录!
您需要 登录 才可以下载或查看,没有帐号?注册
x
想通过matlab拟合得到二元函数y=f(x1,x2)的表达式,但是这里面18组初始的b值是参考别人的模板,因为是用另外的数据生成的,请问怎么用程序实现我这组数据里面的18组初始b值呢?比如用nlinfit,或者cftool,请大佬多多指教! 之前用cftool拟合后得到的poly23多项式函数的函数值域与已知数据相差几个数量级,函数表达式明显是错误的( r7 c/ G0 [1 Q5 D
代码如下:4 {6 @ q6 F" _- ]/ G
clear,clc/ ^# z0 N/ q4 K4 \/ k7 g2 _
x1=[1.2:.1:1.7]';$ c; O6 y6 G i1 g
x2=[1000 5000 7500 10000 15000];
5 Z l- q* R, }1 B5 Z' {6 h& _y=[14964 13166 14235 15550 19200& C! b* U( R8 U' m" m
13479 13479 13090 14235 17014
. D5 L# ]$ E7 z+ Y3 O3 O 13747 13747 12750 13139 149913 h: D2 |# S3 z# e' q
13917 13917 13019 12553 139670 P6 D% F# ]: h. d) j( P/ B# A& U$ z
14065 14065 13386 12821 13188( S* q: L; o! E, f1 S' S$ k4 `
14306 14306 13578 13040 12728];
4 ` G0 K8 z" y/ b. n! W! S0 yn1=length(x1);n2=length(x2);& E* f; R x6 t4 w+ }$ w( V
x1=x1*ones(1,n2);x1=x1( ;
( r- T" M P1 sx2=ones(n1,1)*x2;x2=x2( ;4 ~/ q) D" G9 _% Y% P, k7 j1 W0 t
y=y( ; X=[x1,x2]; n=length(y);
' a8 u' }. b6 A" l1 P; Tstr=num2str([1:n]');) Z) ^, l5 e& D8 r. U
fx1=@(b,x1,x2)(b(1)+b(2)*x2+b(3)*x2.^2+b(4)*x1.*x2+b(5)*x1.^3+b(6)*x2.^3+b(7)*x1.^4+b(8)*x2.^4+b(9)*x1.*x2.^3+b(10)*x1.^2.*x2.^2+b(11)*x1.^5+b(12)*x2.^5)./(1+b(13)*exp(b(14)*x1+b(15)*x2+b(16)*x1.^2+b(17)*x2.^2+b(18)*x1.*x2));5 [* G' Y6 H* T) z2 I [+ n
fx2=@(b,X)(b(1)+b(2)*X(:,2)+b(3)*X(:,2).^2+b(4)*X(:,1).*X(:,2)+b(5)*X(:,1).^3+b(6)*X(:,2).^3+b(7)*X(:,1).^4+b(8)*X(:,2).^4+b(9)*X(:,1).*X(:,2).^3+b(10)*X(:,1).^2.*X(:,2).^2+b(11)*X(:,1).^5+b(12)*X(:,2).^5)./(1+b(13)*exp(b(14)*X(:,1)+b(15)*X(:,2)+b(16)*X(:,1).^2+b(17)*X(:,2).^2+b(18)*X(:,1).*X(:,2)));- b- H. U1 C+ m
b=[345.45 -191.53 41.80 -2.6793 83.3156 -4.4545 -156.1311 0.2333 -0.0086582 0.251804 81.565 -0.0048047 989422 -146.607 4.003 95.35213571 -0.25589 2.752];5 E4 y7 v8 Z7 f9 I8 K2 W R
for l=1:5
- H1 i6 Y" i) K4 }: v: D: O b=lsqcurvefit(fx2,b,X,y);
5 V# g9 @4 \& i8 B( `8 \& o b=nlinfit(X,y,fx2,b);
$ y- i' g# Q* a5 A- A! Wend' N$ o8 _) o% e7 r4 W' A: W
b
9 ?1 I0 N# m* Y8 K; _1 Mfigure(1),clf2 X0 j& o, s" t
plot3(x1,x2,y,'o')
8 F9 `2 l1 `; n2 fstem3(x1,x2,y,'filled')
+ K. Q& }! @ W+ y4 @5 n6 u6 ?text(x1,x2,y+.01,str)
# g% D; Y' z4 Z7 _& z7 Z- K. ^hold on
. I& C; M2 I1 u; Q[x11,x22]=meshgrid(min(x1):range(x1)/80:max(x1),min(x2):range(x2)/80:max(x2));
V9 ~! g$ i( E7 ^yhat=fx1(b,x11,x22);
' E5 B0 r, v* B$ p; a8 GsuRF(x11,x22,yhat)# I3 i3 q& B# M: d/ P) K# d
shading interp
, v8 _/ q0 O) U8 Y7 ealpha(.8)
! O4 n. K4 c: B) T$ A# laxis tight
6 p! C; V$ _6 t! |+ Q- `xlabel('α');ylabel('C'),zlabel('V1')" N/ m/ y X j& v* _ u7 v8 j
SSy=var(y)*(n-1)
/ Q! E( }& B* J- n' ?y1=fx1(b,x1,x2);
4 f& o1 t7 J6 O: l5 VRSS=(y-y1)'*(y-y1)
. J6 Z6 ^2 b8 |2 Drsquare=(SSy-RSS)/SSy- j5 D. [6 K; d
MSe=RSS/(n-length(b))
1 `. e9 o0 l; i9 M! \/ }: B& w: J
|
|