声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 3003|回复: 5

[1stopt] (求助)1stopt非线性约束拟合问题。

[复制链接]
发表于 2009-12-17 00:06 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
各位好,
    小弟有个问题求教。
     我现在用的是1stopt 1.5注册版。我编写的拟合代码如下。我想加进去几个参数的约束条件:a1+a2+a3+a4=100;k1<k2<k3<k4;e1<e2<e3<e4;m1>m2>m3>m4。

    请教一下,我应该怎么修改我的拟合代码,还有就是1stopt 1.5注册版是否可以进行有非线性约束条件的拟合问题


Parameters a(1:4)=[0,100],k(1:4)=[0,],e(1:4)=[0,],m(1:4)=[0,];
Variable x1, x2, x3, y;
Algorithm = GLM;
Function y=a1*exp(-k1*exp(-e1*x1)*x2*x3^m1)+a2*exp(-k2*exp(-e2*x1)*x2*x3^m2)+a3*exp(-k3*exp(-e3*x1)*x2*x3^m3)+a4*exp(-k4*exp(-e4*x1)*x2*x3^m4);
Data;
1.841945596 4000000 10000000 2.104558482
1.841945596 4000000 6666666.667 2.527143784
1.841945596 4000000 5000000 3.451810632
1.841945596 4000000 4000000 3.964651068
1.900142929 4000000 20000000 2.131455708
1.900142929 5000000 20000000 1.198420834
1.900142929 6000000 20000000 0.970093274
1.900142929 4000000 10000000 2.719608376
1.900142929 5000000 10000000 2.383990771
1.900142929 6000000 10000000 2.918647847
1.900142929 4000000 6666666.667 3.802968856
1.900142929 5000000 6666666.667 2.681653402
1.900142929 6000000 6666666.667 3.239322549
1.900142929 4000000 5000000 5.688763236
1.900142929 5000000 5000000 3.298197587
1.900142929 6000000 5000000 3.597653366
1.900142929 4000000 4000000 6.884494347
1.900142929 5000000 4000000 3.75455385
1.900142929 6000000 4000000 4.490043544
1.962137804 4000000 20000000 5.0369538
1.962137804 5000000 20000000 2.810760086
1.962137804 6000000 20000000 2.673584235
1.962137804 7000000 20000000 2.493671681
1.962137804 4000000 10000000 5.70818901
1.962137804 5000000 10000000 3.317623361
1.962137804 6000000 10000000 3.723472611
1.962137804 7000000 10000000 3.102744414
1.962137804 4000000 6666666.667 6.162453266
1.962137804 5000000 6666666.667 4.796671917
1.962137804 6000000 6666666.667 4.907548258
1.962137804 7000000 6666666.667 3.863637043
1.962137804 4000000 5000000 6.600280329
1.962137804 5000000 5000000 6.791549489
3.537619041 6000000 5000000 5.626600759
3.537619041 7000000 5000000 4.22824388
3.537619041 4000000 4000000 7.181559262
3.537619041 5000000 4000000 11.35660641
3.537619041 6000000 4000000 6.116727982
3.537619041 7000000 4000000 4.97359589
3.758720231 4000000 20000000 7.660927596
3.758720231 5000000 20000000 3.451810632
3.758720231 6000000 20000000 3.035202491
3.758720231 4000000 10000000 15.50087117
3.758720231 5000000 10000000 5.441308759
3.758720231 6000000 10000000 5.996288183
3.758720231 4000000 6666666.667 20.39228109
3.758720231 5000000 6666666.667 7.797206873
3.758720231 6000000 6666666.667 6.893161231
3.758720231 4000000 5000000 26.88258166
3.758720231 5000000 5000000 9.56764204
3.758720231 6000000 5000000 8.3115416
3.758720231 5000000 4000000 12.25288174
3.758720231 6000000 4000000 9.261910241
4.00930158 4000000 20000000 21.50582624
4.00930158 5000000 20000000 12.84103441
4.00930158 4000000 10000000 29.81019524
4.00930158 5000000 10000000 21.15167943
4.00930158 4000000 6666666.667 34.19175331
4.00930158 5000000 6666666.667 32.22586497
4.00930158 4000000 5000000 39.11454333
4.00930158 5000000 5000000 36.92421258
4.00930158 4000000 4000000 44.77790363
4.295680264 4000000 20000000 38.78281088
4.295680264 7000000 20000000 17.00890896
4.295680264 4000000 10000000 44.59858879
4.295680264 7000000 10000000 29.80451694
4.295680264 4000000 6666666.667 49.3818121
4.295680264 7000000 6666666.667 33.3594336
4.295680264 4000000 5000000 53.69283966
4.295680264 7000000 5000000 39.9399893
4.295680264 4000000 4000000 60.05552783
回复
分享到:

使用道具 举报

发表于 2009-12-18 16:52 | 显示全部楼层
鉴于时间关系,帮你算了一次,也许还有更好的解。自己验证下是否正确,:

均方差(RMSE): 3.21299918227978
残差平方和(RSS): 732.95882591847
相关系数(R): 0.976603906435484
相关系数之平方(R^2): 0.953755190065048
决定系数(DC): 0.953204477559434
卡方系数(Chi-Square): 28.8026624376944
F统计(F-Statistic): 78.3550749181967
约束条件: a1+a2+a3+a4-(100) = 0
          k2-(k1) = 0.000149765498
          k3-(k2) = 1.231094538E-7
          k4-(k3) = 26.12829726
          m2-(m1) = 0.8322790814
          m3-(m2) = 0.002957005434
          m4-(m3) = 3.387003462
          e2-(e1) = 0.0003348275783
          e3-(e2) = 4.71332591
          e4-(e3) = 7.425811234

参数        最佳估算
----------        -------------
a1                 26.0376395272275
a2                 13.1385093955078
a3                 42.7324353777972
a4                 18.0914156994675
k1                 7.01036479084896E-7
k2                 0.000150466534467939
k3                 0.000150589643921733
k4                 26.1284478495611
e1                 0.280539463839992
e2                 0.280874291418277
e3                 4.99420020162369
e4                 12.4200114356138
m1                 7.10606663318616E-16
m2                 0.832279081387495
m3                 0.835236086821955
m4                 4.22223954927176
 楼主| 发表于 2009-12-19 23:34 | 显示全部楼层
:loveliness:  谢谢dingd大哥,小弟,谢谢了,
 楼主| 发表于 2009-12-20 18:45 | 显示全部楼层
dingd 斑竹,感谢你阅读本贴。你给出的数据,小弟带入之后,发现之前编写的data项把x2和x3的数据列颠倒了。而且我也发现我之前设计的约束条件太简单了。
   请dingd斑竹,帮小弟再跑一次。再次感谢dingd斑竹对我这个菜鸟的帮助!

Parameters a1=[0,100],a2=[0,100],a3=[0,100],a4=[0,100],k1=[0,],k2=[0,],k3=[0,],k4=[0,],e1=[2,25],e2=[2,25],e3=[2,25],e4=[7,25],m1=[0,5],m2=[0,5],m3=[0,5],m4=[0,10];
Variable x1, x2, x3, y;
Function y=a1*exp(-k1*exp(-e1*x1)*x2*x3^m1)+a2*exp(-k2*exp(-e2*x1)*x2*x3^m2)+a3*exp(-k3*exp(-e3*x1)*x2*x3^m3)+a4*exp(-k4*exp(-e4*x1)*x2*x3^m4);
         a1+a2+a3+a4-100=0;
         e2-e1>0;
         e3-e2>0;
         e4-e3>0;
         m1-m2>0;
         m2-m3>0;
         m3-m4>0;
         k1*exp(-e1*x1)*x3^m1-k2*exp(-e2*x1)*x3^m2>0;
         k2*exp(-e2*x1)*x3^m2-k3*exp(-e3*x1)*x3^m3>0;
         k3*exp(-e3*x1)*x3^m3-k4*exp(-e4*x1)*x3^m4>0;
Data;
1.841945596 10000000 4000000 2.104558482
1.841945596 6666666.667 4000000 2.527143784
1.841945596 5000000         4000000 3.451810632
1.841945596 4000000         4000000 3.964651068
1.900142929 20000000 4000000 2.131455708
1.900142929 20000000 5000000 1.198420834
1.900142929 20000000 6000000 0.970093274
1.900142929 10000000 4000000 2.719608376
1.900142929 10000000 5000000 2.383990771
1.900142929 10000000 6000000 2.918647847
1.900142929 6666666.667 4000000 3.802968856
1.900142929 6666666.667 5000000 2.681653402
1.900142929 6666666.667 6000000 3.239322549
1.900142929 5000000         4000000 5.688763236
1.900142929 5000000         5000000 3.298197587
1.900142929 5000000         6000000 3.597653366
1.900142929 4000000         4000000 6.884494347
1.900142929 4000000         5000000 3.75455385
1.900142929 4000000         6000000 4.490043544
1.962137804 20000000 4000000 5.0369538
1.962137804 20000000 5000000 2.810760086
1.962137804 20000000 6000000 2.673584235
1.962137804 20000000 7000000 2.493671681
1.962137804 10000000 4000000 5.70818901
1.962137804 10000000 5000000 3.317623361
1.962137804 10000000 6000000 3.723472611
1.962137804 10000000 7000000 3.102744414
1.962137804 6666666.667 4000000 6.162453266
1.962137804 6666666.667 5000000 4.796671917
1.962137804 6666666.667 6000000 4.907548258
1.962137804 6666666.667 7000000 3.863637043
1.962137804 5000000         4000000 6.600280329
1.962137804 5000000         5000000 6.791549489
3.537619041 5000000         6000000 5.626600759
3.537619041 5000000         7000000 4.22824388
3.537619041 4000000         4000000 7.181559262
3.537619041 4000000         5000000 11.35660641
3.537619041 4000000         6000000 6.116727982
3.537619041 4000000         7000000 4.97359589
3.758720231 20000000 4000000 7.660927596
3.758720231 20000000 5000000 3.451810632
3.758720231 20000000 6000000 3.035202491
3.758720231 10000000 4000000 15.50087117
3.758720231 10000000 5000000 5.441308759
3.758720231 10000000 6000000 5.996288183
3.758720231 6666666.667 4000000 20.39228109
3.758720231 6666666.667 5000000 7.797206873
3.758720231 6666666.667 6000000 6.893161231
3.758720231 5000000         4000000 26.88258166
3.758720231 5000000         5000000 9.56764204
3.758720231 5000000         6000000 8.3115416
3.758720231 4000000         5000000 12.25288174
3.758720231 4000000         6000000 9.261910241
4.00930158 20000000 4000000 21.50582624
4.00930158 20000000 5000000 12.84103441
4.00930158 10000000 4000000 29.81019524
4.00930158 10000000 5000000 21.15167943
4.00930158 6666666.667 4000000 34.19175331
4.00930158 6666666.667 5000000 32.22586497
4.00930158 5000000         4000000 39.11454333
4.00930158 5000000         5000000 36.92421258
4.00930158 4000000         4000000 44.77790363
4.295680264 20000000 4000000 38.78281088
4.295680264 20000000 7000000 17.00890896
4.295680264 10000000 4000000 44.59858879
4.295680264 10000000 7000000 29.80451694
4.295680264 6666666.667 4000000 49.3818121
4.295680264 6666666.667 7000000 33.3594336
4.295680264 5000000         4000000 53.69283966
4.295680264 5000000         7000000 39.9399893
4.295680264 4000000         4000000 60.05552783
发表于 2009-12-20 23:10 | 显示全部楼层
源代码给你,自己可以随时修改计算:

Parameters a1=[0,100],a2=[0,100],a3=[0,100],a4=[0,100],k1=[0,],k2=[0,],k3=[0,],k4=[0,],e1=[2,25],e2=[2,25],e3=[2,25],e4=[7,25],m1=[0,5],m2=[0,5],m3=[0,5],m4=[0,10];
Variable x1, x2, x3, y;
PassParameter d(3);
StartProgram [Pascal];
Procedure MainModel;
var j: integer;
    td1, td2, td3: double;
Begin
    for j := 0 to DataLength - 1 do begin
        y[j]:=a1*exp(-k1*exp(-e1*x1[j])*x2[j]*x3[j]^m1)+a2*exp(-k2*exp(-e2*x1[j])*x2[j]*x3[j]^m2)+a3*exp(-k3*exp(-e3*x1[j])*x2[j]*x3[j]^m3)+a4*exp(-k4*exp(-e4*x1[j])*x2[j]*x3[j]^m4);
        if j = 0 then begin
           td1 := k1*exp(-e1*x1[j])*x3[j]^m1-k2*exp(-e2*x1[j])*x3[j]^m2;
           td2 := k2*exp(-e2*x1[j])*x3[j]^m2-k3*exp(-e3*x1[j])*x3[j]^m3;
           td3 := k3*exp(-e3*x1[j])*x3[j]^m3-k4*exp(-e4*x1[j])*x3[j]^m4;
        end
        else begin
           td1 := min(td1,k1*exp(-e1*x1[j])*x3[j]^m1-k2*exp(-e2*x1[j])*x3[j]^m2);
           td2 := min(td2,k2*exp(-e2*x1[j])*x3[j]^m2-k3*exp(-e3*x1[j])*x3[j]^m3);
           td3 := min(td3,k3*exp(-e3*x1[j])*x3[j]^m3-k4*exp(-e4*x1[j])*x3[j]^m4);
        end;
    end;
    d1 := td1;
    d2 := td2;
    d3 := td3;
    ConstrainedResult := a1+a2+a3+a4-100=0;
    ConstrainedResult := e2-e1>0;
    ConstrainedResult := e3-e2>0;
    ConstrainedResult := e4-e3>0;
    ConstrainedResult := m1-m2>0;
    ConstrainedResult := m2-m3>0;
    ConstrainedResult := m3-m4>0;
    ConstrainedResult := td1>0;
    ConstrainedResult := td2>0;
    ConstrainedResult := td3>0;
End;
EndProgram;
Data;
1.841945596 10000000 4000000 2.104558482
1.841945596 6666666.667 4000000 2.527143784
1.841945596 5000000         4000000 3.451810632
1.841945596 4000000         4000000 3.964651068
。。。。。

一组结果:
均方差(RMSE): 2.45792879014108
残差平方和(RSS): 428.940389555711
相关系数(R): 0.98621286582406
相关系数之平方(R^2): 0.972615816716906
决定系数(DC): 0.972614437652802
F统计(F-Statistic): 133.890501139987

参数                  最佳估算
--------------------        -------------
a1                 11.9706582844648
a2                 28.8527637739238
a3                 53.4020879796602
a4                 5.77177237444129
k1                 1.34505230126704
k2                 56.8870538091767
k3                 6.51913444850881E-12
k4                 16.5798208181613
e1                 3.5871897664056
e2                 3.77098743335754
e3                 4.19805090699037
e4                 21.0165827940815
m1                 4.388965042526
m2                 2.27599478723315
m3                 1.73196980293668
m4                 1.3400953092144

[ 本帖最后由 dingd 于 2009-12-20 23:11 编辑 ]
 楼主| 发表于 2009-12-24 17:32 | 显示全部楼层
谢谢,dingd 大哥,祝各位论坛朋友圣诞快乐,^_^
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-12-29 04:41 , Processed in 0.064088 second(s), 17 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表