From Time Series to Networks in \proglangR with
the \pkgts2net Package

Time series comparison is a fundamental task in data analysis. This comparison can estimated using similarity or distance functions (Deza and Deza, 2009; Esling and Agon, 2012; Wang et al., 2013). Some research domains use the concept of similarity, while others prefer distance functions. Similarity and distance are inverse concepts that have the same usefulness. Here in this paper, for simplicity, I use the distance concept. However, it is important to remember that a similarity function can be generally transformed into distance functions and vice-versa. The distance or metric $d\mathrel{\mathop{\mathchar 58\relax}}\mathcal{X}\times\mathcal{X}\to\mathbb{R}$ is a function that for all time series $X,Y,Z\in\mathcal{X}$, the following properties hold:

• •

  non-negativity: $d(X,Y)\geq 0$,

• •

  symmetry: $d(X,Y)=d(Y,X)$,

• •

  reflexity: $d(X,Y)=0\iff X=Y$,

• •

  triangular difference: $d(X,Y)\leq d(X,Z)+d(Z,Y)$.

Some of the time series comparison measures used in \pkgts2net are not true metrics because they do not follow one or more required properties above (Deza and Deza, 2009). However, for simplicity reasons, this paper uses the term “distance” to generally refer to any of these functions. The rest of Sec. 2.1 reviews the main distance functions implemented in the \pkgts2net package.

A network or a graph (Bondy and Murty, 2008; Barabási and Pósfai, 2016), $G(V,E)$ is composed by a set of $n$ nodes (also called vertices) $V$ and a set of $m$ links (also called edges) $E$ connecting two nodes. The most common representation for a network is the adjacency matrix $A$ whose entry $a_{ij}=1$ if the nodes $i$ and $j$ are connected or $a_{ij}=0$ otherwise. A weighted network can also be represented by $A$. In this case, each entry $a_{ij}=w_{ij}\in\mathbb{R}_{\neq 0}$ corresponds to the weight $w_{ij}$ between nodes $i$ and $j$, or $a_{ij}$ = 0 if there is no link.

A set of time series $\mathcal{X}$ can be transformed into a network by calculating the distance between all pairs of time series. These distances are represented by the distance matrix $D$ whose position $d_{ij}$ stores the distance between the time series $i$ and $j$. The network construction process consists of transform the distance matrix $D$ into an adjacency matrix $A$ where each node represents a time series. The most common building approaches are (Silva and Zhao, 2018; Qiao et al., 2018):

• •

  $k$-NN network: In a $k$ nearest neighbors network ($k$-NN), each node (representing a time series) is connected to the $k$ other ones with the shortest distances. Thus, the construction process requires finding for each row $i$ of $D$ the $k$ other elements with shortest distances, where $k\in\mathbb{W}$ is a parameter defined by the user. Approximation algorithms speed up the construction process (Arya et al., 1998; Chen et al., 2009).

• •

  $\epsilon$-NN network: In a $\epsilon$ nearest neighbors network ($\epsilon$-NN), (Fig. 4-b), each node is connected to all the other ones if their distances are shorter than $\epsilon\in\mathbb{R}_{+}$, defined by the user. An advantage of the $\epsilon$- over $\epsilon$-NN networks is that the $\epsilon$ permits a finer in the threshold choice. Another advantage is the fact that $\epsilon$ can restrict the connections to only strong ones, while the $k$-NN will always connect nodes regardless if the $k$ shortest distances are small or high. Similar to $k$-NN networks, approximations can speed up the construction process (Arya et al., 1998; Chen et al., 2009).

• •

  Weighted network: Instead of limiting the number of links, a weighted network (Fig. 4-c) can be constructed from $D$ by simply connecting all pairs of nodes and using their distances to set weights. In general, link weights have the opposite definition of a distance, which means that the shorter the distance between two time series, the stronger the link connecting the respective nodes. Considering that $D$ is in the interval [0,1], the weighted adjacency matrix can be defined as $A$ = 1 - $D$. Conversely, if $D\in[0,+\infty]$, $D$ can be be first normalized between [0,1] ($D_{\text{norm}}$) and then consider $A$ = 1 - $D_{\text{norm}}$.

• •

  Networks with significant links: In many applications, it is important to minimize the number of spurious links in the network. One way of minimizing these spurious links is by considering statistical tests for the distance functions. In other words, this network construction process connects all pairs of nodes only if their distance is statistically significant. For example, the significance of the Pearson correlation coefficient can be tested using the z-transformation (Eq. 9) described in Sec. 2.1. This network construction process can be seen as a special case of $\epsilon$-NN network where $\epsilon$ is a threshold for significant links. It is important to remember that this building approach minimizes the number of spurious links but does not guarantee that they will not appear.

Networks can be static or temporal (Holme and Saramäki, 2012; Kivelä et al., 2014). A static network $G(V,E)$ uses the whole historical data in a time series to construct a single network. Temporal networks divide the historical period into time windows that are used to construct layers in a temporal (multilayered) network. A temporal network is a sequence of static networks ordered in time. Each layer can be seen as a snapshot of the time window. The temporal network construction process requires the computation of a set of distance matrices $\mathcal{D}=\{D_{1},\cdots,D_{W}\}$, one for each of the $W$ windows. A time window is defined by the user and has a length $\Delta t$ and a step $s$ that defines the number of values that are included and excluded in each window. The temporal network construction consists of transforming each distance matrix $D\in\mathcal{D}$ into a network layer using the same process used for a static network previously described in this section. One of the main advantages of this approach is to make possible the temporal analysis (evolution) of the modeled system.

The network construction process has a combinatorial nature. Considering a set of $n$ time series, the whole construction process requires ${n\choose 2}$ distance calculations in a static network or $W\cdot{n\choose 2}$ in a temporal network, where $W$ is the number of layers. The transformation from the distance matrix to the adjacency matrix requires searches for pairs of shortest nodes, a task that is much faster than the distance calculations. For this reason, the conversion step ($D\to A$) can be neglected in the total time complexity is $O({n\choose 2})$. This means that the number of distance computations can increase very fast for very large time series data sets, making the network modeling approach infeasible.

This section presents methods to transform a single time series into a network (Zou et al., 2019; Silva et al., 2021).

The \pkgts2net package relies on some other important packages. One of the most important ones is \pkgigraph, which is used to model networks (Csardi and Nepusz, 2006). The \pkgigraph package provides many different graph and network analysis tools, which can be used to extract information from the modeled time series. The \pkgigraph objects can be easily transformed or exported to other formats such as graphml, gml and edge list. Besides, other network analysis packages commonly accept these objects. Examples are the \pkgsna (Butts, 2020), \pkgtsna (Bender-deMoll and Morris, 2021), and \pkgmuxViz (De Domenico et al., 2014).

The network construction process normally requires a high number of calculations. For this reason, the most important functions in \pkgts2net are implemented in parallel using multi-core programming. The parallelization was implemented using the \codemclapply() function from the \pkgparallel, which comes with \proglangR. This implementation uses forks to create multiple processes, an operation that is not available in Windows systems. It means that, on Windows, all functions in \pkgts2net with num_cores can only run with one core.

Besides the multicore parallelization, the \pkgts2net was also constructed to be run in supercomputers and computer clusters. The idea is to divide the distance calculations into jobs that could be calculated in parallel and then merged together to construct the distance matrix. The \codets_dist_dir() and \codets_dist_dir() functions make this process easier.

Functions to calculate the distance matrix:

• •
  \code

  ts_dist(): Calculates all pairs of distances and returns a distance matrix $D$. This function receives a list of time series which can either have or not have the same length. It runs serial or in parallel (except in Windows) using \codemclapply() from \pkgparallel package. This function receives a distance function as a parameter. \pkgTs2net provides some distance functions (functions starting with \codetsdist_…()). Other distances can be easily implemented or adapted to be used within this function.

• •
  \code

  ts_dist_part(): Calculate distances between pairs of time series (similarly to \codets_dist()) in part of a list. This function is particularly useful to run in parallel as jobs in supercomputers or computer clusters (HPC). Instead of having a single job calculating all distances (with \codets_dist()), this function permits the user to divide the calculation into different jobs. The result is a data frame with the distances $d_{ij}$ and indexes $i$ and $j$. To construct the distance matrix $D$, the user can merge all the results using the functions \codedist_parts_merge() or \codedist_parts_file_merge().

• •
  \code

  ts_dist_part_file(): This function works similar as \codets_dist_part(). The difference is that it reads the time series from RDS files in a directory. The advantage of this approach is that it does not load all the time series in memory but reads them only when necessary. This means that this function requires much less memory and should be preferred when memory consumption is a concern, e.g., huge data set or very long time series. The disadvantage of this approach is that it requires a high number of file read operations which considerably takes more time during the calculations. The user can merge all the results using the functions \codedist_parts_merge() or \codedist_parts_file_merge().

List of distance functions:

• •
  \code

  tsdist_cor(): Absolute (Eq. 2), positive (Eq. 3) or negative (Eq. 4) correlation distances. This function optionally test the significance for the correlation using the Fisher z-transformation (Eq. 8) (Asuero et al., 2006).

• •
  \code

  tsdist_ccf(): Absolute (Eq. 5), positive (Eq. 6), or negative (Eq. 7) cross-correlation distances.

• •
  \code

  tsdist_dtw(): Dynamic Time warping (DTW) distance (Sec. 2.1.3). This function is a wrapper for the \codedtw() function from the \pkgdtw package (Giorgino, 2009).

• •
  \code

  tsdist_nmi(): Normalized mutual information distance (Eq. 17). The implementation uses functions from the \pkginfotheo package (Meyer, 2014).

• •
  \code

  tsdist_voi(): Variation of information distance (Eq. 18). The implementation uses functions from the \pkginfotheo package (Meyer, 2014).

• •
  \code

  tsdist_mic(): Maximal information coefficient (MIC) distance. This function transforms the \codemine() function from the \pkgminerva package (Albanese et al., 2012) into a distance.

• •
  \code

  tsdist_es(): The events synchronization distance (Eq. 22) proposed by Quian Quiroga et al. (2002). This function can also use the modification proposed by Boers et al. (2019) ($\tau_{\text{max}}$). This function optionally test the significance for the events synchronization by shuffling the events in the time series.

• •
  \code

  tsdist_vr(): Van Rossum distance. This function uses the \codefmetric() function from the \pkgmmpp package (Hino et al., 2015). It optionally test the significance for the events synchronization by shuffling the events in the time series.

Other distances can be easily applied in the \pkgts2net. The only two requirements are that the distance function should accept two arrays and it must return a distance. The packages \pkgTSclust (Montero and Vilar, 2014) and \pkgTSdist (Mori et al., 2016) provide many other distance functions that can be used to create network with \pkgts2net. For example, the euclidean distance (\codediss.EUCL() from \pkgTSclust) can be used to construct the distance matrix: {CodeChunk} {CodeInput} R> install.packages("TSclust") R> library(TSclust) R> D <- ts_dist(ts_list, dist_func=diss.EUCL) It is important to remember that the parameter \codedist_func receives a function, so no parenthesis is necessary (\codedist_func=diss.EUCL() is an error).

Another possibility is to implement specific functions and apply them. For example, the following code snippet implements a distance function that compares the mean values of each time series and applies it to calculate the distance matrix $D$: {CodeChunk} {CodeInput} R> ts_dist_mean <- function(ts1, ts2) abs(mean(ts1) - mean(ts2)) R> D <- ts_dist(ts_list, dist_func=ts_dist_mean)

Multiple time series into a network:

• •
  \code

  net_knn(): Constructs a $k$-NN network from a distance matrix $D$.

• •
  \code

  net_knn_approx(): Creates an approximated $k$-NN network from a distance matrix $D$. This implementation may omit some actual nearest neighbors, but it is considerably faster than \codenet_knn(). It uses the \codekNN() function from the \pkgdbscan package.

• •
  \code

  net_enn(): Constructs a $\epsilon$-NN network from a distance matrix $D$.

• •
  \code

  net_enn_approx(): Creates an approximated $\epsilon$-NN network from a distance matrix $D$. This implementation may omit some actual nearest neighbors, but it is considerably faster than \codenet_enn(). It uses the \codefrNN() function from the \pkgdbscan package.

• •
  \code

  net_weighted(): Creates a full weighted network from a distance matrix $D$. In this case, $D$ should be in the interval $[0,1]$ (see \codedist_matrix_normalize()). The link weights in the resulting network is $w_{ij}=1-d_{ij}$.

• •
  \code

  net_significant_links(): Construct a network with significant links.

A single time series into a network:

• •
  \code

  ts_to_windows(): Extract time windows from a time series. This function can be used to construct a network from time windows. It returns a list of windows that can be used as input for \codets_dist() to construct a distance matrix $D$.

• •
  \code

  tsnet_qn(): Creates transition (quantile) networks.

• •
  \code

  tsnet_vg(): Creates natural and horizontal visibility graphs. This function can be run in parallel.

• •
  \code

  tsnet_rn(): Creates a recurrence network from a time series. This function uses the \coderqa() function from the \pkgnonlinearTseries package to obtain the recurrence matrix.

• •
  \code

  dist_matrix_normalize(): normalizes (min-max) a distance matrix $D$ between zero and one.

• •
  \code

  dist_percentile(): Returns the distance value that corresponds to the desired percentile. This function is useful when the user wants to generate networks with different distance functions but with the same link density.

• •
  \code

  ts_dist_dirs_merge(): \codets_dist_dir() and \codets_dist_dir() calculate parts of the distance matrix $D$. This function merges the results and constructs a distance matrix $D$.

• •
  \code

  tsdist_file_parts_merge(): Merges the distance parts as \codets_dist_dirs_merge(), but reading them from files.

• •
  \code

  events_from_ts(): Function to extract events from a time series. For example, it can return the desired percentile of highest or lowest values. This function is particularly useful for event-based distance functions.

• •
  \code

  dataset_sincos_generate(): Generates a data set with sine and cosine time series with arbitrary noise.

• •
  \code

  random_ets(): Creates random event time series with uniform probability.

• •

  us_cities_temperature_df: A data set containing the temperature in 27 US cities between 2012 and 2017. This data set was adapted from the original (Beniaguev, 2017) for simplicity reasons. Different from the original, this data set comprehends only cities in the US, grouped by month (mean value) temperature and without days with missing data. The data set is stored in a data frame with 61 rows and 28 variables where each column (except date), corresponds to the mean temperature values (°C) for each city.

• •

  us_cities_temperature_list: Similar to us_cities_temperature_df, but in a list, which is the required format for distance calculations (e.g. \codets_dist()).

The \pkgts2net is open-source and free to use according to the MIT License (https://github.com/lnferreira/ts2net/blob/main/LICENSE). All the code is available on Github at https://github.com/lnferreira/ts2net.

This section presents a network modeling example using real data. The data set is composed of monthly temperature data from 27 cities in the US between 2012 and 2017, illustrated in Fig. 7-a. The version presented here was extracted and adapted from Kaggle (Beniaguev, 2017) for simplicity reasons by taking the monthly averages and keeping only US cities.

The following code calculates the distance matrix (Fig. 7-b) using DTW, finds the $\epsilon$ that corresponds to the 30% of the shortest distances, constructs a $\epsilon$-NN network (Fig. 7-c), and find communities using the algorithm proposed by Girvan and Newman (2002):

R> D <- ts_dist(us_cities_temperature_list, dist_func = tsdist_dtw) R> eps <- dist_percentile(D, percentile = 0.3) R> net <- net_enn(D, eps) R> cluster_edge_betweenness(net) {CodeOutput} IGRAPH clustering edge betweenness, groups: 2, mod: 0.33 + groups: $`1`[1]"Portland""Seattle""Albuquerque""Denver"[5]"KansasCity""Minneapolis""SaintLouis""Chicago"[9]"Nashville""Indianapolis""Detroit""Charlotte"[13]"Pittsburgh""Philadelphia""NewYork""Boston"\par$‘2‘ [1] "San Francisco" "Los Angeles" "San Diego" [4] "Las Vegas" "Phoenix" "San Antonio" [7] "Dallas" "Houston" "Atlanta" + … omitted several groups/vertices

Fig. 7-c illustrates the resulting network where nodes are conveniently placed according to their geographical position. Node colors represent the two communities that cluster together colder (blue nodes) and hotter (red nodes) cities.

This section presents how a time series extracted from real data can be modeled as networks. The time series with the monthly atmospheric concentration of carbon dioxide (ppm) measured in the Mauna Loa Observatory between 1959 and 1997 (Keeling et al., 2005), illustrated in Fig. 8a. This data set was mainly chosen by its simplicity and accessibility (available with the base installation of R).

The following code generates a time-window network (net_w) with width 12 and one-value step, a visibility graph (net_vg), a recurrence network (net_rn) with $\epsilon=5$, and a transition network (net_qn) using 100 equally spaced bins using the \pkgts2net package.

R> net_w = as.vector(co2) |> ts_to_windows(width = 12, by = 1) |> ts_dist(cor_type = "+") |> net_enn(eps = 0.25) R> net_vg = as.vector(co2) |> tsnet_vg() R> net_rn = as.vector(co2) |> tsnet_rn(radius = 5) R> net_qn = as.vector(co2) |> tsnet_qn(breaks = 100)

Figures 8b-e illustrate the networks resulting from the transformation of the CO₂ time series. From the topology of the resulting networks, it is possible to extract useful information. For example, the time-window network (Fig. 8b) shows 12 groups representing each month. The visibility graph (Fig. 8c) presents many small-degree nodes that represent the values on valleys while a few hubs connect peaks values of CO₂ after 1990 that connect to many other nodes. These hubs appear due to the increasing trend of CO₂ contraction and some slightly higher peaks during the seasons after 1980. The recurrence (Fig. 8d) and transition (Fig. 8e) networks present some local connections caused by the seasons, and a line shape as a result the increasing trend in the CO₂ concentration over the years.