24小时热门版块排行榜    

查看: 3876  |  回复: 8

cheng1378653

新虫 (初入文坛)

[求助] 用matlab里的nlinfit如何拟合分段函数

Sample Text
t=[0 0.0833 0.25 0.5 1 2];
c=[0 0.0833 1.2491 1.7025 0.2623 0.0757];
函数形式为c=A*(1-exp(-K1*t))       0<=t<=T
                c=A*(1-exp(-k1*T))*exp(-k2*(t-T))     t>=T
如何用脚本式或是函数式M文件拟合,求出A,k1,k2,(如果可以请划一下图形)
回复此楼

» 猜你喜欢

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

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
回帖置顶 ( 共有1个 )

lgycjpcqu

金虫 (正式写手)

【答案】应助回帖

★ ★ ★ ★
cheng1378653: 回帖置顶 2013-03-29 18:09:29
cheng1378653: 金币+2, 有帮助 2013-03-30 22:02:11
fegg7502: 金币+2, 应助指数+1, 鼓励交流 2013-04-02 09:24:40
还是2楼的结果好,用刚才我用1sopt运算的结果作为初始值,用MATLAB进行拟合得到结果和1楼一致,看来非线性最优化还是1sopt强大啊
A =
-0.0367
k1 =
-14.2282
k2 =
3.7015
T =
0.3179
你需要MATLAB 代码的话如下
clear
tic
fun1=inline('beta(1).*(1-exp(-beta(2).*t))','beta','t');
fun2=inline('beta(1).*(1-exp(-beta(2).*beta(4)))*exp(-beta(3)*(t-beta(4)))','beta','t');
fun=@(beta,t)((t<0.5).*fun1(beta,t)+(t>=0.5).*fun2(beta,t));
t=[0 0.0833 0.25 0.5 1 2];
c=[0 0.0833 1.2491 1.7025 0.2623 0.0757];
beta=nlinfit(t,c,fun,[-0.05  -6.8 3.7 0.5]);
A=beta(1);
k1=beta(2);
k2=beta(3);
T=beta(4);
warning off all
c1=beta(1).*(1-exp(-beta(2).*t(1:4)));
c2=beta(1).*(1-exp(-beta(2).*beta(4))).*exp(-beta(3).*(t(4:6)-beta(4)));
cfit(1:4)=c1;
cfit(4:6)=c2;
plot(t,c,'*')
hold on;
plot(t,cfit,'r-')
toc
4楼2013-03-29 10:47:34
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
普通回帖

dingd

铁杆木虫 (职业作家)

【答案】应助回帖


感谢参与,应助指数 +1
fegg7502: 金币+1, 鼓励交流 2013-04-02 09:24:49
1stOpt更简单:
CODE:
Function c=if(t Data;
t=[0 0.0833 0.25 0.5 1 2];
c=[0 0.0833 1.2491 1.7025 0.2623 0.0757];

均方差(RMSE): 0.0282872879644411
残差平方和(SSE): 0.00480102396229929
相关系数(R): 0.999311354215736
相关系数之平方(R^2): 0.998623182664489
决定系数(DC): 0.998184285936899
卡方系数(Chi-Square): 0.0580621164226486
F统计(F-Statistic): 367.165006987895

参数        最佳估算
----------        -------------
t1        0.317859946002189
a        -0.0366738585259589
k1        -14.2282073650289
k2        3.70145325521581

jv1.jpg

2楼2013-03-29 08:50:59
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lgycjpcqu

金虫 (正式写手)

【答案】应助回帖


感谢参与,应助指数 +1
fegg7502: 金币+1, 鼓励交流 2013-04-02 09:24:59
引用回帖:
2楼: Originally posted by dingd at 2013-03-29 08:50:59
1stOpt更简单:

Function c=if(t<T1,A*(1-exp(-K1*t)),A*(1-exp(-k1*T1))*exp(-k2*(t-T1)));
Data;
t=;
c=;


均方差(RMSE): 0.0282872879644411
残差平方和(SSE): 0.00480102396229929
相关系数(R): ...

为啥我的1sopt运算结果不一致呢?
参数        最佳估算
----------        -------------
t1                 0.500000049987681
a                  -0.0568673833778535
k1                 -6.86274581292278
k2                 3.70079073587775
不是说1sopt得到的是唯一解吗?
3楼2013-03-29 10:38:11
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

dingd

铁杆木虫 (职业作家)

【答案】应助回帖


fegg7502: 金币+1, 应助指数+1, 3ks 2013-04-02 09:25:13
引用回帖:
3楼: Originally posted by lgycjpcqu at 2013-03-29 10:38:11
为啥我的1sopt运算结果不一致呢?
参数        最佳估算
----------        -------------
t1                 0.500000049987681
a                  -0.0568673833778535
k1                 -6.86274581292278
k2                 3.70079073587775
...

如果你用的1stOpt老版本(如1.5),可是如下代码:
CODE:
Variable t,c;
Function c=if(t Data;
0        0
0.0833        0.0833
0.25        1.2491
0.5        1.7025
1        0.2623
2        0.0757

5楼2013-03-29 11:15:58
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lgycjpcqu

金虫 (正式写手)

引用回帖:
5楼: Originally posted by dingd at 2013-03-29 11:15:58
如果你用的1stOpt老版本(如1.5),可是如下代码:

Variable t,c;
Function c=if(t<T1,A*(1-exp(-K1*t)),A*(1-exp(-k1*T1))*exp(-k2*(t-T1)));
Data;
0        0
0.0833        0.0833
0.25        1.2491
0.5        1.7025
1        0. ...

这次对了,看来是我的版本太老了
6楼2013-03-29 13:53:09
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

cheng1378653

新虫 (初入文坛)

引用回帖:
4楼: Originally posted by lgycjpcqu at 2013-03-29 10:47:34
还是2楼的结果好,用刚才我用1sopt运算的结果作为初始值,用MATLAB进行拟合得到结果和1楼一致,看来非线性最优化还是1sopt强大啊
A =
-0.0367
k1 =
-14.2282
k2 =
3.7015
T =
0.3179
你需要MATLAB 代码的 ...

你好,你的这个代码是默认t=0.5时取得最大值,而我是想要根据拟合的结果,求出取得最大值得T值,谢谢你的帮助
7楼2013-03-30 22:01:05
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lgycjpcqu

金虫 (正式写手)

【答案】应助回帖


fegg7502: 金币+1, 应助指数+1, 鼓励交流 2013-04-02 09:25:26
引用回帖:
7楼: Originally posted by cheng1378653 at 2013-03-30 22:01:05
你好,你的这个代码是默认t=0.5时取得最大值,而我是想要根据拟合的结果,求出取得最大值得T值,谢谢你的帮助...

fun=@(beta,t)((t<0.5).*fun1(beta,t)+(t>=0.5).*fun2(beta,t)); 这一句中的0.5应该换成beta(4),这个是一个参数,我忘改过来了,不好意思,不过和其它参数一样也必须有一个合理的初始值才能得到好的结果。
8楼2013-03-31 10:34:12
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

cheng1378653

新虫 (初入文坛)

引用回帖:
8楼: Originally posted by lgycjpcqu at 2013-03-31 10:34:12
fun=@(beta,t)((t<0.5).*fun1(beta,t)+(t>=0.5).*fun2(beta,t)); 这一句中的0.5应该换成beta(4),这个是一个参数,我忘改过来了,不好意思,不过和其它参数一样也必须有一个合理的初始值才能得到好的结果。...

谢谢
9楼2013-03-31 12:25:54
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 cheng1378653 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[论文投稿] 售SCI一区T0P文章,我:8.O.55.1.O.5.4,科目全,可+急 +3 jRl3mE6ddZGq 2026-08-03 7/350 2026-08-04 04:50 by k7OM8YghWbkC
[硕博家园] 售SCI一区T0P文章,我:8.O.55.1.O.54,科目齐全,可+急 +3 jRl3mE6ddZGq 2026-08-03 6/300 2026-08-04 04:47 by k7OM8YghWbkC
[找工作] 售SCI一区文章,我:8.O.551.O.5.4,科目全,可伽急 +3 cu9Nq1xK233Z 2026-08-03 8/400 2026-08-04 04:47 by k7OM8YghWbkC
[博后之家] 售SCI一区文章,我:8O5.5.1.O5.4,科目全,可伽急 +3 cu9Nq1xK233Z 2026-08-03 8/400 2026-08-04 04:46 by k7OM8YghWbkC
[公派出国] 售SCI一区文章,我:8.O.551.O.5.4,科目全,可伽急 +3 cu9Nq1xK233Z 2026-08-03 7/350 2026-08-04 04:46 by k7OM8YghWbkC
[教师之家] 售SCI一区文章,我:8.O.551.O.5.4,科目全,可伽急 +4 cu9Nq1xK233Z 2026-08-03 8/400 2026-08-04 04:41 by k7OM8YghWbkC
[有机交流] 一个有机合成实验室都需要哪些设备? 50+3 kf2781974 2026-07-31 9/450 2026-08-04 03:56 by kf2781974
[考研] 售SCI一区文章,我:8.O.55.1.O.54,科目齐全,可伽急 +3 cu9Nq1xK233Z 2026-08-03 8/400 2026-08-04 03:21 by k7OM8YghWbkC
[硕博家园] 售SCI一区T0P文章,我:8.O.55.1.O.54,科目齐全,可+急 +4 cu9Nq1xK233Z 2026-08-03 12/600 2026-08-04 03:17 by k7OM8YghWbkC
[基金申请] 什么时候能放榜呀? +3 Jacob678 2026-08-03 3/150 2026-08-03 16:14 by gltch
[高分子] UV压敏胶开发 +3 ichall 2026-07-30 5/250 2026-08-03 14:40 by Sunrisepay
[基金申请] 面上提前没消息,有中的吗 +14 archvillain 2026-08-02 16/800 2026-08-03 12:23 by fuweiguochen
[基金申请] 面上再次挂了,太难了,躺也躺不了,倦也卷不过,小学校之殇! +19 低垂的野花 2026-07-31 25/1250 2026-08-03 09:25 by gy116024
[基金申请] 2026年国自然面上资助率 +16 布布和一二 2026-07-30 22/1100 2026-08-03 09:20 by gy116024
[教师之家] 基础研究怎么拉横向,学校到款任务越来越多,难以完成 拉横向,都有哪些途径啊 +8 锦衣卫寒战 2026-07-28 12/600 2026-08-03 08:22 by syl200707
[基金申请] 微信指数没变化,科研之友没阅读 +15 wangze12014 2026-07-28 19/950 2026-08-02 20:05 by 蔡棒棒菂
[考博] 2027年申博 50+3 射雕英雄胜 2026-07-30 3/150 2026-08-02 09:36 by lfy8008
[高分子] HXDI做水性聚氨酯乳液,是不是特别容易出渣 15+3 yuyusuv 2026-07-29 3/150 2026-07-31 09:21 by huizingga
[基金申请] 系统今天又提示维护了,估计离放榜不远了 +11 winnerche 2026-07-29 15/750 2026-07-30 23:08 by jnhyjjm
[基金申请] 你们的时间戳变了吗 +3 archvillain 2026-07-30 4/200 2026-07-30 18:53 by levinzhwen
信息提示
请填处理意见