Fractal behavior of the trajectories of the foot centers of pressure during pregnancy

Pregnancy produces significant variations in the gait and postural sway, with increasing falling risk as it progresses. To counteract this risk, pregnant women activate mechanisms in which the forefoot loads are enhanced with respect to the rearfoot ones and this has been indicated as a source of possible musculoskeletal discomfort (Ribeiro et al 2013). The changes in the way the feet rest on the ground produce modifications in the plantar center of pressure (CoP) trajectories and relatively large displacements, with respect to non-pregnant women, were observed in both anteroposterior and medio-lateral directions, mainly in the last trimester of pregnancy (Nyska et al 1997, Ribas and Guirro 2007, Karadag-Saygi et al 2010, Zhang et al 2015, Mei et al 2018).

Plantar CoP trajectories in various situations have been investigated by using instrumented insoles (Winter 1995, Chesnin et al 2000, Martínez-Martí et al 2014, 2016). These devices permit to monitor the movement of both left and right foot plantar CoPs when the subject is either standing or walking. In addition, many others variables related to gait, such as the time support or the position of the foot in swing phase, may be also measured.

The purpose of the present work was to analyze the evolution during pregnancy of the plantar CoP trajectories of pregnant women obtained with instrumented insoles. Instead of the statistical properties of the CoP paths, such as the area enclosing the CoP or the length of the total CoP excursion, already analyzed (Nyska et al 1997, Ribas and Guirro 2007, Karadag-Saygi et al 2010, Ribeiro et al 2013, Zhang et al 2015, Mei et al 2018), we concentrated here on the fractal, non-stationary time characteristics of the CoP trajectories.

To do that, the so-called 'detrended fluctuation analysis' (DFA) technique was used. The DFA method was introduced by Peng et al (1994) and has been considered in many different contexts (Peng et al 1995, Eke et al 2002, Blázquez et al 2009). To characterize the signals analyzed, it uses a scaling exponent that is linked to their fractal behavior.

The procedure we have followed is similar to that considered in previous works to study CoP paths corresponding to the whole body and obtained by means of force platforms (Duarte and Zatsiorsky 2001, Delignières et al 2003, Blázquez et al 2012).

2.1. Instrumented Insoles

In this work, the data corresponding to the plantar CoPs were obtained by using F-Scan^® (TekScan, Boston, USA), a commercial system that includes two instrumented insoles and a datalogger unit that is situated in the lower back of the subject (Tekscan 2018). This system provides dynamic pressure, force and timing information such as the CoP position and trajectories during standing or walking phases, allowing a complete and detailed gait analysis.

The insoles are made of a flexible plastic substrate which contains up to 954 resistive sensing elements that are organized as a set of rows and columns (Tekscan 2018). Their width is 2.5 mm and they are spaced by 5.1 mm, which means a total of 60 rows and 21 columns. The insoles are laminated on both sides with a 0.1 mm flexible protective covering, and the total thickness in sensing region is 0.15 mm. Pressure range goes from 517 to 862 kPa. Insoles are inserted into a sport shoe provided by the research group, and then connected to the datalogger. Data are sent to the computer through a Wi-Fi connection.

The coordinates of the instantaneous position of the plantar CoP were calculated as (Chesnin et al 2000):

$$\begin{eqnarray}&&{x}_{{\rm{COP}}}=\displaystyle \frac{\displaystyle \sum _{i=1}^{n}\,{p}_{i}\,{x}_{i}}{\displaystyle \sum _{i=1}^{n}\,{p}_{i}},\,{y}_{{\rm{COP}}}=\displaystyle \frac{\displaystyle \sum _{i=1}^{n}\,{p}_{i}\,{y}_{i}}{\displaystyle \sum _{i=1}^{n}\,{p}_{i}},\end{eqnarray} \tag{ 1 }$$

where p_i is the pressure measured by the i-th sensor in the insoles and (x_i, y_i) are the coordinates of its center. In our experiment a sampling frequency of 40 Hz was chosen. Given the huge number of sensors per foot in the F-Scan^® system and the definition of the plantar CoP coordinates, we consider that this frequency is enough for this study. Moreover, results can be directly compared with the results of previous works in which data were collected with the same acquisition frequency to study postural sway using a force platform (Duarte and Zatsiorsky 2001, Delignières et al 2003, Blázquez et al 2009, 2010, 2012).

2.2. Characteristics of the subject sample

2.3. The detrended fluctuation analysis method

In the present work, the DFA method has been used to investigate the correlation properties of the plantar CoP trajectories obtained as indicated above. The fractal scaling exponents that characterize the signals in this method are calculated as follows.

Let {ξ} ≡ {ξ_i = ξ(t_i); i = 1, 2, ..., M} be the original series whose data are assumed to be measured with a fixed frequency. The first step is to calculate a new series, {Θ} ≡ {Θ_l; l = 1, 2, ..., M}, whose elements are

$$\begin{eqnarray}&&{{\rm{\Theta }}}_{l}=\displaystyle \sum _{i=1}^{l}\left({\xi }_{i}-\langle \xi \rangle \right).\end{eqnarray} \tag{ 2 }$$

It is simply the accumulated series referred to the average of the M data of the original series, which is given by:

$$\begin{eqnarray}&&\langle \xi \rangle =\displaystyle \frac{1}{M}\,\displaystyle \sum _{i=1}^{M}\,{\xi }_{i}.\end{eqnarray} \tag{ 3 }$$

Now, a set of p subseries labelled as $\{{{\rm{\Theta }}}_{\tau }^{m};\,m=1,2,\ldots ,p\}$ are extracted from the accumulated series {Θ}. These p pieces span the whole series {Θ}, may overlap each other or not, and all of them contain τ consecutive values:

$$\begin{eqnarray}&&\{{{\rm{\Theta }}}_{\tau }^{m}\}\equiv \{{{\rm{\Theta }}}_{{l}_{m}+k};\,k=1,2,\ldots ,\tau \}\,.\end{eqnarray} \tag{ 4 }$$

A polynomial of degree q, ${P}_{q}^{m}$, is fitted to each one the pieces $\{{{\rm{\Theta }}}_{{\tau }}^{{m}}\}$ and the series

$$\begin{eqnarray}&&\{{{\psi }}_{{\tau }}^{{m}}\}\equiv \{{{\psi }}_{{{l}}_{{m}}+{k}};\,{k}=1,2,\ldots ,{\tau }\}\,,\end{eqnarray} \tag{ 5 }$$

is calculated. Here

$$\begin{eqnarray}&&{\psi }_{{l}_{m}+k}\,=\,{{\rm{\Theta }}}_{{l}_{m}+k}\,-\,{P}_{q}^{m}({t}_{{l}_{m}+k})\,,\end{eqnarray} \tag{ 6 }$$

and l_m indicates the initial point of the corresponding subseries that depends on the index m. $\{{\psi }_{\tau }^{m}\}$ is the so-called detrended fluctuation series and is obtained once the local trend of the signal, which is represented by the fitting polynomial, is substracted from the subseries $\{{{\rm{\Theta }}}_{{\tau }}^{{m}}\}$.

Finally, the important quantity is the root mean squared fluctuation of the subseries with size τ:

$$\begin{eqnarray}&&{F}({\tau })\,=\,{\left(\displaystyle \frac{1}{{p}}\displaystyle \sum _{{m}=1}^{{p}}\frac{1}{{\tau }}\sum _{{k}=1}^{{\tau }}{\left[{{\psi }}_{{{l}}_{{m}}+{k}}\right]}^{2}\right)}^{1/2}\,,\end{eqnarray} \tag{ 7 }$$

The original signal is said to show scaling properties if F(τ) ∼ τ^α. In that case, the scaling exponent α determines the fractal characteristics of the signal and its correlation properties. Any completely uncorrelated random series would have α = 0.5, while signals characterized by α-values larger or smaller than 0.5 are positively or negatively correlated, respectively. The first ones are called persistent signals and are characterized by the presence of positive feedback mechanisms in their dynamics. Signals with α < 0.5 are anti-persistent and dynamically ruled by negative control feedback mechanisms (Hurst 1965, Peng et al 1994, Hu et al 2001, Eke et al 2002).

2.4. Specific calculations

According to Hu et al (2001), the value of τ must be conveniently chosen to ensure the feasibility of the DFA method. In the present work we have considered subseries with sizes between τ = 2³ and τ = M/10. We note that the smallest τ value corresponds to 0.2 s, well above the sampling time. On the other hand, we have considered the DFA-1 technique, that is, the fitting polynomial has degree q = 1.

It has been shown that the signals of the CoP of the subject whole body have a multi-scaling character (Duarte and Zatsiorsky 2001, Delignières et al 2003, Blázquez et al 2009, 2010). In this situation the corresponding DFA function F(τ) presents different α-exponents in the small- and large-τ ends, α^S and α^L, respectively. In the present work we aimed at testing if this situation is the same for the plantar CoP trajectories and we calculated these two exponents. To do that the log-log F(τ) versus τ plot was evaluated in 20 τ points, almost equally distributed in $\mathrm{log}\tau$, and straight lines were fitted to the four points closest to each one of the two plot ends. The slope of these fitted lines were the corresponding scaling exponents.

The original signals obtained with the insoles were space-like, that is, they corresponded to the spatial position of the respective plantar CoPs. The scaling exponents of the small- and large-τ ends were labelled ${\pi }_{z}^{{\rm{S}}}$ and ${\pi }_{z}^{{\rm{L}}}$ with z ≡ x or y. As the signal sampling was carried out with constant frequency, it was also possible to obtain velocity-like signals whose elements were the differences between consecutive data of the original series. These signals were also analyzed and the corresponding exponents were labelled as ${\nu }_{z}^{{\rm{S}}}$ and ${\nu }_{z}^{{\rm{L}}}$. The exponents π_z and ν_z are not, however, independent. According to Eke et al (2002) and Hu et al (2001), as the space-like signals are actually the integral of the corresponding velocity-like ones, their scaling exponents must verify the relation:

$$\begin{eqnarray}&&{\pi }_{z}\,=\,{\nu }_{z}\,+\,1\,.\end{eqnarray} \tag{ 8 }$$

In the work by Blázquez et al (2012) it was shown that when using DFA-1 analysis, the feasibility of the estimations of the scaling exponents was ensured if the one uses ${\alpha }_{z}^{{\rm{S}}}\equiv {\nu }_{z}^{{\rm{S}}}$ and ${\alpha }_{z}^{{\rm{L}}}\equiv {\pi }_{z}^{{\rm{L}}}-1$. In our analysis we considered these two exponents.

To have a quantitative estimate of the multi-scaling character of the signals produced by the insoles, we have defined a new quantity,

$$\begin{eqnarray}&&{{\delta }}_{{z}}\,=\,{{\alpha }}_{{z}}^{{\rm{S}}}-{{\alpha }}_{{z}}^{{\rm{L}}}\,.\end{eqnarray} \tag{ 9 }$$

If δ_z ∼ 0 the trajectories present DFA functions with a single slope, then not showing multi-scaling.

To have an overall estimation of the effects investigated, the α^S and α^L exponents obtained with DFA, as well as the differences δ, were averaged over the whole woman sample and the corresponding standard deviations were considered as a measure of the variability observed.

First the values of the aforementioned scaling exponents, averaged over the sample of pregnant women participating in the experiment, were analyzed. Figure 1 shows the results obtained. The first point to be noted is that there are not differences between the values found for the left (solid symbols) and right (open symbols) feet: as we can see, the respective average exponents coincide within the statistical uncertainties, with p-values well above 0.05 in all cases.

Zoom In Zoom Out Reset image size

Also, there is a good agreement between the values found for the x coordinate (represented by the open and solid symbols) and those corresponding to the y direction (shown by the gray rectangles). The only exception is that of the exponents (open and solid squares), in standing conditions (left panel), for the 2nd. and 3rd. trimesters: in these cases the average ${{\alpha }}_{{x}}^{{\rm{S}}}$ values appear to be smaller than those found in the y direction, with p < 0.05, though the results agree within the statistical standard deviation. In all the other cases p > 0.05.

In standing conditions (left panel), the average small-τ end ${\alpha }_{x}^{{\rm{S}}}$ scaling parameters (open and solid squares) coincide, within the uncertainties, with the large-τ end ${\alpha }_{x}^{{\rm{L}}}$ ones (open and solid circles), all of them being around 0.5. The same occurs for the corresponding y average exponents (shown with the gray rectangles). All the comparisons show p-values above 0.05 except in the case of the α^S exponents for the last two trimesters mentioned above.

However, in walking conditions (right panel) a different behavior is observed. Therein we can see that the large-τ end average scaling parameters reduce significantly (p < 0.05) along the pregnancy, reaching values of the order of −0.5 in the 3rd. trimester. On the contrary, the small-τ end ones do not change, remaining at ∼0.5 as in standing conditions, with p-values well above 0.05. The negative values of the scaling parameters ${\alpha }_{z}^{{\rm{L}}}$ indicate that the fluctuations reduce with the scale of observation, the signal becoming smoother for increasing scales.

The results obtained indicate that the foot CoP trajectories show a multi-scaling character, a feature that is accentuated as pregnancy progresses: the differences between the average scaling exponents obtained at the small-τ and the large-τ ends of the DFA scaling functions increase significantly from the 1st. to the 3rd. trimester.

It is also interesting to look at these results in terms of the quantity δ, defined in equation (9). Figure 2 shows the δ values, averaged over the whole woman sample, as a function of the control trimester. The values found in standing conditions (left panel) are practically constant and close to 0, indicating that, in average, the plantar CoP trajectories do not show multi-scaling effects. However, the average δ's obtained for walking women (right panel) grow as pregnancy progresses, pointing out a change of the scaling properties of the signals analyzed between small-τ and large-τ scales.

Zoom In Zoom Out Reset image size

We also investigated the behavior of δ for each individual woman in our sample by calculating, for each of the three control trimesters and for both standing and walking conditions, the quantities δ_x and δ_y for both feet. In general the four values obtained showed up to be rather similar and, to summarize the information, their average and standard deviation were determined. Results are shown in figure 3 for the 42 women for whom the complete experimental information was available. The overall behavior observed in the average δ data of figure 2 can be seen here too. Except in a few cases the values found in standing conditions (open squares) stay close to 0 and show changes much smaller than those corresponding to walking women (solid circles).

Zoom In Zoom Out Reset image size

The results obtained in standing conditions differ from those found in previous works for the CoP of the whole body, for which a multi-scaling behavior is found (Duarte and Zatsiorsky 2001, Blázquez et al 2012). In the latter case, an anti-persistent behavior of the trajectory occurs for large τ values, when the body disequilibrium reaches situations where the falling risk becomes non-negligible. The individual foot signals that we have measured here do not show such behavior in standing conditions but it appears very clearly in the walking outcomes.

The technique described here permits to study the variations in the walking characteristics of pregnants, based on a mathematical analysis of the trajectories of the center of pressure of both feet.

In previous works (Ribas and Guirro 2007, Karadag-Saygi et al 2010, Zhang et al 2015, Mei et al 2018), both anteroposterior and medio-lateral large shifts were observed in the plantar CoP trajectories measured for pregnant women, with respect to those found for non-pregnant ones. This feature became more evident as pregnancy progressed.

The fact that the ${\alpha }_{z}^{{\rm{L}}}$ exponents obtained in our analysis are clearly smaller than 0.5 for the walking CoP trajectories, could be related to these shifts. The small, even negative, values of these exponents indicate that the walking CoP trajectories behave in an anti-persistent way in the large-τ scale. This could be at the base of the mechanisms activated by pregnant women when they try to counteract the risk of falling that increases as pregnancy progresses.

The fractal behavior of the trajectories of the foot centers of pressure in pregnant women was investigated, paying special attention to its evolution as pregnancy develops. Both standing and walking signals for 71 women were collected in the three trimesters of the pregnancy. The analysis of the trajectories was carried out by means of the DFA technique and the scaling exponents corresponding to the small and large ends of the respective scaling functions were extracted and compared.

No differences have been observed between the exponents corresponding to the x and y coordinates and to both feet. In the case of the standing data the exponents found at the two fractal scales analyzed were around 0.5 and remained practically constant through pregnancy. The situation was different for walking data where the large-end exponents reduced considerably with respect to the small-end ones, which maintained roughly the same value as in the standing case. The observed effect increased as pregnancy progressed.

This increase of the anti-persistent behavior of the plantar CoP trajectories may be linked to a displacement in the foot anteroposterior and medio-lateral pressure balances already observed in different previous experiments.

This work has been supported in part by the Junta de Andalucía (FQM0387, FQM0362, P07-FQM03163, P10-TIC5997).

In this appendix, details of the calculations carried out to obtain the scaling exponents are given in order to address some open questions about the procedure we have followed.

Figure 4 shows the whole information corresponding to the y signals, right foot, of subject #31, both in standing (left panels) and walking (right panels) conditions. The original, spatial-like, signals are shown in panels (a) and (b). The length of the standing signal is smaller that that of the walking one. This is a general characteristic of our data set. On the other hand, the walking signal shows the expected pattern because of the periodic loss of pressure when the foot does not rest on the ground. The velocity-like signals, obtained by subtracting consecutive values of the original ones, are shown in panels (c) and (d).

Zoom In Zoom Out Reset image size

Panels (e) and (f) show the log-log plots obtained from the respective signals and used to extract the scaling exponents. In the examples shown we see how the curves behave linearly at these ends. As indicated in subsection 2.4, the linear fits providing such exponents were carried out by considering the four point closest to the small and large ends of the curves. We checked that there are no significant changes if instead three or five points are considered.

The values of the scaling exponents obtained in all cases are also given. In principle, and according to equation (8), it should be verified that ν_y = π_y − 1. However, in the case we are analyzing, this is well satisfied only for the small-τ end of the walking signal. It is also approximately accomplished for the large-τ end of the standing one. But this condition is largely violated in the other two cases. This is the reason why, as explained in subsection 2.4, the exponents considered in this case as those adequately representing the scaling properties of the signal analyzed are ${\alpha }_{y}^{{\rm{S}}}\equiv {\nu }_{z}^{{\rm{S}}}$ and ${\alpha }_{y}^{{\rm{L}}}\equiv {\pi }_{z}^{{\rm{L}}}-1$.

Finally, it is worth to give a look to the role played by the length of the signals. As indicated in table 2, we have analyzed signals with different lengths M. As the size τ of the subseries analyzed ranges between 8 and M/10, the values at the respective large ends would be different and this could affect the ${\alpha }_{z}^{{\rm{L}}}$ scaling parameters. In order to check this we have repeated the analysis producing the figure 1 but neglecting all signals having less than 500 points. Figure 5 show the results obtained. As we can see the differences with figure 1 are small and the overall findings occur in a very similar way.

Zoom In Zoom Out Reset image size