1 Introduction
The analysis of functional data is becoming more and more important in many areas of application such as medicine, economics, or geology (cf. Ullah and Finch (2013), Wang et al. (2016)), where this type of data occurs naturally. In industry, functional data are often a byproduct of continuous monitoring of production processes, yielding great potential for data mining tasks. A common type of functional data are time series, as time series can often be considered as discretized functions over time.
Many researchers publish software implementations of their algorithms, therefore simplifying the access to already established methods. Even though such a readily available, broad range of methods to choose from is desirable in general, it also makes it harder for nonexpert users to decide which method to apply to a problem at hand and to figure out how to optimize their performance. As a result, there is an increasing demand for automated model selection and parameter tuning.
Furthermore, the functionality of available pipeline steps ranges from simple data structures for functional data, to feature extraction methods and packages offering direct modeling procedures for regression and classification. Users are again faced with a multiplicity of software implementations to choose from and, in many instances, combining several implementations may be required. This can be difficult and timeconsuming, since the various implementations utilize a multiplicity of different workflows which the user needs to become familiar with and synchronize in order to correctly carry out the desired analysis.
There is a wide variety of packages for functional data analysis in R R Core Team (2014) available that provide functionality for analyzing functional data. Examples range from the fda (Ramsay et al., 2018) package which includes object types for functional data and allows for smoothing and simple regression, to, e.g., boosted additive regression models for functional data in FDboost (Brockhaus and Ruegamer, 2018). For an extensive overview, see the CRAN task view (Scheipl, 2018).
Many of those packages are designed to provide algorithmic solutions for one specific problem, and each of them requires the user to become familiar with its user interface. Some of the packages, however, such as fda.usc (FebreroBande and Oviedo de la Fuente, 2012) or refund (Goldsmith et al., 2018) are not designed for only one specific analysis task, but combine several approaches. Nevertheless, these packages do not offer unified frameworks or consistent user interfaces for their various methods, and most of the packages can still only be applied separately.
A crucial advantage of providing several algorithms in one package with a unified and principled user interface is that it becomes much easier to compare the provided methods with the intention to find the best solution for a problem at hand. But to determine the best alternative, one still has to be able to compare the methods at their best performance on the considered data, which requires hyperparameter search and, more preferably, efficient tuning methods.
While the different underlying packages are often difficult and sometimes even impossible to extend to new methods, custom implementations and extensions can be easily included in the accompanying software.
We want to stress that the focus of this paper does not lie in proposing new algorithms for functional data analysis. Its added value lies in a large comparison of algorithms while providing a unified and easily accessible interface for combining statistical methods for functional data with the broad range of functions provided by mlr, most importantly benchmarking and tuning. Additionally, the often overlooked possibility of extracting nonfunctional features from functional data is integrated, which enables the user to apply classical machine learning algorithms such as support vector machines (Cortes and Vapnik, 1995) to functional data problems.
In a benchmark study similar to Bagnall et al. (2017a) and Fawaz et al. (2019), we explore the performance of implemented methods, and try to answer the following questions:

Can functional data problems be solved with classical machine learning methods ignoring the functional structure of the data as well as with more elaborate methods designed for this type of data? Bagnall et al. (2017a) measure the performance of some nonfunctionaldataspecific algorithms such as the rotation forest (Rodriguez et al., 2006), but this does not yield a complete picture.

Guidance on the wide range of available algorithms is often hard to obtain. We aim to make some recommendations in order to simplify the choice of learning algorithm.

Do statistical methods explicitly tailored to the analysis of functional data (e.g. FDboost, Brockhaus and Ruegamer, 2018) perform well on classical time series tasks? No benchmark results for these methods, which provide interpretable results, are currently available.

Many methods that represent functional data in a nonfunctional domain have been proposed and are also often applied in practice. Examples for this include either hand crafted features (cf. Stachl and Bühner, 2015), summary statistics (Hyndman et al., 2018), or generally applicable methods such as wavelet decomposition (Aldrich, 2013).

Hyperparameter optimization is a very important step in many machine learning applications. In our benchmark, we aim to quantify the impact of hyperparameter optimization for a set of given algorithms on several data sets.
Contributions
As contributions of this paper, we aim to answer the questions posed above. Additionally, we provide a toolbox for the analysis of functional data. It implements several methods for feature extraction and directly modeling functional data, including a thorough benchmark of those algorithms. This toolbox also allows for full or partial replication of the conducted benchmark comparison.
2 Related Work
In the remainder of the paper, we focus on comparing algorithms from the functional data analysis and the machine learning domain. Functional data analysis traditionally values interpretable results and valid statistical inference over prediction quality. Therefore functional data algorithms are often not compared with respect to their predictive performance in literature. We aim to close this gap. On the other hand, machine learning algorithms often do not yield interpretable results. While we consider both aspects to be important, we want to focus on predictive performance in this paper.
2.1 Feature extraction and classical machine learning methods
In this work, we differentiate between machine learning algorithms that can directly be applied to functional data, and algorithms intended for scalar features, which we call classical machine learning methods.
A popular approach when dealing with functional data is to reduce the problem to a nonfunctional task by extracting relevant nonfunctional features (Ullah and Finch, 2013). Applying classical machine learning methods after extracting meaningful features can then lead to competitive results (cf. Goldsmith and Scheipl, 2014, e.g.)
or at least provide baselines, which are in general not covered by functional data frameworks. In our framework, such functionality is easily available by combining feature extraction, e.g., based on extracting heuristic properties
(cf. tsfeatures; Hyndman et al., 2018) or wavelet coefficients (Mallat, 1989; Aldrich, 2013) and analyzing these derived scalar features with classical machine learning tools provided by mlr.Based on some existing functionality of the listed packages, we adapt different feature extraction methods. Along with different algorithms already proposed in literature, we propose two new custom methods, DTWKernel and MultiResFeatures:

tsfeatures (Hyndman et al., 2018) extracts scalar features, such as autocorrelation functions, entropy and other heuristics from a time series.

fourier
transforms data from the time domain into the frequency domain using the
fast fourier transform
(Brigham and Morrow, 1967). Extracted features are either phase or amplitude coefficients. 
bsignal BSpline representations from package FDboost Brockhaus and Ruegamer (2018)
are used as feature extractors. Given the knots vector and effective degree of freedom, we extract the design matrix for the functional data using
mboost. 
wavelets (Aldrich, 2013) applies a discrete wavelet transform to time series or functional data, e.g., with Haar or Daubechies wavelets. The extracted features are wavelet coefficients at several resolution levels.

PCA
projects the data on their principal component vectors. Only a subset of the principal component scores representing a given proportion of signal variance is retained.

DTWKernel computes the dynamic time warping distances of functional or time series data to (a set of) reference data. We implement dynamic time warping (DTW) based feature extraction. This method computes the dynamic time warping distance of each observed function to a (userspecified) set of reference curves. The distances of each observation to the reference curves is then extracted as a vectorvalued feature. The reference curves can either be supplied by the user, e.g., they could be several typical functions for the respective classes, or they can be obtained from the training data. In order to compute dynamic time warping distances, we use a fast dynamic time warping (Rakthanmanon et al., 2012) implementation from package rucrdtw (BoerschSupan, 2016).

MultiResFeatures extracts features, such as the mean at different levels of resolution (zoomin steps). Inspired by the image pyramid and wavelet methods, we implement a feature extraction method, multiresolution feature extraction where we extract features like mean and variance computed over specified windows of varying widths. Starting from the full sequence, the sequence is repeatedly divided into smaller pieces, where at each resolution level, a scalar value is extracted. All extracted features are concatenated to form the final feature vector.
2.2 Methods for functional data
Without feature extraction, direct functional data modeling (both classification and regression) methods incorporated in our package span two families: The first family of semiparametric approaches includes FGAM (Goldsmith et al., 2018), FDboost (Brockhaus and Ruegamer, 2018), and the functional generalized linear model (FGLM; Ramsay, 2006)
, which are all structured additive models. Those methods use (tensor product) spline basis functions or functional principal components (FPCs)
(Srivastava and Klassen, 2016)to represent effects fitted in a generalized additive model. While FGAM and FGLM use the iterated weighted least square (IWLS) method to generate maximum likelihood estimates, FDboost uses componentwise gradient boosting
(Hothorn et al., 2010) to optimize the parameters. Additionally, the estimated parameters can be penalized. A general formula for this family of methods is , whererepresents a functional of the conditional response distribution (e.g., an expectation or a quantile),
is a vector of (functional) covariates and are partial additive effects of subsets of in basis function representation, cf. Greven and Scheipl (2017) for a general introduction.The second family of methods are nonparametric methods as introduced in Ferraty and Vieu (2006), e.g., based on (semi)metrics which quantify local or global differences or distances across curves. For example, the distance between two instances could be defined by the distance of two curves . Kernel functions are used to average over the training instances and weigh their respective contributions based on the value of their distance semimetric to the predicted instance. Functional nearest neighbors algorithms can also be defined based on such semimetrics. Implementations can be found in packages fda.usc (FebreroBande and Oviedo de la Fuente, 2012) and classiFunc (Maierhofer and Pfisterer, 2018).
2.3 Toolboxes for functional data analysis
The package fda (Ramsay et al., 2018) contains several object types for functional data and allows for smoothing and regression for functional data. Analogously, the Rpackage fda.usc (FebreroBande and Oviedo de la Fuente, 2012) contains several classification algorithms that can be used with functional data. In Python, scikitfda (Grupo de Aprendizaje Automatico  Universidad Autonoma de Madrid, 2019) offers both representation of and (pre)processing methods for functional data, but only a very small set of machine learning methods for classification or regression problems is implemented at the time of writing.
As a byproduct of the TimeSeries Classification Bakeoff Bagnall et al. (2017a), a wide variety of algorithms were implemented and made available. But this implementation emphasizes the benchmark over providing a data analysis toolbox for users, and is therefore not easily usable for inexperienced users.
2.4 Benchmarks
The recently published benchmark analysis TimeSeries Classification Bakeoff by Bagnall et al. (2017a) provides an overview of the performance of stateoftheart algorithms for time series classification. They reimplement (in Java) and compare 18 algorithms designed especially for time series classification on 85 benchmark time series data sets from Bagnall et al. (2017b). In their analysis, they also include results from several standard machine learning algorithms. They note that the rotation forest (Rodriguez et al., 2006) and random forest (Breiman, 2001) are competitive with their time series classification baseline (1nearest neighbor with dynamic time warping distance; Tormene et al., 2008). Their results show that ensemble methods such as collection of transformation ensembles (COTE; Bagnall et al., 2015) perform best, but for the price of considerable runtime.
Deep learning methods applied to time series classification tasks have also shown competitive prediction power. For example, Fawaz et al. (2019) provide a comprehensive review of stateoftheart methods. The authors compared both generative models and discriminative models, including
fully connected neural networks, convolutional neural networks, autoencoders
and echo state networks, whereas only discriminative endtoend approaches were incorporated in the benchmark study.The benchmark study conducted in this work does not aim to replicate or compete with earlier studies like (Bagnall et al., 2017a), but instead tries to extend their results.
3 Functional Data
In contrast to nonfunctional data analysis, where the measurement of a single observation is a vector of scalar components whose entries represent values of the separate multidimensional features, functional data analysis treats and analyses the features themselves as functions over their domain. By learning to represent the underlying function, the carried out analysis is not just restricted to the measured discrete values but it is possible to sample from (and analyze) the entire domain space.
In this work, we focus on pairs of features and corresponding labels for supervised learning. In contrast to nonfunctional data analysis, where the measurement of a single observation is a vector of scalar components, functional features are functionvalued over their domain. The features can thus also be a function, i.e., , . In practice, functional data comes in the form of observed values , where each corresponds to a discrete point on the continuum. Those observed values stem from an underlying function evaluated over a set of points. A frequent type of functional data is time series data, i.e., measurements of a process measured at discrete timepoints.
For example, in some electrical engineering applications, signals are obtained over time at a certain sampling rate, but other domains are possible as well. Spectroscopic data, for example, are functional data recorded over certain parts of the electromagnetic spectrum. One such example is depicted in Figure 1. It shows spectroscopy data of fossil fuels (Fuchs et al., 2015) where the measured signal represents reflected energies in the ultravioletvisible (UVVIS) and the near infrared spectrum (NIR). In the plot, different colors correspond to different instances. This is a typical example of a scalaronfunction regression problem, where the inputs are a collection of spectroscopic curves for a fuel, and the prediction target is the heating value of the fossil fuel.
In Figure 2, we display two functional classification scenarios. The goal in those scenarios is to distinguish the class type of the curve, which can also be understood as a functiononscalar problem. Figure 1(a) shows the vertical position of an actor’s hand while either drawing a toygun and aiming at a target, or just imitating the motion with the blank hand. This position is measured over time. The two different types of classes of the curves can be distinguished by the color scheme.
Figure 1(b) shows a data set built for distinguishing images of beetles from images of flies based on their outlines. While following the outline, the distance to the center of the object is measured which is then used for classification purposes.
The latter data sets are available from Bagnall et al. (2017b).
4 Functional Data Analysis with mlrFDA
Along with the benchmark, we implement the software mlrFDA, which extends the popular machine learning framework mlr for the analysis of functional data. As the implemented functionality is an extension of the mlr package, all of the functionality available in mlr transfers to the newly added methods for functional data analysis. We include a brief overview of the implemented functionality in A.1. A more indetail overview and tutorial on mlr can be found in the mlr tutorial (Schiffner et al., 2016).
mlr provides a unified framework for machine learning methods in R, currently supporting tasks from
main problem types: (multilabel)classification, regression, cluster analysis, and survival analysis. For each problem type, many algorithms (called
learners) are integrated. This yields an extensive set of modeling options with a unified, simple interface. Moreover, advanced techniques such as hyperparameter tuning, preprocessing and feature selection are also part of the package. An additional focus lies on extensibility, allowing the user to integrate their own algorithms, performance measures and preprocessing methods. As
mlrMBO seamlessly integrates into the new software, many different tuning procedures can be readily adapted by the user. Tuning of hyperparameters is usually not integrated in software packages for functional data analysis and thus would require the user to write additional, nontrivial code that handles (nested) resampling, evaluation and optimization methods.mlrFDA contains several functional data algorithms from several R packages, e.g., fda.usc, refund or FDboost. The algorithms’ functionality, however, remains unchanged, only their user interface is standardized for use with mlr. For detailed insights into the respective algorithms, full documentation is available in the respective packages.
Since our toolbox is built on mlr’s extensible class system, our framework is easily extensible to other methods that have not yet been integrated, and the user can include his or her own methods which do not necessarily need to be available as a packaged implementation. Additionally, mlrFDA inherits mlr’s functionality for performance evaluation and benchmarking, along with extensive and advanced (hyperparameter) tuning. This makes our platform very attractive for evaluating which algorithm fits best to a problem at hand, and even allows for large benchmark studies.
5 Benchmark Experiment
In order to enable a comparison of the different approaches, an extensive benchmark study is conducted. This paper does not aim to replicate or reproduce results obtained by Bagnall et al. (2017a) or Fawaz et al. (2019). Instead we focus on providing a benchmark complementary to previous benchmarks. This is done because the experiments require large amounts of computational resources, and
the added value of an exact replication of the experiments (with open source code) is comparatively small. Nonetheless, we aim for results that can be compared, and thus extend the results obtained by
Bagnall et al. (2017a) by staying close to their setup. The experiments were carried out on a high performance computing cluster, supported by the Leibniz Rechenzentrum Munich. Individual runs were allowed up to GB of RAM and hours runtime for each evaluation. We want to stress that this benchmark compares implementations, which does not always necessarily correspond to the performance of the corresponding theoretical algorithm. Additionally, methods for functional data analysis are traditionally more focused on valid statistical inference and interpretable results, which does not necessarily coincide with high predictive performance.5.1 Benchmark Setup
Data sets  51 Data sets, see table 7 

Algorithms  Function (Package) 
Machine Learning:   glmnet (glmnet) 
 rpart (rpart)  
 ksvm (kernlab)  
 ranger (ranger)  
 xgboost (xgboost)  
Functional Data 
 classif.knn (fda.usc) 
 classif.glm (fda.usc)  
 classif.np (fda.usc)  
 classif.kernel(fda.usc)  
 FDboost (FDboost)  
 fgam (refund)  
 knn with dtw (classiFunc)  
Feature Extraction + ML   feature extraction: see table 6 
 in combination with ML algorithms marked with a .  
Measures  mean misclassification error, training time 
Resampling  20fold stratified subsampling; 
class balances and train/test set size as in Bagnall et al. (2017a).  
Tuning  100 iterations of Bayesian optimization (3fold inner CV). 
Corresponding hyperparameterranges can be obtained  
from tables 3 and 5. 
A benchmark experiment is defined by four important characteristics: The data sets algorithms are tested on, the algorithms to be evaluated, the measures used for evaluating predictive performance, and a resampling strategy used for generating train and test splits of the data. A comprehensive overview of the conducted benchmark setup can be obtained from Table 1.
These characteristics are briefly described subsequently before providing and discussing the results. We use a subset of 51 data sets from the popular UCR archive (Bagnall et al., 2017b) in order to enable a comparison of results in Bagnall et al. (2017a) with the additional methods described in this paper. The data sets stem from various application types such as ECG measurements, sensor data, or image outlines, therefore having varying training set sizes or measurement lengths. For more detailed information about the data sets, see Bagnall et al. (2017b).
We selected data using the following criteria: In order to reduce the computational resources we did not run data sets that have multiple versions, exclude data sets with less then 3 examples in each class remove data sets with more than instances or time series longer than
measurements. As some of the classifiers only work with multiclass targets via
1vsall classification, we additionally excluded data sets with more then 40 classes. In essence, we benchmark small and medium sized data sets with a moderate amount of different classes.We add new algorithms and feature extraction methods which can be combined with arbitrary machine learning methods for scalar features (c.f. Table 1). Additionally we test classical machine learning methods, in order to obtain a broader perspective on expected performance if the functional nature of the data is ignored. As we benchmark default settings as well as tuned algorithms, in total different algorithms are evaluated across all data sets. When combining feature extraction and machine learning methods, we fuse the learning algorithm and the preprocessing, thus treating them as a pipeline where data is internally transformed before applying the learner. This allows us to jointly tune the hyperparameters of learning algorithm and preprocessing method. The respective defaults and parameter ranges can be obtained from Table 3 (feature extractors) and Table 5 (learning algorithms). More detailed description of the hyperparameters can be obtained from the respective packages documentation.
In order to generate train/test splits, and thus obtain an unbiased estimate of the algorithm’s performance, we use stratified subsampling. We use 20 different train/test splits for each data set in order to reduce variance and report the average. For tuned models, we use use nested crossvalidation
(Bischl et al., 2012) to ensure unbiased estimates, where the outer loop is again subsampling with 20 splits, and the inner resampling for tuning is a 3fold (stratified) crossvalidation. All compared 80 algorithms are presented exactly the same index sets for the 20 traintest outer subsampling splits.Mean misclassification error (MMCE) is chosen as a measure of predictive performance in order to stay consistent with Bagnall et al. (2017a)
. Other measures, such as area under the curve (AUC) require predicted probabilities and do not trivially extend to multiclass settings.
While Bagnall et al. (2017a) tune all algorithms across a carefully handcrafted grid, we use Bayesian optimization (Snoek et al., 2012). In order to stay comparable, we analogously fix the amount of tuning iterations to 100.
We use mlrMBO (Bischl et al., 2017) in order to perform Bayesian optimization of the hyperparameters of the respective algorithm. Additionally, in order to scale the method to a larger amount of data sets and machines, the Rpackage batchtools (Bischl et al. (2015), Lang et al. (2017)) is used. This enables running benchmark experiments on highperformance clusters. For the benchmark experiment, a job is defined as resampling of a single algorithm (or tuning thereof) on a single version of a data set. This allows for parallelization to an arbitrary number of CPU’s, while at the same time guaranteeing reproducibility. The code for the benchmark is available from https://github.com/compstatlmu/2019_fda_benchmark for reproducibility.
5.2 Results
This Section tries to answer the questions posed in section 1. We evaluate various machine learning algorithms in combination with feature extraction, classical time series classification approaches, the effect of tuning hyperparameters for several methods, and try to give recommendations with respect to which algorithm(s) to choose for new classification problems.
Algorithms evaluated in this benchmark have been divided into three groups: Algorithms specifically tailored to functional data, classical machine learning algorithms without feature extraction and classical machine learning algorithms in combination with feature extraction.
5.2.1 Algorithms for functional data
Performances of algorithms specifically tailored to functional data analysis can be obtained from Figure 3. The knearest neighbors algorithm from package classiFunc (Maierhofer and Pfisterer, 2018) in combination with dynamic time warping (Rakthanmanon et al., 2012) distance seems to perform best across data sets. It is also considered a strong baseline in Bagnall et al. (2017a).
5.2.2 Machine Learning algorithms with feature extraction
Performances of different machine learning algorithms in combination with feature extraction with and without tuning can be obtained from Figure 4.
We conclude that feature extraction using splines (bsignal) and wavelets as well as extracting dynamic time warping distances works well when combined with conventional machine learning algorithms, even at their default hyperparameters. Among the learners, random forests, especially in combination with bsignal show quite advantageous performance. In addition, we find an obvious improvement from hyperparameter tuning for the Fourier feature extraction. In terms of learners, random forest and gradient boosted tree learners (xgboost) perform better than support vector machines.
5.2.3 Machine Learning algorithms without feature extraction
Additionally, we aim to provide some insight with regards to the performance of machine learning algorithms that ignore the functional nature of our data. Figure 5 provides an overview over the performance of different machine learning algorithms that are often used together with traditional tabular data. Performances in this figure are obtained from algorithms directly applied to the functional data without any additional feature extraction. The widely used gradient boosting (xgboost) and random forest (ranger) implementations seem to work reasonably well for functional data even without additional feature extraction.
5.2.4 The effect of tuning hyperparameters
From our experiments, we conclude, that tuning hyperparameters of machine learning algorithms in general has a nonnegligible effect on the performance. Using Bayesian optimization in order to tune algorithm hyperparameters on average yielded an absolute increase in accuracy of 5.4% across data sets and learners.
Figure 6 displays the aggregated time over all data sets, taking into account the time required for hyperparameter tuning. All experiments have been run on equivalent hardware on highperformance computing infrastructure. Due to fluctuations in server load, this does not allow for an exact comparison with respect to computation time, but we hope to achieve comparable results as we repeatedly evaluate on subsamples. Note that we restrict the tuning to algorithms where tuning traditionally leads to higher performances.^{3}^{3}3Additionally, we find significantly improved performance for tuned FDboost in Figure 3
5.2.5 Top 10 Algorithms and recommendations
Table 2 showcases the top algorithms from the benchmark in terms of average rank in predictive accuracy across data sets. With this list, we aim to provide some initial understanding of the performance of different algorithms and feature extraction methods. Note that this list by no means reflects performance on future data sets, but might serve as an indicator, of which algorithms one might want to try first given computational constraints.
Algorithm Setting  Accuracy %  Average Rank 

ranger_wavelet_tuned  0.92  12.90 
xgboost_wavelet_tuned  0.92  14.45 
ranger_bsignal_tuned  0.90  15.02 
knn_dtw_tuned  0.92  15.22 
ranger_none_default  0.90  15.59 
ranger_bsignal_default  0.89  15.71 
ranger_wavelet_default  0.90  16.33 
knn_dtw_default  0.92  16.43 
xgboost_bsignal_tuned  0.90  17.57 
ranger_none_tuned  0.89  18.49 
We observe that wavelet extraction in combination with either ranger or xgboost seems to be very strong. They obtain an average rank of 12.90 and 14.45 (out of 80) respectively. Dynamic time warping distances for knearest neighbors indeed seems to be a strong baseline, even without tuning. Another strong feature extraction method seems to be the extraction of Bspline features. Using the algorithms above allows us to obtain an accuracy within of the maximum on of the data sets.
If the only criterion for model selection is predictive performance, (tuned) machine learning models in combination with feature extraction is a competitive baseline. This class of methods achieves within 95% of the optimal performance on out of data sets, while they include the best performing classifier in cases.
id  type  values  def.  trafo 

bsignal  
bsignal.knots  int  {3,…,500}  10   
bsignal.df  int  {1,…,10}  3   
multires  
res.level  int  {2,…,5}     
shift  num  [0.01,1]     
pca  
rank.  int  {1,…,30}     
wavelets  
filter  chr  d4,d8,d20,la8,la20,bl14,bl20,c6,c24     
boundary  chr  periodic,reflection     
fourier  
trafo.coeff  chr  phase,amplitude     
dtwkernel  
ref.method  chr  random,all  random   
n.refs  num  [0,1]     
dtwwindow  num  [0,1]     
5.2.6 Comparison to classical time series classification
Even though the main purpose of this paper is not a direct comparison with the results from Bagnall et al. (2017a), we can use our results to show that applying functional data approaches and classical machine learning approaches together with feature extraction can still improve classification accuracy compared to current stateoftheart time series classification methods.
In the experiments we conducted, the methods described in this paper improved accuracy on 9 out of the 51 data sets which is displayed in Figure 7. The 9 data sets and the corresponding best learner are displayed in Table 4. For each data set, only the best reached accuracy for both sets of algorithms is displayed.
Name  Algorithm_Setting  Accuracy 

Beef  xgboost_wavelet_tuned  0.83 
ChlorineConcentration  ksvm_none_tuned  0.91 
DistalPhalanxOutlineAgeGroup  ranger_none_default  0.83 
DistalPhalanxOutlineCorrect  ranger_dtwkernel_default  0.83 
DistalPhalanxTW  ranger_bsignal_default  0.76 
Earthquakes  FDboost_none_default  0.80 
Ham  xgboost_wavelet_tuned  0.84 
InsectWingbeatSound  ranger_wavelet_default  0.65 
SonyAIBORobotSurface1  ksvm_wavelet_default  0.94 
Additionally, we evaluate how our learners rank in comparison to the individual bakeoff algorithms from Bagnall et al. (2017a). The algorithm which performs best on a data set obtains the rank 1. The mean rank of the individual learners over all data sets (we take the intersection of the data sets from our benchmark and the ones from Bagnall et al. (2017a)). The average sorted ranks for the top 50% algorithms are displayed in Figure 8. We observe that the ensemble methods get the top ranks, which is no surprise, as for instance the COTE algorithm Bagnall et al. (2015) internally combines several classifiers from 4 different time series domains.
However, compared to the classical time series algorithms from Bagnall et al. (2017a) with the ensemble methods removed, our functional data algorithms obtain an overall good rank in accuracy performance, interleaved with the algorithms from Bagnall et al. (2017a). Note that the benchmarks are not exactly comparable due to minor differences in the benchmark setup, and we instead only include their reported results.
parameter  type  values  default  trafo 
ksvm  
C  num  [15,15]    2^ x 
sigma  num  [15,10]    2^ x 
ranger  
mtry.power  num  [0,1]    
min.node.size  num  [0,0.99]    
sample.fraction  num  [0.1,1]     
xgboost  
nrounds  int  {1,…,5000}  100   
eta  num  [10,0]    2^ x 
subsample  num  [0.1,1]     
booster  chr  gbtree,gblinear     
max_depth  int  {1,…,15}     
min_child_weight  num  [0,7]    2^ x 
colsample_bytree  num  [0,1]     
colsample_bylevel  num  [0,1]     
lambda  num  [10,10]    2^ x 
alpha  num  [10,10]    2^ x 
FDboost  
mstop  int  {1,…,5000}  100   
nu  num  [0,1]  0.01   
df  num  [1,5]  4   
knots  int  {5,…,100}  10   
degree  int  {1,…,4}  3   
6 Summary and Outlook
In this work, we provide a benchmark along with a software implementation that integrates the functionality of a diverse set of Rpackages into a single user interface and API. Both contributions come with a multiplicity of benefits:

The user is not required to learn and deal with the vast complexity of the different interfaces the underlying packages expose.

All of the existing functionality (e.g., preprocessing, resampling, performance measures, tuning, parallelization) of the mlr ecosystem can now be used in conjunction with already existing algorithms for functional data.

We expose functionality that allows us to work with functional data using traditional machine learning methods via feature extraction methods.

Integration of additional preprocessing methods or models is (fairly) trivial and automatically benefits from the full mlr ecosystem.
In order to obtain a broader overview of the performance of the integrated methods, we perform a large benchmark study. This allows users to get an initial overview of potential performances of the different algorithms. Specifically,

We open up new perspectives for time series classification tasks by incorporating methods from functional data analysis, as well as feature transformations combined with conventional machine learning models.

Based on the large scale benchmark, we conclude that many learners have competitive performance (Figure 7) and additionally offer the interpretability of many functional data analysis methods. Our toolbox serves as a strong complement and alternative to other time series classification software.

The presented benchmark study uses stateoftheart Bayesian optimization for hyperparameter optimization, which results in significant improvements over models that are not tuned. This kind of hyperparameter tuning is easy to do with mlrFDA. Tuning, albeit heavily influencing performance is often not investigated. Our benchmark closes this gap in existing literature.

We find that extracting vector valued features and feeding them to a conventional machine learning model can often form competitive learners.

The paretooptimal set in terms of performance on each data set contains different algorithmfeatureextraction combinations. Our toolbox offers the same API for all methods and allows to automatically search over this space, and thus allows users to obtain optimal models without knowing all underlying methods.
Concerning the questions we proposed at the beginning of the paper, we draw the following conclusions:

Tuning only a subset of the presented learners and feature extractions, i.e., the methods listed in Table 2, is sufficient to achieve good performances on almost all data sets in our benchmark.

A simple random forest without any preprocessing can also be a reasonable baseline for time series data. It achieves an average rank of (top 4) in our benchmark.

Most algorithms for functional data (e.g., FDboost) do not perform well in our benchmark study. As those algorithms are fully interpretable and offer statistically valid coefficients, they can still be useful in some applications, and should thus not be ruled out.

Feature extraction techniques, such as bspline representations (bsignal) and wavelet extraction work well in conjunction with machine learning techniques for vector valued features such as xgboost and random forest.

Tuning leads to an average reduction in absolute MMCE of (ranger), (xgboost), (ksvm) (across feature extraction techniques) and (FDboost). This holds for all feature extraction techniques, where improvements range from multires to fourier.
In future work we will continue to expand the available toolbox along with a benchmark of new methods, and provide the R community a wider range of methods that can be used for the analysis of functional data. This includes not only integrating many already available packages, and as a result to enable preprocessing operations such as smoothing (e.g., fda (Ramsay et al., 2018)) and alignment (e.g., fdasrvf (Tucker, 2016) or tidyfun (Scheipl and Goldsmith, 2019)
), but also to explore and integrate advanced imputation methods for functional data. Further work will also extend the current implementation to support data that is measured on unequal or irregular grids. Additionally, we aim to implement some of the current stateofthe art machine learning models from the time series classification bakeoff
(Bagnall et al., 2017a), such as the Collective of TransformationBased Ensembles (COTE) (Bagnall et al., 2015). This enables researchers to use and compare with current stateoftheart methods.References
 Wavelets: a package of functions for computing wavelet filters, wavelet transforms and multiresolution analyses. Note: R package version 0.30 External Links: Link Cited by: §A.2, item 4, §2.1, §2.1.
 Timeseries classification with cote: the collective of transformationbased ensembles. IEEE Transactions on Knowledge and Data Engineering 27 (9), pp. 2522–2535. External Links: Document, ISSN 10414347 Cited by: §2.4, §6.
 The great time series classification bake off: a review and experimental evaluation of recent algorithmic advances. Data Mining and Knowledge Discovery 31 (3), pp. 606–660. Cited by: item 1, §1, §2.3, §2.4, §2.4, Figure 7, Figure 8, §5.1, §5.1, §5.1, §5.2.1, §5.2.6, §5.2.6, §5.2.6, Table 1, §5, §6.
 Timeseries classification with cote: the collective of transformationbased ensembles. IEEE Transactions on Knowledge and Data Engineering 27 (9), pp. 2522–2535. Cited by: §5.2.6.
 The UEA & UCR time series classification repository. Note: www.timeseriesclassification.com Cited by: §2.4, §3, §5.1.

Resampling methods for metamodel validation with recommendations for evolutionary computation
. Evolutionary Computation 20 (2), pp. 249–275. Note: PMID: 22339368 External Links: Document, Link Cited by: §5.1.  BatchJobs and BatchExperiments: abstraction mechanisms for using R in batch environments. Journal of Statistical Software 64 (11), pp. 1–25. External Links: Link Cited by: §5.1.
 mlrMBO: A Modular Framework for ModelBased Optimization of Expensive BlackBox Functions. External Links: Link Cited by: §5.1.
 rucrdtw: fast time series subsequence search in r. The Journal of Open Source Software 1, pp. 1–2. External Links: Document, Link Cited by: §A.2, §2.1.
 Random forests. Mach. Learn. 45 (1), pp. 5–32. External Links: ISSN 08856125, Link, Document Cited by: §2.4.
 The fast fourier transform. IEEE Spectrum 4 (12), pp. 63–70. External Links: Document, ISSN 00189235 Cited by: §2.1.
 FDboost: boosting functional regression models. Cited by: item 3, §1, §2.1, §2.2.
 Supportvector networks. Mach. Learn. 20 (3), pp. 273–297. External Links: ISSN 08856125, Link, Document Cited by: §1.
 Deep learning for time series classification: a review. Data Mining and Knowledge Discovery, pp. 1–47. Cited by: §1, §2.4, §5.
 Statistical computing in functional data analysis: the r package fda.usc. Journal of Statistical Software 51 (4), pp. 1–28. External Links: Link Cited by: §1, §2.2, §2.3.
 Nonparametric functional data analysis: theory and practice. Springer Science & Business Media. Cited by: §2.2.
 Penalized scalaronfunctions regression with interaction term. Computational Statistics & Data Analysis 81, pp. 38–51. External Links: Document Cited by: Figure 1, §3.
 Refund: regression with functional data. Note: R package version 0.117 External Links: Link Cited by: §1, §2.2.
 Estimator selection and combination in scalaronfunction regression. Computational Statistics & Data Analysis 70, pp. 362–372. Cited by: §2.1.
 A general framework for functional regression modelling. Statistical Modelling 17 (12), pp. 1–35. Cited by: §2.2.
 Scikitfda: functional data analysis in python. External Links: Link Cited by: §2.3.
 Classification of time series by shapelet transformation. Data Min. Knowl. Discov. 28 (4), pp. 851–881. External Links: ISSN 13845810, Link, Document Cited by: Figure 2.
 Modelbased boosting 2.0. Journal of Machine Learning Research 11 (Aug), pp. 2109–2113. Cited by: §2.2.
 Tsfeatures: time series feature extraction. Note: R package version 0.1 External Links: Link Cited by: item 4, §2.1, §2.1.
 Driver sleepiness detection system based on eye movements variables. Advances in Mechanical Engineering 5 (), pp. 648431. External Links: Document, Link Cited by: §A.2.
 Introduction to functional data analysis. CRC Press. Cited by: §3.
 Batchtools: tools for r to work on batch systems. The Journal of Open Source Software 2 (10). External Links: Document, Link Cited by: §5.1.
 ClassiFunc: classification of functional data. Note: R package version 0.1.1 External Links: Link Cited by: §2.2, §5.2.1.
 A theory for multiresolution signal decomposition: the wavelet representation. IEEE Trans. Pattern Anal. Mach. Intell. 11, pp. 674–693. Cited by: §2.1.
 EEG feature extraction for classifying emotions using fcm and fkm. In Proceedings of the 7th WSEAS International Conference on Applied Computer and Applied Computational Science, ACACOS’08, Stevens Point, Wisconsin, USA, pp. 299–304. External Links: ISBN 9789606766619, Link Cited by: §A.2.
 R: a language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. External Links: Link Cited by: §1.
 Searching and mining trillions of time series subsequences under dynamic time warping. In Proceedings of the 18th ACM SIGKDD international conference on Knowledge discovery and data mining, pp. 262–270. External Links: Document, Link Cited by: §A.2, §2.1, §5.2.1.
 Fda: functional data analysis. Note: R package version 2.4.8 External Links: Link Cited by: §1, §2.3, §6.
 Functional data analysis. Wiley Online Library. Cited by: §2.2, §3.
 Three myths about dynamic time warping data mining.. In SDM, H. Kargupta, J. Srivastava, C. Kamath, and A. Goodman (Eds.), pp. 506–510. External Links: ISBN 9781611972757 Cited by: Figure 2.
 Rotation forest: a new classifier ensemble method. IEEE Transactions on Pattern Analysis and Machine Intelligence 28 (10), pp. 1619–1630. External Links: Document, ISSN 01628828 Cited by: item 1, §2.4.
 Tidyfun. GitHub. Note: https://github.com/fabians/tidyfun Cited by: §6.
 CRAN task view  functional data analysis. External Links: Link Cited by: §1.
 Mlr tutorial. arXiv preprint arXiv:1609.06146. Cited by: §4.
 Practical bayesian optimization of machine learning algorithms. In Advances in Neural Information Processing Systems 25, F. Pereira, C. J. C. Burges, L. Bottou, and K. Q. Weinberger (Eds.), pp. 2951–2959. External Links: Link Cited by: §5.1.
 On the use of the wavelet decomposition for time series prediction. Neurocomputing 48 (1), pp. 267 – 277. External Links: ISSN 09252312, Document, Link Cited by: §A.2.
 Functional and shape data analysis. Springer. Cited by: §2.2.
 Show me how you drive and i’ll tell you who you are. recognizing gender using automotive driving parameters. Procedia Manufacturing 3, pp. 5587 – 5594. Note: 6th International Conference on Applied Human Factors and Ergonomics (AHFE 2015) and the Affiliated Conferences, AHFE 2015 External Links: ISSN 23519789, Document, Link Cited by: item 4.
 Matching incomplete time series with dynamic time warping: an algorithm and an application to poststroke rehabilitation. Artificial Intelligence in Medicine 45 (1), pp. 11–34. External Links: Document Cited by: §2.4.
 Fdasrvf: elastic functional data analysis. Note: R package version 1.6.0 External Links: Link Cited by: §6.
 Applications of functional data analysis: a systematic review. BMC Medical Research Methodology 13 (1), pp. 43. External Links: ISSN 14712288, Document, Link Cited by: §1, §2.1.
 Functional data analysis. Annual Review of Statistics and Its Application 3 (1), pp. 257–295. External Links: Document, Link Cited by: §1.
Appendix A API overview
For the interested reader, we introduce a brief overview of the API and functionality.
a.1 Representing functional data in mlrFDA
A sketch of the data structure we use to represent functional data can be found in the right part of Figure 9. We assume a data set consists of data for observational units, organized in rows of features, where one row contains all observed features for one observational unit, i.e., each row typically contains several functional and/or scalar covariates. For a classical, nonfunctional data set, the features are single columns (as depicted in the left part of Figure 9). A functional data set, on the other hand, consists of singlecolumn scalar features as well as functional features of different length for each functional covariate, each represented by multiple adjacent and connected columns. Each of these columns contains the evaluations of the functional feature at a certain argument value for all observational units (right part of Figure 9).
As an example which will be used throughout the remainder of this paper, we use the fuelSubset data set from package FDboost, see also Figure 1. It contains a numeric target variable heatan, the fuel’s heating value, a scalar feature h20, the fuel’s water content, and two functional features NIR and UVVIS, measured at 231 and 129 wavelengths, respectively. To start with a clean sheet, we create a data.frame containing all features as separate columns.
R¿ library(mlr) R¿ library(FDboost) R¿ df = data.frame(fuelSubset[c(”heatan”, ”h2o”, ”UVVIS”, ”NIR”)])
The first step when setting up an experiment in any analysis is to make the data accessible for the specific algorithms that will be applied. In mlr, the data itself, and additional information, such as which column corresponds to the target variable is stored as a Task, requiring the input data to be of type data.frame.
The list of column positions of the functional features is then passed as argument fd.features to makeFunctionalData(), which returns an object of type data.frame in which the columns corresponding to each functional feature are combined into matrix columns.
^{4}^{4}4As an alternative, a list of the column names containing the functional features is also valid as argument to fd.features, which is especially useful if columns are already labeled.
R¿ fd.features = list(”UVVIS” = 3:136, ”NIR” = 137:367) R¿ fdf = makeFunctionalData(df, fd.features = fd.features) R¿ str(fdf) ’data.frame’: 129 obs. of 4 variables: $ heatan: num 26.8 27.5 23.8 18.2 17.5 … $ h2o : num 2.3 3 2 1.85 2.39 … $ UVVIS : num [1:129, 1:134] 0.145 1.584 0.814 1.311 1.373 … .. attr(*, ”dimnames”)=List of 2 .. ..$ : NULL .. ..$ : chr ”UVVIS.1” ”UVVIS.2” ”UVVIS.3” ”UVVIS.4” … $ NIR : num [1:129, 1:231] 0.2818 0.2916 0.0042 0.034 0.1804 … .. attr(*, ”dimnames”)=List of 2 .. ..$ : NULL .. ..$ : chr ”NIR.1” ”NIR.2” ”NIR.3” ”NIR.4” …
We additionally specify the name "fuelsubset" and the target variable "heatan". The structure of the functional Task object is rather similar to the nonfunctional Task, with the additional information functionals, which states how many functional features are present in the underlying data.
R¿ tsk1 = makeRegrTask(”fuelsubset”, data = fdf, target = ”heatan”) R¿ print(tsk1) Supervised task: fuelsubset Type: regr Target: heatan Observations: 129 Features: numerics factors ordered functionals 1 0 0 2 Missings: FALSE Has weights: FALSE Has blocking: FALSE Has coordinates: FALSE
After defining the task, a learner is created by calling makeLearner. This contains the algorithm that will be fitted on the data in order to obtain a model. Currently, mlrFDA supports both functional regression and functional classification. A list of supported learners can be found in Table 1
a.2 Machine Learning and Feature Extraction
Classical machine learning algorithms do not take into account the characteristics of functional data and treat the input data as vector valued features. Without additional preprocessing, this typically yields poor performance on, as the models cannot exploit the lower intrinsic dimensionality of the functional covariates nor the fact that they represent observations over a continuum.
In mlrFDA, classical algorithms can be applied to functional data, however, a warning message will be displayed. In our example, we train a partitioning tree on the functional data, while ignoring the functional structure.
R¿ rpart.lrn = makeLearner(”regr.rpart”) R¿ m = train(learner = rpart.lrn, task = tsk1) Functional features have been converted to numerics
For conventional learning algorithms to work well on functional data, informative scalar features need to be extracted from the functional features.
Feature extraction is applied in practice for a multiplicity of reasons, as it often not only reduces the dimensionality of the resulting problem, but also allows researchers to make use of domain knowledge, for example by handcrafting features from measurements of continuous processes. Examples for this include deriving features that allow for sleepiness detection (Jin et al., 2013), or by extracting features from electrocardiogram data in order to detect emotions (Murugappan et al., 2008). The resulting features often have a much lower dimensionality, which often improves fitted models. Other preprocessing methods for functional or time series data include extracting general purpose features such as wavelet coefficients (Aldrich, 2013; Soltani, 2002), principal component scores or Fourier coefficients. The resulting scalar features can then be used with different machine learning methods such as nearest neighbors.
In the following section, we showcase the feature extraction procedure using general purpose features as an example. We want to emphasize that it is also easily possible to write custom feature extraction methods using the makeFeatureExtractionMethod function.
Name  Function  Package 

Discrete Wavelet Transform  extractFDAWavelets()  wavelets 
Fast Fourier Transform  extractFDAFourier()  stats 
Principal Component Analysis  extractFDAPCA()  stats 
BSpline Features  extractFDABsignal()  FDboost 
MultiResolution Feature Extraction  extractFDAMultiResFeatures()   
Time Series Features  extractFDATsfeatures()  tsfeatures 
Dynamic TimeWarping Kernel  extractFDADTWKernel()  rucrdtw 
In our example, we extract the Fourier coefficients from the functional feature UVVIS, and principal component scores from the second functional feature NIR in order to transform the original task with functional data into a conventional task.
R¿ feat.methods = list(”UVVIS” = extractFDAFourier(), ”NIR” = extractFDAPCA()) R¿ extracted = extractFDAFeatures(tsk, feat.methods = feat.methods) R¿ extracted
desc Extraction of features from functional data: Target: heatan Functional Features: 2; Extracted features: 2
As an alternative, the feature extraction can be applied in a wrapper method
makeExtractFDAFeatsWrapper(). In general, a wrapper combines a learner method with another method, thereby creating a new learner that can be handled like any other learner. In our case, a classical machine learning method is combined with the data preprocessing step of feature transformation from functional to nonfunctional data.
R¿ wrapped.lrn = makeExtractFDAFeatsWrapper(”regr.rpart”, feat.methods = feat.methods)
This is suitable for honest crossvalidation of dataadaptive feature extraction methods like principal components. We can now crossvalidate the learner created above using mlr’s resample function with 10fold crossvalidation.
R¿ res = resample(learner = wrapped.lrn, task = tsk1, resampling = cv10)
In the same way, we can train and predict on data, or benchmark multiple learners across multiple data sets. Additionally, we can apply a tuneWrapper to our learner in order to automatically tune hyperparameters of the learner and the preprocessing method during crossvalidation.
Appendix B Data sets used in the Benchmark
Table 7 contains all data sets used in the benchmark along with additional data properties.
Name  Obs.  Classes  Length  Type  Split 

Adiac  781  37  176  IMAGE  0.50 
ArrowHead  211  3  251  IMAGE  0.17 
Beef  60  5  470  SPECTRO  0.50 
BeetleFly  40  2  512  IMAGE  0.50 
BirdChicken  40  2  512  IMAGE  0.50 
Car  120  4  577  SENSOR  0.50 
CBF  930  3  128  SIMULATED  0.03 
ChlorineConcentration  4307  3  166  SIMULATED  0.11 
Coffee  56  2  286  SPECTRO  0.50 
Computers  500  2  720  DEVICE  0.50 
CricketX  780  12  300  MOTION  0.50 
DistalPhalanxOutlineAgeGroup  539  3  80  IMAGE  0.74 
DistalPhalanxOutlineCorrect  876  2  80  IMAGE  0.68 
DistalPhalanxTW  539  6  80  IMAGE  0.74 
Earthquakes  461  2  512  SENSOR  0.70 
ECG200  200  2  96  ECG  0.50 
ECGFiveDays  884  2  136  ECG  0.03 
ElectricDeviceOn  1008  2  360  DEVICE  0.63 
EpilepsyX  275  4  208  HAR  0.61 
FaceAll  2250  14  131  IMAGE  0.25 
FacesUCR  2250  14  131  IMAGE  0.09 
Fish  350  7  463  IMAGE  0.50 
GunPoint  200  2  150  MOTION  0.25 
Ham  214  2  431  SPECTRO  0.51 
Herring  128  2  512  IMAGE  0.50 
InsectWingbeatSound  2200  11  256  SENSOR  0.10 
ItalyPowerDemand  1096  2  24  SENSOR  0.06 
LargeKitchenAppliances  750  3  720  DEVICE  0.50 
Lightning2  121  2  637  SENSOR  0.50 
Lightning7  143  7  319  SENSOR  0.49 
Meat  120  3  448  SPECTRO  0.50 
MedicalImages  1141  10  99  IMAGE  0.33 
MoteStrain  1272  2  84  SENSOR  0.02 
OSULeaf  442  6  427  IMAGE  0.45 
Plane  210  7  144  SENSOR  0.50 
RefrigerationDevices  750  3  720  DEVICE  0.50 
ScreenType  750  3  720  DEVICE  0.50 
ShapeletSim  200  2  500  SIMULATED  0.10 
SmallKitchenAppliances  750  3  720  DEVICE  0.50 
SonyAIBORobotSurface1  621  2  70  SENSOR  0.03 
Strawberry  983  2  235  SPECTRO  0.62 
SwedishLeaf  1125  15  128  IMAGE  0.44 
SyntheticControl  600  6  60  SIMULATED  0.50 
ToeSegmentation1  268  2  277  MOTION  0.15 
Trace  200  4  275  SENSOR  0.50 
TwoLeadECG  1162  2  82  ECG  0.02 
TwoPatterns  5000  4  128  SIMULATED  0.20 
UWaveGestureLibraryX  4478  8  315  MOTION  0.20 
Wafer  7164  2  152  SENSOR  0.14 
Wine  111  2  234  SPECTRO  0.51 
Yoga  3300  2  426  IMAGE  0.09 
Appendix C Failed and missing experiments
Experiments for some algorithm / data set combinations failed due to implementation details or algorithm properties. In order to increase transparency, failed algorithms are listed here, and if available reasons for failure are provided.
At the time of the benchmark, the implementation in the tsfeatures package was not stable enough to be included in the benchmark.

classif.fgam
Data sets: BeetleFly, BirdChicken, Coffee, Computers, DistalPhalanxOutlineCorrect, Earthquakes, ECG200, ECGFiveDays, ElectricDeviceOn, GunPoint, Ham, Herring, ItalyPowerDemand, Lightning2, MoteStrain, ShapeletSim, SonyAIBORobotSurface1, Strawberry, ToeSegmentation1, TwoLeadECG, Wafer, Wine, Yoga
Reason: Too few instances in some classes, such that p ¿ n. 
classif.fdausc.kernel and .np Data sets: ElectricDeviceOn, ShapeletSim

classif.fdausc.knn
Data sets: DistalPhalanxTW, EpilepsyX