### Appendix 1

### Coupling length

The interpretation of coupling length starts with selecting four tracks, two pes and two corresponding manus tracks. Since alternative pairs of manus tracks could have been selected to pair with the given pes tracks, the subsequent analysis begins with one such choice of pairings, which will then be generalized to explore other possible pairings of manus and pes tracks.

Next, we construct two couplers, one for the pair of pes locators and another for the pair of manus locators. In general, a coupler is a line segment representing the pairing of two locators (*x*_{l}, *y*_{l}) and (*x*_{r}, *y*_{r}), and its midpoint (gamma ) is simply given by:

$$ gamma = left( {frac{{x_{l} + x_{r} }}{2},frac{{y_{l} + y_{r} }}{2} } right) $$

(3)

The midpoint (gamma ) of a coupler approximates the center of a limb girdle, therefore the coupling length *CL* is the instantaneous separation between the pes and manus couplers, specifically:

$$ CL = sqrt { left( {x_{gamma ,pes} – x_{gamma ,manus} } right)^{2} + left( {y_{gamma ,pes} – y_{gamma ,manus} } right)^{2} } $$

(4)

where ({x}_{gamma ,pes}) corresponds to the x coordinate of the pes coupler, and so forth. Figure 5 shows the pes and manus couplers as line segments between corresponding pes and manus locators, and the line segment *CL* is shown connecting the midpoints of the two couplers.

Coupling length is subject to measurement uncertainty since it is a computation based on uncertain locator positions. The positional uncertainty of a locator reflects the underlying track measurement error when in support, and that uncertainty is necessarily greater when the limb is stepping between tracks. When the locator is in mid-step, its position is calculated using linear interpolation between the coordinates of the previous and next track centers, using a conventional interpolation variable α (0 ≤ α ≤ 1). Initially (at α = 0), the locator is still at the previous track, and is just finishing support. For intermediate values of α (0 < α < 1) the locator is in transition and its position is interpolated between that of the previous and next tracks. When α reaches 1 the locator is again in support at the next track. Note that α is a relative phase, i.e., relative to when the locator begins transition, which depends on *LP* and *DF* and the particular limb. We will need to also compute the absolute phase (i.e., relative to *P*_{L}) at which to start the interpolation for that locator, which is specific to the limb, the limb phase, and the duty factor. Given that, we can then compute the position of each of the four limbs even if one limb is in transition, at any point within the overall step cycle.

Since we consider only walking gaits (*DF* > 0.5), both hindlimbs are simultaneously in support during two intervals of each step cycle. We measure *CL* only at the end of each such interval (when either pes *P*_{L} or P_{R} just finishes its support phase and begins its next step). This choice of timing will be shown to permit *CL* to converge to *GL*. At those two moments within a cycle, we know the position of the pes coupler ({gamma }_{pes}) since both pes locators are in support, placed on their respective tracks. In the following, we describe only the case for *P*_{L} since the reasoning is similar for *P*_{R}. Since *DF* > 0.5 we also know that one manus locator is also in support, hence we know the position of the track on which it stands, and finally, we can also interpolate the position of the other manus locator. With both manus locators known, we can then solve for the position *γ*_{M} of the manus coupler and thus compute *CL* = (({gamma }_{manus}-{gamma }_{pes})) according to Eq. 4 for any combination of limb phase and duty factor *DF* ≥ 0.5.

### Coupling length for the walk gait

The equation for coupling length will be developed first for the familiar case of a *G2* walk (*LP* = 0.25), then will be generalized to other gaits. Figure 13 shows phase diagrams for the walk for three duty factors (*DF* = 0.5, 0.75, and 0.8). As in Fig. 2, the blue and red intervals indicate support and transition, respectively. In each case, a step cycle begins with the left pes *P*_{L} just beginning support, followed by the left manus *M*_{L} at *t* = 0.25, the right pes *P*_{R} at 0.5, and the right manus *M*_{R} at *t* = 0.75.

Coupling length is defined twice per step cycle (when either the left or right pes completes support and begins transition), but for purposes of development we consider when the left pes *P*_{L} just finishes support and begins to step. The state of the four locators is examined at this instant for three examples of duty factor (Fig. 13). For the very short duty factor of 0.5, note that while three locators are in support (indicated by green dots), *M*_{R} is in transition (red dot), midway between the previous and next track at that instant. Note that for longer duty factors (wherein *P*_{L} remains in support longer), *M*_{R} will complete more of its transition and converge on the location of the next track at that moment of measurement. For sufficiently long duty factors (i.e., *DF* ≥ 0.75), all four locators are in support at that moment.

The left pes *P*_{L} begins support at time *t* = 0.0 and it remains in support until* t* = *DF* whereupon it begins transition. At that moment, the right manus *M*_{R} is still in transition for any duty factor less than 0.75, and its position would have to be interpolated as a fraction α of the distance between the last position for *M*_{R} and its next position. For the case of a *G2* walk (*LP* = 0.25), the interpolation value *α*_{walk} depends only on the duty factor *DF*, and is given by:

$$ alpha_{walk} = frac{0.25}{{1.0 – DF}} $$

(5)

The numerator of (5) is the relative phase between *M*_{R} and *P*_{L}, indicated by the dashed black line in Fig. 13 highlighting the fraction of a step cycle between the end of support for *P*_{L} and the end of support for *M*_{R}. The denominator of (5) is the duration in which the limb is in transition (i.e., it is the relative complement of the duty factor). The quotient *α*_{walk} is therefore the fraction of the way that *M*_{R} has traveled from the previous track to the next track. For *DF* = 0.5, *M*_{R} is half the distance between the previous and next tracks, and for any value *DF* ≥ 0.75 the interpolation would be complete, and *M*_{R} would be at the next track. What are the consequences of interpolation on the coupling length computation *CL*? Since *P*_{L} and *P*_{R} are in support, the pes coupler ({gamma }_{pes}) is fixed, and *M*_{L} is also in support at this moment, so that as *M*_{R} converges on the next track as *DF* increases, *CL* converges on *GA*.

Figure 14 shows this diagrammatically for *DF* = 0.5 (upper trackway diagram) and *DF* = 0.75 (lower). The gait in both cases is a *G2* walk (*LP* = 0.25). Both diagrams show the configuration of the four locators at the instant the left hindlimb is about to initiate another step. For *DF* = 0.5 the right manus is in transition midway between two manus tracks while for *DF* = 0.75 its transition phase is complete, all locators correspond to track positions, and *CL* converges to the conventional *GL* computation (Fig. 1). For values of *DF* > 0.75, *CL* would remain as shown in the lower diagram of Fig. 14.

### Generalization to other gaits

The above discussion, while specific to the *G2* walk, also applies to the other lateral sequence gaits such as *G1* and *G3*. Recall that in the case of the walk, depending upon the duty factor, the right manus might still be in transition when we need to compute the location of the manus coupler. We interpolate the position of *M*_{R} between the previous and next right manus tracks using an interpolation value *α*_{walk}. While this interpolation value depended only upon the duty factor in Eq. 5. The generalization of (5) to apply to other gaits will depend on both *LP* and *DF*.

To illustrate, Fig. 15 diagrams the *G1* gait (lateral sequence lateral couplet, *LP* = 0.125), for three duty factors (*DF* = 0.5, 0.6, and 0.9). Recall that the numerator in (5) was graphically indicated by the dashed black line in Fig. 13. The generalization of the numerator, which will be called the α-phase and given the variable ({phi }_{alpha }) below (in Table 3 and Eq. 6), is graphically indicated by the dashed black line in Fig. 15.

The generalization of the interpolation value *α* is now an equation with two unknowns, the α-phase ϕ_{α} (which is itself a function of *LP* and *DF*), and the duty factor *DF*:

$$ alpha left( {phi_{alpha } , DF} right) = fleft( x right) = left{ { begin{array}{*{20}c} {1 ;left| {if phi_{alpha } ge 1.0 – DF} right.} \ {frac{{phi_{alpha } }}{1.0 – DF};left| {otherwise} right.} \ end{array} } right. $$

(6)

Table 3 shows the values of α-phase and the minimum duty factor *DF*_{Min} required to achieve four-limb support, for each of the eight gaits considered. This generalization extends to diagonal sequence gaits *G5*–*G7*, wherein the three-limb support is provided by *P*_{L}, *P*_{R}, and *M*_{R} and *M*_{L} in transition (and the α-phase is then measured between *M*_{L} and the pes *P*_{L}).

Coupling length *CL* can now be measured as the distance between the pes and manus couplers for any walking gait. In each case, the two pes locators are in support and one or the other manus locator may require linear interpolation, by:

$$ gamma_{manus} = left( {frac{{x_{1} + x_{2, prev} + alpha left( {x_{2, next } – x_{2, prev} } right)}}{2}, frac{{y_{1} + y_{2, prev} + alpha left( {y_{2, next } – y_{2, prev} } right)}}{2} } right) $$

(7)

where (psi 1) is the locator in support, ({(x}_{psi 1}, {y}_{psi 1})) is the position of that locator, and (psi 2) refers to the other locator which is in transition, the position of which is being interpolated.

This formula reduces to the original *GA* formula when α is 0 or 1, i.e., when the duty factor is sufficiently long to permit all four limbs to be in simultaneous ground contact for that given gait. Coupling length is therefore a discretely sampled, periodic value within a trackway that can be measured twice per gait cycle (each time a pes support cycle ends).

Given that the S1 trackway consists of 10 complete step cycles, each trial resulted in 20 coupling length measurements (twice per cycle for 10 cycles). Median *CL* was then computed for each trial for the 20 coupling lengths. Figure 16 plots the results sorted by median *CL* length. No coupling lengths were excluded, resulting in some improbably short-coupled (*CL* = 0.77 m) and some improbably long-coupled (*CL* = 4.7 m) trackmakers for the size of the tracks in this trackway. Each computed value of *CL* also has an associated uncertainty which reflects error propagation from the uncertainty in the track locations. Despite the excellent preservation of this trackway, the small uncertainties in track position result in a large degree of overlap between median coupling length values in neighboring solutions. While trials with shorter *CL* have a greater relative uncertainty than those with longer values of *CL*, the median values by themselves provide very little information. Not finding any distinguishing behavior of median *CL* across any of the 40 distinct gaits, we turned to examining how *CL* varied locally along the length of the trackway for each trial condition, hence we developed a quality-of-fit measure.

The next step is to create a metric by which we can evaluate the variation in coupling length along the trackway for each trial, in order to find a ‘best’ solution.

### Appendix 2

### Evaluating quality of fit

We considered how *CL* varies along the trackway for each condition. We began by focusing on the *G2* walk gait (specifically *G2-1-60* and *G2-2-60*) and plotting their coupling length values versus cycle along the length of the trackway. Figure 17 shows the results for *G2-1-60* and *G2-2-60*, respectively, where the configuration *G2-2* has one additional stride length separation between the pes pair and the manus pair. Each of the following plots has a one-meter range on the ordinate centered on the median coupling length value for that trial. The observed variations in coupling length values along each trial arise from two potential sources. The first is noise, which is inherent in any biological system. The second is a characteristic perturbation to *CL* in response to the track pattern perturbation. It’s this second type of characteristic variation that we will exploit, provided it can be distinguished from the noise.

With reference to Fig. 17, there are clear differences between the responses of the *G2-1-60* and the longer-coupled *G2-2-60* to this trackway. Since we hold constant both the gait and the duty factor, the difference between the two trials is essentially the ‘body length’ (by adding one stride length to the choice of manus tracks to associate with the given pes tracks). We seek to establish parameters derived from the coupling lengths in each trial that can be used to quantify the characteristic variations between trials. These parameters will form a “feature space” (a set of possible combinations of parameter values) that can be used to determine the fitness of each simulation trial as a solution for the trackway. However, as with many complex systems, the characteristic variations are not perfectly orthogonal to the noise. A feature space cannot be created that represents only the characteristic variations within trials. Some amount of noise will be included in each parameter and will have to be accounted for in the analysis using uncertainty and error propagation techniques.

What interpretation should be given to the variations in measured coupling length as it is applied along a trackway? This computation results from the application of fixed specific values for the free parameters to the positions of actual track locations. The irregularity in computed *CL* therefore results from the irregularity in the trackway, and even a trackway as subtly irregular as that in Fig. 7 could result in substantial variation in *CL*. This variation is difficult to explain in biological terms, as that would suggest the axial length of the trackmaker varied quite substantially during its passage. If we reject that conclusion, then it comes down to either having applied the ‘wrong’ choice of free parameters or having incorrectly assumed that the trackmaker walked with a constant gait as it created that trackway. To resolve this, we seek a metric on which to judge the suitability of a given gait as a model for the trackmaker’s movements along an actual trackway. We do not seek to solve for the ‘right’ choice of gait parameters; rather we seek to rule out ‘wrong’ choices, based on a fitness function which we will develop.

### Persistent residuals

Consider a fictional trackway where two simulation trials, *a* and *b*, produce coupling length plots with the same median value but rather different profiles (Fig. 18). We wish to quantify the salience of these differences.

A residual analysis is a common approach to distinguish these cases, wherein the difference between each sample and the expected value is calculated. Each difference is squared to make all values positive and to weigh relatively larger differences more strongly than smaller ones. If we specify the consistent baseline value as the expectation value, we can calculate this residual sum of the squares (*RSS*) as:

$$ RSS = mathop sum limits_{i = 1}^{N} left( {x_{i} – overline{x}} right)^{2} $$

(8)

where (overline{x}) is the baseline value. For the two plots in Fig. 18, however, the *RSS* results are identical, hence a standard residual analysis would conclude that the two trials are equally efficacious solutions for that trackway.

We can improve the residual analysis by considering the physical restrictions on the trackmaker. In trial *a*, the two deviations are isolated and short in duration (consisting of one sample above the baseline in each case) while in trial *b* there is one deviation for twice the duration. Given we are computing coupling length along a trackway, one single step of either a pes or manus that deviates substantially from the median step length could produce an isolated deviation in *CL* from the baseline. A slightly longer (or shorter) than usual pes step at the same moment as a slightly shorter (or longer) manus step could also combine to create a short-isolated deviation in the *CL* computation. But the story is different for trial *b* in Fig. 18, where *CL* deviates from the baseline for two successive measurements. For the trackmaker to sustain a deviation from baseline for two or more successive measurements (one step cycle or more) requires some change in gait parameters, if one is to believe the trackmaker cannot substantially change its axial length during that protracted interval. Taking this into account requires altering the way the residuals are calculated. Our concern is not that a residual exists, but that a residual persists for a longer interval than would reflect an isolated deviation of the trackmaker. We can define this notion of persistence by replacing the square of each residual in the RSS calculation (Eq. 8) with the product of nearest neighbor residuals, i.e.,

$$ PRSS = mathop sum limits_{i = 1}^{N – 1} left| {x_{i} – overline{x}} right|left| {x_{i + 1} – overline{x}} right| $$

(9)

Applied to our examples in Fig. 18, trial *a* has a value of *PRSS* = 0 while trial *b* has a *PRSS* = 0.023. The persistent residual effectively removes noise variations caused by single-sample deviations. With it we can conclude that trial *a* is a more efficacious solution for the trackway. In the example above we specified a reasonable, but somewhat arbitrary, value for the expectation value (underline{x}). For more general simulation trials, *PRSS* is better expressed as:

$$ PRSS = mathop sum limits_{i = 1}^{N – 1} left| {CL_{i} – overline{CL} } right|left| {CL_{i + 1} – overline{CL} } right| $$

(10)

where (overline{CL}) is the median length across the given simulation trial. Applying this *PRSS* formulation (Eq. 10) to the BEB-500 S1 trackway yields the following plot of the persistent residual sum of squares for coupling length (Fig. 19).

It is apparent that there is a scaling issue in these results. We did not observe this in our fictional example because the two trials shared the same expectation (baseline) values. Median coupling values generally differ in actual simulation trials, however. Normalizing the residuals by their median coupling values eliminates the problem:

$$ PRSS = mathop sum limits_{i = 1}^{N – 1} left| {frac{{CL_{i} }}{{overline{CL} }} – 1} right|left| {frac{{CL_{i + 1} }}{{overline{CL} }} – 1} right| $$

(11)

There is another subtle issue caused by boundary conditions at the beginning and end of the trackway. The number of samples, *N*, in a trial depends on the gait for a given trial. For example, manus locators in a long-coupled trial (such as *G2-2*) are further ahead in a trackway than in a short-coupled trial (*G2-1*) at the same moment in a simulation, and therefore will complete fewer samples. This relatively penalizes those simulation trials with more samples, which is resolved by simply dividing the *PRSS* by the number of residuals in the summation:

$$ PRSS = frac{1}{N – 1 }mathop sum limits_{i = 1}^{N – 1} left| {frac{{CL_{i} }}{{overline{CL} }} – 1} right|left| {frac{{CL_{i + 1} }}{{overline{CL} }} – 1} right| $$

(12)

PRSS can now be calculated by Eq. 12 for every simulation trial on a trackway and used to compare one trial to another. For BEB-500 S1 this yields the graph in Fig. 20.

We can now use *PRSS* values to compare the solution efficacy of the two trials *a* and *b* in terms of standard deviations:

$$ Delta PRSS = frac{{|PRSS_{a} – PRSS_{b} |}}{{sqrt {sigma_{a}^{2} } + sigma_{b}^{2} }}, $$

(13)

where ({sigma }_{(a, b)}) are the uncertainties for the *a* and *b* trials.

While pairwise comparisons are useful, a global comparison across all trials would be preferable and allow the solution efficacy of every trial to be ranked at once. We therefore generalize the previous pairwise comparison formula (Eq. 13) to:

$$ Delta PRSS = frac{{|PRSS_{i} – PRSS_{min } |}}{{sqrt {sigma_{i}^{2} + sigma_{min }^{2} } }} $$

(14)

Comparisons are made for each trial relative to the smallest persistent residual sum-of-squares value, *PRSS*_{min}, from the collection of solutions. A value of 0 indicates that there is no observed difference between a particular trial and the ‘best’ solution while values greater than 0 indicate a less efficacious solution than *PRSS*_{min}. This makes *Δ*_{PRSS} suitable as a fitness parameter for persistent residuals (Fig. 21).

### Swing

Consider another fictional trackway where this time two simulation trials *a* and *b* produce the two coupling length *CL* plots in Fig. 22. Each trial has two persistent deviation regions, one large and one small. While the two persistent deviation regions are both positive in trial *a*, in trial *b* they swing from positive to negative.

The value of *PRSS* for these two trials is 0.025 and thus fails to distinguish the differences in the sign of the residuals, which represents a greater range of *CL* values for trial *b* than trial *a*. Another parameter is therefore needed to express the overall range of coupling length values within a trial. We define this parameter, which we will call swing, as the difference between the minimum and maximum coupling length values:

$$ swing = frac{{left| {CL_{max } – CL_{min } } right|}}{{CL_{median} }} $$

(15)

It is normalized by the median coupling length value in order to avoid the scaling issue that arose earlier in the PRSS formulation when comparing trials with different median coupling lengths. For trials *a* and *b* in Fig. 22 the calculated swing values are 0.15 and 0.2, respectively. The smaller swing for trial *a* suggests that *a* is a better solution than trial *b* (i.e., the trackmaker using gait *a* would have had a lesser transient deviation from steady state than if using gait *b*).

But the swing formulation is not yet complete. We previously established the importance of persistence in the formulation of the *PRSS* parameter. Persistence pertains to swing as well. While deviation for the course of one sample can be attributed to noise, deviation for more samples likely indicates more than noise as the underlying cause. To illustrate this, suppose another fictional trackway produced the following two coupling length plots for two trials *a* and *b* (Fig. 23).

In this example, the only difference between trials *a* and *b* is the third coupling length sample, which, in Fig. 23a has zero deviation from the baseline value, while in Fig. 23b that sample has a non-zero deviation from the baseline value. Using the formulation for swing outlined above, these two example trials have the same swing value despite the greater persistence of a deviation in trial *b* than in trial *a*. This is efficiently resolved by introducing a forward-moving average, which is generally defined in the form

$$ overline{x}_{i} = frac{{x_{i} + x_{i + 1} }}{2} $$

(16)

Next, we need to incorporate the uncertainties of the samples into this forward-moving average, resulting in a weighted version:

$$ overline{CL}_{i} = frac{{frac{{CL_{i} }}{{sigma_{i}^{2} }} + frac{{CL_{i + 1} }}{{sigma_{i + 1}^{2} }}}}{{sigma_{i } + sigma_{i + 1} }} $$

(17)

where σ_{i} and σ_{i+1} are the uncertainties in coupling length values at i and i + 1. The weighted averaged coupling length values can then be substituted into the previous swing Eq. (15) as

$$ swing = frac{{left| {overline{CL}_{max } – overline{CL}_{min } } right|}}{{CL_{median} }} $$

(18)

This updated swing formulation now properly distinguishes between the two trials from the previous example with a higher swing value for trial *b* than for trial *a*. Applying this swing formulation (Eq. 18) to BEB-500 S1 yields the plot in Fig. 24.

In the same fashion as the *PRSS* parameter, the swing parameter needs to be converted into a fitness parameter for comparison across trials. We use the same approach for swing as we used for PRSS. That is, the lowest swing value among the trials is used as the highest fitness value, and all other trials are compared to it using the deviation significance calculation:

$$ Delta_{swing} = frac{{left| {swing_{i} – swing_{min } } right|}}{{sqrt {sigma_{i}^{2} + sigma_{min }^{2} } }} $$

(19)

The computed swing deviation is shown in Fig. 25.

### Solution fitness

Now that we have our two fitness parameters Δ_{PRSS} and Δ_{swing} (Eqs. 14 and 19) we combine them by plotting each trial on a scatter plot with the parameters as the values for each axis. To weigh the fitness parameters equally, the values for all of the trials are rescaled to a maximum of 1 (Fig. 26).

The closer a particular trial is to the origin of this fitness parameter space, the more efficacious a solution the trial is for the trackway. Computing the Euclidean distance for each point creates the final fitness parameter ranking for the trials

$$ fitness = sqrt {Delta_{swing}^{2} + Delta_{PRSS}^{2} } $$

(20)

And finally, to more clearly present the quality-of-fitness, the fitness for each trial is subtracted from the poorest fitness to yield the plot in Fig. 11.

## Rights and permissions

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

##### Disclaimer:

This article is autogenerated using RSS feeds and has not been created or edited by OA JF.

Click here for Source link (https://www.springeropen.com/)