Abstract A method is described for computing the Center of Mass (COM) from empirical estimates of the Center of Pressure (COP) obtained by means of a force platform. The method is based on a biomechanical model of sway movements in quiet standing, according to which the horizontal acceleration of the COM is approximately proportional to the COM-COP difference. The equations are solved by approximating the solution with best fitting spline functions. The implications for movement control are discussed.

Introduction. In spite of its apparent simplicity, the nature of the control mechanisms that allow humans to stabilise postural sway is still an object of controversy. Sway is the persistent oscillation of the COM (Center Of Mass), the controlled variable, in the AP (Antero-Posterior) and ML (Medio-Lateral) planes. This variable is not measurable in a direct way, although is frequently confused with the COP (Center Of Pressure) that is typically measured by means of a force platform. The difference between the two trajectories is small but is significant and the study of its structure is an important method for addressing fundamental motor control problems. The most general way to compute the trajectory of the COM is to apply the definition (see, for example, Winter et al. 1998), i.e. measuring the COM of each body segment by means of stereophotogrammetric methods and then adding them up, weighted according to their mass. However, this is very cumbersome, especially if we think of clinical applications, and requires a very critical calibration. More practical methods are limited to force platform data. An empirical approach is based upon the fact that the two trajectories have similar shapes but the COP trajectory appears to be more "noisy"; thus, the computation can be seen as a type of "filtering" (Benda et al. 1994). A more formal approach (Shimba 1984, Zatsiorsky and King, 1998) takes into account that the second time derivative of the horizontal COM trajectory is proportional to the horizontal component of the ground reaction force: this force is measured by means of a fully equipped force platform and integrated twice, with a special technique for recovering the unknown initial conditions. Here we propose a new method that is based on a biomechanical model but does not require to measure the horizontal component of the ground reaction force and thus can operate with a much cheaper force platform. The model. In the following we only consider AP oscillations. However, the model and the algorithm can be applied to ML oscillations as well. Figure 1 shows the human “inverted pendulum” in v which we single out the ground reaction force f = ( f H , f V ) , the force of gravity mg, and the ankle torque τankle due to the muscles. For the foot we can write an equilibrium equation: f H δ + fV u + τ ankle = 0 (1)

1

Computing the COM from the COP

Human Movement Science 18, 759-67, 1999

where δ,u are, respectively, the vertical and horizontal displacements of the COP with respect to the ankle. For the sway movements of the rest of the body the following equations apply f V − mg = m&z& ≈ 0 f H = m&y& (2) &y& d τ ankle + mgy = ( Iθ& ) ≈ Iθ&& ≈ mh 2 k s = mhk s &y& dt h where m is the weight of the body excluding the feet, h is the distance of the COM from the ankle, I is its moment of inertia with respect to the ankle joint, and ks is a shape factor that depends on the distribution of mass in the body. In particular, ks=1 if all the mass were concentrated in the COM (in that case I=mh2) and ks=1.33 if it were distributed in a uniform way along a thin rod (in that case I=4/3mh2). For the human body ks is closer to the latter estimate than the former one. There are three approximations in the model, that are well satisfied for sway in quiet standing: (i) the vertical acceleration of the COM is negligible, (ii) the moment of inertia is constant, and (iii) the angular acceleration is proportional to the horizontal acceleration of the COM. We can now combine Eq. 1 and 2, obtaining the following sway equation: g &y& = ( y − u) (3) (hk s + δ ) or, more simply, &y& =

g ( y − u) he

(4)

having defined an "effective distance" he=ksh+δ. In particular, the equation says that the COM-COP difference is proportional to the horizontal component of the ground reaction force.

Fig. 1. Mechanics of quiet standing. 2

Computing the COM from the COP

Human Movement Science 18, 759-67, 1999

The algorithm. The problem is to integrate Eq. 4. Since we do not know the initial state, we cannot use standard numerical integration techniques that, in any case, have numerical instability problems when applied to unstable systems. However we know that the solution must be a smooth function of time and we can seek it within a suitable class of functions by means of a variational method. This requires to define a functional to be minimised, using the equation as a constraint. We chose a piece-wise polynomial representation of y(t) in terms of B-spline functions (De Boor, 1978) because it is linear_in_the_parameters. Since the system equation (4) is linear then the functional is quadratic and the solution can be provided by the LSE method. In order to have an algorithm whose performance is independent of the duration of the experiment, we adopted a moving-window paradigm, i.e. the LSE is iterated on shifted strings of data of duration ±Tw. The estimation step. Let m be the number of B-splines, n the number of samples in the time window, { p i } the set of unknown coefficients that weight the influence of each B-spline. The B-spline approximation of y(t) can be written as follows: m

y (t ) = ∑ Bi (t ) ⋅ pi

(5)

i =1

It must be m<

&& (t ) ⋅ p &y&(t ) = ∑ B i i

(6)

i =1

If we now consider the whole set of samples within the observation window (t1 … tn), we can write the following linear relationships for y and &y& : y (t1 ) B1 (t1 ) . . Bm (t1 ) p1. . . . ⋅ . = (7) . . . p m y (t n ) B1 (t n ) . . Bm (t n ) && (t ) . . B && (t ) &y&(t1 ) B 1 1 m 1 p1. . . . ⋅ . = . . . p && && (t ) m &y&(t n ) B1 (t n ) . . B m n In matrix notation this can be written in a more compact way: && ⋅ p &y& = B y = [B ] ⋅ p (8) We apply it to Eq. 4 and obtain the following matrix equation && ⋅ p − g B ([B ] ⋅ p − u ) = e (9) he where e is an error-vector, whose norm must be minimised by computing the optimal parameter vector. We introduce a matrix of coefficients [A]: [A] = B&& − g h e [B] (10)

[]

[]

[]

3

Computing the COM from the COP

Human Movement Science 18, 759-67, 1999

and re-write Eq. 9 more simply as [A] ⋅ p + g / he u = e Thus we can solve the LSE problem by using the pseudo-inverse of [A]:

(

)

T pˆ = arg min (e ⋅ e) = − g / he [ A]T ⋅ [ A] p

−1

⋅ [ A]T ⋅ u

and the optimal estimate of the COM profile for the time instant to in the middle of the observation window is computed taking into account Eq. 5: yˆ (to ) = ∑ Bi (to ) pˆ i

(11) (12)

(13)

i

The iteration procedure. The optimization step is performed for each position of the time window, yielding a sequence of estimates of the COM position given by Eq. 13. There are border effects and the simplest way to deal with them is to cut the initial and the final Tw seconds of recorded data. If the breakpoints are uniformly spaced, the pseudo-inverse matrix is a constant and can be computed only once. All the computations were performed using the spline package of Matlab. The parameters of the algorithm are the width of the fitting window (Tw), the number of knots (m), and the biomechanical parameter (g/he). The choice of the first two parameters is not very critical. In qualitative terms we can say that the fitting window must be large enough to include about one main oscillation of the COM: a value of ±1s is appropriate. The number m of knots should be determined according to the frequency band of u(t) which is of the order of 2-3 Hz. Thus a sensible choice is to have one knot every 100-150ms. The last parameter (g/he) must be estimated from anthropometric data. However, it is possible to show that the reconstruction algorithm is weakly dependent upon this parameter. In particular, we may consider the frequency response related to Eq. 4: G=

− g / he Y ( jω ) k = = 2 2 U ( jω ) ( jω ) − g / he ω + k

(14)

where k=g/he. From this it is easy to obtain the following relation dG dk k = (1 − 2 ) (15) G k ω +k Therefore, the percentual error in the COM-COP estimate is always less than the corresponding error in the k parameter. In particular, if we consider a typical case in which he=1m (thus k=9.81 s-2) and take into account that the u-signal is mostly located around 0.5Hz, we have that dG/G=0.501 dk/k; for smaller frequency values the coefficient is reduced even further. The empirical sensitivity analysis supports this evaluation. Experiments and validation. Figure 2 shows an example of application of the algorithm. The trajectory of the COP was measured by means of a force platform developed in the lab, with a surface of 50x50 cm2 and four estensimetric load cells. The COP trajectory was sampled at 50 Hz and then low-pass filtered with a 5 Hz cutoff frequency. The input trajectory of the COP (limited to the AP component) is denoted by the thicker line and the computed trajectory of the COM by the thin line. The he parameter was set by measuring the distance between the iliac spine and the malleolus (for h), the distance between the ankle and the ground (for δ) and choosing a value of 1.2 for ks. The reconstructed COM trajectory is very weakly dependent on the value of he. For example, the trajectories obtained with a ±10% variation of the parameter are graphically coincident with the nominal one except in the peaks: the average difference was of the order of 10-4 mm and the maximum difference less than 0.1 mm.

4

Computing the COM from the COP

Human Movement Science 18, 759-67, 1999

Figure 2. Reconstraction of the COM trajectory. Thick line: measured COP curve; thin line:reconstructed COM curve. The validation of the method was carried out by means of an inverted pendulum of known geometry, stabilised by means of springs. The springs were chose in order to have a resonant frequency of the pendulum similar to the predominant harmonic of human sway. The pendulum was put in motion by gentle pushes in the AP direction that caused a sequence of damped oscillations. The trajectory of the COP was measured by means of the same force platform mentioned above and the corresponding COM trajectory was measured by means of an optical stereophotogrammetric system (MacReflex by Qualisys) with a reflective marker placed as close as possible to the COM. The COP data were then processed according to the described algorithm, yielding the estimated COM trajectory. The discrepancy among the two values (estimated vs. measured COM: the percentual amplitude error) was on average 8.95% of the input range of COP data. The phase error, computed as the relative time shift between corresponding peaks of the two COM curves was undetectable at the sampling rate of 50 Hz. Thus it seems that the proposed method of extracting COM data from a basic force platform is reliable and simple enough to be applied in a clinical environment.

Discussion The proposed method for recovering the COM trajectory from the COP trajectory is based on the fact that during sway movements the horizontal acceleration of the COM is approximately proportional to the COM-COP difference: &y& ≈ g e / h ( y − u ) . We used this equation as a constraint for the least squares fitting problem. For the brain it is the equation of an unstable plant to be controlled, where y is the controlled variable and u is the control variable. In another paper (Morasso et al. 1999) we discuss biologically plausible strategies of control, based on the notion of “internal model” (Wolpert and Kawato 1998). Here we wish to focus on the biological significance of the COM-COP difference, after having established that the stiffness stabilisation of upright posture is not biologically plausible (Morasso and Schieppati 1999). In fact, for any kind of feedback stabilisation mechanism the estimate of y is essential

5

Computing the COM from the COP

Human Movement Science 18, 759-67, 1999

but unfortunately no sensory channel provides a direct measurement of this variable. Thus we must face the problem of figuring out plausible ways for the brain to obtain indirect estimates and Eq. 4 suggests a computational strategy. By using the plantar cutaneous receptors sensitive to the shear forces the brain can recover the COM-COP difference because it is proportional to the horizontal component of the ground reaction. Moreover, the COP position can be detected by other plantar receptors, sensitive to the vertical component of the ground reaction, and thus we can conclude that the COM position is stored in the patterns of activation of plantar receptors and can be recovered by an appropriate process of multisensory fusion. In any case, the COM-COP difference has a direct behavioural significance and thus it might provide useful clinical indicators in the study of balance disorders. The described algorithm is a robust mechanism for obtaining such data from standard posturographic measurements. References Benda, B.J., Riley, P.O., Krebs, D.E. (1994) Biomechanical relationship betwenn the center of gravity and center of pressure during standing. IEEE Trans. on Rehabil. Engineering, 2,3-10. Benson, A.J., Spencer, M.B., Stott, J.R. (1986) Thresholds for detection of the direction of whole body linear movement in the horizontal plane. Aviat. Space Environ. Med., 57, 1088-1096. De Boor, C. (1978) A practical guide to splines. Springer Verlag, New York. Fitzpatrick, R., McCloskey, D.I. (1993) Proprioceptive, visual and vestibular thresholds for the perception of sway during standing in humans. J. Physiol., 478, 173-176. Morasso, P.G., Baratto, L., Capra, R., Spada, G. (1999) Internal models in the control of posture. Neural Networks, 12, no. 7-8. Morasso, P.G., Schieppati, M. (1999) Can muscle stiffness alone stabilize upright stabding? J. Neurophysiol., in press. Nashner, L.M. (1982) Adaptations of human movement to altered environments. Trends neurosci., 5, 358-361. Shimba, T. (1984) An estimation of center of gravity from froce platform data. J. Biomech., 17, 53-60. Winter, D.A., Patla, A.E., Prince, F., Ishac, M., Gielo-Perczak, K. (1998) Stiffness control of balance in quiet standing. J. of Neurophysiol., 80, 1211-1221. Wolpert, D.M., Kawato, M. (1998) Internal models of the cerebellum. Trends in Cognitive Sciences, 2, 338-347. Zatsiorsky, V.M., King, D.L. (1998) An algorithm for determining gravity line location from posturographic recordings. J. Biomech., 31, 161-164.

6