## Abstract

Outside the laboratory, human movement typically involves redundant effector systems. How the nervous system selects among the task-equivalent solutions may provide insights into how movement is controlled. We propose a process model of movement generation that accounts for the kinematics of goal-directed pointing movements performed with a redundant arm. The key element is a neuronal dynamics that generates a virtual joint trajectory. This dynamics receives input from a neuronal timer that paces end-effector motion along its path. Within this dynamics, virtual joint velocity vectors that move the end effector are dynamically decoupled from velocity vectors that do not. Moreover, the sensed real joint configuration is coupled back into this neuronal dynamics, updating the virtual trajectory so that it yields to task-equivalent deviations from the dynamic movement plan. Experimental data from participants who perform in the same task setting as the model are compared in detail to the model predictions. We discover that joint velocities contain a substantial amount of self-motion that does not move the end effector. This is caused by the low impedance of muscle joint systems and by coupling among muscle joint systems due to multiarticulatory muscles. Back-coupling amplifies the induced control errors. We establish a link between the amount of self-motion and how curved the end-effector path is. We show that models in which an inverse dynamics cancels interaction torques predict too little self-motion and too straight end-effector paths.

## 1. Introduction

Understanding how organisms generate voluntary, goal-directed movements is one of the most difficult problems in theoretical neuroscience and remains largely unsolved. At the core of the problem is the fact that movement necessarily involves multiple levels of neural control, which are tightly interrelated and none of which can easily be neglected. Minimally, understanding how the nervous system moves an effector toward a target entails understanding how movements are prepared, timed, and controlled. That movements are prepared, in part at least, before they are initiated is reflected in the fact that the time to movement initiation varies with the amount of information available about an upcoming movement (Rosenbaum, 1980; Ghez et al., 1997). Movement trajectories reflect global properties of the movement task from the very beginning, for instance, in that the end effector moves in the direction of the target from the start of the movement. At the neuronal level, activity in motor, premotor, and parietal cortex precedes movement initiation and depends on movement parameters such as the direction and extent of end-effector motion (Moran & Schwartz, 1999; Georgopoulos, 1995; Sergio & Kalaska, 1998; Cohen & Andersen, 2002). Movement plans may, however, be updated anytime before movement initiation, as well as during the movement (Goodale, Pélisson, & Prablanc, 1986). This occurs involuntarily with a delay of about 100 ms when movement targets are spatially displaced. Even relatively abstract codes for movement goals may be updated at longer delays (Pisella et al., 2000). Movements are timed in the sense that the effector is at the right location at the right time (Schöner, 2002; Warren, 2006). Timing is central to coordination, in which the trajectories of different effectors are kept aligned temporally. Interlimb coordination is important for both rhythmic (Turvey, 1990; Kelso, 1995) and temporally discrete motor acts (Kelso, Southard, & Goodman, 1979, 1984; Gracco & Abbs, 1986). Finally, physically moving an effector entails generating forces and torques, which accelerate and decelerate effectors. Movements are controlled in the sense that muscles and joints are harnessed to generate the desired physical trajectory of the effector system. This entails dealing with constraints internal to the biomechanical system such as inertia, interaction torques, and Coriolis forces, but also to external factors such as gravity or external force fields (Jordan, 1990).

A central difficulty to understanding movement generation comes from the fact that movement preparation, timing, and control are closely and mutually coupled. Changes in a movement task typically imply changes of the values of movement parameters, changes to the timing of the movement, and changes to the torques encountered during execution of the movement. Exposing an effector to an external force field induces adjustments not only at the control level, but also of movement timing and movement parameters (Gribble & Ostry, 2000). No single perturbation of a movement exists that would only affect one of the three aspects of movement generation. A mechanical perturbation, for instance, not only triggers control action but may also reset the timing of the movement (phase resetting; see Yamanishi, Kawato, & Suzuki, 1979) and lead to updated movement parameters, generating corrective actions. Theories of movement generation must therefore address integration across these different levels of motor control.

The integrative nature of movement generation is reflected in the tight coupling of the associated neuronal processes, which are distributed across many neuronal structures. The timing of motor acts is, for instance, affected by feedback loops through the motor cortex and the cerebellum, as well by spinal pattern generators and feedback loops back to the motor cortex (Houk & Wise, 1995; Barto, Fagg, Sitkoff, & Houk, 1999; Graziano, 2006; Desmurget & Turner, 2008). Parietal, premotor, and motor cortex are all involved in extracting movement parameter values, including their updating when sensory information about target locations or current hand position changes (Georgopoulos, Kettner, & Schwartz, 1988, 1993; Schwartz, 1993, 1994). Control in the face of external force fields implicates both the cerebellum and motor cortex (Wolpert, Miall, & Kawato, 1998).

Dynamical systems ideas seem particularly well suited to deal with such rich internal coupling. Dynamical field theory, which is based on the dynamics of neuronal populations in cortical and subcortical structures (Erlhagen, Bastian, Jancke, Riehle, & Schöner, 1999), has provided an account of how movements are prepared and updated (Kopecz & Schöner, 1995; Erlhagen & Schöner, 2002). The same framework has been successfully used to understand the processes governing the initiation and termination of movements (Kopecz, 1995; Trappenberg, Dorris, Munoz, & Klein, 2001; Wilimzig, Schneider, & Schöner, 2006). The timing and coordination of movements can be understood in terms of coupling among the neuronal dynamics controlling individual degrees of freedom (Schöner & Kelso, 1988; Kelso, 1995; Grossberg, Pribe, & Cohen, 1997; Schöner, 2002). The level of control is naturally described in the language of dynamical systems. While much work on control remains at an abstract, computational level (Todorov & Jordan, 2002), the perspective exists to link this language to neuronal structures (Scott, 2004; Shadmehr & Krakauer, 2008). Detailed neuronal grounding has been provided only for the simplest impedance properties of individual joint muscle systems (Feldman, 1966; Mussa-Ivaldi, Hogan, & Bizzi, 1985; Feldman, Adamovich, Ostry, & Flanagan, 1990).

Our strategy in this letter is to remove from movement systems the artificial constraints that induce unique relationships between the different levels of movement generation. The vast majority of experiments on human movement have made use of nonredundant effectors, in which a unique solution exists for the given motor goal. In most real-life movement tasks, however, the human movement system is kinematically redundant, that is, there are more, sometimes many more, mechanical degrees of freedom than task variables. For instance, reaching to a 3D position to point or position the hand requires only three degrees of freedom (or six if the hand is to be oriented in specific ways). But the upper arm alone contains seven degrees of freedom. Natural movement involves the shoulder blade and portions of the upper body, leading to 10 or more degrees of freedom (Yang & Scholz, 2005; Tseng, Scholz, & Hotchkiss, 2003). Similarly, in upright stance, only the Cartesian position of the center of mass is critical to remaining mechanically stable, but many joint angles contribute to setting that position (Hsu, Scholz, Schöner, Jeka, & Kiemel, 2007).

Only recently have movement scientists begun to systematically exploit the inherent redundancy of the movement apparatus to analyze how movement is controlled (Scholz, & Schöner, 1999; Scholz, Schöner, & Latash, 2000; Admiraal, Medendorp, & Gielen, 2002; Tseng, Scholz, Schöner, & Hotchkiss, 2003; Todorov & Jordan, 2002, 2003; Torres & Zipser, 2004). (For precursors to this work, see, e.g., Ivaldi, Morasso, & Zaccaria, 1988; Cruse, Brüwer, & Dean, 1993; Haggard, Hutchinson, & Stein, 1995.) Our strategy in this work is to make use of the inherent redundancy of effector systems to give the central nervous system some freedom for how to realize a given movement task. The idea is to learn from the decisions of the nervous system, from its choices to stabilize some variables and not others, and from the correlations and trade-offs induced by such choices (Schöner, 1995).

The goal of this letter is to develop a key element of an integrative theory of movement generation that addresses how multiple degrees of freedom are harnessed to achieve a motor task. The theory is framed within dynamical systems thinking; it is process oriented and consistent with physiological principles. Choosing the simplest redundant system accessible in both experiment and theory, we report experimental data on pointing movements in a two-dimensional end-effector plane with a four-degree-of-freedom arm. A neuronal dynamics forms the core of our theoretical model, from which we will derive kinematic features, including the end-effector paths and trajectories, the joint velocities, and, in particular, self-motion (the component of the joint velocity vector that does not move the end effector), and compare these to the experimental data. We will selectively manipulate our model to examine the contributions of different theoretical assumptions to the observed movement features. Finally, we will implement a number of alternative theoretical accounts, including models with inverse dynamics, and will show how these alternatives fail to explain the observed kinematic features of redundant arm movement. In a separate work, we will study movement variability and correlation, including the uncontrolled manifold structure of the variance of joint motion (Martin, Schöner, & Scholz, 2007).

## 2. Model

### 2.1. Survey of the Model and Its Neuronal Embedding.

To model the generation of discrete goal-directed movements, we need to specify the processes of movement preparation, movement initiation, timing, virtual joint trajectory formation, and the muscle joint and biomechanical dynamics. Our strategy is to use as much as possible available models for those parts of movement generation that have been studied previously and to thus limit innovation to the core problem of organizing the redundant degrees of freedom. Figure 1 provides a survey of the model, together with its neurophysiological embedding. The neuronal substrate for the processes governing movement generation is highly distributed, and many components of the neural systems have overlapping function. Because the mapping of motor function onto neuronal structures is not one-to-one, the figure provides only a first rough sketch, which we will elaborate on below.

Reaching requires information about the scene, which typically comes from the visual system and includes information about both the identity of objects (ventral stream) and about the pose and position of objects (dorsal stream) (Milner & Goodale, 1995). Pose information must be transformed from visual coordinates into body-related coordinates, a task achieved by neural populations in the parietal cortex (Andersen, Snyder, Bradley, & Xing, 1997). The parietal cortex is also involved in using visual and proprioceptive information to estimate the spatial position of the hand (Gréa et al., 2002). The neural processes underlying these transformations remain outside the scope of the present model (see Simmering, Schutte, & Spencer, 2008, for a dynamic field account of such spatial representations). The outcome of these processes is modeled simply by assuming that estimates of the spatial coordinates of movement targets as well of the initial spatial position of the end effector are continuously available.

Preparing a goal-directed movement involves extracting from such spatial information the values of the movement parameters, most prominently direction and extent. In the presence of multiple possible movement targets, movement preparation also involves selection of one particular movement and suppression of the neuronal activation representing distracter targets. Premotor and motor cortex are involved in both of these aspects of movement preparation (Georgopoulos, Schwartz, & Kettner, 1986; Cisek & Kalaska, 2005), as well as in integrating prior information about upcoming movements (Bastian, Riehle, Erlhagen, & Schöner, 1998, 2003). The model does not address the processes of movement preparation, which has been treated in previous work using the same theoretical framework (Kopecz & Schöner, 1995; Erlhagen & Schöner, 2002). Instead, the values of movement parameters' direction and extent, as well as the desired movement time, are assumed given.

Movement initiation and termination involves a wide array of brain areas, as observed, for instance, when brain activity during discrete movement is contrasted with brain activity during rhythmic movement of the same effector (Schaal, Sternad, Osu, & Kawato, 2004). A comparatively precise understanding of initiation and termination has been achieved for saccadic eye movements in terms of the interaction between neural populations responsible for fixation and populations responsible for generating visually induced saccades (Dorris, Pare, & Munoz, 1997). Goal-directed pointing movements show some of the same signatures, suggesting similar mechanisms (Bekkering, Pratt, & Abrams, 1996). Theoretical accounts for this mechanism have been provided in the language of neuronal dynamical fields in terms of inhibitory interactions among the two relevant neuronal populations (Kopecz, 1995; Trappenberg et al., 2001; Wilimzig et al., 2006). As the processes of initiation and termination of a discrete movement are not the main focus of this letter, we model these processes through a simplified version of these accounts, in which the two populations responsible for fixation and movement are described by two competing dynamical activation variables.

Motor and premotor cortex are involved not only in the specification of movement parameters, but also in the generation of the time courses of goal-directed movements (Graziano, Taylor, Moore, & Cooke, 2002; Hatsopoulos, Xu, & Amit, 2007). The precise manner in which this happens is still largely unknown. What is clear is that a closely coupled ensemble of neural populations in motor cortex and the cerebellum with thalamus as a way station are involved, the basal ganglia playing a regulatory role as well (Houk & Wise, 1995). Motor cortex also interacts with spinal pattern generators (Drew, Kalaska, & Krouchev, 2008; Poppele & Bosco, 2003). The importance of the cerebellum for the timing of voluntary movements is clear from neuropsychological studies (Ivry, 1997). In fact, one account for cerebellar function conceives of the cerebellum as a neuronal mechanisms for the precise measurement of time (Braitenberg, Heck, & Sultan, 1997). Another view is that the cerebellum is a predictor of the sensory consequences of motor commands, a form of forward model (Wolpert et al., 1998). This conception is not necessarily in conflict with the former view, as predicting the time course of movements is what a forward model is essentially about. Much of the evidence about internal models derives from studies of how movement generation adapts to unknown force fields. Because we do not address such adaptation, we can use a simplified model of movement timing. Our functional description of these distributed but closely coupled populations of neurons takes the form of a neuronal oscillator, which can be started and stopped by the initiation system (Schöner, 2002). The critical assumption is that the oscillator generates and predicts the time course of end-effector motion along its movement path. Evidence that movement timing resides at the task level comes, for instance, from studies in which the spatial pattern of coordination between end effectors, not the pattern of joint motions, determines the stability of relative timing (Mechsner, Kerzel, Knoblich, & Prinz, 2001).

The specific joint configuration used to realize a goal-directed movement is reflected in motor cortical activation (Scott & Kalaska, 1995). We believe that a core function of motor cortex is to transduce the neuronal trajectory that predicts task-level motion into joint-level activation patterns, although reciprocal coupling with downstream structures is likely to play an important role as well (Graziano, 2006). The core of our model describes this transduction as a neuronal dynamics of a virtual joint configuration vector, which is driven by the virtual end-effector trajectory generated by the timing system. Two central assumptions structure this system. First, we postulate that within this neuronal dynamics of joint configurations, task-relevant combinations of degrees of freedom that affect the end effector are decoupled from task-irrelevant combinations of degrees of freedom that do not affect the end effector. Second, estimates of the real joint configuration are assumed to couple into this neuronal dynamics within the subspace of task-irrelevant combinations of degrees of freedom. We refer to this input as back-coupling from the effector level. It updates the trajectory of the virtual joint configuration within the null-space of end-effector motion, enabling the neuronal dynamics to yield to changes of joint configuration that do not affect the real end-effector path. The neurophysiological basis for such back-coupling may come from proprioceptive coupling into spinal networks mediated by Ia-interneurons and Renshaw cells (see, e.g., a similar mechanism proposed by Latash, Shim, Smilga, & Zatsiorsky, 2005) as well from transcortical loops that remap the functional relationship between cortical and spinal networks (Graziano, 2006).

Any theory of movement generation must take into account that the ensemble of muscles that converge on a joint, together with local spinal as well as transcortical feedback loops, endows joints with impedance properties that have an impact on how effector systems respond to motor commands (Asatryan & Feldman, 1965; Feldman, 1966; Hogan, 1985; Mussa-Ivaldi et al., 1985; Feldman & Levin, 1995; Ostry & Feldman, 2003). We model the active generation of torques by these distributed neural systems of joint muscle control using a simplified version of an established model (Gribble, Ostry, Sanguineti, & Laboissière, 1998). All muscles acting on a joint are lumped together and controlled by the associated virtual joint position and velocity. Multiarticular muscles are modeled by coupling among joint torques. Finally, the kinetics and kinematics of the arm are modeled using standard techniques.

### 2.2. Task Setting.

Both model and experiment involve the task of moving a pointer tip (end effector) from a start location to a target location at a comfortable speed. Start and target locations lie in a two-dimensional horizontal plane. Configurations of the arm are also restricted to that plane, so that effectively only four degrees of freedom are available. In principle, each associated joint angle (sternoclavicular joint, shoulder, elbow, and wrist) can be moved individually, although these degrees of freedom are effectively coupled both mechanically and at the level of neural control. In this task, the effector system is redundant because only two degrees of freedom are required to achieve a particular location of the end effector in the plane (see Figure 2). Four different end-effector paths defined between three different starting locations and two different target locations sample the work space coarsely. In addition, two of the movements were performed with two different initial arm configurations, leading to six different movement conditions (experimental details are provided in section 3).

### 2.3. Movement Parameters.

Based on this task setting, we assume that estimates of the movement parameters are available. The initial position of the end effector and the target location determine the amplitude of end-effector motion, *U _{i}*, along the two Cartesian directions

*i*= 1, 2. The desired movement time,

*T*, is also assumed given.

### 2.4. Movement Initiation and Termination.

*u*

_{m}) and the resting or fixation state (

*u*

_{r}). These variables evolve as described by a competitive neuronal dynamics: Here,

*h*is a negative resting level of neural activation. Inputs

*I*

_{m}and

*I*

_{r}bias the system toward initiating movement or resting. These inputs depend on the predicted state of the end effector, and this dependence turns the movement state off at the end of the movement (see appendix D for details). Inhibitory coupling between these two activation variables, mediated by a sigmoidal nonlinearity, , makes that only one of the two variables can be activated at the same time.

### 2.5. Timing.

**u**(

*t*) = (

*u*

_{1},

*u*

_{2}), determines the virtual end-effector velocity,

**v**(

*t*) = (

*v*

_{1},

*v*

_{2}), through so that the virtual end-effector velocity tracks the timing signal. The indices refer to the two Cartesian components of the end effector, and β

_{v}is a positive constant. Although it is not neuronally realistic, we use the Hopf normal form (Perko, 1991) as the simplest mathematical representation of a stable limit cycle oscillator that stands for a class of neuronal dynamics that exhibit this type of solution (Schöner, 2002). For each (“excitatory”) component,

**u**, the Hopf equation contains a second (“inhibitory”) component,

**z**: Herein, the Hopf equation, generates a stable limit cycle solution with cycle time,

*T*= 2π/ω

_{h}, relaxation time, 1/2α

_{h}, and amplitude,

*U*. This oscillator is active while the initiation system is in the movement state (). When the resting state is activated (), the timing signal has a stable fixed point at

_{i}**u**= 0.

### 2.6. Neural Dynamics of the Virtual Joint Configuration.

**u**(

*t*), is transformed into a virtual joint trajectory, . This requires inversion of the Jacobian equation, where is the Jacobian matrix built from the partial derivatives of the kinematic model (see appendix A). Because the effector is redundant, the equation is not invertible, however. To understand the implications of redundancy, it is useful to visualize the kinematic model of the effector around a given virtual joint configuration, . Figure 2 illustrates how the same end-effector position can be achieved by multiple joint configurations. The ensemble of joint configurations that leads to the same end-effector position forms a manifold, sketched in Figure 2 on the right. This is the “uncontrolled manifold” (Schöner, 1995), shown in earlier work to structure the variance of multijoint movement (Scholz, & Schöner, 1999; Scholz et al., 2000). Any change of joint configuration along this manifold does not change the position of the end effector, while configuration changes away from the manifold do change the end-effector position.

Instantaneous joint velocity vectors that are tangential to the manifold generate so-called self-motions, that is, internal motions of the joint configuration that leave the end-effector position unchanged. At a given joint configuration, , the linear space spanned by vectors tangent to the manifold is the null-space of the Jacobian, (see Figure 2, right). Basis vectors of the null-space are nontrivial solutions of . In the present case, **v** is two-dimensional and is four-dimensional, so that the null-space is two-dimensional (except at singularities of the effector). Two basis vectors that span the null-space can be determined as the columns of the matrix, **E**, which solves . The null-space and its orthogonal complement, the range space of the Jacobian, divide the space of virtual joint velocities into two subspaces (see Figure 3, left). Note that these subspaces depend on the joint configuration at which the Jacobian is computed.

**s**, is added to the virtual end-effector velocity,

**v**. This augmented Jacobian equation is invertible: where the matrix,

**J**

^{+}, is the Moore-Penrose pseudoinverse (see, e.g., Murray, Li, & Sastry, 1994). This equation decomposes the joint configuration velocity, , into two components: is the range-space and the null-space component (see the left panel of Figure 3).

**E**

^{T}). The same kind of mechanism may occur at the level of joint velocities (second term). The real joint configuration must be sensed and estimated, leading to processing delays (index

*d*; see appendix D for details). This form of back-coupling of the real into the virtual joint configuration dynamics implies both stabilization of the joint configuration within the uncontrolled manifold (through the terms dependent on and ) and driving virtual self-motion (when the terms and are different from zero). The projection of the back-coupling term onto the null-space ensures that the dynamics within the space of self-motion depends on only the components of and within that subspace, so that the range-space and null-space remain decoupled.

**v**, and the self-motion velocity,

**s**, by virtual joint velocities using equations 2.5 and 2.6: To implement the model, the matrices , , , and are computed analytically.

### 2.7. Muscle-Joint Model.

*i*, can then be characterized by a single function, (listed in appendix C), where and are the real joint angle and velocity. At rest and in the absence of external forces, the muscle joint system is at equilibrium at

**T**=

**0**and . Depending on the time course of the virtual joint trajectory, and on the biomechanics of the arm, the realized joint trajectory may deviate significantly from the virtual trajectory. This is why taking into account the nonlinear dependence of muscle force generation on muscle state is important (Gribble et al., 1998).

Approximations consist of neglecting the time-delaying effect on force generation of the muscle calcium kinetics and rewriting the velocity-dependent function (damping), which led us to characterize the muscle active state directly by and following Lussanet, Smeets, & Brenner (2002). Our muscle model accounts for the dependence of active torque on the joint velocity () as in Gribble et al. (1998) and Hogan (1984). Because the damping effect of muscles is proportionally larger at low velocities than at high velocities (Gielen, Houk, Marcus, & Miller, 1984; Houk, Fagg, & Barto, 2002), a nonlinear velocity-dependent function is added. The main constraint for modeling damping is its role in leading to a smooth transition to rest at the end of the movement, which is not critical for our results.

**z**, between the physical torques at all joints,

**T**

_{m}, and the vector of muscle joint torques, , The matrix,

**z**, is listed in appendix C.

### 2.8. Arm Kinematics and Kinetics.

**x**= (

*x*,

*y*), to the joint configuration, : (the upper index

*T*indicates the transpose, so that the joint configuration is a column vector; the equations are listed in appendix A). The model is derived assuming an articulated rigid body with four revolute joints whose axes of rotation are perpendicular to the two-dimensional plane of motion.

**T**

_{m}is the vector of active torques generated at the skeleton joints by muscle forces (all terms are listed in appendix B).

### 2.9. Simulations.

The model was implemented in Matlab version 13 (2002) using the numerical Euler method to solve the differential equation. Appendixes E and C list the parameter values of the model used for all movements anywhere in the work space.

## 3. Experimental Methods and Analysis

Participants in the experiments were three healthy individuals from the University of Delaware community, 21 to 35 years of age. Participants gave informed consent before participation. All participants were right-handed and reached with their right arms to the targets.

### 3.1. Procedure.

Participants sat on chair with a high backrest. Their trunk movements were restrained by a harness about their chest that was attached to the chair. They wore a hand splint with a stylus, the tip of which was aligned with the position of an extended index finger. Reflective spherical markers (1.5 cm in diameter) were placed on bony landmarks of the participants: (1) the sternoclavicular joint, (2) just below the lateral tip of the acromion, (3) the lateral epicondyle of the humerus, (4) the radial styloid process of the wrist, and (5) the tip of stylus fixed to the hand brace. Markers were also placed at the centers of the targets. All markers lay in the same horizontal plane.

There were six conditions involving two target locations, positioned at 90% of each subject's arm length. The right target was located at an angle of 40 degrees to the right of a line passing forward through the right acromion process. The left target was positioned 40 degrees to the left of this line, requiring the participants to reach across their body. Two starting pointer locations were used in which the stylus marker was 7.8 cm anterior to either the sternum marker or the right acromion. The third starting location was selected to be at 20 degrees to the left and at 50% of each subject's arm length.

The table height was set so that the right arm rested on it in the horizontal plane when in the starting position. When the pointer was in the right-most and left-most starting position, two different initial joint configurations were used in separate conditions. Movement time (MT) was kept constant by using a Lafayette Instrument time that provided feedback after every trial. Subjects were asked to reach as quickly and accurately as possible. The movement time was determined during test trials, and that time was maintained for the actual experiments. Trials for which participants deviated by more than 5% from the target MT were repeated. The participants were instructed to reach the target with the tip of the stylus. Emphasis was on both spatial accuracy (“try to touch the center of the target”) and temporal accuracy (“try to reach the target in the designated time”). Additionally, they were instructed to perform the reaching with one continuous movement, with no intentional pauses during the movement.

Before the actual data collection session, participants were given a practice session with verbal feedback provided by the experimenters on spatial as well as temporal accuracies of reaching. In the actual data collection session, kinematic data of the markers were collected by the VICON-370 motion measurement system at a sampling frequency of 120 Hz. The collected data were then filtered in both directions using a fourth-order Butterworth low-pass filter with a cut-off frequency of 5 Hz. Each participant performed 25 reaching movements for each experimental condition.

In order to compute the mean end-effector trajectories and self-motion in experiments, movements must be matched in time from one trial to the next (see also Figure 6). The beginning of the movement was defined as the time when end-effector velocity first reached 1% of its peak value. The end of the movement was defined as the time when end-effector velocity fell below 3% of its peak. Trajectories are time-warped to match the mean movement time in both experiment and simulated data.

### 3.2. Self-Motion Analysis.

**s**, is the component of the joint velocity vector within the null-space of the Jacobian. It was computed by projecting velocity vectors, , onto the basis vectors

**e**

_{1}and

**e**

_{2}(that form the columns of

**E**) at each time sample,

*t*, and for each trial,

*n*: The associated range space component is The amount of self-motion and range-space motion was computed by taking the mean across trials of the lengths of these vectors. No normalization relative to the number of degrees of freedom was needed, as both subspaces have the same dimensionality of two.

## 4. Results

To compare experiment and theory, we look at features of end-effector paths, end-effector trajectories, and joint trajectories with particular emphasis on self-motion. In each case, we present experimental data from all participants to illustrate reproducible patterns. In addition to presenting reasonable fits of these data at one constant plausible set of parameter values of the model, we examine the role of different components of the model by making targeted changes to demonstrate how these affect simulated patterns of movements. These manipulations include, for instance, introducing a component that mimics perfect inverse dynamics or setting high-impedance value for the muscles.

### 4.1. End-Effector Paths.

The end-effector paths observed in experiment (left three panels in the top row of Figure 4) deviate from idealized straight line paths (Morasso, 1981). For most of the movements, the end-effector paths are slightly but consistently curved. This curvature depends on the position in the work space and is qualitatively reproducible across participants. Only movement 4 has an almost straight end-effector path. The variability from trial to trial for each participant is considerable, but does not wash out the consistent pattern, so that end-effector path curvature is a robust feature of these pointing movements. The model (see the top right panel of Figure 4) generates curved paths that qualitatively match the observed pattern in experiment except for movement 6, for which our parameter setting for impedance was too far off (see Figure 10 below, where increased impedance removes the problem). Recall that the simulations use a single reference parameter set (see appendix E), which was determined to be consistent with physiological values and to achieve a rough match of the experimental data. We did not explicitly fit data for every movement and every subject independently and thus neglected that impedance values may vary in space and time (Tsuji et al., 1995; Gomi & Osu, 1998) or across subjects. Note also that the model does not show constant terminal errors, as movement planning is not addressed (but see Erlhagen & Schöner, 2002).

In the model, the virtual end-effector path can be computed from the virtual-joint path, λ_{i}(*t*), on the basis of the kinematic model of the effector. The virtual paths are straight line segments whatever the position in the work space. This demonstrates that the work space dependence of the curvature of the real end-effector path is not a consequence of movement planning in the model. Conversely, the virtual end-effector paths are not isomorphic with the real end-effector paths. In the model, the curved end-effector paths come from the relatively sluggish control at the muscle joint level that results when physiologically realistic parameter values for the muscle model are chosen. This is demonstrated in Figure 10B in which an increase of muscle impedance by a factor of 10 leads to unrealistically straight end-effector paths. We come back to this issue below.

The curved end-effector paths are thus signatures of imperfect control, that is, of a failure to realize the straight virtual end-effector path. Given sluggish muscle joint systems, two factors may contribute to this deviation from the plan: biomechanical coupling among the joints through interaction torques and muscular coupling among joints through multiarticular muscles. Manipulating these factors leads to characteristic changes also at the level of the joint trajectories, so we return to their role later in this section, around Figure 10. The upshot will be that perfect compensation of interaction torques and reduction of muscular coupling leads to straight end-effector paths invariantly within the work space and contrary to the experimental data.

How does redundancy play into the end-effector path? A given end-effector path can be realized by a multitude of task-equivalent joint configurations. The biomechanics of the arm, including the interaction torques, as well as the contributions of biarticular muscles, depend on joint angles, speed, and accelerations. The control problem is thus different for different arm configurations. Given the imperfections of the muscle joint control system, we may expect that such differences may lead to differences in the end-effector paths as well.

This hypothesis was tested in the experiments and in the model by imposing two different sets of initial arm configurations with identical end-effector position. In experiment, the starting joint configurations were chosen so that they could be reproducibly imposed and were comfortable for the participants, which limited the range of joint angle variation. Table 1 lists the mean joint angles for the starting configurations for the two pairs of movements considered. The model made use of these values as initial conditions for the joint angles. Note that the end-effector location was nearly identical for the two initial postures.

Movement Designation . | Joint 1 . | Joint 2 . | Joint 3 . | Joint 4 . |
---|---|---|---|---|

1 | −0.473 | 1.562 | 1.942 | 0.380 |

2 | −0.680 | 1.635 | 2.219 | -0.324 |

4 | −0.659 | 1.642 | 2.204 | -0.317 |

5 | −0.446 | 1.569 | 1.950 | 0.285 |

Movement Designation . | Joint 1 . | Joint 2 . | Joint 3 . | Joint 4 . |
---|---|---|---|---|

1 | −0.473 | 1.562 | 1.942 | 0.380 |

2 | −0.680 | 1.635 | 2.219 | -0.324 |

4 | −0.659 | 1.642 | 2.204 | -0.317 |

5 | −0.446 | 1.569 | 1.950 | 0.285 |

Notes: Each pair of movements had the same starting and target end-effector position (movements 1 and 2, respectively, movements 4 and 5). The joint angles are given in radians.

The bottom row of Figure 4 depicts the mean end-effector paths for the three subjects for the two movements with two initial configurations (1 and 2 are a pair, as are 4 and 5) as well as the associated simulations. In both cases, the end-effector paths are not qualitatively different for the two initial configurations, only small quantitative differences being observable (for movements 4 and 5, differences are in some cases essentially due to a slight shift of the end-effector starting position). This pattern of results is reproduced by the model at the reference parameter set. Thus, the model accounts for the observed motor equivalence of the two initial configurations. Apparently, although muscular control is sluggish in the model, the differences in biomechanical and muscular conditions that can be generated by reasonable variations of joint configuration are not sufficient to induce significant changes in end-effector control. Note that the model does not include compensatory mechanisms for interaction torques and muscle interjoint coupling, so in the model, motor equivalence must ultimately break down. We have verified that end-effector paths start to deviate when the initial joint configurations are made much more different from each other, beyond the range reachable in humans. Thus, the observed motor equivalence is not evidence for inverse dynamics or other forms of compensation for joint coupling. At realistic physiological conditions and within the range of realistic joint configuration variations, sluggish control does not lead to sizable effects on the end-effector path.

### 4.2. Joint and End-Effector Trajectories.

In the model, kinematic invariance of end-effector path under changes of effector configuration comes from the decoupling between joint velocity combinations that move the end-effector and joint velocity combinations that do not (see equation 2.9). A direct illustration of this principle is shown in Figure 5. In the left panel, the four joint trajectories are shown as time series together with the associated virtual joint trajectories, λ_{i}(*t*). Only when the real and virtual joint trajectory differ is muscular activity induced and torque generated. The real joint trajectory therefore lags behind the virtual trajectory. The amount of lag depends on biomechanics and muscle properties. This distance between real and virtual joint trajectory is largest for the most proximal joint 1, which accelerates against the largest inertial moment, while a smaller distance is sufficient to move the distal joint 4 that encounters very little inertia. Thus, the real joint trajectories are not simply time-shifted copies of the virtual joint trajectories.

In the middle panel of Figure 5, an additional acceleration is inserted into the null-space during the movement phase (, **s**(0) = 0 in equation 2.9). This leads to different virtual joint trajectories and consequently also to different real observed joint trajectories. Because the virtual joint trajectories are decoupled across the two subspaces of null- and range-space, the associated virtual end-effector trajectories are exactly identical in both conditions. This is also true for the real end-effector trajectories shown in the right panel of Figure 5. The overlaid traces from the two conditions are nearly indistinguishable. This is true even though dynamical conditions differ when motion in the null-space is induced. In fact, when unrealistically large amounts of self-motion are imposed, motor equivalence breaks down. Within a realistic range of kinematic conditions, however, motor equivalence prevails in the face of kinetic differences.

To examine the time courses of movements more closely, we must first address differences in movement time that occur in experiment. The left panel of Figure 6 shows the absolute value of end-effector velocity in time for different trials from a single participant. End-effector velocity profiles are smooth and bell shaped, invariantly across work space and participants (not shown). Trials differ by the total duration of the movement. The associated velocity profiles have the same shape but are rescaled in time. To compute mean trajectories, such differences in movement duration must be corrected by aligning and warping individual trajectories (see section 3). The right panel illustrates typical end-effector velocity profiles obtained from the model, which are also bell shaped and invariant across work space and match experimental velocity profiles.

The two-dimensional end-effector velocity trajectory reflects the curved end-effector paths. One way to look at that is to decompose the end-effector velocity into a component along the straight line from the initial to the target position and a component perpendicular to that direction. Figure 7 shows these two components of end-effector velocity for the three participants in the experiments and for the model at the reference parameter set for movements 1 and 4. The curved movement 1 has a sizable velocity component perpendicular to the straight path, while the relatively straight movement 4 does not. Both are captured by the model.

### 4.3. Self-Motion.

How do joint motions generate the end-effector trajectory? In a redundant system, joint velocity vectors can be decomposed into two components. Combinations of joint velocities that lie in the range space of the Jacobian move the end effector. Combinations of joint velocities that lie in the null-space of the Jacobian do not move the end effector. Arm motion that involves only such combinations of joint velocities that lie entirely within the null space of the Jacobian does not move the end effector at all. This amounts to a purely internal motion of the redundant effector, also called *self-motion*. The observation of self-motion provides evidence that the redundant degrees of freedom are in fact used during movement tasks that involve the end effector.

To quantitatively assess self-motion, we decompose the observed real joint trajectories into these two components (see section 3) and determine at every moment in time the lengths of these two joint velocity vectors representing the amounts of range-space motion and of self-motion. Because both subspaces are two-dimensional in the system studied here, there is no need to normalize these values to the dimensionalities of the two subspaces. Each row in Figure 8 shows time series of these two components of motion.

The first thing to notice is that there is a considerable amount of self-motion in all cases and for all subjects. The relationship between the amount of self-motion and the amount of range-space motion reaches beyond 50% in some cases and is typically above 30%. This discovery certainly excludes any notion that the total amount of joint velocity would be minimized such as in the pseudo-inverse solution (which is the minimum norm solution and thus generates zero self-motion; see, e.g., Cruse et al., 1993).

The model accounts for the pattern of self-motion observed. The small glitch in the last third of the trajectory comes from the deactivation of the oscillatory drive of the movement and could be smoothed out if a more gradual switching process was used.

One may be tempted to distinguish two fundamental causes of self-motion: self-motion may either be planned or may arise out of the imperfections of the muscle joint control system. The top two panels of Figure 9 contrast the amount of self-motion and range-space motion observed in the model at the level of the real joint trajectory (left) and virtual joint trajectory (right). Clearly, under these conditions, there is “planned” self-motion, although the total amount of real self-motion is larger and the temporal evolution is different.

The distinction between “planned” self-motion and self-motion arising from imperfect control is misleading, however, because planning may be part of the overall feedback loop (Todorov & Jordan, 2002). In our model, the back-coupling from the real joint trajectory into the dynamics of the virtual joint trajectory induces self-motion. This is illustrated by the second row of Figure 9, which comes from a simulation in which this back-coupling term was set to zero (all other parameters as in the reference parameter set). In this case, there is no self-motion at the level of the virtual trajectory at all (right panel), while significant self-motion remains at the level of the real joint trajectory (left panel). Thus, back-coupling contributes to self-motion, but back-coupling is not necessary for self-motion to arise. At model parameter settings that are physiologically plausible and fit a wide range of kinematic characteristics, there is a component of self-motion that is caused by neuromuscular control problems.

In the model, it is easy to explore those contributions by selectively varying model parameters that change the properties of the control system. We vary muscle impedance to establish the role of sluggish muscle joint systems, coupling among muscle joint systems by biarticular muscles and the presence of an inverse dynamics model to establish the role of interaction torques. At the reference parameter setting, the coupling among joints induced by multiarticular muscles is not a major contribution. Setting that coupling to zero barely changed the amount and time course of self-motion (not shown). Increasing the impedance of the muscle joint system by a factor of 10 (see the bottom panels of Figure 9) essentially eliminates all self-motion. Because the real trajectory then closely tracks the virtual trajectory, this also means that the virtual trajectory is self-motion free. In fact, the back-coupling term has no function in this limit case, because the very small deviation of the real from the virtual trajectory sends a very weak signal back to the virtual trajectory dynamics. Self-motion is thus caused in large part by the sluggish control exercised by the muscle-joint system.

### 4.4. Link Between Self-Motion and Curved End-Effector Paths.

We emphasized earlier that sluggish control leads to curved end-effector paths. The relationship between self-motion and curved end-effector paths is explored in Figure 10, which shows the end-effector paths associated with the two manipulations discussed above. Eliminating virtual self-motion by suppressing back-coupling does not affect end-effector paths much (see Figure 10A). Back-coupling is thus not necessary to obtain curved end-effector paths. Making muscle joint systems much less sluggish by increasing their impedance (see Figure 10B), in contrast, does straighten end-effector paths so much that they no longer match the experimentally observed pattern.

Given that we account for self-motion by postulating sluggish control at the muscle joint level, what are the mechanical problems this control system must solve, and how may they contribute to self-motion? To address this, we vary a number of factors of the biomechanical dynamics in Figure 11. The top left shows self- and range-space motions when both the Coriolis and centrifugal components of the interaction torques are eliminated. This is done simply by deleting the corresponding terms from the equation (something that is, of course, physically impossible to do in reality). The amount of self-motion is almost unchanged, which implies that self-motion is not caused primarily by these coupling terms.

To go further, we eliminate all interaction torques, including the nondiagonal elements of the inertial matrix (see equation D.2). This can be thought of as an emulation of a mode of control in which a neuronally computed inverse dynamics solution is used to generate the active torques that cancel the interaction torques. The top-right panel of Figure 11 shows that self-motion is strongly reduced, proving that interactions torques contribute substantially to the control problems that cause self-motion. The remaining level of self-motion is too low in comparison to experiment.

With all interaction torques removed, what is causing the remaining self-motion? The simulation shown in the bottom-left panel of Figure 11 illustrates that coupling among joints through multiarticular muscles plays a role. Such coupling leads to systematic deviations of the real from the virtual joint trajectory, some of which lead to self-motion. On the basis of the reference parameter set, setting all off-diagonal elements of the coupling matrix, **z**, to zero eliminates the contributions of multiarticular muscles, which reduces self-motion notably but not entirely. Muscular coupling among joints thus also contributes to self-motion. Eliminating both multiarticular muscles and canceling interaction torques finally suppresses self-motion almost completely (see the bottom-right panel of Figure 11). The small residual self-motion comes from nonlinearities of the muscle model.

These manipulations of the mechanical properties of the control system confirm the link we uncovered above between self-motion and curved end-effector paths. The bottom panels of Figure 10 display end-effector paths that are obtained for two of these manipulations. When an inverse dynamics is emulated (bottom left), paths become much straighter, to the point of no longer being realistic. The postulate of an inverse dynamics solution for multidegree of freedom control is thus incompatible with both the observed patterns of self-motion and the observed curvature of end-effector paths. Interestingly the curvature for movement 6 fits experiment better than that produced from the reference parameter set (see Figure 4). The opposite is true for movements 1 and 3.

The small amount of curvature that remains is due to the coupling among joints by multiarticular muscles. When that coupling is eliminated while at the same time providing an inverse dynamics cancellation of interaction torques (bottom right panel of Figure 10), end-effector paths become perfectly straight. This is the condition that also eliminates self-motion completely (see the bottom right panel of Figure 11). Overall, the rule seems to be that the less self-motion, the straighter the end-effector paths are. Note that this link between straightness of end-effector paths and the absence of self-motion is a property of the model, not a logical necessity: the very concept of self-motion means that it is does not itself affect the end-effector trajectory and thus also does not affect the end-effector path.

## 5. Discussion

We have analyzed a process model of the control of redundant multidegree-of-freedom arm movement that consists of five components: (1) a neural dynamics to initiate and terminate discrete goal-directed movements, (2) a neural oscillator that generates a timing signal that paces the progress of the end-effector along its path, (3) a neuronal dynamics of equilibrium points of all muscle joint systems of the redundant effector; (4) a muscle model, and (5) the mechanical dynamics of the redundant arm. Only the third component and its reciprocal coupling to the muscle joint system contained new assumptions. These were twofold. First, when locally decomposing the space of virtual joint velocities into the two subspaces in which the virtual end effector either moves (range-space) or does not move (null-space), we assumed that the dynamics within these subspaces are uncoupled. In other words, forces changing virtual joint velocities in the range-space do not influence virtual joint velocities in the null-space, and vice versa. Second, within the null-space, we assumed that the system receives input from the estimated real joint configuration in the form of a back-coupling, which is zero when virtual and real joint configurations and velocities are identical. We postulated that this core module of the model is largely housed in motor cortex but has strong links through mutual coupling into the structures involved in generating the time course of movement, including the cerebellum and the basal ganglia (Houk & Wise, 1995). Neural support for the back-coupling of the sensed effector state exists through that network as well as through reciprocal coupling to the spinal cord (Graziano, 2006).

We examined how this model accounts for the kinematic features of multidegree-of-freedom movement by comparing simulations of the model to results from a behavioral experiment in the modeled two-dimensional task setting with four degrees of freedom. Our strategy was to use parameter values for the muscle joint system that lie within typical physiological estimates of joint impedance and hence are relatively “sluggish.” In the light of criticism that pure servo-control systems without an internal model that anticipates joint torques are unable to generate the torques required to make movements at normal or fast rates (Gomi & Kawato, 1996; Ostry & Feldman, 2003; Kistemaker, Van Soest, & Bobbert, 2006), our assumptions put the theoretical framework of virtual joint trajectories (Feldman, 1966; Feldman & Levin, 1995; Won & Hogan, 1995) to a serious test.

We find that end-effector paths of the redundant experimental system are consistently curved, although the curvature varies with the location of the movement in work space. Our model accounts for curved end-effector paths through the combined effect of the dynamical properties of muscles and the dynamics of the virtual trajectory. Specifically, we find that unrealistically straight trajectories result from assuming higher muscle impedance than physiologically realistic and from assuming that the biomechanical interaction torques are exactly cancelled. Even if we ensure that the virtual end-effector path is straight, the sluggish muscle joint systems and the interaction torques lead to curved end-effector paths.

Our report of curved end-effector paths is in contrast to the classical approximation of nearly straight end-effector paths in nonredundant systems in a similar geometry (Morasso, 1981). The curvature of end-effector paths has been an intense topic of discussion over the years. It has been investigated relative to the question of whether movements are planned and controlled in terms of joint variables or in terms of end-effector variables (see the review in Barreca & Guenther, 2001). Most of these studies employed, however, nonredundant effector systems that, compared to redundant systems, face lesser control problems and thus have lesser potential for control errors to affect the end-effector path. Osu, Uno, Koike, & Kawato, (1997), for instance, showed that participants making two-dimensional movements using two joints could produce straighter end-effector paths if their movements were guided by a template path to do so. Monitoring the EMG of a number of involved muscles, Osu and colleagues excluded increased co-contraction as the control strategy that achieved straighter paths. They conclude that participants planned curved paths. Barreca and Guenther (2001) postulated such curved end-effector paths on other grounds as a way to minimize control effort and avoid joint limits. Note that in our formulation, the virtual joint trajectory is not a fixed plan but can be dynamically updated during the movement. Because back-coupling affects only the part of the virtual joint trajectory that does not affect the end effector, the virtual end-effector path that emerges from this dynamics is straight. Imperfect decoupling of the two subspaces of the virtual joint configuration space may conceivably induce curved virtual end-effector paths as a consequence of control errors.

Wolpert, Ghahramani, and Jordan (1994) proposed perceptual distortions as another potential source of end-effector path curvature by showing that perceptual estimates of straightness showed similar deviations from straightness as end-effector paths. Osu et al. (1997) showed, however, that participants spontaneously generate curved paths even where their perceptual distortions are minimal, such as in the fronto-parallel plane, and that participants can minimize curvature when instructed to do so. This rules out that perceptual distortion is the primary factor producing curved end-effector paths.

Our model predicts motor equivalence in the sense that differences in initial arm configuration do not necessarily lead to differences in end-effector paths. In principle, the active torques needed to produce the same end-effector trajectory differ when initial configuration varies. So one could have expected that a model like ours without inverse dynamics would not lead to equivalent end-effector trajectories. It simply turns out that quantitatively, such dependencies remain insignificant for the joint configurations probed in our experiments.

Finally, we have found consistently across participants and locations in work space that the joint velocity vectors contain a considerable amount of self-motion. Typically at least 30% of the joint velocities do not move the end effector! One might think of self-motion as a potential solution that the nervous system uses to avoid the well-known integrability problem of inverse kinematics in redundant systems (for review, see Mussa-Ivaldi & Hogan, 1991). When a closed path movement of the end-effector is repeatedly performed by a redundant arm using the Moore-Penrose pseudoinverse (which generates no self-motion), then the joint-configuration continues to drift from cycle to cycle, typically without bound. Self-motion makes it possible to compensate for this drift, leading to reproducible configurations at similar end-effector configurations. Relatedly, constraints within the space of degrees of freedom such as Donder's or Listings' law may also give rise to self-motion (Medendorp, Crawford, Henriques, Gisbergen, & Gielen, 2000). Although self-motion was not directly addressed in the latter study, it is plausible that self-motion would be needed in order to fulfill the observed constraints on multidegree-of-freedom movement trajectories.

The model accounts for the observed pattern of self-motion. Within our model, we were able to understand where self-motion comes from. One cause of self-motion is imperfect control of the low-impedance muscle joint system. Control errors are, a priori, as likely to induce self-motion as they are to induce end-effector relevant joint velocities. When we reduced control errors by increasing impedance, self-motion was strongly reduced. The back-coupling of estimated joint configurations and joint velocities into the neuronal dynamics of the virtual joint trajectory means that self-motion induced by control errors affects the virtual joint velocity. This is why in the model we find self-motion even at the level of the virtual joint velocity (although clamping such “planned” self-motion to zero does not eliminate self-motion of the real effector). We showed that an inverse dynamics model that cancels all interaction torques and reduces the amount of control error leads to too little self-motion compared to experiment as well as to end-effector paths that are unrealistically straight. Another contribution to self-motion are multiarticular muscles that couple torque generation at different joints. Eliminating such muscular coupling from the model also reduces self-motion and makes end-effector paths straighter.

The problem of redundant motor control has been addressed in other theoretical and modeling efforts. We share the framework of neuronal dynamics with a series of models (Bullock & Grossberg, 1988; Bullock, Grossberg, & Guenther, 1993, 1996; Guenther, 1994; Guenther & Barreca, 1997). These models learn a one-to-one mapping from task-level variables to the associated joint velocities, in effect selecting a specific solution to the inverse kinematics problem. Feedback at the level of the task variables enables replanning and accounts for how movement goals are reached even if degrees of freedom are blocked. Although a variant of these models solves the nonintegrability problem mentioned above and thus must be producing self-motion (see section 5 in Guenther & Barreca, 1997), no systematic account for self-motion is provided. These models do not address the biomechanical dynamics of redundant effectors and thus do not account for the role played by interaction torques and biarticular muscles in generating the control errors that lead to self-motion and to nonstraight end-effector paths. More recent generalizations (Bullock, Cisek, & Grossberg, 1998; Cisek, Grossberg, & Bullock, 1998) of the VITE class of models (Bullock & Grossberg, 1988) address the neuromuscular level in considerable detail but have not yet been applied to redundant effector systems.

Principles of stochastic optimal control have been used to provide a very different kind of account for redundancy in which control signals are computed by optimizing an effort functional at each point in time during the movement (Todorov & Jordan, 2002, 2003; Todorov, 2004). While the temporal and spatial continuity of the computation is neurophysiologically plausible, the processes through which neuronal networks perform the complex optimality computations in real time remain to be explored (but see Scott, 2004; Shadmehr & Krakauer, 2008). These models have not taken into account real arm kinematics, probably due to the difficulty of solving the optimality conditions for nonlinear geometries, nor have they linked to muscle joint models and biomechanical dynamics (but see Guigon, Baraduc, & Desmurget, 2007, for a first effort in this direction).

One of the limitations of our model is the restriction to two-dimensional end-effector motion and only four joint angles. This is primarily a practical limitation dictated by the desire to make detailed comparisons between theory and experiment, more difficult to achieve when the number of degrees of freedom and of task dimensions increases. The empirical question may arise, however, if there is anything specific about motion constrained to a plane that promotes the signatures of redundancy on which we focused: self-motion and curved end-effector paths. To demonstrate that this is not the case, we provide sample data from a separate experimental study in which participants reached in three spatial dimensions to point to targets using the full set of 10 degrees of freedom from scapular motion to the wrist. Figure 12 shows the end-effector paths in three dimensions for three participants as well as the associated joint velocities in the range-space and the null-space. In this system, the range-space defined in reference to the three-dimensional Cartesian end-effector position is three-dimensional, while the null-space has seven dimensions. Thus, the amount of joint motion in these two subspaces is normalized by dividing the sum of the squared joint velocities by the number of dimensions of each subspace. This down-weights the amount of self-motion. Even so, there is substantial self-motion of typically about a third of the range-space motion, quite similar to the two-dimensional data discussed earlier. The end-effector paths are curved, also quite similar to what we reported in two dimensions. Thus, empirically, motion in three dimensions involving many more degrees of freedom reveals the same signature of redundancy as analyzed and modeled in this letter. We have generalized our model to this problem and found the model compatible with these results. Because this modeling must address a number of new issues, including the effect of gravitational torques, a report on this generalization exceeds the scope of this letter.

At the level of the muscle joint system, compared to much more elaborate muscle models (Chang, Brown, & Loeb, 2000), we have made strongly simplifying assumptions motivated by the practical goal of keeping the model simple and limiting the number of model parameters. The strongest simplifying assumption was neglecting the modulation of stiffness through a modulation of co-contraction during movement. This is a limitation that must be overcome when generalizing to mechanically more challenging conditions such as lifting loads or moving rapidly. We have not included mechanisms that would impose joint limits. This did not become relevant in the comparison to data, as the real movements stayed far from joint limits. The redundancy problem at the muscular level is interesting in its own right, and it seems possible to use ideas similar to the ones discussed here at that level (see Laboissière, Ostry, & Feldman, 1996, for a discussion that has conceptual similarity with the uncontrolled manifold). Addressing this level is, however, beyond the scope of our contribution.

A final major limitation of the model reported here is the absence of an account for variability (see Goodman & Latash, 2006, for a first effort in this direction). We have analyzed the impact of various sources of noise on the movement generated by our model (Martin, 2005). Given the amount of detailed analysis and comparison to experiment this requires, we have decided to discuss the stochastic properties of our theoretical account in a forthcoming article.

At this point, what general conclusions can we draw from our analysis? One insight is certainly that the redundant degrees-of-freedom in a given task are used. The considerable amount of self-motion we observed shows that it is not a priority of the nervous system to select solutions to the degrees-of-freedom problem that are unique to a given task. In the model, this self-motion is not planned but emerges from control problems such as relatively low impedance, interaction torques, and multi-articular muscles. More generally, however, the degrees of freedom that are redundant with respect to one task may very well be used to accommodate another task at the same time. The flexibility offered by allowing self-motion may be the general principle of how the central nervous system accommodates the complex demands on human movement that arise in the real world. This flexibility may be the most important feature of human movement, more important than good control. We have found consistently that a picture with relatively sluggish muscle joint systems and no specific mechanisms to compensate for interaction torques or for muscular interjoint coupling provides an adequate description of multidegree-of-freedom movement.

## Appendix A: Kinematic Model

*l*the length of arm segment

_{i}*i*(

*i*= 1… 4) numbered from proximal to distal segments, and θ

_{i}are the joint angles (

*i*= 1… 4 from the proximal (sterno-clavicular joint) to the distal joints (wrist)), computed in each case relative to the next proximal segment (see Figure 2). The coordinates of the end-effector position, (

*x*,

*y*), are in a Cartesian coordinate system centered at the sterno-clavicular.

## Appendix B: Kinetic Model

The biometric parameters entering the kinematics and kinetics were set based on measurements from one participant and are listed in Table 2. The center of mass and the inertia of the various arm segments are computed following Hanavan (1964). Data for the scapular joint are not available and are estimated to be a quarter of the upper torso biometrics data.

Parameter Name . | Symbol . | Value . | Units . | . |
---|---|---|---|---|

Body mass | M | 55 | kg | |

First segment length | l_{1} | 0.2024 | m | |

Second segment length | l_{2} | 0.3035 | m | |

Third segment length | l_{3} | 0.2586 | m | |

Fourth segment length | l_{4} | 0.1658 | m |

Parameter Name . | Symbol . | Value . | Units . | . |
---|---|---|---|---|

Body mass | M | 55 | kg | |

First segment length | l_{1} | 0.2024 | m | |

Second segment length | l_{2} | 0.3035 | m | |

Third segment length | l_{3} | 0.2586 | m | |

Fourth segment length | l_{4} | 0.1658 | m |

## Appendix C: Muscle Model

_{i}, by a constant amount of co-contraction,

*Co*. The combined torque,

*T*, generated at joint

_{i}*i*by the agonist and antagonist components, is at equilibrium when θ

_{i}= λ

_{i}(in the absence of external torques). Here,

*K*

_{l}is a linear and

*K*

_{nl}is a nonlinear stiffness factor. Two types of viscosity are taken into account: a linear contribution to viscosity (coefficient ) captures the physical properties of muscle tissue, and a nonlinear contribution (coefficient, ) to viscosity summarizes more complex velocity dependent contributions to peripheral control.

Parameter Name . | Symbol . | Value . | Units . |
---|---|---|---|

Co-contraction | Co | π/90 | rad |

Impedance matrix | Z | ||

Linear stiffness | K _{l} | .4 | kg m^{2}s^{−2} |

Nonlinear stiffness gain | K _{nl} | 1 | |

Linear active viscosity | μ_{bl} | .3 | kg m^{2}s^{−1} |

Linear passive viscosity | μ_{rl} | 0.03 | kg m^{2}s^{−1} |

Parameter Name . | Symbol . | Value . | Units . |
---|---|---|---|

Co-contraction | Co | π/90 | rad |

Impedance matrix | Z | ||

Linear stiffness | K _{l} | .4 | kg m^{2}s^{−2} |

Nonlinear stiffness gain | K _{nl} | 1 | |

Linear active viscosity | μ_{bl} | .3 | kg m^{2}s^{−1} |

Linear passive viscosity | μ_{rl} | 0.03 | kg m^{2}s^{−1} |

## Appendix D: Other Model Details

### D.1. Termination Phase of the Movement.

The end of the movement is signaled to the neuronal dynamics, equation 2.1, by a function that probes how close the timing variable is to its resting state: . This function is clamped to zero when the movement is initiated. The beginning of the movement is brought about by setting *I _{m}* to a large positive value for a brief moment, which helps switch the neuronal dynamics in the movement state.

During the resting phase of the system, the virtual joint motion is actively damped (by adding to equation 2.11). Any remaining discrepancies between the desired end-effector location, r_{d} and the estimated end-effector location, **r**_{e} are servo-controlled to zero.

### D.2. State Estimation.

### D.3. Emulation of Inverse Dynamics.

## Appendix E: Parameter Values

There are four classes of parameters for which values must be determined. (1) Biometric parameters are directly estimated from measurements on one participant and listed in Table 2. (2) Parameters of the muscle joint model were determined based on estimates from the empirical literature as discussed above. These parameters are listed in Table 3. (3) A set of model parameters describe the mechanism for initiating and terminating the movement as well as determining the temporal shape of the movement. The values of these parameters have very little impact on the results and are listed in Table 4. (4) A small set of model parameters describes the core neuronal dynamics of virtual joint trajectory formation and are relevant for the model. These four parameters are listed in Table 5.

Parameter Name . | Symbol . | Value . | Units . |
---|---|---|---|

Oscillator limit cycle time | ω_{h} | 16 | Hz |

Oscillator stability | α_{h} | 0.05 | s^{−1} |

Oscillator relaxation constant | β_{o} | 15 | s^{−1} |

Switching dynamic constant | β r, β m | 400 | s^{−1} |

Switching dynamic resting state | h | 1 | |

Switching dynamic constant | Δ | 2 | |

Stopping activity gain (I) _{r} | a | 15 | |

Gain of sigmoid for initiation dynamics | a _{r} | 100 | |

Fine positioning constant | β_{f} | 250 | s^{−1} |

Virtual velocity damping constant | β_{s} | 30 | s^{−1} |

Parameter Name . | Symbol . | Value . | Units . |
---|---|---|---|

Oscillator limit cycle time | ω_{h} | 16 | Hz |

Oscillator stability | α_{h} | 0.05 | s^{−1} |

Oscillator relaxation constant | β_{o} | 15 | s^{−1} |

Switching dynamic constant | β r, β m | 400 | s^{−1} |

Switching dynamic resting state | h | 1 | |

Switching dynamic constant | Δ | 2 | |

Stopping activity gain (I) _{r} | a | 15 | |

Gain of sigmoid for initiation dynamics | a _{r} | 100 | |

Fine positioning constant | β_{f} | 250 | s^{−1} |

Virtual velocity damping constant | β_{s} | 30 | s^{−1} |

Parameter Name . | Symbol . | Value . | Units . |
---|---|---|---|

Virtual relaxation constant | β_{v} | 30 | s^{−1} |

Back-coupling constant position | β_{s1} | 100 | s^{−2} |

Back-coupling constant velocity | β_{s2} | 35 | s^{−1} |

Delay constant | τ | 0.01 | s |

Parameter Name . | Symbol . | Value . | Units . |
---|---|---|---|

Virtual relaxation constant | β_{v} | 30 | s^{−1} |

Back-coupling constant position | β_{s1} | 100 | s^{−2} |

Back-coupling constant velocity | β_{s2} | 35 | s^{−1} |

Delay constant | τ | 0.01 | s |

## Acknowledgments

This work was supported by grant number NS050880 from the National Institute of Neurological Disorders and Stroke. The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institute of Neurological Disorders and Stroke or the National Institutes of Health.

## References

*Science*, 1994,

*263*, 1295–1297