Estimation Procedures Affect the Center of Pressure Frequency Analysis

Even though frequency analysis of body sway is widely applied in clinical studies, the lack of standardized procedures concerning power spectrum estimation may provide unreliable descriptors. Stabilometric tests were applied to 35 subjects (20-51 years, 54-95 kg, 1.6-1.9 m) and the power spectral density function was estimated for the anterior-posterior center of pressure time series. The median frequency was compared between power spectra estimated according to signal partitioning, sampling rate, test duration, and detrending methods. The median frequency reliability for different test durations was assessed using the intraclass correlation coefficient. When increasing number of segments, shortening test duration or applying linear detrending, the median frequency values increased significantly up to 137%. Even the shortest test duration provided reliable estimates as observed with the intraclass coefficient (0.74-0.89 confidence interval for a single 20-s test). Clinical assessment of balance may benefit from a standardized protocol for center of pressure spectral analysis that provides an adequate relationship between resolution and variance. An algorithm to estimate center of pressure power density spectrum is also proposed.


Introduction
Despite being a simple task, maintaining a quiet standing posture requires the interaction of different stimuli, which are coded by specialized sensory receptors, and the coordination of muscle responses to stabilize center of gravity spontaneous sways (1).Since the center of pressure (COP) under the feet and center of gravity displacements are equivalent at low frequencies, the ability to control body sway is usually assessed using a force plate (2).However, the choice of suitable stabilometric parameters is a crucial matter for pathologic classification and diagnosis and for monitoring rehabilitation.
Several studies focused on the identification of reliable and non-redundant parameters to describe COP time series (3)(4)(5).Chiari et al. (6) outlined the influence of five biomechanical factors on force plate measurements, while Lafond et al. (7) determined the number of trials to increase the intraclass correlation coefficient of averaged stabilometric parameters up to 0.90.According to these investigators (7), more than 20 trials are recommended to reduce the within-subject variability of mean and median frequency of COP, although they did not mention the procedures applied to estimate the COP power spectrum.
Frequency domain analysis seems to be useful for evaluating postural disturbances caused by specific diseases or conditions (8)(9)(10)(11).Patients affected by a chronic tension-type headache exhibit increased COP frequencies compared with a group of healthy subjects (9).On the other hand, children with myelomeningocele show an opposite trend, characterized by reduced median frequency values (10).
Although the COP time series result from a stochastic process, most clinical studies reporting frequency analysis T.M.M. Vieira et al. www.bjournal.com.br of COP sway generally neglect estimation errors of the computed power spectral density (PSD) function.Some experimental designs lack a description of the spectral estimation procedures (12), while others do not take the estimate bias and variance into account (3).
The objective of this study was to evaluate the influence of four different procedures on the estimated COP power density spectrum ( COP ) during quiet standing tasks.Two procedures relate to COP signal acquisition, the sampling rate (SR ) and test duration, while the number of segments (K ) and the use of detrending apply to COP time series.Some guidelines are proposed and a script is provided in Appendix A to estimate the actual COP power density spectrum (S COP ).

Material and Methods
Thirty-five healthy subjects (20 males and 15 females), with a mean age of 23.7 years (range, 20-51 years), weight 57.62 kg (range, 54-95 kg) and height 1.64 m (range, 1.6-1.9m), participated in the experiment.The Ethics Committee for medical research of the Federal University of Rio de Janeiro approved the experimental protocol, and all subjects gave written informed consent before enrollment in the study.Subjects presenting neurological or musculoskeletal impairments that affect equilibrium were excluded from the study.
Prior to the 60 s stabilometric test, all subjects were instructed to stand in the upright position with heels 6 cm apart, 10° of feet abduction, arms relaxed along the body and facing forward.Figure 1 briefly depicts the protocol for positioning the foot.To assess reliability according to different test durations, subjects repeated this task five times with a 3-min rest period.COP time series was measured using a force plate (AccusWay PLUSS , AMTI, USA) with a built-in 16 bits A/D converter and 0-100 Hz bandwidth, and thus sampled at 100 Hz and low-pass filtered using a bidirectional second order Butterworth filter with 12.5 Hz cut-off frequency.Spectral analysis was carried out using custom algorithms developed in Matlab 6.5 (The Mathworks, USA).
Applying the fast Fourier transform (FFT) algorithm to a single and short segment of COP time series erroneously presumes an underlying deterministic process, so that the respective spectrum estimates do not converge due to inconsistency and bias effects (13).Some elementary concepts to provide an unbiased and consistent estimator of the PSD of a stabilometric record will be outlined in the following subsections.

Partitioning COP spectrum
Since body sway results from a stochastic process, spectral analysis is affected by statistical properties of the spectral estimator.The estimator based on the discrete Fourier transform ( COP ) applied to a single segment of COP record is inconsistent, since its variance is equal to the quadratic magnitude of S COP (14).Dividing the COP signal into K non-overlapped segments and averaging the resulting spectra is an alternative to reduce the PSD estimator variance known as periodogram or Bartlett method (13): (Equation 1) where (ƒ) is the PSD estimated for each segment of the signal: (Equation 2) Figure 1. Figure 1. Figure 1. Figure 1. Figure 1.Schematics of foot positioning on the force plate and an example of a short epoch of a typical stabilometric record.The trace on the force plate surface corresponds to the centre of pressure (COP) trajectory on the surface, starting at 0 s (square) and lasting 10 s (circle).The bottom panel illustrates the time series of the same COP trace along the anterior-posterior axis in arbitrary units (a.u.).The window size (L) is the nearest integer defined by the total number of samples (N), L = round(N/K ).As the use of a Hamming window attenuates the boundaries of each segment, the normalization factor U assures that (ƒ) total energy is unitary: where z (l ) represents the applied window function.Note in Equation 1 that for a single segment, COP becomes equal to .
For each frequency the representation of I L (ƒ) in the complex Cartesian plane is composed of a real and imaginary part, except for the real valued 0 Hz (mean value) and SR/2 frequencies, and each sample of COP is given by the addition of both I L (ƒ) squared components.As the real and imaginary components are Gaussian random variables and the segments are independent, which was assured by non-overlapping, each value of the estimated COP power spectrum adheres to a χ 2 distribution function with 2K degrees of freedom (13).Thus, a confidence interval for S COP is estimated: (Equation 4) where LL and UL are, respectively, the lower and upper limits for the interval, obtained directly from the table of χ 2 distribution for 2K degrees of freedom and a significance level α.In this study, the confidence interval to S COP was estimated for two, six and twelve degrees of freedom, according to COP data partitioning into 1, 3, and 6 segments, using α = 0.05 and a linear fit to remove each segment trend before applying FFT.

Sampling strategy
A single COP time series results from the multiplication between the continuous signal by an impulse train with 1/ SR period.As in frequency domain multiplication corresponds to the convolution of both signals, S COP becomes periodic at SR and only frequencies within ±SR/2 are represented.The actual signal energy residing at frequencies beyond these limits artificially increases the estimated power of higher frequencies.This phenomenon, known as "aliasing", is prevented when the Nyquist theorem is applied, assuring that SR ≥2 f max , where f max is the maximum frequency of the signal (14).As stabilometric signals are real valued, the estimated power spectra are symmetrical and only the components between 0 Hz and +SR/2 are considered, enclosing a total of N/2 components.

Test duration
Since the discrete PSD is defined at frequencies given by M(SR / N), where M is an integer between 0 and N/2, the length of the COP time series affects COP spectral resolution and the amplitude of harmonics unrepresented in the estimated spectrum.Spectral resolution can be improved by increasing the number of samples while keeping sampling frequency constant, which means a longer stabilometric record.
In frequency domain, a stabilometric record corresponds to spectral convolution between the actual COP power spectrum and a sine(x )/x function, which produces a biased spectrum estimate due to leakage phenomenon.This leakage is proportional to the width of sine(x )/x central lobe, which shortens as the record length increases (see Ref. 14 for a full description).Therefore, a stabilometric test with increased duration improves spectral resolution and reduces bias error of COP .

Detrending flag
Partitioning COP time series into epochs of smaller duration may produce non-stationary epochs, mostly due to the low frequencies in the actual record (i.e., a 0.1-Hz harmonic would be observed as a trend in epochs of 5 s).This effect is suppressed by removing each segment trend (14), in general given by the fitted first order polynomial.A linear detrending also reduces discontinuities at segment boundaries caused by truncation, eliminating spurious components at higher frequencies, at the cost of reduced power in very low frequency components.

Signal processing
Three estimators were applied to evaluate the number of segments and sampling rate effect on COP spectrum, according to the signal division in K segments (1, 3, and 6) and to the length of the digital COP signal (preserved and decimated to 50 and 20 samples/s), respectively.The influence of the test duration was investigated using the initial duration of 60 s (6000 samples), which was then truncated to 40 s (4000 samples), and 20 s (2000 samples), keeping the same 100 Hz sampling rate.Figure 2 shows a brief overview of these procedures.The mean value and the first order polynomial fit were considered for assessing the detrending effects on COP estimated spectra.
To assess K and SR effects, COP spectra were estimated using the Bartlett method described in the previous section, with a fixed number of 6 segments and linear detrending.Considering test duration and trend removal, COP was estimated using the Welch periodogram method (13), dividing the signal into 7 segments with 50% of overlapping and applying a linear detrend on each segment 1 .Although overlapping increases the variance of the averaged power spectrum estimator, it allows an increased frequency resolution due to longer epoch duration.
The median frequency of anterior-posterior COP displacements (F50 ), which divides COP spectrum into two portions of equal power, was the spectral descriptor used to evaluate the influence of acquisition (sampling rate and test duration) and estimation (data segmentation and detrending) parameters on COP .

Statistical analysis
One-way ANOVA was used to test F50 equality between different estimators (15).By means of the t-test for paired samples, the linear detrending effect on F50 was evaluated.As a consistent estimate is expected from COP estimated with the Welch method, F50 reliability was assessed according to the different test durations of 20, 40, and 60 s, using the intraclass correlation coefficient proposed by Shrout and Fleiss (16).This coefficient is the ratio between the variance due to subjects and total variance evaluated for five trials.
The F50 values across subjects decreased significantly when reducing the number of COP signal divisions, increasing test duration or applying mean detrending (Table 1).As the COP spectrum is mainly composed of low frequencies, F50 was invariant with sampling rate.Longer test durations exhibited high intraclass correlation for F50 values (Table 2), which means increased COP reliability.
COP shape and magnitude did not change with sampling rate strategy (Figure 3A), maintaining F50 value almost constant (P = 0.997).To emphasize this effect, all spectra presented in Figure 3A were depicted within 0-1 Hz, although COP upper limit increased as a function of SR.Note that the first harmonic magnitude and position remains the same independent of the SR used.
The confidence interval of S COP estimates reduced with COP time series fragmentation using an increased number of segments (Figure 3B).A reduction of first harmonic power is also evident, leading to a significant in-      www.bjournal.com.brcrease of median frequency (P < 0.001).When K = 1 segment, the estimated COP spectrum (thicker line) was not observed due to a large confidence interval, which corresponded to 2 COP /UL (= 2040 COP ).Shorter record lengths resulted in power attenuation and peak frequency displacement (Figure 4), which corresponds to the first harmonic of each estimated spectrum.As expected, median frequency presented a statistically significant decrease with test duration (P < 0.001), although increasing intraclass correlation and improving F50 reliability (Table 2).
The linear trend removal greatly reduced COP power below 0.3 Hz (Figure 5), providing median frequency values significantly higher compared with mean subtraction (P < 0.001).Moreover, the trend was not removed from the windowed COP segments when estimating its power spectrum using the "psd" function supplied by Matlab 6.5, which provides a biased COP spectrum.By applying t-tests to compare the F50 values estimated with the Matlab "psd" function and the proposed algorithm, a slight increase (P = 0.16) was observed for the median frequency values estimated with the algorithm and using the mean value for detrending.The F50 significantly increased (P < 0.001) for the spectra estimated with the proposed algorithm and removing the linear trend of each segment (Table 3).

Discussion
This study focused on a technical approach to investigate the validity of COP frequency analysis widely applied to assess postural stability in clinical settings.When using direct methods (periodogram) for estimating the PSD of COP displacements (S COP ), the observed median frequency values were dependent on the number of segments, test duration and the trend removal approach.

Spectral invariability with sampling rate
Corroborating other studies (3,6), COP sway during quiet standing was restricted to low frequencies (<1 Hz) independent of the spectral estimation procedure carried out, likely explaining F50 and COP invariability with the time resolution applied.However, contribution of frequencies up to 2 Hz during quiet standing is reported in the Figure 3. Figure 3. Figure 3. Figure 3. Figure 3. Center of pressure (COP) power spectrum estimated for 20, 50, and 100 Hz sampling rates (SR ) by partitioning the record into six non-overlapped epochs (A).B, COP power spectrum (thick line) and its 95% confidence interval (thin line) estimated according to signal partition into 1, 3, and 6 segments (K ).evidence for contributions above 2 Hz, a conservative minimum sampling rate of 5 Hz should be considered, not only to assure COP time series properties, but also the representation of COP higher frequencies.In addition, an anti-aliasing (low-pass, 2 Hz cut-off) filtering of COP signal should precede any A/D conversion or decimation.

Consistent spectrum estimate and segmentation effects
Although the duration effect of a single segment (15, 30, 60, and 120 s) on COP was previously addressed (3), the estimate variability was not.From the theory of periodograms, it is well known that spectrum consistency increases with signal segmentation, as the confidence interval of S COP becomes proportionally narrow (13,14).On the other hand, when using fixed signal duration, the number of segments and their length are inversely related so that partitioning decreases spectral resolution.Power attenuation of lower COP frequencies and leakage to neighboring frequencies caused by signal partition may overestimate the COP higher frequencies contribution.
Since very low frequencies may represent a linear trend for a short COP epoch, shortening segments duration emphasizes the linear trend removal effect, explaining the significant increase of F50 values when increasing K, in addition to the low resolution produced by an increase in frequency spacing.Resolution and consistence must be considered before test application, since both are inversely dependent for the same test duration.

Resolution, detrending and F50 reliability
As test duration increased COP spectrum compressed toward lower frequencies, so that the F50 calculated with 2000 samples segment was 137% greater than the F50 estimated using non-truncated data.The significant increase of F50 obtained for shorter test duration results merely from spectral estimation method related to resolution and detrending effects.While lower spectral resolution caused a larger spacing between harmonics, affecting the median determination, the linear trend removal applied to assure stationarity before DFT calculation produced a major effect on F50 estimations.When reducing test duration, detrending effects were critical, attenuating the power of COP frequencies below 0.3 Hz (Figure 4).
As expected, detrend played a significant role in the median frequency increase.Instead of an energy leakage as observed for shorter segment estimations in Figure 3B, linear detrending affected only lower frequencies, although enough to increase the relative contribution of higher frequencies to F50.On the other hand, the mean removal eliminates only the 0 Hz component preserving COP shape (Figure 5). Figure 5. Figure 5. Figure 5. Figure 5. Figure 5.Effect of the mean (dashed line) and linear (solid line) detrending flag on the center of pressure (COP) power density spectrum averaged for all estimates.Note the difference between COP spectrum estimated with Matlab "psd" function (top) and with the algorithm proposed in Appendix A (bottom), which removes zero frequency value and provides the density spectrum.
literature (17-----19).Isolating the lower limb proprioceptive reflex by strapping subjects to a vertical support and applying a random disturbance of small amplitude, Fitzpatrick et al. (17) observed reflex responses within 1-2 Hz frequency band, suggesting an operating range for lower limb proprioceptive receptors.Furthermore, for COP spectrum estimated using 20 Hz sampling rate and 12.5 Hz cut-off frequency, any energy within 10.0-12.5 Hz would emerge as an "alias" into the 7.5-10.0Hz band and increase F50 values, which was not confirmed.Considering the lack of Table 3. Table 3. Table 3. Table 3. Table 3. Mean ± SD values of F50 estimated with the proposed algorithm (Custom) and the Matlab "psd" function, removing the mean of linear trend of each segment.The periodograms estimated with 50% overlapping allowed a satisfactory relationship between spectral resolution and consistence, with the first harmonic below 0.1 Hz for a 60-s epoch and small variance as well.The lower intraclass coefficient values of F50 mentioned by Carpenter et al. (3), even to longer test duration (0.45-s to a 60-s test duration), may reflect the higher variance of COP spectrum estimated from a single segment.In the present study, intraclass correlation coefficient values ranged from 0.82 to 0.96 for F50 calculated with signal durations ranging from 20 to 60 s, respectively (Table 2).The estimation method based in overlapping segments confers high withinsubject repeatability even for estimates of short duration (20 s), which is adequate in clinical surveys that avoid longer test durations (i.e., fall risk assessment, rehabilitation monitoring).

Standardization
Even though median frequency is usually evaluated in clinical studies, the lack of standardization and description of the procedures for spectral estimation constrains the comparison of outcomes from different studies.Comparing studies that apply spectral analysis in similar samples of 13 young healthy subjects, the values reported for anterior-posterior median frequency in each study were, respectively, 0.06, 0.6, and 0.18 Hz (8,11,12).Such differences likely emerge from anthropometry and base of support features or spectral estimation procedures.However, Chiari et al. (6) observed that F50 is almost independent of individual biomechanical factors and support base configuration, attributing the between-studies divergence in F50 to estimation procedures.The lower F50 values mentioned by Nardone et al. ( 12) may be derived from the longer test duration (51 s) compared with the other two studies (30 s).Unfortunately, except for Davidson et al. (8) who estimated COP spectra using sinusoidal tapers, the others did not provide enough information concerning estimation procedures.
Meyer et al. (20) observed an increase of F50 values when the tactile sensitivity of forefoot was suppressed but not with complete inhibition of plantar region receptors, suggesting that remaining sensorial mechanisms did not compensate for the forefoot sensitivity loss.Although providing information about S COP estimation, the authors did not account for COP variance.COP power spectrum was estimated partitioning a 30-s stabilometric record into 2 segments with 50% overlapping, leading to unreliable F50 values.
The standardization of a protocol that uses spectral analysis to evaluate postural control must take into account the relationship between spectral resolution and estimator variance.Stabilometric tests of 60 s are recommended to reduce COP variance without constraining low frequency analysis by applying overlapped segmentation.When longer tests must be avoided, 20-s duration provides reliable spectra with a low 0.2-Hz spectral resolution.In such cases, the use of modern estimators based on autoregressive modeling could be advantageous, since these methods are adequate for shorter signals (14).
The procedures carried out for spectral analysis, including signal partition, test duration and trend removal, produced major alterations in amplitude distribution of COP spectrum.Anterior-posterior median frequency increased with COP time series segmentation, shorter test durations and linear detrending.These factors must be considered before applying spectral analysis.For general clinical purposes, applying stabilometric tests of 60 s, using a sampling rate of at least 20 Hz and partitioning the whole record into 6 or 7 segments seems to provide consistent estimates for COP spectrum.Tests of 20 s should be restricted to studies focusing on COP higher frequencies, or when longer tests are not advisable, and providing that frequency resolution is not so relevant.
Biol Res 42(7) 2009 Center of pressure frequency analysis www.bjournal.com.brI L (ƒ) is the Fourier transform of each epoch of COP time series, multiplied by a Hamming window with L samples.

Figure 2 .
Figure 2.Figure 2.Figure 2.Figure 2.Figure 2. Schematic representation of the analysis procedures applied to investigate sampling rate (SR ), number of segments (K ) and test duration (T ) effect on the center of pressure (COP) estimated spectrum.A, Sampling each record at 100, 50 and 20 Hz; B, Signal division according to grid lines as 1, 3, and 6 segments; C, COP time series truncation at 20, 40, and 60 s.
Figure 2.Figure 2.Figure 2.Figure 2.Figure 2. Schematic representation of the analysis procedures applied to investigate sampling rate (SR ), number of segments (K ) and test duration (T ) effect on the center of pressure (COP) estimated spectrum.A, Sampling each record at 100, 50 and 20 Hz; B, Signal division according to grid lines as 1, 3, and 6 segments; C, COP time series truncation at 20, 40, and 60 s.
Figure 2.Figure 2.Figure 2.Figure 2.Figure 2. Schematic representation of the analysis procedures applied to investigate sampling rate (SR ), number of segments (K ) and test duration (T ) effect on the center of pressure (COP) estimated spectrum.A, Sampling each record at 100, 50 and 20 Hz; B, Signal division according to grid lines as 1, 3, and 6 segments; C, COP time series truncation at 20, 40, and 60 s.
Figure 2.Figure 2.Figure 2.Figure 2.Figure 2. Schematic representation of the analysis procedures applied to investigate sampling rate (SR ), number of segments (K ) and test duration (T ) effect on the center of pressure (COP) estimated spectrum.A, Sampling each record at 100, 50 and 20 Hz; B, Signal division according to grid lines as 1, 3, and 6 segments; C, COP time series truncation at 20, 40, and 60 s.

Figure 2 .
Figure 2.Figure 2.Figure 2.Figure 2.Figure 2. Schematic representation of the analysis procedures applied to investigate sampling rate (SR ), number of segments (K ) and test duration (T ) effect on the center of pressure (COP) estimated spectrum.A, Sampling each record at 100, 50 and 20 Hz; B, Signal division according to grid lines as 1, 3, and 6 segments; C, COP time series truncation at 20, 40, and 60 s.

Figure 4 .
Figure 4.Figure 4.Figure 4.Figure 4.Figure 4. Spectral estimates for a single COP record with different durations (full line = 60 s, long dashed line = 40 s, and short dashed line = 20 s).COP = center of pressure.

Figure 4 .
Figure 4.Figure 4.Figure 4.Figure 4.Figure 4. Spectral estimates for a single COP record with different durations (full line = 60 s, long dashed line = 40 s, and short dashed line = 20 s).COP = center of pressure.

Figure 4 .
Figure 4.Figure 4.Figure 4.Figure 4.Figure 4. Spectral estimates for a single COP record with different durations (full line = 60 s, long dashed line = 40 s, and short dashed line = 20 s).COP = center of pressure.

Figure 4 .
Figure 4.Figure 4.Figure 4.Figure 4.Figure 4. Spectral estimates for a single COP record with different durations (full line = 60 s, long dashed line = 40 s, and short dashed line = 20 s).COP = center of pressure.

Figure 4 .
Figure 4.Figure 4.Figure 4.Figure 4.Figure 4. Spectral estimates for a single COP record with different durations (full line = 60 s, long dashed line = 40 s, and short dashed line = 20 s).COP = center of pressure.

Table 2 .Table 2 .Table 2 .Table 2 .Table 2 .
Mean values of F50 for each trial and test duration and the intraclass correlation coefficient (ICC) calculated across trials.F50 = median frequency of anterior-posterior center of pressure displacements.ICC and its 95%CI were applied to five repeated test duration measures.T = time in seconds.