24小时热门版块排行榜    

查看: 2427  |  回复: 4

anarky

新虫 (初入文坛)

[求助] 请教matlab非线性拟合结果不唯一的问题 已有2人参与

我用matlab拟合一组实验数据,函数是y=a*(b+x)^c+d*(1-f*exp(-g*x)),现在的问题是发现只要改一次初值就会得到不同的拟合结果,或者改一下上下限,也会得到不同的拟合结果,这么一来,我觉得就不可能找到最优解,虽然上下限可以取得很大,认为是无界,但是初值怎么取?因为我不可能尝试所有的参数初值,因为我要拟合多组实验数据,获得公式中参数a(1)~a(6)的变化规律,如果参数结果不唯一,我怎么确定哪个参数才是最优的参数呢?事实上,如果从物理意义来限定参数取值也可以,但是那样得到的参数并不是最优解,而最优解的参数又不一定符合物理意义。 我该怎么处理呢?哪位能帮忙看一下,非常感谢!
实验数据和代码如下:
主程序
clc
xdata=load('78028-ta.txt');  %read section1
size1=size(xdata);              %section1 行、列数
sizen(1)=size1(1);               %section1数据量
for i=1:size1;
xd(i)=xdata(i,1);
yd(i)=xdata(i,2);
end

set(gca,'Fontsize',15);
plot(xd,yd,'ro');
xlabel('True Strain');
ylabel('True Stress/MPa');

hold on
lb=[-100000, 0,    0, -10000, -10000, -10000];
ub=[ 100000, 0.02, 1,  10000,  10000,  10000];
a1=1898;a2=0.00002;a3=0.136;a4=542;a5=0.784;a6=62.4;
   a0=[a1,a2,a3,a4,a5,a6];
options = optimset('TolFun',1e-6,'LargeScale','on','TolCon',1.0e-6,'MaxFunEvals',2.e10,'MaxIter',2000000,'TolX',1e-6);
exitflag=0;
while (exitflag<=0)                                %反复迭代
    [a,resnorm,residual,exitflag]=lsqnonlin(@SV,a0,lb,ub,options);      %nonlinear fitting, solve coefficient a
    a0=a;
end
R2=1-resnorm/sum((yd-mean(yd)).^2);
standard=sqrt(resnorm/size(residual,2));
stp=1;
out(stp,1)=R2;
out(stp,2)=standard;
out(stp,3)=a(1);
out(stp,4)=a(2);
out(stp,5)=a(3);
out(stp,6)=a(4);
out(stp,7)=a(5);
out(stp,8)=a(6);

hold on
fid=fopen('SV-fit-78028-ta.txt','wt');
for i=1:size1;
     G(i)=a(1)*(a(2)+xd(i))^a(3)+a(4)*(1-a(5)*exp(-a(6)*xd(i)));
   fprintf(fid,'%8.6f   %10.6f\n',xd(i),G(i));
end

fprintf(fid,'     R2        standard       s1             s2           s3         V1         V2         V3\n');
fprintf(fid,'%10.6f   %10.6f   %10.6f   %10.6f   %10.6f   %10.6f   %10.6f   %10.6f\n',out(i,1),out(i,2),out(i,3),out(i,4),out(i,5),out(i,6),out(i,7),out(i,8));
fclose(fid);
plot(xd,G,'b','LineWidth',2)
hold on
legend('EXP','SV');
hold off
clear
函数
function F=SV(a)
xdata=load('78028-ta.txt');
size1=size(xdata);
sizen(1)=size1(1);

for i=1:sizen(1);
xd(i)=xdata(i,1);
yd(i)=xdata(i,2);
end

for i=1:sizen(1);
    F(i)=a(1)*(a(2)+xd(i))^a(3)+a(4)*(1-a(5)*exp(-a(6)*xd(i)));
F(i)=F(i)-yd(i);
end

实验数据
78028-ta.txt
0        223.634
0.0005        261.9
0.001        373.358
0.0015        425.658
0.002        462.61
0.0025        492.699
0.003        513.549
0.0035        524.142
0.004        538.681
0.0045        546.22
0.005        554.575
0.0055        563.233437
0.006        570.27162
0.0065        577.07317
0.007        583.656651
0.0075        590.037892
0.008        596.230537
0.0085        602.246461
0.009        608.09609
0.0095        613.788651
0.01        619.332371
0.0105        624.734633
0.011        630.002102
0.0115        635.14083
0.012        640.156342
0.0125        645.053702
0.013        649.837576
0.0135        654.512279
0.014        659.081814
0.0145        663.549911
0.015        667.920058
0.0155        672.19552
0.016        676.379368
0.0165        680.474496
0.017        684.483635
0.0175        688.409373
0.018        692.254165
0.0185        696.020341
0.019        699.710123
0.0195        703.325629
0.02        706.4514911
0.025        735.7217232
0.03        758.4120509
0.035        778.3482998
0.04        794.7567166
0.045        809.6868333
0.05        822.8422861
0.055        834.8128986
0.06        846.1687589
0.065        856.1567004
0.07        866.0042315
0.075        874.4658077
0.08        882.1951827
0.085        890.0445724
0.09        897.0900675
0.095        903.7413962
0.1        910.0615411
0.105        915.8285085
0.11        921.8803708
0.115        927.2834332
0.12        932.4186891
0.125        937.3558307
0.315        1073.901289
回复此楼

» 收录本帖的淘帖专辑推荐

Matlab

» 猜你喜欢

» 本主题相关价值贴推荐,对您同样有帮助:

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dingd

铁杆木虫 (职业作家)

【答案】应助回帖

★ ★ ★ ★
感谢参与,应助指数 +1
anarky: 金币+4, ★★★很有帮助, 非常感谢! 2015-12-09 21:24:24
用1stOpt,不用初值,结果唯一稳定:
CODE:
Function  y=a*(b+x)^c+d*(1-f*exp(-g*x));
Data;
0        223.634
0.0005        261.9
0.001        373.358
0.0015        425.658
0.002        462.61
0.0025        492.699
。。。。。

均方差(RMSE):6.83618029811499
残差平方和(SSE):2944.20174730514
相关系数(R): 0.99915354759785
相关系数之平方(R^2): 0.998307811677368
确定系数(DC): 0.998307811677368
卡方系数(Chi-Square): 4.49199298531208
F统计(F-Statistic): 6736.83883025038

参数        最佳估算
----------        -------------
a        -1838.97884865162
b        0.00549570600631707
c        -0.0663243340918725
d        3047.85332211066
f        0.0810032510811905
g        817.309082573211

其它软件很难得到上面的结果
请教matlab非线性拟合结果不唯一的问题
c2.jpg

2楼2015-12-08 10:18:31
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

yanshuiting

新虫 (正式写手)

Matlab寻优过程得到的是局部最优解,如果这个解能够满足实际需求,就可以使用,不用太过顾虑全局最优解。当模型复杂时,全局最优解释很难得到的。

发自小木虫Android客户端
3楼2015-12-08 12:36:50
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

anarky

新虫 (初入文坛)

引用回帖:
3楼: Originally posted by yanshuiting at 2015-12-08 12:36:50
Matlab寻优过程得到的是局部最优解,如果这个解能够满足实际需求,就可以使用,不用太过顾虑全局最优解。当模型复杂时,全局最优解释很难得到的。

感谢参与!
4楼2015-12-16 21:03:30
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

ybkooo

至尊木虫 (著名写手)

^^

【答案】应助回帖

算法有局限性,可以考虑 遗传算法、粒子群算法;结果应该稳定
abcd
5楼2015-12-16 21:54:02
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 anarky 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[基金申请] 欢迎发来filecode的Mz6后的代码验证其规律 +55 医学老男孩 2026-08-13 123/6150 2026-08-19 23:44 by stan_et
[基金申请] 今天基金会出结果吗?20260819 +13 kkkl_v 2026-08-19 13/650 2026-08-19 22:37 by 宝贝虫子
[基金申请] filecode,4个jtjc了 +11 ziyangfang 2026-08-19 13/650 2026-08-19 21:28 by aasahr
[教师之家] 为什么余额宝的年化利率越来越低?主要原因有哪些? +5 瞬息宇宙 2026-08-15 5/250 2026-08-19 20:23 by super2002521
[基金申请] 今天放榜没戏了吧 +4 yuleib84 2026-08-19 5/250 2026-08-19 19:52 by hhs666
[基金申请] filecode=后面第一个是大写字母 +10 wangze12014 2026-08-14 12/600 2026-08-19 16:56 by 苦难博士
[基金申请] 2027广东省杰青 +4 奶牛小黑 2026-08-15 10/500 2026-08-19 11:02 by wanfengnew
[基金申请] 朋友圈看到的 +6 wangzilk 2026-08-18 8/400 2026-08-19 10:55 by Haru815
[基金申请] 快农历七夕节了,轻松一下,男人悄悄话,女施主请不要进来。 +6 Tide man 2026-08-14 7/350 2026-08-19 09:56 by ZJTJZ
[基金申请] 什么时候开奖? +6 CrisMessi 2026-08-18 6/300 2026-08-19 08:06 by Equinoxhua
[基金申请] 明天放榜? +5 Shxjjxjkx 2026-08-18 5/250 2026-08-18 18:14 by -大大大大大-
[论文投稿] 投稿咨询 +4 wwm09 2026-08-17 6/300 2026-08-18 15:36 by wwm09
[基金申请] 时间戳又变了8-15 +14 archvillain 2026-08-15 26/1300 2026-08-18 13:37 by phantomgost
[基金申请] 今天维护系统维护 祝所有人 高中 +8 gjjjzhong 2026-08-18 9/450 2026-08-18 13:01 by 家与远方
[基金申请] 哪位老哥知道今年的国自然具体哪一天放榜? +13 Ldrop2023 2026-08-13 16/800 2026-08-18 12:25 by 淀粉搬运工
[基金申请] 93BebMhtakh前后11位开头都是大写 +4 且听虎啸 2026-08-17 5/250 2026-08-18 00:49 by 蔡棒棒菂
[基金申请] 感觉是下周放榜了 +6 angus9576 2026-08-17 11/550 2026-08-17 23:57 by angus9576
[基金申请] 今天系统多次维护,明天很可能放榜! +8 zju2000 2026-08-16 9/450 2026-08-17 12:20 by lmz0216
[精细化工] 招聘 金属平磨液,抛光液研发工程师 +3 小天0311 2026-08-14 3/150 2026-08-16 07:31 by H9PLUS
[基金申请] 各位道友,我要去昆明玩几天,回来见。 +7 Tide man 2026-08-14 8/400 2026-08-15 01:11 by arzu_hma
信息提示
请填处理意见