## Abstract

Smooth, frictionless, kinematic constraints on the motion of a grasped object reduce the motion freedoms at the hand, but add force freedoms, that is, force directions that do not affect the motion of the object. We are studying how subjects make use of these force freedoms in static and dynamic manipulation tasks. In this study, subjects were asked to use their right hand to hold stationary a manipulandum being pulled with constant force along a low-friction linear rail. To accomplish this task, subjects had to apply an equal and opposite force along the rail, but subjects were free to apply a force against the constraint, orthogonal to the pulling force. Although constraint forces increase the magnitude of the total force vector at the hand and have no effect on the task, we found that subjects applied significant constraint forces in a consistent manner dependent on the arm and constraint configurations. We show that these results can be interpreted in terms of an objective function describing how subjects choose a particular hand force from an infinite set of hand forces that accomplish the task. Without assuming any particular form for the objective function, the data show that its level sets are convex and scale invariant (i.e., the level set shapes are independent of the hand-force magnitude). We derive the level sets, or “isocost” contours, of subjects' objective functions directly from the experimental data.

## INTRODUCTION

Consider the task of using a single hand to carry a rigid object from one configuration (position and orientation) to another. The configuration of a rigid object in space has 6 degrees of freedom: 3 translational freedoms and 3 rotational freedoms. Now suppose that the object's motion between the start and goal configuration is subject to smooth frictionless configuration constraints. Such constraints arise, for example, when the object is confined to a linear rail or attached to a low-degree-of-freedom linkage. These constraints limit the object to *m* <6 degrees of freedom, but allow the subject to apply forces against the constraints in the 6 − *m* directions orthogonal to the motion freedoms. To accomplish the constrained carrying task, the subject's motor control system must simultaneously resolve a number of redundancies, including *1*) determining a trajectory of the arm consistent with the object's *m*-dimensional configuration space (where the arm will generally have more than *m* degrees of freedom), *2*) determining the constraint forces applied during the motion, and *3*) determining how joint torques, implied by the trajectory and constraint forces, are distributed across muscle groups.

Although a great deal of work has studied how motion freedoms [redundancy *1*) above] are resolved in unconstrained point-to-point arm motions (e.g., Alexander 1997; Flash and Hogan 1985; Harris and Wolpert 1998; Kuo 1994; Todorov and Jordan 1998; Uno et al. 1989) and how muscle load-sharing [redundancy *3*) above] is resolved in completely constrained single-arm isometric tasks (e.g., Buchanan et al. 1986; Flanders and Soechting 1990; Gomi 2000; van Bolhuis and Gielen 1999), there has been less work on understanding how constraint force freedoms [redundancy *2*) above] are resolved in partially constrained tasks. Because applying forces against constraints has no effect on the task, the manner in which force freedoms are resolved provides powerful clues to the organization of the motor control system. For example, although many hypotheses have been proposed to explain experimental unconstrained arm motion data (minimum Cartesian jerk, minimum rate of torque change, minimum metabolic cost, etc.), these hypotheses' differing predictions of constraint forces can be used to determine whether they are applicable also to the case of constrained manipulation. Ideally an organizing principle would be able to explain both unconstrained and constrained arm motions.

Consider, for example, the following experiment. A subject consistently chooses a trajectory (and associated path) solving a particular point-to-point unconstrained reaching task. Now we place a frictionless guide rail exactly along that chosen path and ask the subject to again perform the reaching task several times. Optimization models of the dynamics may predict that the subject will learn to use the new force freedom by applying forces against the rail to further decrease the objective function. In other words, the presence of the rail allows the subject to “optimize” further (Lynch et al. 2000). This is exactly the effect we are interested in: how subjects naturally take advantage of force freedoms.

There are reports of previous studies on natural interaction with smooth constraints. Perhaps most related to the present paper are studies on how subjects turn a crank (Russell and Hogan 1989; Svinin et al. 2003). Russell and Hogan showed that subjects apply significant radial forces (compressing or extending the crank) even though they are workless. Svinin et al. argued that experimental forces and motions during crank-turning can be described by minimization of a weighted combination of changes in hand force and joint torques. Van der Helm and Veeger (1996) studied the related problem of shoulder muscle activation during wheelchair propulsion, and both their shoulder mechanism model (van der Helm 1994) and experimental results showed that subjects apply forces normal to the rim of the wheel. The efficiency of the motion is also discussed. Gomi (1998) showed that the natural stiffness at the hand during motion is altered in the presence of a guiding constraint. Scheidt et al. (2000) studied persistence of motor adaptation during constrained multijoint arm movements. Their decomposition into kinematic and dynamic criteria influencing disadaptation correspond roughly to our decomposition into trajectory and force freedoms.

This paper reports the results of the simplest possible experiment studying how subjects resolve a constraint force redundancy. Each subject was asked to use the right hand to hold stationary a handle being pulled with constant force along a low-friction linear rail (*m* = 1 motion freedom). The arm was supported in a horizontal plane with the wrist cuffed, so that the arm could be treated as a two-joint shoulder–elbow mechanism. To hold the handle stationary, the subject had to apply an equal and opposite force along the rail, but subjects were also free to apply a force against the constraint, orthogonal to the direction of the pulling force (one force freedom in the horizontal plane because the arm was treated as a two-joint mechanism). Despite the fact that constraint forces increase the magnitude of the total force vector at the hand and have no effect on the task, we found that subjects applied significant constraint forces in a consistent manner dependent on the arm and constraint configuration. We show that the constraint forces can be interpreted in terms of an objective function describing how subjects choose a particular hand force from a family of hand forces that accomplish the task. Without assuming any particular form for the objective function, the data show that its level sets are convex and scale invariant (i.e., the level set shapes are independent of the hand-force magnitude). We derive the level sets, or “isocost” contours, of subjects' objective functions directly from the experimental data. In other words, in contrast to previous work on optimization models that use experimental data to support or invalidate candidate objective functions based on a biomechanical model, we use a new method for directly measuring the level sets of the objective function without assuming any particular form for it. These level sets may be thought of as nonparametric objective functions that act as descriptors and predictors of behavior, independent of any interpretation in terms of biomechanics and neural control. Importantly, the objective functions appear to be independent of the arm configuration when expressed as objective functions on joint torques. We have compared our results to the predictions of several biomechanical models of force generation, and although these results are inconclusive because of uncertainty in subjects' physiological parameters, models based on the sum of muscle tensions and stresses can be effectively ruled out. An objective function that is a simple positive-definite quadratic form on the joint torques appears to fit the data well. If the muscles are springlike (e.g., Mussa-Ivaldi et al. 1985), one interpretation of this objective function is that subjects choose a hand force that satisfies the task while minimizing the potential energy stored in the muscles. Portions of this work have previously appeared in conference form (Pan et al. 2004; Tickel et al. 2002).

## METHODS

### Setup and protocol

Fourteen healthy right-handed male subjects (Table 1) were seated in a custom-made high-backed chair with an adjustable seat, to raise or lower the height of the shoulder plane based on the height of the subject. To fix the shoulder location, subjects were restrained by a 4-point harness. The wrist was immobilized by an over-the-counter commercially available wrist cuff, and the subject grasped a vertical handle on a slider on a horizontal low-friction linear rail. The rail is mounted on a lazy Susan turntable, allowing the rail to be rotated 360° in the plane. The orientation of the rail can be fixed at any angle by clamping the turntable. The handle can spin freely about a vertical axis so that no torques at the hand are involved, and a support plate is attached to the handle to support the weight of the forearm (Fig. 1*A*). This support maintains the arm in a horizontal plane throughout experiments without fatiguing the shoulder.

A 6-axis force-torque sensor (ATI-AI Gamma 15–50) is positioned between the handle and the slider and is used to measure forces against the rail. A cable attached to the slider passes through a pulley system, allowing weights to be suspended under the slider to create a tangential force along the rail.

Each trial consisted of the subject holding the handle while a weight was hung from the cable, causing a tangential pulling force on the slider. The subject then stabilized the position of the handle at the center of the rail. Forces normal to the rail were then recorded for 1 s and averaged. The weight was then removed from the cable.

For 8 subjects (subjects 1 through 8), the slider handle was located at (0, 45 cm) in a frame fixed to the shoulder, as shown in Fig. 1*B.* Sixteen angles of the rail were used, evenly spaced at 22.5° intervals. At each of the 16 test angles, 2 different weights were hung from the cable, 0.858 kg (the light weight) and 1.759 kg (the heavy weight). These resulted in tangential forces of 8.4 and 17.3 N, respectively. For each angle and weight, the experiment was repeated 3 times. Therefore, for each handle position, we collected 16 × 2 × 3 = 96 data points. The ordering of the trials was randomized to minimize any history effect in the results. Subjects were told to suspend the weight as naturally and comfortably as possible, and not to excessively cocontract to stiffen the position of the handle. Fatigue was minimized by the short durations of each experiment, and subjects were permitted to take a break at any time. The total testing time for each subject was about 1 h.

For 6 subjects (subjects 9 through 14), the experimental protocol was similar, except only a single weight of 1.2 kg was used, giving 11.8 N tangential force. Each of these 6 subjects was tested at 5 different positions of the hand, as shown in Fig. 9A.

The protocol was approved by the Northwestern University Institutional Review Board.

### Objective functions and isocost contours

The organizing principle governing how subjects apply constraint forces can be expressed as an objective function describing the “cost” of generating a particular force vector at the hand. This objective function may reflect the “effort” involved in applying a particular force, or it may simply reflect the organization of the motor control system. In either case, the role of the objective function is simply to resolve the freedom in the applied constraint force.

An objective function *g* may be viewed as a mapping from the arm configuration and the force and torque applied at the hand to a nonnegative real number representing the “cost.” For ease of discussion, we will assume that the arm configuration is implied and that the hand force can be written as a 2-vector **f**_{h} = (*f _{hx}, f_{hy}*) in a Cartesian frame aligned with the shoulder frame (Fig. 1B). Then the objective function

*g*(

**f**

_{h}) can be viewed as a function

*g*: ℜ

^{2}→ ℜ. (Any objective function on joint torques or muscle tensions uniquely defines an objective function on the hand forces.) Thus the objective function forms a 2-dimensional surface. A weak assumption on the form of

*g*is that the cost increases monotonically as we move outward along any ray from the origin

**f**

_{h}= (0, 0). For low to moderate hand forces, this seems intuitively correct: when a subject can solve a task with any hand force α

**f**

_{h}for α ≥ 1, the subject is likely to choose approximately the minimum necessary force, α ≈ 1. An objective function satisfying this condition defines a bowl in the hand-force space, and there are no local minima except the global minimum at (0, 0). The level sets (curves) of a bowl-shaped cost function are concentric, closed, and star-shaped (can be written as a function of the polar angle) about the origin. We call these level curves

*isocost contours*in the hand-force space.

Two important subclasses of bowl-shaped objective functions are those whose isocost contours are all convex, and those whose isocost contours are *scale invariant*—each isocost contour is a uniformly scaled version of every other. These properties are shown graphically in Fig. 2. An objective function *g*(**f**_{h}) is convex if its Hessian matrix is positive definite at all *f*_{h}, i.e. (1) (2) An objective function is scale invariant if it satisfies the property for some scaling function *k*(α) satisfying *k*(1) = 1 and monotonically increasing with α. A common example of *k*(α) is a power law α^{p} (*p* > 0). If an objective function is both convex and scale invariant, we refer to it as a CSI objective function. These properties of an objective function generalize immediately to hand forces and torques in more than 2 dimensions.

To use an objective function to predict constraint forces, note that the subject must apply a specific tangential force **f**_{t} along the rail to prevent the slider from moving, but is free to apply any normal force **f**_{n} in the orthogonal constrained direction^{1} This means that the space of hand forces that solves the task is the one-dimensional linear subspace *L* (a line) of the 2-dimensional hand-force space defined as where **f**_{n} is a nonzero force vector in the constrained direction. The objective function predicts that the subject will choose the force in *L* that minimizes the cost. At this point, *L* is tangent to one of the isocost contours. (If the cost function is nonconvex, a tangent point is not necessarily a minimum.) This can be seen graphically in Fig. 3. On the plot of isocost contours, construct the line *L* passing through the point **f**_{t} = (*f _{tx}*,

*f*) and perpendicular to the direction of the rail. This is the linear subspace of forces satisfying the task. The optimal total force

_{ty}**f**

_{h}=

**f**

_{t}+

**f**

_{n}occurs where the line is tangent to an isocost contour.

As shown in Fig. 3, a strange situation can occur if the isocost contours are not convex: multiple optimal hand forces **f**_{h} may be predicted. Equivalently, as the angle of the tangential force **f**_{t} passes smoothly through a critical angle where a nonuniqueness occurs, the optimal applied constraint force **f**_{n} changes discontinuously. An example of a nonconvex objective function is g(**f**_{h}) = |*f _{hx}*|

^{1/2}+ |

*f*|

_{hy}^{1/2}. Models with exponents <1 lead to nonconvexity and apparently have little physical basis.

If the objective function is scale invariant, the direction of the optimal total force **f**_{h} depends only on the direction of **f**_{t}, not on its magnitude. In other words, if the required tangential force **f**_{t} is scaled by α, then the optimal normal force is also scaled by α. This does not hold for more general objective functions.

### Reconstructing isocost contours

The measured force **f**_{h} applied by a subject in a given trial is decomposed into 2 orthogonal components: a component tangential to the rail denoted **f**_{t} with value equal and opposite to the pulling force, and a component perpendicular to the rail, denoted **f**_{n}, which is the object of study. For each rail orientation we define a reference frame such that **f**_{h} = (*f _{n}*,

*f*), where

_{t}*f*and

_{n}*f*are the scalar values of the force against and along the rail, respectively. The +

_{t}*f*-axis is oriented along the rail in the direction the subject must apply a force to prevent motion of the slider, and the +

_{t}*f*-axis is 90° clockwise, as shown in Fig. 1

_{n}*B*.

The results for a single subject at a single slider position can be plotted as shown in Fig. 6 in the results section. There are 2 plots, one corresponding to the light weight and one to the heavy weight. Each plot shows the normal force applied by the subject as a function of the angle of the tangential force **f**_{t} along the rail.

From each plot we can reconstruct an isocost contour of the subject's objective function. That this is possible may not be obvious because the experiments do not keep the subject on the same isocost contour. In fact, there is no way to design the experiments to do so without knowing the objective function in advance. If the objective functions are CSI, however, then it is possible to extract the isocost contours from the data, as described below.

Each point on a normal force plot, as in Fig. 6, indicates a point in the hand force (*f _{hx}*,

*f*) space, at an angle β relative to the +

_{hy}*f*-axis (Fig. 4). At this point, the direction of the normal force

_{hx}*f*is tangent to the isocost contour, as shown in Fig. 3. Therefore, the

_{n}*p*data points of the normal force plot give us a set of angles β

_{i},

*i*= l

*p,*and a tangent direction γ

_{i}associated with each β

_{i}. Choosing a point at an arbitrary radius

*r*

_{1}(say

*r*

_{1}= 1) along a ray at angle β

_{1}from the origin of the (

*f*,

_{hx}*f*) space, follow the tangent angle γ

_{hy}_{1}(i.e., integrate) until β

_{2}is reached. Then using angle γ

_{2}, integrate until β

_{3}is reached, and so on. (More sophisticated interpolating numerical integration could instead be used.) Continue around angularly until the curve reaches β

_{1}again. If the normal force plot comes from a scale-invariant objective function, the curve will close at β

_{1}. The key point is that for scale-invariant objective functions, the tangent direction γ depends only on the angle β of the force

**f**, not the magnitude ‖

_{h}**f**

_{h}‖, and therefore the data do not have to be derived from the same isocost contour to be able to reconstruct an isocost contour.

The procedure outlined above will result in a closed curve only if the normal force data is zero mean—the integral of the normal force curve must be zero (see appendix). All scale-invariant bowl-shaped objective functions imply a normal force plot with zero mean. As we see in the results section, the data are approximately zero mean and support the CSI hypothesis, making the reconstruction possible.

### Transforming to joint space

For the two-joint arm of Fig. 1*B*, isocost contours in the hand-force space are transformed to isocost contours in the joint-torque space by the relation (3) where **f**_{h} = (*f _{hx}*,

*f*)

_{hy}^{T}is a hand-force vector, τ = (τ

_{1}, τ

_{2})

^{T}is a vector of shoulder and elbow torques, and the arm Jacobian

*J*(θ) is

### Ellipse fitting of isocost contours

In the results section we see that isocost contours in joint-torque space appear rather elliptical, so we can fit ellipses to the data. We use the general form *Ax*^{2} + 2*Bxy* + *Cy*^{2} + 2*Dx* + 2*Ey* = 1, where *A*, *B*, and *C* describe the shape of the ellipse, and *D* and *E* describe its offset from the origin. The coefficients *A*, *B*, *C*, *D*, and *E* can be found by a least-squares fit minimizing ∑ (1 − *Ax*^{2} − 2*Bxy* − *Cy*^{2} − 2*Dx* − 2*Ey*)^{2} over the data points. Defining features of the ellipse are its orientation, given by tan^{−1} 2*B*/(*A* − *C*), and its eccentricity, given by where *l _{a}* and

*l*are the half-lengths of the long and short principal axes

_{b}### Biomechanical modeling

Joint torques are caused by a complex set of uniarticular and biarticular muscles crossing both the shoulder and the elbow (An et al. 1981; Meek et al. 1990; Murray et al. 2000; Pigeon et al. 1996; van der Helm 1994; Wood et al. 1989). The torque generated by each muscle is a function of the muscle tension stemming from muscle activation and the joint-angle–dependent moment arms based on the bone attachment points. The maximum tension available from a muscle is roughly a function of the physiological cross-sectional area (*PCSA*) and muscle stretch (and, in nonisometric settings, the rate of lengthening or shortening).

To simplify the model, we follow van Bolhuis and Gielen (1999) and Gomi (2000) and combine the muscles into 6 muscle groups: shoulder extensor and flexor, elbow extensor and flexor, and biarticular extensor and flexor. We define the muscle tension vector φ = (φ_{se}, φ_{sf}, φ_{ee}, φ_{ef}, φ_{be}, φ_{bf})^{T} ∈ ℜ^{6} to capture the tension of each of these groups of muscles. All elements of the vector must be nonnegative, indicating that each muscle group is capable of pulling only. This simplification into muscle groups makes the assumption that all muscles in each group are activated proportionally (van Bolhuis and Gielen 1999). With this model, the joint torques τ are obtained from the muscle tensions φ by (4) where *A*(θ) ∈ ℜ^{2×6} is a matrix of joint-angle–dependent moment arms.

Figure 5 shows a model of the arm with these 6 muscle groups (adapted from Gomi 2000). By *Eq. 3*, torque arising from shoulder monoarticular muscles causes hand forces along the line of the forearm, torque arising from elbow monoarticular muscles causes hand forces along the line through the shoulder, and biarticular muscles with τ_{1} = τ_{2} generate hand forces parallel to the upper arm.

By combining *Eq. 3* and *4*, we get (5) Because **f**_{h} ∈ ℜ^{2} and φ ∈ ℜ^{6}, there is an infinite set of muscle tension vectors φ that generate a specified **f**_{h}. Thus, in addition to the freedom to provide any force normal to the rail, the subject has freedom in how to share the load across muscle groups.

We consider the following minimization models as candidates for interpreting experimental isocost contours. Some of the muscle load-sharing models were considered for isometric force generation by Gomi (2000) and van Bolhuis and Gielen (1999).

HAND Hand-force magnitude ‖

**f**_{h}‖ is minimized. According to this model, the subject applies only forces tangent to the rail. The constraint force is zero.T2 Torque squared, ∑

_{i}τ_{i}^{2}. For a stationary robot arm with identical DC motors at the shoulder and elbow, this solution minimizes the electrical power to the motors.WT2 A quadratic form of the joint torques of the form τ

^{T}*W*τ, where*W*is a positive-definite 2 × 2 symmetric matrix (only 3 unique elements). Isocost contours are ellipses in the joint-torque space, a generalization of the T2 model, where*W*is the identity matrix and the isocost contour is a circle in joint-torque space.MT

*k*Sum of muscle tensions raised to the power of*k*= 1, 2, or 3, ∑_{i}φ_{i}^{k},*i*∈ {*se*,*sf*,*ee*,*ef*,*be*,*bf*}. The model*k*= 1 was proposed by Yeo (1976), and Nelson (1983) and Hogan (1984) suggest that metabolic power consumed by a muscle is proportional to the square of muscle force*k*= 2.MS

*k*Sum of muscle stresses raised to the power of*k*= 1, 2, or 3, ∑_{i}(φ_{i}/*PCSA*)_{i}^{k}where*PCSA*is the physiological cross-sectional area of muscle_{i}*i*(Crowninshield and Brand 1981). This is a measure of the activation of the muscle. There is some evidence that muscle endurance time is inversely proportional to (φ/*PCSA*)^{3}(Prilutsky et al. 1998).

Calculating the predictions of the MT*k* and MS*k* models requires the physiological cross-sectional area *PCSA* and moment arms for each muscle group. Table 2 gives examples of parameters we used, following (Gomi 2000). These values are lumped parameters obtained from data found in Meek et al. (1990). In these parameters, the moment-arm matrix *A* is independent of the joint angles. A definitive test of optimization models would, of course, require a method for obtaining the parameters of Table 2 for each subject.

Each of the 9 models defines a CSI objective function in the hand-force space for a given arm configuration. For the T2 and WT2 models, the isocost contours are ellipses that can be found in closed form. The isocost contours for the HAND model are simply circles centered at the origin. Isocost contours for the other models can be found numerically. The linear models MT1 and MS1 result in convex polygonal isocost contours; the other models have strictly convex isocost contours.

The linear models MT1 and MS1 predict activation of only one of the muscle groups for a given task, whereas higher-order models predict greater sharing of the load across the muscle groups. Each of the models predicts a normal force plot that can be compared to experimental data.

To obtain the prediction for model T2, let α**f**_{n} represent the force applied against the constraint, where **f**_{n} is a unit vector normal to the tangential force **f**_{t} applied by the subject to resist motion. Then α satisfies the equation The prediction for WT2 can be obtained similarly.

To obtain the predictions of models MT*k* and MS*k, k* = 1, 2, 3, we solve for the tension vector φ minimizing the objective function, subject to φ ≥ 0 (all muscles pulling) and requiring that the tangential force be equal to **f**_{t}. This problem is a linear programming problem for *k* = 1 and a nonlinear optimization for *k* = 2, 3. All optimizations were solved using CFSQP (Lawrence et al. 1994), C code implementing sequential quadratic programming.

## RESULTS

The experimental results for subject 1 at the hand position (0, 45 cm) are shown in Fig. 6 as plots of the applied constraint force *f _{n}* as a function of the angle of the +

*f*-axis. Two curves are shown: one for the light weight (

_{t}*f*= 8.4 N) and one for the heavy weight (

_{t}*f*= 17.3 N). The solid line shows the average applied constraint force over the 3 trials. The actual averaged data points are shown as dots, whereas the rest of the curve is a spline interpolation. The shaded region shows the range of normal forces measured over the 3 trials. The dotted line is identical to the solid line of average normal forces, except it has been shifted up or down so that its integral over all test angles is zero (zero mean). Figure 7 shows experimental results for subjects 1 through 8. The solid line shows the average applied constraint force for the light weight and the dotted line shows the average applied constraint force for the heavy weight.

_{t}Simple observation of the data indicates that subjects often apply force normal to the constraint, depending on the angle of the constraint and the direction of the tangent force, even though normal forces are not necessary for the task. In fact, for several of the subjects, the peak value of the normal force is about as large as the required force along the rail. As the subsequent sections show, the experimental data support the hypothesis that the data can be explained by a CSI objective function, allowing us to reconstruct each subject's isocost contours.

### Intertrial consistency

The intertrial variations in the normal force curves over the 3 trials are small, indicating that resolution of the force freedom is systematic rather than random. A glance at Fig. 6 (subject 1) shows that the envelope of normal forces over the 3 trials is narrow relative to the peak normal forces. For subject 1, the average variation in normal forces applied at a particular constraint angle, as a percentage of the range of average normal forces over all angles, was 16.6% for the light weight and 15.3% for the heavy weight. Other subjects displayed similar behavior. For subjects 1 through 8, the average intertrial variation was 22.5 ± 7.3% for the light weight and 19.4 ± 5.0% for the heavy weight.

### Convexity

As described in the methods section, a nonconvex objective function implies discontinuities in the normal force plot as a function of the tangential force angle. Although it is impossible to prove the absence of discontinuities in the underlying normal force curves from sampled data, the data in Fig. 7 appear to be indicative of smooth underlying curves. This supports convexity of the underlying objective function.

### Scale invariance

For our experiment, the scale-invariance hypothesis can be written where *f _{n,heavy}* and

*f*are the normal forces applied by a given subject at a given angle of the tangential force

_{n,light}*f*for the large (17.3 N) and small (8.4 N) tangential forces, respectively. The mean and SD of the residual

_{t}*z*for subjects 1 through 8, in Newtons, are −0.1 ± 2.5, 0.8 ± 1.3, 0.6 ± 2.0, −0.1 ± 1.7, 0.1 ± 2.2, 0.6 ± 1.6, 0.4 ± 1.7, and 0.2 ± 1.2. For the 8 subjects pooled the mean and SD of the residual

*z*is 0.3 ± 1.8. Notice that the mean of

*z*is close to zero.

To further test this hypothesis, we pooled all 8 subjects data and performed a 2-way ANOVA with “constraint angle” and “force magnitude” as experimental factors. The results showed no main effect for the “force magnitude” factor on the data (*F* = 1.143, *P* = 0.286). There was no evidence of system deviation from the scaling hypothesis across subjects for “constraint angle × force magnitude” interaction (*F* = 0.586, *P* = 0.884, adjusted *R*^{2} = 0.62). We then performed 2-way ANOVAs for each subject individually (*n* = 8). The results showed no main effect for the “force magnitude” factor except for one subject. The results also showed that there are some idiosyncratic deviations from the scale-invariance hypothesis within each subject attributed to the “constraint angle × force magnitude” interaction. However, the pooled analysis, along with correlation coefficients for a linear fit to the scaling hypothesis (0.876 ± 0.053), show that the data reasonably support scale invariance for the tangential forces we tested.

### Reconstructing isocost contours

The isocost contour reconstruction procedure outlined in the methods section results in a closed curve only if the normal force data is zero mean. We can see that the experimental normal force curves are in fact approximately zero mean. In general, the amount of uniform shift (raising or lowering) of a normal force curve to achieve zero mean is small relative to the maximum normal force, and in nearly all cases the shifted curve remains within the 3-trial variability (see, for example, Fig. 6). The ratio of the amount of shift |Δ*f*| to the range (from min to max) of the average normal forces for subjects 1 through 8 is 4.3 ± 2.7 percent for the light weight and 3.0 ± 3.0 percent for the heavy weight. This is significant, because although all scale-invariant objective functions predict zero mean normal force curves, this is not true for general non-scale-invariant objective functions. The fact that the experimental data are approximately zero mean is further evidence of the CSI objective function model.

To reconstruct a subject's isocost contour, we first shift the curve of average normal forces by subtracting the mean value of the normal force. This shifts the curve uniformly up or down and produces a zero mean curve. We then apply the numerical integration scheme described in the methods section. The results for subjects 1 through 8 at the hand position (0, 45 cm) in the shoulder frame are shown in Fig. 8. Because of the scale-invariance hypothesis, only the shape of the isocost contours is of interest; their sizes are arbitrary.

We can make a few general observations about the shapes of the isocost contours. As predicted, for each subject the shapes of the isocost contours derived for |*f _{t}*| = 8.4 N and 17.3 N are similar. The isocost contours are stretched in the

*f*direction relative to the

_{hy}*f*direction, indicating that a larger force in the

_{hx}*f*direction has the same “cost” as a smaller force in the

_{hy}*f*direction. This is not surprising for this configuration of the arm, and similar stretching has been observed in experimentally derived stiffness ellipses for the arm (Mussa-Ivaldi et al. 1985). The long axis of the isocost contour for most subjects passes approximately through the shoulder or leans to pass between the shoulder and the elbow. This is another common feature of stiffness ellipses.

_{hx}A local minimum (maximum) of an isocost contour is defined as a point such that no nearby point on the contour is closer to (further from) the origin, in a Euclidean sense. If **f**_{t} lies on a local maximum or minimum of an isocost contour, then the predicted normal force is zero. Therefore, every zero crossing in the normal force data predicts a local extremum in the isocost contour at the angle of the vector **f**_{t}. More specifically, if the slope at the zero crossing is negative (i.e., the normal force changes from positive to negative as the constraint angle increases), it is a local minimum, and otherwise it is a local maximum. Any closed curve containing the origin must attain an equal number of local minima and maxima, with a minimum of one each. Every average normal force curve in our experiments showed 4 zero crossings (Fig. 7), predicting 2 local maxima and 2 local minima in the isocost contours for all subjects.

The key points of the reconstructed isocost contours are that *1*) they are independent of any strong biomechanical modeling assumptions apart from CSI and *2*) they represent how the constraint force freedoms are used by subjects in solving static manipulation tasks.

### Configuration dependency of isocost contours

To investigate the dependency of the isocost contours on arm configuration, we performed similar experiments with subjects 9 through 14 at 5 different hand positions shown in Fig. 9*A*: (0, 45 cm) as before (labeled CENTER), (0, 55 cm) (labeled FAR), (0, 35 cm) (labeled NEAR), (31.8 cm, 31.8 cm) (labeled RIGHT), and (−31.8 cm, 31.8 cm) (labeled LEFT) in the shoulder frame. LEFT and RIGHT are obtained from CENTER by rotating ±45° about the shoulder. The procedure is the same as with the previous experiments, but only one weight (1.2 kg, tangential force of 11.8 N) is used.

The isocost contours for a typical subject are shown in Fig. 9*B*. We can see that the isocost contour at RIGHT and LEFT have shapes similar to the isocost contour at CENTER, rotated approximately ±45°. In contrast, changing both the elbow and shoulder angle (FAR and NEAR positions) causes the isocost contour to both change shape and rotate. For example, the isocost contour at FAR becomes more anisotropic because of the greater extension of the elbow. These results suggest that the objective function may be better expressed in the joint-torque space rather than the hand-force space, as discussed next.

### Invariance in joint space

To investigate the dependency of the isocost contours on arm configuration, we transform isocost contours in the hand-force space to isocost contours in the joint-torque space (*Eq. 3*). Figure 10 shows the isocost contours of Fig. 9 expressed in the joint-torque space. We see that the joint-torque isocost contours are nearly constant over the different arm configurations. These results are typical for all subjects.

Because these joint-torque isocost contours appear rather elliptical, we fit the data to ellipses as described in the methods section. The centers of the fitted ellipses are close to the origin, indicating that flexion and extension torques have similar costs. For each subject (subjects 9 through 14), the mean value and SD of the eccentricities of the fitted ellipses at the 5 hand positions are 0.879 ± 0.061, 0.841 ± 0.071, 0.849 ± 0.052, 0.835 ± 0.064, 0.833 ± 0.128, and 0.818 ± 0.067. The small SDs indicate that the eccentricity is essentially constant as the arm configuration changes. Subjects' ellipse orientations are 39.4 ± 8.0, 40.3 ± 8.6, 35.3 ± 5.1, 50.8 ± 9.8, 47.3 ± 6.6, and 34.5 ± 8.7°, indicating that the orientation of the ellipse changes little with large changes in joint angles. Finally, the similarity of eccentricity and orientation of the fitted ellipses show that the joint-torque isocost contours are approximately the same for all subjects, with eccentricity of 0.843 ± 0.073 and orientation of 41.3° ± 9.4°. Keep in mind that identical joint-torque isocost contours can lead to different hand-force isocost contours, depending on the length of subjects' upper arms and forearms.

Assuming zero offset from the origin, a torque ellipse can be expressed in the form τ* ^{T}W*τ =

*c*, where

*c*> 0 is a constant (i.e., the WT2 model of the methods section). A typical

*W*matrix obtained by a least-squares fit to a joint-torque isocost contour is The diagonal terms of the symmetric 2 × 2 positive-definite matrix

*W*indicate that elbow torques have a somewhat greater cost associated with them than shoulder torques. The negative constants in the off-diagonal entries mean that the total cost is increased if the elbow and shoulder torques have an opposite sign, and decreased if they have the same sign. We believe that the decreased cost associated with having both torques be the same sign can be reasonably attributed to the presence of biarticular muscles, which create torque of the same sign about both joints. Presumably the off-diagonal entries would be zero in the absence of biarticular muscles—the cost would have no dependency on the relative values of shoulder and elbow torques.

In Mussa-Ivaldi et al. (1985), stiffness at the hand is also interpreted in terms of a positive-definite stiffness matrix *K,* invariant to the arm configuration when expressed as a stiffness matrix *R* in joint space. If δθ is the vector of joint angle displacements from the equilibrium point, the potential energy stored in the springlike muscles can be written δθ* ^{T}R*δθ. By the relation τ =

*R*δθ, this potential energy can be rewritten as τ

^{T}R^{−1}τ. Thus the inverse of the stiffness matrix

*R*

^{−1}can be used as the

*W*matrix in the WT2 model with the following interpretation: iso-cost contours are contours of constant potential energy, and subjects choose a hand force that minimizes the stored potential energy in the spring-like muscles. The inverse of the stiffness matrix found in Mussa-Ivaldi et al. (1985), normalized so that the top left element is identical to the

*W*matrix given above, is One difference from our result is that this matrix implies greater cost for shoulder torques than for elbow torques. Nonetheless, as we see in the next section, this

*R*

^{−1}matrix reasonably predicts subjects' behavior, somewhat better than the identity matrix implicit in the T2 model, as a result of the off-diagonal terms.

### Biomechanical modeling

The 9 force-generation models described in the methods section predict hand-force isocost contours as shown in Fig. 11 using the physiological parameters of Table 2. Note the strong anisotropy of the MS*k* isocost contours, attributed to the large *PCSA* of the uniarticular shoulder muscles. For the WT2 isocost contour shown, the positive-definite *W* matrix is the inverse of the stiffness matrix *R*^{−1} of Mussa-Ivaldi et al. (1985), as described above.

For 8 of the force-generation models, we calculated the correlation coefficients between the experimental data of subjects 1 through 8 at the (0, 45 cm) hand position (CENTER) and the predicted data for each subject, using the physiological parameters in Table 2 (Gomi 2000) and the measured upper and forearm lengths for each individual subject. The results are shown in Fig. 12. The HAND model is not included because it predicts zero normal forces for all experiments, and thus can be discarded as a candidate.

The results show that the MT2 (0.81 ± 0.12) and MT3 (0.81 ± 0.11) models fit the experimental data best, followed by the WT2 (0.75 ± 0.15) and T2 (0.66 ± 0.16) models. The MT1, MS1, MS2, and MS3 have low (or even negative) correlation coefficients for most subjects. The MT1 and MS1 models predict polygonal isocost contours, predicting discontinuities in the normal force plots, which are not evident in the data. They can apparently be discarded relative to models with exponents of ≥2.

A weakness of the approach for the MT*k* and MS*k* models is that we are forced to assume physiological parameters, and published parameters demonstrate significant variations. For instance, the parameters used by van Bolhuis and Gielen (1999) are significantly different. Using their parameters, the correlation coefficients are as shown in Fig. 13. With these parameters, the MS3 and MS2 models outperform the MT3 and MT2 models, although the correlation is somewhat less than obtained with the MT2 and MT3 models under Gomi's parameters. One reason the MS*k* models become more competitive is that the *PCSA* values of the shoulder uniarticular muscles given by van Bolhuis and Gielen are closer to those of the other muscles, meaning that the isocost contours of the MS*k* models are more isotropic than those under the Gomi parameters.

For a visual interpretation of the effect of different physiological parameters (muscle *PCSA* and moment arms), Fig. 14 plots the MT*k* and MS*k* (*k* = 2, 3) joint-torque isocost contours using the parameters from Gomi (2000) and van Bolhuis and Gielen (1999). For the MT2 and MT3 models, Gomi's parameters yield isocost contours more closely resembling the experimental data. For the MS2 and MS3 models, the isocost contours obtained using van Bolhuis and Gielen's parameters are superior. This is consistent with the correlation coefficients results.

To measure the sensitivity of the model predictions to the physiological parameters, we randomly generated 16,000 sets (1,000 sets for each of subjects 1 through 8 with both light and heavy weights) of physiological parameters using means and SDs derived from published data (Garner and Pandy 2003; Gomi 2000; Gribble et al. 1998; Lemay and Crago 1996; van Bolhuis and Gielen 1999), as shown in Table 3. For simplicity, in the Monte Carlo tests we assume all physiological parameters are uniformly distributed within 1SD of the mean. The resulting correlation coefficients are 0.28 ± 0.24 for MT1, 0.43 ± 0.25 for MT2, 0.48 ± 0.24 for MT3, 0.07 ± 0.21 for MS1, 0.24 ± 0.21 for MS2, 0.31 ± 0.21 for MS3, 0.66 ± 0.16 for T2, and 0.75 ± 0.15 for WT2. (Note that the T2 and WT2 models do not use the physiological parameters, so their SDs show only the variation between subjects.) The SDs in the MT and MS models are relatively large, meaning that the model predictions change dramatically according to changes in the physiological parameters.

To better understand the large SDs in the Monte Carlo tests, we also calculated a sensitivity index for each parameter, defined as (*P _{j}*/

*C*)(∂

_{i}*C*/∂

_{i}*P*), where

_{j}*C*is the correlation coefficient and

_{i}*P*is the physiological parameter. This sensitivity is a linear estimate of the percentage change in the variable

_{j}*C*caused by a 1% change in the parameter

_{i}*P*The sensitivity index is relatively small for all parameters (the mean value is <1 for most models, and the SD is very small except in the MT1 and MS1 models), which shows that the model predictions are not particularly sensitive to small changes in any particular parameter. The large percentage variations in the

_{j}.*PCSA*parameters reported in the literature, however, when taken together, lead to significant variations in the model predictions.

We also did Monte Carlo tests to compare the predictions with the experimental data for subjects 9 through 14 at 5 different hand positions. The results are similar to the previous results for subjects 1 through 8. Although there is little conclusive that we can say to validate or invalidate particular optimization models (other than the MT1 and MS1 models, which are inferior to their higher-exponent counterparts), the simple WT2 model is competitive for any set of physiological parameters. We return to this in the discussion.

## DISCUSSION

The primary findings of this paper can be summarized by 3 points. *1*) The experimental results show that subjects apply significant constraint forces, even though they are unnecessary for the task. Each subject's applied normal forces follow a consistent pattern, indicating that the constraint force freedom is resolved in a systematic manner. *2*) The data can be interpreted in terms of a convex scale-invariant (CSI) objective function on forces applied at the hand. The level sets of a subject's objective function represent the sets of hand forces with equal “cost” to the subject. These level sets, or isocost contours, can be reconstructed directly from the experimental data without any biomechanical modeling. *3*) The isocost contours, when expressed in the joint-torque space, are approximately invariant to the configuration of the arm. They are also similar across subjects.

Our use of a CSI objective function is as a descriptor and predictor of behavior in constrained tasks. It does not require that we commit to a particular interpretation of it. For instance, it may describe that subjects' natural behavior minimizes some notion of “effort.” Although we find this interpretation appealing, the objective function could also describe properties of the motor control system that resist simple interpretation of behavior as minimization of effort. For example, Krylow and Rymer (1997) argue that smooth motions in motor control may arise largely from intrinsic muscle mechanics.

Our observation that objective functions in the hand-force space are CSI is consistent with most existing models for static muscle load-sharing (sum of muscle group stresses or tensions raised to the power of ≥1). Each of these models defines a CSI objective function in the muscle stress or tension space, and linear mappings (from muscle stress to muscle tension, from muscle tension to joint torques, and from joint torques to hand forces) preserve convexity and scale invariance. Previous research has also found that as the direction of an applied force remains fixed but the magnitude is scaled, the muscle activation and force patterns also scale, further evidence of scale invariance (Buchanan et al. 1986; Flanders and Soechting 1990; Valero-Cuevas et al. 2000). We note that there is recent evidence for an objective function that is the sum of linear and squared terms of the muscle tension (FCT van der Helm, personal communication), which is not scale invariant, but such scaling effects would become more noticeable at larger percentages of maximum force, which were not explored in this study.

### Relationship to stiffness ellipses

The approximately elliptical shape of the joint-torque isocost contours, and their approximate invariance to arm configuration, brings to mind the stiffness ellipse description of natural stiffness at the hand (Gomi and Osu 1998; McIntyre et al. 1996; Milner 2002; Mussa-Ivaldi et al. 1985; Perreault et al. 2001). The joint-torque isocost contours of Figure 10 and the correlation coefficients of Figures 12 and 13 also indicate that the 3 unique entries of the symmetric positive-definite weighting matrix *W* in the WT2 model provide a reasonable low-complexity description of subject behavior. This model better captures the shape of the joint-torque isocost contours (approximately elliptical, not aligned with the τ_{1}–τ_{2} axes) than the joint-torque circle of the T2 model. The *W* matrix used in our analysis is the inverse of the joint stiffness matrix from Mussa-Ivaldi et al. (1985), which leads to the interpretation that subjects minimize the potential energy stored in springlike muscles while stabilizing the manipulandum.

It is important to keep in mind, however, that stiffness ellipses and isocost contours are *not* the same thing; the former express the behavior of the arm in response to brief perturbations, whereas the latter express how subjects actively resolve force freedoms.

### Generalizations and assistive guides

The work described in this paper can be extended in at least 2 ways: by extending to partially constrained reaching tasks, similar to the crank-turning work of Russell and Hogan (1989) and Svinin et al. (2003), where the arm dynamics become significant; and by increasing the number of degrees of freedom of the arm and the number of force freedoms to be resolved. We have begun work toward the former extension by designing and building a planar manipulandum that can implement arbitrary constraint curves in the plane (Worsnopp et al. 2004). This manipulandum is superior to traditional robotic manipulanda at enforcing smooth constraints because the constraint is generated by a steerable wheel rolling on a table. For the latter extension, the notions of orthogonal “tangential” and “normal” (workless) force and torque subspaces must be generalized properly, according to the kinetic energy metric of the manipulandum.

One reason we are interested in these generalizations is that an eventual goal of this work is to design constraint surfaces to assist a human in manipulating a load from one configuration to another. Constraint surfaces are passive and inherently safe to interact with, and a properly designed constraint surface or guide rail can improve the ergonomics of a repetitive material handling task. Before tackling this problem, however, we must understand how humans naturally take advantage of the presence of kinematic constraints. The work described in this paper is a step toward that understanding.

## APPENDIX

To see that the normal force curve must be zero mean to obtain a closed isocost contour by the integration procedure outlined in the results section, consider one step of the integration process in Fig. 4. The angle β and the associated tangent direction γ are derived from the normal force plot, ψ is the constraint angle, *r* is the current radius of the isocost contour, dβ is the increment of β, and d*r* is the increment of *r.* We have (A1) (A2) where *f _{n}* is the normal force and

*f*is the tangential force. Integrating, we can write

_{t}*r*as a function of β (A3) For the isocost contour to close, we must have

*r*(2π) =

*r*(0) implying (A4) From

*Eq. A2*we have and taking the derivative we get So

*Eq. A4*can be written Because we have We know that γ = (π/2) + ψ (for

*f*< 0) or γ = (3π/2) + ψ (for

_{n}*f*> 0), so finally we get as the condition for a closed curve, which means the normal force plot must have zero mean.

_{n}## GRANTS

This work was supported by National Science Foundation Grants IIS-9875469 and IIS-0082957.

## Acknowledgments

We thank the anonymous reviewers for constructive suggestions to improve this work, as well as Drs. F. A. Mussa-Ivaldi, E. Perreault, J. L. Patton, C. Mah, J. P. Dewald, and J. Hidler and the entire robotics lab at the Rehabilitation Institute of Chicago for stimulating discussions on this topic.

## Footnotes

↵1 In generalized force spaces representing both forces and torques, “tangential” forces do work on the manipulandum, whereas “normal” forces are workless. Tangential and normal forces are orthogonal by the inertia metric of the manipulandum rather than the standard Euclidean metric. For our system, the inertia metric is equivalent to the Euclidean metric.

The costs of publication of this article were defrayed in part by the payment of page charges. The article must therefore be hereby marked “

*advertisement*” in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.

- Copyright © 2005 by the American Physiological Society