Calibration of MEMS Triaxial Accelerometers Based on the Maximum Likelihood Estimation Method

As a widely used inertial device, a MEMS triaxial accelerometer has zero-bias error, nonorthogonal error, and scale-factor error due to technical defects. Raw readings without calibration might seriously affect the accuracy of inertial navigation system. )erefore, it is necessary to conduct calibration processing before using a MEMS triaxial accelerometer. )is paper presents a MEMS triaxial accelerometer calibration method based on the maximum likelihood estimation method. )e error of the MEMS triaxial accelerometer comes into question, and the optimal estimation function is established. )e calibration parameters are obtained by the Newton iteration method, which is more efficient and accurate. Compared with the least square method, which estimates the parameters of the suboptimal estimation function established under the condition of assuming that the mean of the random noise is zero, the parameters calibrated by the maximum likelihood estimation method are more accurate and stable. Moreover, the proposed method has low computation, which is more functional. Simulation and experimental results using the consumer low-cost MEMS triaxial accelerometer are presented to support the abovementioned superiorities of the maximum likelihood estimation method. )e proposed method has the potential to be applied to other triaxial inertial sensors.


Introduction
Nowadays, with the gradual rise of microelectromechanical system (MEMS), low-precision inertial sensors, especially low-cost inertial sensors, have been widely used in many fields, such as motion tracking, attitude-controlling system, and unmanned aircraft systems. Of them, low-cost MEMS triaxial accelerometers are widely used in pedestrian navigation, vehicle navigation, unmanned aerial vehicle (UAV), carrier rocket, and other fields [1][2][3][4][5]. However, MEMS triaxial accelerometers have errors due to inaccuracy in manufacturing techniques, such as zero-bias error, scalefactor error, and nonorthogonal error. e existence of these errors will affect the accuracy of the raw readings of the MEMS triaxial accelerometer, thus affecting the accuracy of inertial navigation system. erefore, for low-cost MEMS triaxial accelerometer sensors, calibration is a necessary step.
Under laboratory conditions, the turntable is usually used to calibrate the MEMS triaxial accelerometer [6,7]. However, MEMS triaxial accelerometer calibration is usually carried out without external auxiliary equipment in the field. In order to calibrate the MEMS triaxial accelerometer in the field test where it is inconvenient to use external auxiliary equipment such as a turntable, some methods have been proposed. In references [8][9][10], the 6-position method is used to collect the data of the MEMS triaxial accelerometer under the condition of static position and the error parameters are calibrated by the least square (LS) method. However, the scale-factor error and zero-bias error can be calibrated when the 6-position method is used, but the nonorthogonal error cannot be calibrated. Hence, a complete error model was constructed in references [11][12][13]. e LS method is used to estimate the error parameters, but these methods cannot obtain the optimal estimation because these methods assumed that the mean of the random noises is zero. erefore, in order to improve the accuracy of calibration, the unbiased objective estimation function is constructed by the maximum likelihood (ML) estimation method in references [14][15][16]. In references [14,15], the particle swarm optimization (PSO) algorithm is proposed to compute the objective function established by the ML estimation method, which is used to solve for the parameters. e initial value of iteration is not required, but the algorithm needs a longer time to obtain the calibrated results, so it is not suitable for practical application. In reference [16], the gradient descent method is proposed to estimate the error parameters, this process can be taken as the ML estimation problem. In each iteration, the error parameters are optimized by the gradient descent method. However, this method is relative to the initial rough values of the error parameters. If the initial rough values are not suitable, the estimation cannot get the global optimal solution. e step size and iteration times should also be taken into account when optimizing the objective function with the gradient descent method. If the sampling data are large, it is going to take a long time to compute in the gradient descent method.
Due to the aforementioned drawbacks of the traditional methods, this paper proposes a new calibration method. Firstly, the total error model function is determined, and then the ML estimation method was used for the optimal estimation. In this paper, different from the abovementioned methods, the Newton iteration method is used to optimize the objective function established by the ML estimation method.
e Newton iteration method has the advantages of fast convergence speed and global optimal solution. In practical application, the calibration of the MEMS triaxial accelerometer can be realized without external auxiliary equipment.

Error Analysis and Model Establishment of MEMS Triaxial Accelerometer
e raw readings of the MEMS triaxial accelerometer are affected by errors during manufacturing and installation; the errors are constructed by zero-bias error, scale-factor error and nonorthogonal error, and so on. e zero-bias error is the error caused by the MEMS triaxial accelerometer sensor outputs not being zero when the ideal outputs are zero, which is denoted by b a in this paper. e scale-factor error is mainly caused by the different sensitivities of each sensor's sensitive axis, which is represented by a diagonal matrix C s in (1), where s i (i � x, y, z) represents the scale-factor deviation of each axis [17]: e nonorthogonal error refers to the three sensitive axes that are not orthogonal to each other during the installation process of the MEMS triaxial accelerometer. A schematic diagram of the nonorthogonal error is given in Figure 1, where the (S x -S y -S z ) frame is the orthogonal frame and the (X-Y-Z) frame is the nonorthogonal frame. e nonorthogonal error can be represented by the upper triangular matrix C n in (2), which is calculated from the projection of the orthogonal frame to the nonorthogonal frame.
e MEMS triaxial accelerometer is not only affected by its own error factors but also affected by the measurement noise [18]. In general, the raw readings of the MEMS triaxial accelerometer can be expressed as (3). In (3), y n k represents the local gravity acceleration, which is approximately considered unchanged in a certain area under local navigation frame (n-frame) [19], C b n represents the rotation matrix from the navigation frame (n-frame) to the body frame (b-frame), and ε represents the Gaussian white noise with the variance of σ and the mean value of zero.
Define T a � C S C n and y b k � C b n y n k . e error model (3) can be rewritten as (4), which is the MEMS triaxial accelerometer reading generation model with error: Define R a � T −1 a , the error correction model of MEMS triaxial accelerometer under body frame (b-frame) is obtained as follows:

Maximum Likelihood Estimation Model Construction.
In this section, model (4) is used to construct the ML estimation model for parameters that need to be calibrated. In order to get the likelihood function L(θ), it is necessary to discuss the statistical properties of y m k . e random noise ε in model (4) is Gaussian white noise, which obeys the normal distribution. en, the probability density function of y m k can be obtained as follows: According to the theory of maximum likelihood estimation method, the likelihood function of a continuoustype random variable is where f(y m k ) is the probability density function of sample y m k and θ is the parameter to be estimated, θ � T a , b a , y c k . e idea of an ML estimate is to find the right θ such that the likelihood function is maximized: Finding the maximum value of likelihood function is equivalent to finding the minimum value of ; in order to facilitate the discussion, the constant term 2σ 2 was ignored and then the ML problem is expressed as (9). Equation (10) can be obtained by expanding the constraint conditions in (9) using model (5). Normalizing the MEMS triaxial accelerometer calibrated data, the modular square of the data should be one. Now, due to the noise, the modular square of the calibrated MEMS triaxial accelerometer data is influenced.
According to (11), μ c is not equal to zero; it means that c is a biased noise. If the noise c is taken as the unbiased noise, like the traditional LS method, the estimated results are not optimal. In order to establish a total error estimation model to perform the optimal estimate, a Lagrange constraint term needs to be introduced to eliminate the influence of the biased noises c. e minimization objective function is shown as (12), where λ k denotes the Lagrange multiplier used to constrain the two paradigms of the calibration value y c k . e estimated value at this time is θ , where vec(·) denotes the recomposition of a matrix into a column vector in the order of the columns, and the column vector does not include the three lower triangular elements of the matrix [20].

Newtonian Iteration Solution.
e iteration formula for the Newtonian iteration method is where θ (i+1) represents the (i + 1)-th iteration of parameter θ, θ (i) represents the i-th iteration of parameter θ, ∇f(θ (i) ) and ∇ 2 f(θ (i) ) represent the Jacobian matrix and Hessian matrix of the objective function in the i-th iteration process, respectively. e Jacobian matrix and Hessian matrix are described as follows: In the following formula, define d k � y m k − b a : Mathematical Problems in Engineering , where I denotes the 3 × 3 unit matrix and ⊗ denotes the Kronecker product. e abovementioned formulas are used in the Newton iteration method process. However, the Newton iteration method requires a good initial value to ensure that the iteration converges to a global minimum rather than a local minimum, so good initial estimate parameters are necessary for an accurate iteration result.

Initial Estimate.
In this section, the ellipsoid fitting method of initial estimation is discussed. Since the initial parameters are required to ensure the accuracy before doing the ML estimate, it is important to give accurate initial parameters. Ignoring the noise in (10) and expanding 1 � Bringing the vector y m k � X k Y k Z k T into (16) and then expanding and comparing with (17), we can get the following: Substituting the raw readings into (17), the parameters of the ellipsoid can be obtained by the LS method, and then the initial parameters of calibration are obtained by is the Cholesky decomposition. (3) and (4), the parameters of the simulation are taken to be as (20). e attitude transformation matrix C b n , which is converted from the navigation frame (n-frame) to the body frame (b-frame), is generated by Euler angles. e pitch angle, swing angle, and heading angle are taken to be 120, 120, and 80 degrees, respectively; to simulate the generation of real angle data, the data information was multiplied by random numbers to generate 300 sets of angle data, and then 300 sets of the attitude conversion matrix are converted. e random errors are set as Gaussian white noises with the mean of zero and a standard deviation of 0.1.

Simulation Analysis. According to models
e ideal model of the vectors measured by the MEMS triaxial accelerometer after normalization should be a unit sphere, the center of the sphere is the coordinate origin, and the radius is 1. e model fitted from the measured and calibrated data is represented by the black ellipsoid and the data are represented by the red points in Figures 2(a) and 2(b). It can be observed that the calibrated ellipsoid is a unit sphere with its center at the coordinate origin.
In Figure 3, the MEMS triaxial accelerometer data calibrated by the LS method and the data calibrated by the ML estimation method are significantly improved compared with the uncalibrated MEMS triaxial accelerometer data. e modular square of the MEMS triaxial accelerometer calibrated by the ML estimation method fluctuates less and is closer to one than that of the LS method. Statistically, in Figure 3, the means of modular square are 1.172 (uncalibrated), 0.9945 (calibrated by the LS method), and 1 (calibrated by the ML method), respectively. e standard deviations of modular square are 0.18 (uncalibrated), 0.1052 (calibrated by the LS method), and 0.0403 (calibrated by the ML method), respectively. Both the LS method and the ML method have obvious effects on the calibration, but the results calibrated by the ML estimation method have less fluctuation in the vicinity of modular square of one than the results calibrated by the LS method. e parameters set in this paper and the parameters calibrated by the LS method and the ML estimation method are compared in Table 1.    e experiment results are shown in Figure 4, which shows that the objective function values reach the convergence state after the third iteration. To further highlight the advantages of the ML estimation method over the LS method, the equations in (20) are defined as a metric and 50 Monte Carlo experiments are performed. In (20), T ac and b ac , respectively, represent the parameters calibrated and T a and b a , respectively, represent parameters set in this paper. In Figures 5 and 6, the variation of the objective function specified during the Monte Carlo experiment can be observed. e results of e T and e b , calculated by the parameters calibrated by the LS method, are both represented by the black solid lines. Also, the mean values of the results calibrated by the LS method are both represented by black dashed lines. Similarly, the results of e T and e b , calculated by the parameters calibrated by the ML method, are both represented by the red solid lines and the mean values of the results calibrated by the ML method are both represented by red dashed lines. In Figures 5 and 6, the parameters calibrated by the ML estimation method in the Monte Carlo experiment are closer to the parameters set in this paper than those calibrated by the LS method. In Figure 7, we can see that although the modular square of the data calibrated by the LS method is close to 1, the modular square of the data calibrated by the ML estimation method fluctuates slightly at 1. According to the analysis of Figures 5-7, it can be seen that the ML estimation method is more stable and accurate under the influence of noise than the LS method.

Experimental Analysis.
is section will verify the proposed method through practical experiments. Figures 8(a) and 8(b) show the experimental equipment and sensors of this paper, including the 3D printed frame, MPU9250 module, bubble level, and black desktop. e MPU9250 module is fixed to the 3D printed frame (the angle between the MPU9250 module and the 3D printed frame is not considered). e  bubble level is used to observe the status of the black desktop. e sampling frequency of the MPU9250 module is 100 Hz, and the data are transmitted to the master computer through a serial port and then received and stored in the master computer. e 3D printed frame is placed on the black desktop.
When the x-axis of the 3D frame is up, the 3D printed frame rotates horizontally by 90 degrees at regular intervals, so the static data of the accelerometer can be obtained when the x-axis is up. Similarly, static data can be obtained when the x-axis is down, when the y-axis is up and down, and when the z-axis is up and down. e sensor's static data at 24 positions are recorded as data 1. e 3D printed frame is held in hand for slow movement, and the data obtained during the movement are recorded as data 2.
Calibrating data 1 and data 2 by the LS method and ML estimation method, the modular square results of Figures 9  and 10 are obtained, respectively. e statistical information can be obtained from Table 2. In Table 2, the mean of modular square of the two groups of sensors data calibrated by the LS method is very close to 1, while the mean of modular square calibrated by the ML method is 1, which is closer to the theoretical value 1, so it is more accurate. e standard deviation after calibration by the ML estimation method is significantly reduced compared with that   Before calibration After LS calibration After ML calibration Figure 9: Comparison of modular square before and after calibrating measured data 1.   Mathematical Problems in Engineering calibrated by the LS method. In Figures 9 and 10, the pink curve fluctuates less than the green curve and fluctuates around 1. e two sets of sensors data and calibration data are projected onto the ellipsoid, which is plotted in Figures 11 and 12. By comparing Figures 11(a), 11(b), 12(a), and 12(b), after the calibration, the uncalibrated ellipsoid in black has been calibrated as a unit sphere, and the center of sphere is at the origin of the frame, which verifies the effectiveness of the calibration method proposed in this paper.   X -a x is Y -a x is (b) Figure 12: (a) e ellipsoid fitted by data 2 (black) and the data 2 points (red). (b) e ellipsoid fitted by data 2 calibrated by the ML estimation method (black) and the points of data 2 calibrated by the ML estimation method (red).

Conclusions
In this paper, a complete error model of MEMS triaxial accelerometers is constructed. e influence of biased noises of the traditional objective function model, which is constructed by modulus calculation, is eliminated. To eliminate this influence, a calibration method based on the ML estimation method is proposed in this paper. e total error unbiased objective function is constructed, and the optimal estimation is obtained by the Newton iteration method. To make the Newton iteration method converge to the optimal value, an ellipsoid fitting method is adopted to get the initial parameters for the Newton iteration method. In the simulation, the results of the ellipsoid model and modular square verify the validity of the ML estimation method. Comparing using the Monte Carlo experiment, the results calibrated by the ML estimation method are better and more stable than the LS method. e experiment verifies the feasibility of the ML estimation method under different sets of data collected by the MPU9250 module. e method proposed in this paper is not only limited to calibrating the MEMS triaxial accelerometer but also has potential calibration effects for sensors of the same triaxial structure, such as magnetometers.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.