157
1
INTRODUCTION
Parameterestimationandstatisticaltestingarecrucial
components in any navigational data processing
system. The DIA method for the Detection,
Identification and Adaptation of model
misspecifications combines estimation with testing.
Parameter estimation is conducted to determine
estimatesfortheparametersofinterestandstatistical
testing is conducted to validate the
results with the
aim of removing any unwanted biases that may be
present.Inordertorigorouslycapturetheestimation‐
testing combination, the DIA‐estimator has been
introduced in [1], together with a unifying
probabilisticframework.Thisallowsonetotakeinto
account the intricacies of the combination when
evaluating
the contributions of the decisions and
estimators involved. As it will be shown, this
inclusion of the estimation+testing combination is
crucial to be able to obtain a rigorous quality and
integrity‐risk description of one’s position and
navigational results. Such rigor is needed so as to
avoid too optimistic, and therefore
dangerous,
navigationalintegrityassessments.
2
THEDIA‐METHOD
AbriefreviewisgivenoftheDetection‐Identification‐
Adaptation (DIA) method. Although the method is
applicableforanystatisticalidentificationprocedure,
we here use for illustration purposes the relatively
simple datasnooping procedure for single‐outlier
identification[1].
2.1
StatisticalHypotheses
We first formulate the null‐ and alternative
hypotheses, denoted by H
0 and Hi, respectively. As
alternative hypotheses,we consider those describing
outliers in individual observations. Here we restrict
ourselvestothecaseofoneoutlieratatime.Inthat
casethereareasmanyalternativehypothesesasthere
are observations. Therefore, the observationalmodel
underH
0andHiisgivenas
:() ()
0 yy
E y = A x; D y = QH
(1)
:() ()
iiiyy
E y = A x + c b ; D y = QH
(2)
The DIA-Method for Navigational Integrit
y
P.J.G.Teunissen
1,2
1
DelftUniversityofTechnology,Delft,TheNetherlands
2
CurtinUniversity,Perth,Australia
ABSTRACT:InthiscontributionwepresentareviewoftheDIA‐methodtoensurenavigationalintegrity.The
DIA‐methodrigorouslycombinesparameterestimationandstatisticaltestingfortheDetection,Identification
andAdaptationofmultivariateandmultiplemodelmisspecifications.Wedescribethestatisticalpropertiesof
theso‐
obtainedDIA‐estimatortogetherwithitsprobabilitydensityfunction.Numericalexamplesaregivento
highlightvariousaspectsofthenavigationalDIA‐estimator.
http://www.transnav.eu
the International Journal
on Marine Navigation
and Safety of Sea Transportation
Volume 17
Number 1
March 2023
DOI:10.12716/1001.17.01.16
158
with
()E
the expectation operator,
()
D
the
dispersionoperator,
m
y
thenormallydistributed
randomvectorofobservableslinkedtotheestimable
unknown parameters through the design matrix
m n
A
of
rank( )
A
=n
, and
m m
yy
Q
the
positive‐definite variance matrix of y . The
redundancy of
0
H
is
rank( )r=m- A =m-n
.
i
c
is the
canonical unit vector having one as its i
th
entry and
zeros elsewhere, and b
i is the scalar bias. Note that
[
i
A
c 󰇠
isaknownmatrixoffullrank.Asthenumber
of observations is equal to m, there are also m
alternativehypotheses
i
H
definedin(2); i=1, ,m.
The best linearunbiased estimator (BLUE)of the
unknownparametersxisgivenby
00
::
++
iii
x = A y x = A y
0
ˆˆ
HH
(3)
with
-1 -1 -1
()
+T T
yy yy
A
=AQA AQ
the BLUE‐inverse of A,
-1 -1 -1
()
TT
iiyyiiyy
=AQA AQ
the BLUE‐inverse of
i
ic
A
=P A
and
-1 -1 -1
()
i
TT
cmiiyyiiyy
P=I-ccQc cQ
being an orthogonal
projector that projects onto the orthogonal
complementoftherangespaceofc
i.
2.2
DIA‐datasnoopingprocedure
The DIA‐method has been widely employed in a
varietyofapplications,suchasthequalitycontrolof
positioning and navigation, and the integrity
monitoringofGNSSmodels,seee.g.[2,3].TheDIA
steps are realized using the misclosure vector
r
t
givenas
TT
tt yy
t = B y; Q = B Q B
(4)
wherethe
mr
matrixBisa full‐rankmatrix,with
rank( )
B
=r
,ofwhichtherangespaceisanorthogonal
complement of that of A, i.e.
[
m m
A B 󰇠 is
invertible and
0
T
A
B= . With
()
i
ii yy
y
Ax+c b ,Q~
H
N
for
01i= , ,…,m and
0cb =
00
,themisclosurevectoristhen
distributedas
( ) for 0 1
i
i
T
tiiyy
t
μ
=B cb,Q i= , ,…,m~
H
N
(5)
Astiszero‐meanunder
0
H
andalsoindependent
of
0
x
ˆ
,itprovidesalltheavailableinformationuseful
for validation of
0
H
[1]. Thus, an unambiguous
testing procedure can be established through
assigning the outcomes of t to the statistical
hypotheses
i
H
for 01i= , ,…,m.
The DIA‐datasnooping procedure is specified as
follows:
1.
Detection:Accept
0
H
if
0
t P
with

2
0

tt
r
Q α,r
=t ||t|| kP
(6)
in which
2-1
=() ()
tt
T
Qtt
|| || Q
and
α,r
k
is the
percentage of the central Chi‐square distribution
withrdegreesoffreedom.If
0
H
isaccepted,then
0
x
ˆ
isprovidedastheestimateofx.Otherwise,go
tostep2.
2.
Identification:ComputeBaardaʹsteststatisticforall
alternativehypothesesas[4,5]
-1
-1
i
ii
T
ttt
i
T
tttt
cQt
w=
cQc
(7)
inwhich
i
T
ti
c=Bc
isthei
th
columnofmatrix
T
B
since ci is a canonical unit vector. Select
i 0
H
if
0
i
t P
with
00
{1 }
||=max| |
r
iij
j,…,m
=t w w/PP
(8)
3.
Adaptation: When
i
H
is selected, then
i
x
ˆ
is
providedastheestimateofx.
The partitioning P
i in terms of the (original)
misclosurevectorisintroducedin[1],andanexample
isshowninFigure3in[ibid].
3
THEDIA‐ESTIMATOR
Given the above three steps, estimation and testing
are combined in DIA‐datasnooping. A unifying
framework is presented in [1] to rigorously capture
theprobabilisticpropertiesofthiscombination.
Figure1.EstimationandtestingcombinedleadstotheDIA‐
estimator[6].
As such, the DIA‐estimator
x
was introduced,
which captures the whole estimation‐testing scheme
anditisgivenas(seeFigure1)
0
()
m
jj
j=
x
= x pt
ˆ
(9)
with p
j(t)being theindicator functionof region
j
P
,
i.e.
()=1
j
pt
for
j
t P
and
()=0
j
pt
elsewhere.
Therefore, the DIA‐estimator
x
is a combination of
j
x
ˆ
for
01j= , ,…,m
and the misclosure vector t. The
Probability Density Function (PDF) of
x
under
i
H
reads[1]
0
0
() ( )
()()
j
j
0
j
m
ˆ
xi x,t i
j=
m
ˆ
xjiti
j=
f θ |= fθ, τ | dτ
= f θ +L τ | fτ | dτ
P
P
HH
HH
(10)
wherethesecondequalityisaconsequenceof
159
0
j
j
ˆˆ
x
=x -L t
(11)
with
-1
2
1
j
jtt
+T
j
jt tt
t Q
L= AccQ
|| c ||
(12)
Itisimportanttorealize,althoughweareworking
withlinearmodels(cf.1and2)andlinearestimators
j
x
ˆ
based normally distributed data, that the DIA‐
estimator
x
,whichcapturesthefinaloutcomeofour
combined estimation+testing process, is now not
normally distributed anymore. Its distribution is
givenby(10)anditisthisexpressionthatoneneeds
tousewhenevaluatingthequalityorintegrityrisksof
one’spositioningornavigationalresults.
4
THEDIA‐ESTIMATORINCASEOFASINGLE
ALTERNATIVEHYPOTHESIS
To illustrate thecharacteristicsof the DIA‐estimator,
we take the single alternative hypothesis case as a
simpleexample.Supposethatin(1),thereisonlyone
unknownparameter (n=1)andalsotheredundancy
ofthemodelis
one(r=1),i.e.
x
and
t
.The
canonical form of such a model, applying the
Tienstra‐transformation
T
to the normally
distributedvectorofobservablesy[1],reads
0
0
2
0
2
0
0



















i
,i
i
ˆ
x
ˆ
x
T
t
t
xb
ˆ
x
A
y= y ~ ,
t
b
B
H
TN
(13)
whichisspecifiedfor
{0 }i,a
as
00 0
0
0
: 0 0
:


,
,a a
ˆ
xt
ˆ
ax aata
b,b
bL b,bb
H
H
(14)
for some
{0}
a
b/
, and also
a
L
which
establishesthefollowinglink
0

aa
ˆˆ
x
xL t
(15)
sothat
0
()()
aaa
ˆˆ
Ex| Ex| xHH
.
ThecorrespondingDIA‐datasnoopingprocedureis
thendefinedas:
1.
Detection:Accept
0
H
if
0
t P
with
0
[ , ]
α,1 α,1
=k kP
(16)
Provide
0
x
ˆ
astheestimateofx.
2.
Identification:Select
a
H
if
c
t
0
P
with
c
00
/PP
.
3.
Adaptation: When
a
H
is selected,
a
x
ˆ
is provided
astheestimateofx.
Withtheabovethreesteps,theDIA‐estimatorand
itsPDFunder
i
H
,
{0 }i,a
,aregivenby
00 0
() (1 ())
a
ˆˆ
x
x pt x pt
(17)
and
0
00
0
() ()
()()()




c
ˆ
xix i
ˆˆ
xaixiti
f θ |= fθ |
fθ +L | f θ | f|dτ
P
HH
HHH
(18)
As there is only one alternative hypothesis
a
H
,
therearefoureventstoconsider:CorrectAcceptance
(CA),False Alarm(FA),Missed Detection(MD) and
Correct Detection (CD). Using their probability of
occurrences,the PDFof the DIA‐estimator
()
x
i
f θ |H
canbedecomposedas
0
0
00CAFAFA
MD CD CD
( ) ( ) P ( FA) P
()() P (CD) P
a
a
ˆˆ
xx x|
ˆˆ
xax a x|
f θ |= fθ | + fθ |
f θ |= fθ | + fθ |
HH
HH
(19)
In Figure 2, considering
0
22
05m
ˆ
x
.
,
22
2m
t
and
05
a
L.
,weshowhowthePDFs
0
()
x
f θ |H
[top]
and
()
x
a
f θ |H
[middle and bottom] are formed
according to(19). The solid and dashed blue curves,
respectively, depict
0
0CA
()=(CA)
ˆ
xx|
f θ |fθ |H
and
FA
(FA)
a
ˆ
x|
f θ |
in the top panel, and
0
MD
() (MD)
ˆ
xax|
f θ |= f θ |H
and
CD
(CD)
a
ˆ
x|
f θ |
in the
middle and bottom panels. The black curve shows
0
()
a
ˆ
x
f θ | H
whichisalsoequalto
()
a
ˆ
x
a
f θ | H
.
These results clearly show how the PDF of the
DIA‐estimatordiffersfromthePDFsof
0
x
ˆ
and
a
x
ˆ
.
5
DIACONFIDENCEREGIONINCASEOF
MULTIPLEALTERNATIVEHYPOTHSES
Sofar,we havebeenworkingwith anobservational
model with one unknown parameter and one
redundancy. In this section, we work with the
satellite‐based single point positioning (SPP) model
based on the observations of m satellites with four
unknown parameters
(n=4) and 4
rm
redundancy. As alternative hypotheses, we consider
those given in (2). In that case there are as many
alternativehypothesesasthereareobservations.
Assuming there arem pseudorange observations,
theobservationalmodelunder
i
H
for 01
i,,,m is
givenas
2
: ( ) = [ ]



imiiyypm
x
E y G e c b , Q = σ I
dt
H
(20)
where the
3
m
matrix
1
[]
TTT
m
G = - u , … , - u contains
the receiver‐satellite unit direction vectors u
i as its
rows,ande
misthem–vectorofones,andagainwith
00
0
cb
. The unknown receiver coordinate
componentsandclockerrorare,respectively,denoted
bythe3‐vectorx andscalardt.Thedispersionofthe
observables is characterized through the standard
deviation
p
and the identity matrix
m
I
. At this
stage, in order to simplify our analysis, we do not
consider a satellite elevation‐dependent variance
matrix.
160
Figure2.IllustrationoftheDIA‐estimatorPDF[top]under
null hypothesis
0
H and [middle and bottom] under
alternativehypothesis
a
H
.Thecontributingparametersare
setas
0
22
05m
ˆ
x
.
,
22
2m
t
,
05
a
L
.
.Panels,fromleftto
right, correspond to
3
FA
P10
and
1
FA
P10
. Panels in
themiddleandatthebottomcorrespondtoresp.
CD
P04
.
and
CD
P099 .
.
For integrity purposes, it is common to describe
the quality of an estimator of x by computing its
corresponding confidence region or confidence level
[7]. To construct the confidence region for the
parametersofinterest,onemaydefinebeforehandthe
shapeoftheregion,forexample,spherical,ellipsoidal,
orrectangular,and
thencomputeitssizeforagiven
confidencelevel[8,9],orhavetheconfidenceregion
determined by the contours of the PDF of the
estimator for a given confidence level [10, 11].
Whichever approach is taken, the properties of the
confidenceregionaredeterminedbytheprobabilistic
properties of
the estimator. It is thereforecrucial, in
order to have a realistic confidence region, that one
works with the correct PDF of one’s estimator. This
implies,whenexercisingaDIAprocedure,oneshould
usethePDFoftheDIAestimator
x
ratherthanthe
PDF of its constituent estimators
i
x
ˆ
. To emphasize
this, we compare three different confidence regions
undernull‐hypothesisasfollows:
1.
An ellipsoidal 100(1‐α)% confidence region E(x)
basedonthenormaldistributionof
0
x
ˆ
under
0
H
.
This is the region that one normally uses in
practiceforconstructingaconfidenceregion.Note
however,thatthisregionneglectstheprobabilistic
impactoftesting.
2.
An ellipsoidal 100(1‐α)% confidence region E’(x),
that has the same shape as E(x), but with a size
determined by the nonnormal PDF of the DIA‐
estimatorunderthenull‐hypothesis.
3.
A 100(1‐α)% confidence region C(x) of which the
shapeandsizearedeterminedasthehighestdensity
region[10,11]ofthenonnormalPDFoftheDIA‐
estimatorunderthenull‐hypothesis.
Figure3showstheabovethreeconfidenceregions
fortheSPPmodel(20)correspondingwiththeFigure
4
skyplot.Thelargersizeoftheredellipsecomparing
to the blue one demonstrates that the confidence
region E(x) has a poor coverage in the sense that it
providesatoooptimisticpictureandthatitshouldthus
be made larger in order to contain the required
probabilityofthe
DIAestimator.Itwillbeclearthata
too optimistic assessment is dangerous in case of
safety‐of‐life applications. The green region C(x) in
Fig. 4 is significantly different in shape from the
ellipsoidalconfidenceregionwhichisconventionally
used. Depending on the α value, C(x) can be a
nonconvex
region,whichiscausedbytheshapeofthe
nonnormalPDFoftheDIA‐estimattorunderthenull‐
hypothesis.Thisshowsagaintheinfluencetestinghas
on confidence‐statements and that one therefore
shouldusethePDFoftheDIA‐estimatortoevaluate
the quality of estimated position and navigation
results.
Figure3.Illustrationof100(1‐α)%confidenceregionsforthe
SPPmodel(20)correspondingwiththeskyplotofFig4.The
resultsshowdifferentvaluesofα,withσp=0.5mandPFA=0.1.
In each panel, green area indicates the region C(x), red
ellipseindicatesboundaryofE’(x)andblueellipseindicates
theboundaryofE(x)[6].
Figure4:Skyplotofsixsatellites.
161
6
SUMMARYANDCONCLUSION
Themessageofthiscontributionfindsitsorigininthe
combination of parameter estimation and statistical
testing.Thesetwoactivitiesaretypicallydisconnected
inpractice,whenitcomestodescribingthequalityof
theeventualestimator.Thatis,thedistributionofthe
estimatorunderanidentifiedalternativehypothesis
is
usually used, without regard to the conditioning
process that led to the decision to accept this
hypothesisastheworkingmodel.Weanalysedwhat
the contribution of this negligence is to the actual
integrity risk and showed that for a rigorous
assessmentoneshouldusetheDIA‐estimatorand
its
probabilitydensityfunction.
REFERENCES
[1]TeunissenPJG (2017) Distributional theory for the DIA
method. Journal of Geodesy, Springer, 1–22, DOI
10.1007/s00190‐017‐1045‐7.
[2]Teunissen PJG (1990) An integrity and quality control
procedureforuseinmulti‐sensorintegration.In:Proc.of
IONGPS‐1990,ION,pp513–522
[3]Amiri Simkooei A (2001) Comparisonof reliabilityand
geometrical strength criteria in geodetic networks.
JournalofGeodesy75(4):227–233
[4]Baarda W (1967) Statistical concepts in geodesy.
Netherlands Geodetic Commission, Publ. on geodesy,
Newseries2(4)
[5]Teunissen PJG (2000) Testing theory: an introduction.
DelftUniversityPress,SeriesonMathematicalGeodesy
andPositioning
[6]Zaminpardaz S,
Teunissen PJG (2022): GNSS Detection
andEstimation.EncyclopediaofGeodesy,Springer,MG
Sideris(Ed.),1‐9.
[7]Koch KR (1999): Parameter estimation and hypothesis
testing in linear models. Springer Science & Business
Media.
[8]HooverW(1984):Algorithmsforconfidencecirclesand
ellipses (NOAA Technical Report NOS 107 C&GS 3).
Rockville,
MD: National Oceanic and Atmospheric
Administration.
[9]Hofmann‐Wellenhof B, Legat K, Wieser M (2003):
Navigation: Principles of Positioning and Guidance.
Springer.
[10]HyndmanRJ (1996): Computing and graphing highest
density regions. The American Statistician 50(2): 120‐
126.
[11]TeunissenPJG(2007):Least‐squarespredictioninlinear
modelswithinteger
unknowns.JournalofGeodesy81:
565‐579.