Multi-class classification of biomechanical data: A functional LDA approach based on multi-class penalized functional PLS

Abstract

A functional linear discriminant analysis approach to classify a set of kinematic data (human movement curves of individuals performing different physical activities) is performed. Kinematic data, usually collected in linear acceleration or angular rotation format, can be identified with functions in a continuous domain (time, percentage of gait cycle, etc.). Since kinematic curves are measured in the same sample of individuals performing different activities, they are a clear example of functional data with repeated measures. On the other hand, the sample curves are observed with noise. Then, a roughness penalty might be necessary in order to provide a smooth estimation of the discriminant functions, which would make them more interpretable. Moreover, because of the infinite dimension of functional data, a reduction dimension technique should be considered. To solve these problems, we propose a multi-class approach for penalized functional partial least squares (FPLS) regression. Then linear discriminant analysis (LDA) will be performed on the estimated FPLS components. This methodology is motivated by two case studies. The first study considers the linear acceleration recorded every two seconds in 30 subjects, related to three different activities (walking, climbing stairs and down stairs). The second study works with the triaxial angular rotation, for each joint, in 51 children when they completed a cycle walking under three conditions (walking, carrying a backpack and pulling a trolley). A simulation study is also developed for comparing the performance of the proposed functional LDA with respect to the corresponding multivariate and non-penalized approaches.

Keywords

functional data Linear discriminant analysis multi-class classification PLS regression P-spline penalty

1 Introduction

Biomechanical gait data are commonly used to differentiate between several gait pathologies or different physical activities. Different statistical methodologies for classification of gait data have been developed in the literature with this aim. Gait data are usually collected at different points in a continuous scale, such as time, so that they are generated from an underlying time-varying function. In the applications that motivate this research, acceleration curves are measured at real time every two seconds, and rotation curves are measured in terms of the percentage of gait cycle. In spite of the continuous nature of kinematics curves, in many of the biomechanical studies the statistical analysis is made from a sample of summary measures for each curve such as minima, maxima, area under the curve, angles at heel strike, range of motion or the timing of specific events (Sadeghi et al., 1997; Sadeghi, 2003; Schmidt and Docherty, 2010; Orantes-González et al., 2015).

Currently, data science tools are being developed to build accurate classification and prediction models for gait research by using all the available information. One of them is functional data analysis (FDA) that begins by reconstructing the true functional form of each sample curve from discrete sampling points that could be different for different subjects. This way all the available information about the sample curves is exploited to explain the movement in the statistical analysis. A detailed description of the most common FDA methodologies, and interesting FDA applications with R and Matlab, are described in the pioneer books by Ramsay and Silverman (1997), Ramsay and Silverman (2002) and Ramsay et al.(2009). A complete study of nonparametric FDA techniques can be seen in Ferraty and Vieu (2006). Statistical inference for functional data has been studied in Horvath and Kokoszka (2012). A more recent book is by Kokoszka and Reimherr (2017), where new methods for FDA are collected (functional time series, spatial functional data, between others).

Biomechanical data are usually obtained over a number of discrete time points (snapshots) and assumed to be generated by some underlying smooth function. Because of this, the steps typically used when FDA is applied to the analysis of biomechanical data are: basis expansion representation, curve registration (time normalization or landmark registration) and functional principal component analysis (FPCA) (dimension reduction and variability explanation in terms of uncorrelated scalar variables). The first step consists of obtaining smooth representations of the sample curves in terms of a basis of functions (Fourier, B-splines, wavelets, etc.). In the second step, registration is used to reduce phase variability between curves while preserving the individual curve's shape and amplitude. The effect of different registration methods on cyclical kinematic data is analysed in Crane et al. (2010). In the third step, FPCA is successfully applied to explain the ways of variation in biomechanics of sports injuries and coordination in race-walking, jumping and running (Ryan et al., 2006; Harrison et al., 2007; Donoghue et al., 2008; Donà et al., 2009). On the other hand, multivariate principal component analysis (PCA) on the raw discrete-time observations of the biomechanics curves has also been applied to provide biomechanics interpretation of the principal component scores and to discriminate between healthy and sick subjects (Deluzio et al., 1997; Daffertshofer et al., 2004; Milovanovic and Popovic, 2012).

This work is motivated by two different sets of biomechanical data: the human activity data and the gait data. The first dataset is based on the human activity recognition collected by Anguita et al. (2013). The original experiments were carried out with a group of 30 volunteers within an age bracket of 19–48 years. Each person performed three different activities wearing a smartphone on the waist (walking, walking upstairs and walking downstairs). Using its embedded accelerometer and gyroscope, they captured 3-axial linear acceleration and 3-axial angular velocity at a constant rate of 50Hz every two seconds. The aim is to classify the curves of acceleration according to the categorical variable given by the three mentioned physical activities. The gait dataset comes from a wide experimental study developed in the biomechanics laboratories of the Sport and Health Institute of the University of Granada (iMUDS). In order to collect the data, 26 reflective markers were placed on the skin of a total of 51 children between 8 and 11 years old. The data were recorded by a 3D motion capture system meanwhile each subject completed a cycle walking over the platform in three conditions (walking, carrying a backpack that weighs 20% of the subject's weight and pulling a trolley that weighs 20% of the subject's weight). For each subject, the 3-axial angular rotation were registered for each join (ankle, foot progress, hip, knee, pelvis, thorax), and all conditions, in 101 equidistant points of the complete gait cycle. Now, the aim is to classify the rotation curves according to the three different load types.

Both case studies considered in this article are clear examples of functional data with repeated measures because the functional variables of interest (angular rotation and linear acceleration) are measured repeatedly on the same sample individuals under different experimental conditions. Therefore, the aim of this article is to classify a set of functional data with repeated measures (kinematic curves) according to a categorical variable with more than two categories (measurement conditions). In this article, an extension of linear discriminant analysis (LDA) to the case of functional data with repeated measurements will be considered. Alternative classification approaches based on logit regression were developed and successfully applied in different areas as environment, medicine, marketing and so on (James, 2002; Escabias et al., 2004, 2005, 2007; Delaigle and Hall, 2012a; Aguilera-Morillo and Aguilera, 2013; Escabias et al., 2014). An extension to functional data of nonparametric Bayes classifiers based on simple density ratios was proposed in (Dai et al., 2017) for the case of binary classification.

LDA is a popular and consolidate methodology for classification and dimension reduction (Fisher, 1936). So, LDA provides a sequence of linear combinations of the original predictor variables (linear discriminants) that maximize the between-class variance relative to the within-class variance. An important problem is that LDA overfits the data when we have a large number of highly correlated predictor variables because the within-class covariance matrix is difficult to estimate. This is the case of FDA that works with a large number of discrete-time observations for each sample curve. A way to solve this problem is to introduce some kind of regularization in the estimation of the covariance matrix (Friedman, 1989; Frank and Friedman, 1989; Hastie et al., 1995). A general review of regularized techniques in discriminant analysis was developed in Mkhadri et al. (1997). A Bayesian approach for Fisher's discriminant analysis of stochastic process was introduced in Shin, (2008). An alternative solution, for the case of irregularly sampled curves, is based on using as predictors a spline basis expansion of the discrete-time observations and assuming Gaussian distribution of the vector of basis coefficients with common covariance matrix for all classes, by analogy with LDA (James and Hastie, 2001). This functional LDA approach was used to assessment of embryonic growth (Bottomley et al., 2009). Other way of solution consists of projecting the predictor vectors (multivariate analysis) or curves (FDA) onto a finite dimensional basis so that LDA, or any other discrimination procedures, can then be used taking the resulting basis coefficients as predictors. A common approach is based on the representation of the original predictors by a small set of orthogonal scalar variables given by PCA or partial least squares (PLS) regression. Both approaches are commonly used in chemometrics to discriminate between essential properties of substances from discrete measurements of the near-infrared (NIR) spectra. PCA was applied in Sato (1994) to classify NIR spectral data of vegetable oils. PLS regression was performed in Li et al. (2000) to discriminate edible fats and oils by Fourier transform NIR spectroscopy. Multivariate classification techniques such as LDA, quadratic discriminant analysis (QDA) and logit regression were applied in Indahl et al. (1999). In the functional data context, reduction dimension techniques and regularization can be jointly used to get good classification rates and accurate interpretation of the results. The functional LDA approach presented in this article is inspired in this research line.

Functional PCA and functional PLS regression were introduced as natural extensions of their multivariate counterparts to solve the problems of high dimension and multicollinearity associated with the scalar-on-function linear model (Deville, 1974; Dauxois et al., 1982; Ocaña et al., 1999; 2007; Preda and Saporta, 2005; Aguilera et al., 2016). Both methodologies were compared on different simulated datasets concluding that they have similar forecasting performance, but the estimated parameter function provided by functional PLS regression is more accurate and needs fewer components (Reiss and Ogden, 2007; Aguilera et al., 2010; Delaigle and Hall, 2012b; Delaigle and Hall, 2012c). In addition, the performance of PCA for discrimination may not be optimal, because PCA only identifies gross variability, and then it is not capable to distinguish the between-groups and within-groups variability. In that sense, several authors proposed to use PLS regression for dimension reduction, as for example Barker and Rayens (2003), Pérez-Enciso and Tenenhaus (2003) in the multivariate context and Preda et al. (2007) in the functional context.

From a methodological point of view, the originality of our work lies, on the one hand, in proposing a penalized spline estimation of the PLS components to solve the problem of lack of smoothness, and on the other hand, in considering a multi-class approach for PLS regression based on the split-up variation to solve the problem of repeated measures. Then, classic LDA is carried out on the estimated PLS components. To the best of our knowledge, it is the first approach in the literature for functional linear discriminant analysis (FLDA) with repeated measures based on functional PLS. Then, in order to classify gait curves according to the kind of activity and take into account both, the lack of smoothness in the data and the correlation between repeated measurements, a new multi-class FLDA based on penalized functional PLS regression for repeated measures is proposed in Section 2. The classification results of the considered multi-class functional LDA-PLS approach on simulated data and both biomechanical datasets are presented and discussed in Sections 3 and 4. Conclusions are finally included in Section 5.

2 Multi-class classification of functional data

Let us consider a categorical response variable $Y = {1, \dots, K},$ with $K$ $(K > 2)$ categories and a functional predictor $X = {X (t) : t \in T}$ (continuous and second order stochastic process whose sample paths belong to the Hilbert space of square-integrable functions $L_{2} (T)$ ). Then, the aim is to predict $Y$ from $X (t) .$ That means to classify the sample curves according to the groups defined by the categories of $Y .$

The aim of FLDA is to find linear combinations

Φ_{i} (X) = α_{i} + \int_{T} X (t) β^{i} (t) dt, {i = 1, \dots, K - 1},

(2.1)

so that the between-class variance is maximized with respect to the total variance

\max_{β} \frac{V (E [Φ (X) | Y])}{V (Φ (X))} .

Due to the high dimension of $X (t)$ (the number of observed variables is usually larger than the number of sample observations), the estimation of the discriminant functions is an ill-posed problem because the standard estimate for the within-class covariance matrix is singular so that it is impossible to apply the usual classification rule. This problem is solved in this article by using a dimension reduction approach based on functional PLS regression Preda et al. (2007). Then, the multi-class classification problem from repeated measures is solved by using the following two-step algorithm:

Step 1. Dimension reduction by penalized and non-penalized multi-class functional PLS regression (for repeated measures) of the random vector $\tilde{Y} = ({\tilde{Y}}_{1}, \dots, {\tilde{Y}}_{K - 1})$ on the functional predictor $X (t)$ , with ${{\tilde{Y}}_{i} \in {0, 1} : i = 1, \dots, K - 1}$ being the dummy variables associated with the categorical response Y.

Step 2. Functional LDA of Y on a reduced set of functional PLS components estimated in the above step.

2.1 Multi-class functional PLS regression

In order to reduce the dimension, PLS regression of the vector $\tilde{Y} = ({\tilde{Y}}_{1}, \dots, {\tilde{Y}}_{K - 1})$ on the functional predictor $X (t),$ with ${{\tilde{Y}}_{i} \in {0, 1} : i = 1, \dots, K - 1}$ is performed. In that context, the non-penalized PLS components are generalized linear combinations, of type $t = \int_{T} X (t) w (t) dt,$ estimated by solving the following maximization problem

\begin{matrix} {w, c} & = & argmax {Cov}^{2} (\int_{T} X (t) w (t) dt, \sum_{i = 1}^{K - 1} c_{i} {\tilde{Y}}_{i}) \\ ∥ w ∥ = 1 \\ ∥ c ∥ = 1 \end{matrix}

(2.2)

with $∥ \cdot ∥$ representing the usual norms in the spaces $L_{2} [T]$ and $^{K - 1}$ where the component weights belong to, respectively, and $c$ being the weight vector of $\tilde{Y}$ .

From Proposition 1 in Preda and Saporta (2005), it can be seen that the solution to this maximization problem is reached for $w$ and $c,$ the eigenvectors associated with the largest eigenvalue of $U_{X}$ and $U_{\tilde{Y}},$ respectively, defined by $U_{X} = C_{X \tilde{Y}} \circ C_{\tilde{Y} X}$ and $U_{\tilde{Y}} = C_{\tilde{Y} X} \circ C_{X \tilde{Y}},$ with $C_{\tilde{Y} X}$ and $C_{X \tilde{Y}}$ being the cross-covariance operators.

Then, the first PLS component, $t_{1} = \int_{T} X (t) w_{1} (t) dt,$ is estimated in terms of the associated weight function, $w_{1} (t),$ which is solution to the following problem

C_{X \tilde{Y}} \circ C_{\tilde{Y} X} (w_{1}) = λ_{max} w_{1} .

PLS regression is an iterative algorithm and the first PLS step is completed by ordinary linear regression of $X_{0} (t) = X (t)$ and ${\tilde{Y}}_{0} = \tilde{Y}$ on $t_{1},$ obtaining the residuals, $X_{1} (t)$ and ${\tilde{Y}}_{1},$ necessary to the following step

\begin{matrix} X_{1} (t) & = & X_{0} (t) - p_{1} (t) t_{1} t \in T \\ {\tilde{Y}}_{1, i} & = & {\tilde{Y}}_{0, i} - c_{1, i} t_{1} i = 1, \dots, K - 1, \end{matrix}

where $p_{1} (t) = ((X_{0} (t) t_{1}) / (t_{1}^{2}))$ and $c_{1} = ((Y_{0} t_{1}) / (t_{1}^{2})) .$

In general, let us consider the hth PLS component, $t_{h} = \int_{T} X_{h - 1} (t) w_{h} (t) dt,$ which is estimated by the following problem

C_{X \tilde{Y}}^{h - 1} \circ C_{\tilde{Y} X}^{h - 1} (w_{h}) = λ_{max} w_{h},

where $C_{X \tilde{Y}}^{h - 1}$ and $C_{\tilde{Y} X}^{h - 1}$ are the cross covariance operators of $X_{h - 1} (t)$ and ${\tilde{Y}}_{h - 1}$ , respectively. The hth PLS step is completed by ordinary linear regression of $X_{h - 1} (t)$ and ${\tilde{Y}}_{h - 1}$ on $t_{h},$ obtaining the residuals $X_{h} (t)$ and ${\tilde{Y}}_{h}$

\begin{matrix} X_{h} (t) & = & X_{h - 1} (t) - p_{h} (t) t_{h} t \in T \\ {\tilde{Y}}_{h, i} & = & {\tilde{Y}}_{h - 1, i} - c_{h, i} t_{h} i = 1, \dots, K - 1, \end{matrix}

where $p_{h} (t) = ((X_{h - 1} (t) t_{h}) / (t_{h}^{2}))$ and $c_{h} = ((Y_{h - 1} t_{h}) / (t_{h}^{2})) .$

A basis representation approach

Let us consider now the basis representation of the functional predictor and the weight function, that is,

X (t) = \sum_{j = 1}^{p} α_{j} ϕ_{j} (t) and w (t) = \sum_{j = 1}^{p} w_{j} ϕ_{j} (t) .

As a consequence, the maximization problem in Equation (2.2) can be written as follows:

{w, c} = argmax \frac{w^{'} Ψ Σ_{α \tilde{Y}} c c^{'} Σ_{α \tilde{Y}}^{'} Ψ w}{w^{'} Ψ w + c^{'} c},

with $w = {(w_{1}, \dots, w_{p})}^{'}$ being the vector of basis coefficients of $w (t)$ and $Σ_{α \tilde{Y}} = (σ_{ji})_{p \times K - 1}$ being the cross-covariance matrix between $α$ (the vector of basis coefficients of $X (t)$ ) and $\tilde{Y},$ where $σ_{ji} = E [α_{j} {\tilde{Y}}_{i}],$ with ${j = 1 \dots, p; i = 1, \dots, K - 1}$ .

The cross-covariance operators expressed in terms of the basis expansion of the functional predictor $X (t)$ are given by

\begin{matrix} C_{\tilde{Y} X} : & L^{2} (T) \to^{K - 1} \\ f = \sum_{j = 1}^{p} f_{j} ϕ_{j} (t) \to C_{\tilde{Y} X} (f) = Σ_{α \tilde{Y}}^{'} Ψ f \\ C_{X \tilde{Y}} : & ^{K - 1} \to L^{2} (T) \\ x \to f (t) = \sum_{j = 1}^{K - 1} f_{j} ϕ_{j} (t) : f = Σ_{α \tilde{Y}} x, \end{matrix}

where $f = (f_{1}, \dots, f_{p})^{'}$ is the vector of basis coefficients of the function $f .$

Then, at the first PLS step, the weight function $w_{1} (t)$ is given by the eigenfunction associated with the largest eigenvalue of $U_{X}$ so that

Σ_{α \tilde{Y}} Σ_{α \tilde{Y}}^{'} Ψ w_{1} = λ w_{1},

(2.3)

and the weight vector $c_{1}$ is given by the eigenvector associated with the largest eigenvalue of $U_{\tilde{Y}}$ so that

Σ_{α \tilde{Y}}^{'} Ψ Σ_{α \tilde{Y}} c_{1} = λ c_{1} .

Let us consider now the decomposition $Ψ = (Ψ^{1 / 2}) {(Ψ^{1 / 2})}^{'}$ . Then,

w_{1}^{'} Ψ w_{1} = w_{1}^{'} (Ψ^{1 / 2}) {(Ψ^{1 / 2})}^{'} w_{1} = {\tilde{w}}_{1}^{'} {\tilde{w}}_{1},

with ${\tilde{w}}_{1} = {(Ψ^{1 / 2})}^{'} w_{1}$ $(w_{1} = (Ψ^{- 1 / 2})^{'} {\tilde{w}}_{1}) .$ This way, Equation (2.3) can be expressed as follows

{(Ψ^{1 / 2})}^{'} Σ_{α \tilde{Y}} Σ_{α \tilde{Y}}^{'} Ψ^{1 / 2} \tilde{w_{1}} = λ \tilde{w_{1}} .

(2.4)

At the hth PLS step, the hth PLS component is estimated in terms of the associated weight function by solving to the following problem

{(Ψ^{1 / 2})}^{'} Σ_{α_{h - 1} {\tilde{Y}}_{h - 1}} Σ_{α_{h - 1} {\tilde{Y}}_{h - 1}}^{'} Ψ^{1 / 2} \tilde{w_{h}} = λ \tilde{w_{h}},

(2.5)

where $Σ_{α_{h - 1} {\tilde{Y}}_{h - 1}}$ is the cross-covariance matrix between $α_{h - 1}$ (the vector of basis coefficients of $X_{h - 1} (t)$ ) and ${\tilde{Y}}_{h - 1},$ and ${\tilde{w}}_{h} = {(Ψ^{1 / 2})}^{'} w_{h}$ $(w_{h} = (Ψ^{- 1 / 2})^{'} {\tilde{w}}_{h}) .$

In general, by considering Equations (2.4) and (2.5) and taking into account that ${(Ψ^{1 / 2})}^{'} Σ_{α \tilde{Y}}$ is the cross-covariance matrix between $\tilde{Y}$ and $(Ψ^{1 / 2})^{'} α$ , it can be concluded that FPLS with multiple response is equivalent to ordinary PLS of $\tilde{Y}$ on the random vector $(Ψ^{1 / 2})^{'} α .$

2.2 Penalized multi-class FPLS regression

As the biomechanical data analysed in this article are not smooth enough, a penalized approach for multi-class functional partial least squares (FPLS) regression is proposed.

Let us consider now the roughness penalty function

{Pen}_{d} (w) = \int_{T} [D^{d} (w) (t)]^{2} dt,

with $d = 2$ , as a measure of roughness of $w$ (O'Sullivan, 1986). By considering the basis expansion of $w (t),$ the penalty function can be written as ${Pen}_{d} (w) = w^{'} P_{d} w,$ where $w = {(w_{1}, \dots, w_{p})}^{'}$ is the vector of basis coefficients of $w (t)$ and $P_{d}$ the matrix of the cross inner products of the d-order derivatives of the basis functions. As alternative, working in any context where regression on B-splines is useful, it can be considered $P_{d} = {(▵^{d})}^{'} ▵^{d},$ with $▵^{d}$ the matrix of d-order differences between the adjacent basis coefficients (Eilers and Marx, 1996).

In this section the PLS components are estimated by solving the following problem

\begin{matrix} {w, c} & = & argmax {Cov}^{2} (\int_{T} X (t) w (t) dt, \sum_{i = 1}^{K - 1} c_{i} {\tilde{Y}}_{i}), \\ ∥ w ∥_{λ} = 1 \\ ∥ c ∥ = 1 \end{matrix}

with $∥ \cdot ∥$ being the norm in the space $^{K - 1}$ and $∥ \cdot ∥_{λ}$ being the norm associated with an inner product defined by $⟨ f, g ⟩_{λ} = ⟨ f, g ⟩ + λ (f^{'} P_{d} g),$ with $f = {(f_{1}, \dots, f_{p})}^{'}$ and $g = {(g_{1}, \dots, g_{p})}^{'}$ being the vectors of basis coefficients of $f (t)$ and $g (t),$ respectively.

Then, the hth PLS component ( $h > 1$ ), $t_{h} = \int_{T} X_{h - 1} (t) w_{h} (t) dt,$ is estimated by solving the problem

\begin{matrix} {w, c} & = & argmax {Cov}^{2} (\int_{T} X_{h - 1} (t) w (t) dt, \sum_{i = 1}^{K - 1} c_{i} {\tilde{Y}}_{h - 1, i}) . \\ ∥ w ∥_{λ} = 1 \\ ∥ c ∥ = 1 \end{matrix}

(2.6)

By assuming the basis representation of the functional predictor and the weight function, the problem in (2.6) can be written as follows

{w, c} = argmax \frac{w^{'} Ψ Σ_{α_{h - 1} {\tilde{Y}}_{h - 1}} c c^{'} Σ_{α_{h - 1} {\tilde{Y}}_{h - 1}}^{'} Ψ w}{w^{'} (Ψ + λ P_{d}) w + c^{'} c},

(2.7)

with $w = {(w_{1}, \dots, w_{p})}^{'}$ being the vector of basis coefficients of $w (t),$ $Σ_{α_{h - 1} {\tilde{Y}}_{h - 1}} = (σ_{ji})_{p \times K - 1}$ being the cross-covariance matrix between $α_{h - 1}$ and ${\tilde{Y}}_{h - 1},$ where $σ_{ji} = E [α_{h - 1, j} {\tilde{Y}}_{i}],$ with ${j = 1 \dots, p; i = 1, \dots, K - 1},$ and $λ$ and $P_{d}$ the smoothing parameter and the penalty matrix, respectively.

Assuming the decomposition ${LL}^{'} = Ψ + λ P_{d}$ and defining $\tilde{w} = L^{'} w$ $(w = (L^{- 1})^{'} \tilde{w})$ the problem in (2.7) can be expressed as follows

{w, c} = argmax \frac{{\tilde{w}}^{'} L^{- 1} Ψ Σ_{α \tilde{Y}} c c^{'} Σ_{α \tilde{Y}}^{'} Ψ (L^{- 1})^{'} \tilde{w}}{{\tilde{w}}^{'} \tilde{w} + c^{'} c} .

By analogy with the non-penalized approach, the associated eigenproblem is

L^{- 1} Ψ Σ_{α \tilde{Y}} Ψ^{'} (L^{- 1})^{'} \tilde{w} = λ \tilde{w}, \tilde{w} \in^{p} .

(2.8)

In general, the hth PLS component, $t_{h} = \int_{T} X (t) w_{h} (t) dt,$ is defined by $w_{h} = (L^{- 1})^{'} {\tilde{w}}_{h},$ which is the vector of basis coefficients of $w_{h} (t),$ with ${\tilde{w}}_{h}$ being the eigenvector related to the largest eigenvalue of the problem in (2.8).

Finally, taking into account that $L^{- 1} Ψ Σ_{α \tilde{Y}}$ is the cross-covariance matrix between $\tilde{Y}$ and $L^{- 1} Ψ α,$ the penalized approach is reduced to a classical PLS of $\tilde{Y}$ on the random vector $L^{- 1} Ψ α .$

2.3 The problem of repeated measures

In previous sections, penalized and non-penalized approaches for the multi-class FPLS regression have been proposed. However, these approaches are based on ordinary PLS, and then the between-subject and the within-subject variations are not studied by separate.

As solution to the problem of repeated measures, in this section a multi-class approach for penalized and non-penalized FPLS regression based on the split-up variation is proposed. In the same spirit as in de Nord and Theobald (2005), Westerhuis et al. (2010) and Liquet et al. (2012), this approach first decomposes the variability in the data matrix and then applies the multi-class PLS regression on the within-subject variation matrix.

For a more detailed explanation of the multi-class approach, let us consider ${x_{i} (t) : t \in T, i = 1, \dots, n}$ a sample of the functional variable $X (t)$ and ${y_{1}, y_{2}, \dots, y_{n}}$ a random sample of $Y$ associated with it, with $y_{i} = {1, 2, \dots, K}$ . Henceforth, let $\tilde{Y}$ denote a $N \times (K - 1)$ matrix with zeros and ones, comprising the values associated with the dummy variables obtained from $Y$ .

Then, the sample estimation of the multi-class FPLS described in sections 2.1 (non-penalized) and 2.2 (penalized) is as follows:

Non-penalized multi-class FPLS: multi-class PLS of $\tilde{Y}$ on the matrix $X_{(N \times p)}^{1} = A Ψ^{1 / 2}$

Penalized multi-class FPLS: multi-class PLS of $\tilde{Y}$ on the matrix $X_{(N \times p)}^{2} = A Ψ (L^{- 1})^{'},$ so that ${LL}^{'} = Ψ + λ P_{d},$

with $A$ being the matrix of basis coefficients of the functional variable $X,$ and $Ψ$ and $P_{d}$ as in previous sections.

In practice, functional data are observed in a finite set of points. Because of this, matrix $A$ must be estimated. By assuming that the underline process to biomechanical data is smooth, B-spline basis functions are considered, and the basis coefficients are estimated by the least squares criterion. This approach is called regression splines. Additionally, for data observed with noise, a penalized estimation of the basis coefficients in terms of P-splines can be considered (see Eilers and Marx, 1996; Aguilera and Aguilera-Morillo, 2013, for more details).

Split-up variation

In general, let us consider a data matrix $X_{N \times p}$ containing the information related to $n$ subjects and each of the $K$ possible stimulus ( $N = n \times K$ ), with $p$ being the number of basis functions (or variables in the multivariate setting). In that sense, let $x_{i k}^{j}$ denote the value of the jth column of $X$ for the ith subject under the kth stimulus.

Following the mixed-effect models philosophy, matrix $X$ can be decomposed into an offset term, a between-subject and a within-subject part so that

X = X_{offset} + X_{between - subject} + X_{within - subject} .

But mixed models rely on assumptions such as Gaussian distribution of random effects. As an alternative, a split-up variation approach, which does not require the above-mentioned assumptions, is considered in this article. Exactly, the offset, between-subject and within-subject variation matrices ( $X_{o},$ $X_{b}$ and $X_{w},$ respectively) can be obtained as follows: $X_{o} = 1_{(N \times 1)} x_{\cdot \cdot (1 \times p)}^{'},$ where $x_{\cdot \cdot} = (x_{\cdot \cdot}^{1}, \dots, x_{\cdot \cdot}^{p})^{'}$ and $x_{\cdot \cdot}^{j} = \frac{1}{N} \sum_{k = 1}^{K} \sum_{i = 1}^{n} x_{i k}^{j} .$ $X_{b}$ is obtained by concatenating the matrices $1_{(K \times 1)} x_{b i (1 \times p)}^{'}$ for each subject into $X_{b},$ where $x_{b i} = (x_{i \cdot}^{1} - x_{\cdot \cdot}^{1}, \dots, x_{i \cdot}^{p} - x_{\cdot \cdot}^{p})^{'}$ and $x_{i \cdot}^{j} = \frac{1}{K} \sum_{k = 1}^{K} x_{i k}^{j} .$ $X_{w} = X - X_{i \cdot},$ with $X_{i \cdot}$ obtained by concatenating the matrices $1_{(K \times 1)} x_{i \cdot (1 \times p)}^{'}$ for each subject into $X_{i \cdot},$ with $x_{i \cdot} = (x_{i \cdot}^{1}, \dots, x_{i \cdot}^{p})^{'} .$

Multi-class approach for FPLS

Once the split-up variation has been carried out, penalized and non-penalized multi-class FPLS regression is performed on the within-subject variation matrix. In that sense, the penalized and non-penalized multi-class approaches for FPLS regression are given by

Non-penalized multi-class FPLS: multi-class PLS of $\tilde{Y}$ on the matrix $X_{w}^{1}$

Penalized multi-class FPLS: multi-class PLS of $\tilde{Y}$ on the matrix $X_{w}^{2},$

where $X_{w}^{1}$ and $X_{w}^{2}$ are the within-subject variation matrices from $X^{1}$ and $X^{2},$ respectively.

2.4 Functional LDA on a reduced set of functional PLS components

Let $T$ be a $n \times q$ matrix comprising the columns of the first $q$ PLS component scores (by variability order). By considering the sample estimation of the penalized and non-penalized multi-class approaches for FPLS regression proposed in Section 2.3, the sample estimation of $T$ is as follows

Non-penalized multi-class FPLS: $T = {A Ψ}^{1 / 2} V$

Penalized multi-class FPLS: $T = A Ψ (L^{- 1})^{'} V,$

with $V$ being the matrix comprising the columns of the eigenvectors ${\tilde{w}}_{1}, \dots, {\tilde{w}}_{q}$ associated with the $t_{1}, \dots, t_{q}$ PLS components.

Once the PLS components have been estimated, LDA of the original sample values of $Y$ $(y_{1}, \dots, y_{n})$ on the matrix $T$ is carried out. Then, linear combinations such as $Φ_{i} (X) = α_{i} + T β_{PLS}^{i},$ ${i = 1, \dots, K - 1}$ are estimated, so that the between-class variance is maximized with respect to the within-class variance, with $β_{PLS}^{i} = (β_{{PLS}_{1}}^{i}, \dots, β_{{PLS}_{q}}^{i})^{'}$ being the vector of discriminant coefficients.

But the aim is to estimate the discriminant functions $β^{i} (t),$ $i = 1, \dots, K - 1$ included in the linear combinations given in equation (2.1). To this end, let us consider the basis representation of the sample curves and the discriminant functions $(β^{i} (t) = \sum_{j = 1}^{p} β_{j}^{i} ϕ_{j} (t))$ and, as a consequence, equation (2.1) can be rewritten as

Φ_{i} (X) = α_{i} + \int_{0}^{T} X (t) β^{i} (t) dt = α_{i} + {A Ψ β}^{i},

with $β^{i}$ being the vector of basis coefficients of the ith discriminant function.

Finally, taking into account the sample estimation of $T,$ we can conclude that

Non-penalized multi-class FPLS: $β^{i} = (Ψ^{- 1 / 2})^{'} V β_{PLS}^{i}, i = 1, \dots, K - 1$

Penalized multi-class FPLS: $β^{i} = (L^{- 1})^{'} V β_{PLS}^{i}, i = 1, \dots, K - 1 .$

Once the random vector of discriminant functions $β (t) = (β^{1} (t), \dots, β^{K - 1} (t))^{'}$ is obtained, the classification is then performed in the transformed space based on some distance metric, such as Euclidean distance.

Then given a new sample observation $x_{0} (t),$ it is classified to

\begin{matrix} \underset{k}{argmin} & d (\int_{T} x_{0} (t) β (t), \int_{T} {\overset{̅}{x}}^{k} (t) β (t)), \end{matrix}

(2.9)

where ${\overset{̅}{x}}^{k} (t) = \frac{1}{m_{k}} \sum_{x (t) \in C_{k}} x (t)$ is the mean curve related to the kth class and $m_{k}$ the number of sample observations in the kth class, with $k = 1, \dots, K .$

By assuming the basis representations of $x_{0} (t),$ $β (t)$ and ${\overset{̅}{x}}^{k} (t)$

\begin{matrix} x_{0} (t) = \sum_{j = 1}^{p} a_{0_{j}} ϕ_{j} (t), & β^{i} (t) = \sum_{j = 1}^{p} β_{j}^{i} ϕ_{j} (t), & {\overset{̅}{x}}^{k} (t) = \sum_{j = 1}^{p} {\overset{̅}{a}}_{j}^{k} ϕ_{j} (t), \end{matrix}

problem (2.9) can be rewritten as follows

\begin{matrix} \underset{k}{argmin} & d (a_{0}^{'} Ψ B, {\overset{̅}{a}}^{k^{'}} Ψ B), \end{matrix}

with $a_{0_{p \times 1}}$ being the vector of basis coefficients of $x_{0} (t),$ ${\overset{̅}{a}}_{p \times 1}^{k}$ being the vector of basis coefficients of ${\overset{̅}{x}}^{k} (t),$ $Ψ_{p \times p}$ being the inner product matrix, $B_{p \times (K - 1)}$ the matrix with columns the vectors of basis coefficients of the $K - 1$ discriminant functions (i.e., $β^{1}, \dots, β^{K - 1}$ ). Finally, the Euclidean distance between two vectors $V^{1} = a_{0^{'}} Ψ B$ and $V^{2} = {\overset{̅}{a}}^{k^{'}} Ψ B$ is computed as $d (V^{1}, V^{2}) = \sqrt{\sum_{i = 1}^{K - 1} (V_{i}^{1} - V_{i}^{2})^{2}} .$

The smoothing parameter $λ$ and the number of PLS components used in the LDA approach based on Penalized multi-class FPLS have been selected by leave-one-out cross-validation on a training sample, selecting the combination of both that minimizes the misclassification rate obtained from LDA on the PLS components. In the non-penalized approach, only the number of PLS components is required, and it is chosen by leave-one-out cross-validation on a training sample, minimizing the misclassification rate obtained from functional LDA.

All the analysis in Sections 3 and 4 were carried out in a computer with an Intel Core i5 processor at 3 GHz with 8Gb of RAM, running R-project version 3.5.1 (R Development Core Team, 2008). The two multi-class functional PLS approaches were reduced to multi-class PLS regression (in the multivariate context), and then they were implemented in R using the function spls available in the mixOmics R package developed by Rohart et al. (2017).

Authors are preparing an R package that, among other things, includes the code used in the manuscript. For readers wishing to use the software, it will be available as soon as possible at the R repository website. However, you can contact the corresponding author if you need more information regarding the code.

3 Simulation study

In order to test the performance of the proposed methods, a simulation study has been carried out. Inspired by Cuevas et al.(2004), let us consider $120$ sample curves, $40$ for each of the 3 classes, according to the following process

X_{i}^{k} (t) = m_{k} (t) + ε_{i}^{k} (t), t \in [0, 1], i = 1, \dots, 40, k = 1, 2, 3,

where $m_{i} (t) = t^{(k / 5)} (1 - t)^{(6 - k / 5)} .$ The process $X_{i}^{k} (t)$ was generated in discretized version $X_{i}^{k} (t_{j}),$ for $j = 1, \dots, 101$ , with $t_{1}, \dots, t_{101}$ being equidistant points in $[0, 1] .$ Finally, $ε_{i}^{k} (t_{j})$ are i.i.d. random variables $N (0, σ_{ε} = 0.2) .$

But our work is focused on the presence of repeated measurements in the $k$ classes, and then an additional term must be considered in order to represent the subject effect. In that sense, and inspired by Durbán et al. (2005), a specific-subject term has been included in the above process so that

X_{i}^{k} (t) = m_{k} (t) + a_{i} sin (t π) + ε_{i}^{k} (t), t \in [0, 1], i = 1, \dots, 40, k = 1, 2, 3,

with $a_{i} \sim N (μ_{i}, σ_{s} = 0.02)$ and $μ_{i} \sim Uniform (0, 0.05) .$ The $sin (t π)$ function has been considered in order to simulate curves taking the same value at the beginning and at the end of the observation interval (as in the gait data study). Of the 120 simulated sample curves, 90 (30 for each class) were considered as training sample and the remaining 30 (10 for each class) as test sample.

The sample curves from two different subjects, and for each of the three classes, are shown in Figure 1. As we can see, there is a clear dependency structure between the three observations of the same subject, and the between-subjects variation is not only given by a magnitude effect, but also by a shape variation. This is a realistic simulation of the performance of kinematic data since it makes no sense to assume that the curves of the different subjects are parallel.

Figure 1:

Simulation study. Sample paths related to classes 1, 2 and 3 (solid, dashed and dotted line, respectively) for two subjects A and B, black and grey lines, respectively

In Figure 2, all the sample curves, with and without noise, are displayed highlighting the three classes by different line types. In order to fit the regression splines, a cubic B-spline basis defined on 15 basis knots has been considered.

Figure 2:

Simulation study. Smooth and noisy sample paths (left and right panel, respectively) related to classes 1, 2 and 3 (solid line, dashed line and dotted line, respectively) for 30 subjects

The simulation scheme was run 500 times and the main results are summarized in Figure 3. As we can see, the two functional approaches tend to select a larger number of PLS components than the multivariate version, as covariates in the LDA. But there is not a great difference between them. Regarding the correct classification rates (CCR), there is a clear overfitting in the multivariate approach, which provides CCR close to 1 in the cross-validation. Between the two proposed functional approaches, the one based on penalized multi-class functional PLS regression provides the best performance, with the highest CCR and the minimum variability.

Figure 3:

Simulation study. Box plots showing to the number of PLS components used in the LDA (bottom panel), the correct classification rates from the cross-validation on the training sample and from the classification of the test sample (top-left and top-right panels, respectively). Experiments run 500 times

4 Case studies

In this section, the results from the analysis of two real kinematic datasets are summarized. In general, kinematic data from inertial sensors placed on human body is usually collected in the format of linear acceleration or angular rotation. Linear acceleration or angular rotation represent the trajectory of the human movement, which can be identified with a function in time and frequency domain. In that sense, data from kinematic studies can be considered as a functional dataset, and then functional data techniques can be applied to analyse this type of data. In addition, it is important to remark the presence of repeated measures in kinematic studies. So, both the between-subject and within-subject variations must be taken into account.

4.1 Human activity data

The human activity dataset is part of a wider experiment, focused on the human activity recognition, carried out by Anguita et al. (2013). In this article, let us consider the linear acceleration (metre per second squared), measured on axis X and recorded in 128 equidistant knots at the interval $[0, 2.56]$ seconds, and related to three different stimulus (walking, walking upstairs and walking downstairs). A total of 30 subjects participated in this experiment. Then, a total of 90 sample paths are considered (one sample path per stimulus and subject). Subject 23 was removed from the sample for being an outlier. The remaining sample observations were distributed into training and test sample, with observations (repeated measures) related to 20 and 9 subjects, respectively. The raw data (displayed by stimulus) are shown in Figure 4. In order to visualize the variability of the spectra between stimulus, in Figure 5 the spectra related to the three activities have been overlapped for two subjects A and B, left and right panel, respectively.

Figure 4:

Human activity data. Raw data. Sample paths displayed separately by stimulus: walking, walking upstairs and walking downstairs, from left to right, respectively

The first step in the analysis was to approximate the sample curves by mean of a basis representation and using a cubic B-spline basis defined on 25 equidistant knots. In Figure 6 the raw data together with the regression spline can be seen. Once the sample curves have been approximated, the aim is the human activity pattern recognition. To this end FLDA based on a multi-class approach for functional PLS regression (FLDA-MFPLS) was carried out to classify the sample curves according to the stimulus which produced them. It could be assumed that the raw data, which were collected by a smartphone on the waist of the subjects, are affected by some error or noise. In that sense, and aiming to avoid a possible lack of smoothness in the estimation of the discriminant functions, the penalized version of FLDA-MFPLS (FLDA-PenMFPLS), proposed in section 2.3, was also considered.

The discriminant functions estimated from penalized and non-penalized functional LDA can be seen in Figure 7. The necessity of smoothing in the estimation of the discriminant functions is obvious, with penalized functional LDA providing the smoothest and the most interpretable functions.

Figure 5:

Human activity data. Raw data. Sample paths related to walking (solid line), walking upstairs (dashed line) and walking downstairs (dotted line) for two subjects A and B, left and right panel, respectively

Figure 6:

Human activity data. A sample path for each stimulus (dashed line) together with the regression spline (solid line) estimated with a cubic B-splines basis defined on 25 equidistant knots

In order to evaluate the ability of penalized and non-penalized approaches for human activity pattern recognition, Table 1 summarizes information related to the CCR based on a leave-one-out cross-validation on the training sample ( ${CCR}_{training}$ ), and based on the test sample classification ( ${CCR}_{test}$ ). In this table, a multivariate version of LDA based on multi-class PLS regression (LDA-MPLS) is also considered for comparison purpose.

It can be seen that functional LDA performs better than the multivariate approach. Besides, between the two functional approaches (penalized and non-penalized), the one based on penalized multi-class functional PLS (FLDA-PenMFPLS) achieves the highest CCR in both, training and test samples. A considerable difference between the CCR in training and test samples highlights a possible overfitting mainly in the multivariate and in the non-penalized functional approaches. In the penalized approach this fact is not so remarkable, providing then the best classification rate.

Figure 7:

Human activity data. Discriminant functions estimated by functional LDA based on non-penalized (dashed line) and penalized (solid line) multi-class functional PLS regression, using 4 and 7 PLS components, respectively

Table 1:

Human activity data. Correct classification rate based on a leave-one-out cross-validation on the training sample ( ${CCR}_{training}$ ) and based on the test sample classification ( ${CCR}_{test}$ ), and the number of PLS components ( $# PLS {cp}^{'} s$ ) used in the functional LDA based on penalized and non-penalized multi-class functional PLS regression (FLDA-PenMFPLS and FLDA-MFPLS, respectively). $λ$ denotes the smoothing parameter of the penalized approach. A multivariate version of LDA based on multi-class PLS regression (LDA-MPLS) is also considered for comparison purpose

	LDA-MPLS	FLDA-MFPLS	FLDA-PenMFPLS
${CCR}_{training}$ (%)	83	95	93
${CCR}_{test}$ (%)	56	67	78
$# PLS {cp}^{'} s$	5	4	7
$λ$	–	–	0.41

4.2 Gait data

The gait dataset comes from a wide experimental study developed in the biomechanics laboratories of the iMUDS. A total of 51 participants (25 boys and 28 girls) between 8 and 11 years old were involved in this study. In order to collect the data, 26 reflective markers were placed on the children's skin. The kinematics data was recorded by a 3D motion capture system (Qualisys AB, Göteborg, Sweden). Each subject completed a cycle walking over the platform in three conditions (walking, carrying a backpack that weighs 20% of the subject's weight and pulling a trolley that weighs 20% of the subject's weight). For each subject, the 3-axial angular rotation were registered for each join (ankle, foot progress, hip, knee, pelvis, thorax) in all conditions. Finally, the angular rotation curve represents the observation for the subject's gait cycle for each join, axis direction and experimental condition. Interested readers could request access to the data from Jose M. Heredia-Jiménez and Eva Orantes-González from iMUDS.

In this section, only a part of the above-described experimental dataset has been considered. Exactly, we are interested on the thorax angular position (radians) measured on axis Z. This is a clear example of functional data with repeated measures in the sense that three curves, concerning the three experimental conditions, are available for each subject. In addition, each sample curve was recorded in 101 equidistant points of the gait cycle. So, the interval $[0, 100]$ represents the percentage of gait cycle completed by the subject.

Subject 10 was removed from the sample for being an outlier. The remaining sample observations were distributed into training and test samples, with observations (repeated measures) related to 39 and 11 subjects, respectively. The raw data (displayed by experimental conditions) are shown in Figure 8.

Figure 8:

Gait data. Raw data. Sample paths displayed separately by type of activity: walking, carrying a backpack and pulling a trolley, both weighing 20% of the children's weight, from left to right, respectively

In order to visualize the within subjects and the between subjects variability, in Figure 9 the spectra related to the three gait events have been overlapped for two subjects A and B, left and right panel, respectively. It can be seen that walking and walking with a trolley provide sample paths with a more similar shape than walking with a backpack. Also it is interesting to highlight the between subjects variability (differences in scale and shape).

Figure 9:

Gait data. Raw data. Sample paths related to walking (solid line), carrying a backpack (dashed line) and pulling a trolley (dotted line), both weighing 20% of the children's weight and walking downstairs from two sample individuals A and B, left and right panel, respectively

The sample curves have been approximated by mean of a basis representation using a cubic B-spline basis defined on 20 equidistant knots. In Figure 10 the raw data and the regression splines have been displayed. Once the sample curves have been approximated, the aim is the classification of the thorax angular rotation curves. In that sense, the discriminant functions estimated by penalized and non-penalized functional LDA can be seen in Figure 11. Once again the penalized approach provides the smoothest functions, providing an intuitive interpretation in relation with detecting periods in which the function is positive or negative.

Figure 10:

Gait data. A sample path for each stimulus (dashed line) together with the regression spline (solid line) estimated with a cubic B-splines basis defined on 20 equidistant knots

Figure 11:

Gait data. Discriminant functions estimated by functional LDA based on non-penalized (dashed line) and penalized (solid line) multi-class functional PLS regression, using 5 and 3 PLS components, respectively

The classification performance of the penalized and non-penalized approaches is shown in Table 2. It obvious that functional LDA performs better than the multivariate approach. Between the two functional approaches, the one based on penalized multi-class functional PLS (FLDA-PenMFPLS) achieves the best classification performance (the highest CCR on the test sample).

Table 2:

Gait data. Correct classification rate based on leave-one-out cross-validation on the training sample ( ${CCR}_{training}$ ) and based on the test sample classification ( ${CCR}_{test}$ ), and the number of PLS components ( $# PLS {cp}^{'} s$ ) used in the functional LDA based on penalized and non-penalized multi-class functional PLS regression (FLDA-PenMFPLS and FLDA-MFPLS, respectively). $λ$ denotes the smoothing parameter of the penalized approach. A multivariate version of LDA based on multi-class PLS regression (LDA-MPLS) is also considered for comparison purpose

	LDA-MPLS	FLDA-MFPLS	FLDA-PenMFPLS
${CCR}_{training}$ (%)	83	87	82
${CCR}_{test}$ (%)	33	61	85
$# PLS {cp}^{'} s$	3	5	3
$λ$	–	–	2.87

5 Conclusions

In this work, a methodological solution to the problem of multi-class classification of functional data with repeated measures has been proposed. FLDA has been considered as a classifier, by solving the problem of infinite dimension of the functional data by means of a novel approach of functional PLS regression for repeated measures described in Section 2. This work has been motivated by two real problems related to kinematic data. In both cases, data is affected by some noise, and then some type of penalization must be considered in the estimation of the discriminant functions. To solve this problem in the case of functional data with repeated measures, a multi-class approach for penalized functional PLS, that introduces a P-spline penalty in the definition of the inner product in the PLS algorithm, is proposed in Sections 2.2 and 2.3.

In Sections 3 and 4, the performance of the two functional approaches (penalized and non-penalized, FLDA-PenMFPLS and FLDA-MFPLS, respectively) are compared with a multivariate version of LDA based on multi-class PLS regression (LDA-MPLS). As we can see in Figure 3 and Tables 1 and 2, from a classification point of view, functional approaches perform better than the multivariate version, increasing the CCR on a test sample in a very significant way. Regarding the two functional approaches, the one based on penalized multi-class functional PLS achieves the best classification ability, with a difference in the CCR of more than 10% in both case studies.

The discriminant functions estimated from penalized and non-penalized functional LDA are shown in Figures 7 and 11. It can be seen that the penalized version (FLDA-PenMFPLS) provides smoother functions, which are more interpretable in the sense of detecting periods in which the function is positive or negative.

Finally, it can be concluded that LDA based on a multi-class approach for penalized functional PLS is the most competitive method from both classification and estimation point of view.

Footnotes

Acknowledgments

This research has been supported by the research projects MTM2017-88708-P, Spanish Ministry of Economy and Competitiveness (also supported by the FEDER program) and IJCI-2017-34038, Agencia Estatal de Investigacin, Ministerio de Ciencia, Innovacin y Universidades.

Declaration of conflicting interests

The authors declared no potential conflicts of interest with respect to the research, authorship and/or publication of this article.

Funding

References

Aguilera

Aguilera-Morillo

(2013) Comparative study of different B-spline approaches for functional data. Mathematical and Computer Modelling , 58, 1568–79.

Aguilera

Aguilera-Morillo

Preda

(2016) Penalized versions of functional PLS regression. Chemometrics and Intelligent Laboratory Systems , 154, 80–92.

Aguilera

Escabias

Preda

Saporta

(2010) Using basis expansion for estimating functional PLS regression: Applications with chemometric data. Chemometrics and Intelligent Laboratory Systems , 104, 289–305.

Aguilera-Morillo

Aguilera

(2013) P-spline estimation of functional classification methods for improving the quality in the food industry. Communications in Statistics: Simulation and Computation , 44, 2513–34.

Anguita

Ghio

Oneto

Parra

Reyes-Ortiz

(2013) A Public Domain Dataset for Human Activity Recognition Using Smartphones. In Proceedings of the 21th European Symposium on Artificial Neural Networks, Computational Intelligence and Machine Learning, ESANN 2013, Bruges, Belgium , 437–42.

Barker

Rayens

(2003) Partial least squares for discrimination. Journal of Chemometrics , 17, 166–73.

Bottomley

Daemen

Mukri

Papageorghiou

Kirk

Pexsters

De Moor

Timmerman

Bourne

(2009) Functional linear discriminant analysis: A new longitudinal approach to the assessment of embryonic growth. Human Reproduction , 24, 278–83.

Crane

Cassidy

Rothman

Gerstner

(2010) Effect of registration on cyclical kinematic data. Journal of Biomechanics , 43, 2444–47.

Cuevas

Febrero

Fraiman

(2004) An ANOVA test for functional data. Computational Statistics and Data Analysis , 47, 111–22.

10.

Daffertshofer

Lamoth

CJC

Meijer

Beek

(2004) PCA in studying coordination and variability: A tutorial. Clinical Biomechanics , 19, 415–28.

11.

Dai

Müller

Yao

(2017) Optimal Bayes classifiers for functional data and density ratios. Biometrika , 79, 545–60.

12.

Dauxois

Pousse

Romain

(1982) Asymptotic theory for the principal component analysis of a vector random function: Some applications to statistical inference. Journal of Multivariate Analysis , 12, 136–54.

13.

de Nord

Theobald

(2005) Multilevel component analysis and multilevel PLS of chemical process data. Journal of Chemometrics , 19, 301–07.

14.

Delaigle

Hall

(2012a) Componentwise classification and clustering of functional data. Biometrika , 99, 299–313.

15.

Delaigle

Hall

(2012b) Methodology and theory for partial least squares applied to functional data. Annals of Statistics , 40, 322–52.

16.

Delaigle

Hall

(2012c) Achieving near perfect classification for functional data. Journal of the Royal Statistical Society: Series B , 74, 267–86.

17.

Deluzio

Wyss

Zee

Costigan

Serbie

(1997) Principal component models of knee kinematics and kinetics: Normal vs. pathological gait patterns. Human Movement Science , 16, 201–17.

18.

Deville

(1974) Methodes statistiques et numeriques de l'analyse harmonique. Annales De L'insee , 123, 3–101.

19.

Dona

Preatoni

Claudio Cobelli

Rodano

Harrison

(2009) Application of functional principal component analysis in race walking: An emerging methodology, sports biomechanics. Sports Biomechanics , 8, 284–301.

20.

Donoghue

Harrison

Coffey

Hayes

(2008) Functional data analysis of running kinematics in chronic Achilles tendon injury. Medicine and Science in Sport and Exercise , 40, 1323–35.

21.

Durban

Harezlak

Wand

Carroll

(2005). Simple fitting of subject-specific curves for longitudinal data. Statistics in Medicine , 24, 1153–67.

22.

Eilers

PHC

Marx

(1996) Flexible smoothing with B-splines and penalties. Statistical Science , 11, 89–121.

23.

Escabias

Aguilera

Aguilera-Morillo

(2014) Functional PCA and base-line logit models. Journal of Classification , 31, 296–324.

24.

Escabias

Aguilera

Valderrama

(2004) Principal component estimation of functional logistic regression: Discussion of two different approaches. Nonparametric Statistics , 16, 365–84.

25.

Escabias

Aguilera

Valderrama

(2005) Modeling environmental data by functional principal component logistic regression. Environmetrics , 16, 95–107.

26.

Escabias

Aguilera

Valderrama

(2007) Functional PLS logit regression model. Computational Statistics and Data Analysis , 51, 4891–4902.

27.

Ferraty

Vieu

(2006. Nonparametric Functional Data Analysis. Theory and Practice . Berlin: Springer-Verlag.

28.

Fisher

(1936) The use of multiple measurements in taxonomic problems. Annals of Eugenics , 7, 179–88.

29.

Frank

Friedman

(1989) Classification: Oldtimers and newcomers. Journal of Chemometrics , 3, 463–75.

30.

Friedman

(1989) Regularized discriminant analysis. Journal of the American Statistical Association , 84, 165–75.

31.

Harrison

Ryan

Hayes

(2007) Functional data analysis of joint coordination in the development of vertical jump performance. Sports Biomechanics , 6, 196–211.

32.

Hastie

Buja

Tibshirani

(1995) Penalized discriminant analysis. The Annals of Statistics , 23, 73–102.

33.

Horvath

Kokoszka

(2012. Inference for Functional Data with Applications . Berlin: Springer-Verlag.

34.

Indahl

Sahni

Kirkhus

Naes

(1999) Multivariate strategies for classification based on NIR-spectra–with application to mayonnaise. Chemometrics and Intelligent Laboratory Systems , 49, 19–31.

35.

James

G. M.

(2002). Generalized Linear Models with Functional Predictors. Journal of the Royal Statistical Society. Series B , 64, 411–32.

36.

James

Hastie

(2001) Functional linear discriminant analysis for irregularly sampled curves. Journal of the Royal Statistical Society: Series B , 63, 533–50.

37.

Kokoszka

Reimherr

(2017. Introduction to Functional Data Analysis . London: Chapman and Hall/CRC.

38.

van de Voort

Ismail

Sedman

Cox

Simard

Buijs

(2000) Discrimination of edible oil products and quantitative determination of their iodine value by Fourier transform near-infrared spectroscopy. Journal of the American Oil Chemists’ Society , 77, 29–36.

39.

Liquet

Lê Cao

Hocini

Thiebaut

(2012) A novel approach for biomarker selection and the integration of repeated measures experiments from two assays. BMC Bioinformatics , 13, 325.

40.

Milovanovic

Popovic

(2012) Principal component analysis of gait kinematics data in acute and chronic stroke patients. Computational and Mathematical Methods in Medicine , 2012, 1–8.

41.

Mkhadri

Celeux

Nasroallah

(1997) Regularization in discriminant analysis: An overview. Computational Statistics and Data Analysis , 23, 403–23.

42.

Ocaña FA, Aguilera

Escabias

(2007) Computational considerations in functional principal component analysis. Computational Statistics , 22, 449–65.

43.

Ocaña FA, Aguilera

Valderrama

(1999) Functional principal components analysis by choice of norm. Journal of Multivariate Analysis , 71, 262–76.

44.

Orantes-González E, Heredia-Jimenez J and Soto-Hermoso

(2015) The effect of school trolley load on spatiotemporal gait parameters of children. Gait and Posture , 42, 390–93.

45.

O'Sullivan

(1986) A statistical perspective on ill-posed inverse problems. Statistical Science , 1, 505–27.

46.

Perez-Enciso M and Tenenhaus

(2003) Prediction of clinical outcome with microarray data: A partial least squares discriminant analysis (PLS-DA) approach. Human Genetics , 112, 581–92.

47.

Preda

Saporta

(2005) PLS regression on a stochastic process. Computational Statistics and Data Analysis , 48, 149–58.

48.

Preda

Saporta

Leveder

(2007) PLS classification of functional data. Computational Statistics , 22, 223–35.

49.

R Development Core Team (2008) R: A Language and Environment for Statistical Computing . Vienna: R Foundation for Statistical Computing. URL http://www.R-project.org (last accessed 29 August 2019).

50.

Ramsay

Silverman

(1997. Functional Data Analysis . Berlin: Springer-Verlag.

51.

Ramsay

Silverman

(2002) Applied Functional Data Analysis: Methods and Case Studies . Berlin: Springer-Verlag.

52.

Ramsay

Hooker

Graves

(2009. Functional Data Analysis with R and MATLAB . Berlin: Springer-Verlag.

53.

Reiss

Ogden

(2007) Functional principal component regression and functional partial least squares. Journal of the American Statistical Association , 102, 984–96.

54.

Rohart

Gautier

Singh

Lê Cao

(2017) mixOmics: An R package for ’omics feature selection and multiple data integration. PLoS Computational Biology , 13, 1–19.

55.

Ryan

Harrison

Hayes

(2006) Functional data analysis in biomechanics: A case study of knee joint vertical jump kinematics. Sports Biomechanics , 5, 121–38.

56.

Sadeghi

(2003) Local or global asymmetry in gait of people without impairments. Gait and Posture , 17, 197–204.

57.

Sadeghi

Allard

Duhaime

(1997). Functional gait asymmetry in able-bodied subjects. Human Movement Science , 16, 243–58.

58.

Sato

(1994). Application of principal- component analysis on near-infrared spectroscopic data of vegetable oils for their classification. Journal of the American Oil Chemists’ Society , 71, 293–98.

59.

Schmidt

Docherty

(2010) Comparison of the posture of school children carrying backpacks versus pulling them on trolleys. Clinical Chiropractic , 13, 253–60.

60.

Shin

(2008) An extension of Fisher's discriminant analysis for stochastic processes. Journal of Multivariate Analysis , 99, 1191–1216.

61.

Westerhuis

van Velzen

Hoefsloot

Smilde

(2010) Multivariate paired data analysis: Multilevel PLSDA versus OPLSDA. Metabolomics , 6, 119–28.