立刻注册计量论坛 交流工作中的点滴

您需要 登录 才可以下载或查看,没有帐号?立即注册

x

本帖最后由 史锦顺 于 2017-4-21 12:11 编辑





                              
最小二乘法拟合公式的有关分析  


                                                                                   史锦顺




1 线性拟合公式的应用
       因变量Y与自变量X间,大体有线性关系。

       用实验(实际测量),来确定Y对X的函数关系。测量N次,数据为:

                 X1  X2 …… Xi …… XN

                 Y1  Y2 …… Yi …… YN                                            (1)

       自变量X取准,误差可略,认为无误差。因变量Y是仪器的测得值,有误差。Yi与Xi一一对应。

       拟合得到的线性公式为:

                  Y = a + bX                                                         (2)

       所谓“最小二乘法拟合”,就是依数据列(1)求得公式(2)。或者说是求得直线方程(2)的斜率b和截距a。



       时频界,求频标的漂移率(晶振称老化率)普遍用最小二乘法。

       “自变量”是时间t,ti =iT, T是时间单位(日或半日),误差可略。“因变量”是频标的频率值fi,是仪器测得值,有误差。

       拟合得频率值随时间变化的线性公式为

                  f = fo + Kt                                                          (3)

       将(3)式减去参考频率f标,得频差关系式为:

                  f – f标 = fo – f标 +Kt

                  Δf =Δfo + Kt                                                        (4)

       为方便,通常表示为相对频差

                  Δf /fo =Δfo /fo+ Kt/fo

                  δf = δfo + kt                                                        (5)



       要分清两个不同的问题。

       第一个问题是:已知频率变化率k、初始频差偏差δfo 的条件下,怎样用公式(5)求得总偏差δf。

       (5)式是个简单的代数式。δfo与kt都是偏差,它们二者合成,就是代数运算,二者都是有符号、有数值的量,按“代数计算法则”计算就是了,搞不确定度合成,是看错了对象,画蛇添足,瞎胡算。

       例1 已知晶振甲的老化率k= -3×10-10。2017年4月1日,用铯原子钟测得晶振甲的相对频差是 – 1.0×10-8,问到2017年7月10日,长期加电的晶振甲的频差是多少?

       解答:k= -3×10-8/日,δfo= – 1.0×10-8,t=100日,代入(5)式,得

                δf = – 1.0×10-8 + (-3×10-10/日)×100日

                    = – 1.0×10-8 – 3×10-8

                    = – 4×10-8



       例2 对晶振甲,已知老化率k= -3×10-10,参照检定规程的频率调整法。2017年4月1日,用铯原子钟调整晶振甲的相对频差是 +1.0×10-8,问到2017年7月10日,长期加电的晶振甲的频差是多少?

       解答:k= -3×10-10/日,δfo= +1.0×10-8,t=100日,代入(5)式,得

                δf = +1.0×10-8 + (-3×10-10/日)×100日

                    = +1.0×10-8  – 3×10-8

                    = – 2×10-8



       用最小二乘法,费力气计算得到的公式(5),上两例的应用是正确的。不确定度体系对公式(5)的处理,用什么“不确定度传播律”,是错误的,是抹煞人类的智慧,否定知识。这是“不确定度体系是伪科学”判断的又一个证据。请大家注意一个事实:凡用GUM法处理的计算问题,几乎都是错误的。至于不确定度论说教中的“相关”“不相关”都是没用的假话。初始频差、线性变化率、时间,三者都完全是各自独立的量,该认为是“不相关”;而三者共同构成总的频差的量值,又怎能说不相关?其实初始频差、变化产生的频差,客观上的作用必定是“代数和”,那种关于“相关性”的分析是没有用的。

       如果是估计“范围”,起决定作用的是二项和平方展开式的交叉系数。这里是两项系统误差,算“范围”也必须取“绝对和”,而不能是不确定度体系的“方和根”。“假设不相关”无用;“真不相关”也没用。相关性的判别,与误差合成法无关。



       上边讲的是第一个问题:拟合公式的应用;第二个问题是,测定(5)式中两个量截距与斜率的误差有多大?初始频差的测量误差易得;而求老化率k的误差,是当今世界测量计量界的一大难题。叶德培先生在她的样板评定中,缺失了(或者说弄错了)。

       老史先将老化率k的公式简化;由此方见端倪。简述如下。




2 老化率公式的简化
       简化公式的主要技巧是对称编号。

       取测量次数N为奇数,中间数为0,则上有(N-1)/2个号,下有(N-1)/2个号。

       例如

       1)测量晶振的频率老化率(线性漂移率),测量7天。每天定时(例如9:00)测量。每日一个数据(三个频率的平均值),对称编号就是f(-3)、f(-2) 、f(-1)、f(0)、f(+1)、f(+2)、f(+3)。

       2)测量晶振的频率老化率。前7天。每天定时两次测量(例如9:00、21:00)。每日两个数据(每个数据是三次测量的平均值),第8天9:00得1数据。对称编号就是f(-7)、f(-6) 、f(-5)、f(-4)、f(-3)、f(-2)、f(-1)、f(0)、f(+1) 、f(+2)、f(+3)、f(+4)、f(+5) 、f(+6)、f(+7)。



       要拟合的公式是

                 f = fo + Kt                                                          (3)

       其中,fo是初始频率,K是线性变化率。fo易获知;主要作业是拟合直线斜率K。史锦顺得到的简化公式(参见附件)为:

                K = ∑j=(-n)→n (fj –fr) j / [(N-1)N(N+1)T/12]            (6)

       简化公式的主要技巧是对称编号。  



【晶振日老化率速算法】

(1) 7天,每天1值

       按时间顺序标记数据:中间数标0,由0分界,顺时序标1至3,逆时序标-1至-3。各数据减一常数之后的尾数乘标号,乘积累加,除以28,即得日老化率。

       口诀:
对称编号,去整作差,号乘差累加,除以二十八。  

(2)7周天,每天2点,共15点(第8天测一个点)

       按时间顺序标记数据:中间数标0,由0分界,顺时序标1至7,逆时序标-1至-7。各数据减一常数之后的尾数乘标号,乘积累加,除以140,即得日老化率。

       口诀:
对称编号,去整留零,累加号乘零,除以140。



3 拟合误差
       拟合误差指求截距与斜率的误差。

       斜率K的简化(严格式)表达为:

                  K = ∑j=(-n)→n (fj –fr)j / [(1/12)(N-1)N(N+1)T]            (6)

       (6)式中,fr可以随意取值,因此K的测量误差与测量fj时的系统误差的恒定值无关。测量晶振老化率用原子频标,其系统误差为恒定值(无频率漂移)。因此误差来自原子频标与比对器的随机误差。随机误差取“方和根”。

       令fr=0,有

                  K =∑j=(-n)→n fj j / [(1/12)(N-1)N(N+1)T]                    (7)

       将fj 表达为:

                  fj = fjo + ξ                                                               (8)

       (8)式代入(7)

                  K =∑j=(-n)→n (fjo+ξ)j / [(1/12)(N-1)N(N+1)T]

                     =∑j=(-n)→n fjo j / [(1/12)(N-1)N(N+1)T]

                        +∑j=(-n)→n ξ j / [(1/12)(N-1)N(N+1)T]                  (9)

       由(9)知,K 的误差为

                  ΔK=∑j=(-n)→n ξ j / [(1/12)(N-1)N(N+1)T]

                  Kσ=∑j=(-n)→n σ j / [(1/12)(N-1)N(N+1)T]                     (10)

       测量用同一原子频标与比对器。各次测量,随机误差相同,因此σ可提出来。分子成为j的平方求和,j从1到n,结果乘2.

                  Kσ2=2σ2(12+22+32……+n2)/分母2

                       = 2 (1/6) n(n+1)(2n+1) σ 2 /分母2 (查数学手册得知自然数平方之和)

                       = (1/12)(N-1)N(N+1) σ 2 / 分母2



                  Kσ = σ√[(1/12)(N-1)N(N+1)] / 分母

                      = σ√[(1/12)(N-1)N(N+1)] / [(1/12)(N-1)N(N+1)T]

                  Kσ =(σ/ T) /√[(1/12)(N-1)N(N+1)]                            (11)



       原子频标的随机误差σ到老化率的传递系数,Kσ /σ

       取样间隔时间:日

               N=3    1/√2

               N=5    1/√10

               N=7    1/√28

       取样间隔时间:0.5日(12小时)

               N=15   1/√70





附件 论最小二乘法




论最小二乘法.doc (292 KB, 下载次数: 31)

2017-4-21 11:53 上传

点击文件名下载附件

下载积分: 金币 -1

相关新闻

联系我们

热 线 0755-27784155

电 话180-2695-0976

在线咨询点击这里给我发消息

在线报价

邮件:1981004853@qq.com

工作时间:周一至周五,8:30-18:00

联系微信
联系微信
分享本页
返回顶部