## Abstract

Standing balance is a simple motion task for healthy humans but the actions of the central nervous system (CNS) have not been described by generalized and sufficiently sophisticated control laws. While system identification approaches have been used to extracted models of the CNS, they either focus on short balance motions, leading to task-specific control laws, or assume that the standing balance system is linear. To obtain comprehensive control laws for human standing balance, complex balance motions, long duration tests, and nonlinear controller models are all needed. In this paper, we demonstrate that trajectory optimization with the direct collocation method can achieve these goals to identify complex CNS models for the human standing balance task. We first examined this identification method using synthetic motion data and showed that correct control parameters can be extracted. Then, six types of controllers, from simple linear to complex nonlinear, were identified from 100 s of motion data from randomly perturbed standing. Results showed that multiple time-delay paths and nonlinear properties are both needed in order to fully explain human feedback control of standing balance.

## 1 Introduction

Standing balance is a simple motion task which allows researchers to easily investigate how the human central nervous system (CNS) uses feedback control. Understanding the CNS in human motions can not only help clinical treatments of movement disorders but also can inspire the design of controllers for assistive devices, such as exoskeletons and prostheses, to generate human-like movements. Numerous studies have treated the CNS in human standing balance as a postural feedback controller and system identification methods have been used to find its mathematical model [113]. In the majority of these studies, external mechanical stimuli, such as push/pull forces and standing platform motions, were used to evoke participants' body sway motions at variety statuses. Indirect identification approaches were normally used, in which a closed-loop mathematical model of the human standing balance system is required. It has been reported that the indirect approach (closed-loop system identification) can avoid the bias caused by the open-loop identifications that just rely on the information of CNS inputs (joint motions) and outputs (joint torques) [14,15]. Both nonparametric and parametric system identification methods have been used. However, limitations exist in both study directions which prevented the identification of good generalized models of the postural control system in the CNS.

In nonparametric postural controller identifications, the CNS is described by a frequency response function (FRF) [11,12,16]. The frequency domain identification method was often used to find the FRF of the closed-loop model, including the CNS that can best explain the experimental data [15,1720]. Because the length of data does not matter when transferring to frequency domain through the Fourier transformation [21], this approach can be easily applied on long duration standing balance data that recorded from experiments where the random or multisine external perturbations were applied. In general, longer duration motion data can provide more information about the CNS, which allows identification of more complex, or more generalized models. However, the frequency domain approach requires multiple perturbation sources to identify a multi-input and multi-output system [21], which is typical for human standing balance task since at least the ankle and hip strategies were used. Developing the hardware to provide perturbations across multiple body segments is difficult for both research and clinical applications. In addition, this approach treats the system as linear, which eliminates the ability to identify nonlinear properties of the CNS. Nonlinear effects, such as threshold and saturation, are known to exist in sensory and spinal circuits. Furthermore, the human body has many nonlinear components, e.g., nonlinear multibody dynamics and nonlinear muscle properties. Studies have shown that healthy humans use different feedback control gains to control posture under different amplitudes of external perturbations [4,7], indicating that the generalized standing balance controller is nonlinear.

Parametric identification has been used to find gains and other parameters in a predefined control structure that best explain the experimental data. It has the advantage that a single stimulus or perturbation source can be used [22]. Trajectory optimization can be used to solve this identification problem [4], to find controller parameters that best explain the measured output of the system. In this approach, both the plant dynamics and the CNS controller can be nonlinear since identification is performed in the time-domain. However, only simple linear controllers have been identified, using data from short duration balance motions that were perturbed by short ramp displacements [4,7,23]. One difficulty is that the shooting method [24] used in these studies have difficulty with long duration trajectory optimizations, as the standing balance is naturally an unstable system. Eigenvalue constraints have been used which increased the chance of getting stable motions [4], but these are not necessary conditions for the stability of nonlinear systems. Shooting methods also have a tendency to find local optima [25] when data have detailed features in the time domain. Identifications on long duration data are essential to get generalized motion controllers that can explain human responses in different circumstances, such as under multiple ramp perturbations that have different amplitudes. Only recently, parametric identification has been done on long duration randomly perturbed balance data [22]. However, the objective was to compare the FRFs in the frequency domain, which still requires the assumption that the system is linear.

The goal of this study is to show that trajectory optimization with the direct collocation method can identify generalized and nonlinear standing balance controllers from long duration motion data. Direct collocation has recently been used in predictive simulations of human movement [2628], performing better than shooting methods, but to our knowledge, has not been used for identification of human feedback control. Therefore, in this paper we first validate this method by doing identification on simulated data where the feedback control parameters are known (simulation study). In the simulation study, we also investigated the effect of data length on the accuracy of identified control gains. Then, we identify six types of feedback controller, from simple linear to complex nonlinear, on 100 s of measured human standing balance motions (experimental study).

## 2 Methods

An indirect identification approach was used in this study [3,4,17,22], in which a mathematical model was created to represent the standing balance system and an optimization method was used to fit to experimental data by optimizing the model parameters. The mathematical model was a closed-loop system, where a body dynamics (plant) model and a feedback controller model were both included. Body dynamics was simplified as a double-link pendulum without the knee joint, since ankle and hip strategies are mostly used for standing balance [29,30]. Equations of the body dynamics model can be found in Appendix  A. The feedback control parameters P were optimized to produce the best fit between the output of the simulated system and the experimental data (Fig. 1).

Fig. 1
Fig. 1
Close modal

### 2.1 Standing Balance Experiments.

Experiments were performed on eight able-bodied participants (seven males, one female, age 18–35 years) with approval from the Institutional Review Board of Cleveland State University (study No. IRB-FY2018-40). Participants were asked to stand on an R-Mill instrumented treadmill (Forcelink, Amsterdam, the Netherlands) with their arms crossed in front of their chest and instructed to keep balance without taking a step. The perturbation consisted of random anterior–posterior platform motions. The perturbation signal was designed using random square pulses with five position values ([−5, −2.5, 0, 2.5, 5] cm), and six durations ([0.25, 0.5, 0.75, 1.0, 1.25, 1.5] s). Amplitudes and durations were randomly combined to generate a 300 s perturbation signal. This random signal was sent to the treadmill control module in the d-flow software (Version 3.26.0, Motek, Amsterdam, the Netherlands) and executed through the treadmill's “sway” actuation. Twenty-seven reflective markers [31] were placed on each participant to record their reactions using a 10-camera motion capture system (Osprey, Motion Analysis Corp., Rohnert Park, CA). Hip and ankle joint angles were calculated from these recorded marker data after applying a double second-order Butterworth filter with the cutoff frequency of 16 Hz. Dynamic parameters of the human body model were generated for each participant, in which the segment lengths were calculated based on the recorded marker data and mass properties were calculated based on the body weight and the anthropometry table of Winter [32]. Five markers were placed on the rigid frame of the treadmill to record the executed sway motion. After applying the 16 Hz low pass filter, acceleration was calculated using a three-point finite difference on the averaged coordinates of these five markers [33]. All data were sampled at 100 Hz.

The experiment consisted of a 300 s quiet standing trial, followed by a 300 s perturbation trial. Then participants were asked to sit down and rest for 300 s. After the rest period, a 300 s perturbation trial, using the same perturbation as the previous trial, was recorded. Finally, participants were asked to do another trial of 300 s quiet standing. General information of the participants are shown in Table 1. The detailed description of the experiment, collected data, data processing, as well as its preliminary analysis were included in the Zenodo project.2

Table 1

Information of the eight participants in the order of data collection date

IdGenderAge (yr)Height (m)Mass (kg)
1Male221.6074.29 ± 0.26
2Female48.37 ± 0.21
3Male181.8079.12 ± 0.20
4Male271.7863.10 ± 0.16
5Male321.7970.56 ± 0.19
6Male351.6558.24 ± 0.27
7Male281.7568.75 ± 0.17
8Male271.6360.33 ± 0.19
IdGenderAge (yr)Height (m)Mass (kg)
1Male221.6074.29 ± 0.26
2Female48.37 ± 0.21
3Male181.8079.12 ± 0.20
4Male271.7863.10 ± 0.16
5Male321.7970.56 ± 0.19
6Male351.6558.24 ± 0.27
7Male281.7568.75 ± 0.17
8Male271.6360.33 ± 0.19

### 2.2 Simulated Motion Data.

The identification method was first validated through a simulation study. The closed-loop mathematical model in the indirect identification platform was used to generate simulated data. The dynamic parameters of participant 4 were used for the double-link pendulum model. The perturbation input was the measured platform acceleration in the first perturbation trial of participant 4. The postural feedback controller, a full state feedback proportional-derivative (PD) controller with time delay, was the same as in the identification study by Goodworth and Peterka [22]. This controller is specified in Eq. (3). Pink noise (similar to Goodworth's study) was added at the motion feedback loop of each joint as sensor noise. One hundred seconds of simulation data was generated using the midpoint Euler method [34] and sampled at 100 Hz. Ten different realizations of sensor noise (with different random number seeds) were used to simulate ten different trials. Figure 2 shows one trial of simulated motion data with one realization of sensor noise under the external perturbation.

Fig. 2
Fig. 2
Close modal

### 2.3 Controller Identification on Simulated Data.

Controller identifications on the simulated motion data were done through the described indirect approach. Five lengths (10, 30, 50, 70, and 90 s) of the simulation data (Fig. 2) were selected for the identification to check the effects of data length on the identification results. Lower and upper bounds of the identified control parameters were set as 0 and 2 times their true values. Ten optimizations were performed on each controller identification problem (controller identification on each period of data) with random initial guesses within the bounds to increase the chance of finding global optimum. The identification result with the lowest objective function among ten optimizations was selected as the best result for the identification problem and was shown in the results section.

### 2.4 Controller Identification on Experimental Data.

One hundred seconds of experimental data was used from the middle of each perturbation trial to avoid the adaptation period at the beginning, while also minimizing the effect of fatigue at the end. The quiet standing trials were not used. Data from the last six participants (Table 1) were used, since the first two participants were in the preliminary testing of the experiment protocol. In total, 12 periods of experimental data (six participants and two perturbation trials for each participant) were used to identify the standing balance controllers.

Six controller types, from simple linear to complex nonlinear, were identified on each selected data trial. Four of them are linear: the PD controller, the full-state proportional-derivative (FPD) controller, the full-state proportional-derivative feedback with time delay controller (FPDTD) (which was used by Goodworth [22]), and the linear state combination with time delay (LSCTD) controller. The other two are nonlinear: the neural network (NN) controller and the neural network with time delay (NNTD) controller.

The linear controllers were mathematically described as:

Proportional-derivative controller (six parameters)
$[T1(t)T2(t)]=[K11B110000K22B22][r1−θ1(t)0−θ˙1(t)r2−θ2(t)0−θ˙2(t)]$
(1)
Full-state proportional-derivative controller (ten parameters)
$[T1(t)T2(t)]=[K11B11K21B21K12B12K22B22][r1−θ1(t)0−θ˙1(t)r2−θ2(t)0−θ˙2(t)]$
(2)
Full-state proportional-derivative feedback with time delay controller (14 parameters)
$[T1(t+td1)T2(t+td2)]=[Kpas100Kpas2][r1−θ1(t+td1)r2−θ2(t+td2)]+[K11B11K21B21K12B12K22B22][r1−θ1(t)0−θ˙1(t)r2−θ2(t)0−θ˙2(t)]$
(3)
Linear state combination with time delay controller (34 parameters)
$[T1(t)T2(t)]=∑i=0D([K11iB11iK21iB21iK12iB12iK22iB22i][r1−θ1(t−i·Δt)0−θ˙1(t−i·Δt)r2−θ2(t−i·Δt)0−θ˙2(t−i·Δt)])$
(4)

where $T1(t)$ and $T2(t)$ are the ankle and hip joint torques; $θ1(t)$ and $θ2(t)$ are the ankle and hip joint angles; r1 and r2 are the reference angles of the ankle and hip joints; K11 and B11 are the proportional and derivative gains for the ankle joint; K22 and B22 are the proportional and derivative gains for the hip joint; K21 and B21 are the proportional and derivative gains for the ankle joint from the hip joint angle and angular velocity feedback; K12 and B12 are the proportional and derivative gains for the hip joint from the ankle joint angle and angular velocity feedback; $Kpas1$ and $Kpas2$ are the passive proportional gains for the ankle and hip joints, respectively; $td1$ and $td2$ are the delay time for the ankle and hip joints in the FPDTD controller; $θ1(t−i·Δt)$ and $θ2(t−i·Δt)$ are the inputs of the LSCTD controller, representing the ankle and hip joint angles at $i·Δt$ time prior to the current time t; Ki and Bi are the proportional and derivative gains that multiply with the state at $i·Δt$ time prior to the current time t; and D is the number of delayed state feedback paths. The number of time delays D was 3 and they were multiples of $Δt$ = 60 ms.

The NN controller was defined as regular feed-forward neural network [35] with one hidden layer and four hidden nodes. It is a nonlinear controller (policy), since its activation function is the smoothed leaky-ReLU function: $f(x)=x+0.7(x−x2+0.00012)$. Inputs of the NN controller were the states of the closed-loop model (two joint angles and two angular velocities) and a constant node (unit input). The hidden layer included a constant node (unit input) also. Outputs of the NN controller were the two joint torques. This controller had 30 parameters (weights) that were optimized.

The NNTD controller also used the regular feed-forward neural network but with one hidden layer and eight hidden nodes. Inputs of the NNTD controller were the current state information and time-delayed state feedback from 60, 120, and 180 ms prior, as in the LSCTD controller. Outputs of the NNTD controller were the two joint torques. The activation function used was the same as in the NN controller. This controller had 154 parameters (weights) that were optimized.

### 2.5 Trajectory Optimization With Direct Collocation.

The indirect identification approach mentioned above transforms the identification problems into trajectory optimizations where the state trajectories and control parameters were optimized at the same time. The direct collocation method [36,37], rather than the shooting method [38], was used to perform the trajectory optimizations. The trajectories were represented by values at N discrete time points and the Midpoint Euler approximation for state derivative $x˙$ was used to convert the system dynamics into nonlinear constraints. The controller identification problem can then be written as a large-scale nonlinear program (NLP)
$Find:state trajectory: X={x1,x2,…,xN}and control parameters: Pto minimize the objective function: F=∑n=1N(θ1n,m−θ1n)2+(θ2n,m−θ2n)22·Nsubject to:system dynamics:f(xi+1+xi2,xi+1−xih,P,ai+1+ai2)=0(i=1…N−1)and bounds on statet rajectories:xlow≤xi≤xupp (i=1…N)and bounds on control parameters:Plow≤P≤Pupp$
(5)

where h is the time-step (20 ms); xi is the state at node i, comprised of ankle and hip joint angles and angular velocities: $x=[θ1,θ˙1,θ2,θ˙2]$; $θ1n,m,θ2n,m$ are the measured joint angles at time point n; $f(x,x˙,P,a)=0$ represents the dynamics of the closed-loop model, with platform acceleration a(t). Details are provided in Appendix  A. The objective function, objective gradient, dynamic constraints, and constraint gradients were coded in Python (3.6.0). An interior point optimizer (Ipopt) [39] was used to solve the large-scale nonlinear program. Optimizations of the linear controllers were performed on a local laptop. The other optimizations were performed on a Linux cluster in the Ohio Supercomputer Center [40]. The optimization code can be found in this GitHub repository.3

For each identification problem, ten optimizations were performed with random initial guesses of control parameters. The motion data were used as initial guess for the state trajectory. An identification problem was defined as the identification of one type of controller on one trial of experimental data. The result with the lowest objective function was selected as the solution for the identification problem. For the NNTD controller type, only one initial guess was used for each identification problem due to the long solution time. In total, 72 identification problems (12 periods of data and 6 types of controllers) were solved, resulting in 612 optimizations.

For the PD and FPD controllers, stochastic identification [41] was used to help identify practically stable controllers by enforcing the closed-loop model generates practically stable motions with multiple episodes of process noise. Two and three episodes were used in the stochastic identifications for the PD and FPD controllers, respectively.

## 3 Results

### 3.1 Identification Results on Simulated Data.

Figure 3 shows the means and standard deviations of the identified control gains among ten realizations of sensor noise (identified gains in each realization of sensor noise were from the best fit result in ten random initial guess optimizations). The leftmost bar in each subplot indicates the “true” value of the corresponding control gain that was used to generate the simulated motion data. The other five bars from left to right in each subplot indicate the identified control gains from five lengths of simulated data (10, 30, 50, 70, and 90 s). In general, the identified control gains were close to the true values, except for the two passive proportional gains. However, the differences between the true values and identified results were small when adding the passive and active proportional control gains together (bottom two subplots). Another significant result was that with the increase of the data length, the standard deviations of identified gains among different realizations of sensor noise were decreasing.

Fig. 3
Fig. 3
Close modal

Bias errors and variability across ten noise realizations, for all identified control parameters, are shown in Fig. 4. This result was from the identification of 50 s of simulated data. The bias errors of all identified parameters were less than 2.5% of the true values, except the two passive proportional gains. The variations in identified control gains among ten realizations of noise were below 10% of the averaged values, except for the two passive proportional gains.

Fig. 4
Fig. 4
Close modal

### 3.2 Identification Results on Experimental Data.

Figure 5 shows the root-mean-square (RMS) of the measured motions (first two box plots) and the fits (difference) between data and the outputs from the identified closed-loop models with six controller types (other box plots). With increasing controller complexity, the RMS tended to decrease, indicating a better fit. The RMS value of the PD controller was about half of the variation in the data, showing that the PD controller can only explain about 50% of the variance in standing posture. The RMS value of the NNTD controller was about one fourth of the variation in experimental data, showing that the NNTD controller can explain about 75% amplitude of the standing balance motion. Considering there was only one optimization in the NNTD controller identification, whereas the PD controller result was the best of ten optimizations, the RMS of NNTD controller could be even lower. Fits of ankle and hip joint motions of participant 3 are shown in Appendix  B.

Fig. 5
Fig. 5
Close modal

The identified control gains of the PD, FPD, and FPDTD controllers are shown in Figs. 68. Before averaging, gains were normalized to a typical male with weight of 70 kg and height of 1.75 m [42]. Joint torques were normalized by the product of body mass and height; joint angular velocities were normalized by the square root of gravity divided by the body height. In general, the identified feedback control gains were mostly consistent among subjects and two perturbation trials. In all three types of controllers, proportional control gains K11 ($K11+Kpas1$) of the ankle joint had larger values than the hip joint's proportional gains K22 ($K22+Kpas2$). The $K21,B21,K12,B12$ in the FPD and FPDTD controller types were significantly different from zero. In the FPDTD controller, the time delay parameter of the ankle joint was smaller than the hip joint's, which is not as expected, as the distance for the nervous signal translation is longer for ankle joint so that larger delay time should exist. The proportional feedback gains K22 ($K22+Kpas2$) of the hip joint had higher values in the second perturbation trial than the first perturbation trial. Parameter values for the NN and NNTD controllers are not shown because no meaningful interpretation is possible.

Fig. 6
Fig. 6
Close modal
Fig. 7
Fig. 7
Close modal
Fig. 8
Fig. 8
Close modal
Fig. 9
Fig. 9
Close modal

## 4 Discussion

### 4.1 Identifications on Simulated Data.

Identifications on simulated data showed that the identification method (trajectory optimization with direct collocation) was able to find correct control parameters. More importantly, the identified control gains had less bias error and variation comparing to Goodworth and Peterka's identification study [22]. One possible reason is that the external stimulus in this study was much larger. The peak value of the stimulus in their study was only 5 N·m. We used a stimulus which the peak value of equivalent torque perturbation reached around 30 N·m for the hip joint and 50 N·m for the ankle joint. In general, large perturbation generates a low noise-to-signal ratio (NSR) which helps find true control parameters. This indicates that the external stimulus used in this study was large enough to collect useful human reaction data from experiments. This is consistent with recent results of Schut et al. [43] who showed that the 10 cm peak-to-peak perturbation amplitude is sufficient to collect low NSR data in a standing balance experiment.

Another possible reason that Goodworth and Peterka's results have larger variation is their frequency-domain approach to parameter optimization. Optimizations were performed to fit the FRF of the closed-loop model to the FRF of the experimental data. It is hard to calculate the analytical derivatives of the FRF with respect to the identified parameters, which limited the optimization performance. In our work, the gradient of the objective function and the Jacobian matrix of dynamic constraints were provided analytically helped the optimizer quickly find good solutions.

### 4.2 Identifications on Experimental Data.

The RMS fit errors for the six controller types suggest that both multiple time delay paths and nonlinear properties are needed to fully explain the human standing balance control under random external perturbations. The FPD controller had lower RMS than the PD controller indicating that cross-joint feedback is important for human standing balance. However, the FPDTD controller did not have a significant lower RMS than the FPD controller, suggesting that the additional time delay component was not sufficient to better explain the CNS. The LSCTD controller had a more generalized time delay structure, which generated better fit with experimental data than the FPDTD controller with one time delay for each torque. The NN controller also had lower RMS than the FPD controller, showing the importance of nonlinearity in the controller. Finally, the NNTD controller had both generalized time delay and nonlinearity, which resulted the lowest RMS among all six types of controllers. This helps explain previous studies [4,44] which found that control gains depended on the amplitude of the ramp perturbations, suggesting nonlinearity.

Feedback control gains identified in this study were similar to those found by Park et al. [4]. The proportional control gains were larger for ankle torque than for the hip torque in all shown controller types, possibly because both need to balance the large trunk mass, which is farther away from the ankle. The nonzero cross-joint control gains confirmed that humans use cross-joint feedback to control standing balance. Cross-joint feedback is likely neural, rather than due to elastic properties of two-joint muscles because there are no muscles that cross both hip and ankle.

Identification results from the second perturbation trial were generally more consistent between participants than those from the first perturbation trial. We observed that participants had a relatively smaller body sway motion in second trial, indicating that participants had adapted and no longer needed to explore to optimize their performance. This learning effect may also be responsible for the larger hip feedback gains in the second perturbation trial.

Due to the long computation time, only one optimization was applied in the NNTD controller identification, possibly resulting in a local minimum. With a more thorough optimization, the fit might even have been better. Neural network controllers have a general nonlinear property which is good for testing the potential of nonlinear controllers. However, it may not be a good controller in practice because the control performance is not predictable beyond the observed range of joint motions. Controller performance could be evaluated by simulating the identified system with other perturbation inputs that were not used for the identification and comparing to measured motions. This was not included in this study because the NNTD controller was just used to demonstrate the ability of the proposed identification method. For practical nonlinear controllers, one might also use deterministic nonlinear structures, such as polynomial functions, in which the identified controller parameters have a physical meaning and their performance is more predictable when operating in different conditions.

The identified cross-joint control gains K12 in the FPD controller type were mostly negative, which suggesting a positive cross-joint feedback. However, the cross-joint feedback gains were all positive in the identified FPDTD controllers. This might be because the FPD controller is too simple to represent the human CNS. We would suggest future studies can examine to what extent simple models can extract unbiased information about the control system.

The identified controllers in the experimental study were all from male participants (the last six). Thus, the control parameters showed in this paper may not represent the overall properties of the standing balance CNS of all humans due to the potential differences across gender. However, this paper is mainly focusing on proposing a new identification method. The standing balance controllers of females can be directly identified through the proposed identification method once the motion data of female participants was collected.

In the experimental study, only the middle period (150–250 s) of the recorded balance data was used for the controller identification. This was done because each trial started with a slow drift in the posture, possible due to an adaptation process. This was needed because we assume constant controller parameters. The initial period, however, may contain important information that is clinically relevant. Proper identification would require modeling the dynamics of the adaptation process, which is conceptually possible, but beyond the scope of this study and left for future work.

The identification approach used in this study does not require assuming a linear model and can extract human postural controllers from a long duration experimental data. As far as we aware, it is the first time that complex nonlinear controllers were identified from 100 s of standing balance data under random perturbations. This method can be applied to more complex controllers and motion tasks, such as identification of balance control during human walking. We would like to encourage others who are trying to understand the human motion control using this identification method to explore complex control structures from complex human motions.

## 5 Conclusions

In this work, we showed that trajectory optimization with the direct collocation method can correctly identify control parameters from synthetic motion data. Using 100 s of standing balance kinematic data, under random perturbations, six types of postural controllers, from simple linear to complex nonlinear, were identified. Identification results suggested that the human CNS, even in the simple standing balance task, has nonlinear properties and multiple time delay paths.

## Acknowledgment

Thanks go to Dr. Jason Moore also, for the fundamental work of the postural controller identification and toolboxes that he developed and contributed: Opty, Sympy, and PyDy.

## Funding Data

• National Science Foundation (Grants Nos. 1344954 and 1544702; Funder ID: 10.13039/100000001).

## Footnotes

### Dynamic Model of the Human Standing Balance

The system model of human standing balance includes two components: human body dynamic model and the postural feedback controller (Fig. 1). The human body was simplified to a two-link pendulum. The standing plate is a movable base where the displacement perturbation was applied. The dynamic equation of the human model is
$[(IL+IT+dL2⋅mL+mT⋅(dT2+2dTlL cos(θ2)+lL2))(IT+dT⋅mT⋅(dT+lT cos(θ2)))(IT+dT mT⋅(dT+lT⋅ cos(θ2)))IT+dT2 ⋅mT][θ¨1θ¨2]+[dT⋅lL⋅mT⋅θ˙12 sin(θ2)−dTlT⋅mT⋅(θ˙1+θ˙2)2⋅ sin(θ2)dT⋅lL⋅mT⋅θ˙12 sin(θ2)]+[−dL⋅mT sin(θ1)−dT⋅mT⋅ sin(θ1+θ2)−lL⋅mT⋅ sin(θ1)−dT⋅mT⋅ sin(θ1+θ2)]⋅g+[−dL⋅mT cos(θ1)−dT⋅mT⋅ cos(θ1+θ2)−lL⋅mT⋅ cos(θ1)−dT⋅mT⋅ cos(θ1+θ2)]⋅a=[τ1τ2]$
(A1)

where θ1 and θ2 represent the ankle and hip joint angles, respectively; lL represents the length of leg; mL and mT represent the masses of leg and trunk; dL and dT represent the center of mass location in leg and trunk; $τ1$ and $τ2$ represent the joint torques at ankle and hip joints; g is the gravity; and a is the acceleration of the standing plate due to external perturbation, defined positive in the posterior direction.

In the identification study, the joint torque $τ1$ and $τ2$ were generated by state feedback controllers
$[τ1τ2]=[f1(θ1,θ2,θ1,θ2,P1)f2(θ1,θ2,θ1,θ2,P2)]$
(A2)

where f1 and f2 are the control equations of the ankle and hip joints; P1and P2 represent the control gains inside the ankle and hip controllers.

By putting the state feedback controller into equation (A1), the closed-loop dynamic model of the human standing balance system can be described in this format
$M(q)⋅q¨+C(q,q˙)+G(q)+D(q,a)−F(q,q˙,P)=0$
(A3)

where q = [θ1; θ2] includes the joint angles of the system; M(q) is the mass matrix; $C(q,q˙)$ represents the Coriolis and centrifugal forces; G(q) is the gravity term; D(q;a) is the joint torque term that caused by mechanical perturbation; and $F(q,q˙,P)$ is the state feedback control term.

In trajectory optimizations, this dynamic equation was implemented as a equality constraint. A more generalized format of the system dynamics can be wrote as
$f(x,x˙,P,a)=0$
(A4)

where $x=[q,q˙]$ represents the system state.

### Motion Fit of the Identified Six Types of Controllers

Motion fits of the identified six types of controllers of participant 3 are shown here (Fig. 9). Top two subplots are the external stimulus signal (standing platform translation). Left side subplots are the ankle joint motion fits. Right side subplots are the hip joint motion fits. It is clear that the fits got better with the increase of the controller complexity.

## References

1.
Kuo
,
A. D.
,
1995
, “
An Optimal Control Model for Analyzing Human Postural Balance
,”
IEEE Trans. Biomed. Eng.
,
42
(
1
), pp.
87
101
.10.1109/10.362914
2.
Oie
,
K. S.
,
Kiemel
,
T.
, and
Jeka
,
J. J.
,
2002
, “
Multisensory Fusion: Simultaneous Re-Weighting of Vision and Touch for the Control of Human Posture
,”
Cognit. Brain Res.
,
14
(
1
), pp.
164
176
.10.1016/S0926-6410(02)00071-X
3.
Peterka
,
R.
,
2002
, “
Sensorimotor Integration in Human Postural Control
,”
J. Neurophysiol.
,
88
(
3
), pp.
1097
1118
.10.1152/jn.2002.88.3.1097
4.
Park
,
S.
,
Horak
,
F. B.
, and
Kuo
,
A. D.
,
2004
, “
Postural Feedback Responses Scale With Biomechanical Constraints in Human Standing
,”
Exp. Brain Res.
,
154
(
4
), pp.
417
427
.10.1007/s00221-003-1674-3
5.
Alexandrov
,
A. V.
,
Frolov
,
A. A.
,
Horak
,
F.
,
Carlson-Kuhta
,
P.
, and
Park
,
S.
,
2005
, “
Feedback Equilibrium Control During Human Standing
,”
Biol. Cybern.
,
93
(
5
), pp.
309
322
.10.1007/s00422-005-0004-1
6.
Maurer
,
C.
,
Mergner
,
T.
, and
Peterka
,
R.
,
2006
, “
Multisensory Control of Human Upright Stance
,”
Exp. Brain Res.
,
171
(
2
), pp.
231
250
.10.1007/s00221-005-0256-y
7.
Torres-Oviedo
,
G.
, and
Ting
,
L. H.
,
2007
, “
Muscle Synergies Characterizing Human Postural Responses
,”
J. Neurophysiol.
,
98
(
4
), pp.
2144
2156
.10.1152/jn.01360.2006
8.
Kiemel
,
T.
,
Elahi
,
A. J.
, and
Jeka
,
J. J.
,
2008
, “
Identification of the Plant for Upright Stance in Humans: Multiple Movement Patterns From a Single Neural Strategy
,”
J. Neurophysiol.
,
100
(
6
), pp.
3394
3406
.10.1152/jn.01272.2007
9.
Mergner
,
T.
,
2010
, “
A Neurological View on Reactive Human Stance Control
,”
Annu. Rev. Control
,
34
(
2
), pp.
177
198
.10.1016/j.arcontrol.2010.08.001
10.
Bingham
,
J. T.
,
Choi
,
J. T.
, and
Ting
,
L. H.
,
2011
, “
Stability in a Frontal Plane Model of Balance Requires Coupled Changes to Postural Configuration and Neural Feedback Control
,”
J. Neurophysiol.
,
106
(
1
), pp.
437
448
.10.1152/jn.00010.2011
11.
Goodworth
,
A. D.
, and
Peterka
,
R. J.
,
2012
, “
Sensorimotor Integration for Multisegmental Frontal Plane Balance Control in Humans
,”
J. Neurophysiol.
,
107
(
1
), pp.
12
28
.10.1152/jn.00670.2010
12.
Boonstra
,
T. A.
,
Schouten
,
A. C.
, and
Van der Kooij
,
H.
,
2013
, “
Identification of the Contribution of the Ankle and Hip Joints to Multi-Segmental Balance Control
,”
J. Neuroeng. Rehabil.
,
10
(
1
), p.
23
.10.1186/1743-0003-10-23
13.
Engelhart
,
D.
,
Schouten
,
A. C.
,
Aarts
,
R. G.
, and
van der Kooij
,
H.
,
2015
, “
Assessment of Multi-Joint Coordination and Adaptation in Standing Balance: A Novel Device and System Identification Technique
,”
IEEE Trans. Neural Syst. Rehabil. Eng.
,
23
(
6
), pp.
973
982
.10.1109/TNSRE.2014.2372172
14.
van der Kooij
,
H.
,
van Asseldonk
,
E.
, and
van der Helm
,
F. C.
,
2005
, “
Comparison of Different Methods to Identify and Quantify Balance Control
,”
J. Neurosci. Methods
,
145
(
1–2
), pp.
175
203
.10.1016/j.jneumeth.2005.01.003
15.
Engelhart
,
D.
,
Boonstra
,
T. A.
,
Aarts
,
R. G.
,
Schouten
,
A. C.
, and
van der Kooij
,
H.
,
2016
, “
Comparison of Closed-Loop System Identification Techniques to Quantify Multi-Joint Human Balance Control
,”
Annu. Rev. Control
,
41
, pp.
58
70
.10.1016/j.arcontrol.2016.04.010
16.
Van Der Kooij
,
H.
, and
Peterka
,
R. J.
,
2011
, “
Non-Linear Stimulus-Response Behavior of the Human Stance Control System is Predicted by Optimization of a System With Sensory and Motor Noise
,”
J. Comput. Neurosci.
,
30
(
3
), pp.
759
778
.10.1007/s10827-010-0291-y
17.
Van Der Kooij
,
H.
, and
De Vlugt
,
E.
,
2007
, “
Postural Responses Evoked by Platform Pertubations Are Dominated by Continuous Feedback
,”
J. Neurophysiol.
,
98
(
2
), pp.
730
743
.10.1152/jn.00457.2006
18.
Pasma
,
J. H.
,
Engelhart
,
D.
,
Maier
,
A. B.
,
Aarts
,
R. G.
,
Van Gerven
,
J. M.
,
Arendzen
,
J. H.
,
Schouten
,
A. C.
,
Meskers
,
C. G.
, and
van der Kooij
,
H.
,
2016
, “
Reliability of System Identification Techniques to Assess Standing Balance in Healthy Elderly
,”
PLoS One
,
11
(
3
), p.
e0151012
.10.1371/journal.pone.0151012
19.
Pasma
,
J. H.
,
Engelhart
,
D.
,
Maier
,
A. B.
,
Schouten
,
A. C.
,
van der Kooij
,
H.
, and
Meskers
,
C. G.
,
2015
, “
Changes in Sensory Reweighting of Proprioceptive Information During Standing Balance With Age and Disease
,”
J. Neurophysiol.
,
114
(
6
), pp.
3220
3233
.10.1152/jn.00414.2015
20.
Pasma
,
J. H.
,
Engelhart
,
D.
,
Maier
,
A. B.
,
Meskers
,
C. G.
,
Aarts
,
R. G.
,
Schouten
,
A. C.
, and
van der Kooij
,
H.
,
2015
, “
Assessing Standing Balance Using Mimo Closed Loop System Identification Techniques
,”
IFAC-PapersOnLine
,
48
(
28
), pp.
1381
1385
.10.1016/j.ifacol.2015.12.325
21.
Pintelon
,
R.
, and
Schoukens
,
J.
,
2012
,
System Identification: A Frequency Domain Approach
,
Wiley
,
New York
.
22.
Goodworth
,
A. D.
, and
Peterka
,
R. J.
,
2018
, “
Identifying Mechanisms of Stance Control: A Single Stimulus Multiple Output Model-Fit Approach
,”
J. Neurosci. Methods
,
296
, pp.
44
56
.10.1016/j.jneumeth.2017.12.015
23.
Torres-Oviedo
,
G.
, and
Ting
,
L. H.
,
2010
, “
Subject-Specific Muscle Synergies in Human Balance Control Are Consistent Across Different Biomechanical Contexts
,”
J. Neurophysiol.
,
103
(
6
), pp.
3084
3098
.10.1152/jn.00960.2009
24.
Press
,
W. H.
,
Teukolsky
,
S. A.
,
Vetterling
,
W. T.
, and
Flannery
,
B. P.
,
2007
,
Numerical Recipes 3rd Edition: The Art of Scientific Computing
,
Cambridge University Press
,
Cambridge, England
.
25.
Vyasarayani
,
C. P.
,
Uchida
,
T.
,
Carvalho
,
A.
, and
McPhee
,
J.
,
2011
, “
Parameter Identification in Dynamic Systems Using the Homotopy Optimization Approach
,”
Multibody Syst. Dyn.
,
26
(
4
), pp.
411
424
.10.1007/s11044-011-9260-0
26.
Ackermann
,
M.
, and
van den Bogert
,
A. J.
,
2010
, “
Optimality Principles for Model-Based Prediction of Human Gait
,”
J. Biomech.
,
43
(
6
), pp.
1055
1060
.10.1016/j.jbiomech.2009.12.012
27.
Koelewijn
,
A. D.
, and
van den Bogert
,
A. J.
,
2016
, “
Joint Contact Forces Can Be Reduced by Improving Joint Moment Symmetry in Below-Knee Amputee Gait Simulations
,”
Gait Posture
,
49
, pp.
219
225
.10.1016/j.gaitpost.2016.07.007
28.
Dembia
,
C. L.
,
Bianco
,
N. A.
,
Falisse
,
A.
,
Hicks
,
J. L.
, and
Delp
,
S. L.
,
2019
, “
Opensim Moco: Musculoskeletal Optimal Control
,”
BioRxiv
.10.1101/839381
29.
Horak
,
F. B.
,
1987
, “
Clinical Measurement of Postural Control in Adults
,”
Phys. Therapy
,
67
(
12
), pp.
1881
1885
.10.1093/ptj/67.12.1881
30.
Moore
,
J.
, and
van den Begert
,
A.
,
2015
, “
Human Standing Controller Parameter Identification With Direct Collocation
,”
15th International Symposium on Computer Simulation in Biomechanics (ISB)
, Edinburgh, UK, July 9–11.
31.
Koelewijn
,
A. D.
,
Heinrich
,
D.
, and
van den Bogert
,
A. J.
,
2019
, “
Metabolic Cost Calculations of Gait Using Musculoskeletal Energy Models, a Comparison Study
,”
PLoS One
,
14
(
9
), p.
e0222037
.10.1371/journal.pone.0222037
32.
Winter
,
D. A.
,
2009
,
Biomechanics and Motor Control of Human Movement
,
Wiley
,
New York
.
33.
Wang
,
H.
, and
van den Bogert
,
A.
,
2020
, “
Standing Balance Experiment With Long Duration Random Pulses Perturbation
,”
Zenodo
.10.5281/zenodo.3631958
34.
Crank
,
J.
, and
Nicolson
,
P.
,
1947
, “
A Practical Method for Numerical Evaluation of Solutions of Partial Differential Equations of the Heat-Conduction Type
,”
,
6
, pp.
207
226
.
35.
Hagan
,
M. T.
,
Demuth
,
H. B.
, and
Beale
,
M.
,
1996
,
Neural Network Design
,
PWS Publishing
,
Boston, MA
.
36.
Hargraves
,
C. R.
, and
Paris
,
S. W.
,
1987
, “
Direct Trajectory Optimization Using Nonlinear Programming and Collocation
,”
J. Guid., Control, Dyn.
,
10
(
4
), pp.
338
342
.10.2514/3.20223
37.
Kelly
,
M.
,
2017
, “
An Introduction to Trajectory Optimization: How to Do Your Own Direct Collocation
,”
SIAM Rev.
,
59
(
4
), pp.
849
904
.10.1137/16M1062569
38.
Betts
,
J. T.
,
1998
, “
Survey of Numerical Methods for Trajectory Optimization
,”
J. Guid., Control, Dyn.
,
21
(
2
), pp.
193
207
.10.2514/2.4231
39.
Wächter
,
A.
, and
Biegler
,
L. T.
,
2006
, “
On the Implementation of an Interior-Point Filter Line-Search Algorithm for Large-Scale Nonlinear Programming
,”
Math. Program.
,
106
(
1
), pp.
25
57
.10.1007/s10107-004-0559-y
40.
Ohio Supercomputer Center,
1987
, “Ohio Supercomputer Center,” Ohio Supercomputer Center, Columbus, OH, http://osc.edu/ark:/19495/f5s1ph73
41.
Wang
,
H.
, and
van den Bogert
,
A. J.
,
2020
, “
Identification of the Human Postural Control System Through Stochastic Trajectory Optimization
,”
J. Neurosci. Methods
,
334
, p.
108580
.10.1016/j.jneumeth.2020.108580
42.
Fryar
,
C. D.
,
Kruszan-Moran
,
D.
,
Gu
,
Q.
, and
Ogden
,
C. L.
,
2018
, “
Mean Body Weight, Weight, Waist Circumference, and Body Mass Index Among Adults: United States, 1999–2000 Through 2015–2016
,” National Health Statistics Reports, No. 122.
43.
Schut
,
I.
,
Pasma
,
J.
,
Mestdagh
,
J. D. V.
,
Van Der Kooij
,
H.
, and
Schouten
,
A.
,
2019
, “
Effect of Amplitude and Number of Repetitions of the Perturbation on System Identification of Human Balance Control During Stance
,”
IEEE Trans. Neural Syst. Rehabil. Eng.
,
27
(
12
), pp.
2336
2343
.10.1109/TNSRE.2019.2943206
44.
Welch
,
T. D.
, and
Ting
,
L. H.
,
2009
, “
A Feedback Model Explains the Differential Scaling of Human Postural Responses to Perturbation Acceleration and Velocity
,”
J. Neurophysiol.
,
101
(
6
), pp.
3294
3309
.10.1152/jn.90775.2008