|
|
EDA365欢迎您登录!
您需要 登录 才可以下载或查看,没有帐号?注册
x
想通过matlab拟合得到二元函数y=f(x1,x2)的表达式,但是这里面18组初始的b值是参考别人的模板,因为是用另外的数据生成的,请问怎么用程序实现我这组数据里面的18组初始b值呢?比如用nlinfit,或者cftool,请大佬多多指教! 之前用cftool拟合后得到的poly23多项式函数的函数值域与已知数据相差几个数量级,函数表达式明显是错误的- @* P5 \7 N. _4 s1 }
代码如下:) T5 e" n$ `& N7 ?- q# M
clear,clc
: l, \( `" ~0 j4 n% zx1=[1.2:.1:1.7]';1 O3 ?7 B/ I) s% c$ B% z$ {/ v
x2=[1000 5000 7500 10000 15000];
$ i- D5 l H8 s; p. t0 z9 by=[14964 13166 14235 15550 192007 i; C! q3 E C+ k- o7 O9 P
13479 13479 13090 14235 17014* y% O- |- A/ z$ V4 M
13747 13747 12750 13139 14991
; Z. n5 H* b* F9 Z0 Y3 X# }. O3 F6 P 13917 13917 13019 12553 13967: W7 e' h4 z( x1 x3 K
14065 14065 13386 12821 13188
) V$ ~1 r; [7 X* D6 O1 \ 14306 14306 13578 13040 12728];5 |& m F+ z9 p! T% F) Q
n1=length(x1);n2=length(x2);
& R5 ~. [6 @9 l ]" L% wx1=x1*ones(1,n2);x1=x1( ;% J: t6 i& a! c$ Y( D7 ^
x2=ones(n1,1)*x2;x2=x2( ;
/ e9 K. L8 }6 c t$ E' ly=y( ; X=[x1,x2]; n=length(y);
- M0 [7 s8 m4 [; j ]str=num2str([1:n]');
! L% E+ F) h' ~: w( l8 Rfx1=@(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));7 f- h( j; }% z4 R
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)));' o7 R, i$ U0 o% S$ M% F
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];
+ q9 A' g- V9 M; z. Sfor l=1:5! V( O3 ]* ?$ {3 U
b=lsqcurvefit(fx2,b,X,y);) I% U# c5 A* v/ G# ]! z/ L! M
b=nlinfit(X,y,fx2,b);$ j& X0 s& J0 t2 a d
end
* N$ ]5 m* b' Rb
/ t9 P# p: Q9 [* \" g* O( O. ]figure(1),clf
0 b3 z5 R6 F+ wplot3(x1,x2,y,'o')
7 S3 O' g$ F* Cstem3(x1,x2,y,'filled')
% y6 [$ U# w1 ?- \* p8 ftext(x1,x2,y+.01,str)- t) c8 Q4 r: U
hold on
; x: y X5 S9 G# B7 m* A2 p7 N% N- s[x11,x22]=meshgrid(min(x1):range(x1)/80:max(x1),min(x2):range(x2)/80:max(x2));
# n% g4 n5 |0 p" p! iyhat=fx1(b,x11,x22); T' e, X( w) x' I1 b6 A
suRF(x11,x22,yhat)
- K% T/ e w' f4 v/ Wshading interp- E6 ^9 K! Q" b' C5 ^3 B
alpha(.8)
3 W3 c" d6 {" v& j R4 vaxis tight) C9 r/ q- e9 r6 c* I" O
xlabel('α');ylabel('C'),zlabel('V1'): l$ s9 |8 e$ G; i
SSy=var(y)*(n-1)( r4 c7 a) [$ I8 Y3 G H/ b
y1=fx1(b,x1,x2);
; n; s; z6 g) Y* D" Z9 }( cRSS=(y-y1)'*(y-y1)0 X( j; w1 @* [/ ?* g( x
rsquare=(SSy-RSS)/SSy Z7 A% {% t& E& J7 O7 {$ x1 y
MSe=RSS/(n-length(b))9 R, T: Z u4 T; c! \/ h
$ y- O) s- |7 n6 ^5 {* C6 F G: P" ? |
|