Adjustment of Measurements with Multiplicative Errors: Error Analysis, Estimates of the Variance of Unit Weight, and Effect on Volume Estimation from LiDAR-Type Digital Elevation Models

Modern observation technology has verified that measurement errors can be proportional to the true values of measurements such as GPS, VLBI baselines and LiDAR. Observational models of this type are called multiplicative error models. This paper is to extend the work of Xu and Shimada published in 2000 on multiplicative error models to analytical error analysis of quantities of practical interest and estimates of the variance of unit weight. We analytically derive the variance-covariance matrices of the three least squares (LS) adjustments, the adjusted measurements and the corrections of measurements in multiplicative error models. For quality evaluation, we construct five estimators for the variance of unit weight in association of the three LS adjustment methods. Although LiDAR measurements are contaminated with multiplicative random errors, LiDAR-based digital elevation models (DEM) have been constructed as if they were of additive random errors. We will simulate a model landslide, which is assumed to be surveyed with LiDAR, and investigate the effect of LiDAR-type multiplicative error measurements on DEM construction and its effect on the estimate of landslide mass volume from the constructed DEM.


Introduction
Theory and methods of adjustment have been developed, both with the advance of measurement technology and with our deepened understanding of measurement errors. Although manufacturers of surveying instruments would always provide accuracy specifications [1], these specified numbers may not correctly reflect the nature of the errors of collected measurements. To understand the performance of a surveying instrument in practice, or equivalently, to obtain the accuracy of the instrument for a specific set of measurements, Airy [2] and Helmert [3] invented the theory and methods of variance component estimation, which has since become one of the important topics in statistics, geodesy and beyond (see e.g., [4][5][6]). If measurements are supposed to contain outliers, outlier detection and robust statistics are developed (see e.g., [7][8][9][10][11][12][13]). If prior information is available, Bayesian statistics and least squares collocation are accordingly borne naturally. To precisely determine the best trajectory of an object in motion, the theory and methods of optimal filtering, smoothing and control are incepted with the celebrated publication of Kalman [14]. Recently, with the advance of space observation technology such as global positioning system (GPS) and interferometric synthetic aperture radar (InSAR) (see e.g., [15][16][17]) and their wide applications, we have to handle observational models that contain real-valued and integer unknowns, which have been coined as mixed integer (linear) models by Xu [18].
Almost all theory and methods of adjustment have been developed on the basis of the following mathematical or functional model: where is a vector of measurements, () f  is the linear or nonlinear functional,  is the vector of unknown parameters to be estimated. In the case of GPS and InSAR, the vector  may contain real-valued and integer unknowns. is the vector of random errors of the measurements . An important feature of the observational equation (1) is that the random errors will disturb the true values () f  of the measurements in the additive manner, implicitly indicating that the sizes of the random errors will not change with the sizes of measurements. In other words, has nothing to do with  or () f  . However, in geodetic practice, we realize that some random errors of measurements change with the sizes of measurements themselves. For example, the accuracy of electronic distance measurements (EDM), GPS and/or VLBI baselines can be represented by using the following well-known formula: (see e.g., [1,19,20] where stands for the random error of the baseline L, and are the additive and multiplicative random errors, which are assumed to be statistically independent, with zero means and variances a 2 and b 2 , respectively. in (3) is said to be a multiplicative error, since it will be propagated into the error of the baseline L through the baseline itself. As a result, this part of the error is proportional to the length of the baseline. The longer the baseline, the less accurate it is. However, the effect of the additive error remains constant, irrelevant to the length of the baseline itself. As a second example, a SAR observational equation can be written as follows: (see e.g., [21][22][23][24][25][26] (3) and (4) are obviously fundamentally different from (1). From this point of view, adjustment theory and methods, as developed on the basis of the model (1) (with only additive errors), cannot serve as a solid theoretical foundation to process measurements that are collected on the model (3) and/or (4).
The observational model (4) is called a multiplicative error model, which is widely known as a generalized linear model in statistics (see e.g., [27,28]). In engineering, it is also called a multiplicative speckle model (see e.g., [21,24]). The most widely used method in statistics to handle such a model is the quasi-likelihood method, which was invented by Wedderburn [27]. Quasi-likelihood has since become a standard method for parameter estimation in the model (4) and been widely applied in practice (see e.g., [29][30][31][32]). Mathematically, the method directly solves a set of nonlinear equations for the unknown parameters without any link to an objective function. As a result, Wedderburn [27] showed that the method is equivalent to maximum likelihood, if and only if the joint probability density function (pdf) of measurements is of exponential class. Unlike Wedderburn [27], Xu and Shimada [33] did not assume any pdf for the measurements with multiplicative errors. Instead, they directly started with the least squares (LS) method and derived the bias-corrected LS method for parameter estimation in the model (4). A digital elevation model (DEM) is a numerical or digital representation of the topography of the Earth. A number of mathematical interpolators have been proposed for the construction of a DEM on the (implicit) assumption that the measurements are contaminated with additive random errors (see e.g., [34][35][36]). Recently, robust methods have also been proposed to construct a DEM by allowing the measurements to contain a certain percent of erroneous data (outliers) (see e.g., [37,38]). DEMs have found important applications in hazard assessment for disaster prevention and/or mitigation, for example, in landslide analysis [39,40], in construction of hazard maps [41] and to help reconstruct the history of hydrological and glacial events (see e.g., [42]). Practically, a DEM can be constructed by using remote sensing techniques such as (In)SAR and light detection and ranging (LiDAR). Although (In)SAR and LiDAR measurements have been well known to be of speckle (or multiplicative) noise, they have been used to construct a DEM as if these types of measurements were of additive random errors (see e.g., [37,39,[43][44][45]). One of the purposes of this paper is to extend the methods of DEM construction to the case in which measurements are contaminated with multiplicative random errors and investigate their effect on the estimate of volume computed from LiDAR-type DEM, which can be important in practical hazard evaluation of landslides.
Although the errors of GPS, EDM and VLBI baselines and measurements of InSAR and LiDAR have been shown to be of multiplicative nature, almost nothing can be found in the geodetic and DEM literature, except for Xu et al. [46]. In this paper, we will limit ourselves to the model of multiplicative errors and substantially extend the work of Xu and Shimada [33]. We will not assume any pdf for the measurements with multiplicative errors. The paper is organized as follows: Section 2 will briefly outline three LS-based methods for parameter estimation. We will focus on analytical error analysis of all adjusted quantities of interest in Section 3, which further supplements the work of Xu and Shimada [33] in the sense that they mainly discussed parameter estimation and stochastic/numerical simulations. The error analysis of quantities other than model parameters is not covered by Xu et al. [46] either. In Section 4, we will derive the estimators of variance of unit weight in association with the three LS-based methods. Finally, we will simulate an example of DEM construction from LiDAR-type data contaminated with multiplicative random errors and investigate their effect on the estimate of volume from the constructed DEM.

Representation of Models with Multiplicative Errors
As in the case of adjustment with additive random errors, given a set of measurements with multiplicative errors, we have the corresponding adjustment model as follows: (see e.g., [21,33]), or equivalently in vector form: where y i is a measurement, is the functional model of the measurement y i , is the random error with mean zero. If all these scalar quantities are represented in vector form, the corresponding vectors of the measurements y i , the functional models , and the random errors can be denoted by , () f  and , respectively. Here  is the t-dimensional vector of unknown parameters and stands for the Hadamard product of matrices and/or vectors. In the observational equations (6), 1 stands for a vector with all its elements being equal to unity. Without loss of generality and for convenience of derivations, we will further assume that the elements of the random vector are all independent, namely, the variance-covariance matrix of is assumed to be equal to , where is the identity matrix. If the functional models are linear with respect to ,  namely, , is a t-dimensional vector, in this case, the model (6) with multiplicative errors can be rewritten as follows: where . It is obvious from the model (7) that the errors of the measurements are proportional to the strengths of the signals. The stronger the signals  , the larger the errors of .
Mathematically, the variance-covariance matrix of the measurements is the function of the unknown parameters . In the following, we will briefly outline three LS methods for parameter estimation, namely, the ordinary LS method, the weighted LS method and the bias-corrected weighted LS method.

The Three LS-Based Methods for Parameter Estimation
When the ordinary LS method is applied to the linear model (7) with multiplicative errors, we have the following minimization problem: It is well known that the optimal solution to (8) is the LS estimate of the parameters , i.e.
It is trivial to prove that the ordinary LS estimate LŜ  is an unbiased estimate of  .
If the random errors in the model (7) are assumed to be statistically independent with the same variance, the measurements are still statistically independent but are no longer of the same variance or weight. In this case, the weight matrix of the measurements now depends on the unknown parameters . By denoting , with standing for a diagonal matrix, we have the variance-covariance matrix of the measurements in (7) as follows: Since the elements of D 1 are the functions of the parameters  , those of the variance-covariance matrix of the measurements y are also the functions of  . If one replaces  in (10) with its approximate value, one can obtain an approximate variance-covariance matrix and then go on to estimate  with it. If an approximate value is sufficiently precise, the approximate estimate of  should be almost the same as the optimal estimate with a true variance-covariance matrix of the measurements. However, if approximate knowledge is not sufficient to obtain a precise approximation of y Σ , the difference between the approximate and optimal estimates cannot be neglected, as can be seen in Xu et al. [46].
As a result of (10), we can apply the weighted LS method to (7) and obtain the following minimization problem: Differentiating (11) with respect to the unknown parameters  and letting it equal zero, we have: where WLŜ  stands for the weighted LS estimate of  , and is the estimate of . According to the rule of matrix differentiation (see e.g., [47]), we have: Inserting (13) into (12) and after some arrangement, we can obtain the weighted LS estimate of , which is denoted by WLŜ  and given as follows: The weighted LS estimate WLŜ  of (14) is obviously nonlinear with respect to the measurements , even though the functional models are linear in the unknown parameters  , since the matrix and the elements of the second term on the right hand side of (14) are all the nonlinear functions of the weighted LS estimate itself. In general, we can only use numerical methods to solve for . Statistically, Xu and Shimada [33] proved that the weighted LS estimate is biased, even though are linear; this is fundamentally different from the weighted LS estimation in a linear model with additive random errors. They also derived the bias of , and further proved that the bias of the weighted LS estimate is solely due to the fact that the variance-covariance matrix of the measurements is the function of the unknown parameters  . In other words, the nonzero partial derivatives in (13) or (14) should be totally responsible for the bias of the weighted LS estimate [33,46]. In the case of linear models with additive random errors, since the variance-covariance matrix of measurements is independent of the unknown parameters, all are equal to zero. Therefore, they suggested removing the second term on the right hand of (14) and constructed an unbiased estimate as follows: where bĉ  is the bias-corrected weighted LS estimate.
We should note that formula (15) is the same as that derived by using the quasi-likelihood method, which nevertheless requires the class of exponential distributions in order for quasi-likelihood to be equivalent to maximum likelihood [27]. However, in this paper, we do not assume any distribution for the measurements y. Actually, the estimate (15) is obtained solely from applying the weighted LS principle, indicating that more complicated adjustment with multiplicative errors can also be solved by using the conventional LS estimation principle.

Error Analysis of Quantities in the Model with Multiplicative Errors
Since Xu and Shimada [33] only focused on the estimation of parameters and numerical simulations to compute the biases and accuracy of the estimated parameters, we will complement their study by providing an analytical error analysis of the adjusted quantities. We should point out that error analysis is well documented in linear models with additive random errors, but nothing along the same line has ever been available in association with measurements contaminated by multiplicative errors, at least, to our best knowledge. Alternatively, one may carry out an error analysis with an approximate variance-covariance matrix of the measurements, as can be seen in geodetic adjustment of distance networks; but such an analysis cannot serve as a solid theoretical foundation of statistical inference in models with multiplicative errors, the extent of approximation will depend on the model itself, the stochastic model of measurements and prior knowledge on the unknown parameters. In this section, we will only assume the first and second central moments for the measurements y and will derive the accuracy formulae of the estimated quantities for the three LS-based estimation methods. We will limit ourselves to the first order approximation. We should note that some of the theoretical formulae to be given below may contain the unknown parameters, which should be substituted by their approximate values or estimates in practical computation.

Accuracy of the Estimated Parameters
Applying the error propagation law to the ordinary LS estimate (9) and after taking the variance-covariance matrix (10) into account, we can readily obtain the variance-covariance matrix of LŜ  as follows: Unlike the ordinary LS estimate LŜ  , (14) and (15) have clearly shown that both the weighted LS estimate and the bias-corrected weighted LS estimate bĉ  are nonlinear with respect to the measurements . In order to derive their accuracy formulae, we will have to represent and bĉ  in terms of the random errors ε and then follow the definition of variance to derive the variancecovariance matrices of WLŜ  and bĉ  . In this paper, we will limit ourselves to the first order approximation for accuracy assessment. Since the weighted LS estimate is not an unbiased estimate, we will also have to take its bias into account. In other words, we will have to compute the mean squared error matrix of WLŜ  for proper error analysis. In order to derive the bias of the weighted LS estimate WLŜ  , Xu and Shimada [33] Taylorexpanded with respect to the random errors ε and then truncated the expansion up to the second order term, which is directly written as follows: where: and is a second order term of ε, and its mathematical expectation is the bias of the weighted LS estimate.
We now apply the error propagation law to (17a) up to the first order approximation of ε and as a result, obtain the first order approximation of the variance-covariance matrix of the weighted LS estimate as follows: Denoting the bias of WLŜ  by E( ) b  , we can finally obtain the mean squared error matrix as follows: where . Since the essential difference between the weighted LS estimate WLŜ  and the bias-corrected weighted LS estimate bĉ  is that bĉ  is unbiased up to the second order approximation. In other words, their first order representations in terms of are identical. Therefore, we can simply write the variance-covariance matrix of bĉ  below:

Accuracy of the Adjusted Measurements
By using the ordinary LS estimate LŜ  in (9), we can readily compute the corresponding LS adjusted values of the measurements: According to the error propagation law, we can then obtain the variance-covariance matrix of the adjusted measurements as follows: In the similar manner, based on the weighted LS estimate and the bias-corrected weighted LS estimate and by using the formulae (19) and (20), we can directly obtain the adjusted measurements computed from the weighted LS estimate and the bias-corrected weighted LS estimate, and their mean squared error (MSE) and variance-covariance matrices, respectively, as follows:

Accuracy of the Corrections of Measurements
As is well known, the corrections of measurements are important quantities in practical data processing and quality assessment/control. We denote the corrections of measurements as follows:

V y y 
where V is the correction vector of the measurements y. By substituting (21a), (22a) and (23a) into the above formula, we can readily obtain the corrections of measurements with the ordinary LS, weighted LS and bias-corrected weighted LS estimates, respectively, as follows: By linearizing (24a), (25a) and (26a) with respect to the random errors ε, and after applying the error propagation law to the linearized corrections of measurements, we can then derive their corresponding variance-covariance matrices:

The Covariances of the Adjusted Quantities
In practical data processing, we are also interested in computing the cross-covariance matrices of the original measurements, the estimated parameters, the adjusted measurements and the corrections of measurements. With the three LS-based methods in hand, if we follow the above lines of thought for the variance-covariance matrices of these quantities, we can then easily derive and obtain their covariances. If we further limit ourselves to the linear approximation in terms of the random errors  , then the weighted LS estimate and the bias-corrected weighted LS estimate should lead to the same covariance representations. Thus, by omitting the details of derivations, we simply list the crosscovariance matrices among the measurements, the estimated parameters, the adjusted measurements and the corrections of measurements for the ordinary LS and the bias-corrected weighted LS estimates in Table 1.

The Estimates of Variance of Unit Weight with the Three LS-Based Methods
The variance of unit weight is one of the most important quantities in statistical quality evaluation and hypothesis testing. Since the measurements are not equally weighted, the conventional estimator of the variance of unit weight by using the ordinary LS residuals of measurements and the redundant number of (n-t) is not unbiased [48]. To estimate the variance of unit weight, we have to start with the naï ve or weighted sum of square of the corrections of measurements, find its mathematical expectation in terms of the variance of unit weight 2  and then estimate 2  . In what follows, we will derive the estimates of the variance of unit weight in association with the three LS-based methods. From the formula (24b), we can readily derive the mathematical expectation of the sum of square of the ordinary LS corrections of measurements: from which we can then construct the following estimate of the variance of unit weight: where LS1 r is the equivalent redundant number of measurements and given by If we use the weighted sum of square of the corrections of measurements in (28a), then the corresponding mathematical expectation can be rewritten as follows:

X( X X ) X I D X( X X ) X
where the second redundant number LS2 r of measurements is given by For the weighted LS and bias-corrected weighted LS methods, if the bias of the weighted LS estimate is not significant, both methods should lead to almost the same estimates of the variance of unit weight. In the similar manner to the derivation of the estimate of the variance of unit weight in association with the ordinary LS method, we can readily use (25b) and (26b) to obtain the estimates of the variance of unit weight by using the weighted LS and the bias-corrected weighted LS estimates, which are simply listed, respectively, as follows:  stand for the estimates of the variance of unit weight by using the weighted LS and bias-corrected weighted LS estimates, respectively. If we remove the biases of the corrections of measurements due to the bias in the weighted LS estimate from the corrections of measurements, we can then obtain a new estimate of the variance of unit weight as follows: where  b is the bias of the weighted LS estimate.

Numerical Examples and Practical Effect on the Estimate of Volume from LiDAR-Type DEM
DEM has been playing an increasing role in hazard assessment (see e.g., [39][40][41]). A precise DEM will reduce the uncertainty of hazard mapping and could help make the operation and management of disaster rescue and support more focused. LiDAR has been widely used to construct DEM (see e.g., [37,39,40,[43][44][45]). Although LiDAR measurements have been proved, both theoretically and experimentally, to be of multiplicative random errors (see e.g., [49][50][51]), they have been treated as if they were of additive error nature in DEM construction. Thus, this section is to serve two major purposes through numerical simulations: (i) to investigate the effect of LiDAR-type measurements on DEM construction; and (ii) to collectively use the error analysis and the estimates of the variance of unit weight given in the previous two sections to estimate the errors of the volume of landslide mass, since a precise estimate of the volume of a landslide from the constructed DEM can be important in practical hazard evaluation.
In this section, we will follow the website http://geology.com/usgs/landslides/ with information on landslides and choose a rotational landslide model to simulate our numerical examples. A landslide of rotational type may create a curved surface. If the reader is interested in the worldwide damaging landslides in the 20th and 21st centuries, she/he may refer to the website http://landslides.usgs.gov/learning/majorls.php of the United States Geological Survey (USGS) for more information, with a landslide ranging from a few to many thousands millions cubic meter. The largest landslide in record seems to be the 1920 Haiyuan landslide, with an area of 50,000 square km and causing 100,000 plus deaths, due to the Haiyuan earthquake of magnitude M8.5, according to the above USGS website.
To start with the simulation, we first assume a simple mountain of trapezoidal prism, with a slope of 50°, a height of 700 m and the length of 500 m. We then assume that a landslide occurs along the slope of the mountain, say due to a large earthquake. The design of landslide examples may be illustrated in Figure 1, with the left panel showing the original trapezoidal mountain and the right panel corresponding to the same mountain after the landslide. The volume of landslide mass can be computed by the following integration: x, y S2 x, y dxdy (33) where D stand for the integration domain of the landslide surface S2(x, y) projected on the xy plane, the plane function S1(x, y) and the curved surface S2(x, y) are respectively represented as follows: In this simulated example, the volume of the landslide mass is about 24.5686 million m 3 . We now assume that DEMs are constructed by using airborne LiDAR before and after the landslide occurs, respectively. The height of flight is fixed at 900 m and the trajectory of the flight is assumed to be free of errors. For simplicity of simulation, we may also assume that the slope plane S1(x, y) before the landslide is sufficiently precise and treated as error-free. Since the surface function S2(x, y) has only six unknown parameters to be estimated for the complete reconstruction of the DEMs after the landslide, we use the uniform distribution to generate 30 LiDAR points for S2(x, y); this number of points is sufficiently large to estimate the variance of unit weight reliably. The standard deviation of multiplicative errors is equal to 0.5 per cent, i.e.,  = 0.005, which is transformed into an equivalent standard deviation of unit weight of 1.268 m for the LiDAR ranging measurements.
With the simulated LiDAR data at hand, we can then estimate the surface S2(x, y) by using each of the three estimation methods given in Section 2 and further compute the estimated volume of the mass due to the landslide, which can be symbolically written as follows: where j sm V is the estimate of the volume sm V of the landslide, the superscript j stands for each of the three estimation methods, namely, the ordinary LS, bias-corrected weighted LS and weighted LS methods. V 0 and the coefficients kl b are constant, independent of the measurements and the estimated parameters but calculated by integrating the corresponding given functions of S1(x, y) and S2(x, y) without the unknown parameters  2 in the domain D , and: ( , , ,  Figure 1. Illustrated model of landslide of rotational type. The left and right panels stand for the model mountains before and after the landslide, respectively. Table 2 are the true and estimated volumes of the simulated landslide mass by the three methods. It is clear from Table 2 that the bias-corrected weighted LS method leads to a difference of 38,664 m 3 from the true value and performs better than the ordinary and weighted LS methods. The difference between the true volume and that estimated from the ordinary LS method is 80,790 m 3 , about the double of the difference by using the bias-corrected weighted LS method; a landslide with a size of this difference by itself is sufficiently large to cause severe damaging, when compared to the large landslides of the 20th and 21st centuries listed in the USGS website http://landslides.usgs.gov/learning/majorls.php. Nevertheless, the example has indicated that all the three methods can lead to good results, if the measurements are sufficiently precise. The abbreviations OLS, BCLS and WLS stand for the ordinary LS, bias-corrected and weighted LS methods, Volume for the volume of the landslide mass, and RMSE for the root mean squared errors of the estimated volumes of the landslide mass. σ is the standard deviation of unit weight.

Listed in
In order to evaluate the errors of the estimated volume by each of the three methods, we apply the results of error analysis in Section 3 and the estimates of variance of unit weight in Section 4 to the estimate of volume (35) and readily obtain the error estimates for the estimated volumes of landslide mass. By the ordinary LS method one usually means to use ( ) More precisely, to estimate the accuracy of the estimated volume of the landslide, we first estimate the variance of unit weight by using (28b-c) for the ordinary LS method, (30) and (31) for the weighted LS and the bias-corrected weighted LS methods, respectively, and then compute the variance-covariance matrix of the estimated parameters in the case of the ordinary LS and bias-corrected weighted LS methods and the MSE matrix in the case of the weighted LS method. With these results at hand, we can finally apply the error propagation law to (35) and obtain the accuracy of the estimated volume of the landslide by each of the three methods. The MSE roots (RMSE) of the estimated volumes are given in row RMSE of Table 2. Obviously, the RMSE value by the bias-corrected weighted LS method is the smallest, which is nevertheless only better than the ordinary LS method by about 3.5 per cent in this specific example. The results of this example indicate that if LiDAR measurements are sufficiently precise, the three methods are able to precisely estimate the volume of landslide mass with sufficiently good accuracy. One may wonder whether the improvement of the bias-corrected weighted LS method over the ordinary LS method is marginally small. To answer this question, we simulate 500 independent examples with 30 different sampling points and repeat the above procedure to compute the RMSE values. The accuracy improvement of volume estimation with the bias-corrected weighted LS method over the ordinary LS method is shown in Figure 2. Depending on the locations of the measurements, the maximum improvement reaches about 19 per cent, which can have a significant consequence in hazard evaluation. Since the LiDAR measurements are assumed to be very precise, the biases of the estimated parameters from the weighted LS method are very small, as can be seen from the RMSE value in Table 2 in this example.

Concluding Remarks
Adjustment has been founded on the basis of observational models with additive random errors to process data in geoscience. The most important feature of an observational model with additive errors is that random errors do not change with the size of a signal or a measurement. However, with the advance of modern space observation technology, we realize that some of random errors change with the sizes of measurements, which may be attributed to the random nature of physical media along the path of observation. For example, GPS, VLBI and EDM baselines, InSAR and LiDAR measurements all show the feature of multiplicative random errors in the sense that their accuracy always contains one term that is proportional to the length of the baseline or the strength of signal. Theoretically speaking, the theory and methods developed on the basis of observational models with additive errors cannot meet the need of models with multiplicative errors. New theory and methods have to be developed accordingly. Xu and Shimada [33] demonstrated that the LS principle can be used to adjust observational models with multiplicative errors, unless the induced bias in the weighted LS estimate can be removed. This paper has substantially extended the work by Xu and Shimada [33] by supplementing a complete error analysis for the quantities of interest. We have derived the cross-covariance matrices among such quantities. Since the variance of unit weight has been one of the most important quantities in statistical analysis of measurements and hypothesis testing, we have also constructed five estimators that correspond to the ordinary LS, the weighted LS and the bias-corrected weighted LS methods. Although LiDAR measurements have been known, theoretically and experimentally, to be of multiplicative error nature (see e.g., [37,39,[43][44][45]), they have been used to reconstruct DEM as if they were of additive random errors. We have simulated an example of reconstructing DEM from LiDAR-type measurements and investigated the effect on volume estimation from the reconstructed DEM. The results have shown that all the three methods can well reconstruct the DEM, if measurements are sufficiently precise. The bias-corrected LS method can perform much better than the ordinary LS method, with a maximum improvement of about 19 per cent from 500 independent repetitions of the example experiment. Such an improvement could be important in hazard evaluation of landslides, for example.