# A New Fuzzy Logic-Based Similarity Measure Applied to Large Gap Imputation for Uncorrelated Multivariate Time Series

Authors: Thi-Thu-Hong Phan; André Bigand; Émilie Poisson Caillault
Journal: Applied Computational Intelligence and Soft Computing
 (2018)

---

## Abstract

The completion of missing values is a prevalent problem in many domains of pattern recognition and signal processing. Analyzing data with incompleteness may lead to a loss of power and unreliable results, especially for large missing subsequence(s). Therefore, this paper aims to introduce a new approach for filling successive missing values in low/uncorrelated multivariate time series which allows managing a high level of uncertainty. In this way, we propose using a novel fuzzy weighting-based similarity measure. The proposed method involves three main steps. Firstly, for each incomplete signal, the data before a gap and the data after this gap are considered as two separated reference time series with their respective query windowsQ
b and Q
a. We then find the most similar subsequence (Q
b
s) to the subsequence before this gap Q
b and the most similar one (Q
a
s) to the subsequence after the gap Q
a. To find these similar windows, we build a new similarity measure based on fuzzy grades of basic similarity measures and on fuzzy logic rules. Finally, we fill in the gap with average values of the window following Q
b
s and the one preceding Q
a
s. The experimental results have demonstrated that the proposed approach outperforms the state-of-the-art methods in case of multivariate time series having low/noncorrelated data but effective information on each signal.


---

## Body

## 1. Introduction

Nowadays huge time series can be considered due to the availability of effective low-cost sensors, the wide deployment of remote sensing systems, internet based measure networks, etc. However, collected data are often incomplete for various reasons such as sensor errors, transmission problems, incorrect measurements, bad weather conditions (outdoor sensors), for manualmaintenance, etc. This is particularly the case for marine samples [1] that we consider in this paper. For example, the MAREL-Carnot database characterizes sea water in the eastern English Channel, in France [2]. The data contain nineteen time series that are measured by sensors every 20 minutes as nitrate, fluorescence, phosphate, pH, and so on. The analysis of these data with remarkable size and shape allows sea biologists to reveal events such as algal blooms, understand phytoplankton processes [3] in detail, or detect sea pollution and so on. But the data have a lot of missing values: 62.2% for phosphate, 59.9% for nitrate, 27.22% for pH, etc., and the size of missing data varies from one-third hour to several months.Most proposed models for multivariate time series analysis often have difficulties processing incomplete datasets, despite their powerful techniques. They usually require complete data. Then the question is how can missing values be dealt with? Ignoring or deleting is a simple way to solve this drawback. But serious problems regularly arise when applying this solution. This is prominent in time series data where the considered values depend on the previous ones. Furthermore, an analysis based on the systematic differences between observed and unobserved data leads to biased and unreliable results [4]. Thus, it is important to propose a new technique to estimate the missing values. The imputation technique is a conventional method to handle incompleteness problems [5].Considering imputation methods for multivariate time series, taking advantage of the correlations between variables is commonly applied to predict lacking data [6–11]. This means that relations permit using the values of available features to estimate the missing values of other features. However, considering multivariate datasets having low/noncorrelations (for instance the MAREL-Carnot dataset), the observed values of full variables cannot be utilized to complete attributes containing missing values. To handle missing data in this case, we must employ the observed values of the unique variable with the missing data to compute the incomplete values. Therefore the proposed method has to manage the high level of uncertainty of this kind of signal.Particularly, imperfect time series can be modelled using fuzzy sets. The fuzzy approach makes it possible to handle incomplete data, vague, and imprecise circumstances [12], which provide a high uncertainty environment to make decision. This property enables modelling and short-term forecasting of traffic flow in urban arterial networks using multivariate traffic data [13, 14]. Recent works to urban traffic flow prediction [15] and to lane-changes prediction [16] have been proposed with success. Furthermore, the successful use of fuzzy-based similarity measure in pattern recognition [17], in retrieval systems [12], and in recommendation systems [18] leads us to study its ability to complete missing values in uncorrelated multivariate time series. Wang et al. [19] proposed using information granules and fuzzy clustering for time series long-term forecasting. But according to our knowledge, there is no application devoted to complete large gap(s) in uncorrelated multivariate time series using a fuzzy-weighted similarity measure.Thus, this paper aims to propose a new approach, named FSMUMI, to fill large missing values in low/uncorrelated multivariate time series by developing a new similarity measure based on fuzzy logic. However, estimating the distribution of missing values and whole signals is very difficult, so our approach makes an assumption of effective patterns (or recurrent data) on each signal.The rest of this paper is organized as follows. In Section2, related works to imputation methods and fuzzy similarity measure are reviewed. Section 3 introduces our approach for completing large missing subsequences in low/uncorrelated multivariate time series. Next, Section 4 demonstrates our experimental protocol for the imputation task. Section 5 presents results and discussion. Conclusions are drawn and future work is presented in the last section.


## 2. Related Works

This section presents, first, related work about multivariate imputation methods, followed by a review on the fuzzy similarity measure and its applications.


### 2.1. Classical Multivariate Imputation Methods

Up to now, numerous successful researches have been devoted to complete missing data in multivariate time series imputation such as [10, 11, 20–28]. Imputation techniques can be categorized in different perspectives: model-based or machine learning-based and clustering-based imputation techniques.In view of the model-based imputation, two main methods were proposed. The first method was introduced by Schafer [20]. With the hypothesis that all variables follow a multivariate normal distribution, this approach is based on the multivariate normal (MVN) model to determine completion values. And, the second method, namely, MICE, was developed by van Buuren et al. [21] and Raghunathan et al. [22]. This method uses chained equations to fill in incomplete data: for each variable with missing values, MICE computes the imputation data by exploiting the relations between all other variables.According to the concept of machine learning-based imputation, many studies focus on completion of missing data in multivariate time series. Stekhoven and Bühlmann [6] implemented missForest based on the Random Forest (RF) method for multivariate imputation. P.Bonissone et al. [29] proposed a fuzzy version of RF that they named fuzzy random forest FRF. At the moment FRF is only devoted to classification and in our case FRF may be only interesting to separate correlated and uncorrelated variables in multivariate time series if necessary. In [25], Shah et al. investigated a variant of MICE which fills in each variable using the estimation generated from RF. The results showed that the combination of MICE and RF was more efficient than original methods for multivariate imputation. K-Nearest Neighbors (k-NN)-based imputation is also a popular method for completing missing values such as [11, 26, 27, 30–32]. This approach identifies k most similar patterns in the space of available features to impute missing data.Besides these principal techniques, clustering-based imputation approaches are considered as power tools for completing missing values thanks to their ability to detect similar patterns. The objective of these techniques is to separate the data into several clusters when satisfying the following conditions: maximizing the intercluster similarity and minimizing intracluster dissimilarity. Li et al. [33] proposed the k-means clustering imputation technique that estimates missing values using the final cluster information. The fuzzy c-means (FcM) clustering is a common extension of k-means. The squared-norm is applied to measure the similarity between cluster centers and data points. Different applications based on FcM are investigated for the imputation task as [7–9, 34–38]. Wang et al. [19] used FcM based on DTW to successfully predict time series in long-term forecasting.In general, most of the imputation algorithms for multivariate time series take advantage of dependencies between attributes to predict missing values.


### 2.2. Methods Based on Fuzzy Similarity Measure

Indeed similarity-based approaches are a promising tool for time series analysis. However, many of these techniques rely on parameter tuning, and they may have shortcomings due to dependencies between variables. The objective of this study is to fill large missing values inuncorrelated multivariate time series. Thus, we have to deal with a high level of uncertainty. Mikalsen et al. [39] proposed using GMM (Gaussian mixture models) and cluster kernel to deal with uncertainty. Their method needs ensemble learning with numerous learning datasets that are not available in our case at the moment (marine data). So we have chosen to model this global uncertainty using fuzzy sets (FS) introduced by Zadeh [40]. These techniques consider that measurements have inherent vagueness rather than randomness.Uncertainty is classically presented using three conceptually distinctive characteristics: fuzziness, randomness and incompleteness. This classification is interesting for many applications, like sensor management (image processing, speech processing, and time series processing) and practical decision-making. This paper focuses on (sensor) measurements treatment but is also relevant for other applications.Incompleteness often affects time series prediction (time series obtained from marine data such as salinity and temperature). So it seems natural to use fuzzy similarity between subsequences of time series to deal with these three kinds of uncertainties (fuzziness, randomness, and incompleteness). Fuzzy sets are now well known and we only need to remind the basic definition of “FS.” Considering the universe X, a fuzzy set A
∈
X is characterized using a fuzzy membership function μ
A:(1)
μ
A
:
X
→
0,1
,whereμ
A
(
x
) represents the membership of x to A and is associated to the uncertainty of x. In our case, we will consider similarity values between the subsequences as defined in the following. One solution to deal with uncertainty brought by multivariate time series is to use the concept of fuzzy time series [41]. In this framework, the variable observations are considered as fuzzy numbers instead of real numbers. In our case the same modelling is used considering distance measures between subsequences and then we compute the fuzzy similarity between these subsequences to find similar windows in order to estimate the missing values in observations.Fuzzy similarity is a generalization of the classical concept of equivalence and defines the resemblance between two objects (here subsequences of time series). Similarity measures of fuzzy values have been compared in [42] and have been extended in [43]. In [42], Pappis and Karacapilidis presented three main kinds of similarity measures of fuzzy values, including(i)
measures based on the operations of union and intersection,
(ii)
measures based on the maximum difference,
(iii)
measures based on the difference and the sum of membership grades.In [44, 45], the authors used these definitions to propose a distance metric for a space of linguistic summaries based on fuzzy protoforms. Almeida et al. extended this work to put forward linguistic summaries of categorical time series [46]. The introduced similarity measure takes into account not only the linguistic meaning of the summaries but also the numerical characteristic attached to them. In the same way, Gupta et al. [12] introduced this approach to create a hybrid similarity measure based on fuzzy logic. The approach is used to retrieve relevant documents. In the other research, Al-shamri and Al-Ashwal presented fuzzy weightings of popular similarity measures for memory-based collaborative recommend systems [18].Concerning the similarity between two subsequences of time series, we can use the DTW cost as a similarity measure. However, to deal with the high level of uncertainty of the processed signals, numerous similarity measures can be used to compute similarity like the cosine similarity, Euclidean distance, Pearson correlation coefficient. Moreover, a fuzzy-weighted combination of scores generated from different similarity measures could comparatively achieve better retrieval results than the use of a single similarity measure [12, 18].Based on the same concepts, we propose using a fuzzy rules interpolation scheme between grades of membership of fuzzy values. This method makes it possible to build a new hybrid similarity measure for finding similar values between subsequences of time series.


## 2.1. Classical Multivariate Imputation Methods

Up to now, numerous successful researches have been devoted to complete missing data in multivariate time series imputation such as [10, 11, 20–28]. Imputation techniques can be categorized in different perspectives: model-based or machine learning-based and clustering-based imputation techniques.In view of the model-based imputation, two main methods were proposed. The first method was introduced by Schafer [20]. With the hypothesis that all variables follow a multivariate normal distribution, this approach is based on the multivariate normal (MVN) model to determine completion values. And, the second method, namely, MICE, was developed by van Buuren et al. [21] and Raghunathan et al. [22]. This method uses chained equations to fill in incomplete data: for each variable with missing values, MICE computes the imputation data by exploiting the relations between all other variables.According to the concept of machine learning-based imputation, many studies focus on completion of missing data in multivariate time series. Stekhoven and Bühlmann [6] implemented missForest based on the Random Forest (RF) method for multivariate imputation. P.Bonissone et al. [29] proposed a fuzzy version of RF that they named fuzzy random forest FRF. At the moment FRF is only devoted to classification and in our case FRF may be only interesting to separate correlated and uncorrelated variables in multivariate time series if necessary. In [25], Shah et al. investigated a variant of MICE which fills in each variable using the estimation generated from RF. The results showed that the combination of MICE and RF was more efficient than original methods for multivariate imputation. K-Nearest Neighbors (k-NN)-based imputation is also a popular method for completing missing values such as [11, 26, 27, 30–32]. This approach identifies k most similar patterns in the space of available features to impute missing data.Besides these principal techniques, clustering-based imputation approaches are considered as power tools for completing missing values thanks to their ability to detect similar patterns. The objective of these techniques is to separate the data into several clusters when satisfying the following conditions: maximizing the intercluster similarity and minimizing intracluster dissimilarity. Li et al. [33] proposed the k-means clustering imputation technique that estimates missing values using the final cluster information. The fuzzy c-means (FcM) clustering is a common extension of k-means. The squared-norm is applied to measure the similarity between cluster centers and data points. Different applications based on FcM are investigated for the imputation task as [7–9, 34–38]. Wang et al. [19] used FcM based on DTW to successfully predict time series in long-term forecasting.In general, most of the imputation algorithms for multivariate time series take advantage of dependencies between attributes to predict missing values.


## 2.2. Methods Based on Fuzzy Similarity Measure

Indeed similarity-based approaches are a promising tool for time series analysis. However, many of these techniques rely on parameter tuning, and they may have shortcomings due to dependencies between variables. The objective of this study is to fill large missing values inuncorrelated multivariate time series. Thus, we have to deal with a high level of uncertainty. Mikalsen et al. [39] proposed using GMM (Gaussian mixture models) and cluster kernel to deal with uncertainty. Their method needs ensemble learning with numerous learning datasets that are not available in our case at the moment (marine data). So we have chosen to model this global uncertainty using fuzzy sets (FS) introduced by Zadeh [40]. These techniques consider that measurements have inherent vagueness rather than randomness.Uncertainty is classically presented using three conceptually distinctive characteristics: fuzziness, randomness and incompleteness. This classification is interesting for many applications, like sensor management (image processing, speech processing, and time series processing) and practical decision-making. This paper focuses on (sensor) measurements treatment but is also relevant for other applications.Incompleteness often affects time series prediction (time series obtained from marine data such as salinity and temperature). So it seems natural to use fuzzy similarity between subsequences of time series to deal with these three kinds of uncertainties (fuzziness, randomness, and incompleteness). Fuzzy sets are now well known and we only need to remind the basic definition of “FS.” Considering the universe X, a fuzzy set A
∈
X is characterized using a fuzzy membership function μ
A:(1)
μ
A
:
X
→
0,1
,whereμ
A
(
x
) represents the membership of x to A and is associated to the uncertainty of x. In our case, we will consider similarity values between the subsequences as defined in the following. One solution to deal with uncertainty brought by multivariate time series is to use the concept of fuzzy time series [41]. In this framework, the variable observations are considered as fuzzy numbers instead of real numbers. In our case the same modelling is used considering distance measures between subsequences and then we compute the fuzzy similarity between these subsequences to find similar windows in order to estimate the missing values in observations.Fuzzy similarity is a generalization of the classical concept of equivalence and defines the resemblance between two objects (here subsequences of time series). Similarity measures of fuzzy values have been compared in [42] and have been extended in [43]. In [42], Pappis and Karacapilidis presented three main kinds of similarity measures of fuzzy values, including(i)
measures based on the operations of union and intersection,
(ii)
measures based on the maximum difference,
(iii)
measures based on the difference and the sum of membership grades.In [44, 45], the authors used these definitions to propose a distance metric for a space of linguistic summaries based on fuzzy protoforms. Almeida et al. extended this work to put forward linguistic summaries of categorical time series [46]. The introduced similarity measure takes into account not only the linguistic meaning of the summaries but also the numerical characteristic attached to them. In the same way, Gupta et al. [12] introduced this approach to create a hybrid similarity measure based on fuzzy logic. The approach is used to retrieve relevant documents. In the other research, Al-shamri and Al-Ashwal presented fuzzy weightings of popular similarity measures for memory-based collaborative recommend systems [18].Concerning the similarity between two subsequences of time series, we can use the DTW cost as a similarity measure. However, to deal with the high level of uncertainty of the processed signals, numerous similarity measures can be used to compute similarity like the cosine similarity, Euclidean distance, Pearson correlation coefficient. Moreover, a fuzzy-weighted combination of scores generated from different similarity measures could comparatively achieve better retrieval results than the use of a single similarity measure [12, 18].Based on the same concepts, we propose using a fuzzy rules interpolation scheme between grades of membership of fuzzy values. This method makes it possible to build a new hybrid similarity measure for finding similar values between subsequences of time series.


## 3. Proposed Approach

The proposed imputation method is based on the retrieval and the similarity comparison of available subsequences. In order to compare the subsequences, we create a new similarity measure applying a multiple fuzzy rules interpolation. This section is divided into two parts. Firstly, we focus on the way to compute a new similarity measure between subsequences. Then, we provide details of the proposed approach (namely, Fuzzy Similarity Measure Based Uncorrelated Multivariate Imputation, FSMUMI) to impute the successive missing values of low/uncorrelated multivariate time series.


### 3.1. Fuzzy-Weighted Similarity Measure between Subsequences

To introduce a new similarity measure using multiple fuzzy rules interpolation to solve the missing problem, we have to define an information granule, as introduced by Pedrycz [47]. The principle of justifiable granularity of experimental data is based on two conditions: (i) the numeric evidence accumulated within the bounds of numeric data has to be as high as possible and, (ii) at the same time, the information granule should be as specific as possible [19].To answer the first condition, we take into account 3 different distance measures between two subsequencesQ (Q
=
q
i
,
i
=
1
,
…
,
T) and R (R
=
r
i
,
i
=
1
,
…
,
T) including Cosine distance, Euclidean distance (these two measures are widely used in the literature), and Similarity distance (this one was presented in our previous study [48]). These three measures are defined as follows:(i)
Cosine distance is computed by (2). This coefficient presents the cosine of the angle between Q and R
(2)
C
o
s
i
n
e
Q
,
R
=
∑
i
=
1
T
q
i
·
r
i
∑
i
=
1
T
q
i
2
·
∑
i
=
1
T
r
i
2
(ii)
Euclidean distance is calculated by(3)
E
D
∗
Q
,
R
=
∑
i
=
1
T
q
i
-
r
i
2
To satisfy the input condition of fuzzy logic rules, we normalize this distance to[
0,1
] by this function E
D
=
1
/
(
1
+
E
D
∗
(
q
,
r
)
).
(iii)
Similarity measure is defined by the function (4). This measure indicates the similarity percentage between Q and R
(4)
S
i
m
Q
,
R
=
1
T
∑
i
=
1
T
1
1
+
q
i
-
r
i
/
max
⁡
Q
-
min
⁡
QTo answer the second condition, we use these 3 distance measures (or attributes) to generate 4 fuzzy similarities (see Figure2), then applied to a fuzzy inference system (see Figure 1) using the cylindrical extension of the 3 attributes which provides 3 coefficients to calculate a new similarity measure. The universe of discourse of each distance measure is normalized to the value 1.Figure 1
Computing scheme of the new similarity measure.Figure 2
Membership function of fuzzy similarity values.And, finally, the new similarity measure is determined by(5)
F
B
S
M
=
w
1
∗
C
o
s
i
n
e
Q
,
R
+
w
2
∗
E
D
Q
,
R
+
w
3
∗
S
i
m
Q
,
Rwhere w
1, w
2, and w
3 are the weights of the Cosine, ED, and Sim measures, respectively. Thus uncertainty modelled using FS is kept during the similarity computation and makes it possible to deal with a high level of uncertainty as shown in the sequel. The coefficients w
i are generated from the fuzzy interpolation system (Figure 1). We use FuzzyR R-package [49] to develop this system. All input and output variables are expressed by 4 linguistic terms as low, medium, medium-high, and high. A trapezoidal membership function is handled in this case to match input and output spaces to a degree of membership (Figure 2). The multiple rules interpolation is applied to create the fuzzy rules base. So, 64 fuzzy rules are introduced. Each fuzzy rule is presented in the following form:Rule R:IF (C
o
s
i
n
e is l
v
1) and (E
D is l
v
2) and (S
i
m is l
v
3)THEN (w
1 is l
w
1) and (w
2 is l
w
2) and (w
3 is l
w
3) in which l
v
i
,
l
w
i
∈
low
,
medium
,
medium-high
,
high, and i
=
1,2
,
3.


### 3.2. FSMUBI Approach

Let us consider some notations about multivariate time series and the concept of large gap. A multivariate time series is represented as a matrixX
N
×
M with M collected signals ofsize
N. x
(
t
,
i
) is the value of the ith signal at time t. x
t
=
x
t
,
i
,
i
=
1
,
…
,
M is the feature vector at the t-th observation of all variables. X is called an incomplete time series when it contains missing values. We define the term gap of T-size at position t as a portion of X where at least one signal of X between t and t
+
T
-
1 contains consecutive missing values ∃
i
∣
∀
t
∈
t
,
t
+
T
-
1
,
x
t
,
i
=
N
A.Here, we deal with large missing values in low/uncorrelated multivariate time series. For isolated missing values (T
=
1) or small T-gap, conventional techniques can be applied such as the mean or the median of available values [50, 51]. A T-gap is large when the duration T is longer than known change process. For instance, in phytoplankton study, T is equal to one hour to characterize Langmuir cells and one day for algal bloom processes [52]. For small time series (N
<
10,000) without prior knowledge of an application and its change process, we set a large gap when T
≥
5
%
N.The mechanism of FSMUMI approach is demonstrated in Figure3. Without loss of generality, in this figure, we consider a multivariate time series including 3 variables whose correlations are low. The proposed approach involves three major stages. The first stage is to build two queries Q
a and Q
b. The second stage is devoted to find the most similar windows to the queries. This stage includes two minor steps, comparing sliding windows to queries by using the new similarity measure and selecting the similar windows Q
b
s and Q
a
s. Finally, the imputation values are computed by averaging values of the window following Q
b
s and the one preceding Q
a
s to complete the gap.Figure 3
Scheme of the completion process:(
1
) building queries, (
2
) comparing sliding windows, (
3
) selecting the most similar windows, and (
4
) completing gap.This method concentrates on filling missing values in low/uncorrelated multivariate time series. For this type of data, we cannot take advantage of the relations between features to estimate missing values. So we must base our approach on observed values on each signal to complete missing data on itself. This means that we can complete missing data on each variable, one by one. Further, an important point of our approach is that each incomplete signal is processed as two separated time series, one time series before the considered gap and one time series after this gap. This allows increasing the search space for similar values. Moreover, applying the proposed process (one by one), FSMUMI makes it possible to handle the problem of wholly missing variables (missing data at the same index in the all variables).The proposed model is described in Algorithm1 and is mainly divided into three phases:(i)
The first phase: Building queries (cf. 1 in Figure2)
For each incomplete signal and eachT-gap, two referenced databases are extracted from the original time series and two query windows are built to retrieve similar windows. The data before the gap (noted D
b) and the data after this gap (denoted D
a) are considered as two separated time series. We noted Q
b is the subsequence before the gap and Q
a is the respective subsequence after the gap. These query windows have the same size T as the gap.
(ii)
The second phase: Finding the most similar windows (cf. 2 and 3 in Figure2)
For theD
b database, we build sliding reference windows (noted R) of size T. From these R windows, we retrieve the most similar window (Q
b
s) to the Q
b query using the new similarity measure f
b
s
m as previously defined in Section 3.1. Details are in the following:
We first find the threshold, which allows considering two windows to be similar. For each increments
t
e
p
_
t
h
r
e
s
h
o
l
d, we compute a f
b
s
m similarity measure between a sliding window R and the query Q
b. The t
h
r
e
s
h
o
l
d is the maximum value obtained from the all f
b
s
m calculated (Step a: in Algorithm 1).
We then find the most similar window to the queryQ
b. For each increment similar window s
t
e
p
_
s
i
m
_
w
i
n, a f
b
s
m of a R sliding reference and the query Q
b is estimated. We then compare this f
b
s
m to the t
h
r
e
s
h
o
l
d to determine if this R reference is similar to the query Q
b. We finally choose the most similar window Q
b
s with the maximum f
b
s
m of all the similar windows (Step b: in Algorithm 1).
The same process is performed to find the most similar windowQ
a
s in D
a data.
In the proposed approach, the dynamics and the shape of data before and after a gap are a key-point of our method. This means we take into account both queriesQ
a (after the gap) and Q
b (before the gap). This makes it possible to find out windows that have the most similar dynamics and shape to the queries.
(iii)
The third phase (cf. 4 in Figure2)
When results from both referenced time series are available, we fill in the gap by averaging values of the window precedingQ
a
s and the one following Q
b
s. The average values are used in our approach because model averaging makes the final results more stable and unbiased [53].Algorithm 1:FSMUMI algorithm.
Input:
X
=
x
1
,
x
2
,
…
,
x
M: incomplete uncorrelated multivariate time series
N: size of time series.
t: index of a gap (position of the first missing of the gap)
T: size of the gap
s
t
e
p
_
t
h
r
e
s
h
o
l
d: increment for finding a threshold
s
t
e
p
_
s
i
m
_
w
i
n: increment for finding a similar window
Output:  Y - completed (imputed) time series
(
1
)
foreach incomplete signal x
j
∈
X
d
o
(
2
)
foreach gap at t index in x
j
d
o
(
3
)
Divide x
j into two separated time series D
a, D
b: D
a
=
x
j
t
+
T
:
N
,
D
b
=
x
j
[
1
:
t
-
1
]
(
4
)
Completing all lines containing missing parameter on D
a
,
D
b by a max trapezoid function
(
5
)
Construct queries Q
a, Q
b-temporal windows after and before the gap Q
a
=
D
a
[
1
:
T
], Q
b
=
D
b
[
t
-
T
+
1
:
t
-
1
]
(
6
)
for
D
b datad
o
(
7
)
Step a: Find the threshold in the D
b database
(
8
)
i
←
1; F
S
M
←
N
U
L
L
(
9
)
while
i
≤
length
(
D
b
)
d
o
(
10
)
k
←
i
+
T
-
1
(
11
)
Create a reference window: R
(
i
)
=
D
b
[
i
:
k
]
(12)Calculate a fuzzy-based similarity measure between Q
b and R
(
i
): f
b
s
m
(13)Save the f
b
s
m to F
M
S
(14)
i
←
i
+
s
t
e
p
_
t
h
r
e
s
h
o
l
d
(15)
end while
(16)
return
t
h
r
e
s
h
o
l
d
=
max
⁡
{
F
B
M
S
}
(17)
Step b: Find similar windows in the D
b database
(18)
i
←
1; L
o
p
b
←
N
U
L
L
(19)
while
i
≤
length
(
D
b
)
d
o
(20)
k
←
i
+
T
-
1
(21)Create a reference window: R
(
i
)
=
D
b
[
i
:
k
]
(22)Calculate a fuzzy-based similarity measure between Q
b and R
(
i
): f
b
s
m
(23)
if
f
b
s
m
≥
t
h
r
e
s
h
o
l
d
then
(24)Save position of R
(
i
) to L
o
p
b
(25)
end if
(26)
i
←
i
+
s
t
e
p
_
s
i
m
_
w
i
n
(27)
end while
(28)
return  position of Q
b
s - the most similar window to Q
b having the maximum fuzzy similarity measure in the
L
o
p
b list.
(29)
end for
(30)
for
D
a datad
o
(31)Perform Step a and Step b for D
a data
(32)
return position of Q
a
s - the most similar window to Q
a
(33)
end for
(34)Replace the missing values at the position t by average vector of the window after Q
b
s and the one previous Q
a
s
(35)
end for
(36)
end for
(37)
return
Y - imputed time series


## 3.1. Fuzzy-Weighted Similarity Measure between Subsequences

To introduce a new similarity measure using multiple fuzzy rules interpolation to solve the missing problem, we have to define an information granule, as introduced by Pedrycz [47]. The principle of justifiable granularity of experimental data is based on two conditions: (i) the numeric evidence accumulated within the bounds of numeric data has to be as high as possible and, (ii) at the same time, the information granule should be as specific as possible [19].To answer the first condition, we take into account 3 different distance measures between two subsequencesQ (Q
=
q
i
,
i
=
1
,
…
,
T) and R (R
=
r
i
,
i
=
1
,
…
,
T) including Cosine distance, Euclidean distance (these two measures are widely used in the literature), and Similarity distance (this one was presented in our previous study [48]). These three measures are defined as follows:(i)
Cosine distance is computed by (2). This coefficient presents the cosine of the angle between Q and R
(2)
C
o
s
i
n
e
Q
,
R
=
∑
i
=
1
T
q
i
·
r
i
∑
i
=
1
T
q
i
2
·
∑
i
=
1
T
r
i
2
(ii)
Euclidean distance is calculated by(3)
E
D
∗
Q
,
R
=
∑
i
=
1
T
q
i
-
r
i
2
To satisfy the input condition of fuzzy logic rules, we normalize this distance to[
0,1
] by this function E
D
=
1
/
(
1
+
E
D
∗
(
q
,
r
)
).
(iii)
Similarity measure is defined by the function (4). This measure indicates the similarity percentage between Q and R
(4)
S
i
m
Q
,
R
=
1
T
∑
i
=
1
T
1
1
+
q
i
-
r
i
/
max
⁡
Q
-
min
⁡
QTo answer the second condition, we use these 3 distance measures (or attributes) to generate 4 fuzzy similarities (see Figure2), then applied to a fuzzy inference system (see Figure 1) using the cylindrical extension of the 3 attributes which provides 3 coefficients to calculate a new similarity measure. The universe of discourse of each distance measure is normalized to the value 1.Figure 1
Computing scheme of the new similarity measure.Figure 2
Membership function of fuzzy similarity values.And, finally, the new similarity measure is determined by(5)
F
B
S
M
=
w
1
∗
C
o
s
i
n
e
Q
,
R
+
w
2
∗
E
D
Q
,
R
+
w
3
∗
S
i
m
Q
,
Rwhere w
1, w
2, and w
3 are the weights of the Cosine, ED, and Sim measures, respectively. Thus uncertainty modelled using FS is kept during the similarity computation and makes it possible to deal with a high level of uncertainty as shown in the sequel. The coefficients w
i are generated from the fuzzy interpolation system (Figure 1). We use FuzzyR R-package [49] to develop this system. All input and output variables are expressed by 4 linguistic terms as low, medium, medium-high, and high. A trapezoidal membership function is handled in this case to match input and output spaces to a degree of membership (Figure 2). The multiple rules interpolation is applied to create the fuzzy rules base. So, 64 fuzzy rules are introduced. Each fuzzy rule is presented in the following form:Rule R:IF (C
o
s
i
n
e is l
v
1) and (E
D is l
v
2) and (S
i
m is l
v
3)THEN (w
1 is l
w
1) and (w
2 is l
w
2) and (w
3 is l
w
3) in which l
v
i
,
l
w
i
∈
low
,
medium
,
medium-high
,
high, and i
=
1,2
,
3.


## 3.2. FSMUBI Approach

Let us consider some notations about multivariate time series and the concept of large gap. A multivariate time series is represented as a matrixX
N
×
M with M collected signals ofsize
N. x
(
t
,
i
) is the value of the ith signal at time t. x
t
=
x
t
,
i
,
i
=
1
,
…
,
M is the feature vector at the t-th observation of all variables. X is called an incomplete time series when it contains missing values. We define the term gap of T-size at position t as a portion of X where at least one signal of X between t and t
+
T
-
1 contains consecutive missing values ∃
i
∣
∀
t
∈
t
,
t
+
T
-
1
,
x
t
,
i
=
N
A.Here, we deal with large missing values in low/uncorrelated multivariate time series. For isolated missing values (T
=
1) or small T-gap, conventional techniques can be applied such as the mean or the median of available values [50, 51]. A T-gap is large when the duration T is longer than known change process. For instance, in phytoplankton study, T is equal to one hour to characterize Langmuir cells and one day for algal bloom processes [52]. For small time series (N
<
10,000) without prior knowledge of an application and its change process, we set a large gap when T
≥
5
%
N.The mechanism of FSMUMI approach is demonstrated in Figure3. Without loss of generality, in this figure, we consider a multivariate time series including 3 variables whose correlations are low. The proposed approach involves three major stages. The first stage is to build two queries Q
a and Q
b. The second stage is devoted to find the most similar windows to the queries. This stage includes two minor steps, comparing sliding windows to queries by using the new similarity measure and selecting the similar windows Q
b
s and Q
a
s. Finally, the imputation values are computed by averaging values of the window following Q
b
s and the one preceding Q
a
s to complete the gap.Figure 3
Scheme of the completion process:(
1
) building queries, (
2
) comparing sliding windows, (
3
) selecting the most similar windows, and (
4
) completing gap.This method concentrates on filling missing values in low/uncorrelated multivariate time series. For this type of data, we cannot take advantage of the relations between features to estimate missing values. So we must base our approach on observed values on each signal to complete missing data on itself. This means that we can complete missing data on each variable, one by one. Further, an important point of our approach is that each incomplete signal is processed as two separated time series, one time series before the considered gap and one time series after this gap. This allows increasing the search space for similar values. Moreover, applying the proposed process (one by one), FSMUMI makes it possible to handle the problem of wholly missing variables (missing data at the same index in the all variables).The proposed model is described in Algorithm1 and is mainly divided into three phases:(i)
The first phase: Building queries (cf. 1 in Figure2)
For each incomplete signal and eachT-gap, two referenced databases are extracted from the original time series and two query windows are built to retrieve similar windows. The data before the gap (noted D
b) and the data after this gap (denoted D
a) are considered as two separated time series. We noted Q
b is the subsequence before the gap and Q
a is the respective subsequence after the gap. These query windows have the same size T as the gap.
(ii)
The second phase: Finding the most similar windows (cf. 2 and 3 in Figure2)
For theD
b database, we build sliding reference windows (noted R) of size T. From these R windows, we retrieve the most similar window (Q
b
s) to the Q
b query using the new similarity measure f
b
s
m as previously defined in Section 3.1. Details are in the following:
We first find the threshold, which allows considering two windows to be similar. For each increments
t
e
p
_
t
h
r
e
s
h
o
l
d, we compute a f
b
s
m similarity measure between a sliding window R and the query Q
b. The t
h
r
e
s
h
o
l
d is the maximum value obtained from the all f
b
s
m calculated (Step a: in Algorithm 1).
We then find the most similar window to the queryQ
b. For each increment similar window s
t
e
p
_
s
i
m
_
w
i
n, a f
b
s
m of a R sliding reference and the query Q
b is estimated. We then compare this f
b
s
m to the t
h
r
e
s
h
o
l
d to determine if this R reference is similar to the query Q
b. We finally choose the most similar window Q
b
s with the maximum f
b
s
m of all the similar windows (Step b: in Algorithm 1).
The same process is performed to find the most similar windowQ
a
s in D
a data.
In the proposed approach, the dynamics and the shape of data before and after a gap are a key-point of our method. This means we take into account both queriesQ
a (after the gap) and Q
b (before the gap). This makes it possible to find out windows that have the most similar dynamics and shape to the queries.
(iii)
The third phase (cf. 4 in Figure2)
When results from both referenced time series are available, we fill in the gap by averaging values of the window precedingQ
a
s and the one following Q
b
s. The average values are used in our approach because model averaging makes the final results more stable and unbiased [53].Algorithm 1:FSMUMI algorithm.
Input:
X
=
x
1
,
x
2
,
…
,
x
M: incomplete uncorrelated multivariate time series
N: size of time series.
t: index of a gap (position of the first missing of the gap)
T: size of the gap
s
t
e
p
_
t
h
r
e
s
h
o
l
d: increment for finding a threshold
s
t
e
p
_
s
i
m
_
w
i
n: increment for finding a similar window
Output:  Y - completed (imputed) time series
(
1
)
foreach incomplete signal x
j
∈
X
d
o
(
2
)
foreach gap at t index in x
j
d
o
(
3
)
Divide x
j into two separated time series D
a, D
b: D
a
=
x
j
t
+
T
:
N
,
D
b
=
x
j
[
1
:
t
-
1
]
(
4
)
Completing all lines containing missing parameter on D
a
,
D
b by a max trapezoid function
(
5
)
Construct queries Q
a, Q
b-temporal windows after and before the gap Q
a
=
D
a
[
1
:
T
], Q
b
=
D
b
[
t
-
T
+
1
:
t
-
1
]
(
6
)
for
D
b datad
o
(
7
)
Step a: Find the threshold in the D
b database
(
8
)
i
←
1; F
S
M
←
N
U
L
L
(
9
)
while
i
≤
length
(
D
b
)
d
o
(
10
)
k
←
i
+
T
-
1
(
11
)
Create a reference window: R
(
i
)
=
D
b
[
i
:
k
]
(12)Calculate a fuzzy-based similarity measure between Q
b and R
(
i
): f
b
s
m
(13)Save the f
b
s
m to F
M
S
(14)
i
←
i
+
s
t
e
p
_
t
h
r
e
s
h
o
l
d
(15)
end while
(16)
return
t
h
r
e
s
h
o
l
d
=
max
⁡
{
F
B
M
S
}
(17)
Step b: Find similar windows in the D
b database
(18)
i
←
1; L
o
p
b
←
N
U
L
L
(19)
while
i
≤
length
(
D
b
)
d
o
(20)
k
←
i
+
T
-
1
(21)Create a reference window: R
(
i
)
=
D
b
[
i
:
k
]
(22)Calculate a fuzzy-based similarity measure between Q
b and R
(
i
): f
b
s
m
(23)
if
f
b
s
m
≥
t
h
r
e
s
h
o
l
d
then
(24)Save position of R
(
i
) to L
o
p
b
(25)
end if
(26)
i
←
i
+
s
t
e
p
_
s
i
m
_
w
i
n
(27)
end while
(28)
return  position of Q
b
s - the most similar window to Q
b having the maximum fuzzy similarity measure in the
L
o
p
b list.
(29)
end for
(30)
for
D
a datad
o
(31)Perform Step a and Step b for D
a data
(32)
return position of Q
a
s - the most similar window to Q
a
(33)
end for
(34)Replace the missing values at the position t by average vector of the window after Q
b
s and the one previous Q
a
s
(35)
end for
(36)
end for
(37)
return
Y - imputed time series


## 4. Experiment Protocol

The experiments are performed on three multivariate time series with the same experiment process and the same gaps, described in detail below.


### 4.1. Datasets Description

For the assessment of the proposed approach and the comparison of its performance to several published algorithms, we use 3 multivariate time series, one from UCI Machine Learning repository, one simulated dataset (this allows us to handle the correlations between variables and percentage of missing values), and finally a real time series hourly sampled by IFREMER (France) in the eastern English Channel.(i)
Synthetic dataset [54]: The data are synthetic time series, including 10 features, 100,000 sampled points. All data points are in the range -0.5 to +0.5. The data appear highly periodic but never exactly repeated. They have structure at different resolutions. Each of the 10 features is generated by independent invocations of the function: (6)
y
=
∑
i
=
3
7
1
2
i
sin
⁡
2
π
2
2
+
i
+
rand
⁡
2
i
t
;
0
≤
t
≤
1
whererand
⁡
(
x
) produces a random integer between 0 and x.
These data are very large so we choose only a subset of 3 signals for performing experiments.
(ii)
Simulated dataset: In the second experiment, a simulated dataset including 3 signals is produced as follows: for the first variable, we use 5 sine functions that have different frequencies and amplitudes F
=
f
1
,
f
2
,
f
3
,
f
4
,
f
5. Next, 3 various noise levels are added to data F, S
=
F
,
F
+
n
o
i
s
e
1
,
F
+
n
o
i
s
e
2
,
F
+
n
o
i
s
e
3. We then repeat S 4 times (this dataset has 32,000 sampled points). In this study, we treat with missing data in low/uncorrelated multivariate time series. So to satisfy this condition, the two remaining signals are generated based on the first signal with the correlations between these signals are low (≤
0.1
%). We apply the Corgen function of ecodist R-package [55] to create the second and the third variables.
(iii)
MAREL-Carnot dataset [2]: The third experiment is conducted on MAREL-Carnot dataset. This dataset consists of nineteen series such as phosphate, salinity, turbidity, water temperature, fluorescence, and water level that characterize sea water. These signals were collected from the 1
s
t January 2005 to the 9
t
h February 2009 at a 20 minute frequency. Here they were hourly sampled, so they have 35,334 time samples. But the data include many missing values, the size of missing data varying on each signal. To assess the performance of the proposed method and compare it with other approaches, we choose a subgroup including fluorescence, water level, and water temperature (the water level and the fluorescence signals are completed data, while water temperature contains isolated missing values and many gaps). We selected these signals because their correlations are low.
After completing missing values, completion data will be compared with the actual values in the completed series to evaluate the ability of different imputation methods. Therefore, it is necessary to fill missing values in the water temperature. To ensure the fairness of all algorithms, filling in the water temperature series is performed by using the na.interp method ([56]).


### 4.2. Multivariate Imputation Approaches

In the present study, we perform a comparison of the proposed algorithm with 7 other approaches (comprising Amelia II, FcM, MI, MICE, missForest, na.approx, and DTWUMI) for the imputation of multivariate time series. We use R language to execute all these algorithms.(1)
Amelia II (Amelia II R-package) [57]: The algorithm uses the familiar expectation-maximization algorithm on multiple bootstrapped samples of the original incomplete data to draw values of the complete data parameters. The algorithm then draws imputed values from each set of bootstrapped parameters, replacing the missing values with the drawn values.
(2)
FcM-Fuzzy
c
-means based imputation: This approach involves 2 steps. The first step is to group the whole data into k clusters using fuzzy-c means technique. A cluster membership for each sample and a cluster center are generated for each feature. The second step is to fill in the incomplete data by using the membership degree and the center centroids [33]. We base on the principles of [33] and use the c-means function [58] to develop this approach.
(3)
MI: Multiple Imputation (MI R-package) [59]: This method uses predictive mean matching to estimate missing values of continuous variables. For each missing value, its imputation value is randomly selected from a set of observed values that are the closest predicted mean to the variable with the missing value.
(4)
MICE: Multivariate Imputation via Chained Equations (MICE R-package) [60]: For each incomplete variable under the assumption of MAR (missing at random), the algorithm performs a completion by full conditional specification of predictive models. The same process is implemented with other variables having missing data.
(5)
missForest (missForest R-package) [6]: This algorithm uses random forest method to complete missing values. For each variable containing missing data, missForest builds a random forest model on the available data. To estimate missing data this model is applied in the variable, repeating the procedure until it meets a stopping condition.
(6)
Linear interpolation: na.approx (zoo R-package) [61]: This method is based on an interpolation function to predict each missing point.
(7)
DTWUMI[62]: For each gap, this approach finds the most similar window to the subsequence after (resp. before) the gap based on the combination of shape-features extraction and Dynamic Time Warping algorithms. Then, the previous (resp. following) window of the most similar one in the incomplete signal is used to complete the gap.


### 4.3. Imputation Performance Measurements

In order to estimate the quantitative performance of imputation approaches, six usual criteria in the literature are used as follows:(1)
Similarity evaluates the similar percent between the estimated values (y) and the respective real values (x). This index is defined by(7)
S
i
m
y
,
x
=
1
T
∑
i
=
1
T
1
1
+
y
i
-
x
i
/
max
⁡
x
-
min
⁡
x
where T is the number of missing values. The similarity tends to 1 when the two curves are identical and tends to 0 when the amplitudes are strongly different.
(2)
R
2 score is determined as the square of correlation coefficient between two variables y and x. This indicator makes it possible to assess the quality of an imputation model. A method presents better performance when its score is higher (R
2
∈
[
0,1
])
(3)
RMSE (Root Mean Square Error) is computed as the average squared difference betweeny and x. This is an appreciate coefficient to measure global ability of a completion method. In general, a lower RMSE highlights a better imputation performance.(8)
R
M
S
E
y
,
x
=
1
T
∑
i
=
1
T
y
i
-
x
i
2
It is now well admitted that good imputation performance does not lead automatically to good estimation performance. It is why other indices like FSD, FA2, and FB (that enable evaluating the shape of the two signals) are used in this study.
(4)
FSD (Fraction of Standard Deviation) is defined as(9)
F
S
D
y
,
x
=
2
∗
S
D
y
-
S
D
x
S
D
y
+
S
D
x
This fraction points out whether a method is acceptable or not. Applying to the imputation task, when FSD value approaches 0, an imputation method is impeccable.
(5)
FB: Fractional Bias: determines the rate of predicted valuesy are overestimated or underestimated relative to observed values x. This indicator is given by (10). An imputation model is considered ideal as its FB equals 0.(10)
F
B
y
,
x
=
2
∗
m
e
a
n
y
-
m
e
a
n
x
m
e
a
n
y
+
m
e
a
n
x
(6)
FA2 defines the percentage of outlier between two variablesy and x. It is described by (11)
F
A
2
y
,
x
=
length
0.5
≤
y
/
x
≤
2
length
x
When FA2 value is close to 1, a model is considered perfect.


### 4.4. Experimental Process

Indeed, evaluating the ability of imputation methods cannot be done because the actual values are lacking. So we must produce artificial missing data on completed time series in order to compare the performance of imputation approaches. We use a technique based on three steps to assess the results detailed in the following:(i)
The first step: Generate simulated missing values by removing data values from full time series.
(ii)
The second step: Apply the imputation methods to fill in missing data.
(iii)
The third step: Evaluate the ability of proposed approach and compare with state-of-the-art methods using different performance indices above-mentioned.In this paper, we perform experiments with seven missing data levels on three large datasets. On each signal, we create simulated gaps with different rates ranging from 1%, 2%, 3%, 4%, 5%, 7.5%, to 10% of the data in the complete signal (here the biggest gap of MAREL-Carnot data is 3,533 missing values corresponding to 5 months of hourly sampled). For every missing ratio, the approaches are run 5 times by randomly choosing the positions of missing in the data. We then perform35 iterations for each dataset.


## 4.1. Datasets Description

For the assessment of the proposed approach and the comparison of its performance to several published algorithms, we use 3 multivariate time series, one from UCI Machine Learning repository, one simulated dataset (this allows us to handle the correlations between variables and percentage of missing values), and finally a real time series hourly sampled by IFREMER (France) in the eastern English Channel.(i)
Synthetic dataset [54]: The data are synthetic time series, including 10 features, 100,000 sampled points. All data points are in the range -0.5 to +0.5. The data appear highly periodic but never exactly repeated. They have structure at different resolutions. Each of the 10 features is generated by independent invocations of the function: (6)
y
=
∑
i
=
3
7
1
2
i
sin
⁡
2
π
2
2
+
i
+
rand
⁡
2
i
t
;
0
≤
t
≤
1
whererand
⁡
(
x
) produces a random integer between 0 and x.
These data are very large so we choose only a subset of 3 signals for performing experiments.
(ii)
Simulated dataset: In the second experiment, a simulated dataset including 3 signals is produced as follows: for the first variable, we use 5 sine functions that have different frequencies and amplitudes F
=
f
1
,
f
2
,
f
3
,
f
4
,
f
5. Next, 3 various noise levels are added to data F, S
=
F
,
F
+
n
o
i
s
e
1
,
F
+
n
o
i
s
e
2
,
F
+
n
o
i
s
e
3. We then repeat S 4 times (this dataset has 32,000 sampled points). In this study, we treat with missing data in low/uncorrelated multivariate time series. So to satisfy this condition, the two remaining signals are generated based on the first signal with the correlations between these signals are low (≤
0.1
%). We apply the Corgen function of ecodist R-package [55] to create the second and the third variables.
(iii)
MAREL-Carnot dataset [2]: The third experiment is conducted on MAREL-Carnot dataset. This dataset consists of nineteen series such as phosphate, salinity, turbidity, water temperature, fluorescence, and water level that characterize sea water. These signals were collected from the 1
s
t January 2005 to the 9
t
h February 2009 at a 20 minute frequency. Here they were hourly sampled, so they have 35,334 time samples. But the data include many missing values, the size of missing data varying on each signal. To assess the performance of the proposed method and compare it with other approaches, we choose a subgroup including fluorescence, water level, and water temperature (the water level and the fluorescence signals are completed data, while water temperature contains isolated missing values and many gaps). We selected these signals because their correlations are low.
After completing missing values, completion data will be compared with the actual values in the completed series to evaluate the ability of different imputation methods. Therefore, it is necessary to fill missing values in the water temperature. To ensure the fairness of all algorithms, filling in the water temperature series is performed by using the na.interp method ([56]).


## 4.2. Multivariate Imputation Approaches

In the present study, we perform a comparison of the proposed algorithm with 7 other approaches (comprising Amelia II, FcM, MI, MICE, missForest, na.approx, and DTWUMI) for the imputation of multivariate time series. We use R language to execute all these algorithms.(1)
Amelia II (Amelia II R-package) [57]: The algorithm uses the familiar expectation-maximization algorithm on multiple bootstrapped samples of the original incomplete data to draw values of the complete data parameters. The algorithm then draws imputed values from each set of bootstrapped parameters, replacing the missing values with the drawn values.
(2)
FcM-Fuzzy
c
-means based imputation: This approach involves 2 steps. The first step is to group the whole data into k clusters using fuzzy-c means technique. A cluster membership for each sample and a cluster center are generated for each feature. The second step is to fill in the incomplete data by using the membership degree and the center centroids [33]. We base on the principles of [33] and use the c-means function [58] to develop this approach.
(3)
MI: Multiple Imputation (MI R-package) [59]: This method uses predictive mean matching to estimate missing values of continuous variables. For each missing value, its imputation value is randomly selected from a set of observed values that are the closest predicted mean to the variable with the missing value.
(4)
MICE: Multivariate Imputation via Chained Equations (MICE R-package) [60]: For each incomplete variable under the assumption of MAR (missing at random), the algorithm performs a completion by full conditional specification of predictive models. The same process is implemented with other variables having missing data.
(5)
missForest (missForest R-package) [6]: This algorithm uses random forest method to complete missing values. For each variable containing missing data, missForest builds a random forest model on the available data. To estimate missing data this model is applied in the variable, repeating the procedure until it meets a stopping condition.
(6)
Linear interpolation: na.approx (zoo R-package) [61]: This method is based on an interpolation function to predict each missing point.
(7)
DTWUMI[62]: For each gap, this approach finds the most similar window to the subsequence after (resp. before) the gap based on the combination of shape-features extraction and Dynamic Time Warping algorithms. Then, the previous (resp. following) window of the most similar one in the incomplete signal is used to complete the gap.


## 4.3. Imputation Performance Measurements

In order to estimate the quantitative performance of imputation approaches, six usual criteria in the literature are used as follows:(1)
Similarity evaluates the similar percent between the estimated values (y) and the respective real values (x). This index is defined by(7)
S
i
m
y
,
x
=
1
T
∑
i
=
1
T
1
1
+
y
i
-
x
i
/
max
⁡
x
-
min
⁡
x
where T is the number of missing values. The similarity tends to 1 when the two curves are identical and tends to 0 when the amplitudes are strongly different.
(2)
R
2 score is determined as the square of correlation coefficient between two variables y and x. This indicator makes it possible to assess the quality of an imputation model. A method presents better performance when its score is higher (R
2
∈
[
0,1
])
(3)
RMSE (Root Mean Square Error) is computed as the average squared difference betweeny and x. This is an appreciate coefficient to measure global ability of a completion method. In general, a lower RMSE highlights a better imputation performance.(8)
R
M
S
E
y
,
x
=
1
T
∑
i
=
1
T
y
i
-
x
i
2
It is now well admitted that good imputation performance does not lead automatically to good estimation performance. It is why other indices like FSD, FA2, and FB (that enable evaluating the shape of the two signals) are used in this study.
(4)
FSD (Fraction of Standard Deviation) is defined as(9)
F
S
D
y
,
x
=
2
∗
S
D
y
-
S
D
x
S
D
y
+
S
D
x
This fraction points out whether a method is acceptable or not. Applying to the imputation task, when FSD value approaches 0, an imputation method is impeccable.
(5)
FB: Fractional Bias: determines the rate of predicted valuesy are overestimated or underestimated relative to observed values x. This indicator is given by (10). An imputation model is considered ideal as its FB equals 0.(10)
F
B
y
,
x
=
2
∗
m
e
a
n
y
-
m
e
a
n
x
m
e
a
n
y
+
m
e
a
n
x
(6)
FA2 defines the percentage of outlier between two variablesy and x. It is described by (11)
F
A
2
y
,
x
=
length
0.5
≤
y
/
x
≤
2
length
x
When FA2 value is close to 1, a model is considered perfect.


## 4.4. Experimental Process

Indeed, evaluating the ability of imputation methods cannot be done because the actual values are lacking. So we must produce artificial missing data on completed time series in order to compare the performance of imputation approaches. We use a technique based on three steps to assess the results detailed in the following:(i)
The first step: Generate simulated missing values by removing data values from full time series.
(ii)
The second step: Apply the imputation methods to fill in missing data.
(iii)
The third step: Evaluate the ability of proposed approach and compare with state-of-the-art methods using different performance indices above-mentioned.In this paper, we perform experiments with seven missing data levels on three large datasets. On each signal, we create simulated gaps with different rates ranging from 1%, 2%, 3%, 4%, 5%, 7.5%, to 10% of the data in the complete signal (here the biggest gap of MAREL-Carnot data is 3,533 missing values corresponding to 5 months of hourly sampled). For every missing ratio, the approaches are run 5 times by randomly choosing the positions of missing in the data. We then perform35 iterations for each dataset.


## 5. Results and Discussion

This section provides experiment results obtained from the proposed approach and compares its ability with the seven published approaches. Results are discussed in three parts, i.e., quantitative performance, visual performance, and execution times.


### 5.1. Quantitative Performance Comparison

Tables1, 2, and 3 illustrate the average ability of various imputation methods for synthetic, simulated, and MAREL-Carnot time series using 6 measurements as previously defined. For each missing level, the best results are highlighted in bold. These results demonstrate the improved performance of FSMUMI to complete missing data in low/uncorrelated multivariate time series.Table 1
Average imputation performance indices of various imputation algorithms on synthetic dataset (100,000 collected points).
Gap size
Method
Accuracy indices
Shape indices
1-Sim
1-R
2
RMSE
FSD
FB
1-FA2
1%
FSMUMI
0.136
0.261
0.051
0.358
3.253
0.364
Amelia
0.275
0.999
0.143
0.409
2.252
0.773
FcM
0.231
0.722
0.096
1.889
2.208
0.996
MI
0.275
0.999
0.142
0.421
2.091
0.773
MICE
0.258
0.944
0.13
0.406
2.452
0.72
missForest
0.248
0.915
0.122
0.389
3.976
0.744
na.approx
0.052
0.066
0.019
0.054
0.29
0.074
DTWUMI
0.257
0.713
0.88
0.725
0.405
0.69
2%
FSMUMI
0.1
0.295
0.046
0.155
0.395
0.337
Amelia
0.259
0.998
0.147
0.275
2.005
0.803
FcM
0.208
0.686
0.104
1.863
2.289
0.987
MI
0.259
0.998
0.147
0.268
2.11
0.81
MICE
0.244
0.968
0.14
0.255
7.616
0.759
missForest
0.239
0.968
0.133
0.279
3.156
0.792
na.approx
0.104
0.278
0.047
0.224
0.398
0.347
DTWUMI
0.237
0.775
0.867
0.509
8.449
0.646
3%
FSMUMI
0.113
0.341
0.056
0.219
0.852
0.322
Amelia
0.218
0.911
0.127
0.133
6.128
0.76
FcM
0.214
0.601
0.1
1.832
1.759
0.989
MI
0.253
0.993
0.141
0.236
2.295
0.775
MICE
0.21
0.873
0.118
0.208
5.118
0.703
missForest
0.188
0.796
0.102
0.215
1.846
0.627
na.approx
0.148
0.43
0.072
0.372
2.382
0.577
DTWUMI
0.231
0.799
0.874
0.332
27.952
0.69
4%
FSMUMI
0.06
0.146
0.037
0.099
0.738
0.299
Amelia
0.208
1
0.14
0.213
2.171
0.807
FcM
0.155
0.759
0.095
1.85
2.09
0.986
MI
0.208
0.999
0.14
0.196
2.302
0.807
MICE
0.209
0.987
0.138
0.22
3.748
0.801
missForest
0.196
0.968
0.127
0.216
3.94
0.827
na.approx
0.145
0.721
0.092
0.252
5.251
0.689
DTWUMI
0.148
0.586
0.918
0.185
12.688
0.719
5%
FSMUMI
0.055
0.132
0.032
0.058
0.098
0.201
Amelia
0.214
0.997
0.15
0.147
2.238
0.79
FcM
0.179
0.715
0.108
1.818
2.194
0.993
MI
0.231
0.996
0.167
0.206
3.094
0.808
MICE
0.221
0.968
0.152
0.222
2.3
0.79
missForest
0.212
0.944
0.143
0.315
4.547
0.819
na.approx
0.16
0.8
0.118
0.352
18.217
0.622
DTWUMI
0.186
0.885
0.88
0.213
0.723
0.694
7.5%
FSMUMI
0.049
0.071
0.027
0.069
0.505
0.184
Amelia
0.197
0.998
0.147
0.045
1.305
0.792
FcM
0.158
0.809
0.104
1.813
1.866
0.991
MI
0.2
0.992
0.15
0.038
1.645
0.797
MICE
0.205
0.988
0.15
0.057
10.744
0.799
missForest
0.188
0.97
0.136
0.284
4.396
0.812
na.approx
0.192
0.971
0.142
0.669
2.163
0.712
DTWUMI
0.133
0.653
0.908
0.064
1.113
0.571
10%
FSMUMI
0.061
0.181
0.043
0.114
0.511
0.26
Amelia
0.202
0.999
0.147
0.034
4.062
0.788
FcM
0.164
0.872
0.104
1.837
2.201
0.992
MI
0.21
0.997
0.155
0.12
2.954
0.785
MICE
0.209
0.996
0.15
0.055
3.994
0.779
missForest
0.194
0.97
0.135
0.308
3.024
0.811
na.approx
0.183
0.997
0.129
0.372
1.455
0.719
DTWUMI
0.155
0.782
0.893
0.026
1.182
0.626Table 2
Average imputation performance indices of various imputation algorithms on simulated dataset (32,000 collected points).
Gap size
Method
Accuracy indices
Shape indices
1-Sim
1-R
2
RMSE
FSD
FB
1-FA2
1%
FSMUMI
0.083
0.515
1.033
0.159
2.51
0.574
Amelia
0.157
1
2.206
0.232
3.619
0.794
FcM
0.118
0.998
1.483
1.98
2.015
0.998
MI
0.16
0.999
2.241
0.2
0.915
0.799
MICE
0.159
0.998
2.201
0.214
1.449
0.801
missForest
0.127
0.998
1.608
0.836
12.034
0.861
na.approx
0.146
0.992
1.901
0.393
18.997
0.777
DTWUMI
0.09
0.552
1.156
0.007
6.022
0.562
2%
FSMUMI
0.068
0.487
1.166
0.194
1.971
0.611
Amelia
0.12
0.998
2.312
0.107
2.191
0.794
FcM
0.093
0.999
1.672
1.985
1.96
0.998
MI
0.12
1
2.307
0.123
3.949
0.789
MICE
0.119
0.999
2.282
0.114
8.881
0.789
missForest
0.096
1
1.769
0.941
2.777
0.858
na.approx
0.118
1
2.261
0.721
2.059
0.786
DTWUMI
0.074
0.523
1.545
0.008
3.686
0.583
3%
FSMUMI
0.068
0.453
1.053
0.076
10.649
0.582
Amelia
0.13
0.999
2.212
0.062
3.779
0.794
FcM
0.098
0.999
1.526
1.984
2.22
0.997
MI
0.13
0.999
2.197
0.078
9.374
0.795
MICE
0.129
1
2.19
0.067
1.938
0.792
missForest
0.102
0.999
1.626
0.855
2.407
0.851
na.approx
0.116
0.997
1.938
0.518
1.974
0.818
DTWUMI
0.073
0.526
1.189
0.01
8.725
0.567
4%
FSMUMI
0.064
0.412
1.067
0.061
1.374
0.568
Amelia
0.122
1
2.305
0.032
2.446
0.764
FcM
0.096
1
1.607
1.982
2.325
0.997
MI
0.125
1
2.261
0.043
2.391
0.792
MICE
0.124
0.999
2.233
0.045
42.495
0.791
missForest
0.101
1
1.726
0.876
2.901
0.854
na.approx
0.109
1
1.99
0.475
1.94
0.811
DTWUMI
0.066
0.465
1.172
0.004
2.079
0.547
5%
FSMUMI
0.063
0.404
1.062
0.062
4.508
0.577
Amelia
0.122
1
2.273
0.028
4.109
0.798
FcM
0.092
1
1.619
1.984
2.192
0.998
MI
0.123
1
2.287
0.024
5.582
0.797
MICE
0.121
1
2.267
0.044
2.326
0.792
missForest
0.097
0.999
1.731
0.923
2.473
0.859
na.approx
0.114
1
1.988
0.567
2.247
0.809
DTWUMI
0.063
0.454
1.166
0.003
1.594
0.545
7.5%
FSMUMI
0.06
0.408
1.063
0.049
4.843
0.566
Amelia
0.117
1
2.232
0.034
3.306
0.792
FcM
0.09
1
1.605
1.981
3.562
0.998
MI
0.119
0.999
2.259
0.025
1.946
0.793
MICE
0.118
1
2.238
0.032
9.359
0.794
missForest
0.094
0.999
1.695
0.907
1.259
0.858
na.approx
0.108
1
1.958
0.461
3.089
0.816
DTWUMI
0.065
0.477
1.19
0.004
3.851
0.566
10%
FSMUMI
0.061
0.4226
1.086
0.051
5.558
0.572
Amelia
0.117
1
2.269
0.021
3.074
0.793
FcM
0.089
1
1.607
1.981
2.683
0.997
MI
0.118
0.9996
2.233
0.02
2.05
0.793
MICE
0.118
0.9998
2.254
0.018
3.424
0.793
missForest
0.094
0.9999
1.702
0.909
1.87
0.857
na.approx
0.11
1
1.958
0.541
2.006
0.798
DTWUMI
0.067
0.5371
1.293
0.012
3.093
0.577Table 3
Average imputation performance indices of various imputation algorithms on MAREL-Carnot dataset (35,334 collected points).
Gap size
Method
Accuracy indices
Shape indices
1-Sim
1-R
2
RMSE
FSD
FB
1-FA2
1%
FSMUMI
0.051
0.156
1.532
0.044
0.081
0.191
Amelia
0.187
0.544
5.132
0.378
0.354
0.482
FcM
0.156
0.342
4.037
0.4
0.347
0.338
MI
0.192
0.561
5.282
0.396
0.365
0.497
MICE
0.166
0.608
5.596
0.423
0.35
0.436
missForest
0.165
0.472
4.422
0.385
0.355
0.381
na.approx
0.061
0.171
1.748
0.067
0.06
0.161
DTWUMI
0.084
0.181
2.466
0.214
0.149
0.198
2%
FSMUMI
0.045
0.037
1.446
0.053
0.083
0.182
Amelia
0.146
0.369
4.743
0.211
0.222
0.429
FcM
0.116
0.06
3.418
0.415
0.237
0.231
MI
0.146
0.364
4.72
0.218
0.228
0.435
MICE
0.129
0.369
4.711
0.197
0.21
0.413
missForest
0.116
0.155
3.575
0.33
0.193
0.258
na.approx
0.06
0.07
2.012
0.045
0.094
0.214
DTWUMI
0.042
0.018
1.095
0.029
0.066
0.154
3%
FSMUMI
0.053
0.11
1.294
0.134
0.08
0.166
Amelia
0.176
0.503
4.694
0.426
0.224
0.478
FcM
0.139
0.251
3.35
0.441
0.237
0.314
MI
0.17
0.531
4.474
0.354
0.221
0.476
MICE
0.157
0.552
4.905
0.34
0.184
0.429
missForest
0.139
0.345
3.556
0.422
0.184
0.346
na.approx
0.068
0.224
1.79
0.062
0.056
0.169
DTWUMI
0.096
0.216
2.587
0.329
0.136
0.223
4%
FSMUMI
0.059
0.058
1.466
0.094
0.101
0.183
Amelia
0.171
0.44
4.389
0.287
0.2
0.456
FcM
0.126
0.152
2.779
0.285
0.203
0.727
MI
0.166
0.41
4.234
0.277
0.204
0.444
MICE
0.15
0.379
4.15
0.268
0.19
0.411
missForest
0.129
0.234
3.134
0.23
0.187
0.303
na.approx
0.077
0.13
2.006
0.068
0.135
0.268
DTWUMI
0.07
0.105
1.77
0.15
0.12
0.138
5%
FSMUMI
0.051
0.22
2.025
0.227
0.152
0.167
Amelia
0.151
0.551
4.924
0.303
0.189
0.461
FcM
0.113
0.337
3.606
0.301
0.199
0.254
MI
0.143
0.567
4.612
0.249
0.123
0.448
MICE
0.131
0.523
4.75
0.274
0.188
0.419
missForest
0.104
0.371
3.443
0.229
0.147
0.274
na.approx
0.065
0.213
2.071
0.175
0.038
0.233
DTWUMI
0.067
0.275
2.363
0.22
0.157
0.242
7.5%
FSMUMI
0.043
0.056
1.52
0.075
0.039
0.189
Amelia
0.14
0.42
4.546
0.191
0.197
0.437
FcM
0.104
0.123
3.12
0.328
0.198
0.23
MI
0.142
0.427
4.624
0.222
0.222
0.443
MICE
0.126
0.38
4.375
0.206
0.208
0.437
missForest
0.112
0.202
3.587
0.329
0.228
0.288
na.approx
0.073
0.081
2.043
0.092
0.107
0.243
DTWUMI
0.06
0.102
1.999
0.071
0.074
0.215
10%
FSMUMI
0.053
0.098
1.642
0.083
0.055
0.191
Amelia
0.14
0.3
4.294
0.24
0.142
0.442
FcM
0.1
0.098
3.68
0.136
0.101
0.303
MI
0.14
0.112
4.294
0.24
0.142
0.442
MICE
0.12
0.42
4.066
0.152
0.077
0.383
missForest
0.097
0.461
3.049
0.104
0.117
0.255
na.approx
0.071
0.529
1.873
0.098
0.094
0.253
DTWUMI
0.081
0.381
3.293
0.119
0.124
0.224Synthetic Dataset. Table 1 presents a comparison of 8 imputation methods on synthetic dataset that contains 7 missing data levels (1-10%). The results clearly show that when a gap size is greater than 2%, the proposed method yields the highest similarity, R
2, FA2, and the lowest RMSE, FB. With this dataset, na.approx gives the best performance at the smallest missing data level for all indices and is ranked second for other ratios of missing values (2-5%) for similarity and FA2, RMSE (2-4%), and R
2 (the 1
s
t rank at 2% missing rate, the 2
n
d at 3%, 5%). The results can explain that the synthetic data are generated by a function (6). na.approx method applies the interpolation function to estimate missing values. So it is easy to find a function to generate values that are approximate real values when missing data rates are small. But this work is more difficult when the missing sample size rises; that is why the ability of na.approx decreases as missing data levels increase, especially at 7.5% and 10% rates. Although this dataset never exactly repeats itself and our approach is proposed under the assumption of recurrent data the FSMUMI approach proves its performance for the imputation task even if the missing size increases.Among the considered methods, the FcM-based approach is less accurate at lower missing rates but it provides better results at larger missing ratios as regards the accuracy indices.Simulated Dataset. Table 2 illustrates the evaluation results of various imputation algorithms on the simulated dataset. The best values for each missing level are highlighted in bold. Our proposed method outperforms other methods for the imputation task on accuracy indices: the highest similarity, R
2, and the lowest RMSE at every missing ratio. However, when considering other indices such as FA2, FSD, and FB, FSMUMI no longer shows its performance. It gains only at a 4% rate for the FB index and at 10% ratio for FA2. In contrast to FSMUMI, DTWUMI provides the best results for FSD indicator at all missing levels and FA2 at the first 5 missing ratios (from 1% to 5%).Different from the synthetic dataset, on the simulated dataset, the FcM-based method is always ranked the third at all missing rates for similarity and RMSE indicators. Following FcM is missForest algorithm for the both indices.Although, in the second experiment, data are built by various functions but they are quite complex so that na.approx does not provide good results.MAREL-Carnot Dataset. Once again, as reported in Table 3, our algorithm demonstrates its capability for the imputation task. FSMUMI method generates the best results as regarding accuracy indices for almost missing ratios (excluding at 2% missing level on all indices, and at 5% missing rate on R
2 score). But when considering shape indicators, FSMUMI only provides the highest FA2 values at several missing levels (3%, 5%-10%). In particular, our method illustrates the ability to fill in incomplete data with large missing rates (7.5% and 10%): the highest similarity, R
2, FA2, and the lowest RMSE, FSD (excluding at 7.5%), and FB. These gaps correspond to 110.4 and 147.2 days sampled at hourly frequency.In contrast to the two datasets above, on the MAREL-Carnot data, na.approx indicates quite good results: the permanent second or third rank for the accuracy indices (the1
s
t order at 5% missing rate on R
2 score), the lowest FSD (from 3% to 5% missing rates), and FB at some other levels of missing data. But when looking at the shape of imputation values generated from this method, it absolutely gives the worst results (Figure 6).Other approaches (including FcM-based imputation, MI, MICE, Amelia, and missForest) exploit the relations between attributes to estimate missing values. However, three considered datasets have low correlations between variables (roundly 0.2 for MAREL-Carnot data,≤
0.1 for simulated and synthetic datasets). So these methods do not demonstrate their performance for completing missing values in low/uncorrelated multivariate time series. Otherwise, our algorithm shows its ability and stability when applying to the imputation task for this kind of data.DTWUMI approach was proposed to fill large missing values in low/uncorrelated multivariate time series. However, this method is not as powerful as the FSMUMI method. DTWUMI only produces the best results at 2% missing level on the MAREL-Carnot dataset and is always at the second or the third rank at all the remaining missing rates on the MAREL-Carnot and the simulated datasets. That is because the DTWUMI method only finds the most similar window to a query either before a gap or after this gap, and it uses only one similarity measure, the DTW cost, to retrieve the most similar window. In addition, another reason may be that DTWUMI has directly used data from the window following or preceding the most similar window to completing the gap.


### 5.2. Visual Performance Comparison

In this paper, we also compare the visualization performance of completion values yielded by various algorithms. Figures4 and 5 illustrate the form of imputed values generated from different approaches on the synthetic series at two missing ratios 1% and 5%.Figure 4
Visual comparison of completion data of different imputation approaches with real data on the1
s
t signal of synthetic series with the gap size of 1000.Figure 5
Visual comparison of completion data of different imputation approaches with real data on the1
s
t signal of synthetic series with the gap size of 5000.Figure 6
Visual comparison of completion data of different imputation approaches with real data on the2
n
d signal of MAREL-Carnot dataset with the gap size of 353.At a 1% missing rate, the shape of imputation values produced by na.approx method is closer to the one of true values than the form of completion values given by our approach. However, at a 5% level of missing data, this method no longer shows the performance (Figure5). In this case, the proposed method proves its relevance for the imputation task. The shape of FSMUMI’s imputation data is almost similar to the form of true values (Figure 5).Looking at Figure6, FSMUMI one more time proves its capability for uncorrelated multivariate time series imputation: completion values yielded by FSMUMI are virtually identical to the real data on the MAREL-Carnot dataset. When comparing DTWUMI with FSMUMI, it is clear that FSMUMI gives improved results (Figures 4, 5, and 6).


### 5.3. Computation Time

Besides, we perform a comparison of the computational time of each method on the synthetic series (in second - s). Table4 indicates that na.approx method requires the shortest running time and DTWUMI approach takes the longest computing time. The proposed method, FSMUMI, demands more execution time as missing rates increase. However, considering the quantitative and visual performance of FSMUMI for the imputation task (Table 1, Figures 5 and 6), the required time of the proposed approach is fully acceptable.Table 4
Computational time of different methods on the synthetic series in second (s).
Method
Gaps size
1%
2%
3%
4%
5%
7.5%
10%
FSMUMI
353.9
427.5
701.9
1037.8
1423.6
2525.5
3556.8
Amelia
3.2
3.4
5.2
3.2
3.2
3.2
3.2
FcM
40.9
39.8
40.0
41.1
41.2
46.7
45.6
MI
844.1
714.0
739.1
723.3
724.5
719.7
726.5
MICE
7021.1
9187.7
21909.6
13041.9
14833.9
19417.7
23812.6
missForest
26833.8
24143.8
22969.9
32056.6
36485.8
42424.1
28521.1
na.approx
0.11
0.089
0.167
0.09
0.088
0.088
0.094
DTWUMI
5002.67
15714.8
37645.82
64669.71
86435.38
180887.78
273879


## 5.1. Quantitative Performance Comparison

Tables1, 2, and 3 illustrate the average ability of various imputation methods for synthetic, simulated, and MAREL-Carnot time series using 6 measurements as previously defined. For each missing level, the best results are highlighted in bold. These results demonstrate the improved performance of FSMUMI to complete missing data in low/uncorrelated multivariate time series.Table 1
Average imputation performance indices of various imputation algorithms on synthetic dataset (100,000 collected points).
Gap size
Method
Accuracy indices
Shape indices
1-Sim
1-R
2
RMSE
FSD
FB
1-FA2
1%
FSMUMI
0.136
0.261
0.051
0.358
3.253
0.364
Amelia
0.275
0.999
0.143
0.409
2.252
0.773
FcM
0.231
0.722
0.096
1.889
2.208
0.996
MI
0.275
0.999
0.142
0.421
2.091
0.773
MICE
0.258
0.944
0.13
0.406
2.452
0.72
missForest
0.248
0.915
0.122
0.389
3.976
0.744
na.approx
0.052
0.066
0.019
0.054
0.29
0.074
DTWUMI
0.257
0.713
0.88
0.725
0.405
0.69
2%
FSMUMI
0.1
0.295
0.046
0.155
0.395
0.337
Amelia
0.259
0.998
0.147
0.275
2.005
0.803
FcM
0.208
0.686
0.104
1.863
2.289
0.987
MI
0.259
0.998
0.147
0.268
2.11
0.81
MICE
0.244
0.968
0.14
0.255
7.616
0.759
missForest
0.239
0.968
0.133
0.279
3.156
0.792
na.approx
0.104
0.278
0.047
0.224
0.398
0.347
DTWUMI
0.237
0.775
0.867
0.509
8.449
0.646
3%
FSMUMI
0.113
0.341
0.056
0.219
0.852
0.322
Amelia
0.218
0.911
0.127
0.133
6.128
0.76
FcM
0.214
0.601
0.1
1.832
1.759
0.989
MI
0.253
0.993
0.141
0.236
2.295
0.775
MICE
0.21
0.873
0.118
0.208
5.118
0.703
missForest
0.188
0.796
0.102
0.215
1.846
0.627
na.approx
0.148
0.43
0.072
0.372
2.382
0.577
DTWUMI
0.231
0.799
0.874
0.332
27.952
0.69
4%
FSMUMI
0.06
0.146
0.037
0.099
0.738
0.299
Amelia
0.208
1
0.14
0.213
2.171
0.807
FcM
0.155
0.759
0.095
1.85
2.09
0.986
MI
0.208
0.999
0.14
0.196
2.302
0.807
MICE
0.209
0.987
0.138
0.22
3.748
0.801
missForest
0.196
0.968
0.127
0.216
3.94
0.827
na.approx
0.145
0.721
0.092
0.252
5.251
0.689
DTWUMI
0.148
0.586
0.918
0.185
12.688
0.719
5%
FSMUMI
0.055
0.132
0.032
0.058
0.098
0.201
Amelia
0.214
0.997
0.15
0.147
2.238
0.79
FcM
0.179
0.715
0.108
1.818
2.194
0.993
MI
0.231
0.996
0.167
0.206
3.094
0.808
MICE
0.221
0.968
0.152
0.222
2.3
0.79
missForest
0.212
0.944
0.143
0.315
4.547
0.819
na.approx
0.16
0.8
0.118
0.352
18.217
0.622
DTWUMI
0.186
0.885
0.88
0.213
0.723
0.694
7.5%
FSMUMI
0.049
0.071
0.027
0.069
0.505
0.184
Amelia
0.197
0.998
0.147
0.045
1.305
0.792
FcM
0.158
0.809
0.104
1.813
1.866
0.991
MI
0.2
0.992
0.15
0.038
1.645
0.797
MICE
0.205
0.988
0.15
0.057
10.744
0.799
missForest
0.188
0.97
0.136
0.284
4.396
0.812
na.approx
0.192
0.971
0.142
0.669
2.163
0.712
DTWUMI
0.133
0.653
0.908
0.064
1.113
0.571
10%
FSMUMI
0.061
0.181
0.043
0.114
0.511
0.26
Amelia
0.202
0.999
0.147
0.034
4.062
0.788
FcM
0.164
0.872
0.104
1.837
2.201
0.992
MI
0.21
0.997
0.155
0.12
2.954
0.785
MICE
0.209
0.996
0.15
0.055
3.994
0.779
missForest
0.194
0.97
0.135
0.308
3.024
0.811
na.approx
0.183
0.997
0.129
0.372
1.455
0.719
DTWUMI
0.155
0.782
0.893
0.026
1.182
0.626Table 2
Average imputation performance indices of various imputation algorithms on simulated dataset (32,000 collected points).
Gap size
Method
Accuracy indices
Shape indices
1-Sim
1-R
2
RMSE
FSD
FB
1-FA2
1%
FSMUMI
0.083
0.515
1.033
0.159
2.51
0.574
Amelia
0.157
1
2.206
0.232
3.619
0.794
FcM
0.118
0.998
1.483
1.98
2.015
0.998
MI
0.16
0.999
2.241
0.2
0.915
0.799
MICE
0.159
0.998
2.201
0.214
1.449
0.801
missForest
0.127
0.998
1.608
0.836
12.034
0.861
na.approx
0.146
0.992
1.901
0.393
18.997
0.777
DTWUMI
0.09
0.552
1.156
0.007
6.022
0.562
2%
FSMUMI
0.068
0.487
1.166
0.194
1.971
0.611
Amelia
0.12
0.998
2.312
0.107
2.191
0.794
FcM
0.093
0.999
1.672
1.985
1.96
0.998
MI
0.12
1
2.307
0.123
3.949
0.789
MICE
0.119
0.999
2.282
0.114
8.881
0.789
missForest
0.096
1
1.769
0.941
2.777
0.858
na.approx
0.118
1
2.261
0.721
2.059
0.786
DTWUMI
0.074
0.523
1.545
0.008
3.686
0.583
3%
FSMUMI
0.068
0.453
1.053
0.076
10.649
0.582
Amelia
0.13
0.999
2.212
0.062
3.779
0.794
FcM
0.098
0.999
1.526
1.984
2.22
0.997
MI
0.13
0.999
2.197
0.078
9.374
0.795
MICE
0.129
1
2.19
0.067
1.938
0.792
missForest
0.102
0.999
1.626
0.855
2.407
0.851
na.approx
0.116
0.997
1.938
0.518
1.974
0.818
DTWUMI
0.073
0.526
1.189
0.01
8.725
0.567
4%
FSMUMI
0.064
0.412
1.067
0.061
1.374
0.568
Amelia
0.122
1
2.305
0.032
2.446
0.764
FcM
0.096
1
1.607
1.982
2.325
0.997
MI
0.125
1
2.261
0.043
2.391
0.792
MICE
0.124
0.999
2.233
0.045
42.495
0.791
missForest
0.101
1
1.726
0.876
2.901
0.854
na.approx
0.109
1
1.99
0.475
1.94
0.811
DTWUMI
0.066
0.465
1.172
0.004
2.079
0.547
5%
FSMUMI
0.063
0.404
1.062
0.062
4.508
0.577
Amelia
0.122
1
2.273
0.028
4.109
0.798
FcM
0.092
1
1.619
1.984
2.192
0.998
MI
0.123
1
2.287
0.024
5.582
0.797
MICE
0.121
1
2.267
0.044
2.326
0.792
missForest
0.097
0.999
1.731
0.923
2.473
0.859
na.approx
0.114
1
1.988
0.567
2.247
0.809
DTWUMI
0.063
0.454
1.166
0.003
1.594
0.545
7.5%
FSMUMI
0.06
0.408
1.063
0.049
4.843
0.566
Amelia
0.117
1
2.232
0.034
3.306
0.792
FcM
0.09
1
1.605
1.981
3.562
0.998
MI
0.119
0.999
2.259
0.025
1.946
0.793
MICE
0.118
1
2.238
0.032
9.359
0.794
missForest
0.094
0.999
1.695
0.907
1.259
0.858
na.approx
0.108
1
1.958
0.461
3.089
0.816
DTWUMI
0.065
0.477
1.19
0.004
3.851
0.566
10%
FSMUMI
0.061
0.4226
1.086
0.051
5.558
0.572
Amelia
0.117
1
2.269
0.021
3.074
0.793
FcM
0.089
1
1.607
1.981
2.683
0.997
MI
0.118
0.9996
2.233
0.02
2.05
0.793
MICE
0.118
0.9998
2.254
0.018
3.424
0.793
missForest
0.094
0.9999
1.702
0.909
1.87
0.857
na.approx
0.11
1
1.958
0.541
2.006
0.798
DTWUMI
0.067
0.5371
1.293
0.012
3.093
0.577Table 3
Average imputation performance indices of various imputation algorithms on MAREL-Carnot dataset (35,334 collected points).
Gap size
Method
Accuracy indices
Shape indices
1-Sim
1-R
2
RMSE
FSD
FB
1-FA2
1%
FSMUMI
0.051
0.156
1.532
0.044
0.081
0.191
Amelia
0.187
0.544
5.132
0.378
0.354
0.482
FcM
0.156
0.342
4.037
0.4
0.347
0.338
MI
0.192
0.561
5.282
0.396
0.365
0.497
MICE
0.166
0.608
5.596
0.423
0.35
0.436
missForest
0.165
0.472
4.422
0.385
0.355
0.381
na.approx
0.061
0.171
1.748
0.067
0.06
0.161
DTWUMI
0.084
0.181
2.466
0.214
0.149
0.198
2%
FSMUMI
0.045
0.037
1.446
0.053
0.083
0.182
Amelia
0.146
0.369
4.743
0.211
0.222
0.429
FcM
0.116
0.06
3.418
0.415
0.237
0.231
MI
0.146
0.364
4.72
0.218
0.228
0.435
MICE
0.129
0.369
4.711
0.197
0.21
0.413
missForest
0.116
0.155
3.575
0.33
0.193
0.258
na.approx
0.06
0.07
2.012
0.045
0.094
0.214
DTWUMI
0.042
0.018
1.095
0.029
0.066
0.154
3%
FSMUMI
0.053
0.11
1.294
0.134
0.08
0.166
Amelia
0.176
0.503
4.694
0.426
0.224
0.478
FcM
0.139
0.251
3.35
0.441
0.237
0.314
MI
0.17
0.531
4.474
0.354
0.221
0.476
MICE
0.157
0.552
4.905
0.34
0.184
0.429
missForest
0.139
0.345
3.556
0.422
0.184
0.346
na.approx
0.068
0.224
1.79
0.062
0.056
0.169
DTWUMI
0.096
0.216
2.587
0.329
0.136
0.223
4%
FSMUMI
0.059
0.058
1.466
0.094
0.101
0.183
Amelia
0.171
0.44
4.389
0.287
0.2
0.456
FcM
0.126
0.152
2.779
0.285
0.203
0.727
MI
0.166
0.41
4.234
0.277
0.204
0.444
MICE
0.15
0.379
4.15
0.268
0.19
0.411
missForest
0.129
0.234
3.134
0.23
0.187
0.303
na.approx
0.077
0.13
2.006
0.068
0.135
0.268
DTWUMI
0.07
0.105
1.77
0.15
0.12
0.138
5%
FSMUMI
0.051
0.22
2.025
0.227
0.152
0.167
Amelia
0.151
0.551
4.924
0.303
0.189
0.461
FcM
0.113
0.337
3.606
0.301
0.199
0.254
MI
0.143
0.567
4.612
0.249
0.123
0.448
MICE
0.131
0.523
4.75
0.274
0.188
0.419
missForest
0.104
0.371
3.443
0.229
0.147
0.274
na.approx
0.065
0.213
2.071
0.175
0.038
0.233
DTWUMI
0.067
0.275
2.363
0.22
0.157
0.242
7.5%
FSMUMI
0.043
0.056
1.52
0.075
0.039
0.189
Amelia
0.14
0.42
4.546
0.191
0.197
0.437
FcM
0.104
0.123
3.12
0.328
0.198
0.23
MI
0.142
0.427
4.624
0.222
0.222
0.443
MICE
0.126
0.38
4.375
0.206
0.208
0.437
missForest
0.112
0.202
3.587
0.329
0.228
0.288
na.approx
0.073
0.081
2.043
0.092
0.107
0.243
DTWUMI
0.06
0.102
1.999
0.071
0.074
0.215
10%
FSMUMI
0.053
0.098
1.642
0.083
0.055
0.191
Amelia
0.14
0.3
4.294
0.24
0.142
0.442
FcM
0.1
0.098
3.68
0.136
0.101
0.303
MI
0.14
0.112
4.294
0.24
0.142
0.442
MICE
0.12
0.42
4.066
0.152
0.077
0.383
missForest
0.097
0.461
3.049
0.104
0.117
0.255
na.approx
0.071
0.529
1.873
0.098
0.094
0.253
DTWUMI
0.081
0.381
3.293
0.119
0.124
0.224Synthetic Dataset. Table 1 presents a comparison of 8 imputation methods on synthetic dataset that contains 7 missing data levels (1-10%). The results clearly show that when a gap size is greater than 2%, the proposed method yields the highest similarity, R
2, FA2, and the lowest RMSE, FB. With this dataset, na.approx gives the best performance at the smallest missing data level for all indices and is ranked second for other ratios of missing values (2-5%) for similarity and FA2, RMSE (2-4%), and R
2 (the 1
s
t rank at 2% missing rate, the 2
n
d at 3%, 5%). The results can explain that the synthetic data are generated by a function (6). na.approx method applies the interpolation function to estimate missing values. So it is easy to find a function to generate values that are approximate real values when missing data rates are small. But this work is more difficult when the missing sample size rises; that is why the ability of na.approx decreases as missing data levels increase, especially at 7.5% and 10% rates. Although this dataset never exactly repeats itself and our approach is proposed under the assumption of recurrent data the FSMUMI approach proves its performance for the imputation task even if the missing size increases.Among the considered methods, the FcM-based approach is less accurate at lower missing rates but it provides better results at larger missing ratios as regards the accuracy indices.Simulated Dataset. Table 2 illustrates the evaluation results of various imputation algorithms on the simulated dataset. The best values for each missing level are highlighted in bold. Our proposed method outperforms other methods for the imputation task on accuracy indices: the highest similarity, R
2, and the lowest RMSE at every missing ratio. However, when considering other indices such as FA2, FSD, and FB, FSMUMI no longer shows its performance. It gains only at a 4% rate for the FB index and at 10% ratio for FA2. In contrast to FSMUMI, DTWUMI provides the best results for FSD indicator at all missing levels and FA2 at the first 5 missing ratios (from 1% to 5%).Different from the synthetic dataset, on the simulated dataset, the FcM-based method is always ranked the third at all missing rates for similarity and RMSE indicators. Following FcM is missForest algorithm for the both indices.Although, in the second experiment, data are built by various functions but they are quite complex so that na.approx does not provide good results.MAREL-Carnot Dataset. Once again, as reported in Table 3, our algorithm demonstrates its capability for the imputation task. FSMUMI method generates the best results as regarding accuracy indices for almost missing ratios (excluding at 2% missing level on all indices, and at 5% missing rate on R
2 score). But when considering shape indicators, FSMUMI only provides the highest FA2 values at several missing levels (3%, 5%-10%). In particular, our method illustrates the ability to fill in incomplete data with large missing rates (7.5% and 10%): the highest similarity, R
2, FA2, and the lowest RMSE, FSD (excluding at 7.5%), and FB. These gaps correspond to 110.4 and 147.2 days sampled at hourly frequency.In contrast to the two datasets above, on the MAREL-Carnot data, na.approx indicates quite good results: the permanent second or third rank for the accuracy indices (the1
s
t order at 5% missing rate on R
2 score), the lowest FSD (from 3% to 5% missing rates), and FB at some other levels of missing data. But when looking at the shape of imputation values generated from this method, it absolutely gives the worst results (Figure 6).Other approaches (including FcM-based imputation, MI, MICE, Amelia, and missForest) exploit the relations between attributes to estimate missing values. However, three considered datasets have low correlations between variables (roundly 0.2 for MAREL-Carnot data,≤
0.1 for simulated and synthetic datasets). So these methods do not demonstrate their performance for completing missing values in low/uncorrelated multivariate time series. Otherwise, our algorithm shows its ability and stability when applying to the imputation task for this kind of data.DTWUMI approach was proposed to fill large missing values in low/uncorrelated multivariate time series. However, this method is not as powerful as the FSMUMI method. DTWUMI only produces the best results at 2% missing level on the MAREL-Carnot dataset and is always at the second or the third rank at all the remaining missing rates on the MAREL-Carnot and the simulated datasets. That is because the DTWUMI method only finds the most similar window to a query either before a gap or after this gap, and it uses only one similarity measure, the DTW cost, to retrieve the most similar window. In addition, another reason may be that DTWUMI has directly used data from the window following or preceding the most similar window to completing the gap.


## 5.2. Visual Performance Comparison

In this paper, we also compare the visualization performance of completion values yielded by various algorithms. Figures4 and 5 illustrate the form of imputed values generated from different approaches on the synthetic series at two missing ratios 1% and 5%.Figure 4
Visual comparison of completion data of different imputation approaches with real data on the1
s
t signal of synthetic series with the gap size of 1000.Figure 5
Visual comparison of completion data of different imputation approaches with real data on the1
s
t signal of synthetic series with the gap size of 5000.Figure 6
Visual comparison of completion data of different imputation approaches with real data on the2
n
d signal of MAREL-Carnot dataset with the gap size of 353.At a 1% missing rate, the shape of imputation values produced by na.approx method is closer to the one of true values than the form of completion values given by our approach. However, at a 5% level of missing data, this method no longer shows the performance (Figure5). In this case, the proposed method proves its relevance for the imputation task. The shape of FSMUMI’s imputation data is almost similar to the form of true values (Figure 5).Looking at Figure6, FSMUMI one more time proves its capability for uncorrelated multivariate time series imputation: completion values yielded by FSMUMI are virtually identical to the real data on the MAREL-Carnot dataset. When comparing DTWUMI with FSMUMI, it is clear that FSMUMI gives improved results (Figures 4, 5, and 6).


## 5.3. Computation Time

Besides, we perform a comparison of the computational time of each method on the synthetic series (in second - s). Table4 indicates that na.approx method requires the shortest running time and DTWUMI approach takes the longest computing time. The proposed method, FSMUMI, demands more execution time as missing rates increase. However, considering the quantitative and visual performance of FSMUMI for the imputation task (Table 1, Figures 5 and 6), the required time of the proposed approach is fully acceptable.Table 4
Computational time of different methods on the synthetic series in second (s).
Method
Gaps size
1%
2%
3%
4%
5%
7.5%
10%
FSMUMI
353.9
427.5
701.9
1037.8
1423.6
2525.5
3556.8
Amelia
3.2
3.4
5.2
3.2
3.2
3.2
3.2
FcM
40.9
39.8
40.0
41.1
41.2
46.7
45.6
MI
844.1
714.0
739.1
723.3
724.5
719.7
726.5
MICE
7021.1
9187.7
21909.6
13041.9
14833.9
19417.7
23812.6
missForest
26833.8
24143.8
22969.9
32056.6
36485.8
42424.1
28521.1
na.approx
0.11
0.089
0.167
0.09
0.088
0.088
0.094
DTWUMI
5002.67
15714.8
37645.82
64669.71
86435.38
180887.78
273879


## 6. Conclusion

This paper proposes a novel approach for uncorrelated multivariate time series imputation using a fuzzy logic-based similarity measure, namely FSMUMI. This method makes it possible to manage uncertainty with the comprehensibility of linguistic variables. FSMUMI has been tested on different datasets and compared with published algorithms (Amelia II, FcM, MI, MICE, missForest, na.approx, and DTWUMI) on accuracy and shape criteria. The visual ability of these approaches is also investigated. The experimental results definitely highlight that the proposed approach yielded improved performance in accuracy over previous methods in the case of multivariate time series having large gaps and low or non-correlation between variables. However, it is necessary to make an assumption of recurrent data and sufficiently large dataset to apply the algorithm. This means that our approach needs patterns (in our case the two queries (before and after the considered gap)) existing somewhere in the database. This enables us to predict missing values if the patterns occur in the past or in the following data from the considered position. Thus a satisfactory and sufficient dataset (large dataset) is required.In future work, we plan to (i) combine FSMUMI method with other algorithms such as Random Forest or Deep learning in order to efficiently fill incomplete values in any type of multivariate time series; (ii) investigate this approach applied to short-term/long-term forecasts in multivariate time series. We could also investigate complex fuzzy sets ([63]) instead of ordinary fuzzy sets that have given good results using an adaptive scheme in the case of the bivariate time series with small dataset.



---

*Source: 9095683-2018-08-09.xml*