stanbul 2004 
wy! exists; 
es of p in the 
p, containing 
ans of some 
such a way to 
(3) 
EF Ni] 
coefficients. 
ive problem: 
1red points. 
lly arises: 
(4) 
models: they 
estimation of 
(5) 
neighbours). 
value can be 
he stochastic 
global effect 
surrounding 
OUTLIERS 
th dependent 
aches. 
omputational 
ures, requires 
ur method, 
rminants and 
n. For this 
ata has to be 
; each. 
al maximum 
ion processes 
difference 
1d the global 
their failure. 
d: 
In order to deal with the last requirement, in the following 
section, the so-called Maximum Likelihood (ML) method for 
estimate unknown parameters within SAR model is presented. 
3.1 Maximum Likelihood computations 
For our purposes, a maximum likelihood approach for the 
estimate of unknown SAR parameters has been chosen. Let us 
start from the general likelihood function: 
((6.0,0?) » 2no^) "|I - pW| expt - A9) z(z- ao) (6) 
where: 
Z-(1-pW)' (1- pW) 
is the weight matrix, symmetric and positive definite; 
unfortunately and differently with respect to ordinary 
estimations, here the weight matrix contains an unknown term 
as p. It is then necessary to maximize (6) not only with respect 
to 0 and c? , but also with respect to p. 
It can be performed in stages (Pace, Barry and Sirmans, 1998) 
by selecting a vector of length fof values over [0,1] labelled as: 
pv=l1 P2 Pr] 
and considering the log-likelihood function: 
n 
L(8,5,0?) - nfl pW] (5) In(eo eg —2peg eg *p^eq eq) 
where: 
e eg are the residuals from an Ordinary Least Squares (OLS) 
regression of zon A, 
e ey are the residuals from an OLS regression of Wz on A. 
Thus, to maximize (6), the following m terms are evaluated: 
18,007) | | I-pW Ineg eo —2pieq eo +pieq ea) 
L(6,p,,0°) = InI-p,W 4 In(egleg —2p5e, e +phea"e4) (7) 
Loo) |Hll-prW|] ^ [Ine "es -2prea eo *prea ea) 
and the value py giving the maximum log-likelihood value L 
among those in (7) is assumed as the ML estimation p of p. 
The use of a finite set of p will cause some small granularity in 
the chosen values pm, , but it should not prove difficult to 
make the granularity small relative to the statistical precision of 
estimated py, . While this approach may suffer a small loss of 
precision relative to non-linear maximization, the evaluation of 
log-likelihood function over an interval ensures robustness, the 
main property of this approach (Griffith and Layne, 1999). 
Once so obtained p, it is possible to ML estimate SAR 
unknowns and also a new weight (optimised) matrix $ 
0-(ATEA)!ATzz (8.1) 
6? - n^! (z - A8)! Z(z - A0) (8.2) 
E-(I-9w)! a - 8W) (8.3) 
     
  
  
  
  
  
   
   
  
  
   
     
   
   
    
   
   
    
  
   
  
    
   
    
   
   
  
   
   
    
   
    
  
  
    
   
    
   
    
   
   
   
    
   
   
   
   
  
  
    
International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Vol XXXV, Part B3. Istanbul 2004 
3.2 Spatial outliers searching 
A spatial outlier is defined as *an observation which is unusual 
with respect to its neighbouring values" (Haining, 1990). For 
our purposes, the way to assess spatial outlyiness is to compute 
individual departure from the fitted polynomial trend surface. 
To accomplish this goal, starting from (4), the vector e = cle 
of standardised residuals is estimated as: 
e-6^!(I- 6W)(z - A8) (9) 
where p, 6,6 are the unknown parameters of the SAR model 
just simultaneously estimated by means of (7) and (8). 
Afterwards, its n components are inferentially evaluated to find 
which measures do not fit the estimated surface: vector e 
defines in fact the lack of fit statistic ele. 
Standardised residuals e over residuals £ have been thus 
preferred, since they allow a robust spatial autocorrelation 
estimation, which we believe is a sensible property for the 
purpose of detecting spatial outliers. 
From the methodological point of view, the main property of 
our algorithm allowing to detect LIDAR outliers is to perform 
estimations (7), (8) and (9) on different subset of the whole 
dataset. In particular, we start from an initial subset of LIDAR 
data up to include all the dataset of the sub-zone to be filtered. 
4. IMPLEMENTATION OF AN ITERATIVE 
SEARCHING PROCEDURE 
4.1 Block Forward Search for SAR models 
An interesting algorithm to perform iterative SAR estimations 
on increasing datasets is the so-called “Block Forward Search” 
(BFS) proposed by Atkinson and Riani (2000) and Cerioli and 
Riani (2003) for econometric purposes. 
It makes possible to proceed to the joint robust estimation p 
and 6 at each step of the search, starting from a partition of 
datasets in n blocks of contiguous spatial location, and 
considers these blocks as its elementary unit. In the general case 
of grid data, each block is a set of cells, while handling raw data 
is difficult to univocally create the blocks and so the block 
dimension is merely unitary (UFS, Unitary Forward Search). 
The basic idea of the BFS approach is to repeatedly fit the 
postulated model p — A8 to increasing subsets size, selecting 
for any new iteration the observations z best fitting the previous 
subset, that is having the minimum standardised residuals e 
computed by (9). It must to be stressed as in equation (9), p 
and 6 are estimated on the subset outlier free only, while z, A, 
W and c terms are referred to the whole dataset with outliers. 
Thanks to the strategy of block growing, the outliers present 
into z values are included only at the end of the BFS procedure. 
4.2 SFS implemented algorithm 
The proposed algorithm, called simply SFS (Spatial Forward 
Search), implements the Forward Search approach, but, since 
raw LIDAR data are irregularly located, unitary blocks to 
increase the subset size were chosen: in other words, only one 
measured point enlarges the subset at each iteration. 
The SFS algorithm has been implemented as a software tool 
using Matlab® language. Its main steps are (see Figure 1):