Is Intermittent Control the Source of the Non-Linear Oscillatory Component (0.2–2Hz) in Human Balance Control?

Objective: To explain the 0.2–2Hz oscillation in human balance. Motivation: Oscillation (0.2–2 Hz) in the control signal (ankle moment) is sustained independently of external disturbances and exaggerated in Parkinson's disease. Does resonance or limit cycles in the neurophysiological feedback loop cause this oscillation? We investigate two linear (non-predictive, predictive) and one non-linear (intermittent-predictive) control model (NPC, PC, IPC). Methods: Fourteen healthy participants, strapped to an actuated single segment robot with dynamics of upright standing, used natural haptic-visual feedback and myoelectric control signals from lower leg muscles to maintain balance. An input disturbance applied stepwise changes in external force. A linear time invariant model (ARX) extracted the delayed component of the control signal related linearly to the disturbance, leaving the remaining, larger, oscillatory non-linear component. We optimized model parameters and noise (observation, motor) to replicate concurrently (i) estimated-delay, (ii) time-series of the linear component, and (iii) magnitude-frequency spectrum and transient magnitude response of the non-linear component. Results (mean±S.D., p<0.05): NPC produced estimated delays (0.116±0.03s) significantly lower than experiment (0.145±0.04s). Overall fit (i)–(iii) was (79±7%, 83±7%, 84±6% for NPC, PC, IPC). IPC required little or no noise. Mean frequency of experimental oscillation (0.99±0.16 Hz) correlated trial by trial with closed loop resonant frequency (fres), not limit cycles, nor sampling rate. NPC (0.36±0.08Hz) and PC (0.86±0.4Hz) showed fres significantly lower than IPC (0.98±0.2Hz). Conclusion: Human balance control requires short-term prediction. Significance: IPC mechanisms (prediction error, threshold related sampling, sequential re-initialization of open-loop predictive control) explain resonant gain without uncontrolled oscillation for healthy balance.

H UMAN balance requires external forces to be matched using internally generated muscle forces. In healthy people a well-practiced neurophysiological control system integrates all sensory feedback to regulate activity in the muscles required to keep the center of mass (CoM) within the base of support of their unstable mechanical structure. When standing on a force plate, the control signal maintaining position of the CoM is summarized by the moving point of application of the ground reaction force [1], [2]. Control of sagittal CoM position is generated predominantly by the calf and tibialis anterior (TA) muscles [1]- [3]. Oscillation in the control signal at 0.2-2Hz can be observed, and this oscillation is exaggerated in Parkinson's disease [4], [5]. This oscillation, has been associated with "overgeneration" of corrective control signal and when exaggerated has been explained as an abnormal resonance in the feedback control loop [4], [5]. Our objective is to provide a feedback control explanation of this oscillation (0.2-2Hz) and to provide insight into the neurophysiological processes of human balance.

A. Question -What Is the Source of the 0.2-2 Hz Oscillation?
Postural sway, and also manually controlled systems with dynamics similar to an upright standing human both show oscillation which sustains independently of any external disturbance [6], [7]. Within a feedback control loop there are several possible explanations of this oscillation. Sensory or motor noise provides input which can be colored by a closed loop system to provide a spectrum of sustained oscillation [6], [7]. Resonance within a closed feedback loop can produce amplitude peaks at certain frequencies. That resonance can arise from a combination of delays and poorly tuned parameters, particularly if the controller is not optimal. Also, a non-linearity within the feedback loop such as a threshold, a switch, or event triggered open loop control can cause sustained oscillation without requiring noise as an input [8], [9]. A periodic return to the same state without external input (limit cycle) can occur in some circumstances [8], [9]. Our recent analysis of manually controlled systems with dynamics equivalent to a standing adult, shows that intermittent predictive control (IPC) with aperiodic sampling can explain linear power and non-linear remnant without addition of sensory or motor noise [7]. However, balance is different from manual control, and explaining human balance is challenging. Any convincing model has to reproduce concurrently and adequately the linear response to an external disturbance with accurate physiological delays and timing/phase characteristic, and also the oscillation which is not related linearly to the disturbance.
What are the main candidate models? The best validated, model of human balance is time delayed, continuous state feedback [1], [2], [6], [11]. We represent this currently unsurpassed model of human balance by generic non-predictive control (NPC, Fig. 2(a)). The cerebellum is relevant to balance and provides short term prediction [12]. Hence we consider also predictive state feedback control (PC, Fig. 2(b)). Predictive control requires a mechanism to accommodate prediction error: furthermore, cerebellar function is associated intimately with switching function in linked basal ganglia circuits [12], [13]. Thus we consider intermittent predictive control (IPC, Fig. 2(c)) as a logical extension of short-term predictive control [10],  [7], [10]. (a) Continuous non-predictive control (NPC): uses delayed linear observer state feedback. (b) Continuous predictive control(PC): uses standard linear observer, predictor state feedback. (c) Event-driven intermittent predictive control (IPC): the generalized hold implements a continuous linear observer predictor state feedback controller matched to the existing, underlying system-observer-predictor-state feedback. The hold is operated open-loop and the initial state is reset intermittently. The reset is triggered by the predicted hold state x h deviating from the observed state x o by more than a threshold θ. Message: Compared with PC, IPC uses an additional discrete control loop to sequentially re-initialize a continuous open-loop controller. [14], [15]: IPC represents functionality associated with central cerebellar-basal ganglia networks [16].

B. Overview of Approach and Hypotheses
We implemented a task replicating the essence of human balance while allowing precise measurement of balance control and disturbance rejection (Fig. 1) [17]. We investigated three models NPC, PC and IPC for their potential to reproduce concurrently the delayed linear, and non-linear components of human balance (Fig. 2).
We tested three hypotheses (Fig. 3

A. Apparatus and Balance Task
Standing on a stable surface, participants were strapped to a single degree of freedom actuated device, named Whole Body Mover (WBM) with programmable dynamics replicating standing. Participants used visual-haptic feedback to control forward-backward movement of their own body using natural muscle activity from the calf and TA muscles ( Fig. 1(a)) [17].
The WBM ( Fig. 1(a)) comprised a vertical board rotating around a joint collinear with the ankles, connected to a direct drive linear actuator (XTA3810S, Servotube Actuator, Copley Motion, U.K.) at approximately 1m above the axis of rotation. An incremental position encoder is located in the linear actuator. Using a proportional-integral-derivative (PID) controller, the actuated position of the WBM was controlled to follow the output of a real time simulated dynamic system. To replicate postural balance, the system used the equation of motion Jθ = mgh(1 − c)θ − Bθ + u e + d for an inverted pendulum where θ is forward angle, u e is the experimental control signal generated by the participant, d is the external disturbance, mgh is the gravitational toppling moment per unit angle, c is the passive stiffness relative to mgh, and B is the passive ankle joint viscosity. We use a mass m of 70 kg, gravitational acceleration g of 9.81 ms −2 , and a center of mass height, h of 0.92 m giving mgh = 632 N rad −1 [3]. The moment of inertia J given by kmh 2 where k is a shape factor of 1.3 was 77 kg m 2 . The passive ankle stiffness c was 0.85 and passive ankle viscosity B was 2.9 Nm rad −1 s [3], [18]. Following [19] we used the transfer function 6.9722/ (s 2 + 0.03721 s -1.231) where denominator coefficients determine passive stability and system time constants and the numerator coefficient represents coupling between experimental control signal u e and position θ. For comparison, and to include the circumstance where passive ankle stiffness matches the toppling torque due to gravity (c = 1) we tested also the related marginally stable system 6.9722/ (s 2 +0.03721 s).
Position control of the WBM eliminated the influence of the mechanical structure of the WBM [17]. During these tasks, power in the position signal is confined to low frequencies and using cross correlation the delay between simulated output and measured position of the WBM (4±3 ms, mean ± S.D.) was negligible with respect to physiological delays. If the WBM exceeded a range of motion of ±10°the WBM was deactivated and the trial terminated.
The task was implemented using Simulink, compiled using Real-Time Workshop and executed on a PC using Real-Time Windows Target within MATLAB (Math Works, Natick, MA, USA) with a step time of 1 ms. Hardware signals were interfaced via a data-acquisition card (DAQ card 6036E, National Instruments, USA) at a sample rate of 1kHz to 16-bit precision. All signals were saved at 100 Hz.
The experimental control signal u e was generated using a myoelectric interface recording from the ankle plantar and dorsiflexion muscles. This control signal represents physiologically delayed neuromotor output solely, and excludes passive ankle impedance. Surface electrical activity (sEMG) was recorded from the calf muscles (intersection of gastrocnemius medialis and soleus) and TA of both legs (Trigno, Delsys, USA). The Trigno system contributes a fixed delay of 48 ms to the output of sEMG signals. At the beginning of each experimental session, we recorded sEMG (background electrical noise) in all muscles at rest to set dead-zone values to subtract noise from input to the myoelectric control signal. Throughout the task, sEMG signals were processed in real-time through a high-pass filter (f cutoff = 50 Hz), then rectified and then low pass filtered using a second order filter 1/(1+τ s) 2 of time constant τ = 100 ms. For each participant, the output of each electrode was scaled to a common external force measured using a load cell in series with the linear actuator. With the WBM locked, the participant strapped into the WBM, relaxed one leg at a time and with the other leg increased push (calf) or pull (TA) against the board to match a pre-set visual target equivalent to ∼10 Nm ankle moment. The net plantarflexion (backward acting) control signal was generated by the sum of the two scaled calf muscles minus the sum of the two scaled tibialis anterior muscles. During task familiarization the gain applied to output of the myoelectric interface was adjusted to ensure that it was neither too high nor too low for the participant. We find participants adjust easily to halving or doubling the 'myoelectric gain' with little effect on performance [17].
With the WBM unlocked and using an absolute position potentiometer mounted on the rotational axis, the WBM was set to a reference position of 2°forward with respect to the vertical line, to approximate physiological standing. A small, (∼5 Nm) constant, forwards acting external disturbance was applied requiring a low-level tonic plantar flexion contraction to maintain the reference position (Fig. 1). Using this tonic bias disturbance, the task felt exceptionally natural and similar to normal standing. In the absence of a disturbance challenging balance, participants maintained their position with minimal sway and participants reported the task felt as though they were not doing anything.
To enable identification of the human controller, an external disturbance providing random direction, constant amplitude (∼4 Nm), step changes in input was applied ( Fig. 1(b)). We use a small disturbance relative to normal fluctuation of the control signal so as to study natural balance rather than responses to large perturbations. Thirty-two paired steps were implemented in random sequence from (eight inter-step intervals (0.1, 0.2, 0.3, 0.5, 0.8, 1.4, 2.4, 4.0 s), each interval used four instances, including two unidirectional and two bidirectional pairs. Unidirectional pairs were followed by a double size third step returning to zero. A random recovery of 4-5 s was used after each pair, and after each 3rd return step [20]. The trial duration was 250 s.
The experimental control signal represents muscle activity in units of volts. While a force is associated with active muscle contraction, the forces generated have no effect on motion of the WBM. However, to aid understanding of the task, we provide ( Fig. 1(c)) an estimate of the ankle moment derived from the system equation of motion and the measured relationship between u e in volts and output in degrees. This shows the tonic forward bias (∼5 Nm), and shows the normal step disturbance (∼4 Nm) is less than the oscillation associated with regulating balance. The rectified EMG signals ( Fig. 1(c)) confirm the tonic calf activity and the bursts of TA activity when the control signal reverses sign. At these low levels of muscle contraction, rectified EMG provides a good estimate of the neural drive with minimal non-linearity arising from amplitude cancellation [21].
Participants received natural visual, vestibular and ankle-foot related haptic feedback and also contact sensation from the vertical board to which they were strapped. The position of the board and the reference position were also displayed on a screen mounted at eye level a couple of meters in front of the participant ( Fig. 1(a)).
The experiments reported in this study, conducted at Manchester Metropolitan University (MMU), were approved by the Academic Ethics Committee of the Faculty of Science and Engineering, and conform to the Declaration of Helsinki. Participants gave written, informed consent to the experiment.

B. Participants and Procedures
Fourteen healthy participants (10 male, age 32 ± 12 years, mean ± S.D.), attempted five 250s trials in randomized order including eyes open unstable (EO US), eyes closed unstable (EC US), eyes open stable (EO S), eyes closed stable (EC S), lower amplitude disturbance eyes open unstable (L EO US). Participants were first prepared for sEMG recording, then baseline noise thresholds were recorded, then myoelectric signals scaled as above. Participants were strapped to the WBM and given a familiarization, until they were comfortable with the task, the perturbations and the experimental conditions. If necessary participants returned a second day for the actual experiment. Participants were instructed only to keep the WBM within the range of motion (± 10°). + u e nonlin ). Together these linear and non-linear components comprise 100% of the signal (Fig. 4 Left). The linear portion where the control signal is linearly coherent with the disturbance tends to be largest at low frequency, and at these low frequencies larger than the non-linear portion. The non-linear portion tends to be largest at mid-frequencies of 0.2-2 Hz and increase transiently following changes in input disturbance (Fig. 4 Right). We fitted a non-parametric (high order) linear time-series model to generate the linear portion of the experimental control signal.

C. Overview of Analysis
Stage 2 Model estimation (Fig. 5). We fit parametric control models concurrently to the linear and non-linear portions. We simulated a control signal u e sim , and separate into linear u e sim_lin and nonlinear components u e sim_nonlin as above. We calculated a normalized root mean square error (nrmse) for the linear (nrmse lin ) and non-linear (nrmse nonlin ) components and minimized their sum nrmse lin + nrmse nonlin. Normalized root mean square error calculates the ratio of rms error to rms signal. We used the formula nrmse(u e sim, u e ) = u e sim -u e / u e -mean(u e ) where indicates the 2-norm.

D. Practical Details Stage 1: Extraction of Linear and Non-Linear Components
Analysis was applied to experimental and simulation times series down sampled to 10 ms timestep. A time-series model ARX (timestep 10 ms, autoregressive in u e with exogeneous input disturbance d including a dead-time) was used to extract the component of the human control signal related linearly to the disturbance ( Fig. 4(a)). Using Akaike's Information Criterion (AIC), and mindful of the Trigno delay of 48ms delay and minimal lower limb peripheral spinal feedback delay of 40 ms, the deadtime (n k timesteps) was selected from range 0.09, 0.10 … to 0.4 s using an 8th order model (n a = n b = 8 coefficients). The data was then split into equal training and validation halves. Using 'AIC' and this estimated deadtime, model order was selected from 40 or less. Using this chosen deadtime and model order (n a = n b ), L 2 regularization, a constraint for stability and the full time-series, the ARX model was estimated using prediction error method (one time-step ahead) and used to simulate the linear component. Subtraction of the linear component from the experimental control signal gave the non-linear component (Fig. 4). For data simulated by the parametric control model (below), we used the same procedure to separate into linear and non-linear time series. We used MATLAB 2020a for all analysis: functions here were 'delayest', 'selstruc', 'arx', 'arxRegul'.

E. Practical Details Stage 2: Identification of Model Parameters
Using non-linear optimization of parameters, we fitted three models (NPC, PC, IPC, Fig. 2 For NPC and PC models we optimized 10 parameters (delay t d , 4 state feedback gains k, 5 linear quadratic weightings for observer gain Q o , c.f. Appendix). Since gains k were fitted freely rather than designed using system parameters, mismatch gain g ac was redundant and set to 1. We also optimized 4 amplitudes of gaussian observation noise v y and one amplitude of gaussian motor noise v u ). Our previous analysis of remnant noise spectra showed motor noise and position observation noise exhibit a 1/f rather than gaussian pattern [7], thus we integrated the motor noise v u and observation noise applied to the position state prior to input to the simulation.

F. Practical details: Objective Error Function
From simulated timeseries u e sim , using the same procedure as for experimental timeseries (Section D), we calculated the estimated delay n k_sim , the linear timeseries u e sim_lin , the fractional error in delay %Δ lin = (n k_sim -n k )/n k and the normalized error nrmse lin = nrmse (u e lin, u e sim_lin ). Using the same procedure for experimental and simulated data, we calculated also the nonlinear time series u e sim_nonlin , its absolute value |u e sim_nonlin | and its analytic envelope env sim_nonlin using a Hilbert filter (2s). From u e sim_nonlin we calculated magnitude frequency response |u e sim_nonlin(f) | at 0.01, 0.02 … 50 Hz (Fig. 5) and its normalized error nrmse(|u e sim_nonlin(f) |, |u e nonlin(f) |). From |u e sim_nonlin | and using absolute value of first differential of the disturbance |d(t) -d(t-1)| as input, we calculated a non-linear deadtime n k_sim_nonlin , a fractional error in non-linear deadtime %Δ nonlin = (n k_sim_nonlin -n k_nonlin )/n k_nonlin , transient impulse response |u e sim_nonlin | IRF (Fig. 5) and its normalized error nrmse(|u e sim_nonlin | IRF, |u e nonlin | IRF ). For the analytic envelope env _sim_nonlin we calculate its autocorrelation function (env ACF ) and normalized error nrmse(env ACFsim_nonlin , env ACF_nonlin ). For the system output y sim (i.e., body position) we calculated magnitude frequency spectra and its normalized error nrmse(y sim(f) , y (f) ). We used MATLAB functions 'ar', 'bode', 'impulse', 'envelope', 'xcorr'.
Optimization iteration one minimized nrmse lin + nrmse(|u e sim_nonlin(f) |, |u e nonlin(f) |) alone. Some simulations, showed drift in system output y unlike experimental data. So iteration 2, incorporated nrmse(y sim(f) , y (f) ) within the cost function. Analysis following iteration 2, showed amplitude of non-linear oscillation in the simulated data decaying too much rather than sustaining following intermittent corrections. So iteration 3 included nrmse(env sim_nonlin(f) , env nonlin(f) in the cost function. Analysis following iteration 3 showed partial improvement to the temporal variation in amplitude of nonlinear oscillation, showed experimental amplitude of oscillation was related transiently to step changes in disturbance (

G. Analysis of Resonance and Limit Cycles in Estimated Models
For each optimized model (npc, pc, ic, icn), we report the time constant of the least damped conjugate complex pole, and from its inverse, the resonant frequency of the closed-loop system f res [9] and the closed-loop frequency response [22]. For IPC models (ic, icn) we conducted a sensitivity analysis [9]. This sensitivity analysis reports the maximum eigenvalue of the closed-loop system when sampled discretely at each possible open-loop interval (Fig. 9). Eigenvalues greater than one indicate unstable, sustained oscillation. Eigenvalues of unity indicate periodic limit cycles. We tested frequencies of resonance, periodic limit cycles and also rates of intermittent sampling for correlation with the frequency of peak amplitude of the non-linear component. Finally, for experimental and simulated control signals, we calculated sample entropy (m = 2, r = 0.2 S.D., Chebyshev distance) as a non-linear measures of complexity [23] not included in the optimization cost function.

III. RESULTS
Thirteen participants completed all five 250 s trials including eyes open (EO) v eyes closed (EC), unstable (US) v stable (S), and the reduced amplitude disturbance trial (L). One participant failed to complete two 250 s trials. The myoelectric control signal, delays and sway statistics approximated typical values for natural standing and single segment constrained balance for the sensory conditions studied (Table I) [4], [6]. The control signal deadtime was 145±36 ms overall and higher for EO than EC (F 1,19.6 = 7.2, p = 0.014). Position sway was larger (F 1,14.0 = 84.0, p<<0.00001) and faster (F 1,15.6 = 15.0, p = 0.001) for EC than EO. Linearity (linear fit to the control signal) was 28%±10% overall and smaller for low amplitude disturbance (F 1,13.9 = 30.1, p = 0.00008).

A. Presence of Non-Linear Oscillation in All Conditions (H1)
Combining all participants, the non-linear component of the control signal show a broad amplitude peak in the range 0.2-2 Hz, referred to as the non-linear oscillation (NLO), for all the  conditions studied (Fig. 6). Non-linear power at low frequencies (≤0.

B. Model Fit to Linear and Non-Linear Components (H2)
Fitting scores are summarized in Table II. Cost Function: Overall fit (cost function) of the NPC model was significantly lower than PC and IPC (ic, icn). Main effect of model (F 3,22.4 = 11.5, p = 0.00009), post hoc NPC v PC (F 1,26 = 22.6, p<0.001). There was no significant difference between PC or IPC (ic or icn). Delay: Estimated delays using NPC simulation were significantly lower than experiment (Table I). Main effect of model F 4,35.6 = 6.4, p = 0.0005), post hoc NPC v Experiment (F 1,16 = 24.1, p = 0.001). There was no significant difference in delay between any of PC, IPC (ic or icn) and experiment. The inability to reproduce the experimental delays while fitting concurrently the linear and non-linear response, rules NPC out as a viable explanation of this balance task.
Considering only PC and IPC (icn) as leading, viable candidate models: Linear: Fit to the linear component is significantly higher in IPC (icn) than PC (F 1,14.1 = 5.6, p = 0.03), (Table II). Non-Linear: there was no significant difference in fit to the non-linear component (F 1,12.7 = 1.8, p = 0.20). Delay fit _non-linear : Fit to the delay before onset of transient increase in non-linear amplitude was significantly higher in IPC (icn) than PC (F 1,15.4 = 8.5, p = 0.01), (Table II) (Tables II, I).
IPC use aperiodic sampling with a distribution of open loop intervals ranging from the model delay (Table IV) to 362, 599 and 891ms (25%, 50%, 75% percentiles), (Fig. 8)  intervals above 1.5 seconds, giving limit cycles with periods above 1.5s (Fig. 9 left). The marginally stable external systems had eigenvalues less than unity and thus no limit cycles. For both unstable and marginally stable external systems, IPC shows eigenvalues approaching unity, i.e., an instability boundary, at open-loop intervals less than 1s.
In the closed-loop transfer function T, PC models show resonance, and IPC models showed larger resonance at around 1 Hz (Fig. 9 Middle). NPC showed resonance around 0.36 Hz. The underlying continuous controllers used for the open-loop hold in IPC models (Fig. 9 Middle) T shows an additional resonance at 10 Hz that is not present in the complete IPC model. For IPC, frequency of the resonance in T increased for EO v EC conditions (Fig. 9 right).

IV. DISCUSSION
Using a bespoke task, providing accurate measurement of human balance, we seek a feedback control explanation of the sustained oscillation (0.2-2 Hz) in the control signal. Do limit cycles or resonance in the neurophysiological feedback loop cause this oscillation? We tested two linear continuous models, non-predictive and predictive control (NPC, PC) and one nonlinear model intermittent predictive control (IPC) for their ability to simulate concurrently the delays, the delayed linear response and the non-linear remnant.

A. Key Results
We have clear answers to our three hypotheses (Fig. 3):-H1: The human control signal is mostly non-linear remnant (linearity 28%±10%, Table I). The salient feature of the remnant is the 0.2-2 Hz amplitude peak referred to as the non-linear oscillation (NLO). This NLO is present in all conditions with eyes open or closed, when controlling stable or unstable external systems and when the amplitude of disturbance is small or large (Fig. 6). This NLO is an inherent part of the balance process. H2: PC and IPC can reproduce concurrently the delays, remnant NLO and linear response (Table II, Fig. 7). PC and IPC provide possible explanations of human balance. Non-predictive control (NPC) is not a viable explanation as it cannot reproduce the estimated delays concurrently with linear and non-linear components. H3: Limit cycles are rejected as a general explanation of the NLO, since the NLO was present for stable as well as unstable external systems (Fig. 9). Closed loop resonance, combined with continuous predictive (PC) or intermittent predictive IPC control is a possible explanation of the NLO (Fig. 8, Table III).

B. Do These Results Generalize to Natural Human Balance?
We use a bespoke task to test balance control. Estimation of the human control system is informed by precise measurement of the disturbance, control signal and system output (position) of a known external system, and a known neuromuscular system converting EMG into force. In natural standing the control signal for a multi-segment system is hard to define, the neuromuscular and mechanical system are also hard to define precisely and system output (whole body CoM) is difficult to measure precisely. In natural balance, separation of neuromotor from passive contributions to the control signal is imprecise.
Do our findings apply to natural balance? In our task, balance is natural. Participants use their own muscles and their own natural senses to control movement of their own bodies which are strapped to the board. The actuated board becomes part of their body. The difference from normal standing is that movement of their body is constrained to one degree of freedom namely forward and backward motion around the ankle joints. The disturbance steps are small (∼4 Nm) in relation to oscillation of the control signal associated with balance and the sway statistics, e.g., rms sway, rms speed and frequency of non-linear oscillation (NLO) ( Table I) are similar to natural standing [4], [6].
Both this task and normal postural balance, require participants to engage in the same processes to estimate and generate the muscular forces required to balance the effect of external forces and control body position within a finite range. So the neurophysiological control processes studied here should be representative of natural balance.

C. Short Term Prediction Is Required to Explain Balance and the Non-Linear Oscillation (0.2-2Hz)
NPC is not a viable explanation of this balance task as it could not reproduce the experimental delays concurrently with linear and non-linear components. This result is consistent with known function of the cerebellum as a short term predictor and known role of the cerebellum for balance [12], [24].
Discrimination between predictive and non-predictive control requires a) experimental stimulus-response data accurate and precise to high frequencies, b) a neuromotor response signal uncontaminated by passive non-motor components, c) accurate, precise knowledge of the neuromuscular and external system being controlled, d) control models representing NPC and PC, and e) convincing fit to experimental delays concurrent with reproducing the linear response and all features of the non-linear remnant.
This study used stepwise changes in force disturbance giving high frequency range to the independent stimulus. The measured response (myoelectric signal) was a direct output of the human neuromotor system, and was the signal actually used to control balance. Quality of this stimulus-response data is shown by statistically significant differences in experimental delay (EO v EC) from only 14 participants. We propose the ability to distinguish experimental delays (EO v EC) as one benchmark of data quality for discriminating PC from NPC. We propose ability to fit individual experimental trials without averaging as a second benchmark for data quality. We propose ability to fit delays, linear response and all features of the non-linear remnant as a third benchmark.
These previous datasets of constrained and natural balance may lack stimulus-response data to the accuracy, precision and frequencies required, and lack knowledge of the neuromuscular system and plant to the accuracy required to discriminate NPC from PC. Briefly, these studies compared model simulations to (i) the average (8 subjects) complex frequency response (constrained body position/support surface angle) and remnant (stochastic) magnitude frequency response, at 0.017-1.3 Hz [6]; (ii) the average (18 subjects) complex frequency response alone (e.g., composite EMG signals/estimated leg-trunk segment angles), at 0.025-5 Hz [11]; (iii) to the complex frequency response (ankle torque/platform acceleration) and remnant magnitude frequency response at 0.06 -4.4 Hz [1]. (iv) In the time domain, a composite EMG signal was fitted to a delayed, linear combination of derivatives (0, 1, 2) of estimated CoM position averaged over all push-pull platform accelerations [2] A linear fit and estimated delay alone does not discriminate NPC from PC. Discriminating NPC from PC requires models with accurate values for the external system and accurate values of the position (CoM), control and disturbance signals.

D. Limit Cycles Are Eliminated As a General Explanation of the Non-Linear Oscillation (0.2-2Hz)
Limit cycles are a property of non-linear systems, and have been considered a possible feature of human balance [8], [26], [27]. Limit cycles are a periodic return of the system to the same state, without external input [8]. IPC (Fig. 2(c)) can produce limit cycles when the generalized open-loop hold is mismatched to the underlying closed loop system [9].
When the generalized hold is system matched, the IPC system is stable no matter what the open loop interval [28]. However when the hold is based upon an inaccurate model of the external system (g ac ࣔ 1, limit cycles can occur at certain open loop intervals [9], but only for unstable external systems. Two facts rule out limit cycles as an explanation of the observed NLO. First the NLO occurs in stable as well as unstable systems. Second, sensitivity analysis of the eigenvalues v. open loop interval, showed that even for the unstable systems studied, the period of limit cycles (>1.5 s), is too long to account for the mean <f> and peak f pk power frequencies 0.99±0.2 Hz and 0.87±0.3 Hz (Table I, Fig. 9) of the NLO. The NLO is not a limit cycle.

E. Resonance Explains the Non-Linear Oscillation (0.2-2Hz)
For predictive control models (PC, IPC) the frequency of closed-loop resonance f res correlated trial by trial with the mean power frequency <f> of NLO Fig. 8, Table III). For IPC, frequency of sampling did not correlate with <f> of NLO. Thus closed-loop resonance provides a general, linear or nonlinear explanation of the NLO. For IPC the frequency of resonance f res equals the mean frequency <f> of NLO with no significant systematic error, though for PC f res was significantly lower than <f> (Fig. 8, Tables I, III).
A resonant circuit requires an input to excite an oscillation. For linear models, sensorimotor noise provides a stochastic input. We used observation noise consistent with noisy state-estimation [1], [6], [11]. For IPC, excitation of the NLO occurs deterministically (Fig. 5) with or without sensorimotor noise [7]. IPC includes a sequential process of event triggered sampling, discrete reset of the hold states and continuous open loop implementation of the reinitialized hold [10], [15]. The event trigger, sampling and hold processes (Fig. 2(c)) are not present in linear control

TABLE III RESONANCE
Shows Frequency of Largest Resonance (mean ± SD) for All Models (npc, Pc, Ic, icn) and Conditions Including Eyes Open (EO), Eyes Closed (EC), Marginally Stable (S), Unstable (US) System, and Low Amplitude (L) Disturbance and All Combined (All).
( Fig. 2(a), (b)). Inaccurate prediction, resulting from unpredicted disturbance, inaccurate model of the system or added noise, will trigger sampling when the prediction error exceeds a threshold [7], [15]. Discontinuity in the control signal injects energy into the closed loop at each iteration of event triggered sampling even when IPC is wholly deterministic and noise free. Sequential injection of energy, at the frequency of sampling provides one possible, but rejected, explanation of the frequency of NLO.
Resonance shown by a peak in closed-loop gain T (Fig. 9) indicates a large response at certain timescales in relation to an unpredicted disturbance. Falling over takes approximately ∼1s. A large response on this timescale is valuable for survival. This large response is unhelpful only if it turns into oscillation which cannot be turned off. In neurological conditions (essential tremor, Parkinson's, dystonia) uncontrollable tremor illustrates the undesirability of oscillation that cannot be stopped.
A linear controller contains no mechanism for interrupting resonant oscillation. In linear controllers, closed loop resonance arising from delays and poorly tuned parameters, may indicate that control is not optimal.
IPC provides mechanisms within the main feedback loop for resetting the state of the controller (Fig. 2(c)). When a reset is triggered by prediction error exceeding a threshold, the current state of the hold (magnitude and phase) has no influence on the new state value. The open-loop hold implements a fast underlying continuous controller from the new re-initialized hold state. Evolution of the hold state is not influenced by measured or observed feedback. See how the underlying continuous controller icn ucc has higher closed loop gain at higher frequencies, and additional resonance at 10 Hz, compared with PC, and with NPC (Fig. 9).
We propose, that reinitialization and open loop operation of fast continuous controllers within a main slow discrete control loop, benefits robust control of balance [29]- [31].

F. Neurophysiological Interpretation of PC and IPC Models
Valid explanations must reproduce the observed experimental delay (145±40 ms, Table I). The experimental delay includes a precise "Trigno" delay of 48 ms, which conveniently equals the physiological delay (49.7±7 ms) between onset of EMG and onset of force during voluntary (as opposed to electrically stimulated) contraction [32]. Second order neuromuscular dynamics combined with noise or non-linear variation explain latency longer than the model delay to observe statistically significant responses.
NPC model delays were 85±20 ms (Table IV). However, since NPC could not reproduce experimental delays NPC is eliminated as a valid model.
Using model delays of 108±40 ms (Table IV), PC reproduced the observed delay (145±40 ms, Table I). The PC delay for all unpredicted disturbance is fixed whereas physiological reaction times show a distribution. Furthermore, 108±40 ms is inconsistent with previous, estimates of the delay from other authors (∼88±7 ms) [2] and has the unphysiological attribute of excluding lower limb spinal feedback. However 108±40 ms is reasonable as a mean delay representing balance mediated by trans-cortical feedback loops [30], [33].
IPC reproduced observed delays using model delays of 88±20 ms (Table IV). The IPC delay represents the minimum value of a distribution which includes additional variable time for each unpredicted disturbance to cause prediction error to exceed a threshold. This latency (88±20 ms) and the distribution of open loop intervals (Fig. 8) is consistent with a minimum delay defined by spinal feedback and a main contribution from transcortical and central pathways [30].
While there is no significant difference in cost function fit, IPC fits features of the data significantly better than PC (Linear, Delay fit _non-linear , Sample Entropy, Table II; mean power  frequency Table III, Fig. 8). PC uses noise to reproduce the remnant. PC sample entropy is significantly higher than experiment, thus PC is unphysiological. IPC produces the remnant mechanistically and sample entropy is marginally lower than experiment (Tables I, II). Additional sources of variability (e.g., time varying parameters) are possible, thus IPC is compatible with physiology.
IPC confers a functional benefit. Unlike linear control, IPC can use a resonant response without the disadvantage of uncontrollable oscillation. Open loop implementation of an underlying is the maximum eigenvalue of the periodically sampled discrete system including a hold matching the underlying continuous control system, and a series mismatch gain g ac [9]  continuous controller, within a main discrete feedback loop that can reset the states, allows IPC to use shorter latency, higher gain underlying continuous control than the equivalent linear controller (Fig. 9, icn ucc v PC, NPC).This combination of fast continuous feedback loops within a main slow discrete feedback loop is observed in vertebrate neurophysiology [29], [30]. In all vertebrates, the basal ganglia provides a metabolically costly, and thus important, main function of generalized tonic inhibition of fast trans-cortical and subcortical sensorimotor loops. The basal ganglia loop also provides discrete selection (disinhibition) and reinforcement of beneficial responses. Functionally, the basal ganglia provides a slow pathway within the main sensorimotor feedback loop. Working together the cerebellum and basal ganglia can sequentially switch on or off, and change the gain of transcortical and spinal reflex loops, can potentially shortcut, or not, the main slow loop, and can learn from and model, these reflex loops [12], [13], [24], [30]. In principle, vertebrates contain the neurophysiological machinery for IPC.
The exaggeration and emergence of tremor in neurological conditions (essential tremor, Parkinson's, dystonia), adds weight to the idea that processes of IPC are associated with distributed basal ganglia, cerebellar loops [12], [13], [30]. These IPC processes not present in linear PC or NPC, include short term prediction, detection of prediction error, sequential event triggered sampling, modelling of systems within a generalized hold, reselection of initial hold states, and gating of feedback to allow open loop implementation of the continuous hold. This processes are worth investigating for their potential to explain uncontrolled tremor and other deficits of balance in neurological conditions.

G. Current Limitations and Possible Future Extensions
Our analysis used a single, main feedback loop to represent balance control and tested three models (NPC, PC, IPC). This approach lumps spinal, transcortical and more voluntary feedback into a single feedback controller. Our results showed that experimental sample entropy lies between that predicted by PC and IPC and was closer to IPC. IPC appears close to neurophysiological reality with some amendment required to incorporate continuous control more fully.
For example, spinal, transcortical and voluntary feedback could be combined in parallel. Following c.f. Fig. 1 [16], the external system and neuromuscular system augmented by fast continuous, spinal/transcortical feedback could be modelled as an augmented system controlled by an intermittent controller. IPC feedback could add to fast continuous spinal/transcortical feedback (parallel control) or could provide a setpoint to spinal/transcortical feedback (cascade control).
Evidence suggests short and long latency reflexes are limited in magnitude and are pulsatile (short duration) in nature [34], [35]. This evidence supports a serial concept in which a short duration reflex response is the first response to prediction error arising at short latency relative to the original discrete disturbance. That first response is followed by sequential responses at longer latencies relative to the original disturbance. The IPC model illustrates the serial concept, and already combines continuous with intermittent control (Fig. 2(b)). For IPC, the observer generating estimated states, the predictor generating future states one closed loop delay ahead (e.g., 100 ms), and the open loop hold generating time evolving states from an initial state are all continuous. The hold models the continuous closed loop system that it is matching. In the context of this serial model, spinal/transcortical reflexes are therefore understood as the initial triggered response to an unpredicted disturbance rather than as the response of a separate feedback system.

V. CONCLUSION
We report a whole body balance task, providing data suitable for discriminating non-predictive, predictive, and intermittent predictive models of human balance control. We provide evidence that: 1) the non-linear oscillation (NLO) (0.2-2Hz) present in healthy balance control is explained by closed-loop resonance, not limit cycles 2) short-term prediction is required for human balance. The standard state-estimation, state feedback model of human balance did not replicate concurrently the experimental delays, linear response and oscillatory non-linear remnant.
3) intermittent predictive control (IPC) is a viable explanation of human balance, fitting and explaining the mean frequency of NLO, sample entropy, the linear response, and non-linear delay better than continuous predictive control (PC). 4) processes of IPC (prediction error, threshold related sampling, sequential re-initialization of generalized hold, continuous open-loop implementation of predictive control) enable high gain fast reflexes (resonant closed loop gain) without uncontrolled oscillation. Significance: IPC provides new model based concepts to investigate balance in healthy and neurological conditions.

APPENDIX
Description of models [10] Definition of external and neuro-muscular system. The human participant controls an external 2 nd order system using their net myoelectric signal. The explanatory model controls a dynamic system including the external system and a neuromuscular system representing the generation of the net myoelectric signal from state feedback. The external and neuromuscular system are connected in series (Fig. 2). Fig. 2 shows the external 2nd order system, labelled "system" in Fig. 2(a), (b), and a 2nd order linear approximation of the neuro-muscular system, with a time-constant of 100ms, labelled "neuromusc. system". The system is augmented by a disturbance observer with integral action to compensate for any constant disturbances [10], [15].
Linear controllers (Fig. 2(a), (b)): A standard continuous-time state-space controller [15], [25] containing an optimal state observer together with a state-predictor which compensates the model delay t d and state feedback (not optimal), is used as a linear continuous-time predictive controller (PC) modelling the human operator ( Fig. 2(b)). The linear, continuous nonpredictive controller (NPC), Fig. 2(a)) omits the state predictor. Four measured system outputs, y o (position, velocity, and neuromuscular states), are taken as observer inputs, together with the control signal u. The state observer is designed using standard steady-state linear quadratic methods, which involve minimizing a quadratic cost function of the weighted control signals, system states and output signals [7]. The state feedback gains are adjusted freely by model fitting.
Intermittent predictive controller (Fig. 2(c)): The intermittent controller (IPC) [15] is based on the same structure as the linear continuous time predictive control (PC), but instead of continuous feedback of the observer state, feedback is only used at discrete time points, t i (indicated by the dashed line). The sampled observer state, x o (t i ), is fed to a predictor and subjected to a computational delay, t d which represents the predicted physiological delay, resulting in the predictor state x p (t i ). The predictor state is used as the initial condition for a system matched hold element (labelled as "hold") with dynamics which correspond to those of the equivalent continuous time predictive control loop [15]. Thus, in the absence of disturbances, the hold state x h follows the observer state x o . A gain g ac applied between state feedback and the neuromuscular system simulates mismatch between the assumed-predicted system and actual system. When disturbances or uncertainties affect the loop, x h will diverge from x o , resulting in a non-zero prediction error, e p (t) = x h (t)−x o (t). A quadratic switching function of the form e p T Q t e p > 1, with Q t a positive semi-definite matrix, is defined as an event trigger to reset the hold state x h to the observer state x o [10], [15]. All elements of e p corresponding to the velocity (e p vel ), position (e p pos ), two neuromuscular states, and estimated disturbance are considered, and Q t is a diagonal matrix with five positive elements, θ 1-5 ), forming the axes of an elliptic switching surface. The time between trigger events is the intermittent (open-loop) interval, Δ ol i = t i −t i−1 . A new trigger event can only occur if Δ ol exceeds a minimal open loop interval, Δ ol min > t d + Δ s .

ACKNOWLEDGMENT
With appreciation we acknowledge the contribution of Des Richards to the construction of the apparatus, to all our participants who willingly gave their time, to Martin Lakie for his contribution to this inquiry and also to the reviewers for their substantial contribution to improve the manuscript.