Multi-quadric collocation model of horizontal crustal movement

To establish the horizontal crustal movement velocity field of the Chinese mainland, a Hardy multi-quadric fitting model and collocation are usually used. However, the kernel function, nodes, and smoothing factor are difficult to determine in the Hardy function interpolation. Furthermore, the covariance function of the stochastic signal must be carefully constructed in the collocation model, which is not trivial. In this paper, a new combined estimation method for establishing the velocity field, based on collocation and multi-quadric equation interpolation, is presented. The crustal movement estimation simultaneously takes into consideration an Euler vector as the crustal movement trend and the local distortions as the stochastic signals, and a kernel function of the multi-quadric fitting model substitutes for the covariance function of collocation. The velocities of a set of 1070 reference stations were obtained from the Crustal Movement Observation Network of China, and the corresponding velocity field was established using the new combined estimation method. A total of 85 reference stations were used as checkpoints, and the precision in the north and east component was 1.25 and 0.80 mm yr−1, respectively. The result obtained by the new method corresponds with the collocation method and multi-quadric interpolation without requiring the covariance equation for the signals.


Introduction
Horizontal movement velocity fields provide the main basic data for Earth science research.Because measuring a station's velocity by repeated observation is always limited, with many points that cannot be directly measured, a mathe-matical velocity field model is always used to obtain the velocity field.How to obtain reliable horizontal velocity fields from measured points has been the focus of a lot of research all over the world (Argus and Gordon, 1991;Huang et al., 1993;Liu et al., 2001Liu et al., , 2002;;Chai et al., 2009;Jiang and Liu, 2010;Hu and Wang, 2012;Zeng et al., 2012Zeng et al., , 2013)).As we know, the horizontal velocity of a ground point mainly consists of two aspects: the overall movement and the local deformation reflected ground motion.The motion model for horizontal velocity is often used in geophysical models and the statistical fitting method.The geophysical model is often implemented using the Euler vector method (Argus and Gordon, 1991).Generally, the condition for using the Euler vector method is to divide the block reliably and treat each division as a rigid body, but many areas do not satisfy the requirements of a rigid body, which limits the usefulness of the method.The statistical fitting method mainly includes multiquadric functions (Huang et al., 1993;Liu et al., 2001;Zeng et al., 2013) and collocation (Liu et al., 2002;Chai et al., 2009;Jiang and Liu, 2010;Zeng et al., 2012).The multiquadric functions (Hardy, 1978) can be used for fitting the parameters of measured points and estimating the parameters of unmeasured points.The mathematical methods are used in this case, while their physical meaning is not clearly considered.The key issues and difficult problems in their application are the choice of kernel function, smoothing factor, and node.Some scholars have systematically researched their application in horizontal velocity field models in the Chinese mainland (Huang et al., 1993;Nie et al., 2007;Zeng et al., 2013).In order to consider the changing information of velocity fields in different local regions, the least squares collocation method, whose key problem is to determine the co-variance function of signal vectors, can be adopted (Yang, 1992).To balance the contribution of the estimation results from the signal covariance matrix and observation noise, we can use the variance components to estimate the collocation solutions (Yang and Liu, 2002;Yang et al., 2008) or adaptive collocation solutions (Yang and Zeng, 2009;Yang et al., 2009Yang et al., , 2011)), but both of them require iterative calculations.Because the collocation covariance function is very difficult to build, some scholars have analyzed and compared the relationship between the multi-quadric function and covariance collocation (Tao et al., 2002); some scholars, using the compensation principle of least squares, have pointed out that choosing the appropriate regularization parameters after the semi-parametric model may include the collocation model (Tao and Yao, 2003); and some scholars have established a parameter estimation model combining the multiquadric function and configuration models and used the distance function as the signal covariance when estimating the gravity field to obtain results close to the multi-quadric function model (Wang and Ou, 2004).
In view of the above-mentioned facts, in this paper, we take Euler vectors as the long-term overall movement trend of the function model and regard the local variations of horizontal movements as stochastic parameters.At the same time, the stochastic characters are taken into account.We utilize a multi-quadric kernel function to obtain its normalized matrix and combine the characters of the collocation and multi-quadric functions, which not only avoids having to establish the covariance function but also achieves an integrated effect using both collocation and multi-quadric functions (Tao and Yao, 2003).On the basis of the method proposed above, we established a Chinese mainland horizontal movement velocity field by using the velocities of a set of 1041 reference stations, obtained from the Crustal Movement Observation Network of China.In this paper, a new combined estimation method using the kernel function of a multi-quadric fitting model to replace the covariance function of collocation was proposed and tested.

Collocation model of horizontal movement
It is well known that the horizontal velocity of a ground point can be described in two parts.One is the holistic long-term trend of horizontal movement expressed by the Euler vector, which is related to the block where the sites are.The other part is the local movement change, modeled as stochastic signals, which is mainly affected by local surface displacement (Freymueller et al., 2013).Therefore, the function model of point horizontal velocity is where A represents the designed matrix of size 2n × 3, B denotes the coefficient matrix of size n × u, where n and u are the station number and estimated parameter number, re-spectively; ω represents the estimation of the Euler vector because of plate movement, and L is the observation vector of horizontal velocity on a station; the covariance matrix is e = σ 2 0 P −1 e , where P e denotes its weight matrix.Ŝ is the signal estimation of observed points, and its covariance matrix is S = σ 2 0 P −1 S , where P S denotes its weight matrix.The relation between the long-term movement trend of horizontal velocity on a ground station and the Euler vector of plate movement can be expressed as follows (Jin et al., 2006;Yang and Zeng, 2009): where ν e and ν n denote the east and north components of the horizontal velocity of ground stations, respectively, in a topocentric coordinate system.R represents the radius of earth, and λ , φ are the latitude and longitude of the ground point.Thus, the matrix A in Eq. ( 1) is Generally, the measured point signal is the local motion variation without the plate rigid motion, which is the horizontal residual velocity in measured points; so B = I (unit matrix) in Eq. ( 1).
If we take the local movement variation signal Ŝ of a nonobserved point into account and assume that there is a covariance matrix SS = T SS = 0 between the non-observed point's signal Ŝ and the observed point's signal Ŝ then, if the covariance matrix is known, the collocation solution considering the non-observed point's signal simultaneously (Yang and Zeng, 2009;Zeng et al., 2012) is where P L = B S B T + e −1 .The Euler vector estimation ω of plate movement and velocity estimation Ŝ of local movement in the collocation model above depend not only on the covariance matrix e of observed horizontal velocity L but also on the covariance matrix S of the local velocity signal Ŝ.Before data processing, the covariance matrices S and SS must be known, which is the key point of the problem.This is the most difficult aspect of applying the collocation method.In order to overcome the difficulty, the elements in the covariance matrix are calculated according to the covariance function.There is a wide variety of stochastic signal covariance functions, such as the Gauss exponent function and the Hirvonen function.
The principles for determining the covariance function are also diverse.Different covariance functions can be based on different data and principles (Zeng et al., 2012).The covariance functions are usually dependent on the following three aspects: (1) the observed data for fitting, (2) the fitting principles, and (3) the parameter setting in the fitting process.As a result, different researchers can obtain different signal covariance functions even though the used data and principles are the same; that is to say, different physical applications take a lot of effort in determining the covariance function, which greatly limits the collocation model applications.

Multi-quadric function of horizontal movement
In practice, the most difficult aspect of adopting the collocation method is that the signal's a priori variance cannot be determined accurately.The signal covariance functions are generally built by adopting observed data based on certain principles (Zeng et al., 2012).Thus, the priori covariance matrix determined in this way has a strong correlation with the current measured data.It not only affects the optimality of collocation, but also reduces its range of application.Allowing for that, the requirements of precision and density of observed data are relatively high when estimating covariance functions.At the same time, the user should be well informed about the physical field that it is used in.Thus, it is difficult to obtain a general application and influence the optimality of the collocation model (Tao and Yao, 2003).We noticed that the signal covariance function of collocation is a kind of function about distance in general (Zeng et al., 2012).Therefore, in this paper, we do not look upon the local deformation parameters as random signals, like the collocation, but as non-random variables, taking into account their random nature (Tao and Yao, 2003), reflected by the distance function, to construct the covariance matrix of the local signals.
The covariance matrix of the local deformation is where g is the value obtained by the distance function.
We still use the Eq. ( 1) as the observation equation.Because of this, the distance function does not explicitly consider the physical nature of the local deformation.Using the compensation least squares method estimation criterion, Then, according to Tao and Yao (2003) and Wang and Ou (2004), where the matrix A and B are the same as the matrix A and B in Eq. ( 1).
Considering the collocation problem of the non-observed point's signal parameter S , we use the distance function to contact the measured point's signal parameter S and the unmeasured point's signal parameter S ; so the SS in Eq. ( 4) is transformed to Then the estimation of the unmeasured point's signal parameter S is where matrix M S in Eq. ( 9) denotes the multi-quadric kernel function matrix, which here is the same as G S S in Eq. ( 8).
We have revised it in the next version.Matrix D in Eq. ( 8) represents the variance matrix of observation error.Thus, the key to the above problem is still to determine the covariance matrix of the local deformation parameters; in this paper we will transform the local deformation parameter covariance matrix to the covariance function to determine it.
Multi-quadric functions were first proposed in 1978 by Hardy (Hardy, 1978).The basic idea is that any smooth surface can be approached at any precision by using a series of finite and regular mathematical functions, and the nonobserved points can be estimated by making use of the observed points.Because it can be designed flexibly and the controllability is strong, the method has been widely used in the interpolation problems involved in geoscience since the approach was proposed.Tao and Yao (2003) discussed the relationship between the multi-quadric function and collocation in detail and thought that the multi-quadric function is a kind of special covariance function.Based on this idea, we introduce a multi-quadric kernel function to determine the covariance matrix of local deformation.Thus, the covariance matrix ( 5) and (8) in Eq. ( 4) are completely determined by the multi-quadric kernel function.Currently, the most commonly used multi-quadric kernel functions in addition to the conical surface are positive double curved surface, inverted double curved surface, positive thrice curved surface, and inverted thrice curved surface.
1. Positive double curved surface: The symbol δ represents the smoothing factor.Obviously, if δ = 0, the positive double curved surface is a conical surface.
2. Inverted double curved surface: Again, δ represents the smoothing factor, and if δ = 0, the inverted double curved surface is a conical surface.
3. Positive thrice curved surface: 4. Inverted thrice curved surface: As can be seen, the actual covariance matrix S of the signal is replaced by the matrix G S determined in the multi-quadric function, which is exactly the same as the standard configuration model in the form; namely, Eqs. ( 4) and ( 6) are exactly identical.Although the two methods are similar in form, their theoretical basis is different.According to the rule of maximum posterior estimation, the standard configuration model obtains the optimal estimate under the condition of the signal a priori variance being known.However, the latter considers the randomness of the parameters, gives a certain compensation and balance, and still considers it as a parameter estimation of non-random parameters, and the optimization is better than the previous method (Wang and Ou, 2004).In this way, we not only overcome the problem that it is difficult to obtain the covariance function of the classical collocation model, but also keep the formula of the collocation model stable.Because of its simplicity and rationality, the application range of the model may be greatly expanded (Tao et al., 2002;Tao and Yao, 2003).

Calculation and analysis
We used the same observation data as in Yang and Zeng (2009), that is, 1041 repeat observation stations of the Crustal Movement Observation Network project, whose precision of horizontal velocity is superior to fitting point coordinates' velocity by 3 mm yr −1 .These data include 85 highprecision points (29 consecutive observation stations and 56 campaign-mode observation stations) as the external inspection points to assess the accuracy of the constructed velocity field model, and a net of the remaining 985 regular observation points after 56 irregular observation basis points are subtracted from the 1041 repeat observation stations to establish the velocity field model.In order to evaluate the internal precision of the model, the calculation formula for the root-mean-square error (RMSE) of the horizontal residual velocity of the 985 observation stations is calculated: where V is residuum calculated by different solutions, that is, the station horizontal residual velocity.
The RMSE of the external checkpoint horizontal velocity can be defined by where Ẋm is the measured horizontal velocity of the 85 checkpoints, and Ẋc is the horizontal velocity of the 85 stations, calculated by different schemes.

Calculation results for different kernel functions
Since the multi-quadric function is used to determine the signal priori random information in the multi-quadric collocation model, different kernel functions have different results.
We use inverted and positive, double curved, twice and thrice curved surface to obtain results.Table 1 shows statistics for velocity residuals for the different kernel functions and Table 2 shows statistics for the 85 points of external inspection accuracy of the different kernel functions.It can be seen that the different kernel functions have different effects on the accuracy of the results.
1. From the 85 external checkpoint precision statistics, the RMSE for east and north of the positive thrice curved surface function is 1.25 and 0.89 mm yr −1 , respectively, but from the velocity residuals statistics, it is 1.29 and 1.14 mm yr −1 , respectively.The accuracy of the external statistics is basically the same as that of the internal statistics, which shows that the fitting and prediction results are stable.The external precision of the inverted double curved surface can also reach 1.38 and 1.18 mm yr −1 , but its residual velocity is only 0.56 and 0.51 mm yr −1 , and the RMSE is obviously better than that of the positive thrice curved surface.
2. From the internal and external statistical accuracy, the positive surface function and the inverse surface function of the statistical results have obvious differences.
The 85-point accuracy of the inverted double curved surface and inverted twice and thrice curved surface in the north direction are, respectively, 1.18, 1.45, and 1.77 mm yr −1 ; in the east direction they are 1.38, 1.67, and 1.93 mm yr −1 .All of them are superior to 2 mm yr −1 .The variation of the different inverse curve functions is small, and this may be due to the fact that inverse curve performs as the reciprocal of distance.The closer the distance, the greater the correlation between points, and the physical characteristics of the velocity field are basically identical.However, the 85-point accuracy of the positive double curved surface and positive twice and thrice curved surface in the north and east directions is, respectively, 1.27, 5.02, and 0.89 mm yr −1 , and 2.61, 5.02, and 1.25 mm yr −1 .The difference between different positive surface functions is obvious which may be due to the positive surface function and the distance of the positive correlation.The greater the distance between the points, the larger the correlation, and this is contrary to the physical characteristics of the velocity field, which therefore has not been well described by the physical properties of the velocity field.

Comparison with other common methods
As we all know, there are three classical and common methods for establishing velocity field, i.e., Euler vector method, the collocation method, and the multi-quadric function method.In order to examine the effective of the proposed multi-quadric collocation model in calculating the horizontal velocity field, the above-mentioned three methods were used for comparison.Scheme 1: with the physical significance of the plate tectonics Euler vector model being exact, we utilize the least squares principle to solve the unified Euler vector and the local horizontal velocity changes of the GPS stations in the Chinese mainland.
Scheme 2: we use the collocation method to estimate the long-term overall movement trend (Euler vector) and the local velocity variance of the GPS stations.In the process of Schemes 1-4 represent the models using Euler vector method, the collocation method, the multi-quadric function method, and the proposed multi-quadric collocation model.Table 3. Euler vectors of Chinese mainland.
Scheme ω x (rad Myr −1 ) ω y (rad Myr −1 ) ω z (rad Scheme 3: the multi-quadric function method, which is the academic mathematics method and is used to establish the velocity field, has been systematically researched in other research (Zeng et al., 2013).Meanwhile, the method does not have the exact physical significant, which could not be utilized to solve the Euler vector, and only is used to estimate the coefficient of the node.Because the key issue for establishing a horizontal velocity field model using the multiquadric function method is how to ensure the kernel function, smooth factor, and node number, we take advantage of the kernel function which defined the hyperbolic function, the smooth factor which ensured 10, and the node which consisted of the 56 based points.Scheme 4: the multi-quadric collocation model, being used to solve the points' the horizontal velocity, has a specific physical significance and is utilized to estimate the longterm overall movement trend (Euler vector) and the velocity changes of local points, with using Eqs.( 7), ( 9), and the inverted twice curved surface, in process of the calculation.We determined the best value of smoothing factor in the inverted twice curved surface function based on the actual distribution of data.
Table 3 lists the Euler vectors of the Chinese mainland calculated by different schemes.Table 4 lists the external precision of the Euler vectors of the Chinese mainland calculated by different schemes.Table 5 lists the velocity residuals statistics for different schemes.Table 6 shows the precision statistics of 85 external checkpoints.Figure 1 shows the velocity residual statistical results of 85 external checkpoints.It reveals that for the multi-quadric collocation model proposed in this study, (1) the numbers of external checkpoints whose velocity residuals were between −0.5 and 0 mm yr −1 in east and north components are 24 and 26, respectively, and (2) the numbers of external checkpoints whose velocity residuals were between 0 and 0.5 mm yr −1 in east and north components are 16 and 17, respectively.These numbers were obviously higher than the results derived from Euler vector method, collocation method, and multi-quadric function method, which demonstrated that the proposed multi-quadric collocation method outperformed the other three methods.
A comprehensive analysis of the results shows the following.
1.As can be seen from Table 3, the Euler vector, calculated by Euler vector estimation method (Scheme 1), collocation (Scheme 2), and multi-quadric collocation models (Scheme 4) are equivalent in magnitude and trends, but as a result of the difference of the function model and the signal covariance function between collocation and the multi-quadric collocation model, they are slightly different in values and in the estimated remaining horizontal velocity.
2. The least squares model based on the Euler vector only takes into account the entire horizontal movement of the Chinese mainland, so its precision is poor.From the perspective of internal precision (Table 4), horizontal residuals for the east and north directions are up to 30.17 and 27.69 mm yr −1 maximum, respectively, and statistical precision for east and north is 4.63 and 6.19 mm yr −1 , respectively.Viewed from the external checkpoint, the external precision (Table 6) of the Euler vector least squares solutions for east and north are 5.89 and 7.06 mm yr −1 , respectively.This is because the Chinese mainland is not a single block; the Euler vector least squares solution can only determine the overall orientation of the Chinese mainland, and does not take the systemic change of the local area into account.
3. Solving using the collocation method not only takes into account the overall orientation, but also determines the local area random effects.It makes a more reasonable distribution of the residuals (Table 5).The internal precision is 1.60 for east and 1.43 mm yr −1 1 for north.The external check precision also improved, to 2.14 for east and mm yr −1 for north.
www.solid-earth.net/7/817/2016/Solid Earth, 7, 817-825, 2016 5. Synthesizing the collocation method and multi-quadric function method not only obviously improved internal precision, with east and north precision reaching 0.78 and 0.73 mm yr −1 , respectively, which is obviously superior to the collocation method and the multi-quadric function method; but external precision has also obviously improved, with east and north direction precision reaching 1.67 and 1.45 mm yr −1 , respectively, which is slightly superior to the collocation method and the multi-quadric function method.This approach integrates the characteristics of the two methods.It avoids the establishment of the covariance function, and achieves a combined effect.
6.The Euler vector, collocation, and multi-quadric collocation models have physical meaning, and can work out the Euler vectors of Chinese mainland blocks.However, the Euler vector method has larger differences from the other two methods, which, in particular, have smaller differences in latitude and longitude components and rotating angle speed aspects.The Euler vectors determined by the collocation method and the multi-quadric function configuration model are very close.
According to the comparison and analysis for the calculation, we have established a 1 • × 1 • gridded velocity model (for details see Fig. 2) and a 1 • × 1 • velocity field model deducting the Chinese continental crustal deformation background field (for details see Fig. 3).
In Fig. 2, the Chinese continental crust shows eastward movement as a whole.Meanwhile, the China continental re-gion possesses significantly clockwise rotation, but velocity values show difference in some parts of China, of the horizontal velocity field, removing the continental movement background in China, in the Fig. 3.There is an eastwardsoutheastward-southwestward movement model with the feature of higher velocity values in the west by 104 longitude degrees in the Sichuan-Yunnan region.In the northwest region, the continental crust displays an NNW motion model.In the Northeast Plain, there was westward movement.Conversely, South China represents an SSE motion model with small velocity values relative to the others in the China.

Conclusions
Multi-quadric functions and the collocation method are the most commonly used methods of fitting the observed points and estimating the non-observed points.However, it is difficult to select the kernel function, smoothing factor, and node when we use the multi-quadric function method.The key to applying the collocation model is to establish a reliable covariance function.Usually, it tends to adopt an empirical formula in a gravitational field, where the majority of applications are established, by using a measured data covariance model, which not only affects the optimality of the collocation method but also reduces its usable area.
Taking into account the establishment of the covariance function is a relatively high requirement for data, and generally, it is difficult to achieve.Bearing in mind that the covariance function is a function of the distance, the local deformations covariance matrix is often chosen as a simple function of distance in multi-quadric functions.The estimating effects of this integrated approach are similar to the collocation method.Because the solution is simple and reasonable, it can greatly expand the scope of application of the collocation model.

Figure 1 .
Figure 1.The velocity residuals statistics of 85 external checkpoints.The horizontal axes represent the horizontal velocity residuals (N: north; E: east; m yr −1 : meters per year) of a checkpoint.The vertical axes indicate the numbers of checkpoints (i.e., Point num) with corresponding residuals.Schemes 1-4 represent the models using Euler vector method, the collocation method, the multi-quadric function method, and the proposed multi-quadric collocation model.

Table 2 .
External inspection accuracy of different kernel functions (mm yr −1 ).