A Functional approach for Two Way Dimension Reduction in Time Series

This paper aims to develop an approach to map and represent the multiple time series variables (features) in a compact and efficient dimension reduced manner by leveraging the temporal dependencies within and between the involved variables. We begin by introducing useful notations, discussing the prior arts and problem definition, before proceeding to give the details of Bi-Functional Autoencoders (BFAE). Let us assume that for the $i^{th}$ ($i\in\{1,2,...,N\}$) independent subject, we have $R$ features, that are continuously recorded within a compact time interval $\mathcal{S}\subseteq\mathbf{R}$. In particular, the observed timepoints of the $r^{th}$ feature for the $i^{th}$ subject are given by a $M^{(i,r)}_{s}$-dimensional vector $\mathbf{S}^{(i,r)}=[S^{(i,r)}_{1},...,S^{(i,r)}_{j},...,S^{(i,r)}_{M^{(i,r)}_{s}}]^{T}$, with $M^{(i,r)}_{s}$ representing the number of timepoints in the time series and $S^{(i,r)}_{j}\in\mathcal{S}$ for $i=1,...,n;r=1,...,R;j=1,...,M^{(i,r)}_{s}$. The corresponding time series are denoted as $\mathbf{X}^{(i,r)}=[X^{(i,r)}_{1},...,X^{(i,r)}_{j},...,X^{(i,r)}_{M^{(i,r)}_{s}}]^{T}$. The subscript $M^{(i,r)}_{s}$ reflects the fact that the measuring timestamps may vary across features and subjects. For the applicability of prior arts like PCA and AE, the data needs to be regular time series, therefore, we assume the number of timepoints across all features and subject is $M$. The observed data can be denoted as $\{\mathbf{X}^{(i,1)},...,\mathbf{X}^{(i,R)}\}_{i=1}^{N}$. We then concatenate the $R$ temporal features as $\mathbf{{X}}^{(i)}=[\mathbf{{X}}^{{(i,1)}^{T}},...,\mathbf{{X}}^{{(i,R)}^{T}}]^{T}$. Given samples $\{\mathbf{{X}}^{(i)}\}_{i=1}^{N}$, the problem definition can be stated as an unsupervised nonlinear multi-dimensional functional representation learning problem, formally defined as follows:

where $\mathbf{Z}^{(i)}=\{\mathbf{Z}^{(i,1)},...,\mathbf{Z}^{(i,R^{\prime})}\}_{i=1}^{N}$ is a latent representation of the time series $\mathbf{{X}}^{(i)}$ with $R^{\prime}<R$.

This problem formulation has a few disadvantages. Biases may get introduced when performing data pre-processing for having the same number of timepoints for each time series and these widely used prior arts have their own limitations in solving the dimension reduction problems in time series data when they assume them to be scalar values. These multivariate approaches are unable to capture the temporal dependencies among $\mathbf{{X}}^{(i,r)}$, $r=1,...,R$.

In the paper so far, we have seen the dimension reduction problem from a multivariate time series point of view. Let us look at an alternative problem formulation using functional data analysis (FDA). In FDA, we learn from random functions which are dynamically varying data over a continuum. There is a lot of application for FDA in a wide variety of fields, like time series, sensor data, image, and spatial data [14, 9]. We can use this ability of FDA for dealing with continuous underlying curves $X^{(i,r)}(s),s\in\mathcal{S}$ that is observed at discrete timepoints along a time interval. In functional analysis, the input or output has to be functional curves if not both. For given functional features $\{X^{(i,1)}(s),...,X^{(i,R)}(s),s\in\mathcal{S}\}_{i=1}^{N}$, our problem of finding a latent representation of the data can be given as follows:

where $\mathbf{Z}^{(i)}(t)=\{\mathbf{Z}^{(i,1)}(t),...,\mathbf{Z}^{(i,R^{\prime})}(t),t\in\mathcal{S}\}_{i=1}^{N}$ is a latent representation of the time series $\mathbf{{X}}^{(i)(s)}$ with $R^{\prime}<R$ and each time series is observed at $M^{\prime}$ timepoints ($M^{\prime}\leq M$).

The popularity of FDA in time series application is evident with the boom of recent interest in this area of research. The assumption of smoothness can deviate in FDA as shown in [28], allowing us to directly apply different FDA techniques to time series data. Also, compared to the common sequential models (RNNs, LSTMs) which learn fixed parameters over time, FDA enables us to learn feature effects that change over the interval $\mathcal{S}$ and build more efficient models [25]. In the functional field, we have FPCA [6, 1] which can encode the information in a low-dimensional latent space but the mapping is linear. This can be a limiting factor because, in many real-world applications, the multiple functional features can have complex relations not just among them but across different timepoints as well.

The growing interest in deep learning and FDA has resulted in multiple successive deep functional models [19, 26, 16] for the different analytical tasks. FAE [8] is a deep learning approach for dimension reduction, a generalization of the autoencoder from the vector space to the functional setting. In this network, only the first and last layer can accommodate functions with the help of functional neurons given by [19, 20], and the rest of the network consists of traditional neurons. Because of this, FAE can only give a latent representation of functional data in the scalar form, which is limiting. It is very common for functional features to have complex relations that can’t be embedded in a scalar form. Thus, the scalar representation of functional data is inadequate.

In the next section, we propose a non-linear function-on-function autoencoder that leverages the power of fully connected Neural Networks while keeping the functional identity of the time series.

The main goal of this paper is to learn a mapping $\mathbf{f}(\cdot)$ for the time series data to a low dimensional latent space as shown in Figure 1. To achieve this goal we make use of the idea of an autoencoder and Functional Neural Networks [19, 16, 17]. An autoencoder (AE) is a popular dimension reduction technique which can learn non-linear latent presentation in a $R^{\prime}$-dimensional space $\mathbf{R}^{R^{\prime}}$ from a $R$-dimensional vector-valued input space $\mathbf{R}^{R}$ (where $R^{\prime}<R$) using the help of an encoder. The decoder part of the AE helps to map the learned latent representation back to the $R$-dimensional vector-valued input space $\mathbf{R}^{R}$. The input and output of an AE are the same and we learn by measuring the reconstruction error of the output compared to the input. We want to use this idea but in a functional setting where we have an encoder that is mapping from a $R$-dimensional functional space to a $R^{\prime}$-dimensional functional space and the decoder maps $R^{\prime}$-dimensional functional space back to the $R$-dimensional functional space. We make use of continuous neurons defined in [16, 17] to enable us to achieve this task in a functional setting. The continuous neurons take in functional inputs and produce a functional output, thus leading to the preservation of the functional nature of the data.

The framework for our Bi-Functional Autoencoder (BFAE) Model, as seen in Figure 1, involves development of a continuous mapping from layer to layer by using multiple continuous hidden layers consisting of multiple continuous neurons as defined in [16, 17]. In Figure 1, we have three types of layers: an input layer, continuous hidden layers, and a continuous output layer. The input layer is the first layer that takes in the functional features, the continuous hidden layers consist of multiple continuous neurons, and the continuous output layer provides the reconstructed input functional features using the continuous neurons. The model learns with the help of functional weights that are continuous over time (or some other continuum). These functional weights consists of both univariate and bivariate functions, that is, we learn functions as well as surfaces. The key idea here is that we want to preserve the functional nature of the data throughout the whole reconstruction process across the network. This will give a richer structure to transform the data into low dimensional latent representation while exploiting the domain information and continuity of the functional features. We define the $l^{th}$ continuous hidden layer and its $r^{th}$ continuous neuron as:

where $l=1,2,3,...,L$, $H^{(i,r)}_{(0)}(s)=X^{(i,r)}(s)$, $H^{(i,r)}_{(L)}(s)=\widehat{X^{(i,r)}(s)}$, $b^{(r)}_{(l)}(s)\in\mathcal{L}^{2}(\mathcal{S})$ is the unknown univariate intercept function, $w^{(r,j)}_{(l)}\in\mathcal{L}^{2}(\mathcal{S}\times\mathcal{S})$ is the bivariate parameter function for the $r^{th}$ continuous neuron in the $l^{th}$ continuous hidden layer coming from the $j^{th}$ continuous neuron of the $(l-1)^{th}$ continuous hidden layer and $\sigma(\cdot)$ is a standard activation function. The value of $J$, that can be considered as the number of features in each continuous hidden layer, is fixed for the input layer and the continuous output layer as $R$. We can tune $J$ to any scalar value for the rest of the continuous hidden layers. Similarly, the timepoints at which the output of the continuous neurons is observed is fixed for the input layer and the continuous output layer and can be changed for all the other continuous neurons in the remaining continuous hidden layers.

We will now look into the working of our network. We use the continuous neurons defined by Equation 3, to get from the input layer to a dimension reduced form of the input at an intermediate continuous hidden layer (say, it is the $l^{\prime th}$ layer). This is the encoder part of the BFAE approach called as the functional encoder. We can get the low dimensional representation of the inputs in the $l^{\prime th}$ continuous hidden layers as $\mathbf{Z}^{(i)}(s)=H^{(i)}_{(l^{\prime})}(s)$ where $\mathbf{Z}^{(i)}(s)=\{\mathbf{Z}^{(i,1)}(s),...,\mathbf{Z}^{(i,R^{\prime})}(s),s\in\mathcal{S}\}_{i=1}^{N}$ and $H^{(i)}_{(l^{\prime})}(s)=\{H^{(i,1)}_{(l^{\prime})}(s),...,H^{(i,R^{\prime})}_{(l^{\prime})}(s),s\in\mathcal{S}\}_{i=1}^{N}$, $R^{\prime}<R$ and the function $\mathbf{Z}^{(i)}$ is observed at $M^{\prime}\leq M$ timepoints. We move from this latent representation of the input towards the continuous output layer using the continuous neurons in the continuous hidden layers to get the reconstructed functional input values. This part of our approach is called the functional decoder. Therefore, we have developed a mapping from a $R$-dimensional functional space to $R^{\prime}$-dimensional functional space using the functional encoder, and the functional decoder maps $R^{\prime}$-dimensional functional space back to the $R$-dimensional functional space. This latent representation $\mathbf{Z}^{(i)}(s)$ is very flexible in nature, where we can set the number of functional features as $R^{\prime}$ and the number of timepoints observed as $M^{\prime}$, according to the task in hand.

For simplicity, let us assume that Figure 1 has $R$ functional features in the input layer, $L$ continuous hidden layers, each with $J$ incoming and $R$ outgoing connections. BFAE is similar to an autoencoder in nature, but for functional data, we have to define the number of continuous hidden layers, the number of continuous neurons in each of the continuous hidden layers, and the number of timepoints at which these functions are observed for each continuous neuron. The activation function used in the continuous neurons is ReLU, tanh, sigmoid/logistic, or linear (continuous output layer). The forward propagation of the network is straightforward with the help of Equation 3. In the backpropagation phase, we learn using the functional surfaces (bivariate) and the functional intercepts (univariate) directly. We need the functional gradients to learn the functional parameters of the BFAE network. These functional gradients measure the change in a functional to a change in a function on which the functional depends.

We use the functional gradients given below to optimize our network’s functional parameters. The optimization approach we use is the traditional gradient descent. The necessary assumptions and mathematical tools are discussed in [19, 26, 12]. We use Fréchet derivatives, from the calculus of variation, for computing the functional gradients. The partial derivatives needed for the backpropagation are as follows:

where $l=1,2,3,...,L$, $H^{(i,r)}_{(0)}(s)=X^{(i,r)}(s)$, $H^{(i,r)}_{(L)}(s)=\widehat{X^{(i,r)}(s)}$ and $\sigma^{\prime}(\cdot)$ represents the first derivative of $\sigma(\cdot)$. In the backpropagation phase, we pass through the network backward from the continuous output layer to the input layer and calculate the partial derivatives of the continuous neurons with respect to the bivariate functional weight and the functional intercept. The loss function that we minimize is $\mathcal{L}(\theta)$, where $\theta$ is the collection of all functional parameters in the BFAE network. The loss function is given as follows:

The number of continuous hidden layers, the number of continuous neurons in each of the continuous hidden layers, and the number of observed timepoints for each continuous neuron can be considered as hyperparameters. We can also adjust our approach to accommodate irregular functional data or scalar variables. All these offer a lot of flexibility for our approach to work with different kinds of problems and adjust things according to the downstream process for different analytical tasks like prediction, classification, clustering, forecasting, and more.

Some of the current dimension reduction methods can be represented as a special case of BFAE. As we know, a single layer autoencoder (AE) with a linear activation function is analogous to PCA. In the same way, we can adjust the parameters of our model to act like FPCA. The FPCA is a special case of BFAE when our model has the following specifications: linear activation functions, a single continuous hidden layer, and the functional weights are constrained to be orthonormal. The Functional Autoencoder (FAE) [8] paper discusses how AE is a special case of FAE as FAE replaces the scalar weights and inner products of the AE with functional weights and inner products. In the same manner, FAE is a special case of our approach, BFAE. If we specify the functional encoder in our model to learn a scalar ($M^{\prime}=1$) latent representation of the functional features, our approach essentially acts like FAE. Therefore, our approach is a generalization of most of the existing methods like PCA, FPCA, AE, and FAE.

In the first experiment, we consider an individual predictor function (i.e. $R=1$) which is observed on a dense and regular grid. We generate $N=100,1000$ iid random curves $\left\{X^{(1)}(t),\cdots,X^{(n)}(t)\right\}$ from a Gaussian process with mean $0$ and covariance given as follows:

where this is the Matérn covariance function and $K_{\nu}$ is the modified Bessel function of the second kind. The value of $\rho=0.5,$ $\nu=5/2$, and $\sigma^{2}=1.$ These curves are realized at equally-spaced timepoint on the interval $[0,1]$ where $M=50,250$ timepoints. We have specifically chosen a case where $N<M$ to check if it causes any issues with any of the approaches.

After generating the curves using a Gaussian process, we add some random noise ($\epsilon_{i}(t)$) to them. The results are averaged over 100 simulations, and we divide the 100 samples into 80 for training and 20 for testing. The same ratio split is used for the case of 1000 samples. We report the Root Mean Square Error (RMSE) for the reconstruction of the functional features given by Equation (7) and in the case of the scalar approaches, we just report the standard RMSE ignoring the temporal aspect. For our approach in the first experiment, we consider multiple combinations of continuous hidden layers ($L$= 1, 3) and only one continuous neuron ($J=1$). The grid (s) value for latent representation in the continuous neurons is either $M$ or $M^{\prime}=M/5$. We keep the network of AE similar to BFAE, and for PCA and FPCA we follow the $99\%$ variance of explained rule. Note that we can select $M^{\prime}$ to be any scalar value in the range of ($1,M$) but we select this particular value ($M^{\prime}=M/5$) to show the flexibility of our approach.

We observe from Table I that PCA and AE perform similarly on the temporal curves. There is some gain in performance for FPCA but our approach performs the best irrespective of the sample size ($N$) and the number of timepoints ($M$). The difference in performance increases when the sample size and the number of timepoints increases, as our approach is a deep learning model, it learns better with more information. When we set the number of observed timepoints as $M^{\prime}$, our approach loses some performance because we are restricting the true data to a smaller latent space. But, even under this scenario of representing the true curves using 10 or 50 timepoints rather than 50 or 250, we are still performing better than all the other approaches.

Figure 2 shows us the reconstructed curves for a random sample from different approaches for the simulation setting of $N=1000$, $R=1$, and $M=50$. We see a similar story compared to Table I where PCA is performing poorly, and the shape is very different from the truth. PCA is a linear approach and can’t capture the temporal information in the data. AE and FPCA perform better than PCA by getting the general shape correct but are still far off from the truth. Our approach, BFAE, irrespective of the number of timepoints observed, is able to reconstruct the truth well and captures the shape of the curve correctly over the whole time interval.

We extend the first experiment to include more functional features. In the second experiment, we increase the number of features to $R=10$. We follow the same procedure to generate these 10 curves for the same set of values of $N$ and $M$. The results discussed below are an average of 100 simulations after adding random noise and using the same splits for training and testing. For our approach, we consider multiple combinations of continuous hidden layers ($L$= $1$, 3), continuous neurons ($J$= 1, 2, 4) and we set $R^{\prime}=4$. The grid (s) value for latent representation in the continuous neurons is either $M$ or $M^{\prime}=M/5$.

Table II shows that the increase in the functional features results in a deterioration in performance for PCA. While the other two approaches perform better than PCA, our approach performs the best irrespective of the sample size ($N$) and the number of timepoints ($M$). The difference in performance again increases when the sample size and number of timepoints increases. We observe that the errors have increased marginally compared to the previous table as we have added more features. The behavior of decreasing the number of timepoints to $M^{\prime}$ is similar as well.

We consider two real data sets. The first is the speech recognition data from TIMIT (available at http://statweb.stanford.edu/tibs/Ele-mStatLearn/) where we have speech signals for different phonemes. Audio information is available abundantly, but they are measured at a very high frequency requiring high storage capacity. A low dimensional latent representation of this data would be very useful in the industry. The other example deals with the relation between electricity demand and temperature in the city of Adelaide, Australia. This data is important because of the high costs related to the storage of electricity and understanding the effect of temperature on demand can lead to information based operational steps.

For speech recognition data, we set up the experiment similar to [10, 16], where we have voice signals ($R$=1) for two phonemes (response) transcribed as follows: “aa” as the vowel in “dark” and “ao” as the first vowel in “water”. We build a classifier model with the help of FLM after applying dimension reduction to the voice signals. Figure 3 shows that it is difficult to separate the two groups where we see the two phoneme curves computed by a log-periodogram of 150 ($M=150$). Each phonemes information is recorded over 150 points, reduction of such data is beneficial, especially with the potential to expand this low dimension representation to much higher voice signals. We have 800 functional samples, we split these samples into 640 samples for training and 160 samples for testing.

The Adelaide data (available in $fds$ package in R-software) has the temperature and electricity demand records from 7/6/1997 to 3/31/2007 (508 weeks) measured daily at half-hourly rate ($M=48$). We consider each day of the week as a feature ($R$=7) and try to map the relation between temperature and demand using Functional Linear model (FLM) [4, 9, 14]. Before we use FDA to model this relation, we use BFAE to represent the temperature curves in a compact low dimensional form. We split this data into 400 samples for training and 108 samples for testing. Figure 4 shows the half-hourly temperature and electricity demand (Megawatts) for the whole 508 weeks. We can see from the image that the temperature curve for each day of the week follows a similar pattern and learning that pattern can lead to a reduction in dimension with minimum loss of information.

We can observe from Table III that PCA again is having difficulties in reducing the information for both the real data sets. While AE and FPCA are performing better than PCA, BFAE performs the best. For phonemes data, we even reduced the points to $M^{\prime}=30$ and still BFAE performs better than other approaches. The difference between BFAE with timepoints as $M$ and $M^{\prime}$ is negligible, indicating that the low dimension representation is very rich in information. While for Adelaide data, we reduce the number of features to 4 ($R^{\prime}=4$) and the timepoints to $M^{\prime}=12$ and produce the best results using our approach with timepoints as $M$. Our intuition to reduce the temperature feature is certainly validated as seen from Figure 4 and Table III.

Modeling results can be seen in Table IV where we build a classifier for the phonemes data, to predict if the voice signal is for “aa” or “ao”. For the Adelaide data, we map the temperature information for the 7 days of the week to predict the electricity demand for the 7 days of that same week. Both these models are built with the help of the original data and the low dimension latent representation learned using BFAE. We ignore the other approaches as the goal here is to demonstrate that our approach gives modeling results at least as good as the original data by capturing the signals. We can see from Table IV that our approach not only does a good job of modeling but also performs as well as the original data and has a lower tendency to over fit as it contains less noise. The error increases in both cases as we reduce the timepoints to $M^{\prime}$ but the performance is still competitive.

Overall, our approach does the best job in reducing the information into a low dimensional latent space while still maintaining a competitive performance at modeling different tasks compared to the original data. We are able to represent the information in a compact data rich manner and are also able to reconstruct it back to the original scale without much loss of information.