1 Motivation

Several approaches based on machine learning techniques have been proposed to increase MEMS inertial sensors’ performances. Generally, such approaches have assumed that errors cannot be explained by linear models, as well as the observed error probability distribution is not pure Gaussian. However, according to manufacturers the error non-linearity observed in MEMS sensors can be considered in most cases negligible. Therefore, it seems necessary to evaluate linear models in the context of MEMS inertial sensor errors. In particular the application of Time delayed Multiple Linear Regression models. (TD-MLR)

2 Main Hypotesis

Simpler Machine learning linear models are capable of compensating errors observed in MEMS sensors.

3 Models considered

Three models are considered:

  1. Time Delayed - Multi Linear Regression (TD-MLR) : A simple linear model
  2. Multi Layer Perceptron (MLP), a non linear model
  3. Moving Average (MA): the common strategy used in navigation systems for dealing with noise

A brief description of each one of the models is presented in the following subsections.

3.1 Multi Linear Regression Mode (MLR).

A TD-MLR model for representing the relationship between an explanatory variable \(x\) and a response variable \(y\) can be written as,

\[y_t = \alpha + \beta_t L^q( x_{t}) + \epsilon_{t}\]

where, the subscript \(t\) refers to the \(t^{th}\) time unit in the population of size \(M \in \mathbb{N}\). \(y_t\) denotes the observed response for experimental time unit \(t\). \(\epsilon_t\) and \(\alpha\) are respectively the residual error and the intercept terms for the same time unit \(t\). The expression \(L^q( x_{t})\) represents a tapped delay line for the \(x_t\) variable. A tapped delay line is a well-known method to model non-stationary time series implemented as a short-term memory with \(q\) unit delay operators

3.2 Multi Layer Perceptron (MLP)

(From Wikipedia) A multilateral perceptron (MLP) is a feed forward artificial neural network model that maps sets of input data onto a set of appropriate outputs. An MLP consists of multiple layers of nodes in a directed graph, with each layer fully connected to the next one. Except for the input nodes, each node is a neuron (or processing element) with a nonlinear activation function. MLP utilizes a supervised learning technique called back-propagation for training the network.[1][2] MLP is a modification of the standard linear perceptron and can distinguish data that are not linearly separable.

3.3 Moving Average (MA)

The usual technique for dealing with correlated noise consists of applying a Moving Average (MA) filter. The MA filter operates by averaging a predefined number of points from the input signal to produce each point in the output signal. As its name implies, a MA is an average that moves along the input signal values. An old value is dropped as a new value comes available. The number of values to keep will depend of the size of the time window considered. More formally, given an input signal \(x\) of size \(M\) defined as \({x_t}_{(t=1)}^{M}\), a MA with a time window of size \(q\) will output a new signal \({y_t}_{(t=1)}^{M-q}\) defined from the \(x_i\) by taking the arithmetic mean of sub-sequences of \(q\) terms,

\[y_t = \frac{1}{q} \sum_{s=0}^{M-q+1} x_{t-s}\]

4 Experiments

4.1 Dataset Description

The performance of TD-MLR models is evaluated on the well-known Melbourne dataset , which was kindly provided by research members from The Ohio State University (OSU) and the University of Melbourne. Since its publication in 2011, the Melbourne dataset has become a reference in the research field of MEMS inertial sensors.

The dataset contains measures provided by a mobile platform that includes various grades of IMUs, several GPS receivers and data-logging software. The IMUs included in the mobile platform varies from low-quality MEMS sensors to high-quality, navigation-grade sensors. The information provided by IMUs were logged during a predefined trajectory performed by a ground vehicle. As shown in Fig. , the trajectory has two well-defined stretches.

  trajectory1=read.csv("/home/harpo/Dropbox/shared/MEMS-ANN/datasets/dgps_T1.txt",header=F)
  trajectory2=read.csv("/home/harpo/Dropbox/shared/MEMS-ANN/datasets/dgps_T2.txt",header=F)
  par(mfrow=c(1,2))
  plot(trajectory1$V1,trajectory1$V2,type="l",lwd=4,ylab="Latitude [degree]",xlab="Longitude [degree]",main="Stretch 1 (T1)",col="orange",cex=1.2)
  plot(trajectory2$V1,trajectory2$V2,type="l",lwd=4,ylab="Latitude [degree]",xlab="Longitude [degree]",main="Stretch 2 (T2)",col="orange",cex=1.2)

In the first stretch (T1, for short) the vehicle performs several 8-pattern loops at various velocities and accelerations in a parking area at OSU for about 16 minutes. Then, in the second stretch (T2) the vehicle is driven on an internal OSU road for an 8 minutes trip. Notice that T1 has more rich dynamics than T2. The complete dataset covers about 24 minutes.

In the present work, four MEMS IMUs are considered:

  1. Gladiator Landmark 10 (Gladiator).
  2. Xbow IMU400CD (Crossbow)
  3. Xsens (xsns1)
  4. Xtal (xtal)

These four IMUs represent MEMS inertial sensors of different grades, Gladiator is a low quality IMU very sensible to noise, while the remaining are mid-range quality IMUS with variable sensor characteristics . A fifth Non-MEMS IMU is also considered as reference. Namely, the Honeywell H764G-1 (Honeywell), a well-established navigation-grade IMU in use on most military aircraft

4.2 Finding the optimal models

TD-MLR and MLP models are built for each one of the available sensors (i.e AccX, AccY, GyroZ). The models are built using the response observations (i.e. \(y_t\)) from nav-grade Honeywell IMU while the explanatory variables (i.e. \(X_t\)) will consist of time delayed values coming from the corresponding low and mid quality IMUs sensors (i.e Gladiator, Xbow, Xtal and Xsens).

Either for TD-MLR and for MA models, the optimal size of \(q\) should be determined. The strategy followed in this work consists of starting with a very short \(q\) value and incrementing it by one period. For each new \(q\) value the Root Square Mean Error (RMSE) of the resulting model is calculated. For improving RSME estimate, a k-fold cross validation with \(k=10\) is used to evaluate the models for every new lag term. The process is repeated until \(q\) reaches a predefined and sufficiently large value. Then, the optimal value for \(q\) is determined by the model with the lower RSME.

In the case of MLP, not only \(q\) but also the different number of hidden neurons \(h\) are evaluated. In this case, a similar approach is followed for a subset of the of possible network typologies (i.e. \(h\)), again a k-fold cross validation with \(k=10\) is used to evaluate the performance of every \(q\) for each one of the selected \(h\). Similarly to previous models, the model with the lower RMSE is used for determining the optimal values of \(q\) and \(h\)

Following standard machine learning methodology, the models will be adjusted and evaluated using data from Stretch T1. Then, after the optimal lag length (\(q\)) for each model had been selected, an independent evaluation on stretch T2 is conducted in order to get a better estimation of the error on unseen data.

4.2.1 Results for the X1 (AccX) sensor

The following table show the optimal \(q\) value for each model (MLR,MA and MLP) on the four selected IMUS. The \(q\) values was obtained according to the strategy described in the previous section.

best_models %>% filter(sensor=="X1") %>% arrange(imuid)
##    imuid model sensor timedelay  rsme_mean      rsme_sd
## 1      5 mlp10     X1        80 0.18425508 0.0027877182
## 2      5   mlr     X1        60 0.19033079 0.0017188523
## 3      5    ma     X1        40 0.21240464 0.0017426428
## 4      6  mlp5     X1        30 0.09815762 0.0020489154
## 5      6   mlr     X1        25 0.09632335 0.0005943079
## 6      6    ma     X1        15 0.12967777 0.0009625942
## 7      8 mlp30     X1        40 0.07542331 0.0005309780
## 8      8   mlr     X1        10 0.07592635 0.0003954836
## 9      8    ma     X1         1 0.08199104 0.0007452200
## 10     9 mlp10     X1        10 0.10235128 0.0017000333
## 11     9   mlr     X1        15 0.09939805 0.0008778526
## 12     9    ma     X1         1 0.12979918 0.0005084080

The plot shows the RSME distributions for sensor X1 provided by TD-MLR, MA and MLP considering different values for \(h\) (i.e. topologies)

bw_data=inner_join(best_models %>% filter(sensor=="X1"),results, by=c("timedelay","imuid","model"))
## Warning in inner_join_impl(x, y, by$x, by$y, suffix$x, suffix$y): joining
## factors with different levels, coercing to character vector
bwplot(X1~model|imu,groups=timedelay,data=bw_data,scale='free'
       ,auto.key = F
       ,panel = panel.superpose
       ,panel.groups=function(x,y,group.number,...){
         xt <- x[x==min(x)] # find latest year
         yt <- y[y==min(y)] # find value at latest year 
         panel.text(x,yt,labels=levels(factor(bw_data$timedelay))[group.number],pos=1,...)
         panel.bwplot(x,y,...)
             },ylab="X1"
         )

4.2.2 Results for the X2 (AccY) sensor

The following table show the optimal \(q\) value for each model (MLR,MA and MLP) on the four selected IMUS. The \(q\) values was obtained according to the strategy described in the previous section.

best_models %>% filter(sensor=="X2") %>% arrange(imuid)
##    imuid model sensor timedelay  rsme_mean      rsme_sd
## 1      5  mlp5     X2        35 0.22548026 0.0041683144
## 2      5   mlr     X2        40 0.22510823 0.0017614691
## 3      5    ma     X2        30 0.60081515 0.0018267406
## 4      6  mlp5     X2        10 0.11772926 0.0019768169
## 5      6   mlr     X2        25 0.10840339 0.0008605284
## 6      6    ma     X2        15 0.15593408 0.0010565909
## 7      8 mlp40     X2         5 0.08914980 0.0013538570
## 8      8   mlr     X2        10 0.08961391 0.0006216814
## 9      8    ma     X2         1 0.11238620 0.0010835257
## 10     9 mlp10     X2         5 0.13556195 0.0030066884
## 11     9   mlr     X2        20 0.12320711 0.0010895449
## 12     9    ma     X2         1 0.15872687 0.0016744040

The plot shows the RSME distributions for sensor X2 provided by TD-MLR, MA and MLP considering different values for \(h\) (i.e. topologies)

bw_data=inner_join(best_models %>% filter(sensor=="X2"),results, by=c("timedelay","imuid","model"))
## Warning in inner_join_impl(x, y, by$x, by$y, suffix$x, suffix$y): joining
## factors with different levels, coercing to character vector
bwplot(X2~model|imu,groups=timedelay,data=bw_data,scale='free'
       ,auto.key = F
       ,panel = panel.superpose
       ,panel.groups=function(x,y,group.number,...){
         xt <- x[x==min(x)] # find latest year
         yt <- y[y==min(y)] # find value at latest year 
         panel.text(x,yt,labels=levels(factor(bw_data$timedelay))[group.number],pos=1,...)
         panel.bwplot(x,y,...)
             },ylab="X2"
         )

4.2.3 Results for the X6 (GiroZ) sensor

The following table show the optimal \(q\) value for each model (MLR,MA and MLP) on the four selected IMUS. The \(q\) values was obtained according to the strategy described in the previous section.

best_models %>% filter(sensor=="X6") %>% arrange(imuid)
##    imuid  model sensor timedelay   rsme_mean      rsme_sd
## 1      5   mlp5     X6       100 0.015165424 1.015794e-04
## 2      5    mlr     X6        90 0.015824808 6.944545e-05
## 3      5     ma     X6        30 0.017012517 1.255337e-04
## 4      6 mlp100     X6        40 0.002399567 2.757186e-05
## 5      6    mlr     X6        55 0.002056533 1.294075e-05
## 6      6     ma     X6         5 0.004202308 1.969213e-05
## 7      8  mlp80     X6        50 0.004528022 3.006731e-05
## 8      8    mlr     X6        60 0.004247821 2.011378e-05
## 9      8     ma     X6        10 0.010260791 1.898733e-05
## 10     9 mlp100     X6        40 0.003815108 2.781940e-05
## 11     9    mlr     X6        45 0.003341583 2.151424e-05
## 12     9     ma     X6         1 0.005300210 3.498887e-05

The plot shows the RSME distributions for sensor X6 provided by TD-MLR, MA and MLP considering different values for \(h\) (i.e. topologies)

bw_data=inner_join(best_models %>% filter(sensor=="X6"),results, by=c("timedelay","imuid","model"))
## Warning in inner_join_impl(x, y, by$x, by$y, suffix$x, suffix$y): joining
## factors with different levels, coercing to character vector
bwplot(X6~model|imu,groups=timedelay,data=bw_data,scale='free'
       ,auto.key = F
       ,panel = panel.superpose
       ,panel.groups=function(x,y,group.number,...){
         xt <- x[x==min(x)] # find latest year
         yt <- y[y==min(y)] # find value at latest year 
         panel.text(x,yt,labels=levels(factor(bw_data$timedelay))[group.number],pos=1,...)
         panel.bwplot(x,y,...)

             },ylab="X6"
      
         )

The plot shows the RSME distributions for sensor X1 provided by TD-MLR, MA and MLP considering different values for \(h\) (i.e. topologies)

## reshape results for selecting best models for ma and MLP
results_reduced=results %>% select(imuid,model,timedelay,X1,X2,X6)
results_reduced=melt(results_reduced,id.vars=c("imuid","model","timedelay"),variable.name="sensor",value.name="value")
best_taps=group_by(results_reduced,model,imuid,sensor,timedelay) %>% summarise(mean=mean(value),sd=sd(value)) %>% filter(mean==min(mean)) %>% arrange(desc(imuid)) %>% select(imuid,sensor,timedelay)
## Adding missing grouping variables: `model`
best_taps=reshape2::dcast(best_taps,imuid+model~sensor,value.var = 'timedelay')

4.2.4 Final Models

In summary, after a 10-Fold Cross-validation the models with the lower RMSE are:

4.2.4.1 For MLP

best_models %>% filter(grepl("mlp",model) & sensor=='X1')
##   imuid model sensor timedelay  rsme_mean     rsme_sd
## 1     5 mlp10     X1        80 0.18425508 0.002787718
## 2     6  mlp5     X1        30 0.09815762 0.002048915
## 3     8 mlp30     X1        40 0.07542331 0.000530978
## 4     9 mlp10     X1        10 0.10235128 0.001700033
best_models %>% filter(grepl("mlp",model) & sensor=='X2')
##   imuid model sensor timedelay rsme_mean     rsme_sd
## 1     5  mlp5     X2        35 0.2254803 0.004168314
## 2     6  mlp5     X2        10 0.1177293 0.001976817
## 3     8 mlp40     X2         5 0.0891498 0.001353857
## 4     9 mlp10     X2         5 0.1355619 0.003006688
best_models %>% filter(grepl("mlp",model) & sensor=='X6')
##   imuid  model sensor timedelay   rsme_mean      rsme_sd
## 1     5   mlp5     X6       100 0.015165424 1.015794e-04
## 2     6 mlp100     X6        40 0.002399567 2.757186e-05
## 3     8  mlp80     X6        50 0.004528022 3.006731e-05
## 4     9 mlp100     X6        40 0.003815108 2.781940e-05

4.2.4.2 For TD-MLR

best_models %>% filter(grepl("mlr",model) & sensor=='X1')
##   imuid model sensor timedelay  rsme_mean      rsme_sd
## 1     5   mlr     X1        60 0.19033079 0.0017188523
## 2     6   mlr     X1        25 0.09632335 0.0005943079
## 3     8   mlr     X1        10 0.07592635 0.0003954836
## 4     9   mlr     X1        15 0.09939805 0.0008778526
best_models %>% filter(grepl("mlr",model) & sensor=='X2')
##   imuid model sensor timedelay  rsme_mean      rsme_sd
## 1     5   mlr     X2        40 0.22510823 0.0017614691
## 2     6   mlr     X2        25 0.10840339 0.0008605284
## 3     8   mlr     X2        10 0.08961391 0.0006216814
## 4     9   mlr     X2        20 0.12320711 0.0010895449
best_models %>% filter(grepl("mlr",model) & sensor=='X6')
##   imuid model sensor timedelay   rsme_mean      rsme_sd
## 1     5   mlr     X6        90 0.015824808 6.944545e-05
## 2     6   mlr     X6        55 0.002056533 1.294075e-05
## 3     8   mlr     X6        60 0.004247821 2.011378e-05
## 4     9   mlr     X6        45 0.003341583 2.151424e-05

4.2.4.3 For MA

best_models %>% filter(grepl("ma",model) & sensor=='X1')
##   imuid model sensor timedelay  rsme_mean      rsme_sd
## 1     5    ma     X1        40 0.21240464 0.0017426428
## 2     6    ma     X1        15 0.12967777 0.0009625942
## 3     8    ma     X1         1 0.08199104 0.0007452200
## 4     9    ma     X1         1 0.12979918 0.0005084080
best_models %>% filter(grepl("ma",model) & sensor=='X2')
##   imuid model sensor timedelay rsme_mean     rsme_sd
## 1     5    ma     X2        30 0.6008151 0.001826741
## 2     6    ma     X2        15 0.1559341 0.001056591
## 3     8    ma     X2         1 0.1123862 0.001083526
## 4     9    ma     X2         1 0.1587269 0.001674404
best_models %>% filter(grepl("ma",model) & sensor=='X6')
##   imuid model sensor timedelay   rsme_mean      rsme_sd
## 1     5    ma     X6        30 0.017012517 1.255337e-04
## 2     6    ma     X6         5 0.004202308 1.969213e-05
## 3     8    ma     X6        10 0.010260791 1.898733e-05
## 4     9    ma     X6         1 0.005300210 3.498887e-05

4.3 Evaluation on T2 Stretch

The 100% of the training dataset (T1) is used for building the models and then evaluated on unseen data. For performing the evaluation the model on unseen data, a sampling method is used on the Stretch T2. In this case, a fixed window encompasing the 10% of the Stretch T2 is randomly and uniformly sampled with replacement. Such process is repeteated 40 times.

load("/home/harpo/Dropbox/shared/MEMS-ANN/results/results_samples_mlr_own_method.Rda")
# Selecting best models for MA -----------------------------------------
ma_imus=list()
ma_imus[[5]]=(best_models %>% filter(imuid==5 & model=='ma'))$timedelay
ma_imus[[6]]=(best_models %>% filter(imuid==6 & model=='ma'))$timedelay
ma_imus[[8]]=(best_models %>% filter(imuid==8 & model=='ma'))$timedelay
ma_imus[[9]]=(best_models %>% filter(imuid==9 & model=='ma'))$timedelay

# Selecting best models for MLR 
mlr_imus=list()
mlr_imus[[5]]=(best_models %>% filter(imuid==5 & model=='mlr'))$timedelay
mlr_imus[[6]]=(best_models %>% filter(imuid==6 & model=='mlr'))$timedelay
mlr_imus[[8]]=(best_models %>% filter(imuid==8 & model=='mlr'))$timedelay
mlr_imus[[9]]=(best_models %>% filter(imuid==9 & model=='mlr'))$timedelay


# Selecting best models por MLP  using ANOVA
mlp_imus=list()

mlp_imus[[5]]=(best_models %>% filter(imuid==5 & grepl('mlp',model)))$timedelay
mlp_imus[[6]]=(best_models %>% filter(imuid==6 & grepl('mlp',model)))$timedelay
mlp_imus[[8]]=(best_models %>% filter(imuid==8 & grepl('mlp',model)))$timedelay
mlp_imus[[9]]=(best_models %>% filter(imuid==9 & grepl('mlp',model)))$timedelay


# Selecting the best MLP models
mlp_models=list()
mlp_models[[5]]=(best_models %>% filter(imuid==5 & grepl('mlp',model)))$model
mlp_models[[6]]=(best_models %>% filter(imuid==6 & grepl('mlp',model)))$model
mlp_models[[8]]=(best_models %>% filter(imuid==8 & grepl('mlp',model)))$model
mlp_models[[9]]=(best_models %>% filter(imuid==9 & grepl('mlp',model)))$model

4.3.1 Comparing MLR vs Raw

This section analyze the results in terms of RMSE for the TD-MLR compared with the original RMSE.

4.3.1.1 Evaluation on Gladiator

glad_test_results=test_results %>% filter(imuid==5 &  model=='mlr') %>%select(imuid,X1,X1vsRaw,X2,X2vsRaw,X6,X6vsRaw)
names(glad_test_results)=c("imuid","X1_mlr","X1_raw","X2_mlr","X2_raw","X6_mlr","X6_raw")
glad_test_results=tidyr::gather(glad_test_results,"sensor","value",2:7)
glad_test_results=tidyr::separate(glad_test_results,sensor,c("sensor","source"),sep="_") 
px1=wilcox.test(value~source,data=glad_test_results %>% filter(sensor=="X1" & (source=="mlr"|source=="raw")))
px2=wilcox.test(value~source,data=glad_test_results %>% filter(sensor=="X2" & (source=="mlr"|source=="raw")))
px6=wilcox.test(value~source,data=glad_test_results %>% filter(sensor=="X6" & (source=="mlr"|source=="raw")))

The Figure below show the error distribution for the original Gladiator and the TD-MLR on the three sensors

A man-Whitney-Wilcox test is conducted for the three sensors and results are presented in the following table As can be observed the differences are statistically significant for the three sensors

Sensor Name X1 X2 X6
P Value 1.860340410^{-23} 1.860340410^{-23} 1.860340410^{-23}

4.3.1.2 Evaluation on Xbow

xbow_test_results=test_results %>% filter(imuid==6 &  model=='mlr') %>%select(imuid,X1,X1vsRaw,X2,X2vsRaw,X6,X6vsRaw)
names(xbow_test_results)=c("imuid","X1_mlr","X1_raw","X2_mlr","X2_raw","X6_mlr","X6_raw")
xbow_test_results=tidyr::gather(xbow_test_results,"sensor","value",2:7)
xbow_test_results=tidyr::separate(xbow_test_results,sensor,c("sensor","source"),sep="_") 
px1=wilcox.test(value~source,data=xbow_test_results %>% filter(sensor=="X1" & (source=="mlr"|source=="raw")))
px2=wilcox.test(value~source,data=xbow_test_results %>% filter(sensor=="X2" & (source=="mlr"|source=="raw")))
px6=wilcox.test(value~source,data=xbow_test_results %>% filter(sensor=="X6" & (source=="mlr"|source=="raw")))

The Figure below show the error distribution for the original Xbow and the TD-MLR on the three sensors

bwplot(value~source|sensor,data=xbow_test_results,groups=source,layout=c(3,1),scales=list(relation='free'))

A man-Whitney-Wilcox test is conducted for the three sensors and results are presented in the following table As can be observed the differences are statistically significant for the three sensors

Sensor Name X1 X2 X6
P Value 6.935339110^{-6} 5.002387410^{-8} 1.860340410^{-23}

4.3.1.3 Evaluation on Xsens

xsens_test_results=test_results %>% filter(imuid==8 &  model=='mlr') %>%select(imuid,X1,X1vsRaw,X2,X2vsRaw,X6,X6vsRaw)
names(xsens_test_results)=c("imuid","X1_mlr","X1_raw","X2_mlr","X2_raw","X6_mlr","X6_raw")
xsens_test_results=tidyr::gather(xsens_test_results,"sensor","value",2:7)
xsens_test_results=tidyr::separate(xsens_test_results,sensor,c("sensor","source"),sep="_") 
px1=wilcox.test(value~source,data=xsens_test_results %>% filter(sensor=="X1" & (source=="mlr"|source=="raw")))
px2=wilcox.test(value~source,data=xsens_test_results %>% filter(sensor=="X2" & (source=="mlr"|source=="raw")))
px6=wilcox.test(value~source,data=xsens_test_results %>% filter(sensor=="X6" & (source=="mlr"|source=="raw")))

The Figure below show the error distribution for the original Xsens and the TD-MLR on the three sensors

bwplot(value~source|sensor,data=xsens_test_results,groups=source,layout=c(3,1),scales=list(relation='free'))

A man-Whitney-Wilcox test is conducted for the three sensors and results are presented in the following table As can be observed the differences are statistically significant for the three sensors

Sensor Name X1 X2 X6
P Value 0.005236 2.556541410^{-6} 1.860340410^{-23}

4.3.1.4 Evaluation on Xtal

xtal_test_results=test_results %>% filter(imuid==9 &  model=='mlr') %>%select(imuid,X1,X1vsRaw,X2,X2vsRaw,X6,X6vsRaw)
names(xtal_test_results)=c("imuid","X1_mlr","X1_raw","X2_mlr","X2_raw","X6_mlr","X6_raw")
xtal_test_results=tidyr::gather(xtal_test_results,"sensor","value",2:7)
xtal_test_results=tidyr::separate(xtal_test_results,sensor,c("sensor","source"),sep="_") 
px1=wilcox.test(value~source,data=xtal_test_results %>% filter(sensor=="X1" & (source=="mlr"|source=="raw")))
px2=wilcox.test(value~source,data=xtal_test_results %>% filter(sensor=="X2" & (source=="mlr"|source=="raw")))
px6=wilcox.test(value~source,data=xtal_test_results %>% filter(sensor=="X6" & (source=="mlr"|source=="raw")))

The Figure below show the error distribution for the original Xtal and the TD-MLR on the three sensors

bwplot(value~source|sensor,data=xtal_test_results,groups=source,layout=c(3,1),scales=list(relation='free'))

A man-Whitney-Wilcox test is conducted for the three sensors and results are presented in the following table As can be observed the differences are statistically significant for the three sensors

Sensor Name X1 X2 X6
P Value 8.51802710^{-5} 1.117076710^{-6} 1.860340410^{-23}

4.3.2 Comparing MLR vs MLP

The proposed method TD-MLR is now compared with an standard MLP on the four IMUS. Notice that model parameters result from Cross-validation on Stretch T1. For TD-MLR and MLP, the models with the Lower RSME are selected.

4.3.2.1 Evaluation on Gladiator

The Figure below show the error distribution for the TD-MLR and MLP for the three sensors on Gladiator

bwplot(value~model|sensor,data=glad_test_results %>% filter((model=="mlp"|model=="mlr")),groups=factor(model),layout=c(3,1),scales=list(relation='free'))

px1=wilcox.test(value~model,data=glad_test_results %>% filter(sensor=="X1" & (model=="mlp"|model=="mlr")))
px2=wilcox.test(value~model,data=glad_test_results %>% filter(sensor=="X2" & (model=="mlp"|model=="mlr")))
px6=wilcox.test(value~model,data=glad_test_results %>% filter(sensor=="X6" & (model=="mlp"|model=="mlr")))

A man-Whitney-Wilcox test is conducted for the three sensors and the results are presented in the following table As can be observed the differences are NOT statistically significant for the three sensors

Sensor Name X1 X2 X6
P Value 0.8972211 0.5495696 0.9351801

4.3.2.2 Evaluation on Xbow

xbow_test_results=test_results %>% select(model,imuid,sample,X1,X2,X6) %>% filter(imuid==6 & 
  (   model=='mlr' | 
      model==as.character(mlp_models[[6]][1]) | 
      model==as.character(mlp_models[[6]][2]) |
      model==as.character(mlp_models[[6]][3]) 
  ))

xbow_test_results=melt(xbow_test_results,id.vars=c("imuid","model","sample"),variable.name="sensor",value.name="value")

xbow_test_results=xbow_test_results %>% filter(model=='mlr' | 
  ( model==as.character(mlp_models[[6]][1]) & sensor=="X1") | 
  ( model==as.character(mlp_models[[6]][2]) & sensor=="X2") | 
  ( model==as.character(mlp_models[[6]][3]) & sensor=="X6")  
    )
xbow_test_results$model=factor(xbow_test_results$model, levels=c(levels(xbow_test_results$model), "mlp"))
# rename mlpxxx to mlp
xbow_test_results[grepl("mlp",xbow_test_results$model),]$model='mlp'

The Figure below show the error distribution for the TD-MLR and MLP for the three sensors on Xbow

bwplot(value~model|sensor,data=xbow_test_results %>% filter( model=="mlp"|model=="mlr") ,groups=factor(model),layout=c(3,1),scales=list(relation='free'))

px1=wilcox.test(value~model,data=xbow_test_results %>% filter(sensor=="X1" & (model=="mlp"|model=="mlr")))
px2=wilcox.test(value~model,data=xbow_test_results %>% filter(sensor=="X2" & (model=="mlp"|model=="mlr")))
px6=wilcox.test(value~model,data=xbow_test_results %>% filter(sensor=="X6" & (model=="mlp"|model=="mlr")))

A man-Whitney-Wilcox test is conducted for the three sensors and the results are presented in the following table As can be observed the differences are NOT statistically significant for the X1 and X2 sensors. However in this case the difference for the X6 sensor are significant with a p.value <0.05

Sensor Name X1 X2 X6
P Value 0.9123817 0.0773742 0.0350024

4.3.2.3 Evaluation on XSENS

xsens_test_results=test_results %>% select(model,imuid,sample,X1,X2,X6) %>% filter(imuid==5 & 
  (   model=='mlr' | 
      model==as.character(mlp_models[[8]][1]) | 
      model==as.character(mlp_models[[8]][2]) |
      model==as.character(mlp_models[[8]][3]) 
  ))

xsens_test_results=melt(xsens_test_results,id.vars=c("imuid","model","sample"),variable.name="sensor",value.name="value")

xsens_test_results=xsens_test_results %>% filter(model=='mlr' | 
  ( model==as.character(mlp_models[[8]][1]) & sensor=="X1") | 
  ( model==as.character(mlp_models[[8]][2]) & sensor=="X2") | 
  ( model==as.character(mlp_models[[8]][3]) & sensor=="X6") 
    )
xsens_test_results$model=factor(xsens_test_results$model, levels=c(levels(xsens_test_results$model), "mlp"))
# rename mlpxxx to mlp
xsens_test_results[grepl("mlp",xsens_test_results$model),]$model='mlp'

The Figure below show the error distribution for the TD-MLR and MLP for the three sensors on Xsens

bwplot(value~model|sensor,data=xsens_test_results %>% filter( model=="mlp"|model=="mlr") ,groups=factor(model),layout=c(3,1),scales=list(relation='free'))

px1=wilcox.test(value~model,data=xsens_test_results %>% filter(sensor=="X1" & (model=="mlp"|model=="mlr")))
px2=wilcox.test(value~model,data=xsens_test_results %>% filter(sensor=="X2" & (model=="mlp"|model=="mlr")))
px6=wilcox.test(value~model,data=xsens_test_results %>% filter(sensor=="X6" & (model=="mlp"|model=="mlr")))

A man-Whitney-Wilcox test is conducted for the three sensors and results are presented in the following table As can be observed the differences are NOT statistically significant for the three sensors

Sensor Name X1 X2 X6
P Value 0.9275738 0.5431922 0.7630755

4.3.2.4 Evaluation on Xtal

xtal_test_results=test_results %>% select(model,imuid,sample,X1,X2,X6) %>% filter(imuid==5 & 
  (   model=='mlr' | 
      model==as.character(mlp_models[[9]][1]) | 
      model==as.character(mlp_models[[9]][2]) |
      model==as.character(mlp_models[[9]][3]) 
  ))

xtal_test_results=melt(xtal_test_results,id.vars=c("imuid","model","sample"),variable.name="sensor",value.name="value")

xtal_test_results=xtal_test_results %>% filter(model=='mlr' | 
  ( model==as.character(mlp_models[[9]][1]) & sensor=="X1") | 
  ( model==as.character(mlp_models[[9]][2]) & sensor=="X2") | 
  ( model==as.character(mlp_models[[9]][3]) & sensor=="X6") 
    )
xtal_test_results$model=factor(xtal_test_results$model, levels=c(levels(xtal_test_results$model), "mlp"))
# rename mlpxxx to mlp
xtal_test_results[grepl("mlp",xtal_test_results$model),]$model='mlp'

The Figure below show the error distribution for the TD-MLR and MLP for the three sensors on Xtal

bwplot(value~model|sensor,data=xtal_test_results %>% filter( model=="mlp"|model=="mlr") ,groups=factor(model),layout=c(3,1),scales=list(relation='free'))

px1=wilcox.test(value~model,data=xtal_test_results %>% filter(sensor=="X1" & (model=="mlp"|model=="mlr")))
px2=wilcox.test(value~model,data=xtal_test_results %>% filter(sensor=="X2" & (model=="mlp"|model=="mlr")))
px3=wilcox.test(value~model,data=xtal_test_results %>% filter(sensor=="X6" & (model=="mlp"|model=="mlr")))

A man-Whitney-Wilcox test is conducted for the three sensors and results are presented in the following table As can be observed the differences are NOT statistically significant for the three sensors

Sensor Name X1 X2 X6
P Value 0.8972211 0.7196787 0.7630755

4.3.2.5 Discussion

After evaluating MLR and MLP, it is possible to conclude that both model performances are comparable in terms of RMSE for the four IMUs evaluated. MLR has shown an equivalent performance in most of the the cases yet outperforming in the case of sensor X6 (GyroZ) for the Crossbow(xbow).

4.3.3 Comparing MLR vs MA

In this section the proposed MLR model is compared with the standar MA method for filtering noise.

4.3.3.1 Evaluation on Gladiator

glad_test_results=test_results %>% select(model,imuid,sample,X1,X2,X6) %>% filter(imuid==5 & (model=='ma' | model=='mlr'))
glad_test_results=melt(glad_test_results,id.vars=c("imuid","model","sample"),variable.name="sensor",value.name="value")

The Figure below show the error distribution for the TD-MLR and MA for the three sensors on Gladiator

bwplot(value~model|sensor,data=glad_test_results %>% filter((model=="ma"|model=="mlr")),groups=factor(model),layout=c(3,1),scales=list(relation='free'))

px1=wilcox.test(value~model,data=glad_test_results %>% filter(sensor=="X1" & (model=="ma"|model=="mlr")))
px2=wilcox.test(value~model,data=glad_test_results %>% filter(sensor=="X2" & (model=="ma"|model=="mlr")))
px6=wilcox.test(value~model,data=glad_test_results %>% filter(sensor=="X6" & (model=="ma"|model=="mlr")))

A man-Whitney-Wilcox test is conducted for the three sensors and results are presented in the following table As can be observed the differences are NOT statistically significant for the two out of the three sensors

Sensor Name X1 X2 X6
P Value 0.0790219 1.860340410^{-23} 0.5885823

4.3.3.2 Evaluation on Xbow

xbow_test_results=test_results %>% select(model,imuid,sample,X1,X2,X6) %>% filter(imuid==6 & (model=='ma' | model=='mlr' | model=='mlp80'| model=='mlp5'))
xbow_test_results=melt(xbow_test_results,id.vars=c("imuid","model","sample"),variable.name="sensor",value.name="value")
xbow_test_results=xbow_test_results %>% filter(model=='mlr' | model=='ma' | ( model=='mlp80' & (sensor=='X1' | sensor=='X6')) | ( model=='mlp5' & sensor=='X2'))

The Figure below show the error distribution for the TD-MLR and MA for the three sensors on Xbow

bwplot(value~model|sensor,data=xbow_test_results %>% filter( model=="ma"|model=="mlr") ,groups=factor(model),layout=c(3,1),scales=list(relation='free'))

px1=wilcox.test(value~model,data=xbow_test_results %>% filter(sensor=="X1" & (model=="ma"|model=="mlr")))
px2=wilcox.test(value~model,data=xbow_test_results %>% filter(sensor=="X2" & (model=="ma"|model=="mlr")))
px6=wilcox.test(value~model,data=xbow_test_results %>% filter(sensor=="X6" & (model=="ma"|model=="mlr")))

A man-Whitney-Wilcox test is conducted for the three sensors and results are presented in the following table As can be observed the differences are statistically significant for the three sensors

Sensor Name X1 X2 X6
P Value 4.355415510^{-5} 1.310551310^{-7} 5.303148110^{-14}

4.3.3.3 Evaluation on XSENS

xsens_test_results=test_results %>% select(model,imuid,sample,X1,X2,X6) %>% filter(imuid==8 & (model=='ma' | model=='mlr' | model=='mlp80'| model=='mlp100'))
xsens_test_results=melt(xsens_test_results,id.vars=c("imuid","model","sample"),variable.name="sensor",value.name="value")
xsens_test_results=xsens_test_results %>% filter(model=='mlr' | model=='ma' | ( model=='mlp100' & (sensor=='X2' | sensor=='X6')) | ( model=='mlp80' & sensor=='X1'))

The Figure below show the error distribution for the TD-MLR and MA for the three sensors on Xsens

bwplot(value~model|sensor,data=xsens_test_results %>% filter( model=="ma"|model=="mlr") ,groups=factor(model),layout=c(3,1),scales=list(relation='free'))

px1=wilcox.test(value~model,data=xsens_test_results %>% filter(sensor=="X1" & (model=="ma"|model=="mlr")))
px2=wilcox.test(value~model,data=xsens_test_results %>% filter(sensor=="X2" & (model=="ma"|model=="mlr")))
px6=wilcox.test(value~model,data=xsens_test_results %>% filter(sensor=="X6" & (model=="ma"|model=="mlr")))

A man-Whitney-Wilcox test is conducted for the three sensors and results are presented in the following table As can be observed the differences are statistically significant for the three sensors

Sensor Name X1 X2 X6
P Value 0.0057428 2.556541410^{-6} 1.860340410^{-23}

4.3.3.4 Evaluation on Xtal

xtal_test_results=test_results %>% select(model,imuid,sample,X1,X2,X6) %>% filter(imuid==9 & (model=='ma' | model=='mlr' | model=='mlp80'| model=='mlp100'| model=='mlp5'))
xtal_test_results=melt(xtal_test_results,id.vars=c("imuid","model","sample"),variable.name="sensor",value.name="value")
xtal_test_results=xtal_test_results %>% filter(model=='mlr' | model=='ma' | ( model=='mlp80' & sensor=='X1') | ( model=='mlp100' & sensor=='X2')| ( model=='mlp100' & sensor=='X6'))

The Figure below show the error distribution for the TD-MLR and MA for the three sensors on Xtal

bwplot(value~model|sensor,data=xtal_test_results %>% filter( model=="ma"|model=="mlr") ,groups=factor(model),layout=c(3,1),scales=list(relation='free'))

px1=wilcox.test(value~model,data=xtal_test_results %>% filter(sensor=="X1" & (model=="ma"|model=="mlr")))
px2=wilcox.test(value~model,data=xtal_test_results %>% filter(sensor=="X2" & (model=="ma"|model=="mlr")))
px6=wilcox.test(value~model,data=xtal_test_results %>% filter(sensor=="X6" & (model=="ma"|model=="mlr")))

A man-Whitney-Wilcox test is conducted for the three sensors and results are presented in the following table As can be observed the differences are statistically significant for the three sensors

Sensor Name X1 X2 X6
P Value 8.51802710^{-5} 1.055711710^{-6} 1.860340410^{-23}

4.3.3.5 Discussion

In general the TD-MLR outperforms in three out of four IMUS. Namely, Xbow, Xsens and Xtal, while in the case of Gladiator, the performance of TD-MLR has not show a statistically significant difference. Clearly, given the simplicity of the MA method, it is fair to say that for the Gladiator IMU the application of a TD-MLR model is not justified. However, the MLR has shown a significant difference in the remaining IMUs.

5 Concluding Remarks