Fractal Analyses Reveal Independent Complexity and Predictability of Gait

TL;DR

This study applies advanced non-linear analysis techniques to gait time series, revealing that complexity and predictability indexes can distinguish different walking conditions and may improve clinical assessments.

Contribution

It introduces the combined use of Hurst exponent and Minkowski fractal dimension to analyze gait, providing more sensitive indexes for detecting locomotion abnormalities.

Findings

Walking forward shows highest complexity and adaptability.

Perturbations decrease gait complexity.

Gait predictability varies with walking direction and vestibular stimulation.

Abstract

Locomotion is a natural task that has been assessed since decades and used as a proxy to highlight impairments of various origins. Most studies adopted classical linear analyses of spatio-temporal gait parameters. Here, we use more advanced, yet not less practical, non-linear techniques to analyse gait time series of healthy subjects. We aimed at finding more sensitive indexes related to spatio-temporal gait parameters than those previously used, with the hope to better identify abnormal locomotion. We analysed large-scale stride interval time series and mean step width in 34 participants while altering walking direction (forward vs. backward walking) and with or without galvanic vestibular stimulation. The Hurst exponent $\alpha$ and the Minkowski fractal dimension $D$ were computed and interpreted as indexes expressing predictability and complexity of stride interval time series, respectively. We show that $\alpha$ and $D$ accurately capture stride interval changes in function of the experimental condition. Walking forward exhibited maximal complexity ($D$) and hence, adaptability. In contrast, any perturbation (walking backward and/or stimulation of the vestibular system) decreased it. Furthermore, walking backward increased predictability ($\alpha$) through a more stereotyped pattern of the stride interval and galvanic vestibular stimulation reduced predictability. The present study demonstrates the complementary power of the Hurst exponent and the fractal dimension to improve walking classification. These holistic indexes can easily be interpreted in the framework of optimal movement complexity. Our developments may have immediate applications in rehabilitation, diagnosis, and classification procedures.

Tables

Table 1. Table 1: Mean ± plus-or-minus \pm standard deviations values for the parameters computed in the four experimental conditions. FW (BW) stands for forward (backward) walking and S0 (S + {\rm+} ) for no GVS (GVS on).

	FWS0	FWS+	BWS0	BWS+
$T$ (s)	$1.12\pm 0.07$	$1.11\pm 0.17$	$1.51\pm 0.16$	$1.55\pm 0.17$
$C{V}_T$	$0.020\pm 0.007$	$0.017\pm 0.006$	$0.069\pm 0.018$	$0.063\pm 0.020$
$w$ (s)	$0.104\pm 0.029$	$0.098\pm 0.027$	$0.245\pm 0.065$	$0.244\pm 0.065$
${\theta }_0$ (o)	$21.3\pm 1.5$	$21.2\pm 1.4$	$14.2\pm 1.9$	$14.5\pm 1.8$
$\alpha$	$0.725\pm 0.062$	$0.726\pm 0.048$	$0.757\pm 0.060$	$0.716\pm 0.064$
$\gamma$	$0.507\pm 0.216$	$0.500\pm 0.167$	$0.655\pm 0.194$	$0.617\pm 0.202$
$D$	$1.72\pm 0.12$	$1.66\pm 0.11$	$1.62\pm 0.10$	$1.65\pm 0.11$

Table 2. Table 2: Significance level ( p − limit-from 𝑝 p- values) of the post hoc (Holm-Sidak) pairwise multiple comparisons between the different experimental conditions for the parameters analyzed. Statistically significant differences are written in bold font.

	FWS0 vs FWS+	BWS0 vs BWS+	FWS0 vs BWS0	FWS+ vs BWS+
$T$ (s)	.529
$C{V}_T$	.308	$.012$
$w$		.644		$.005$
${\theta }_0$	.176
$\alpha$	.857	.420	$.024$	.157
$D$	$.042$	.143		.592

Table 3. Table 3: Median scores (in %) and interquartile range (IQR) results of the MSAQ. The first (Q1) and third quartiles (Q3) are shown between square brackets.

	FWS+	BWS+
Central	27 $\left[13-51\right]$	22 $\left[11-40\right]$
Gastrointestinal	14 $\left[11-19\right]$	14 $\left[11-17\right]$
Peripheral	11 $\left[10-15\right]$	11 $\left[10-19\right]$
Sopite-related	19 $\left[14-25\right]$	17 $\left[11-22\right]$

Equations

$$\alpha =\frac{1+\gamma }{2}$$

Full text

Fractal Analyses Reveal Independent Complexity and Predictability of Gait.

Dierick Frédéric1,2, Nivard Anne-Laure1, White Olivier3,4, Buisseret Fabien1,5,*

1 Forme et Fonctionnement Humain Research Unit, Department of Physical Therapy, Haute Ecole Louvain en Hainaut (HELHa), Rue Trieu Kaisin, 136, 6061 Montignies-sur-Sambre, Belgium

2 Faculté des Sciences de la Motricité, Université catholique de Louvain, Louvain-la-Neuve, Belgium

3 Université de Bourgogne INSERM-U1093 Cognition, Action, and Sensorimotor Plasticity, Campus Universitaire, BP 27877, 21078 Dijon, France.

4 Acquired Brain Injury Rehabilitation Alliance, School of Health Sciences, University of East Anglia, Norwich, Norfolk, UK.

5 Service de Physique Nucléaire et Subnucléaire, Université de Mons, UMONS Research Institute for Complex Systems, Place du Parc 20, 7000 Mons, Belgium

• buisseretf@helha.be

Abstract

Locomotion is a natural task that has been assessed since decades and used as a proxy to highlight impairments of various origins. Most studies adopted classical linear analyses of spatio-temporal gait parameters. Here, we use more advanced, yet not less practical, non-linear techniques to analyse gait time series of healthy subjects. We aimed at finding more sensitive indexes related to spatio-temporal gait parameters than those previously used, with the hope to better identify abnormal locomotion. We analysed large-scale stride interval time series and mean step width in 34 participants while altering walking direction (forward vs. backward walking) and with or without galvanic vestibular stimulation. The Hurst exponent $\alpha$ and the Minkowski fractal dimension $D$ were computed and interpreted as indexes expressing predictability and complexity of stride interval time series, respectively. We show that $\alpha$ and $D$ accurately capture stride interval changes in function of the experimental condition. Walking forward exhibited maximal complexity ($D$) and hence, adaptability. In contrast, any perturbation (walking backward and/or stimulation of the vestibular system) decreased it. Furthermore, walking backward increased predictability ($\alpha$) through a more stereotyped pattern of the stride interval and galvanic vestibular stimulation reduced predictability. The present study demonstrates the complementary power of the Hurst exponent and the fractal dimension to improve walking classification. These holistic indexes can easily be interpreted in the framework of optimal movement complexity. Our developments may have immediate applications in rehabilitation, diagnosis, and classification procedures.

Introduction

The stride interval of normal human walking is the time period between consecutive heel strikes of the same foot [1]. For more than two decades, a line of research focused on the understanding of the nature of the subtle variations observed in stride intervals and the origin of typical long-range structures in these variations. Today, these investigations are of paramount importance since they could provide a better understanding of the physiological mechanisms involved in normal human walking and in alterations observed in clinical practice. The nature of these stride interval variations could arise either from noisy neural processes that result in errors in the motor output or from alterations in the motor command that account for balance instabilities [2].

Normal gait is characterized by the presence of autocorrelations in the stride interval when considering walking on a sufficiently long time scale [3, 1]. The origin of these autocorrelations may be attributed to neural central pattern generators (CPGs) [3, 1] or a super CPG coupled to a forced Van der Pol oscillator [4], and/or to the biomechanics of walking [5, 6]. For many years, gait analysis has been studied with classical methods adopting biomechanical models in which variability was not of interest. More recent techniques derived from chaos theory are well adapted to analyse time series that exhibit long-range autocorrelation. Importantly, they treat variability as a meaningful interpretable signal. Since the pioneering works of Hausdorff et al. [3, 1], long-range autocorrelations in time series are estimated by the Hurst or fractal exponent ($\alpha$). A fractal, introduced in 1975 by the French mathematician Benoît Mandelbrot (1924–2010) [7], is defined as a geometrical structure that has a regular or an uneven shape repeated over all scales of measurement. It is characterized by a fractal dimension ($D$) greater than the spatial dimension of the structure [8]. A famous example of such object is a snow flake. Objects that are statistically self-similar – parts of it show the same statistical properties at many scales – exhibit strong autocorrelation. The Hurst exponent $\alpha$ is a statistical measure of long-term memory of time series (see e.g. [9] for a review) and is usually associated to fractal-like behaviour. In particular, the peculiar behavior of the stride interval may be referred to as “fractal behavior” [3].

The theoretical model of optimal movement complexity [10] is based on the complementary concepts of predictability and complexity. Nature let us find optimal behavior in terms of skills and variability through evolution. The optimal state of a biological system is characterized by chaotic temporal variations in the steady state output that correspond to maximal predictability. Too few practice results in high disorder (randomness, no predictability) and excessive practice leads to high order (periodic signal, maximal predictability). Adaptation of a system to external stimuli is maximal only at an intermediate state of predictability. Furthermore, a signal from a dynamical system also holds some inherent complexity. A decrease of complexity of a physiological system results from either a reduction in the number of structural components or an alteration in the coupling function between these components. For instance, a joint can become rigid with senescence, hence decreasing the degree of freedom of the system and consequently, its complexity. A holistic approach to study these mechanisms requires to associate specific measurements to these two concepts. The Hurst exponent $\alpha$ captures part of the story and is well suited to reflect predictability. The Minkowski fractal dimension $D$ provides good measurability of the “apparent rugosity” of fractals [11] and reflects complexity. Here, we use these parameters to complement the usual quantification of autocorrelation $\alpha$ in unusual and perturbed gait conditions in an attempt to probe adaptability in the framework of the model of optimal movement complexity [10].

As of today, the vast majority of studies explored autocorrelation in the stride interval during natural forward walking. In one notable exception however, Bollens et al. [12] also tested backward walking in a small sample of young healthy subjects. The authors did not find significant differences in long-range autocorrelation between both walking directions. However, backward walking measures revealed to be more sensitive than forward walking measures to classify elderly fallers compare to non-fallers [13]. The study of backward walking under the perspective of fractal analyses is therefore promising to provide more reliable predictive index of fallers, as previously proposed for forward walking [14]. Backward walking is also frequently used in sports and in rehabilitation settings, and a better understanding of the variability of stride interval in this condition is needed since it is believed that backward walking is at least partly controlled by specialized neural circuits [15].

The vestibular system provides an essential sensory contribution to the maintenance of balance during human walking [16]. Individuals with vestibular disorders show a decreased walking stability accompanied by an increased risk to fall [17]. Therefore, perturbing the vestibular system of healthy subjects with galvanic vestibular stimulation (GVS) is a well targeted mean to probe gait: it is standardized, well tolerated by subjects, generated by currently affordable electrostimulators, and easy to implement when a large number of stride intervals are recorded with an instrumented treadmill. The use of GVS is also an increasingly common clinical intervention on locomotion [18, 19, 20].

Previously, autocorrelations in stride interval time series have been identified not only in healthy young adults [3] but also in children [21] and elderly [22], and even – although significantly modified – in several neurodegenerative conditions. In particular, the cases of Huntington’s disease [22], amyotrophic lateral sclerosis [23], and Parkinson’s disease have been studied [24, 25], with a hope of connecting the observed modifications of fractal behavior to some relevant evaluation of the risk of falling [14]. Here, we hypothesize that the combined effects of walking direction (WD) and the application of GVS on long-range autocorrelations in the stride interval could enhance the sensitivity of fractal analysis to identify impaired gait. We measured $\alpha$ and $D$ during forward and backward walking, with and without the application of binaural and monaural GVS. We speculate that these two indexes should be able to capture differences between experimental conditions and therefore provide better indexes to classify patients.

Material and methods

Participants

Thirty-four undergraduate and graduate healthy students (18 males, 16 females) in physiotherapy took part to this study and were recruited at Haute Ecole Louvain en Hainaut (Charleroi, Belgium). Mean age was 23 years (standard deviation, $SD$=2), height was 173 cm ($SD$=9), mass was 69 kg ($SD$=10), body mass index was 23 kg m*-2* ($SD$=3), and lower limb length ($L$), measured in standing position as the distance between the floor and great trochanter, was 88 cm ($SD$=5).

Subjects were not medicated and did not exhibit any neuromusculoskeletal, orthopaedic, respiratory, or cardiovascular disorders that could influence their gait. Exclusion criteria included vestibular disorders in addition to specific GVS exclusion criteria: presence of a heart pacemaker, pregnancy, metallic brain implants, epilepsy, and skin damage behind the ears or forehead. Eligible participants were required to be able to respond to verbal questions, comprehend questionnaires, and understand instructions during the procedures of the study. Prior to participating, subjects read and signed an informed consent form. The study was approved by the ethics committee of Grand Hôpital de Charleroi and conducted in accordance with the declaration of Helsinki.

Experimental procedure

Subjects walked on an instrumented treadmill (70 cm wide, 185 cm long) with an integrated force plate and an overhead safety frame (N-Mill, Motekforce Link, The Netherlands). They wore comfortable running shoes, a safety harness, and were asked to keep their eyes fixed straight ahead. Four walking conditions were studied during two measurement sessions on two different days: forward walking without GVS (FWS0) and with GVS (FWS+) and backward walking without GVS (BWS0) and with GVS (BWS+). Session 1 included FWS0 and FWS+ conditions and session 2 included BWS0 and BWS+ conditions. During each condition, subjects walked on the treadmill for 15 minutes. Before each session, subjects were given five minutes to familiarize themselves with the treadmill and the conditions.

Subjects walked at their comfortable speed that was determined during the familiarization procedure by the same experimenter (NAL) by tuning the speed of the belt while the subject was walking without GVS. The same speed was then imposed when GVS was applied. Vertical ground reaction force ($F_{v}$) and centre of pressure ($CoP$) of each foot was recorded at a sampling rate of 500 Hz using the manufacturer’s software (CueFors 2, Motekforce Link, The Netherlands). Time series stride interval were computed from heel strikes (during forward walking) or toe strikes (during backward walking) of the right foot identified on $F_{v}$-time histories and time series step width from maximal medio-lateral displacement of $CoP$ of two consecutive steps. At the completion of both sessions, four time series containing the values of the stride intervals in the different conditions were obtained for each subject. Typical plots are displayed in Fig. 1.

Bipolar GVS was applied with a regulated, direct-current device (Compex 3 Professional, Compex Medical SA, Switzerland) with a maximum output current of 20 mA by steps of 0.125 mA. The carbon electrodes (20 cm2) were covered with a saline-soaked sponge held in place over the mastoids or forehead with a strap. The 34 subjects were randomly exposed to one of the three different transcranial stimulation conditions: binaural ($n=12$), unilateral left ($n=11$), and unilateral right ($n=11$). For the binaural stimulation, electrodes were randomly placed over the mastoids with a cathode-left anode-right or cathode-right anode-left montage. For the monoaural stimulations, the cathode was randomly placed over the right mastoid and the anode on the right part of the forehead (right stimulation) or the cathode was placed over the left mastoid and the anode on the left part of the forehead (left stimulation). Fig. 2 shows a schematic representation of the location of the electrodes over the head for the 3 stimulation conditions. The intensity was set at the highest sensory tolerance threshold, that was determined by increasing the current intensity slowly by 0.125-mA steps. That intensity was maintained constant throughout the walking period. Mean current density for subjects was 0.07 mA cm*-2* (range: 0.04–0.08 mA cm*-2*). The intensity and duration of the GVS adhered to the safety criteria for transcranial direct current stimulation [27]. After FWS+ and BWS+ conditions, each subject completed a home-made French translation of the Motion Sickness Assessment Questionnaire (MSAQ) [26], that consists of 16 questions, allowing to differentiate motion sickness symptoms along the gastrointestinal, central, peripheral, and sopite-related dimensions.

Data analysis

The treadmill software directly computed the mean stride interval, $T$, and the mean step width, $w$, for each subject in each condition. The stride amplitude ($\theta_{0}$), i.e. the angle between the leg and the vertical at heel strike, has then been computed from the relation $v\,\frac{T}{4}=L\tan\theta_{0}$, displayed in Fig. 3, where $v$ is the walking speed and $L$ is the lower limb length of the subject.

The temporal analysis of our experimental data has to go a step beyond mean values to study the information contained in stride interval variability. Let $\bm{T}=\{T_{i}:i=1,\dots,n_{0}\}$ be a time series, where $T_{i}$ is the stride interval of cycle $i$ and where $n_{0}$ is the number of cycles recorded during 15 minutes. The first indicator of variability is the coefficient of variation $CV_{T}=SD(\bm{T})/T$. Because $CV_{T}$ provides no information on the dynamics of the stride interval fluctuations, several indexes characterizing the dynamical structure of the time series have been computed such as the Hurst exponent, $\alpha$, the spectral exponent, $\gamma$, and the Minkowski fractal dimension, $D$.

We assess the presence of long-range autocorrelations with the Hurst exponent $\alpha$ and the spectral exponent $\gamma$. These parameters quantify the “predictability” of the time series. As pointed out in e.g. [28, 29], using a single parameter may not be generally enough to assess the presence or not of such autocorrelations. The Hurst exponent, computed by using the Detrended Fluctuation Analysis (DFA) with a linear detrending [30], provides a diagnostic on the long-range trend of the time series. DFA consists in several steps. First, one has to compute the shifted time series $\bm{T}_{\tau}=\{T_{i+\tau}:i=1,\dots,n_{0}-\tau\}$ and the cumulated time series ${\bm{Z}}=\{Z_{i}=\sum^{i}_{j=1}(T_{j}-T):i=1,\dots,n_{0}\}$. One has then to divide the cumulated time series ${\bm{Z}}$ into windows of length $l$, leading to the samples ${\bm{Z}}^{(m)}(l)$, $m$ labelling the window. For each window, a local least squares linear fit is calculated, leading to the fitted values ${\bm{\bar{}}Z}^{(m)}(l)$. Second, one computes the fluctuation function $F(l)=\sqrt{\frac{1}{n_{0}}\sum^{n_{0}}_{i=1}\left[Z_{i}^{(m)}(l)-\bar{Z}_{i}^{(m)}(l)\right]^{2}}$. The Hurst exponent is then defined as the scaling exponent of $F(l)$, i.e. $F(l)\propto l^{\alpha}$. Stationary times series originating from long-range (anti)correlated processes correspond to $0.5<\alpha\lesssim 1$ ($0\leq\alpha<0.5$), respectively. When $\alpha=0.5$, the process is random. Values larger than 1 correspond to unbounded, unstable, processes [9]. In our point of view, a strongly autocorrelated signal can be denoted “predictable”: Its value at a given step is strongly dependent of the system’s previous state.

The spectral exponent $\gamma$ can be extracted from the low-frequency behaviour of the power spectral density $P(f)$ of $\bm{T}$, $f$ being the frequency: $P(f)\propto f^{-\gamma}$. Actually $P(f)$ is the Fourier transform of the autocorrelation function $C(\tau)=E(\bm{T}\bm{T}_{\tau})$, where $E$ denotes the average value. The parameter $\gamma$ is expected to take values between 0 and 1 for long-range autocorrelated processes. For large enough time series, the asymptotic relation

[TABLE]

should relate $\alpha$ and $\gamma$ [9].

We finally compute the Minkowski fractal dimension $D$ of the time series, defined through the box-counting method stating that, if $N(\epsilon)$ is the number of square boxes of size $\epsilon$ needed to fully cover the time series once plotted, then one has $N(\epsilon)\propto\epsilon^{-D}$ for small $\epsilon$. For time series such as the ones we deal with, $D$ will typically lie between 1 (differentiable curve) and 2 (surface with differentiable boundary). Even if $\bm{T}$ does not define a fractal in a rigorous mathematical way, $D$ can be thought of as a relevant estimator of the “apparent roughness” of the corresponding curve [31]. In our opinion, assessing the roughness of the stride interval time series – i.e. the variability of fluctuations from one stride to another – may be associated to the “complexity” of the process [32].

It is worth mentioning at this stage that $D$ and $\alpha$ may be seen as independent variables characterizing a time series. The link $\alpha=2-D$ is actually valid only for some widely studied random walks, but processes with arbitrary values of $\alpha$ and $D$ can be built and could be more representative of realistic time series [33].

For completeness we show in Fig. 4 typical plots of the necessary computations performed in view of computing the parameters $\alpha$, $\gamma$ and $D$. All these quantities are actually the slopes of the different linear regressions performed.

Statistical analysis

All data were checked for normality (Shapiro-Wilk) and equal variance tests. A two-way (WD $\times$ GVS) repeated measures ANOVA with post hoc Holm-Sidak method for pairwise multiple comparisons has been performed and used to examine the effects of walking direction ($FW$ or $BW$), GVS ($S+$ or $S0$), and their interaction on the computed parameters. The significance level has been set at $p=.05$ for all analyses and post hoc statistical power has been calculated. The correlation between parameters are provided as Spearman’s rank correlation coefficient $\rho$, Pearson’s correlation coefficient $r$, and prinicipal component analysis (PCA).

Finally, the median ($Mdn$) scores and interquartile ranges ($IQR$) related to the four symptom dimensions assessed by the MSAQ questionnaire have been computed from the 16 items scores according to [26].

All statistical procedures were performed with SigmaPlot software version 11.0 (Systat Software, San Jose, CA). Indexes have been computed using R free software environment (v. 3.2.2) [34].

Results

Subjects adopted mean comfortable forward speed of 4.4 km h*-1* ($SD$=0.4) and backward speed of 2.2 km h*-1* ($SD$=0.5). The rather small SDs indicate that our sample was quite homogeneous with respect to that parameter. Results of all parameters analysed according to the four experimental conditions are reported in Table 1. Some parameters particularly relevant for the Discussion are also shown as box plots in Fig. 5. In any of the four experimental conditions, stride interval time series are such that $0<\gamma<1$ and $0.5<\alpha<1$. These values confirm that time series characterise long-range autocorrelated processes with memory [28]. Note that our value of $\alpha$ in the FWS0 condition ($95\%$ CI of $[0.78-0.90]$) are fully compatible with the interval reported in the meta-analysis [35] in healthy subjects ($95\%$ CI of $[0.85-0.97]$). Moreover, Pearson’s correlation coefficient between $\alpha$ and $\gamma$ (computed over the complete data set) is equal to $r=0.287$, $p<.01$. Parameters $\alpha$ and $\gamma$ can then be seen as linearly correlated as expected from Eq. (1). We decided to keep $\alpha$ since it is the most widely used autocorrelation-related parameter in the literature [35]. The parameter $\gamma$ will not be considered in the following anymore since it only confirms the presence of long-range autocorrelations. In contrast, Pearson’s and Spearman’s correlation coefficients between $\alpha$ and $D$ are equal to $r=-0.081$, $p=.350$ and $\rho=-0.031$, $p=.719$ respectively. Moreover, a PCA of our complete data set, including the parameters listed in Table 1, shows that nearly 66$\%$ of the total variance is carried by the first two dimensions, the angle between $\alpha$ and $D$ being equal to 131o while that between $alpha$ and $\gamma$ being equal to 14o. $\alpha$ and $D$ can then be considered as independent parameters and are both retained in our analysis.

We first report significant differences found in the two-way RM ANOVA related to WD and GVS factors. WD induces significant differences in all variables (see also Fig. 5): $T$ ($F$=216, $p<.001$, partial effect size $\eta_{p}^{2}$= .991), $CV_{T}$ ($F$=324, $p<.001$, $\eta_{p}^{2}$= .999), $w$ ($F$=13.6, $p<.001$, $\eta_{p}^{2}$= .851), $\theta_{0}$ ($F$=429, $p<.001$ , $\eta_{p}^{2}$= .997), $\alpha$ ($F$=5.594, $p=.024$, $\eta_{p}^{2}$= .222) and $D$ ($F$=12.1, $p=.001$, $\eta_{p}^{2}$= .276). GVS has a less important influence but has nevertheless a significant impact on $T$ ($F$=5.58, $p=.024$, $\eta_{p}^{2}$= .135), $CV_{T}$ ($F$=6.24, $p=.018$, $\eta_{p}^{2}$= .872) and $w$ ($F$=11.0, $p=.002$, $\eta_{p}^{2}$= .190). The two-way RM ANOVA also reported significant interaction between WD and GVS for $T$ ($F$=1.2, $p=.003$, $\eta_{p}^{2}$= .236), $w$ ($F$=4.75, $p=.036$, $\eta_{p}^{2}$= .126), $\theta_{0}$ ($F$=11.4, $p=.002$, $\eta_{p}^{2}$= .257) and $D$ ($F$=7.69, $p=.009$, $\eta_{p}^{2}$= .189).

Post hoc pairwise multiple comparisons results are listed in Table 2. Only $w$ and $D$ are significantly modified when the GVS was active during forward walking. The GVS had a stronger impact in backward walking, with a significant modification of $T$, $CV_{T}$, $\theta_{0}$, $\alpha$. The factor having the major impact on the variables is the WD. In FWS0 vs BWS0 comparison, all parameters were significantly different. In FWS+ vs BWS+ comparison, all parameters were significantly different, except $\alpha$ and $D$.

MSAQ results are summarized in Table 3 according to the four dimensions listed in [26]. The highest median scores were reached for the items related to the central nervous system dimension. It is worth noticing that beyond these four scores, the third item of the questionnaire [26] – “I felt annoyed/irritated” scored highest in the FWS+ condition ($Mdn$=4, IQR=$[2-6]$) and in the BWS+ condition ($Mdn$=3, IQR=$[1-5]$).

Finally, one can ask the question as to how basic kinematic parameters relate to high level indexes such as $\alpha$ and $D$? Spearman’s coefficient between $\alpha$ and $\theta_{0}$, calculated from each subject’s data, is equal to $\rho=-0.207$, $p=.017$. It is therefore relevant to consider that both parameters are correlated. Figure 6 depicts the evolution of $\alpha$ vs $\theta_{0}$ in the FWS0 and BWS0 conditions and compares it to that of data presented in Ahn and Hogan’s model (full circles). Figure 7 reports the evolution of the fractal dimension versus the Hurst exponent in the four experimental conditions. It shows that changing walking direction with respect to the control condition $FW\textsubscript{S0}$ leads to an increase of $\alpha$ and a decrease of $D$ (independently of GVS). However, turning on GVS during forward walking only leads to a decrease of $D$. There is some decoupling between $\alpha$ and $D$.

Discussion

Locomotion has been used as a tool to identify and characterise diverse impairments. Here, we set out to use techniques beyond classical linear analyses of spatio-temporal gait parameters in order to define more sensitive indexes. We asked participants to walk under perturbed conditions induced either by reversing the direction of walking or perturbing the vestibular system, and measured proxies of “complexity” and “predictability” through the Hurst exponent and the fractal dimension.

Walking direction and galvanic vestibular stimulation

Walking direction has a major impact on our results. Stride interval increased in backward walking compared to forward walking which is in line with previous work e.g. [36]. This fact is coherent with the lower walking speed spontaneously chosen by the subjects when walking backward. The kinematic parameters $CV_{T}$ and $w$ are larger in the backward walking conditions than in forward walking. Fluctuations in the stride interval during backward walking are indeed larger than in forward walking. The increase of $CV_{T}$ has been reported previously [12, 13] in young and elderly adults while, to our knowledge, the corresponding increase in the step width has not been reported elsewhere.

While the existence of autocorrelations in backward and forward walking has been acknowledged in [12], that study did not find a significant difference in the Hurst exponent, presumably because of the smaller sample size (12 subjects). The Hurst exponent appears to be larger in backward compared to forward walking. While neurological diseases generally decrease $\alpha$ [35] (more random motion), an increase in $\alpha$ has been reported in children up to typically 14 years old [21]. Both backward walking in adults and forward walking in children can be related to learning processes, with stereotyped, more predictable, motion.

In their recent paper, Ahn and Hogan have shown that long-range autocorrelations may emerge from the dynamics of a particular pendular model of walking described in detail in [37, 6]. It is, to our knowledge, the only model linking kinematic variables and autocorrelations indexes. Here, we compare experimental data with that model for the first time. One of the key plots of [6] shows the variation of the Hurst exponent computed on a range of 500 strides for realistic values of the different parameters versus stride amplitude $\theta_{0}$, which is fixed in their model. A comparison of this plot to our results is shown in Fig. 6. Results from the FWS0 condition are in agreement with the model. The BWS0 condition can be compared too since the dynamical equations presented in [37] are invariant with respect to time reversal. Our results match quite well with the trend of the model. Moreover, GVS has no significant influence on $\alpha$. This is in favour of a mostly mechanical origin of long-range stride interval autocorrelation.

Galvanic vestibular stimulation is a known procedure to electrically stimulate vestibular afferents [39, 40, 38]. Here, we used continuous GVS with an average current of 1.4 mA, which was well tolerated by all subjects (low MSAQ median scores).

Previously, vestibular inputs were thought to be primarily required for stabilizing the head to ensure stable gaze control during gait and for spatial orientation in navigational tasks [42, 41, 38]. More recently, it has also been suggested that vestibular inputs play a role in maintaining dynamic walking stability since they generate phase-dependent influences on lower body control during walking by fine tuning the timing and magnitude of foot displacement [43, 16, 44]. In agreement with those recent findings, our results show that GVS significantly modifies $T$, $CV_{T}$ and $w$. The magnitude of $CV_{T}$ is regularly associated to the risk of fall [45]. We observe that turning on GVS significantly decreases $CV_{T}$ in backward walking; hence training in BWS+ may be relevant to decrease the risk of fall. Previous research has demonstrated that GVS mostly affects stability in the medio-lateral direction [39, 40, 38, 46]. This is in accordance with our results that show an increased $w$ during forward walking.

Interestingly, we found that GVS induced larger $D$ in FWS0 compared to FWS+ condition. This indicates a less complex stride interval time series in the non-stimulated condition of forward walking. This result is in line with [43] showing that planning of the foot placement at heel contact is modulated by vestibular information. Here, we provide another evidence of the influence of GVS on walking variability. It is known that the vestibular system is essential to the maintenance of balance throughout the stepping cycle, with phase-velocity/cadence-dependent modulation on the activity of hip, knee and ankle muscles [16]. Vestibular-muscle coupling is specific for each muscle, probably organised according to each muscle’s functional role in whole-body stabilization during walking. Our analyses suggest that the less complex nature of the stride interval time series reflects the disruption of dynamic balance evoked by GVS.

The optimal complexity model

Stergiou and Decker [10] proposed that time series originating from human motion could be classified by using two indices. The first catches signal complexity and the other measures its predictability. In this context, a healthy motion should be chaotic, characterized by a maximal complexity reflecting the adaptability of the subject to exterior perturbation, and an intermediate predictability. Pathological motion should be characterized by a lower complexity (fewer adaptability) and a predictability that could be either lower (random, “drunken-sailor-like”, motion) or higher (robotic motion) than the healthy motion.

Following on that line, we interpret $D$ as a measure of the complexity of gait time series. Indeed, a large fractal dimension is associated to an apparently rough time series, with abrupt relative changes of values stepwise. A complex time series may be the signature of an adaptable behavior. The more a subject is able to change his/her stride interval from one cycle to another, the more s/he should be able to modify his/her pattern. Therefore, $D$ could be a good indicator of complexity during walking. Moreover, we think that the Hurst exponent – that was independent of $D$ – could be a relevant predictability index. Indeed, $\alpha$ can discriminate between a random motion ($\alpha=0.5$) and a far more predictable, strongly autocorrelated, time series ($0.5<\alpha\lesssim 1$). So $\alpha$ may provide an answer to the question as to how much a stride interval depends on history? This is exactly what predictability stands for. Healthy subject should be characterised by a maximal value of $D$ (high complexity, good adaptation skills) and an optimal value of $\alpha$ (good but not too high predictability). Any significant deviations from these values could indicate pathological motion linked to any or both dimensions (predictability or complexity).

Our results are displayed in a $(\alpha,D)$-plane in Fig. 7. It clearly appears that the FWS0 condition – the healthy motion – has the higher complexity and an intermediate predictability as argued in [10]. The other conditions, non-standard but not pathological either, have lower complexities. Walking backward without GVS leads to a larger value of $\alpha$, that is a more stereotyped, more predictable walking. GVS slightly decreased $\alpha$ in backward walking. In that condition, walking gets closer to a random process, presumably because of the perturbation of the vestibular system. It is worth noticing that at least one of the two parameters $(\alpha,D)$ is significantly modified when going from one condition to another. A two-dimensional representation is necessary to classify all the experimental conditions we study. Hence it can be conjectured that the only study of $\alpha$ in pathological cases may be too restrictive to discriminate between the pathologies and that more non-linear indexes are worth to be added.

Previous studies only computed the Hurst exponent and implicitly considered that the fractal dimension and the exponent were related. Here, we computed $D$ beside $\alpha$ and show for the first time that these two parameters are actually decoupled in some conditions. As already pointed out, children walking forward also have a larger $\alpha$ than healthy young adults walking forward. We have calculated from a freely available dataset [21] the average $D$ for the 50 children having participated to the study. We found values equal to $1.45\pm 0.21$, a significantly lower value than in our FWS0 condition ($t=-6.77$, $p<.001$), as expected. In Parkinson’s disease, $\alpha$ are smaller than in young healthy adults. It has been shown [25] that $\alpha$ decreases with disease’s severity. We have computed $D=1.69\pm 0.10$ from the data of [25] (20 patients with Parkinson’s disease walking for 10 min, with $\alpha=0.70\pm 0.09$). Although lower than our maximal, FWS0, the difference between both values is not significant ($t=-0.941$, $p=.351$). Similarly, neurodegenerative pathologies have actually been shown to generally decrease $\alpha$ with respect to its optimal value [35]. As can be deduced from the above discussion, the parameters $D$ and $\alpha$ are good candidates to disentangle and characterise the main long-term features of walking. The Hurst exponent $\alpha$ is a widely used indicator of long-term autocorrelations, and adding $D$ opens new classification perspectives. Our findings may have immediate applications in rehabilitation, diagnosis, and classification procedures.

We also think that this field could benefit in a near future from new techniques such as a representation of the stride interval time series in terms of complex networks (visibility graphs) [47, 48]. This technique has already proven to be efficient to distinguish healthy from epileptic EEG signals [49], hence it can reasonably be assumed that visibility graphs could provide relevant information on the structure of stride interval time series. New classification schemes allowed nowadays by machine learning could also shed new light on walking dynamics. Algorithms like random forests could help to find better indices to disentangle the different experimental conditions [50]. There is hope that such new techniques could better classify the stride interval time series, but with less common indices, either less intuitive or less easily compared to the literature in walking analysis. Such a research program is beyond the scope of the present study and we leave it for future investigations.

Conclusion

Our findings show that stride interval dynamics behave as a chaotic system exhibiting long-range autocorrelations independently of walking direction. The Hurst exponent is increased when walking backward, suggesting that the more predictable fluctuations of the stride interval reflect more stereotyped motion adopted by subjects in response to this nonstandard condition. The magnitude of these fluctuations are however larger in backward walking, due to the weaker stability of the subjects. The Minkowski fractal dimension complements the characterisation of stride interval variability by considering complexity, or, more intuitively the adaptive capacities of the subject in motion. Any nonstandard condition reduced complexity. The present study thus opens new avenues as to how more accurately classify healthy or pathological walking according to the complexity and predictability of stride interval time series.

Acknowledgments

The authors thank C Detrembleur and T Warlop for stimulating discussions at early stages of this work and for having provided us the data of [25], and M Scohier for useful comments.