EDA365电子论坛网
标题:
matlab拟合得到二元函数
[打印本页]
作者:
zzz.dan
时间:
2020-12-3 13:50
标题:
matlab拟合得到二元函数
想通过matlab拟合得到二元函数y=f(x1,x2)的表达式,但是这里面18组初始的b值是参考别人的模板,因为是用另外的数据生成的,请问怎么用程序实现我这组数据里面的18组初始b值呢?比如用nlinfit,或者cftool,请大佬多多指教!
之前用cftool拟合后得到的poly23多项式函数的函数值域与已知数据相差几个数量级,函数表达式明显是错误的
`1 P: O% ] I' |3 k- [% _
代码如下:
# l# y! _& G$ B# A( U5 m$ o' k( A- i2 I
clear,clc
+ c( V8 B& W1 N/ E$ H6 ^3 o5 v$ r
x1=[1.2:.1:1.7]';
4 }8 L$ e- O" S0 _7 S
x2=[1000 5000 7500 10000 15000];
. z- J/ ^! B3 p/ `, R+ s, t4 a
y=[14964 13166 14235 15550 19200
7 S% {2 M+ S. V: t' B- |4 ?4 T
13479 13479 13090 14235 17014
5 d4 u, W% X+ a8 F+ S$ \
13747 13747 12750 13139 14991
; x y$ w) W7 ^) ^ O* t+ m* M3 M
13917 13917 13019 12553 13967
( H7 b" v. O* R: v2 t
14065 14065 13386 12821 13188
; J* |! D' ?, ~7 b$ A* _/ Q4 f
14306 14306 13578 13040 12728];
$ a: ?! m* k1 g) J/ Y2 a
n1=length(x1);n2=length(x2);
- @2 u7 A7 _6 l
x1=x1*ones(1,n2);x1=x1(
;
, O$ _$ m' C+ g/ g8 ^" Y3 K& Y
x2=ones(n1,1)*x2;x2=x2(
;
- k9 k( ~% T' g0 r# L4 G
y=y(
; X=[x1,x2]; n=length(y);
/ y# o* I H$ N+ L9 |
str=num2str([1:n]');
4 p& {( v! l& q8 z
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));
4 F+ n- ^9 g, C3 m1 |: \% ?3 X4 ^
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)));
9 Q J! m" h! {' K0 V1 t/ t. H
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];
; t% Q, E9 ~: G7 D3 ]4 a1 \
for l=1:5
, m/ o" ?5 C9 y9 A% [
b=lsqcurvefit(fx2,b,X,y);
- F1 h" v7 y( U
b=nlinfit(X,y,fx2,b);
$ A2 s2 w9 e1 s0 d; ~
end
/ `$ m: q$ \# h# W, H
b
$ Q# p/ I. ~3 A( X u
figure(1),clf
) a& {- L) W5 U, V! G
plot3(x1,x2,y,'o')
4 Q. G, ^ U$ {* }- b, A2 z5 ?
stem3(x1,x2,y,'filled')
% s9 l" l' _; u0 ~# J3 m6 t8 j
text(x1,x2,y+.01,str)
8 a$ X. v1 n; o0 j% }
hold on
9 q+ J$ m6 [8 t7 s. d" L
[x11,x22]=meshgrid(min(x1):range(x1)/80:max(x1),min(x2):range(x2)/80:max(x2));
& c! l5 u @" b6 Q9 H2 Q
yhat=fx1(b,x11,x22);
7 D9 L& }% g2 N; ~( d$ R
surf(x11,x22,yhat)
8 F+ q9 g+ j( `: M( v
shading interp
6 O' ~+ r2 E7 b! C' }9 R
alpha(.8)
5 S1 h1 I2 c5 L o/ m& J5 ]( H
axis tight
0 ]% S" p7 B) g- }+ d7 p
xlabel('α');ylabel('C'),zlabel('V1')
K8 `, l5 f: t1 p( U& `
SSy=var(y)*(n-1)
; P) h) D. G& }/ Q5 J5 ]
y1=fx1(b,x1,x2);
p$ D! D z$ Y8 h5 v; m) A
RSS=(y-y1)'*(y-y1)
, A0 E1 B, i. E. m; r; q9 z
rsquare=(SSy-RSS)/SSy
7 |* L& i& F' |/ e: X5 C
MSe=RSS/(n-length(b))
7 _/ ^3 G4 X. _- y3 P* u
+ R$ w8 N8 _) G" D3 G- u+ x: T6 O
作者:
Uifhjvv
时间:
2020-12-3 14:26
帮你顶一下
作者:
nkkopd
时间:
2020-12-3 14:31
数据不同,不能用此前的模型和初值。另,只有少量的数据是不能用很复杂(包含很多效应项)的模型的。
作者:
小小鲁班
时间:
2020-12-3 16:42
欢迎光临 EDA365电子论坛网 (https://bbs.eda365.com/)
Powered by Discuz! X3.2