Assessment of the underlying systems involved in standing balance: the additional value of electromyography in system identification and parameter estimation

Background Closed loop system identification (CLSIT) is a method to disentangle the contribution of underlying systems in standing balance. We investigated whether taking into account lower leg muscle activation in CLSIT could improve the reliability and accuracy of estimated parameters identifying the underlying systems. Methods Standing balance behaviour of 20 healthy young participants was measured using continuous rotations of the support surface (SS). The dynamic balance behaviour obtained with CLSIT was expressed by sensitivity functions of the ankle torque, body sway and muscle activation of the lower legs to the SS rotation. Balance control models, 1) without activation dynamics, 2) with activation dynamics and 3) with activation dynamics and acceleration feedback, were fitted on the data of all possible combinations of the 3 sensitivity functions. The reliability of the estimated model parameters was represented by the mean relative standard errors of the mean (mSEM) of the estimated parameters, expressed for the basic parameters, the activation dynamics parameters and the acceleration feedback parameter. To investigate the accuracy, a model validation study was performed using simulated data obtained with a comprehensive balance control model. The accuracy of the estimated model parameters was described by the mean relative difference (mDIFF) between the estimated parameters and original parameters. Results The experimental data showed a low mSEM of the basic parameters, activation dynamics parameters and acceleration feedback parameter by adding muscle activation in combination with activation dynamics and acceleration feedback to the fitted model. From the simulated data, the mDIFF of the basic parameters varied from 22.2–22.4% when estimated using the torque and body sway sensitivity functions. Adding the activation dynamics, acceleration feedback and muscle activation improved mDIFF to 13.1–15.1%. Conclusions Adding the muscle activation in combination with the activation dynamics and acceleration feedback to CLSIT improves the accuracy and reliability of the estimated parameters and gives the possibility to separate the neural time delay, electromechanical delay and the intrinsic and reflexive dynamics. To diagnose impaired balance more specifically, it is recommended to add electromyography (EMG) to body sway (with or without torque) measurements in the assessment of the underlying systems. Electronic supplementary material The online version of this article (10.1186/s12984-017-0299-x) contains supplementary material, which is available to authorized users.

Background Impaired balance is a common complaint in elderly and patients with specific diseases like vestibular disorders, stroke or Parkinson's disease [1][2][3][4][5][6]. To maintain standing balance, several underlying systems interact, such as the nervous system, sensory systems and motor system. With age, specific diseases and medication use, these systems deteriorate and compensate for each other's deteriorations [7][8][9][10], which makes it difficult to detect the underlying cause of impaired balance. To diagnose and intervene impaired balance with targeted interventions it is important to detect the underlying cause of impaired balance. Closed loop system identification (CLSIT) combined with perturbations is a method to distinguish the contribution of underlying systems in standing balance by taken into account the interrelation between the underlying systems and to describe the underlying systems with physiologically meaningful parameters. This gives the possibility to identify the underlying changes in standing balance and therefore to diagnose impaired balance more specifically [11,12].
Several studies used CLSIT to assess standing balance in a variety of patient groups, such as elderly [13][14][15][16], vestibular loss patients [17], Parkinson's disease patients [18][19][20] and stroke patients [21], by describing underlying changes in standing balance with physiologically meaningful parameters. To 'open' the closed loop, sensory and/or mechanical perturbations were applied to disentangle cause and effect and to assess the contribution of the sensory systems (i.e. proprioception, vision and vestibular system), sensory reweighting [13,17], the dynamics of the human body [22,23] or the (inter)limb stabilization [14,[18][19][20][21]. All these studies used measurements of the kinetics and/or kinematics by motion capture systems and force plates to estimate the sensitivity functions to mechanical and/or sensory perturbations, the control dynamics and/or the human body dynamics [24]. Only two of these studies [22,23] used muscle activation as outcome measure. They used muscle activation in combination with sensory perturbations to identify the human body dynamics including activation dynamics (i.e. the mapping from muscle activation to body segment angles). However, adding muscle activation to CLSIT in combination with sensory perturbations could result in more insight in the underlying mechanisms involved in standing balance, as explained below.
Adding muscle activation to CLSIT makes it possible to separate the neural time delay of the neural pathways, i.e. the transmission and synaptic delay, from the time delay due to the activation dynamics, i.e. the electromechanical delay. Previous studies showed that the lumped time delay in standing balance, consisting of both the neural time delay and the electromechanical delay, increases with age and diseases [13,[15][16][17].
However, it is unclear whether this is due to changes in the activation dynamics or to changes in the time delay of the other neural pathways. Furthermore, including muscle activation in CLSIT will give more and better insight in the contribution of the intrinsic (i.e. passive) and reflexive (i.e. active) dynamics in standing balance [25], as the muscle activation only represents the contribution of the reflexive dynamics and does not include the contribution of the intrinsic dynamics. Previous studies identifying other human motion control systems, like stabilization of the trunk and of the ankle, wrist and shoulder joint, already used muscle activation and showed the possibility to distinguish the contribution of the intrinsic and reflexive dynamics in the stabilization of joints [26][27][28][29].
Adding the muscle activation to CLSIT requires changes to the balance control model, which is used for parameter estimation to describe the underlying systems with physiologically meaningful parameters. A commonly used balance control model is the independent channel (IC) model [17,30]. This model consists of sensory pathways for each sensory system, the intrinsic and reflexive dynamics, modelled by a PD controller, and an inverted pendulum. To separate the contribution of the electromechanical delay and the neural time delay, the activation dynamics has to be added to the model, which represent the mapping from muscle activation to corrective torques. Furthermore, acceleration feedback might be required in the balance control model to describe the muscle activation response accurately, as previous studies showed that the response of the muscle activation is explained by position, velocity and acceleration feedback [31,32].
In this study we investigated the additional value of adding muscle activation of the lower legs measured with electromyography (EMG) to the assessment of the underlying systems in standing balance, especially in the assessment of the intrinsic and reflexive contributions and the assessment of the activation dynamics. Furthermore, we investigated whether it is needed to measure both ankle torque and body sway to estimate reliable and accurate parameters, which will indicate whether the experimental set up can be simplified or not. Experimental data of healthy participants were used to investigate the additional value of EMG in identifying the underlying control mechanisms represented by model parameters. To check the accuracy of the estimated parameters, a model validation study was performed.

Participants
Twenty healthy young participants (10 men, age 23.6 years (SD 2.9 years)) participated in the study. The participants gave written informed consent prior to participation. The protocol was approved by the medical ethics committee of The Medical Spectrum, Enschede, the Netherlands and was in accordance with the Declaration of Helsinki.
Apparatus and recording Figure 1 shows the experimental set up of this study. A Bilateral Ankle Perturbator (BAP) (Forcelink B.V., Culemborg, the Netherlands) was used to disturb the proprioceptive information of both ankles by applying support surface (SS) rotations around the ankle axis ( Fig. 1a) [33]. The actual angles of rotation and the applied torques to the left and right SS of the BAP were measured with a sample frequency of 1000 Hz and were stored for further analysis (Fig. 1b).
The body kinematics of the lower and upper body were measured in anterior-posterior direction using two drawwire potentiometers (Sentech SP2, Celesco, Chatsworth, CA, United States) by connecting them to the participant's trunk and hip. Together with the motor angles and motor torques, the body dynamics were measured using a Matlab interface with a sample frequency of 1000 Hz and were stored for further analysis (Fig. 1b).
Muscle activity of three lower leg muscles on each side (i.e. the M. gastrocnemius medialis, M. soleus and M. tibialis anterior of both legs) was measured with a wireless EMG system (Cometa Systems, Bareggio, Italy) using pairs of self-adhesive Ag-AgCl surface electrodes, which were placed approximately 2 cm apart and longitudinally on the muscle belly, according to the Seniam guidelines [34]. Together with the motor angles, motor torques and the body dynamics, the EMG was synchronously recorded using Vicon software (Vicon Motion Systems, Oxford, UK) with a sample frequency of 2000 Hz, subsequently down sampled to 1000 Hz and stored for further analysis.

Procedure
During all experiments the participants stood on the BAP with stocking feet. The participants were instructed to stand with their arms in front of the chest and to keep both feet on the support surface. The support surface rotated following a continuous perturbation signal with an increasing perturbation amplitude over the trials, 0.5, 1, 2, 4 and 8 degrees peak-to-peak, both with eyes open and eyes closed, resulting in 10 trials. Each trial lasted 2 min in which the perturbation signal was 6 times repeated. Before each trial the participant was given sufficient time to get accustomed to the perturbation. The participants wore a safety harness to prevent falling, which did not constrain normal body sway and did not provide support or orientation information.

Perturbation signal
A pseudorandom ternary sequence (PRTS) of numbers was designed and used as SS angular velocity. Integration of the velocity signal provided an unpredictable perturbation signal of the SS rotation with a wide spectral bandwidth ( Fig. 1c) [35]. An 80 state PRTS signal with a time increment of 0.25 s was generated, resulting in a signal with a period time of 20 s. The trial consisted of 6 complete cycles of the perturbation signal resulting in a trial of 2 min.

Preprocessing
Data analysis was performed with Matlab (The Mathworks, Natick, MA, United States). Leg and hip angles were calculated from potentiometer data resulting in the segment angle of the leg relative to the vertical and the joint angle of the trunk relative to the leg. The body sway was represented by the angle of the Center of Mass (CoM) relative to the vertical, which was calculated using the leg and hip angles and body geometry of individual segments [36]. The data of the motor angles were used to measure the real SS rotation. The ankle torques were obtained from the recorded motor torques.
The EMG signals were combined according to the method described by Kiemel et al. (2008) to obtain one summed EMG signal representing the muscle activation, in which the EMG signals were weighted to optimize the coherence between the perturbation and EMG signal [22]. From each time series of the 6 lower leg muscles the mean was subtracted. The time series were normalized by dividing by their standard deviation computed from all trials for the given subject. Subsequently, each time series was high pass filtered with a bidirectional first-order Butterworth filter with a cut-off frequency of 16 Hz and rectified. To achieve one summed weighted EMG signal for both legs, each time series was multiplied with a weighting factor and summed. The weighting factors were obtained with an optimization function (Matlab function: fmincon) in order to optimize the average coherence between the SS rotation and the summed weighted EMG signal. Average coherence was calculated by averaging the complex coherence (γ SS,EMG ) across 10 conditions and averaging the squared norm of this coherence (|γ SS,EMG | 2 ) across the frequencies in the perturbation signal. The complex coherence was calculated using eq. 1 In which Φ SS,EMG represents the Cross Spectral Density (CSD) of the SS rotation and the EMG of a specific muscle, Φ SS,SS the Power Spectral Density (PSD) of the SS rotation, and Φ EMG,EMG the PSD of the EMG of a specific muscle.
The weighting factors were constraint to be positive for dorsal flexors (i.e. left and right M. tibialis anterior) and negative for plantar flexors (i.e. left and right M. soleus and M. gastrocnemius medialis). In addition, the sum of the absolute weighting factors was constraint to be unity.
The time series of the body sway, SS rotation, ankle torque and muscle activation were segmented into six data blocks of 20 s (i.e. the length of the perturbation signal). Figure 1c shows a typical example of the mean and standard deviation of the SS rotation, body sway, ankle torque and muscle activation. Figure 2 shows the flowchart of the data analysis to identify the mechanisms underlying balance behaviour and to investigate the additional value of the muscle activation in the assessment of the underlying systems in standing balance. First, the data obtained from the human experiment were transformed to sensitivity functions using CLSIT. Next, model parameters describing the underlying systems were estimated by fitting a model on the sensitivity functions, in which three models were used, 1) without activation dynamics, 2) with activation dynamics and 3) with activation dynamics and acceleration feedback. The first model was fitted on the sensitivity functions of the body sway, ankle torque or a combination, resulting in 3 parameter sets. The second model was fitted on all possible combinations of sensitivity functions (i.e. muscle activation, body sway and ankle torque), resulting in 6 parameter sets. The third model was fitted on all combination of sensitivity functions including at least the muscle activation, resulting in 3 parameter sets. The parameter sets were used to find the combination of model and sensitivity functions with the most reliable estimated parameters.

Sensitivity functions
To describe and obtain a non-parametrical description of the human balance control, the sensitivity function of the output of the human balance control (e.g. body sway, ankle torque, muscle activation) to the perturbation was obtained by estimating Frequency Response Functions (FRFs). Therefore, the perturbation, ankle torque, body sway, and muscle activation were transformed to the frequency domain. The PSDs and CSDs were computed to calculate the FRFs [25]. Only the excited frequencies were analysed. The FRFs were estimated using the indirect approach [25]: In which Φ SS,x represents the CSD of the SS rotation and x, which represents the ankle torque (T), body sway (BS) or muscle activation (MA), and Φ SS,SS the PSD of the SS rotation. The FRF magnitude and the FRF phase represent the amplitude ratio and the relative delay, respectively, between the SS rotation and the ankle torque, body sway or muscle activation.
The individual FRFs of the experimental data were averaged across 20 participants resulting in one FRF per condition. The mean FRFs were used for further analysis.

Parameter estimation
To give physiological meaning to the sensitivity functions, a human balance control model ( Fig. 3) was used to describe the sensitivity functions. The balance control model was based on a previously described balance control model [13,17] and consisted of a single inverted pendulum controlled by a neuromuscular controller using position and velocity feedback of the sensory systems, consisting of a time delay and neural controller. Dependent on the model that was used for fitting, the activation dynamics and the acceleration feedback were added to the model. Equation 3 shows the theoretical description of the sensitivity functions, i.e. the transfer functions.
This model consists of several parameters describing the behaviour of the system. The parameters describing the sensitivity functions were estimated using the mathematical transfer functions of the balance control model (see Appendix). In this case, we are searching for model parameters such that the behaviour of the model matches with the experimentally measured behaviour. Table 1 summarizes which model parameters were estimated and fixed during fitting. The mass and length were measured and used to calculate the CoM height defined as the length between the ankle axis and the CoM. The moment of inertia was calculated using the mass and length according to Winter (1990) [36]. The mass, CoM height and moment of inertia were used as TIME SERIES . All combinations of the sensitivity functions were used in model fits to obtain estimated parameters. Three models were fitted, 1) without activation dynamics, 2) with activation dynamics and 3) with activation dynamics and acceleration feedback. The resulting 12 parameter sets were used to investigate the accuracy and the reliability fixed parameters. The parameters of the force feedback were fixed at a constant value to reduce interaction with the other estimated parameters.
The model was fitted on the mean experimental FRFs (averaged across 20 participants) using a nonlinear leastsquare fit (Matlab function: lsqnonlin) by minimizing the sum squared error (E): In which γ SS,x 2 represents the coherence between SS rotation and ankle torque, body sway or muscle activation, S exp (f ) the experimental or simulated sensitivity function and S est (f,p) the estimated sensitivity function based on the estimated model parameters and N the number of estimated data points.
The coherence varies between 0 and 1, in which a coherence close to one reflects a good signal to noise ratio. By fitting the model on the experimental data, the gain of the EMG, the activation dynamics and the intrinsic dynamics were kept constant over conditions, i.e. they were condition-independent. All other parameters, i.e. reflexive dynamics, time delay and proprioceptive weight, varied over conditions and were therefore condition-dependent.

Reliability
To evaluate the goodness of the model fits, first the Goodness of Fit (GOF) in the frequency domain was calculated for each sensitivity function using eq. 5.
In which S est (f k ,p) represents the estimated sensitivity function per frequency and parameter set and S exp (f k ) the experimental sensitivity function per frequency.
Second, the Akaike information criterion (AIC) was calculated using eq. 6, which represents the relative quality of the different models and therefore the consequence of adding data and parameters to the parameter estimation. In other words, it shows whether adding parameters or data points to the analysis has additional value by increasing the quality of the model. The model  Fig. 3 Model of the human balance control. The human body is modelled by an inverted pendulum controlled by the neuromuscular controller, consisting of a passive and active part. The passive part consists of the intrinsic dynamics representing the intrinsic properties of the muscles and produces a passive torque (T p ). The active part consists of the neural controller with a time delay and activation dynamics producing an active torque (T a ). The neural controller receives feedback from the sensory systems by force feedback (i.e. Golgi tendon organs) and position, velocity and acceleration feedback (i.e. proprioception, vision and vestibular system). The activation dynamics are added to the model simulation and in two of the three models used for parameter estimation does not necessarily explain more variance of the experimental data. A lower AIC indicates a higher quality of the model.
In which E is the summed squared error, N the number of estimated data points (i.e. the number of frequency points) and d the number of estimated parameters.
Third, the standard error of the mean (SEM) was calculated per parameter to describe the reliability of the parameters. A low SEM means a more precise estimation and therefore a low minimal detectable change. This indicates that the parameters are more sensitive for detecting differences between groups or changes within groups. The SEM was calculated using the diagonal of the estimated covariance matrixP [37]: In which J is the Jacobian (matrix of partial derivatives of the prediction error (ɛ, equation 4) to each parameter) and E the summed squared error. A high partial derivative corresponds with a low SEM; a change in the estimation results in a high change of the prediction error, which allows a more precise estimation.
The mean SEM (mSEM) was calculated by averaging the relative SEM values of the parameters, in which the relative SEM is obtained by dividing the SEM by the estimated value. To investigate the effect of adding the muscle activation and the acceleration feedback to the model on mSEM, mSEM was calculated for the activation dynamics parameters (mSEM ACT; relative damping and eigenfrequency), acceleration feedback parameter (mSEM K a ) and for the other parameters (mSEM basic; intrinsic dynamics, reflexive dynamics, time delay and proprioceptive weight) separately.

Model validation
A model validation study was performed to investigate the accuracy of the estimated parameters for all combinations of sensitivity functions and models by applying the method on simulated data using estimated model parameters obtained in the experiments. Therefore, one dataset was simulated using the same balance control model as used for the parameter estimation including the activation dynamics and acceleration feedback, which was implemented in Matlab (The Mathworks, Natick, USA) ( Fig. 3 and Appendix). Pink noise was added according to Van der Kooij and Peterka [38] to mimic sensory and motor noise and simulated with Simulink. The parameters were set to specific values found in the human experiment of the condition with 2 degrees peak-to-peak perturbation amplitude and eyes open.
The model was perturbed with a sensory perturbation of the proprioceptive information by a continuous support surface rotation of 2 degrees peak-to-peak amplitude as described in detail in 'Perturbation signal'. The simulation time was 2 min in which the perturbation signal was 6 times repeated. Time series of the support surface rotation, ankle torque, body sway and muscle activation were sampled at 1000 Hz and used for further analysis as explained in 'Data analysis' resulting in sensitivity functions and estimated model parameters.

Accuracy and reliability
A comparison was made between the estimated values and the values used in the simulation represented by the relative difference per parameter. These differences were averaged representing the mean relative difference (mDIFF) indicating the accuracy. mDIFF was calculated for the activation dynamics parameters (mDIFF ACT), the acceleration feedback parameter (mDIFF K a ) and for the other parameters (mDIFF basic), separately.
To assess the reliability, the same measures were calculated as for the experimental data (i.e. the GOF, AIC, mSEM basic, mSEM ACT and mSEM K a ).

Results
A human experiment was performed to show the effect on the estimated parameters by adding muscle activation and therefore activation dynamics and acceleration feedback to CLSIT. Next, a model validation study was performed to investigate the accuracy. Below, the results of both studies are described. Figure 4 gives an overview of the goodness of parameter estimation on experimental data for all 10 conditions (0.5, 1, 2, 4 and 8 degrees peak-to-peak amplitude with eyes open and eyes closed) in terms of mSEM, AIC and GOF. Adding the activation dynamics to the fitted model results in a higher mSEM of the basic parameters compared with a model without activation dynamics, especially when the torque sensitivity function is used in combination or not with the body sway sensitivity function. Adding also the sensitivity function of the muscle activation to the fitting procedure results in a lower mSEM of the basic parameters, except when only the body sway sensitivity function is used. Adding the acceleration feedback in the fitted model resulted in a low mSEM of the basic parameters, which is comparable with the mSEM without activation dynamics. The mSEM of the activation dynamics parameters is only low for the torque sensitivity function or combined with the body sway sensitivity function. In case also the acceleration feedback is added, the mSEM of the activation dynamics is only low when the body sway sensitivity function or combined with the torque sensitivity function is used. The relative SEM of the acceleration feedback itself was also low in case the body sway sensitivity function was used or in combination with the torque sensitivity function (< 0.14).

Human experiment
Adding the activation dynamics did not affect the GOF compared with parameter estimation without activation dynamics. Adding the sensitivity function of the muscle activation affected the GOF of the sensitivity function of the body sway and of the ankle torque. Also a low GOF of the muscle activation was found. Adding the acceleration feedback results again in a high GOF for all sensitivity functions. These results are also shown in Fig. 5; the fit on the sensitivity functions are comparable for all combinations of sensitivity functions and models, with the best fit in case the acceleration feedback is added to the model.
Adding the activation dynamics did not influence the AIC of the fits substantially. Adding also the muscle activation resulted in a lower AIC, which decreased even more with adding the acceleration feedback. Figure 6 shows the results of the parameter estimation of one condition (2 degrees peak-to-peak amplitude with eyes open) with adding the activation dynamics, the muscle activation and the acceleration feedback to CLSIT. Additional files 1 and 2 show the results of all conditions. Adding the activation dynamics in the fitted model results in comparable estimates as the estimates without the activation dynamics. However, the time delay, the proprioceptive weight and intrinsic stiffness show a high SEM with the addition of activation dynamics.
Adding the muscle activation results in a lower estimate of the intrinsic dynamics and time delay and a higher estimate of the proprioceptive weight. Furthermore, a high SEM of the activation dynamics was found. The SEM of the time delay and proprioceptive weight was lower compared with the situation without use of the muscle activation sensitivity function, except when the body sway sensitivity function was combined with the muscle activation sensitivity function. Adding the acceleration feedback shows low SEMs for all estimated parameters, except when only the torque sensitivity  Table 2 shows the goodness of the parameter estimation on the simulated data with each combination of sensitivity functions and the fitted model including the activation dynamics and acceleration feedback or not. A model fit without activation dynamics and without the use of the muscle activation shows a high mDIFF (22.2-22.4%). Adding the activation dynamics to the fitted model increased the mean difference of the basic parameters (23.3-28.0%) and resulted in a high mDIFF of the activation dynamics parameters (i.e. the relative damping and eigenfrequency) (47.6-59.4%).

Model validation
Adding also the muscle activation to the fitting procedure resulted in a higher mDIFF of the basic parameters (47.4-48.2%) and a smaller mDIFF of the activation dynamics parameters (30.9-31.0%). When also the acceleration feedback was added to the fitted model, a decrease of the mDIFF of both the basic parameters (13.1-15.1%) and the activation dynamics parameters (0.7-0.8%) was found. Also, the mDIFF of the acceleration feedback parameter was small (13.9-15.4%).
Adding the activation dynamics to the fitted model resulted in a higher mSEM of the basic parameters (0.46-1.76) compared with a fitted model without the activation dynamics (0.08-0.11). Adding also the muscle activation to the fitting procedure resulted in an even higher mSEM of the basic parameters (1.94 e 12 -2.19 e 12 ). When also the acceleration feedback was added to the fitted model, again a low mSEM (0.07-0.08) was found, which is comparable with the low mSEM as found before, without muscle activation and activation dynamics. The mSEM of the activation dynamics decreased with adding the muscle activation from 0.26-0.44 to 0.22-0.28 and decreased even more with adding the acceleration feedback to 0.08-0.09. The mSEM of the acceleration feedback parameter was low (0.07).
Adding the activation dynamics in the fitted model did not affect the GOF of the sensitivity functions of the body sway and the ankle torque substantially. However, adding the muscle activation to the fitting procedure resulted in a decrease of the GOF. Also adding the acceleration feedback to the fitted model resulted again in a high GOF of all sensitivity functions.
The AIC of the parameter estimation increased with adding the muscle activation, i.e. adding data points. The AIC decreased again with adding also the acceleration feedback to the fitted model, i.e. adding parameters. Figure 7 shows the comparison between the estimated parameters and the parameters used in the simulation with the corresponding SEM. Without activation dynamics in the fitted model and muscle  activation in the fitting procedure, the intrinsic stiffness and time delay are estimated higher compared to the original value. The reflexive dynamics (i.e. the reflexive stiffness and reflexive damping) and the proprioceptive weight are estimated lower, while the intrinsic damping is comparable with the original parameter. Adding the activation dynamics to the fitted model resulted in even more differences between the estimated and original values, i.e. the time delay, intrinsic damping and activation dynamics parameters were not comparable. Furthermore, the intrinsic dynamics, proprioceptive weight and time delay show a high SEM. Adding also the muscle activation to the fitting procedure results in less difference between the estimates of the activation dynamics and reflexive damping. However, the intrinsic damping and time delay were much lower than the original value. Furthermore, a lower SEM was found for the proprioceptive weight, intrinsic stiffness and time delay. Adding the acceleration feedback to the fitted model resulted in more comparable estimates and also a low SEM for all parameters. The effective stiffness (i.e. the sum of the intrinsic and reflexive stiffness) is comparable with the original effective stiffness as used in the simulation in all situations.

Discussion
In this study we investigated the additional value of adding EMG of the lower legs (representing the muscle activation) and the activation dynamics and acceleration feedback to CLSIT and parameter estimation to assess the contribution of the underlying systems involved in standing balance. With experimental data we showed that both the activation dynamics and acceleration feedback must be added to the fitted model in combination with the muscle activation to the parameter estimation to obtain a reliable estimate of the parameters describing the underlying systems in standing balance. A model validation study confirmed these findings and showed that this resulted in more accurate estimated parameters describing the underlying systems in standing balance.
(See figure on previous page.) Fig. 6 Overview of estimated parameters from the experimental data for each combination of sensitivity functions by adding activation dynamics, muscle activation and acceleration feedback to closed loop system identification. Parameter values are given with standard error of the mean (SEM) of one condition (2 degrees peak-to-peak amplitude with eyes open (EO)). BS: body sway, T: ankle torque, MA: muscle activation, ACT: activation dynamics, K a : acceleration feedback Table 2 Overview of mean relative differences (mDIFF), mean standard error of the mean (mSEM), goodness of fit (GOF) and Akaike information criterium (AIC) for each combination of sensitivity functions and models fitted on simulated data mDIFF (%) mSEM GOF (%) AIC  Adding the activation dynamics to the fitted model is only of additional value, when this addition is accompanied by addition of the muscle activation to the parameter estimation and the acceleration feedback to the fitted model. Adding only the activation dynamics result in a higher mSEM due to the less reliable estimation of the intrinsic stiffness, time delay and proprioceptive weight. This indicates an over parameterization of the model, as the goodness of fit was not influenced. When also the muscle activation was added in the parameter estimation, the mSEM decreases again, which indicates a more reliable estimation of the intrinsic stiffness, time delay and proprioceptive weight. However, a low and physiologically unexplainable value of the time delay was found. As also the GOF shows low values for the sensitivity function of the ankle torque and the muscle activation, which indicates that the model does not describe the experimental data well, we added an acceleration feedback to the fitted model. Previous studies already showed that the muscle activation must be explained by acceleration feedback, next to position and velocity feedback [31,32]. Adding this parameter resulted in a higher GOF, a more realistic time delay and a lower AIC, which indicates that the quality of the model increased.
Adding the activation dynamics and acceleration feedback did not influence the estimates of the proprioceptive weight, which indicates that the proprioceptive weight can reliably be identified with addition of these parameters. The estimated proprioceptive weight is comparable with previous studies [13,17].

Distinction between neural time delay and electromechanical delay
Including the activation dynamics to the fitted model used for parameter estimation, makes it possible to separate the neural time delay and the electromechanical delay due to the activation dynamics. Adding only the activation dynamics to the fitted model resulted in less reliable estimated time delay indicated by a high SEM. Adding also the muscle activation to the fitting procedure resulted in more reliable estimated parameters but an unrealistically low time delay (<0.001 s) compared with previous studies [15][16][17]. Previous studies showed that the effective time delay is around 200 ms [17], of which the time delay of the muscle activation is around 12 ms [39]. As in the current study both the neural time delay and electromechanical time delay were separated by including the activation dynamics, we expected to find a neural time delay around 190 ms. By adding the acceleration feedback in the fitted model, more physiologically realistic and reliable values of the time delay were found.
The activation dynamics estimated with addition of the acceleration feedback are comparable with previous studies in which the activation dynamics is modelled as a second order system [22,32,40].

Distinction between intrinsic and reflexive contributions
It was expected that adding the activation dynamics to the fitted model and the muscle activation to the parameter estimation would result in more accurate and reliable estimates of the intrinsic and reflexive dynamics. The results of the human experiment indeed shows reliable estimates of the intrinsic and reflexive dynamics in all situations, except in case only the activation dynamics are added. However, the model validation study shows that these parameters are only accurately estimated when adding the activation dynamics and acceleration feedback to the fitted model in combination with the muscle activation.
The sum of the estimated reflexive and intrinsic stiffness (i.e. the effective stiffness, Fig. 6) was comparable for all models and combinations of sensitivity functions. This indicates that by using these models and sensitivity functions it is possible to estimate the sum of the reflexive and intrinsic stiffness. However, for separating these parameters, muscle activation needs to be added to the parameter estimation, and acceleration feedback and activation dynamics should be added to the fitted model.
A low value was found for the estimation of the intrinsic stiffness and damping of all combinations of sensitivity functions without activation dynamics. Including the activation dynamics into the fitted model resulted in physiologically more realistic intrinsic dynamics. Previous studies estimating the balance behaviour, which did not add the muscle activation, showed that it is difficult to separate the intrinsic and reflexive dynamics using parameter estimation. Peterka (2002) showed the results of a model with and without the intrinsic properties. He indicated that including the intrinsic properties did not always result in good fits [17]. Furthermore, Pasma et al. (2015) and Engelhart et al. (2016) did not add the intrinsic properties in their models, due to unreliable estimates of the intrinsic properties [13,14]. Also Kiemel et al. (2008), who did add the muscle activation, indicated that the intrinsic damping might to be too small to detect [22]. However, adding the activation dynamics in combination with the muscle activation and the acceleration feedback results in more accurate and physiologically meaningful parameters.

Effect of increasing perturbation amplitude and closing the eyes
Previous studies showed that with increasing perturbation amplitude and closing the eyes sensory reweighting, which is an important aspect of standing balance, can be investigated. By sensory reweighting the contribution of sensory systems changes based on the accuracy of the information that these systems receive. Additional files 1 and 2 show the estimated parameters, obtained with all combinations of sensitivity functions and models, of all conditions with increasing perturbation amplitude performed with eyes open and eyes closed.
The results show that the proprioceptive weight decreased with increasing perturbation amplitude (i.e. sensory reweighting), which is in accordance with previous studies [13,17]. This phenomenon is clearly visible and is not influenced by the addition of the activation dynamics and acceleration feedback to the fitted model and the muscle activation to the parameter estimation. This result suggests that with all used models and combinations of sensitivity functions, sensory reweighting can be quantified.
Also the time delay decreases with increasing perturbation amplitude, which is in accordance with previous studies [17]. An explanation for this might be the difference in time delay between each sensory pathway. With increasing perturbation amplitude, proprioception will contribute less while the other sensory systems contribute more. This will result in a change of the transporting time of the sensory information and therefore the time delay [17].
The intrinsic dynamics were kept constant over conditions, resulting in condition-independent parameters. Keeping these parameters constant resulted in smaller number of estimated parameters. The AIC showed that the quality of the fitted model remains the same with fewer parameters (data not shown).

Experimental set up
In this study the parameter estimation was performed with several combinations of sensitivity functions. The results showed that including only the ankle torque or body sway or a combination of these influences the reliability of the estimated parameters. In case the activation dynamics is combined with the muscle activation, the use of only the body sway shows an unreliable estimate of all parameters, while adding the acceleration feedback makes the use of only the torque less reliable for estimating all parameters. These results indicate that the reliability depends on the fitted model and the included time series used in CLSIT and parameter estimation.

Model validation
The aforementioned results are confirmed by the model validation study. The results of the model validation study show that both the neural time delay and activation dynamics can be estimated reliably and accurately in case both the activation dynamics and acceleration feedback are added to the fitted model in combination with the muscle activation. This also allows for a reliable and accurate distinction between the intrinsic and reflexive dynamics.
The results of the model validation study did not confirm the observation that it depends on the fitted model which sensitivity functions must be added to get reliable estimated parameters. The parameters are estimated with the same reliability and accuracy with all combinations of sensitivity functions.

Recommendations
We recommend to add activation dynamics, acceleration feedback and muscle activation to the fitting procedure to accurately and reliably separate the time delay due to the activation dynamics and the neural time delay, and the intrinsic and reflexive dynamics. We therefore expect that using such a model in clinical populations could help to distinguish the different underlying causes of impaired balance found with age and disease and therefore to specifically diagnose impaired balance and develop and apply targeted interventions to improve standing balance [13,41].
It is recommended to combine EMG measurements with a fitted model including acceleration feedback in the neural controller. This is in contrast with previous studies using system identification and parameter estimation. Welch et al. (2008) showed that the muscle activation during stance must be explained by a position, velocity and acceleration feedback gain. Ignoring the acceleration feedback gain results in a low and a not physiologically realistic time delay [31], which is in agreement with our results.
When no EMG is available it is recommended to exclude the activation dynamics and acceleration feedback from the parameter estimation or give them fixed values obtained from literature, to minimize the effect on the accuracy and reliability of the estimation of the activation dynamics, time delay and proprioceptive weight. On the other hand, it is not recommended to use only muscle activation in CLSIT as this results in unreliable and less accurate estimations of all parameters (data not shown).
The results show low values of the intrinsic properties, which is against expectation [42]. This might be due to the used perturbation signal. The intrinsic dynamics are modelled as a stiffness and damping without a time delay and mainly affects the low frequencies. As the low frequencies are underrepresented, a perturbation signal containing more low frequencies might result in a better estimation of the intrinsic dynamics. Therefore, we recommend to investigate the effect of the perturbation signal in general on the parameter estimation.
We recommend to always consider carefully which model and data you will use to describe human balance behaviour, as the results show that it depends on the fitted model whether the muscle activation must be combined with the body sway, ankle torque or both to obtain the most reliable estimated parameters. To get the most accurate and reliable estimated parameters, it is recommended to combine muscle activation with at least body sway measurements. This means that a motion capture system or another method to measure the body sway must be combined with EMG and ideally with force plates.

Conclusions
In conclusion, adding the activation dynamics and acceleration feedback to the fitted model and the muscle activation to CLSIT improves the accuracy and reliability of the estimated parameters describing the underlying systems in standing balance. This shows that it is possible to separate the muscle activation dynamics from the neural time delay and the intrinsic dynamics from the reflexive dynamics, which gives more insight in the contribution of the underlying systems in standing balance. More detailed information about the underlying systems and therefore the underlying changes with age and disease gives the opportunity to diagnose impaired balance more specifically and improve standing balance with targeted interventions. Therefore, it is recommended to measure EMG in combination with body sway (with or without ankle torque) and add activation dynamics and acceleration feedback to the fitted model used for parameter estimation to assess the underlying systems involved in standing balance and to improve diagnosis of impaired balance.