Tensor Data Imputation by PARAFAC with Updated Chaotic Biases by Adam Optimizer

The big data pattern analysis suffers from incorrect responses due to missing data entries in the real world. Data collected for digital movie platforms like Netflix and intelligent transportation systems is Spatio-temporal data. Extracting the latent and explicit features from this data is a challenge. We present the high dimensional data imputation problem as a higher-order tensor decomposition. The regularized and biased PARAFAC decomposition is proposed to generate the missing data entries. The biases are created and updated by a chaotic exponential factor in Adam's optimization, which reduces the imputation error. This chaotic perturbed exponentially update in the learning rate replaces the fixed learning rate in the bias update by Adam optimization. The idea has experimented with Netflix and traffic datasets from Guangzhou, China.


INTRODUCTION
The low quality and homogenous data require the verification and validation process, which is more timeconsuming. Data mining approaches are used to manage a large amount of data in a particular context. Data mining schemes increases the quality of data with a large amount of data set. The meaningful pattern and rules are discovering in a large amount of data by exploration and analyzing them. In data mining approaches, various numbers of models and methods are considered to filter large data. The data mining methods are divided into different categories based on objective evaluation. The main categories are classification, estimation, prediction, clustering, association, and profiling. As we studied in the literature [1][2][3][4][5][6][7][8][9][10][11], various data mining models are proposed for Big Data analytics. The missing data problem is a setback for big data analytics. In the modern era, an intelligent transportation system (ITS) is used to collect large traffic data. Various mobile and stationary sensors are placed in the ITS configuration. The correct information collection in the transportation system represents the state of the system for the development of correct analysis. The collected data is not directly connected to the transportation system due to the presence of traditional data (traffic volume) and a new kind of data (Bluetooth and GPS data). The missing data problem has commonly arrived in a large traffic dataset. The missing data can degrade the performance of the intelligent transportation system. To improve the reliability of the traffic data collection and record the data mining approaches are used in state of the art. Manuscript  Previously various data imputation methods are proposed to solve the missing traffic data problem. The traffic data imputation problems are categorized into two sections: traffic volume imputation and traffic speed data estimation. Both traffic volume and traffic speed data have similar characteristics, like temporal and spatial correlation. The high dimensional data is converted into the tensor as it inherits the high algebraic information than the matrix. The purpose of the selection of traffic and Netflix data is, these contain the temporal as well as spatial resolutions. The pattern is learned by tensor matrix factorization for a part of the data, and that learning is used to impute the missing entries. This paper aims to propose the advanced tensor factorization for data imputation.

Contribution
Our contribution to the work is two-fold: regularized tensor factorization and Adam optimization with dynamic learning rate. The data used in our work is a tensor matrix that is a multidimensional matrix. Several tensor factorization schemes are available like CANDECOMP/PARFAC, tucker, SVD, etc. [20,[22][23][24][25][26][27][28], but the PARAFAC method provides more stable tensor decomposition with faster convergence. The tensor is decomposed to rank one matrix. Few researchers used tensor decomposition for the high dimensional data imputations too [29] [30]. The PARAFAC decomposition still doesn't deal with the unexpected differences in the data values, which are counted as noise in the data. We propose here modified regularized PARAFAC, which decomposes the tensor into three factors [20] as well as into bias matrices for those tensor factors which avoid the fall into local optima. We used the dynamic learning rate to update the bias values during Adam's optimization. The fixed learning rate of 0.001 is changed to a logistic mapped chaotic exponential decaying learning rate with every iteration.

Paper Organization
Further, in this article, we have discussed previous related work in section 2. Section 3 discusses the problem in data imputation and PARAFC tensor factorization. The proposed extended PARAFC tensor factorization for imputation has been discussed in section 4, which is followed by simulation results in the next section. We conclude the results in section 5.

II. LITERATURE REVIEW
State of the art schemes is utilized for the data mining process to achieve processed data. A high-quality data mining approach was presented in [1] by using the K-nearest neighbor algorithm. A predictive data mining framework was presented with the help of various up to date algorithms schemes like machine learning [2].
The machine learning algorithms improves the computational efficiency in removal of unwanted noise from the massive dataset. The Ghost algorithm was used to reconstruct the missing data from a large dataset that recover off period segments of missing data [3]. The ghost algorithm searched in a sequential manner dataset to achieve the data segment. A caching approach was also used to minimize the search space and improves the computational complexity to linear. In [4], various deep generative models were compared in terms of data mining. A Graph Neural Network (GNN) model was enhanced for the missing data imputation (MDI) application. With the help of GNN, a Graph convolutional autoencoder was designed to reconstruct the complete dataset [5]. The adversarial loss and global information were included in the dataset during the reconstruction phase. A Distributed Neural Network (DNN) was designed for the imputation of missing value in terms of Big data context [6]. The outcomes of the DNN data imputation method were compared with the K-nearest neighbor approach. They provided a better response. A selfdeveloped competitive neural network Adaptive Response Theory 2 (ART2) was proposed for the data imputation [7]. The ensemble approach provided accurate results on intra cluster non-missing value elements. The ART2 approach also improved imputation accuracy. The Big Data K-mean and Big Data Fuzzy imputation cluster-based methods were presented in [8]. Similarly, in [9], two algorithms were used to remove the noisy data from Big Data set. The two methods were Heterogeneous and Homogenous ensemble filters. These algorithms were formed with a combination of Deep learning classification models. An MLClean approach was presented for the Big Data preprocessing to trained the Deep Learning model [10]. It provided data cleaning, unfair mitigation, and data sanitization. The MLClean approach was commonly used in Big Data analytics cases. The ROSEFW-RF algorithm was developed based on the Mapreduce function for imbalanced data mining [11]. The medical field data is imbalanced due to its number of parameters. To achieve balanced data, the relevant features, classes, category, and values are classified with the proposed ROSEFW-RF approach. In a study [12], the author provided the Bayesian CP factorization to handle the noisy tensor data. A non-deterministic model was developed using Bayesian CP factorization, which estimated the rank of CP automatically. A robust Bayesian generative model developed in [13] for tensor factorization. It provided the robustness to the boundary conditions and reduced the Gaussian noise, which minimizes the overfitting problem. The low-rank decomposition using a Bayesian model for multiway tensors data with missing observations was presented in [14]. The single and coupled tensor factorization models were tested by the full Bayesian Inference using VB [15]. The full Bayesian Inference method was easily implemented to any large scale model. A fully Bayesian Probabilistic Matrix Factorization (PMF) model was presented in [16]. The overall model parameters and hyperparameters controlled the model capacity. A factor-based algorithm was presented for the classification application known as tensor factorization. The fully Bayesian Probabilistic Tensor Factorization (PTF) was used to analyze the large scale dataset [17]. Some real-world problems of classification were analyzed by the Bayesian PTF model with accurate outcomes. A Bayesian Gaussian CANDECOMP/PARAFAC (BGCP) decomposition model was used to impute the multidimensional imbalance traffic data [18]. In tensor representation with BGCP, the accurate outcomes were found in the third-order case for both random missing and fiber missing traffic data scenarios. The multiple metrics and tensors combined joints were presented by the Bayesian Multi-tensor Factorization model [19]. The Bayesian MTF model analyzed the total number of data collection factors. A multidimensional EEG dataset was analyzed by the Bayesian Tensor Factorization model [20]. The noninformative EEG signals were removed from the BTF method, which improved the classification response in the medical research field. A fully Bayesian system was presented for the automatically learning parameters of missing values traffic data model with Variational Bayes (VB) [21]. The data was collected from the urban traffic speed data set collected in Guangzhou, China. The two methods were tested on the traffic data imputation. The Bayesian Augmented Tensor Factorization (BATF) provided better accuracy in terms of missing values imputations. A three procedure framework was proposed for the missing traffic data recovery with traffic patterns recognition. The proposed STD approach provided higher accuracy in terms of missing traffic data [22]. A matrix factorization technique based on a K-nearest neighbor method was proposed for the Netflix missing data recovery [23]. The matrix factorization based KNN model provided higher accuracy, and it can be easily implemented on realworld datasets. The collaborative filtering approaches were applied for the data filter in [24,25]. The unique properties of implicit feedback datasets were recognized by a factor model [25] III.

PROBLEM STATEMENT
Due to vast applications of tensor, its decomposition for multidimensional data has been the researcher's interest. Their work can be broadly classified into two categories: decomposition by CANDECOMP/PARAFAC (Canonical Decomposition/ Parallel Factor), also referred to as CP and tucker factorization. Both are higher-order principal component analysis (PCA). The convex hull analysis in [36] proves the CP decomposition is better than tucker decomposition on the criteria of good fit analysis. The PARAFAC itself takes a long time in reaching the local minima solution and also has the probability to stuck in the local minima. To avoid it with local optima, we regularized the PARAFAC Tikhonov regularization factor. Still, this form of PARAFAC only considers the latent features, but the data is Spatio-temporal. So, we need to add the biases for all the three dimensions of the data, which can extract the temporal features too. Besides the decomposition part, the rank selection is also uncertain in the tensor decomposition. The Bayesian approach has also been suggested by a few researchers for rank selection [12][13][14][15][16][17][18][19][20][21]. But these focused either on prediction problems or factorization with missing data. Whereas we couldn't find it suitable for data imputation problem as Bayesian optimization works preferably for automatic Rank detection for higher-order low-rank factorization. In the imputation application of CP, the normalized imputation error increases with the increase in rank. We have tested the scheme of regularized augmented PARAFACfor various ranks from 1-10. A plot in figure 1 demonstrates that rank 1 is suitable for data imputation scheme and imputed data for this rank by PARAFAC, and regularized augmented PARAFAC is shown in figure 1(b). It supports our convention that PARAFAC lacks in providing a global solution for missing data generation. The matrices ∈ ℝ × , ∈ ℝ × , ∈ ℝ × , and ∈ ℝ × are generated by the random entries drawn from the standard normal distribution. A new matrix ∈ ℝ × × is generated by a random setting of % of the entries of to be missing. So = * , here binary tensor is represented by the which has an equal dimension of tensor and , , = 0 for the set , , to missing. Figure 2 shows the tensor with random missing entries.

Fig. 2 Arrangement of incomplete tensor by random
indexing.

IV. PROPOSED TENSOR FACTORIZATION A. Regularized-PARAFAC
As discussed above, we factorize the tensor matrix by the PARAFAC method. For the data imputation in our work, the data is of third-order tensor ℝ × × and as for large rank ( > 10), the PARFAC decomposition can be shown as in equation 1 [29] [33].
≈ , , Equation 1 holds for the sufficiently large rank , but since in our case, we are interested in the decomposition up to = 1 [13]. The tensor matrix can be decomposed into N factors of rank 1 as: Where ∘ represents the cross product. For the third-order tensor, it is ∘ ∘ . The conviction of decomposition of equation 1 to the first-order tensor is: The || . || 2 Represents the Frobenius norm. This one is the non-convex optimization problem and can be solved iteratively as The equation 4 has an issue of local optimum solution due to non-uniform scaling of factors [33], so Paatero [28] suggested the addition of a Tikhonov regularization factor in equation 4 which resolved the issue of many local solutions in equation 3 and updated it as ℒ ; , , = Where is a regulation parameter with > 1. In our work, we are using the data set from China Traffic Data [13]. This kind of tensor data has three factors for road segments, the number of time slots, and recorded time intervals. Equation 5 maps the combined rating for users and content in the latent space, but it is not always the case. The temporal dependency on the data is not examined in that, due to which the imputation might be less correct. This biasing for user and content must be included in equation 5 [29]. The augmented equation 5 is written in equation 6 2 ) (7) Where ∅ , , are the bias factors for the three tensor dimensions, is the average of all three factors [31]. The biases are the deviation of each factor from the average . The proposed tensor is shown in figure 3. The factor matrix of the tensor can also be represented as , , = + ∅ + + + [ , , ].

Fig. 3 Proposed augmented PARAFAC factorization for 3 rd order tensor
The proposed modified-PARAFAC is an iterative process, and equation 7 is minimized by gradient descent optimization [13], Alternating least square method [31] recently. But these methods are easily stuck in local minima and less rate of convergence. The advancement in deep learning is dependent upon the fast and global optimized Adam optimization [35]. We factorize the tensor by Adam optimization. The bias is also updated iteratively in each stochastic step. These are initialized at some random value between 0 and 0.1. Though the initial feeding point for these biases can be random, however for our application, we set 0,0.1 , where represents the biases. Biases are updated in each iteration of Adam's optimization. Equation  6 is the objective function of adam optimization, and the error thus calculated is used to update the bias values. We multiply the error from equation 6 by an exponentially decreasing regularization parameter, which reduces with iteration to minimize the error. To add the chaotic perturbation in error, we multiply it with a less random factor, which is calculated by a logistic map. The biases are updated as: Here − , , is the deviation in the factors, is a constant whose values if fixed to 20 in our case. and represent the current iteration and maximum iteration in the adam optimization.
is the perturbed parameter by logistic mapping which is calculated as = 1 − , = 0,1,2,3 ….. The is a system parameter 0,4 0,1 . This way and are also updated. The pseudo-code for the update of biases is shown in algorithm1. The Netflix dataset provides a testset containing 2,817,131 movie pairs with ratings. The latest movie pairs ratings are selected from a recent subset of users in the training dataset.

C. Evaluation
We have followed the random missing entries scheme for missing values generation. After creating the missing values for 10%, 20%, and 50%, we use rank 1 from figure 1 for the tensor data completion. The proposed decomposition solution with bias update is compared with the state of the art schemes. We compared the imputation with Tucker decomposition [30], Bayesian inference [29] for the random missing data. The Netflix dataset is also checked for only random missing entries. The proposed scheme is also evaluated for gradient descent (GD) and Adam optimization for the optimal solution of equation 6. The evaluation of results is statistically done based on mean absolute error (MAE), mean absolute percentage error (MAPE), and root means square error (RMSE). If the and are estimated and true missing values, then these performance measures are defined as: Netflix dataset and traffic datasets are evaluated by the proposed scheme and compared with states of the art SVD combined tucker decomposition [30], BATF [29], BGCP [20], BCPF [21], and modified PARAFAC variants. The schemes in [20], [21], [29] have common Bayesian Gaussian tensor factorization. The CANDECOMP/PARAFAC factors are trained by Bayesian optimization. The [21] presented the hyperparameter of automatic rank determination. In contrast, the solution doesn't differ significantly from other variants in [20] and [29] as the researchers only used the higher number of hyperparameters for optimization. However, the Bayesian optimization of factorization parameters performed well over Tucker-SVD decomposition [30]. Table 1 shows the performance comparison of the proposed data imputation scheme for random missing entries for traffic datasets. It can be observed from table 1 that performance is somewhat good at lower missing ratio; however, the change is not significant with increasing missing values. The dynamic, chaotic learning rate and modified regularized tensor decomposition lead to superior performance than state of the art methods. Figure 4 represents the imputed tensor for different random missing ratios for the traffic dataset by the proposed scheme. The proposed scheme doesn't get affected by the change in the missing rate like other schemes.  Fig. 4 The imputed tensor with random missing entries percentage (a) 50% (b) 30 % and (c) 10% The bias to generate the missing values is updated by gradient descent (GD) optimization in [30]; however, as discussed in the previous section, the Adam optimizer is better than GD. We have followed the update as in algorithm 2. The SVD-Tucker decomposition proposed the constant learning rate to update the bias. Still, it has been proven recently that for the sparse dataset, the adaptive learning rate is better than a constant learning rate. The learning rate can be adaptive dependent upon time, step size, or exponentially changing with iteration. Other Deep learning toolboxes like Keras also provide the feature to use the adaptive learning rate during training. We have perturbed the exponentially decaying learning rate by logistic mapping. It is a linear mapping with variable = 1 − , = 0,1,2,3 ….. The is a system parameter 0,4 0,1 . The logistic mapping shows different behavior for different values of . Figure 5 shows the stationary, periodic, and complete bifurcation diagram for the logistic map. Figure 5(a) and 5(b) shows the behavior of logistic mapping if the system parameter is less than 1 and 3 ≤ < 1 + 6 respectively. For 0 ≤ < 1, updated biases with this logistic map value can lead to premature convergence of the training. Conclusively, data imputation won't be similar to extracted latent features of the rest of the data. The non-diminishing oscillatory behavior ( figure 5(b)) also doesn't converge at all. The chaotic behavior starts beyond 3.56994, and figure 5(c) shows the logistic mapping for a complete range of ∈ 0,4 . The area between ∈ 3.54409,4 is the stable oscillations area, and convergence in optimization can be achieved in this area.
The solid lines in figure 5(c) point to the stable solution. The more uniform is the bias, the less deviation is observed in variance for a dimension. Just to recall, the three biases are for the three dimensions of the Spatio-temporal data. The size of a bias matrix is equal to the size of the decomposed factor matrices. The first row in figure 6 shows the uniformity comparison for the bias matrices for the final selected solution for a 50% missing ratio. We have used the ternary plot since we have 3-dimensional tensor. The final updated bias for the minimum gradient in equation 7 is recorded for plotting. The vertices of the ternary plot represent the updated bias by Adam optimized modified regulated PARAFAC with a dynamic learning rate, with a fixed learning rate and proposed PARAFAC segmentation by GD optimization. These three plots are for each dimension of the data. The dots inside the plots are the bias values that arrange themselves in an elliptical shape with alignment towards the proposed scheme. The elliptical shape is formed due to uniformity in the bias values; the least is the standard deviation; more is the uniformity. The alignment towards proposed scheme vertex indicates the clustering of biases is greatly affected by larger bias magnitudes generated by the proposed scheme. However, the uniformity of bias < ∅ & , as indicated in figure  6(b) for bias for the number of days. The second row in figure 6 has the plots for imputed values and original entries. We have plotted these curves for a few samples. The sparse tensor's entries are also plotted in 6(d) whose entries are all zeros for these samples. (b) for = .  For the Netflix dataset, we have compared the proposed scheme with its three other variants only due to unavailability of state of the art schemes results on this data. The figure shows the performance plot comparison for the Netflix dataset. This dataset has primarily tested for the Netflix prize competition for recommending the most suitable movies to the users. The modified PARAFAC with the proposed dynamic learning rate has minimum MAPE and RMSE. The MAPE has shown no significant change with the increase in the missing entries, but with 80% missing data, all variants have degraded. However, no such behavior is recorded for MAE and RMSE curves. The proposed scheme has the least MAPE and RMSE values at various missing ratios. This dataset also strengthens the theory that Adam optimization outperforms the PARAFAC with gradient descent optimization. The bias for each dimension of the Netflix data is shown in the second row of figure 7.

VI. CONCLUSION
The data imputation is the first step in the data denoising and to process the data with deep learning. In this research article, we have profoundly compared the performance of imputed tensor data with state of the art schemes. The Spatio-temporal real-world tensor data can be decomposed into the low-rank matrix for imputation using two kinds of decomposition primarily: CD/PARAFAC and Tucker decomposition. We compared the proposed modified PARAFAC with a dynamic learning rate with the latest work of Tucker-SVD decomposition on the Traffic dataset. It has been observed that our method is able to reduce the mean absolute error by 20% for the 50% missing entries in the data. Few other researchers have focused on the Bayesian learning training for the tensor decomposition. It has also been suggested previously that automatic rank determination leads to a better generation of missing entries. Our dynamic Adam optimization scheme has powered the modified PARAFAC with a 0.48% reduced root mean square error for 50 % missing data. We have tested the algorithm on Netflix movies dataset too. The size of the dataset has backed the performance due to our hardware constraint, but it has been noticed that dynamically learned modified PARAFAC had won the competition. In the next part of our research, we will test the performance on Hadoop and Spark as the large size of the data created the havoc in the testing. We will also cluster the imputed tensor and semantically compared the imputation.