TE D

Neuroinformatics Journal Copyright ©Humana Press Inc. All rights of any nature whatsoever are reserved. ISSN 1539-2791/03/????????/$25.00

Original Article

EC

Event Identification in Movement Recordings by Means of Qualitative Patterns

RR

Eric Fimbel,1 Anne Sophie Dubarry,2 Maxime Philibert,3 and Anne Beuter4 Département de génie électrique, École de technologie supérieure, Montréal, and Centre de recherche, Institut universitaire de gériatrie de Montréal; 2Département de génie électrique, École de technologie supérieure, Montréal; 3Département de génie électrique, École de technologie supérieure, Montréal, and Centre de neuroscience de la cognition, Université du Québec à Montréal; 4Institut de biologie, laboratoire de physiologie, Université de Montpellier 1, and Centre de neuroscience de la cognition, Université du Québec à Montréal.

CO

1

Introduction

Abstract

DR AF

T-

UN

We present a pattern-matching technique for detecting events in movement recordings. The events are defined as sequences of qualitative changes in the speed and/or the higher order derivatives (e.g., in a speed peak, the acceleration changes from positive to negative). The technique uses qualitative patterns that are sequences of qualitative states (e.g., negative, infinitesimal, positive...) of the speed and the higher order derivatives. A fast pattern-matching algorithm is presented. Its sensitivity can be tuned by means of a filtering parameter, and a multiscale analysis method is proposed for detecting events of different amplitudes and durations. An application to the assessment of the irregularity of rapid movement in Parkinson’s disease is presented. Index Entries: Movement analysis; event detection; qualitative patterns; pattern matching.

Movement recordings by means of threedimensional position trackers provide sampled data representing the position of the extremity of a limb, its segments, and/or its joints. Twodimensional systems, such as digitizing pads, provide similar information, excepting that the position of the hand is assimilated to that of the pointing device. Movement recordings are mostly used in biomechanics and motor control, but also in ergonomics (Walker, 1993) and for analyzing motor disorders like Parkinson’s disease (PD), Huntington’s disease (HD) or cerebellar dysfunctions. An important issue is the quantification of the smoothness (or the irregularity) of a movement. Several criteria are used to determine whether a movement is optimally smooth. Among these criteria, we have the minimality

*Address to which all correspondence and reprint requests should be sent. E-mail: [email protected] 1

RR

EC

TE D

2__________________________________________________________________________________Fimbel et al.

UN

CO

Fig. 1. Motor events during a fast planar arm movement (normal subject). (A) position in the X-Y plane; (B) speed, (C) acceleration; (D) jerk; (E) snap. (F) qualitative changes of the variables represented by means of three values: negative, infinitesimal, positive; from bottom to top: speed, acceleration, jerk, and snap.

DR AF

T-

of the total jerk (derivative of acceleration) (Flash and Hogan, 1985) and of the total mechanical action (Lebedev et al., 2001) as well as the relation between the instant speed and the local curvature, which should follow a power law (Lacquaniti et al., 1983). All these criteria are roughly equivalent (e.g., Wann, 1988; Flash and Viviani, 1995) and more or less accurately verified depending on the type of movement (e.g., Nagasaki, 1987). However, the metrics based upon these criteria give little insight regarding the detail of the irregularities or about their origin, central (Todorov and Jordan, 1998) or biomechanical (Gribble and Ostry, 1996; Schaal and Sternad, 2001) Frequency-domain analysis has been used for periodic irregularities, e.g., static or dynamical tremor (as in Beuter et al., 1999), periodic irregularities in slow finger movements (Wessberg and Vallbo, 1996), and for discrim-

inating tremor from myoclonus (Salazar et al., 2000). However, it is of little help when the anomalies do not correspond to a stationary distribution, for instance in the case of an asymmetrical velocity profile (Nagasaki, 1989) or for few, non-periodic irregularities. Non-periodic irregularities have been interpreted as the result of the segmentation of a movement into short sub-movements (Soechting and Terzuolo, 1987; Milner, 1992), which occurs with normal subjects even in the simplest tasks (pointing, Novak, 2002; moving at constant speed, Doeringer and Hogan, 1998). In motor disorders like HD or PD, the movements appear to be excessively segmented in a variety of conditions. Whatever the origin of this hyper-segmentation (see for instance Hallett and Koshbin, 1980; Sheridan and Flowers, 1990; Godeaux et al., 1992; Phillips, 1994; Berardelli et al., 2001; Pfann et al. 2001),

Neuroinformatics_________________________________________________________________ Volume 1, 2003

Event Identification with Qualitative Patterns _____________________________________________________3

RR

EC

TE D

tion (e.g., second order derivatives: Hildreth, 1983; Canny, 1986; derivatives of any order: Chen, 1992) or contour recognition (e.g., Bengtsson, 1991). The present technique generalizes and extends these works by handling sequences of (rather than single) singularities and by the use of up to the fourth derivative (rather than by that of a unique derivative in all the situations). Another novel aspect of our work is the representation of discontinuities or steep variations by means of qualitative values (q-values): unbounded negative, finite negative, infinitesimal (i.e., null), finite positive, and unbounded positive (derived from Fimbel, 2002). For instance, a step in the linear speed corresponds to an unbounded acceleration. The present qualitative pattern-matching technique allows identifying predefined sequences of singularities called qualitative patterns (q-patterns). The q-patterns are sequences of qualitative states (q-states), defined by the qvalues of the speed and the higher order derivatives. Each q-state may have a variable duration (within some limits) and terminates when a singularity, i.e., a qualitative change in some derivative(s), occurs. During the pattern matching, the signal is iteratively filtered and differentiated up to the fourth derivative. The instantaneous qualitative values of the derivatives are calculated and compared to the qstate. Amismatch indicates the end of a q-state. The new q-values must match the next q-state, otherwise the q-pattern is not recognized. Because the q-patterns are symbolic descriptions, qualitative pattern-matching differs from numerical pattern-matching in several aspects. First, it focuses on the singularities, thus avoiding a complete comparison of the pattern and the signal. Second, the time intervals between the singularities (i.e., the durations of the qstates that separate the singularities) are weakly constrained and independent from one another. Therefore time warping is not required. Third, the q-patterns do not allow

{FIGURE 1 CALLOUT}

DR AF

T-

UN

CO

the resulting movement can be described as a sequence of simple, regular sub-units separated by discrete singularities termed motor events. A variety of indicators have been used for quantifying motor events: local extrema (peaks) of speed or acceleration (Eichhorn et al., 1996; Cirstea and Levin, 2000), zero-crossings of acceleration (Phillips, 1994), number of halts, zero-crossings of speed and jerk (Cobbah and Fairhurst, 2000), and zero-crossings of snap, i.e., the derivative of jerk (Novak, 2002). In spite of capturing only part of the motor events, these indicators may reflect more accurately the subject’s motor condition than global smoothness metrics. For instance, the number of speed peaks provides a better assessment of post-stroke recovery than the total jerk (Rohrer, 2002). Although simple statistics like the number or the density of events are indicators of irregularity, the timing of the events is required for a better understanding of the deficits in motor control. For instance, a higher density of events towards the end of a movement may indicate abnormal corrections in HD (Smith, 2000). Also, identification techniques may be required for classifying events of different nature and/or origin (e.g., motor commands or biomechanical constraints like joint stiffness). The identification technique presented here assumes that the events are sequences of qualitative changes (also termed zero-crossings, singularities, or sign inversions) in the speed and/or the higher derivatives of the movement: acceleration, jerk, and snap. For instance, a speed peak is a change from a positive to a negative acceleration. Although the generality of this definition may be debated, it is sufficient for representing all the aforementioned types of events. Figure 1 illustrates how small irregularities in the linear speed of a rapid planar movement generate qualitative changes in the jerk and the snap. The use of the signal derivatives is wellestablished in image processing, for edge detec-

Neuroinformatics_________________________________________________________________ Volume 1, 2003

4__________________________________________________________________________________Fimbel et al.

TE D

Table 1 Glossary of Terms, Acronyms, and Notations position , i.e., length of the trajectory up to the current time. speed, i.e., derivative of P acceleration, i.e., derivative of S jerk, i.e., derivative of A snap, i.e., derivative of J

q-value q-state q-pattern

qualitative value. qualitative state. A state where data and derivatives have constant q-values qualitative pattern. A sequence of q-states

PD HD SD EMG

Parkinson’s disease Huntington’s disease atypical Parkinsonian syndrome surface electromyogram.

CO

RR

EC

P S A J N

DR AF

T-

UN

quantifying degrees of resemblance (e.g., correlation between pattern and signal), therefore they favor binary pattern-matching, i.e., the signal and the pattern match or not. Qualitative pattern-matching uses thresholds (for calculating the q-values), as do all the aforementioned event-based approaches. In general terms, the thresholds are set so that the noise cannot cause artifacts. Thus they indirectly depend on the filtering parameters. However, in movement analysis, the measurement systems generally remove the noise, and therefore little importance is given to thresholding and filtering in the literature. As a consequence, it is often unclear whether the same results would have been obtained with different thresholds and filtering parameters. Indeed, the thresholds and the filtering control the scale of analysis, i.e., the amplitude and/or the duration of the detectable events. We posit here that in order to provide significant results, the scale of an event-based analysis should be adjustable and that the results should be reproduced at different scales to ensure that they are not accidental (Lowe,

1985). Multiscale analysis can be easily performed with the present technique. First, the thresholds are calibrated. The simple method proposed here is to use the results of human judges to calibrate the algorithm (however more sophisticated methods are equally possible). Then, the sensitivity (i.e., the scale of analysis) can be adjusted by means of a single filtering parameter. The technique is presented in the context of an analysis of the irregularity of rapid movements in Parkinson’s disease. “Method” presents this analysis in more detail, specifically the q-patterns, the pattern-matching algorithm, and the multiscale analysis. “Experimental Study of Rapid Movements in PD” are discussed later. The results of this study are explored in “Discussion.”

Method Qualitative Patterns The q-patterns are sequences of qualitative states defining the q-value of the position, speed, acceleration, jerk, and snap. The possi-

Neuroinformatics_________________________________________________________________ Volume 1, 2003

RR

EC

TE D

Event Identification with Qualitative Patterns _____________________________________________________5

CO

Fig. 2. Example of q-pattern: a sharp speed peak. Left: definition.The q-pattern has 3 q-states and 2 transitions. Notice that the second q-state (peak) may be absent (duration ≈). Right: example of event identified as a sharp speed peak. S: speed; A: acceleration.

DR AF

T-

UN

ble q-values are unbounded negative (—), finite negative (–),infinitesimal (≈), finite positive (+), unbounded positive (++), and the unconstrained q-value (?), i.e., the union of all the qvalues. These values are part of an extension of the sign-based algebra (negative, null, positive) (e.g., de Kleer, 1984) to infinitesimal and infinite (unbounded) values, following the rules of calculus of nonstandard analysis (e.g., Keisler, 1994). The duration of each q-state is defined separately. It characterizes the number of consecutive samples where the signal and the q-state match, given as a range [min samples, max samples]. For convenience, the duration can be defined in four different ways: 1) [a, b], i.e., direct definition of the range, 2) (≈), i.e., infinitesimal duration, equivalent to [0, 1] (the qstate appears on a single sample or between two samples); 3) (+), i.e., finite duration, equivalent to [1,+∞]; 4) (?), i.e., any duration, equiv-

alent to [0, +∞]. Figure 2 gives an example of q-pattern. During the pattern matching, each q-state qi of a q-pattern behaves like a finite-state automaton defined by two variables: a counter cnti of the number of samples matching the q-state (ending with the current sample; cnti is zero if the current sample does not match) and a flag reci indicating whether qi is recognized or not (reci is true if cnti is within the range of duration of qi). The pattern matching examines each sample of the signal, and updates the active q-states. The first q-state of a q-pattern is always active, and when a q-state is recognized (i.e., reci is true), it activates the next qstate. A q-state becomes inactive when a mismatch occurs (i.e., cnti is set to zero). An occurrence of a q-pattern is found when its last q-state is recognized. The occurrences of the q-patterns may be totally or partially overlapped in the signal.

Neuroinformatics_________________________________________________________________ Volume 1, 2003

Fig 2

CO

RR

EC

TE D

6__________________________________________________________________________________Fimbel et al.

For instance, Fig. 3 shows two patterns that overlap completely. This occurs when some q-state contains the unconstrained q-value (?), or because the durations are ambiguous, e.g., (+) vs (≈). It is also possible that a complex q-pattern overlaps several simple q-patterns. For instance, a double peak (“M-shape”) overlaps two single peaks. In this case, the total number of occurrences of q-patterns is no longer significant. In practice, it may not always be possible to define a set of q-patterns that is orthogonal, i.e., such that a given signal matches one and only one of these q-patterns, therefore the existence

DR AF

{FIGURE 3 CALLOUT}

T-

UN

Fig. 3.Weighting overlapped occurrences. Left: non-orthogonal q-patterns (sharp and flat speed peaks; the difference is the duration of the 2nd. q-state). Center: an event corresponding to overlapped occurrences. S: speed A: acceleration. Right: weights of the samples. Double boarder = extension (the same for both occurrences); transitions = 3–4 and 4–5. From left to right: sample, q-value of speed s and acceleration a, number of occurrences containing the sample coI, number of adjacent transitions atI, weight wI=atI/2coI. Total weight=2, i.e., 2 significant transitions in the signal.

of overlapping occurrences of q-patterns should be considered the general case. Here, we propose to count the transitions rather than the q-patterns. However, the transitions themselves may not be clearly separated (e.g., in the case of infinitesimal durations). The solution adopted here is to weight them, as seen below. We define the location of a transition as the two adjacent samples at which the transition takes place. The extension of an occurrence of q-pattern is then defined as the union of the locations of its transitions. The extension is generally discontinuous because it does not con-

Neuroinformatics_________________________________________________________________ Volume 1, 2003

Event Identification with Qualitative Patterns _____________________________________________________7

TE D

2. It adjusts the thresholds of the qualitative values (infinitesimal and unbounded) separately for each derivative. This is necessary because the derivatives have very different maximal amplitudes. The thresholds θ0 and θ∞ initially given for the speed s are scaled by the maximal amplitude of each higher order derivative v as follows:θv0=θ0×|vmax|/ |smax|, θv∞=θ∞×|vmax|/|smax|.

RR

EC

3. For all the samples of the recording, it calculates the q-values of the derivatives (using the foregoing thresholds). The q-value of a derivative v is unbounded positive iff. its absolute value is above θv∞, finite positive between θv0 and θv∞, infinitesimal between –θv0 and θv0, finite negative between –θv∞ and –θv0 , and unbounded negative otherwise. 4. It tries to match all the q-patterns for the N samples. As seen in the previous section, the q-states qk of a q-pattern behave like automata whose state variables (cntk and reck) are stored for the N samples. For each sample, the active q-states of each q-pattern are examined and matched against the current sample. When the last q-state of a q-pattern is recognized, an occurrence of the q-pattern is found, and stored in a two dimensional array (rows=samples, columns=q-patterns). Then the weights of the samples of its extension are updated.

UN

CO

tain the samples between two consecutive transitions. We compute a weight wi representing the number of transitions adjacent to the sample i as follows. Let noi equal the number of occurrences whose extension contains i, and ati equal the number of transitions adjacent to i (between i–1 and i or between i and i+1). We define wi=ati/2noi (or zero if noi is null). Notice that for overlapping transitions, wi is no longer an integer. The sum of the weights represents the number of significant singularities (i.e., the transitions belonging to some q-pattern), notwithstanding the existence of overlapping q-patterns or of q-patterns of variable length (see Fig. 3). The q-patterns are symbolic representations, independent of the amplitude and other numerical features of the data, and relatively independent of the time scale. Therefore, each q-pattern matches a broad class of signals, which is to say that it contains little information. This allows for computationally efficient pattern-matching algorithms such as the one presented in the next subheading.

Qualitative Pattern-Matching Algorithm

T-

The input of the algorithm is a movement recording of N samples, composed of a planar position x–y, and a time stamp t. The algorithm matches M q-patterns and determines all their occurrences. It returns the total weight of the samples as an index of irregularity. We now give a brief description of the steps of the algorithm. The detailed algorithm is presented in Figs. 4 and 6 (a program is available at: http://www.ele.etsmtl.ca/profs/efimbel).

DR AF

{FIGURE 4 CALLOUT}

1. For all the samples of the recording, the algorithm calculates the linear displacement p and the derivatives (speed s, acceleration a, jerk j, snap n). Each derivative is smoothed (lowpass filtered) before the next one is calculated. This iterative smoothing prevents the amplification of high-frequency noise by the derivation operator.

5. The algorithm returns the sum of the weights of the samples. The variables and parameters used by the algorithm are now presented. In the first place, the input variables are as follows. • x [1..N], y [1..N], t [1..N]: instant planar coordinates and time stamps The signal is represented by the following variables, calculated from the input. • rawP[1..N], filteredP[1..N]: position, i.e., linear displacement, raw and smoothed. • rawS[1..N], filteredS[1..N], qS[1..N]: speed (raw, filtered and q-value) • rawA[1..N], filteredA[1..N], qA[1..N]: acceleration (raw, filtered and q-value)

Volume 1, 2003_________________________________________________________________ Neuroinformatics

TE D EC RR CO UN TDR AF Fig. 4. Qualitative pattern-matching algorithm. Comments are in italics.

EC

TE D

Event Identification with Qualitative Patterns _____________________________________________________9

eters. • τ: time constant for smoothing the variables (acts as a sensitivity parameter). • θ 0, θ ∞ : thresholds of the infinitesimal and unbounded q-values (for the speed s). • θs0, θs∞, θa0, θa∞, θj0, θj∞, θn0, θn∞,: specific thresholds calculated from θ0, θ∞, for each derivative.

DR AF

T-

UN

CO

• rawJ[1..N], filteredJ[1..N], qJ[1..N]: jerk (raw, filtered and q-value) • rawN[1..N], filteredN[1..N], qN[1..N]: snap (raw, filtered and q-value) The q-patterns are noted qpj, j=1..M and the q-states of the q-pattern qpj are q1..qkj. Each q-state qk is represented by the following variables. • minDk, maxDk: duration, i.e., minimal and maximal number of samples. • qSk, qAk, qJk, qNk: desired q-values of speed, acceleration, jerk and snap. The state variables used in the pattern matching are represented as follows. • cntk[1..N], reck[1..N]: counter and flag (recognized) of q-state qk for samples 1..N. • actk: flag (active) for q-state qk. The algorithm uses the following variables for calculating the weights. • no[1..N]: number of occurrences whose extension contains a sample. • at[1..N]: number of transitions adjacent to a sample. The output variables of the algorithm are the following. • occ[1..N,1..M]: occurrences of q-patterns. occ[i,k]=1 is an occurrence of qpk ends at i, 0 otherwise. • weight[1..N]: weight of the samples. The algorithm also uses the following param-

RR

Fig. 5. Convolution kernel used for smoothing the variables. Vertical: coefficients applied to the samples during the smoothing. Horizontal: samples (n=current sample). Left τ=3.5. Center: τ=0.5. Right: τ=0.

Remark

For simplicity, we suppose that the boundaries of any array var[1..N] are handled as follows. var[I]= var[1] for any I ≈ 0 and var[I] = var[N] for any I ≈ N.

Initial State of the Variables The arrays x, y, and t contain the movement recording and any other variable is set to zero. This algorithm uses four external functions: qValue, match, smooth, and markOccurrence that are now presented.

qValue(x, θs0, θs∞: real): q-value This function calculates the q-value corresponding to x as seen before.

match(qS1, qA1, qJ1, qN1, qS2, qA2, qJ2, qN2: q-values): boolean This function checks that qV1 = qV2 or qV2 =?, for V=(S, A, J, N). When this occurs, the func-

Volume 1, 2003_________________________________________________________________ Neuroinformatics

UN

CO

RR

EC

TE D

10 ________________________________________________________________________________Fimbel et al.

Fig. 6.Algorithm for marking an occurrence of q-pattern. Comments are in italics.

T-

tion returns true. Otherwise there is a mismatch and the function returns false.

smooth(v[1..N]: real, i: integer, τ: real): real

DR AF

{FIGURE 5 CALLOUT}

This function calculates a filtered value for sample I by means of a symmetrical trapezoidal convolution (Fig. 5). The integer part int(τ) defines the size of a rectangular window centered on i. The total width of this window is 1+2×int(τ) and all its samples have the weight 1. For instance with τ=1, the window is [i–1, i, i+1], and with τ=0, the window is reduced to [i] i.e., there is no smoothing. The decimal part dec(τ) defines the weight of the two samples surrounding the window. For instance, with τ=1.5 the window

is [i–1, i, i+1] and the convolution is 0.5×v[I–2] +(v[i–1]+v[i]+v[i+1]) + 0.5×v[i+1]/ (0.5+3+0.5). The complete formula is as follows. (∑j=i–int(τ)..i+int(τ)v[j]+dec(τ)×(v[i–τ–1]+v [i+τ+1]))/(2τ+1) The smoothing function is a low-pass filter whose cut-off frequency decreases with τ (e.g., τ=2: 9 Hz; τ=4: 5 Hz). Because the convolution kernel is symmetrical, the smoothing does not shift the derivatives (which would impair detecting qualitative transitions involving several derivatives). Moreover, it is easy to calibrate because its response varies slowly with τ.

markOccurrence(no[1..N], at[1..N], i, j: integer)

Neuroinformatics_________________________________________________________________ Volume 1, 2003

Event Identification with Qualitative Patterns ____________________________________________________11

TE D

EC

Multiscale Analysis

Multiscale movement analysis may be done as follows. 1) Determine the thresholds so that the results of the algorithm correspond to those of human judges on a sample of data. 2) Execute the algorithm on the sample of data with different values of the sensitivity parameter and determine the intervals of stability, in which the results are nearly constant. 3) Analyze the complete data with several values of the sensitivity parameter, one for each interval of stability.

T-

UN

CO

The algorithm works in linear time: 1) the computational time is proportional to the number of samples, and 2) in the worst case, the time increases linearly with the total size of the q-patterns. Formally, its time complexity is O(L×N), where L is the total number of q-states of the q-patterns and N is the number of samples. Linear algorithms are very robust (scalable), because their computational time increases smoothly with the number of patterns and/or the size of the data. Therefore the technique may afford high numbers of q-patterns, high sampling rates, and long duration recordings. Moreover, real-time analysis on downsized computing devices such as microcontrollers or digital signal processors (DSP) will presumably be easy.

body requires a constant number of operations, and it is executed for all the q-states and all the samples, i.e., L×N times. Consider now the function markOccurrence (line 12). Its internal loop (line 15) requires a constant number of operations, and is executed for each q-state of the q-patterns that have been recognized, thus, at most L×N times. Therefore, the overall execution time is in O(L×N).

RR

Computational Complexity

DR AF

{FIGURECALLOUT}

This function marks the extension of the occurrence of q-pattern qpj that has been detected at sample i. It updates the arrays no[1..N] (number of occurrences containing each sample) and at[1..N] (number of transitions adjacent to each sample) for all the samples of the extension. This function examines all the qstates of qpj starting from the last one. For each q-state qk, it determines the sample where the q-state started to match (using the counter cntk[I] to move backwards) and it marks the two samples surrounding the transition. The current sample is set just before the transition and the process is iterated up to the first qstate. The corresponding algorithm is presented in Fig. 6.

Determination of the Thresholds θ0, θ∞

A representative set of movement recordings r1..rk is chosen. The parameter τ is set to a small value τ0 (i.e., a high sensitivity), the filtered speed profiles are printed and the number of singularity of each recording ns1..nsk are determined by human judges. Then θ0 and θ∞ are determined by means of a trial-and-error process minimizing ε=[∑I=1..k(wi–nsi)]1/2, i.e., the least-squared difference between the indices of irregularity computed by the algorithm w1..wk and ns1..nsk.

Sketch of Proof

Exploration of the Sensitivity Parameter τ

Notice that the convolution used in the function smooth is updated for each sample with a constant number of operations. Therefore, the instructions that are executed more often are the body of the triple loop of lines 4, 5, and 6 (samples, q-patterns, q-states). Consider first the internal loop on the q-states (line 6). Its

Once θ0and θ∞ are fixed, the algorithm is executed on the sample recordings r1..rk, and the total index of irregularity σ=∑I=1..kwi is plotted as a function of τ. The intervals of stability [τl, τu] where σ is nearly constant are determined. Each of them corresponds presumably to a range of scales (i.e., a combina-

Neuroinformatics_________________________________________________________________ Volume 1, 2003

12 ________________________________________________________________________________Fimbel et al.

DR AF

T-

UN

CO

RR

EC

Data Analysis One value τi is chosen within each interval of stability I and the whole set of recordings is processed for all these values τ1.. τn. As said before, there are few (or no) events whose scales are within the intervals of stability. Therefore, the scales of the events are more likely to be found in the ranges between consecutive intervals. Notice that the number of events whose scale is between τI and τi+1 (or to be more precise, the events that are detected with sensitivity τi but not with sensitivity τi+1) can be obtained by the subtraction σ(τi)–σ(τi+1). The foregoing method is not always applicable. For instance, when the amplitude of the noise is too high (e.g., the signal-to-noise ratio is lower than 1), no filtering method permits a visual detection of the events. Also, the algorithm must be able to provide the same results as the visual detection by changing θ0 and θ∞. Whether or not this is the case must be determined empirically. Finally, when the scale of the events is scattered rather that grouped into well-delimited clusters, there will be no intervals of stability. For instance, drawing movements under visual feedback presumably contain corrections of all the amplitudes (i.e., events of all the scales) and the results of the algorithm may be threshold-dependent (thus not very informative). However, these problems do not occur in the application presented in the next subheading.

ments (Flash et al., 1992), unconstrained, remembered-aiming movements (Poizner et al., 1998), rapid reach-and-grasp movements (Alberts, 2000), sequences of pointing movements (Phillips, 1994), handwriting (Cobbah and Fairhurst, 2000), ballistic arm movements with indirect visual cues (Hallett et Koshbin, 1980), rapid elbow flexion (Pfann et al., 2001), and slow arm movements (Isenberg, 1994). The present study is aimed to assess irregularity in the conditions where patients with PD are at their best, i.e., rapid externally triggered movements (for avoiding akynesia, Siegert et al., 2002), with no visual feedback (for avoiding the reliance on visual cues, Flash et al., 1992), no target (therefore a compensatory strategy for imprecise movements is not required, Sheridan and Flowers, 1990), and a relatively small amplitude (therefore there is no need for complementary EMG activity, even if the initial EMG burst is abnormal, Hallett and Koshbin, 1980; Berardelli et al., 2001). A preliminary study showed that even in these conditions, patients with PD seem to present irregular movements (Fimbel et al., 1999).

TE D

tion of amplitude and duration) where there are few events. Otherwise σ would vary strongly within the interval.

Experimental Study of Rapid Movements in PD Aim of the Study

The lack of smoothness of the movements in PD has been consistently reported in a variety of conditions, e.g., rapid-aiming move-

Subjects Participants included 18 patients (9 males) diagnosed with idiopathic Parkinson’s disease (PD) or atypical Parkinsonian syndrome (PS), at stages 2, or 3 on the Hoehn et Yahr scale, under normal l-dopa medication, aged between 42 and 80 years (mean 61.8, SD 10.5) and 18 control subjects (9 males) aged between 40 and 80 years (mean 59.9, SD 12.05) with no history of neurological problems. Preceding the test, patients were rated on the United Parkinson’s Disease Rating Scale. Prior informed consent was obtained from all subjects. All subjects were right-handed (as determined by means of the Edinburgh handedness inventory) and used their right hand to perform the task.

Neuroinformatics_________________________________________________________________ Volume 1, 2003

EC

TE D

Event Identification with Qualitative Patterns ____________________________________________________13

Task

RR

Fig. 7. q-patterns used for detecting inversions of the jerk.

Qualitative Patterns

Four q-patterns were used: flat and sharp speed peaks, left and right s-shapes (Fig. 7). {FIGURE 7 These patterns allow detection of all the pos- CALLOUT} sible singularities generated by inversions of the jerk (from positive to negative and vice versa) during ballistic movements (in which the EMG only presents an initial triphasic burst, Hallett et Koshbin, 1980). An inversion of the jerk may produce different patterns depending on its timing and its amplitude (left s-shapes or peaks during the acceleration phase, right s-shapes or peaks during the deceleration phase, peaks during the two phases).

DR AF

T-

UN

CO

Subjects were instructed to draw a diagonal line as fast as possible on a digitizing pad (Wacom A4, X-Y spatial precision 0.025 mm; sampling rate 100 Hz). The movements corresponded to planar forearm extensions (average amplitude 20 cm, angular variation: 30±5° depending on the forearm’s length), with the elbow resting on the plane and remaining partly flexed at the end of the movement. The movement had to start after a signal (visual and auditory: 800 Hz, approximate level 80 db). The signal was emitted after a random delay (generated from an exponential distribution that prevents any anticipation). Subjects were instructed to maintain their fixed gaze on the signal device placed in front of their head in order to avoid visual feedback. In case of overshooting, the subject was verbally instructed to reduce the amplitude; undershooting was left uncorrected. Previous training allowed the subjects to remain within the pad without visual feedback. After the training, subjects performed between 15 and 30 movements depending on their motor condition and their degree of fatigue.

Determination of the Thresholds Eight movements were selected, 6 from a patient with PD and 2 from a control (Fig. 8). The sensitivity parameter τ was fixed to 0.05, sufficient for removing the remaining digitizing noise. The filtered speed profiles were printed and given to four human judges who had to count the patterns of Fig. 7 present in each profile (results shown in Table 2). The threshold θ∞=1000 was chosen in order to avoid unbounded q-values (not used in the

Neuroinformatics_________________________________________________________________ Volume 1, 2003

{FIGURE 8 CALLOUT}

{TABLE 2 CALLOUT}

14 ________________________________________________________________________________Fimbel et al. Table 2—Number of Irregularities Detected by the Judges in the Sample Recordings (Fig. 9) judge 1

judge 2

judge 3

judge 4

average number of irregularities

standard deviation

1

3

3

3

4

3.25

0.50

2

13

13

14

14

13.50

0.58

3

2

2

2

2

2.00

0.00

4

3

3

3

4

3.25

0.50

5

23

23

23

23

23.00

0.00

6

18

17

16

17

17.00

0.82

7

9

9

8

9

8.75

0.50

8

5

5

5

5

5.00

0.00

UN

CO

RR

EC

TE D

Curve

{FIGURE 9 CALLOUT}

DR AF

T-

Fig. 8. Sample movement recordings used for determining the thresholds. Horizontal: time. Vertical: speed. Each speed profile has a different scale so that to normalize the height.

Fig. 9. Difference between the results of the algorithm and of the human judges. Difference ε(θ0)=(Σi=1..k(wi–nsi))1/2 in percentage as a function of θ0(horizontal:logarithmic scale) for τ=0.05 i.e.,a minimal filtering (wi:index of irregularity;nsi:average number of singularities detected by the judges).

patterns). An exploration of the threshold θ0 was performed (θ0=10–6 to 1, i.e., θ0=10–k for k=–6 to 0 by steps of 0.5), and the least squared error between the results of the algorithm and those of the judges was plotted as a function of θ0 (Fig. 9). The error is almost constant (1.7%) for θ0<=10–3, hence θ0=10–4 has been chosen.

Determination of the Intervals of Stability The algorithm was executed for values of τ in the range [0, 20] (by steps of 0.25 between 0 and 2; by steps of 0.5 between 2 and 5; by steps of 1 between 5 and 20). The sum of the indices of irregularity returned by the algorithm (σ=Σi=1..8 {FIGURE 10 CALLOUT} wi) was plotted as a function of τ (Fig. 10).

Neuroinformatics_________________________________________________________________ Volume 1, 2003

EC

TE D

Event Identification with Qualitative Patterns ____________________________________________________15

UN

CO

RR

Fig. 10. Intervals of stability. Vertical: σ(τ), i.e., sum of the indices of irregularity for the 8 sample recordings for θ0=10–4 as a function of τ. Left: τ in the range [0..5], Right: τ in the range [0..20]. Thick horizontal bars: intervals of quasi-stability (where σ(τ) varies of at most 5%).

DR AF

T-

Four intervals of stability can be observed: [0, 1], [1.75, 4], [4.5, 6], [8,..]. The intervals [4.5, 6] and [8,..] were discarded because the sensitivity was too low (only the main speed peaks and a few singularities were detected), and one value was chosen in each of the remaining intervals, i.e., τ =0.25 and τ =2.

Data Analysis

The reaction time (RT), the length of the trajectory (LT), the movement time (MT) and the maximal speed V max were computed. Incomplete (i.e., final speed above 20% of the maximum speed Vmax) and anticipated movements (i.e., 5% of Vmax reached before the start signal) were discarded. The main speed peak

Fig. 11.Average index of irregularity per movement for each subject. Left: controls, Right: patients with PD or PS.

was determined (threshold: 1% of Vmax) and the rest of the recording was truncated. The index of irregularity was computed with two values of sensitivity: τ =0.25 (AE, “all events”) and τ =2 (LE, “large events”). The number of “small events” (SE) was obtained by subtraction (SE=AE–LE). Analysis of variance (ANOVA) between controls and PD (and PS) were performed for RT, LT, MT, Vmax, AL, LE, and SE.

Results The reaction time (RT) and movement length (ML) were similar in the two groups. However, patients with PD and PS reached lower maximum speeds Vmax, F(1,34)=15.47, p=0.00039,

Neuroinformatics_________________________________________________________________ Volume 1, 2003

16 ________________________________________________________________________________Fimbel et al.

RR

EC

TE D

Indeed, a mere detection of jerk singularities (e.g., Cobbah and Fairhurst, 2000) is enough for quantifying the irregularity of a movement. However, it does not allow for precise analysis of the anomalies. For instance, the main speed peak (which is a mere consequence of the cinematic movement, hence uninformative) would be counted as a jerk singularity, whereas it would be easy to identify by means of a q-pattern. Moreover, the significance of the jerk has been questioned, and it has been argued that the snap is a better indicator (Wiegner and Werblinka, 1992). The present techniques do not rely on any theory or model of human movement, and the choice of an indicator (or a combination of indicators) is left to the modeler. The definition of the q-patterns is clearly the bottleneck of the technique, because they must account for both the underlying theory (if any) and for the variety of resulting motor events. It is worth noting that a unique phenomenon, like a jerk inversion, may take several forms, as illustrated by the three q-patterns of Fig. 7. It is thus predictable that many q-patterns will be necessary for describing all the events that may occur in complex movements. Although these q-patterns may be defined incrementally (because they are not necessarily orthogonal, i.e., “mutually exclusive”), the definition of large numbers of q-patterns may require unreasonable efforts. Automatic construction of q-patterns (e.g., automatic extraction from numerical recordings) may therefore help extend the application of the technique with complex movements. Nonparametric classification techniques like artificial neural networks (in the case of signal processing see Fa-Long, 1999) are also possible, although their results may be difficult to interpret in terms of movement analysis. To our knowledge, none of these techniques has been used yet for movement analysis.

CO

Discussion

T-

UN

The experiment presented here shows that significant results on rapid movements can be obtained by means of an event-based, qualitative, pattern-matching technique. The motor events present in these movements are easily detected and classified under the form of qpatterns. The technique is presumably applicable to any type of segmented movements, i.e., movements presenting discrete, nonperiodic singularities. For instance, it is currently used for analyzing alternating movements (unpublished observations). On the other hand, frequency-based analysis methods are less informative than eventbased approaches for this kind of nonstationary signal (although a frequency-based detection of singularities is still possible, e.g., Pawlak, 1992). Optimality criteria like the jerk or the mechanical action, or the two-thirds power law may be useful for quantifying a degree of optimality, but, as mentioned in the introduction, they do not provide information about the movement anomalies.

DR AF

{FIGURE 11 CALLOUT}

and consequently showed longer execution times (ET) than controls, F(1,34)=27.05, p=0.000009. The index of irregularity (LE) differed significantly between the two groups (4.74 for PD and PS vs 1.61 singularities(trajectory, F(1,34)=41.12, p=0.0000003). Figure 11 shows the average number of large events per movement for all the subjects. We verified that the effect was independent of movement duration by computing the number of events per second. This indicator also showed a highly significant difference between the two groups (7.30 for PD and PS vs 4.71 for controls, F(1,34)=23.43 , p=0.000028). The difference between patients and controls was also significant for the total number of events (AE) (10.23 vs 4.83, F(1,34)=32.43, p=0.0000021), and for the number of small events (SE) (3.32 vs 5.48, F(1,34)=22.84, p=0.000033).

Neuroinformatics_________________________________________________________________ Volume 1, 2003

Event Identification with Qualitative Patterns ____________________________________________________17

RR

EC

TE D

approximated from a sampled signal. In some cases, their numerical values may be unstable, i.e., varying with the sampling frequency (e.g., stiff systems, Petzold, 1998). The effect of this possible numerical instability on the q-values of the derivatives should be studied in more detail. In spite of the foregoing limitations, the qualitative patterns-matching technique presented here provides a simple, flexible tool for multiscale analysis of cinematic variables. This technique could also be useful for other types of signals (inasmuch as their features are adequate), and for other types of applications, such as biofeedback or direct interfaces.

Acknowledgments

This work was supported by the ETS grant FIR 1-8346-000. We thank Jean-Marc Beaulieu, Christopher Fuhrman, and Christian Gargour for their useful comments and suggestions.

DR AF

T-

UN

CO

The calibration of the thresholds also deserves discussion. We proposed to calibrate them by comparing the algorithm with human judges on a sample of data. This simple method has proven efficient in the case of rapid movement: the results of the judges are consistent (Table 2) and several intervals of stability can easily be determined. However, this method has obvious limitations. The most important is that it cannot work with high levels of noise (in which case the visual recognition of events may be impaired). Moreover, external factors (i.e., the scale of the printed curves, the resolution of the printer, etc.) should be carefully controlled. An automatic threshold calibration method (whether general or domain-dependent) may thus extend the field of application of the present technique. Finally, the use of qualitative (rather than numerical) patterns has to be justified. We posit that qualitative patterns contain minimal symbolic information, thus requiring less work from the modeler, and that this information is sufficient for efficient event identification. The algorithm presented here works in linear time. This is to be compared with the high computational cost of numerical pattern-matching techniques like time warping, correlation coefficients compensating for noise, offset and gain factors, or automatic fitting (e.g., for postsynaptic events identification, Clements and Bekker, 1997). Moreover, qualitative, symbolic descriptions like the q-patterns are presumably easy to understand and recognize, as the similarity of performance between the algorithm and the human judges seems to indicate. On the other hand, the present qualitative pattern-matching technique is limited when the signal-to-noise ratio is too low (because artifacts may appear), but also when the motor events present a continuum of amplitudes (because the results will vary with the sensitivity parameter). Another caveat is the systematic use of derivatives that must be

References Alberts J.L., Saling M., Adler C.H. (2000) Disruptions in the reach-to-grasp actions of Parkinson’s patients. Exp. Brain Res. 134, 353–362. Bengtsson A., Eklundh J-O. (1991) Shape representation by multiscale contour approximation. IEEE Trans. on Pattern Recognition and Machine Intelligence 13, 85–93. Berardelli A., Rothwell JC., Thompson P.D., and Hallett M. (2001) Pathophysiology of bradykynesia in Parkinson’s disease. Brain 124, 2131–2146 Beuter A., Edwards R., de Geoffroy A., Mergler D., and Hudnell, K. (1999) Quantification of neuromotor function for detection of the effects of manganese. Neurotoxicology 20, 355–366. Canny J. (1986) A computational approach to edge detection. IEEE Trans. on Pattern Recognition and Machine Intelligence 8, 679–698. Chen J.G., Ren X.H., Tian Q. (1992) One Dimensional Multiple Scale Edge Detection. in Proceedings of 11th IAPR International Conference on Pattern Recognition III (C), pp.435–438.

Neuroinformatics_________________________________________________________________ Volume 1, 2003

18 ________________________________________________________________________________Fimbel et al.

RR

EC

TE D

Hallett M. and Khoshbin S. (1980) A physiological mechanism for bradykinesia. Brain 103, 301–314. Hildreth E.C. (1983) The detection of intensity changes by computer and biological vision systems. Computer Vision, Graphics and Image Processing 22, 1–27. Isenberg C., Conrad B. (1994) Kinematic properties of slow arm movements in Parkinson’s disease. J. Neurol. 241, 323–330. Keisler H.J. (1994) The hyppereal line. in Real Numbers, Generalizations of the Reals, and Theories of Continua, (Erlich, P., ed), Kluwer Academic Publishers, Boston, MA pp. 207–237. Lacquaniti F., Terzuolo C. and Viviani P. (1983) The law relating the kinematic and figural aspects of drawing movements. Acta Psychologica 54, 115–130. Lebedev S., Tsui W.H., and Van Gelder, P. (2001) Drawing movements as an outcome of the principle of least action. J. Math Psychol. 45, 43–52. Lowe D.G. (1985) Perceptual Organization and Visual Recognition. Kluwer Academic Publishers, Boston, MA. Milner T.E. (1992) A model for the generation of movements requiring endpoint precision. Neuroscience 49, 487–496. Novak K.E., Miller L.E., and Houk J.C. (2002) The use of overlapping submovements in the control of rapid hand movements. Exp. Brain Res. 144, 351–364. Pawlak M. (1992) On the detection and measurement of discontinuities. in Proceedings of Pattern Recognition Conference 1992, Vol.III. Conf. C., pp. 378–381. Petzold L.M., Ascher U.R. (1998) Computer Methods for Ordinary Differential Equations and Differential-Algebraic Equations. SIAM editions. Pfann K.D., Buchman A.S., Comella C.L., and Corcos D.M. (2001) Control of movement distance in parkinson’s disease. Movement Disorders, 16, 1048–1065. Phillips J.G., Martin K.E., Bradshaw J.L., and Jansek, R. (1994) Could bradykinesia in Parkinson’s Disease simply be compensation? J. Neurol. 241, 439–447. Poizner H., Fookson O. I., Berkinblit M.B., Hening W., Fledman G., and Adamovich, S. (1998) Pointing to remembered targets in 3-D space in Parkinson’s disease. Motor Control 2, 251–277. Rohrer B.R., Fasoli S., Igo Krebs H., et al. (2002) Movement smoothness changes during stroke recovery. J. Neurosci. 22, 8297–8304.

DR AF

T-

UN

CO

Cirstea M.C. and Levin M.F. (2000) Compensatory strategies for reaching in strokes. Brain 123, 940–953. Clements J. D. and Bekkers J.M. (1997) Detection of spontaneous synaptic events with an optimally scaled template. Biophysics J. 73, 220–229. De Kleer J., Brown S. (1984) A qualitative physics based on confluences. Artificial Intelligence 24, 7–83. Doeringer J.A., Hogan N. (1998) Intermittency in preplanned elbow movements persists in the absence of visual feedback. J. Neurophysiol. 80,1787–1799. Eichhorn T.E., Gasser T., Mai N., Marquardt C., Arnold G., Schwartz J., and Oertel W.H. (1996) Computational analysis of open loop handwriting movements in Parkinson’s Disease: A rapid method to detect dopamimetic effects. Movement Disorders 11, 289–297. Fa-Long L., Unbehauen R., and Luo F-L (1999) Applied Neural Networks for Signal Processing. Cambridge University Press, Cambridge, MA. Fimbel E., Guerfali W., Joanette Y. Lecours A.R., Gilbert B., and Plamondon, R. (1999) Estudio del perfil de velocidad durante movimientos rápidos en sujetos Parkinsonianos. in Proceedings of sexto congreso latinoamericano de neurosicología, SLAN99 Varadero, Cuba, pp. 67. Fimbel E. (2002) An algebra for qualitative modeling. Technical report RT020206, Electrical engineering department, École de technologie supérieure, Montréal. Available at: http://www.ele.etsmtl.ca/profs/efimbel. Flash T and Hogan N. (1985) The coordination of arm movements: an experimentally confirmed mathematical model. J Neurosci. 5, 1688–1703. Flash T., Inzelberg R., Schechtman E., and Korczyn A. (1992) Kinematic Analysis of Upper Limb Trajectories in Parkinson’s Disease. Experimental Neurology 118, 215–226. Godeaux E., Koulischer D. and Jacquy J. (1992) Parkinsonian bradykynesia is due to depression in the rate of rise of muscle activity. Ann Neurol. 31, 93–100. Gonzalez, R.C., Woods, R.E. (2002) Digital Image Processing, Second edition. Prentice Hall, Upper Saddle River, NJ. Gribble PL, Ostry DJ. (1996) Origins of the power law relation between movement velocity and curvature: modeling the effects of muscle mechanics and limb dynamics. J Neurophysiol. 76, 2853–2860.

Neuroinformatics_________________________________________________________________ Volume 1, 2003

Event Identification with Qualitative Patterns ____________________________________________________19

RR

EC

TE D

dicts the speed profiles of complex arm movements. J. Neurophysiol. 80, 696–714. Viviani P., Flash T. (1995) Minimum-jerk, two-thirds power law, and isochrony: converging approaches to movement planning. J Exp Psychol Hum Percept Perform. 21, 32–53. Wann J., Nimmo-Smith I., and Wing A.M. (1988) Relation between velocity and curvature in movement: equivalence and divergence between a power law and a minimum-jerk model. J Exp Psychol Hum Percept Perform. 14, 622–637. Walker N., Meyer D.E., Smelcer J.B. (1993) Spatial and temporal characteristics of rapid cursor-positioning movements with electromechanical mice in human-computer interaction. Human Factors 35, 431–458. Wiegner A.W., Wierzbicka M.M. (1992) Kinematic models and human elbow flexion movements: quantitative analysis. Exp Brain Res. 88, 665–673.

DR AF

T-

UN

CO

Schaal S., Sternad D. (2001) Origins and violations of the 2/3 power law in rhythmic three-dimensional arm movements. Exp. Brain Res. 136, 60–72. Sheridan M.R., Flowers K.A. (1990) Movement Variability and BradyKinesia in Parkinson’s Disease. Brain 113, 1149–1161. Siegert R.J., Harper D.N., Cameron F.B., and Abernethy D. (2002) Self-initiated versus externally cued reaction times in Parkinson’s Disease. J. Clin. Exp. Neuropsy. 24, 146–153. Smith M.A., Brandt J., and Shadmehr R. (2000) Motor disorder in Huntington’s Disease begins as a dysfunction in error feedback control. Nature 403, 544–549. Soechting J.F. and Terzuolo CA. (1987) Organization of arm movements. Motion is segmented. Neuroscience 23, 39–51. Todorov E., Jordan M.I. (1998) Smoothness maximization along a predefined path accurately pre-

Neuroinformatics_________________________________________________________________ Volume 1, 2003