# An Improved Moving Least Squares Method for Curve and Surface Fitting

Authors: Lei Zhang; Tianqi Gu; Ji Zhao; Shijun Ji; Ming Hu; Xiangbo Li
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

The moving least squares (MLS) method has been developed for the fitting of measured data contaminated with random error. The local approximants of MLS method only take the error of dependent variable into account, whereas the independent variable of measured data always contains random error. Considering the errors of all variables, this paper presents an improved moving least squares (IMLS) method to generate curve and surface for the measured data. In IMLS method, total least squares (TLS) with a parameterλ based on singular value decomposition is introduced to the local approximants. A procedure is developed to determine the parameter λ. Numerical examples for curve and surface fitting are given to prove the performance of IMLS method.


---

## Body

## 1. Introduction

Methods of approximation or interpolation of measured data are considered by many authors [1–6]. The moving least squares (MLS) as an approximation method has been introduced by Shepard [7] in the lowest order case and generalized to higher degree by Lancaster and Salkauskas [8]. The object of those works is to provide an alternative to classic interpolation for approximating a function from its values given at irregularly spaced points by using weighted least squares approximations. This method is to start with a weighted least squares formulation for an arbitrary fixed point and then move this point over the entire parameter domain, where a weighted least squares fit is calculated and evaluated for each measured point individually. The MLS approximation has been largely documented in the literature and used by many scholars for optimization problems and widely used in many types of application in the engineering literatures [9, 10]. Among them, Breitkopf et al. [11] and Naceur et al. [12] have presented an extended approach of pattern search algorithms with a fixed pattern panned and zoomed in a continuous manner across the design space.But MLS method determines local approximants in the sense of ordinary least squares (OLS). When errors occur in all of the variables, it makes more sense to determine the local approximants in the sense of total least squares (TLS) [13]. TLS approach is tightly related to the maximum likelihood principal component analysis method introduced in chemometrics by Schuermans et al. [14]. TLS, as same as OLS, could not well apply to the complicated curve and surface fitting. But compared with OLS method, TLS method considers the errors of all variables.In this paper, an improved moving least squares (IMLS) method is proposed, in which TLS with a parameterλ based on singular value decomposition for considering the varying errors of all variables is introduced to the local approximants of MLS method for determining the local coefficients. In Section 2, a brief description is given for MLS method. IMLS method is presented in detail in Section 3. And curve and surface fitting examples are given in Section 4 for comparing IMLS method with MLS method. Conclusions are drawn in Section 5.


## 2. MLS Method

In the following, only the main aspects of MLS method will be emphasized. In the MLS approximation, the trial function can be expressed as(1)u(x)=∑i=1mpi(x)ai(x)=pT(x)a(x),
where pT(x)=[p1(x),p2(x),…,pm(x)] is a primary function of mth-degree and  a(x)=[a1(x),a2(x),…,am(x)]T represents the corresponding undetermined matrix of coefficients, in which these coefficients are a function of the coordinates of point x.For the local approximants atx, a scattered normal form is defined as
(2)J(x)=∑I=1nw(x-xI)[∑i=1mpi(xI)ai(x)-uI]2,
where xI is the discrete point to be fitted within the influence domain of x. And w(x-xI)  is a weight function which defines the influence domain of  x and attributes a weight to each discrete point depending on its position relative to x. The weight function is positive and its value increases with the decrease of the distance ∥x-xI∥ between x and xI. Various ways of choosing such functions can be found in the literature [15–17].


## 3. IMLS Method

### 3.1. TLS Method

The TLS problem [18, 19]
(3)minΔA,ΔB,X∥[ΔAΔB]∥F2s.t.(A+ΔA)X=B+ΔB
is a generalization of the OLS method when an errors-in-variables (EIV) model
(4)A=A0+A~,B=B0+B~,A0X0=B0
is considered. Here A~, B~ are measurement errors and A0,  B0 are true values. The TLS estimate of X0, that is, the solution of (3), corresponding to X, is proven to provide a consistent estimator of  X0, when the elements of A~ and B~ are zero mean i.i.d.Let(5)C:=[AB]=UΣVT,whereΣ=diag(σ1,σ2,…,σn+d),
be a singular value decomposition of C, let  σ1≥⋯≥σn+d, be the singular values of C, and define the partitionings
(6)V:=[V11V12V21V22],Σ:=[Σ100Σ2].A TLS solution exists if and only ifV22 is nonsingular. In addition, it is unique if and only if σn≠σn+1. In this case, TLS solution is given by
(7)X∧tls=-V12V22-1.


### 3.2. IMLS Method

In this paper, a close attention is paid to the local approximants of MLS method. Considering the errors of all variables, TLS method with a parameterλ set between A and B is introduced to the local approximants of MLS method, forming an improved moving least squares (IMLS) method.In IMLS method, defineCx for local approximants at x as
(8)Cx:=Wx[DAxTBx]=UxΣxVxT,
where D=diag(λ,…,λ)n×n and T=diag(1-λ,…,1-λ)n×n are matrices related to λ and Wx=diag(w(x-x1),…,w(x-xn)) is a weight matrix with compact support. In this paper, the following weight function given by [20] is used:
(9)w(x-xI)={e-(∥x-xI∥/αd)2-e-1/α21-e-1/α2if∥x-xI∥≤d,0else,
where the parameter α allows controlling the weight function (here α is fixed to 1/2) and d defines the influence domain of x.In the influence domain atx, the coefficients of local approximants are solved by
(10)a(x)=[a1(x),a2(x),…,am(x)]=-Vx12Vx22-1.In this paper, the linear approximation of local approximants is considered only.


### 3.3. The Determination of Parameterλ

As known, it is a difficult task to analyze the fitting error for complicated curve. Unlike the straight line fitting, the data is obtained by measurement device from a real object or mechanical part, and the profile of object has not an explicit analytical expression [21].The following procedure is carried out in this paper to determine the parameterλ in IMLS method.Step 1.
Add the known random error(σx,σy,σz) to the measured data (xim,yim,zim) for getting tested data (xit,yit,zit).Step 2.
λ is set to λ0.Step 3.
Fit the tested data(xit,yit,zit) using IMLS for getting fitting value (xif,yift,zif).Step 4.
According to [22], calculate the fitting error s of the measured value zim and fitting value zif(11)s=∑i=1n|zim-zif|and record the value ofs.Step 5.
Varyλ with an increment of certain step h.Step 6.
Repeat Steps3–5 until s reaches the minimum value.Step 7.
Recordλ corresponding to the minimum value of s.Step 8.
Repeat Steps1–7  r  times.Step 9.
Average the recorded values ofλ and take it as the final value for  λ  in IMLS method.


## 3.1. TLS Method

The TLS problem [18, 19]
(3)minΔA,ΔB,X∥[ΔAΔB]∥F2s.t.(A+ΔA)X=B+ΔB
is a generalization of the OLS method when an errors-in-variables (EIV) model
(4)A=A0+A~,B=B0+B~,A0X0=B0
is considered. Here A~, B~ are measurement errors and A0,  B0 are true values. The TLS estimate of X0, that is, the solution of (3), corresponding to X, is proven to provide a consistent estimator of  X0, when the elements of A~ and B~ are zero mean i.i.d.Let(5)C:=[AB]=UΣVT,whereΣ=diag(σ1,σ2,…,σn+d),
be a singular value decomposition of C, let  σ1≥⋯≥σn+d, be the singular values of C, and define the partitionings
(6)V:=[V11V12V21V22],Σ:=[Σ100Σ2].A TLS solution exists if and only ifV22 is nonsingular. In addition, it is unique if and only if σn≠σn+1. In this case, TLS solution is given by
(7)X∧tls=-V12V22-1.


## 3.2. IMLS Method

In this paper, a close attention is paid to the local approximants of MLS method. Considering the errors of all variables, TLS method with a parameterλ set between A and B is introduced to the local approximants of MLS method, forming an improved moving least squares (IMLS) method.In IMLS method, defineCx for local approximants at x as
(8)Cx:=Wx[DAxTBx]=UxΣxVxT,
where D=diag(λ,…,λ)n×n and T=diag(1-λ,…,1-λ)n×n are matrices related to λ and Wx=diag(w(x-x1),…,w(x-xn)) is a weight matrix with compact support. In this paper, the following weight function given by [20] is used:
(9)w(x-xI)={e-(∥x-xI∥/αd)2-e-1/α21-e-1/α2if∥x-xI∥≤d,0else,
where the parameter α allows controlling the weight function (here α is fixed to 1/2) and d defines the influence domain of x.In the influence domain atx, the coefficients of local approximants are solved by
(10)a(x)=[a1(x),a2(x),…,am(x)]=-Vx12Vx22-1.In this paper, the linear approximation of local approximants is considered only.


## 3.3. The Determination of Parameterλ

As known, it is a difficult task to analyze the fitting error for complicated curve. Unlike the straight line fitting, the data is obtained by measurement device from a real object or mechanical part, and the profile of object has not an explicit analytical expression [21].The following procedure is carried out in this paper to determine the parameterλ in IMLS method.Step 1.
Add the known random error(σx,σy,σz) to the measured data (xim,yim,zim) for getting tested data (xit,yit,zit).Step 2.
λ is set to λ0.Step 3.
Fit the tested data(xit,yit,zit) using IMLS for getting fitting value (xif,yift,zif).Step 4.
According to [22], calculate the fitting error s of the measured value zim and fitting value zif(11)s=∑i=1n|zim-zif|and record the value ofs.Step 5.
Varyλ with an increment of certain step h.Step 6.
Repeat Steps3–5 until s reaches the minimum value.Step 7.
Recordλ corresponding to the minimum value of s.Step 8.
Repeat Steps1–7  r  times.Step 9.
Average the recorded values ofλ and take it as the final value for  λ  in IMLS method.


## 4. Numerical Examples and Analysis

In this section, three examples are given to verify the performance of IMLS method. MLS method is also applied in the examples to make a comparison.Example 1.
Consider the function(12)y=e-x/3sin(3x)
and choose a uniformly distributed set of points (xi,yi), i=1,2,…,n, determined by (12)   and add σx and σy, which are normally distributed random errors with mean value of zero, to xi and yi, respectively, forming a set of measured points. And λ in the IMLS method is determined by the above-mentioned procedure. MLS and IMLS are both applied to fit the measured points. The fitting performance is characterized by the fitting error of the real value and fitting value
(13)s0=∑i=1n|yi-yif|,
where yi is the real value and yif here is the fitting value for measured data.The fitting results are summarized in Table1. The fitting curve for Example 1 using IMLS method is illustrated in Figure 1.Table 1
Thes0 of two methods for Example 1.
σ
x
σ
y
s
0
MLS
IMLS
λ
0.000001
0.001
1.320980
0.175673
0.49310417
0.00001
0.001
1.323738
0.172833
0.49310563
0.0001
0.001
1.316369
0.172238
0.49312047
0.001
0.001
1.320659
0.177658
0.49313993
0.001
0.0001
1.317848
0.157265
0.49312350
0.001
0.00001
1.318954
0.146745
0.49311082
0.001
0.000001
1.315052
0.145775
0.49309692Figure 1
The fitting curve of IMLS method for Example1.Example 2.
Consider the function(14)y=cx21+1-(1+k)c2x2,
where c=1/1083  is the reciprocal of the curvature radius of the base vertex and k=-1 is the constant of the quadric surface. The measured points are generated in the same way as Example 1 and fitted by MLS and IMLS, respectively. And fitting results are listed in Table 2. The fitting curve for Example 2 using IMLS method is illustrated in Figure 2.Table 2
Thes0 of two methods for Example 2.
σ
x
σ
y
s
0
MLS
IMLS
λ
0.000001
0.001
0.050892
0.041262
0.50019680
0.00001
0.001
0.049675
0.042708
0.50018100
0.0001
0.001
0.049482
0.040264
0.50015440
0.001
0.001
0.049789
0.042724
0.50010860
0.001
0.0001
0.041764
0.013429
0.50013496
0.001
0.00001
0.041718
0.013170
0.50013892
0.001
0.000001
0.041467
0.013092
0.50014048Figure 2
The fitting curve of IMLS method for Example2.Example 3.
Consider the function(15)z=(x2-y2)2
defined on the region Ω=[-1,1]×[-1,1]. Take a random uniformly distributed set of points in Ω. The measured points are generated in the same way as Examples 1 and 2 and fitted by MLS and IMLS methods, respectively. The fitting performance is characterized by the fitting error of the real value and fitting value
(16)s0=∑i=1n|zi-zif|,
where zi is the real value and zif here is the fitting value for measured data. And the fitting results are listed in Table 3. The fitting surface for Example 3 using IMLS method is illustrated in Figure 3.Table 3
Thes0 of two methods for Example 3.
σ
x
,
y
σ
z
s
0
MLS
IMLS
λ
0.000001
0.001
2.525123
0.592784
0.500148
0.00001
0.001
2.531644
0.594762
0.500158
0.0001
0.001
2.525867
0.593739
0.500166
0.001
0.001
2.578295
0.754762
0.500534
0.001
0.0001
2.555807
0.664834
0.500680
0.001
0.00001
2.553651
0.667688
0.500666
0.001
0.000001
2.558244
0.666255
0.500656Figure 3
The fitting surface of IMLS method for Example3.From the above examples, it can be seen that the fitting results of IMLS are always better than those of MLS method. The variations of parameterλ with the errors of variables for three examples based on the above procedure are shown in Figure 4. The results of three examples show that parameter λ may be the same value when the errors of independent variable and dependent variable are totally different. Moreover, parameter λ illustrates a nonstrict symmetry and the valley (peak) value is not necessarily determined at the point lg(σx/σy)=0  (or lg(σx,y/σz)=0).The variation of parameterλ with the errors of variables.
(a)
Example1
(b)
Example2
(c)
Example3The computation for case study is implemented on a personal computer with MATLAB. The extra time is required for determining the parameterλ in IMLS. In Example 1, when m=2,  n=160,  d=0.3/10,  h=0.00001,  r=100, and λ varies from 0.49 to 0.51, the time required for determining λ is about 1.5 hours.


## 5. Conclusions

The advantage of the MLS approximation is to obtain the shape function with higher order continuity and consistency by employing the basis functions with lower order and choosing a suitable weight function with compact support. The local approximants of MLS only consider the error of dependent variable. However, all variables in the measurement are always perturbed. The paper presents a new curve and surface fitting method called improved moving least squares (IMLS) method which considers the errors of all variables. Different from MLS method, TLS with a parameterλ based on singular value decomposition is introduced to the local approximants in IMLS method. A procedure is given to determine the parameter λ. The curve and surface fitting results by IMLS method for the discrete points generated by numerical simulation are compared with those by MLS method on the same condition. IMLS method gives better performance than MLS method for the varying errors of variables, which confirms the validity of the proposed IMLS method in this paper.



---

*Source: 159694-2013-11-20.xml*