当前位置:网站首页>[GNSS] robust estimation (robust estimation) principle and program implementation
[GNSS] robust estimation (robust estimation) principle and program implementation
2022-07-06 12:53:00 【Proletarians】
Robust estimation
The principle of robust estimation
Robust estimation is a category of modern measurement adjustment , Also known as robust estimation (robust estimate), According to academician Yang, the CAS system likes to call it robust estimation , Wu Da likes to call it robust estimation . Our measurements are random variables , According to the normal distribution , If there is gross error (gross error) Words , We will make the result deviate from the true value when we apply the least squares or Kalman filter ( Filter divergence ) The phenomenon of . When we solve gross errors or systematic errors , It can be understood from two aspects , Mean shift or variance inflation , Robust estimation belongs to variance inflation model , That is, the mean value remains unchanged , The phenomenon of variance change . We can reduce the weight of observations with gross errors .
I compare lazy , I don't want to write a formula , You should have seen a similar derivation , In fact, the formula is the same as the minimum second difference , The place of change is the power matrix , Replaced by equivalent weight , Then what affects the equivalent weight is the equivalent weight function , Common are huber、IGG III etc. , I recommend IGG III, The code is easy to implement .
When we do it , This power matrix is very important , The weight matrix is the inverse of what we call the variance matrix , It is divided into independent and non independent , That is, whether the observed values are related , Whether it is relevant also affects the form of the equivalent weight function . Take the simplest example of independent observation distance
Examples of robust estimation
Let's assume that for a distance 10 Sub independent observation ,10 The observed values are :5.09、5.10、5.13、5.09、5.12、5.08、5.46、7.81、5.10、5.11( Unit is m), I made up the data casually .
From the above paragraph , We can know the number of measurements n=10, Necessary observation m=1, The redundant observation is n-m=9. Independent observation, then the weight matrix is diagonal , Prior weight matrix is unit diagonal matrix . And what is visible to the naked eye is 7.81 This observation is a gross error ,5.46 Well, I'm not sure if .
Robust program implementation
matlab Realized
%% function 3- Robust estimation
% % Suppose a length is observed ,5.09 5.10 5.13 5.09 5.12 5.08 5.46 7.81 5.10 5.11
close all
clear all
clc
k0=1.0;k1=2.5;k=1.0;
B=ones(10,1);% Design matrix
l=[5.09 5.10 5.13 5.09 5.12 5.08 6.46 7.81 5.10 5.11]';
% P=diag(ones(10,1));% A priori weight matrix 1, Equivalent weight
P=diag(1./l);% A priori weight matrix 2, Weight according to the reciprocal of length
x_prev=0;
x0=5.10;% Assign initial value to
l=l-x0;
while(1)
W=B'*P*l;
N=B'*P*B;
x=inv(N)*W;% Parameter vector to be evaluated
v=B*x-l;% Residual vector
sigma0=sqrt(v'*P*v/(10-1));% Mean square error of unit weight
disp([' Parameter vector to be evaluated :',num2str(x),' Mean square error of unit weight :',num2str(sigma0)]);
Q=inv(P);
Qvv=Q-B*inv(N)*B';
for i=1:10
v_=v(i)/(sigma0*sqrt(Qvv(i,i)));
if abs(v_)<=k0
k=1.0;
elseif abs(v_)>k1
k=1e-8;
else
k=(k0/abs(v_))*((k1-abs(v_))/(k1-k0))^2;
end
P(i,i)=P(i,i)*k;
end
if x_prev==0
x_prev=x;
continue;
elseif abs(x_prev-x)<0.01
break;
end
x_prev=x;
end
x=x0+x_prev;% Find the final length initial value + Correct number
disp(['length is ',num2str(x)])
The results are as follows :
The result of using equivalent weight is :
Parameter vector to be evaluated :0.309 Mean square error of unit weight :0.8512
Parameter vector to be evaluated :0.042222 Mean square error of unit weight :0.11331
Parameter vector to be evaluated :0.0025 Mean square error of unit weight :0.01472
Parameter vector to be evaluated :-0.00031553 Mean square error of unit weight :0.0084525
length is 5.1025
Using the reciprocal of length to determine the weight, the result is :
Parameter vector to be evaluated :0.2218 Mean square error of unit weight :0.31128
Parameter vector to be evaluated :0.03985 Mean square error of unit weight :0.048702
Parameter vector to be evaluated :0.0024523 Mean square error of unit weight :0.0065128
Parameter vector to be evaluated :-0.0003396 Mean square error of unit weight :0.0037381
length is 5.1025
The first line in the figure is the result of our least squares ,10 Only one of the results is close to the value obtained by the least square , Obviously, it is inconsistent with our measured results , Therefore, the existence of a gross error will have a devastating effect on the result of the least squares . So if we apply robust estimation , after 3 The final result and unit weight mean square error are obtained in the second iteration .
If I'm on the seventh value , Making a change , by 5.46 become 6.46. The results are shown below :
The result of using the unit right is :
Parameter vector to be evaluated :0.409 Mean square error of unit weight :0.91426
Parameter vector to be evaluated :0.1257 Mean square error of unit weight :0.38585
Parameter vector to be evaluated :0.0025 Mean square error of unit weight :0.01472
Parameter vector to be evaluated :-0.00031551 Mean square error of unit weight :0.0084526
length is 5.1025
Using the reciprocal of length to determine the weight, the result is :
Parameter vector to be evaluated :0.2218 Mean square error of unit weight :0.31128
Parameter vector to be evaluated :0.03985 Mean square error of unit weight :0.048702
Parameter vector to be evaluated :0.0024523 Mean square error of unit weight :0.0065128
Parameter vector to be evaluated :-0.0003396 Mean square error of unit weight :0.0037381
length is 5.1025
Conclusion
The accuracy of parameter estimation is improved by using robust estimation , In addition, the results of different prior weight matrices are inconsistent in the iterative process , For this example, the reciprocal of length is used to determine the weight , Whether or not the minimum two difference is applied , The results are close to the real value , The mean square error of unit weight is also small .
Robust estimation of correlation between observations
The main thing is the three-step iteration :
1、 Standardized residuals , Mean square error of unit weight
2、 Variance amplification factor
3、 Equivalent variance
Welcome to exchange ~
Updated on 2020/11/10 14:20:00
In the above example , Use standardized residuals , Thank you for reminding me by Mr. Xiao of China Academy of testing , I haven't seen classic adjustment for a long time , The covariance matrix propagation formula is really forgotten , Ha ha ha .
Expand : In fact, we are gnss spp in , There is also no correlation between the observations , So you see , Find the corresponding parameters and apply them directly .
边栏推荐
- 2021.11.10 compilation examination
- FairyGUI增益BUFF數值改變的顯示
- [算法] 剑指offer2 golang 面试题1:整数除法
- (5) Introduction to R language bioinformatics -- ORF and sequence analysis
- [Chongqing Guangdong education] reference materials for regional analysis and planning of Pingdingshan University
- (课设第一套)1-5 317号子任务 (100 分)(Dijkstra:重边自环)
- 使用rtknavi进行RT-PPP测试
- [offer78] merge multiple ordered linked lists
- 闇の連鎖(LCA+树上差分)
- 1041 be unique (20 points (s)) (hash: find the first number that occurs once)
猜你喜欢
Affichage du changement de valeur du Buff de gain de l'interface graphique de défaillance
FairyGUI循环列表
抗差估计在rtklib的pntpos函数(标准单点定位spp)中的c代码实现
[Nodejs] 20. Koa2 onion ring model ----- code demonstration
[algorithme] swordfinger offer2 golang question d'entrevue 2: addition binaire
[算法] 剑指offer2 golang 面试题4:只出现一次的数字
341. Flatten nested list iterator
【无标题】
Liste des boucles de l'interface graphique de défaillance
[algorithm] sword finger offer2 golang interview question 10: subarray with sum K
随机推荐
FairyGUI循环列表
Special palindromes of daily practice of Blue Bridge Cup
地球围绕太阳转
Unity场景跳转及退出
NovAtel 板卡OEM617D配置步骤记录
FairyGUI簡單背包的制作
[algorithm] sword finger offer2 golang interview question 12: the sum of the left and right sub arrays is equal
Database course design: college educational administration management system (including code)
[算法] 剑指offer2 golang 面试题2:二进制加法
Mysql database index
Unity scene jump and exit
Mixed use of fairygui button dynamics
1041 be unique (20 points (s)) (hash: find the first number that occurs once)
Introduction to the daily practice column of the Blue Bridge Cup
JUC forkjoin and completable future
Fabrication of fairygui simple Backpack
Unity3D基础入门之粒子系统(属性介绍+火焰粒子系统案例制作)
Acwing-116 pilot brother
rtklib单点定位spp使用抗差估计遇到的问题及解决
[leetcode15] sum of three numbers