找回密码
 欢迎注册
查看: 31206|回复: 17

[讨论] 一道很棘手的函数拟合问题

[复制链接]
发表于 2010-2-4 18:52:32 | 显示全部楼层 |阅读模式

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

您需要 登录 才可以下载或查看,没有账号?欢迎注册

×
某人根据理论推导得知两个物理量存在的关系为:
$y=\frac{1}{p_1+p_2*x^{p_3}}+p_4*x^{p_5}$
接着,他对这两个量进行实验测量,得到下面的一堆的数据,可在拟合这五个参数的值时犯难了,很多软件效果都不理想,特来问问在座的各位。

附:Mathematica拟合函数在默认选项下的表现
Untitled-1.png

data.txt

2.54 KB, 下载次数: 11, 下载积分: 金币 -1 枚, 经验 1 点, 下载 1 次

毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2010-2-4 21:11:31 | 显示全部楼层
要不然可以如此试验一下:
只使用那些x很小的数据,可以忽略掉最后一项,然后使用公式
$1/y-p_1=p_2*x^(p_3)$
可以先大概估计出$p_1$,比如取所有$1/y$的下界,然后代入$p_1$后两边取对数线性拟合可以得到$p_2,p_3$,数据代回去后可以再次矫正$p_1$,然后继续重新计算$p_2,p_3$,计算几轮后可以得到$p_1,p_2,p_3$比较精确的值。
然后再次使用所有数据,可以算出$p_4$
然后可以继续来回迭代几次
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2010-2-4 21:16:23 | 显示全部楼层
附件:data.txt下载后怎么是乱码,文件名为attachment.php
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2010-2-4 21:23:51 | 显示全部楼层
若拟合度始终很低,只能说明建立的模型不准确,或者缺少一些重要的参数(即与结论相关的物理量)
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2010-2-4 21:26:50 | 显示全部楼层
我猜测软件使用的算法大概是按照梯度跌代改进。

这个方法对于比较平滑的参数来说,效果是很不错的。

问题是你这个式子有两个参数都在x的指数上,这样会导致跌代法失效。

如果初值设置得太随意的话,跑不了几步就掉进坑里了。

即使初值设置得比较好,也很难爬到山顶,因为这个地形太奇怪了。

$mathe$的方法很不错,楼主可以试试。

$data.txt$很正常啊。

数据如下:

160.73        6266.7
159.82        6151.9
158.84        6035.1
157.86        5920.9
156.87        5812.6
155.88        5702.2
154.89        5594.9
153.96        5491.3
152.97        5385
151.98        5282.2
150.99        5181.3
150.06        5084.8
149.08        4988.8
148.09        4892.2
147.1        4796.9
146.17        4701
145.18        4608
144.2        4515.2
143.2        4429.6
142.21        4342.9
141.25        4255.6
140.2        4167.1
139.14        4077.6
138.05        3987.9
136.96        3898.9
135.94        3808.5
134.84        3717.7
133.74        3628.9
132.65        3543
131.57        3456.3
130.55        3372.5
129.47        3292.8
128.4        3212.7
127.33        3133.6
126.34        3056.6
125.29        2985.5
124.26        2912.5
123.23        2842.9
122.21        2774.1
121.21        2708.3
120.27        2642.1
119.27        2580.2
118.29        2518.7
117.32        2459.1
116.42        2401.1
115.48        2344.3
114.55        2290.9
113.62        2237.5
112.7        2189
111.85        2138.8
110.94        2089.4
110.03        2042.4
109.13        1998.1
108.28        1953.6
107.38        1906.6
106.48        1867.8
105.6        1824.5
104.72        1784.3
103.91        1745
103.05        1704.5
102.2        1668.7
100.51        1590.4
99.739        1552.5
98.913        1514.7
98.103        1476.4
97.308        1444.3
96.513        1411.4
95.78        1378.5
95.002        1344.8
94.239        1307.8
93.482        1276.1
92.776        1243.5
92.039        1212.9
91.314        1178.7
90.604        1148.4
89.942        1115.9
89.244        1084.5
88.559        1051.5
87.889        1029.7
87.226        996.16
86.569        965.86
85.963        937.72
85.323        907.87
84.694        877.58
84.081        838.17
83.473        819.48
82.876        797.76
82.287        768.54
81.811        749.96
81.178        724.39
80.614        697.24
80.118        674.67
79.574        649.49
79.011        629.83
78.478        614.6
78.012        591.87
77.494        573.43
76.927        558.94
76.512        539.45
75.962        526.99
75.472        514.14
75.014        504.11
74.566        484.4
74.123        473.23
73.608        468.93
73.183        453.77
72.774        448.58
72.369        447.73
71.897        431.79
71.503        432.45
71.116        432.15
70.741        420.71
70.3        427.26
69.935        419.76
69.572        407.28
69.148        408.04
68.796        393.71
68.448        403.74
68.114        408.8
67.717        401.26
67.374        400.81
67.037        401.89
66.741        408.68
66.416        398.49
66.015        414.14
65.373        419.78
64.769        426.82
64.109        418.42
63.44        446.32
62.772        451.55
62.111        473.27
61.508        499.69
60.908        523.66
60.219        551.47
59.699        593.53
59.119        608.69
58.547        658.08
57.992        712.27
57.483        769.4
56.969        826.48
56.472        896.05
55.989        957.57
55.513        1065.1
55.088        1114.1
54.651        1195
54.237        1271.5
53.836        1355.6
53.318        1483.2
52.701        1690
52.08        2245.9
51.431        2470.4
50.877        2719.1
50.298        2957.5
49.74        3155.2
49.2        3279.4
48.702        3546.4
48.182        3741
47.681        4021
47.213        4015.1
46.768        4304.7
46.368        4127.9
45.956        4530.9
45.55        4802.9
45.157        5047.4
44.799        4804.3
44.43        5164.1
44.078        4781
43.727        5175.5
43.384        5708.6
43.079        5679.6
42.899        5161.8
42.719        5399.1
42.547        5483
42.253        4839.4
41.962        5360.7
41.691        5622
40.602        5772.3

不知道复制粘贴的时候有没有操作失误,请楼主核对一下。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2010-2-4 22:48:57 | 显示全部楼层
哈哈,我找到了一种方法:
由于
x的取值范围是40.602-----160.73
y的取值范围是 393.71------6266.7
所以,我可以把x除以40,y除以1000,来代替原来的x,y
这样拟合后效果非常不错:
Untitled-1.png
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2010-2-4 23:01:21 | 显示全部楼层
此题来源于网站(在页面的后面):
http://www.7d-soft.com/Regression_Test.htm

该网站是广告打得很响的1stOpt软件的官网。

七维高科是由留学欧美日归国人员创立的位于北京的高科技公司。公司组成人员的四分之一拥有国外博士学位,平均年龄不到三十岁。公司的宗旨是开发创建拥有自主知识品牌,自主核心技术,领先世界的数值分析,优化计算软件。

美国国家标准与技术研究院(NIST: National Institute of Standards and Technology)提供有一套27道非线性拟合测试题,世界上几乎所有著名的数据分析软件包都以能通过该套测试题集为验证标准。经对比测试,1stOpt是目前唯一不依赖使用NIST提供的初始值,而能以任意随机初始值就可求得全部最优解的软件包(如果使用NIST提供的初始值,则更可轻易求得最优解)。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2010-2-4 23:49:57 | 显示全部楼层
呵,方法很简洁而有效,即对于非线性数据拟合,有必要先需要将x,y均归一化处理(通过线性变换使x,y取值均在0--1之间),然后拟合。。。。这样可能得到的答案更有效
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2010-2-6 10:43:54 | 显示全部楼层
本帖最后由 KeyTo9_Fans 于 2010-2-6 10:46 编辑

wayne太欺负人了!

我像蜗牛一样爬了两天才爬到山顶

没想到这座山居然那么高,高得无法想象。

难怪拟合软件做出来的效果不理想。

我设$y=p_1/(p_2+x^{p_3})+p_4*x^{p_5}$,

和楼主的式子本质上是一样的。

然后随机1亿个初值计算方差,从方差最小的初值开始,以固定的幅度随机抖动5个参数:

如果新的一组参数计算出来的方差更小,则保留新参数,否则把参数还原。

结果$p_3$、$p_4$、$p_5$都比较稳定,而$p_1$和$p_2$在不停地往上滚。

$p_1$初值是100万左右。

我以为往上滚一个数量级就差不多了,所以把步长设置为$-16384$~$16383$的一个随机数。

没想到连续滚了一个晚上还没有停下来。

仔细观察了一下,才发现原来$p_1$和$p_2$的数量级已经上了几个档次。

所以原来的步长显得太小了,计算出来的方差几乎没有什么变化。

于是把步长调大两个数量级(100倍)。

于是方差显著地减小。

说明之前运行一个晚上(500分钟)的效果现在看来只是运行了5分钟而已。

所以我对现在这个速度感到很满意。

然后就这样又跑了一天,发现$p_1$和$p_2$的数量级又上了几个档次。

所以原来的步长又显得太小了,计算出来的方差几乎又没有什么变化了。

于是把步长调大三个数量级(1000倍)。

于是方差又显著地减小了。

说明之前运行一天(1000分钟)的效果现在看来只是运行了1分钟而已。

所以我对现在这个速度感到很满意。

然后就这样又跑了一个晚上,发现$p_1$和$p_2$的数量级又上了几个档次。

所以原来的步长又显得太小了,计算出来的方差几乎又没有什么变化了。

然后终于才反应过来:

步长不应该设置成固定的大小,应该设置成固定的比例!

于是把步长设置为当前参数的千分之一。

然后看到$p_1$和$p_2$指数级递增。

而方差指数级递减,很快就稳定下来了。

最后把步长设置为当前参数的百万分之一、十亿分之一,进行最后的微调。

最终无论如何改动参数,方差都不再减小了。

说明此时我们已经爬到山顶了。

得到的5个参数如下:

$p_1=1.3738627e+32$
$p_2=2.4369393e+28$
$p_3=16.719432$
$p_4=0.0012054476$
$p_5=3.0453779$

换回楼主的式子,5个参数如下:

$p_1=0.00017737867$
$p_2=7.2787479e-33$
$p_3=16.719432$
$p_4=0.0012054476$
$p_5=3.0453779$

函数图像如下:

DataFit.PNG

看起来和6#的很像。

标准差为$104.8$

相关系数为$0.998350$

wayne能否把坐标换回去,看看我们的参数是否一致?

#####

明明几分钟就能完成的工作,被我拖了两天才完成,真是太不应该了。

经过这次拟合试验,我吸取教训了。

下次绝对不会再犯这种低级错误了。

评分

参与人数 1威望 +6 鲜花 +8 收起 理由
wayne + 6 + 8 感谢!!

查看全部评分

毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2010-2-6 12:55:06 | 显示全部楼层
,完全一致,你是怎么做出来的,用啥软件了,可否分享一下经验


wayne KeyTo9_Fans
p1 0.0001773786744870752` 0.00017737867
p2 7.27870095803893*10^-33 7.2787479 *10^-33
p3 16.719433196419864` 16.719432
p4 0.0012054477426615047` 0.0012054476
p5 3.045377916319412` 3.0453779


毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

小黑屋|手机版|数学研发网 ( 苏ICP备07505100号 )

GMT+8, 2024-5-13 03:23 , Processed in 0.071689 second(s), 20 queries .

Powered by Discuz! X3.5

© 2001-2024 Discuz! Team.

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