Abstract
An idealized energy fluctuation model of a bird's body undergoing horizontal flapping flight is developed, focusing on the biomechanical power discernible to a bodymounted accelerometer. Expressions for flight body power constructed from root mean square dynamic body accelerations and wingstroke frequency are derived from first principles and presented in dimensionally appropriate units. As wingstroke frequency increases, the model generally predicts a gradual transition in power from a linear to an asymptotically cubic relationship. However, the onset of this transition and the degree to which this occurs depends upon whether and how forward vibrations are exploited for temporary energy storage and retrieval. While this may vary considerably between species and individual birds, it is found that a quadrature phase arrangement is generally advantageous during level flight. Gravityaligned vertical acceleration always enters into the calculation of body power, but, whenever forward acceleration becomes relevant, its contribution is subtractive. Several novel kinematic measures descriptive of flapping flight are postulated, offering fresh insights into the processes involved in airborne locomotion. The limitations of the model are briefly discussed, and departures from its predictions during ascending and descending flight evaluated. These findings highlight how bodymounted accelerometers can offer a valuable, insightful and noninvasive technique for investigating the flight of freeranging birds and bats.
1. Introduction
Birds flap their wings in order to achieve weight support and locomotion [1–3]. Experiments using highframerate video footage to monitor wing and body motions of birds or bats flying in wind tunnels have combined the findings with aerodynamic results and accelerometry to estimate overall energy expenditure during flight [4–6]. Such approaches have been experimentally valuable and theoretically illuminating, enabling the refinement of aerodynamic theory, but in studies involving freeranging animals where trailing wires and heavy equipment cannot be used, ambulatory recording of body acceleration offers a viable and practical alternative. Accelerometry was initially restricted to wind tunnel work [7], but has now been miniaturized and demands relatively little electrical power. Commercially available microelectromechanical transducers are now capable of faithfully recording highfrequency vibrations, offering a new means of studying the characteristics, kinematics and energetics of freeranging avian flight [8,9] and, indeed, animal locomotion in general [10]. Because the longterm study of birds in the wild is becoming increasingly feasible, there is new scope to assess some of the difficult choices birds face during longrange migrations [11,12].
Traditional techniques for monitoring the metabolic rate or power input of freeranging vertebrates include doubly labelled water [13,14] and heart rate, f_{h}, derived from electrocardiography [15,16]. The latter approach offers good temporal resolution but has historically necessitated calibration of f_{h} against measurements of oxygen consumption, . However, the direct translation of f_{h} to may now be possible for endotherms undergoing primary mode locomotion if augmented by knowledge of heart and body mass [17]. Accelerometry has a similar potential to monitor instantaneous biomechanical power output [18] during locomotion, and strong correlative relationships between body acceleration and have been found in animals running on treadmills [18,19]. As one might expect, body accelerations during flight are generally elevated compared with other forms of locomotion [14,20–22]; however, a theoretical understanding of how body acceleration relates to the biomechanical power of flapping flight has not yet been elucidated. This study aims to address this by setting out a mathematical model that assists the interpretation of accelerometry data captured from birds during flight. Novel measures descriptive of flight kinematics, integral to this modelling, shall also be derived.
Instruments that log acceleration are, for practical reasons, generally attached to the torso of a flying animal. As wings are coupled to the body, this offers much promise as a noninvasive tool that can help estimate the biomechanical power (and indirectly or proportionally, the metabolic costs) associated with flight [18,19,23], with the potential to augment or replace existing methods [15,24]. Activityrelated accelerations can be decomposed into the sum of dynamic and static accelerations which can be respectively derived by high and lowpass filtering either in the time or the frequency domain. To date, biologists have found the dynamic component most informative with regard to correlations with energy expenditure. Two timeaveraged measures of dynamic body acceleration (DBA) have been used when studying the energetics of animals [18,25]. Overall dynamic body acceleration (ODBA) is a running average of the L^{1}norm of the dynamic acceleration [10,18]. The variant using the L^{2}norm, which accurately encapsulates vectorial length, is known as vectorial dynamic body acceleration (VeDBA) [26,27]. An immediate difficulty with using any acceleration measure as a proxy for estimating biomechanical power in the absence of empirical calibration is that the fundamental units of acceleration, namely LT^{−2}, are different both to those of power, ML^{2}T^{−3}, and those of massspecific power, L^{2}T^{−3}. Overcoming this inevitably requires the development of some theoretical model descriptive of the biomechanics of bird flight with respect to body acceleration.
Under the hypothesis that DBA is closely related to overall metabolic costs, accelerometry has been successfully applied to a wide variety of other animals [10]. Consequently, it is not unreasonable to expect that a correlation could also exist between DBA and the biomechanical power directly discernible using a bodymounted accelerometer (body power) during steady horizontal flight and the metabolic rates of birds during flight (subject to the additional uncertainties of the value for the mechanochemical conversion efficiency of muscle [28]). Therefore, the ansatz is adopted here that the kinetics of the body should reflect the kinetics of the wings, thereby ultimately allowing the biomechanical costs incurred during avian locomotion to be estimated. This work primarily focuses on the interpretation of data from accelerometers attached to the body of a flying bird in the absence of additional information, a constraint demanding the development of some mathematical model to theoretically bridge the divide between body vibrations and overall biomechanical costs. The relationship between decomposed vertical and horizontal dynamic accelerations, and the various components of the energy associated with the body are investigated, lateral components being neglected due to the symmetrical beating of the wings. Birds must find ways of contending with the weight of gravity when airborne and flapping flight demands significant energy expenditure [1,3]. Motions of the body on the vertical axis differ from motions within the horizontal plane as they involve changes in gravitational potential. Thus, if the biomechanical power during flight is to be estimated solely from accelerometry, then it is essential to pay heed to the direction of gravity and hence also the absolute orientation of the transducer. The ultimate aim of this work is to derive estimates of accelerationbased biomechanical body power during flapping flight from first principles, potentially leading to the future interpretation of accelerometry from flying animals without the need for direct calibration. Novel statistical measures derived here may also be informative of flight kinematics, pertaining to energetically significant transitions in the wingstroke frequency, such as the relative phase and amplitude of forward and vertical body oscillations.
2. Developing a model
2.1. Preliminaries
Due to the pulsatile character of avian flight associated with the periodic contraction of powerful wing muscles, the energy associated with each wingstroke is delivered sporadically. Efficient flight confers evolutionary advantages, so the effort required by a bird to flap its wings is likely to achieve useful goals such as forward propulsion against aerodynamic drag, the countering of gravity, changes in velocity, ascent/descent and general manoeuvring. The mechanical energy of a bird will fluctuate in time and in still air there is a metabolic cost when the total energy increases. When it decreases, the dissipation of mechanical energy into the surrounding air is used to accomplish these various flight goals [29]. In the absence of gravity and an atmospheric medium, the mean mechanical power required by a vibrating but dissipationless mechanical system is precisely zero, because the total energy of the system is constant at all times. However, energy will necessarily be transformed or exchanged between different elements of the system within individual vibration cycles. Birds, on the other hand, remain airborne and sustain forward momentum despite air resistance by doing mechanical work which they never recover. Notwithstanding this, birds may be able to temporarily store and retrieve energy within individual wingbeat cycles in a similar manner to an idealized dissipationless system, so this should be considered when developing the model, as such a mechanism may provide significant flexibility to execute flight more efficiently under certain circumstances.
We assume a triaxial accelerometer of negligible mass securely affixed to the body of a bird undergoing horizontal flight at a steady air speed. It has been known for several decades that, provided adequate consideration is given to harness design, accelerometers can be externally mounted to birds with negligible oscillation relative to the body [30]. Raw data from an accelerometer can be reoriented using mathematical transformations, an issue returned to later, so this analysis proceeds by taking the zaxis to be oppositely aligned to gravity and the yaxis to correspond to the direction of forward motion. An inertial frame of constant velocity comoving with the bird is considered. The majority of the power in the measured accelerations resides at the fundamental wingstroke frequency f, the angular frequency of the wingstrokes being ω = 2πf. This is a crucial assumption, backed up by the acceleration power spectrum obtained via accelerometry experimentally collected from a flying pigeon, as presented in figure 1. Because power is concentrated at the fundamental frequency, the oscillations of the body in the vertical axis can be accurately approximated by a sinusoid with maximum excursion B relative to the mean altitude. The vertical displacement z and vertical velocity therefore vary as follows: 2.1and 2.2
It is assumed that no energy is stored either elastically or as rotational kinetic energy, leaving the gravitational potential energy U = m_{b}gz, the kinetic energy associated with the vertical axis and the kinetic energy associated with the forward axis to be considered. Here, m_{b} is the mass of the body of the bird, and g is the gravitational acceleration. The beating of both wings is assumed to always be in phase with one another. Because of the bilateral symmetry in birds, we neglect any kinetic energy associated with the xaxis.
The wings also have gravitational and kinetic energies (both vertical and horizontal). However, these are not directly estimable from a bodymounted accelerometer, so terms for these quantities are not explicitly included in the analysis. Instead, it is assumed that the activity detected by the accelerometer will, to a first approximation, be an attenuated reflection of the total biomechanical output of the bird. It is reasonable to think that energy associated with the wings can be subsumed into existing terms by constructive or destructive superposition. This follows from the fact that the addition of two arbitrarily scaled sinusoids of the same frequency, but different phase results in a rescaled sinusoid of the same frequency: 2.3where and H is the Heaviside step function in which H(ξ) = 1 if ξ ≥ 0, otherwise H(ξ) = 0. It can be seen that the resultant sinusoid will generally have an intermediate phase shift.
2.2. Relative phase between vertical and horizontal oscillations of the body
Vibrations on both axes will exhibit simple harmonic motion with both K_{y} and K_{z} varying at the same frequency. However, before an expression for the variation in y can be written as a sinusoid, there is a need to carefully consider its relative phase θ with respect to z. Even if the bird has considerable freedom to adjust the phase of its motion on the horizontal axis, arrangements that minimize the power required in order to sustain flight confer evolutionary advantages. Introducing a horizontal vibration amplitude A, the counterpart to B for the vertical axis, one can write 2.4and 2.5
The convention is adopted that neither A nor B can be negative. The total energy associated with the bird's body, E_{b}(t), is 2.6and 2.7
During each wingstroke, the flight muscles of the bird must supply energy when E_{b} is increasing. An important aspect of the present model is the assumption that when E_{b} decreases, energy is irrecoverably lost to the environment and exploited so as to achieve propulsion and weight support. This can be compared with aerodynamic models which assume the kinetic energy of the wings is never recovered [31,32]. The mean power associated with the body, 〈P_{b}〉, can then be determined according to 2.8
Consider, for now, the case in which A = 0. One then has K_{y} = 0 and E_{b} = U + K_{z}. Gravitational potential energy attains a minimum at t = π/ω when z = −B and a maximum at t = 0 when z = B. At both these times, K_{z} reaches its minimum of zero. Because sin^{2}ωt = (1 − cos 2ωt)/2, the vertical kinetic energy varies sinusoidally at twice the frequency of the gravitational energy. Thus, E_{b}(t = π/ω) is always a minimum. Although U is always a maximum at t = 0, it is possible that E_{b} is not a maximum at that time if K_{z} thereafter increases more rapidly than U decreases. This possibility is apparent upon inspection of the time derivatives of E_{b}, the roots of which correspond to stationary points: 2.9and 2.10Owing to the sin(ωt) term, stationary points exist at t = 0 and t = π/ω. At t = π/ω, one has , so this stationary point is always a minimum. However, for the stationary point at t = 0, which corresponds to a maximum when ω^{2} < g/B and a minimum when ω^{2} > g/B. If both these are minima, additional stationary points occur when cos(ωt) = g/Bω^{2}. Using cos(2ωt) = 2g^{2}/B^{2}ω^{4} − 1, one then finds that at those times , confirming that these points correspond to maxima of E_{b}.
Therefore, when ω^{2} < g/B, the mean body power is simply 〈P_{b}〉 = (ω/2π) (E_{b}^{max} − E_{b}^{min}), where and . However, when ω > g/B, then due to the existence of a new maximum in E_{b} at t ≠ 0, will exceed E_{b}(0), and the mean power will inevitably rise. During intense flight, there will be a high wingstroke frequency, and for each wingstroke the variation in kinetic energy will increase and the variation in gravitational energy will decrease, allowing the variations in kinetic energy to become dominant. However, if it were possible to temporarily store some of the vertical kinetic energy as horizontal kinetic energy, and retrieve it later in the wingstroke cycle, then this elevated maximum in E_{b} could be avoided, and the concomitant increase in power is eliminated.
The transition between the two regimes occurs at g = Bω^{2}, corresponding to the peak gravitational energy U^{max} = m_{b}gB being equal to twice the peak vertical kinetic energy . Clearly, if the total kinetic energy K_{yz} = K_{y} + K_{z} did not fluctuate at all, then K_{z} could be arbitrarily large without incurring any additional increase in mean power. The total kinetic energy is 2.11and 2.12If K_{yz} is constant, then its derivative 2.13must be zero at all times. It is apparent from this expression and (2.3) that K_{yz} varies sinusoidally. The amplitude of these fluctuations vanishes when 2.14and 2.15Because θ is constant, this would demand that t is also constant. However, K_{yz} can be zero, if the sin 2ωt and the cos 2ωt terms are simultaneously zero. One then finds that the cos 2ωt term vanishes if sin 2θ = 0, which is satisfied when . The sin 2ωt term is zero when A^{2} cos 2θ + B^{2} = 0, which yields θ = arccos(−B^{2}/A^{2})/2. When θ = nπ, one has cos 2θ = 1 which must be rejected as it predicts A^{2} = −B^{2}. However, θ = (n ± 1/2)π yields A^{2} = B^{2}, which is acceptable. Therefore, for K_{yz} to remain constant requires A = B and θ = ±π/2. This can be seen in the following 2.16and 2.17A relative phase shift of ±π/2 between the horizontal and vertical axes corresponds to what is commonly termed a quadrature phase arrangement. Only in this circumstance does it hold that K_{y} is a maximum when K_{z} is a minimum, and vice versa. This maximizes the potential for shuttling energy back and forth between the two axes, a useful property that the bird might be able to exploit to decrease its mean power. Although K_{y} will be maximal at t = 0 and t = π/ω, because these maxima are equal, they have no effect on the difference E_{b}(t = 0) − E_{b}(t = π/ω). In fact, this holds for any value of θ because sin^{2}x = sin^{2}(x + π), but K_{yz} is constant only when K_{y} and K_{z} are in antiphase. As the kinetic energies vary at double the fundamental wingstroke frequency, this occurs when θ =±π/2.
While variations in K_{yz} can be completely eliminated, it may not be necessary for the bird to do so because, as will be shown, the bird can in some cases also minimize power when A/B < 1 with excessive fluctuations in K_{z} being completely tamed by smaller fluctuations in K_{y}. This may be preferable as it reduces extraneous energy losses and ameliorates the vibrations transmitted to the head of the bird, which might otherwise make flight an unnecessarily uncomfortable experience compromising visual acuity [33]. Note also that if A/B > 1, then fluctuations in K_{y} may not be adequately absorbed by fluctuations in K_{z}.
Because birds are expected to have a maximum forward velocity at the end of the downbeat, the phase that makes maximal is chosen (θ = −π/2). Hence, θ can be eliminated from the expressions for y and by writing 2.18and 2.19The upbeat commences at t = 0 when is maximal. The body is then at its maximum height above the ground. The body and the wings are in antiphase on the vertical axis, and the same should also be approximately true of the forward axis. To summarize the findings of this section, power reduction is possible only when g/Bω^{2} < 1 and is best achieved by a quadrature phase arrangement.
2.3. Quadrature phase flight
Unless expressly stated otherwise, the analysis now proceeds by assuming quadrature phase applies. In order to assess the biomechanical power in the body, one is interested in determining the maxima and minima (stationary points) of the body energy, necessitating looking for roots of the firsttime derivative and inspecting their signs by taking the secondtime derivative. Both K_{yz} and U vary sinusoidally but because they do not vary at the same frequency their sum is not a simple sinusoid, demanding that calculus be used. Body energy now reads 2.20and the first derivative is 2.21and 2.22Stationary points exist when either cos(ωt) = gB/ω^{2}(B^{2} − A^{2}) or sin(ωt) = 0. The first condition has real solutions only if ω^{2} ≥ gB/B^{2} − A^{2}. When ω^{2} > g/B, the smallest value of A that prevents a maximum in E_{b} from arising anywhere but at t = 0 occurs when cos(ωt) = 1, from which the smallest value of A that minimizes the power can be determined. If this is an overriding consideration with the need to maintain visual acuity a secondary concern, this value of A would seem to be optimal: 2.23
Note that A_{opt} must be real (because ω^{2} > g/B) and because at all times, birds have no need to fly with A > B. Moreover, for given values of A and B, the optimum angular frequency comes to which exceeds the critical frequency above which power can be reduced if A ≠ 0. The optimum ratio of A/B ensures fluctuations in K_{y} are sufficiently large to avoid excessive fluctuations in K_{z}, its value depending only on the ratio ω/ω_{c}: 2.24
This is plotted in figure 2. The second derivative of the body energy is 2.25
The condition that at t = 0 yields A = A_{opt}, showing that the maximum at t = 0 is then also a stationary point of inflection. Now, consider the stationary points that arise when cos(ωt) = gB/ω^{2}(B^{2} − A^{2}). Using the identity cos(2ωt) = 2cos^{2}(ωt) − 1, one finds that 2.26and 2.27If A^{2} < B^{2}, the sign of must be negative if gB < ω^{2}(B^{2} − A^{2}). This satisfies the condition required for the stationary points to exist, establishing that those corresponding to cos(ωt) = gB/ω^{2}(B^{2} − A^{2}) < 1 must be maxima if A^{2} < B^{2}. Considering the alternative situation in which A^{2} > B^{2}, the sign of must be positive if gB < ω^{2}(A^{2} − B^{2}). Hence, the stationary points at which cos(ωt) = gB/ω^{2}(B^{2} − A^{2}) < 1 must then correspond to minima.
Note that the stationary point in E_{b} at t = π/ω has second derivative 2.28
This is negative if gB < ω^{2}(A^{2} − B^{2}) so that when A^{2} > B^{2}, which gives rise to stationary points at cos(ωt) = −gB/ω^{2}(A^{2} − B^{2}), the stationary point at t = π/ω becomes a maximum. This is a point of inflection when gB = ω^{2}(A^{2} − B^{2}) corresponding to the maximum value of A at which power is minimized: 2.29
This has a similar form to the earlier expression . The bird has minimized its mean body power if A_{opt} ≤ A ≤ A_{max} or equivalently if ω^{2} < gB/B^{2} − A^{2}. When this condition is satisfied and ω^{2} > g/B, flapping flight demands less power than when forward axis vibrations are absent altogether (A = 0). For a fixed value of B, this is illustrated in figure 3.
3. Estimation of body power
The calculation of mean body power here involves the integration of E_{b} only as it increases according to (2.8). This can be easily accomplished by subtracting the minima of E_{b} from successive maxima of E_{b}. If t = 0 at t = t_{0}, t = π/ω at t = t_{1} and t = ±arccos[gB/ω^{2}(B^{2} − A^{2})]/ω at t = t* then only these times need be considered in order to determine the minima and maxima of the body energy variation. With reference to figure 4, it can be seen that there are three distinct cases to consider with stationary points located at the following times:

case 1: t = t_{0} and t = t_{1}: {gB > ω^{2}B^{2}−A^{2}},

case 2: t = t_{0}, t = t_{1} and t = t* {gB < ω^{2}(B^{2}−A^{2}); A^{2} < B^{2}; E_{b}(t*) > E_{b}(t_{0})} and

case 3: t = t_{0}, t = t_{1} and t = t* {gB < ω^{2}(A^{2}−B^{2}); A^{2} > B^{2}; E_{b}(t*) < E_{b}(t_{1})}.
The body energy at t_{0} and t_{1} is simply 3.1and 3.2When t = t*, sin^{2}(ωt*) is obtained as 1 − cos^{2}(ωt*): 3.3
After some algebra, the total energy of the body when t = t* is found to be 3.4
3.1. Case 1
This case corresponds to the shaded region in figure 3 and the central traces of figure 4 for which the calculation of mean power is particularly simple: 3.5
3.2. Case 2
This case corresponds to the zone below the shaded region in figure 3 (see also the uppermost trace in figure 4). The stationary point at t_{0} is now a minimum and two new maxima arise at t = t*. Recalling that in this case, B^{2} > A^{2}, the power is 3.6and 3.7
3.3. Case 3
This case corresponds to the zone above the shaded region in figure 3 (see also the lowermost trace in figure 4). The stationary point at t_{1} is now a maximum and new minima arise at t = t*. One has A^{2} > B^{2}, and the power now comes to 3.8and 3.9The poles appearing in the expressions for 〈P_{b2}〉 and 〈P_{b3}〉 as A → B are avoided because when ω^{2} > gB/B^{2} − A^{2} one has A^{2} ≠ B^{2}. It is apparent that the results for cases 2 and 3 are equivalent but for the reversal of sign in the B^{2} − A^{2} terms. Therefore, it would be acceptable to take the modulus of either expression without expressly checking whether A^{2} > B^{2}. In the limiting cases where gB = ω^{2}B^{2} − A^{2}, the prediction of case 1 coincides with that of case 2 or 3. For instance, when then B^{2} − A^{2} = gB/ω^{2} and it can be seen that 〈P_{b2}〉 reduces to 〈P_{b1}〉: 3.10
The importance of correctly distinguishing between case 1 and cases 2 and 3 is stressed, because 〈P_{b2}〉 and 〈P_{b3}〉 overpredict the true power when gB > ω^{2}B^{2} − A^{2} and 〈P_{b1}〉 underpredicts the true power when gB < ω^{2}B^{2} − A^{2}.
3.4. Deriving body power using data from accelerometers
The root mean square (or r.m.s.) value of a quantity is a measure commonly used in physics and engineering and can lend itself to the description of DBA (table 1). As such, it offers an alternative to ODBA and VeDBA. The r.m.s. value of a discrete set of N dynamic single axis accelerations a_{dyn} is , and that of a simple sinusoid such as ψ = β sin(t) is 3.11
Static and dynamic accelerations can be respectively derived from raw acceleration data using lowpass and highpass filtering techniques. While the static acceleration is useful in determining the vertical, gravityaligned axis, there are many circumstances where the direction of the forward axis is more ambiguous. However, birds will generally adjust the roll of their bodies during flight so that the static acceleration vector remains dorsally aligned. It is therefore likely that the r.m.s. value of the dynamic sway, , will be appreciably smaller than the r.m.s. value of the dynamic surge, . In principle, this allows for reorientation of accelerometry data during postprocessing by application of a rotation matrix whose components can be inferred by analysis of the data. A method for achieving reorientation is outlined in appendix A. Estimates of the dynamic surge, , and the dynamic heave, , are readily obtained after highpass filtering of the reoriented acceleration components. These relate to ω, A and B as follows: 3.12and 3.13Hence, body power can be directly evaluated from r.m.s. heave and surge: 3.14and 3.15These expressions conveniently obviate the need to double integrate acceleration data in order to obtain the values of A and B directly, which is generally challenging due to the baseline drift introduced when integrating. However, A and B each feature in the trueorfalse test gB < ω^{2}B^{2} − A^{2} that determines which expression for power is valid. Therefore, a reformulation of the discriminant is also desirable, and because it follows that 3.16
It is therefore possible to determine, without knowledge of either A or B, nor indeed ω, the regime in which the bird is flying. 〈P_{b1}〉 should be used when , otherwise 〈P_{b2,3}〉 is applicable.
The critical frequency ratio, ω/ω_{c}, can similarly be translated using (3.13) and the knowledge that 3.17
It is also possible to express ω/ω_{opt} using only r.m.s. accelerations: 3.18
3.5. Relative body power
According to this model, body power retains linearity with wingstroke frequency until, and as illustrated in figure 5, the angular frequency exceeds 3.19
This can also be conveniently gauged from r.m.s. accelerations alone because 3.20
The expressions for 〈P_{b2}〉 and 〈P_{b3}〉 can be recast as 3.21and 3.22It is now apparent that when , the response becomes asymptotically cubic at higher wingstroke frequencies, i.e. 〈P_{b2,3}〉 ∝ω^{3}. However, at lower wingstroke frequencies, flight proceeds within the linear regime of (3.5). For case 1, it is possible to write , allowing the ratio of 〈P_{b1}〉/〈P_{b2,3}〉 to be expressed in a particularly simple form 3.23
The response is plotted in figure 6 and can be compared with figure 5 where the wingstroke frequency is expressed in units of ω/ω_{c}. If ω < ω_{0}, then the expression for 〈P_{b2,3}〉 is not physically meaningful, and the mean power is always given by 〈P_{b1}〉. When ω > ω_{0}, power can always be reduced, because the ratio 〈P_{b2,3}〉/〈P_{b1}〉 then exceeds unity. Furthermore, it grows without limit as ω → ∞.
Power is always minimized when A = B but it is interesting to know how the ratio 〈P_{b2,3}〉/〈P_{b1}〉 grows for other values of A/B when ω > ω_{0}. First, note that ω_{0} is related to ω_{c} according to 3.24
The relative flight cost ratio 〈P_{b2,3}〉/〈P_{b1}〉 can be expressed either in terms of ω/ω_{0} or A/B and ω/ω_{c}: 3.25and 3.26A colourcoded contour plot of this function is presented in figure 7 (online).
4. Nonquadrature phase
We reiterate that this modelling pertains to steady, horizontal flight. The possibility exists that departures from quadrature phase may be advantageous during ascent or descent, but the mathematics in such cases is considerably more involved. Nevertheless, these situations were numerically investigated under the assumption that A = B with an additional power component representing the change in gravitational energy with time, whose mean value can be estimated using an altimeter or GPS device [34,35]. The results are presented in figure 8. It can be seen that predictions of the quadrature phase model are still accurately upheld in most circumstances involving realistic rates of ascent or descent.
Although the optimal arrangement is never far from quadrature phase, the error is sensitive to the rate of ascent, and is asymmetrical in that it grows faster with descent than ascent. This is evident from the plot presented in figure 9. It has been reported that the power requirements of moderate ascending and descending flight in pigeons can be accurately estimated by summing the power required for level flight with the rate of change of gravitational potential [36]. Interestingly, the same study found that a discrepancy did arise for high descent rates but not for high climb rates, descent being clearly more expensive than anticipated. When flying steeply downwards at a descent angle of −60°, the pigeons flew at a horizontal velocity of 3.6 m s^{−1} and a vertical rate of descent exceeding 3 m s^{−1}. Hence, one possible explanation for this could be that flight costs during rapid descent are so minor and so rarely encountered that there is little need for birds to acquire biomechanical flexibility that would allow significant deviations from quadrature phase. However, this descent rate is more rapid than that which would be expected if the bird were simply gliding (typically no more than 2 m s^{−1} for a pigeon).
If more accurate estimates of biomechanical power are required then direct integration is an option, obtaining relative velocities and displacements for the forward and vertical axes from which the individual energy terms can be derived, combining this information with altitude data (if available) to obtain body energy, and calculating power by numerically evaluating the mean rate of energy increase with time, being careful to ignore periods when the total energy of the bird is decreasing. This approach may offer improved precision, particularly when power is not all concentrated at the fundamental wingstroke frequency or when large departures from quadrature phase are expected.
5. Discussion
Despite the considerable complexities involved in flapping flight, by focusing on the consequences for the mechanical motions of the body, this oscillatory energy fluctuation model provides a useful initial step in theoretically underpinning the use of bodymounted accelerometers to estimate the relative costs of horizontal flapping flight in birds or bats. Accelerationbased proxies for the power detectable in the body of an animal during flight have been derived from first principles ((3.5), (3.7) and (3.9)) and in the appropriate units of ML^{2}T^{−3}. Within these equations, the DBA formalism is encapsulated very naturally via r.m.s. acceleration ((3.14) and (3.15)), or aRMS. These expressions represent a substantial improvement over attempting to estimate the biomechanical energy expenditure or metabolic rate of flying animals from ODBA or VeDBA alone. The magnitude of ODBA varies with orientation, and different rescaling factors must be applied to single axis projections of ODBA and VeDBA if they are to be used to estimate or (appendix B). We point out that aRMS, which is more closely related to VeDBA than ODBA, is just as straightforward to calculate. Furthermore, aRMS may be more universally applicable in future mathematical and empirical studies of animal locomotion. DBA measures lack the units of power and so will always require a direct calibration against the rate of energy turnover. Their use in estimating the biomechanical costs of flapping flight has generally lacked a firm theoretical basis, pays no heed to sensor orientation relative to gravity, ascribes undue significance to accelerations in the horizontal plane and overlooks the significance of wingstroke frequency. As can be seen in (3.22), body power is rather sensitive at times to wingstroke frequency, so may in itself provide a valuable means of gauging flight effort independently of body power estimation.
This analysis furnishes several novel parameters that may provide insights into the kinematics of flight. In particular, it predicts that there may be two regimes of flight with regard to power production and wingstroke frequency and that there may be a transition from slow to moderate intensity flight, when equation (3.23) and figures 6 and 7 suggest power savings may sometimes be possible. For steady horizontal flight involving sinusoidal vibrations on the forward and vertical axes, it has been shown that a quadrature phase arrangement is potentially advantageous. If flight strategies could exploit this, temporarily storing and retrieving energy during each wingbeat cycle, it might be possible for some species to avoid or curb what would otherwise be a thirdorder sensitivity in mean body power to wingstroke frequency. Due to the ubiquity of predators, the need to catch airborne prey and the general requirement for economy of locomotion, there may have been considerable evolutionary pressure for birds to waste very little energy at wingstroke frequencies exceeding ω_{0}. Nevertheless, bodymounted accelerometry is not privy to the subtleties of wing flexion, angle of attack and feathering, so this model leaves open the possibility that birds have considerable scope to adjust their flight style without necessarily compromising efficiency.
When birds are flying in the linear regime then equation (3.14), which reflects the costs of combating gravity, shows that only the vertical, gravityaligned component of the acceleration should enter into the calculation of body power. Although forward accelerations become relevant to power estimation for the asymptotically cubic regime, equation (3.15) shows that their contribution to body power is subtractive. This somewhat counterintuitive need to subtract from in (3.7) stems from the fact that variations of gravitational energy do not always mask variations in kinetic energy. However, quadrature oscillations in K_{y} tend to erase fluctuations in K_{yz}, not reinforce them.
While the model directly considers vibrations only on the body of the bird, owing to mechanical coupling, the same kind of oscillations and tradeoffs should also apply to the horizontal and kinetic components of the wings. However, because the centre of mass of the wings must travel a much greater distance during each wingstroke than the centre of mass of the body, for any given wingstroke frequency, fluctuations in wing kinetic energy grow quadratically with wingstroke excursion, but variations in wing gravitational potential grow only linearly. Therefore, the onset of the transition between the linear and asymptotically cubic flight power regimes might be generally expected to occur at a lower frequency for the wings than for the body. The kinematics of the wings, which cannot be directly measured by the accelerometer, should generally dominate the overall biomechanical costs of flight. Indeed, kinetic energy fluctuations might easily become a more important consideration than compensating for the gravitational energy losses of each wingbeat. If the gravity terms are neglected, then (3.7) and (3.9) simplify somewhat and predict . This could be especially true of energetic highspeed flight, even though there may then be significant wing retraction and supination to avoid undue aerodynamic drag [37,38]. Because the forces involved in wing retraction tend to cancel on the body due to bilateral symmetry, the costs are hidden from bodymounted accelerometers. Hence, the ratio of perceived body power to true total biomechanical power might be somewhat reduced during intense flight for some species, which may well require the tailoring of flight models to each species in the future, following empirical observations and extensions of the modelling.
A number of other original summary statistics could prove useful to the interpretation of bodymounted accelerometry data obtained from flying birds. The ratio A/B may help to characterize the mode of flight performance and perhaps also evaluate the skill and dexterity of individual birds. With the possible exception of hovering flight when birds may be able to recoup some of the kinetic energy stored in the air during the previous halfstroke, flying animals generally have no means of recovering energy lost to the environment. However, vibrations on the forward axis offer convenient energy storage which may also be exploited to reduce pitching of the body. Many birds use various reflex mechanisms during flight to subdue head vibrations and thereby avoid vision impairment [33]. Body accelerations have a total amplitude which for small A/B can be approximated by ω^{2}B(1 + A^{2}/2B^{2}). Hence, the additional loss of visual acuity due to a small A/B ratio would be relatively imperceptible, implying little need for A/B → 0.
The ratio A/B may also be particularly sensitive to effort in realistic situations, correlating with wingstroke frequency and increasing at higher forward velocities. The relative phase lag θ between the vertical and forward axes may act as a marker of ascending or descending flight, or reflect efforts to synchronize wingstroke frequency with other birds during Vformation flight. We also expect the ratios ω/ω_{c}, ω/ω_{opt} and ω/ω_{0} to be informative regarding flight intensity and flight efficiency. Collectively, these measures may also crudely encode hints as to the altitude at which the bird is flying. Due to the complexities of wing kinematics and anatomical constraints, it is very likely that no simple unifying pattern will adequately summarize all species, but departures from normality are often the most interesting aspects of biological research and so additional parameters can prove very useful in highlighting departures from nonconformity. Therefore, these flight variables may be particularly valuable in helping to unpick the challenges involved in flying efficiently. The static acceleration also offers a potentially illuminating variable for flying animals which has been largely ignored to date. Birds can sustain prolonged banking when circling or jostling for position within a cluster flock [39], and any drift in the mean direction of the momentum vector induces a nongravitational contribution to the static acceleration. Thus, when the static acceleration deviates appreciably from gravity, it would suggest that the bird is not undergoing steady horizontal flight. However, the converse is not true, because one also expects the static acceleration to tally with gravity during steady nonbanking ascending or descending flight. Therefore, the distribution and time variability of the static acceleration can be informative.
Dimensional considerations may allow the results obtained here to be extrapolated to some degree, particularly regarding the estimation of biomechanical power from accelerometry for aquatic animals. Due to the buoyancy afforded by water, the estimation of inertial costs for aquatic species during swimming is not encumbered by gravitational considerations [2]. This invulnerability to gravity suggests that the cost of swimming should correlate with the product of body mass, the period of the swimming stroke and some function of the square of the decomposed r.m.s. accelerations, dependent upon the anatomy of the species under consideration. Locomotion costs in terrestrial animals are likely to be more complex: weight support can either be provided continuously or episodically by the ground.
Naturally, there are limitations to what a bodymounted accelerometer alone can glean about flight. During freeranging flights, there could be circumstances where basic inferences may be misleading, particularly if the rate of ascent or descent is unknown. In addition, special care may be needed when attempting to disentangle the static and dynamic accelerations for birds using intermittent modes of flight such as flapgliding or flapbounding. Accelerometers cannot infer absolute velocities in any direction, and many birds exploit the assistance of thermals, following winds and airflow over uneven terrain, all of which are capable of drastically altering the power requirements of flight. Nevertheless, the present model offers a practical and noninvasive method of extracting from accelerometry a variety of parameters that could be informative concerning flight style and performance, while also providing an explicit procedure for determining biomechanical body power in freeflying birds which may be generally proportional to overall flight power. In the complex processes that transform the biochemical energy of birds into atmospheric vortices, turbulence and heat, aerodynamic costs lie downstream of the biomechanical costs. While it might eventually be possible to incorporate them within an extended model, the formidable challenges of contending with complicated wake patterns, vortex interactions and chaotic flow patterns continue to plague theoretical models, and quasistatic approximations to the Navier–Stokes equations commonly used in aerodynamic analysis inherently limit their accuracy and usefulness [40]. A more realistic nearterm goal would be to broaden the present modelling to include wing kinematics and morphology. This will inevitably necessitate the input of anatomical information allowing the body power relationship to accommodate allometric differences between species and also address wingpropelled locomotion in the media of differing densities. Due to the intrinsic complexities, we anticipate that experimental data collected from a variety of species will be required. Birds are graceful aeronauts, skilfully adjusting their posture and technique in ways we have only started to perceive [41]. However, a battery of physiological, biomechanical and aerodynamic techniques can augment and refine one another when quantifying flight costs. These complementary approaches to studying avian energetics hold much promise in arriving at a more unified understanding of the compromises involved when animals fly—whether they are foraging for food, migrating, chasing airborne quarry or evading predators.
Funding statement
This work was supported by Bangor University and the Biotechnology and Biological Research Council, UK (BBSRC) grant BB/F015615/1 awarded to Dr C. M. Bishop and Prof P. J. Butler.
Appendix A. Reorientation during postprocessing
When flying in a straight line, the static acceleration , where represents the timeaveraged mean of the kth accelerometer channel, should equal the gravitational acceleration, g. When this does not hold, it is possible that the bird is turning, undulating, bounding or flying through zones containing vertical air currents. Centripetal acceleration combines with gravity during turns, the direction of the resultant vector determining the degree of banking necessary in order that the bird experiences the net acceleration dorsally so that the forces on the wings are symmetrically balanced. At such times, the increase in the static acceleration relative to gravity is a useful guideline as to the departure from linear motion. For uniform rectilinear motion, the mode of the distribution of the static acceleration magnitude can be identified with gravity. To some degree, this allows for the selfcalibration of accelerometry data.
When a_{s} ≈ g, a condition that is easily checked, it is straightforward to calculate the vertically aligned component of the acceleration a_{vert} using the scalar projection of a onto a_{s}: A1
This vertical acceleration is the combination of both static and dynamic components, so the dynamic vertical acceleration a_{z} can be obtained by subtracting the timeaveraged value of a_{vert} using A2
Projections of a orthogonal to the vertical axis will then lie in the horizontal plane: A3and A4As before, the dynamic acceleration is obtained by subtracting the static acceleration ( and ). A means of determining a_{x} and a_{y} from a_{h1} and a_{h2} is then required. By symmetry, one expects the dynamic acceleration of the forward (y) axis to exceed that on the lateral (x) axis. One way to proceed would be to first determine the fundamental wingstroke frequency using a_{z} then bandpass filter a_{h1} and a_{h2} using fast Fourier transforms so that only frequencies near the fundamental wingstroke frequency are retained. Following this, one could compute the angles ϕ = arctan(a_{h}_{2}/a_{h}_{1}), compile a circular histogram of the results and use the angle corresponding to the peak in the histogram Φ, to reorient a_{h1} and a_{h2} as follows: A5and A6There is a potential ambiguity in this result concerning the polarities of a_{x} and a_{y}. It arises due to the fact that one expects two peaks in the angular histogram separated by π. Therefore, one should also evaluate a_{x} and a_{y} using Φ ← Φ + π. In practice, there may be no need to calculate a_{x} because it is unlikely to represent interesting information, but the polarity of a_{y} is potentially important if one is eager to know, for example, how the phase shift between the forward and vertical axes varies with time. This model expects that the displacement on the yaxis will always lag behind that of the zaxis, and hence a_{z} should always lead a_{y}, resolving the ambiguity. One can again apply bandpass filtering to a_{y} and a_{z} around the detected wingstroke frequency in order to test which value of Φ is appropriate. The phase shift can be accurately measured in the recovered time domain after bandpass filtering using linear interpolation between samples in the vicinity of the positive/negative going zero crossings. If the orientation of the accelerometer is fixed with respect to the body of the bird, then this process need only be performed once, and the value of Φ can then be reused without recalculation.
Appendix B. Relationship between overall dynamic body acceleration, vectorial dynamic body acceleration and r.m.s. acceleration
Consider the instantaneous dynamic acceleration a_{d} = (a_{x}, a_{y}, a_{z}) derived from a triaxial accelerometer. The generalized mean or L^{p}norm of the components of a_{d} is defined as B1
Therefore, the magnitude of a_{d} according to VeDBA is whereas its magnitude according to ODBA is a_{x} + a_{y} + a_{z}. Converting to spherical coordinates (r, θ, ϕ), one has B2 B3 B4The magnitude of a_{d} using VeDBA is thus B5
This confirms the standard expectation of Euclidean trigonometry and the Pythagorean theorem. However, the magnitude of a_{d} in the case of ODBA generally disagrees with this because B6and B7Due to the second term, ODBA will exceed VeDBA unless at least two of the three acceleration components are zero. The two will agree only when there is alignment of the acceleration with one of the three measurement axes. Hence, the response of ODBA varies according to orientation.
To find the maximum error in ODBA, let a_{y} = a_{x} + α and a_{z} = a_{x} + β where, without loss of generality, it can be assumed that neither a_{x}, a_{y} nor a_{z} are negative. Now let ζ = (ODBA^{2}/VeDBA^{2} − 1)/2 so that the stationary points of ζ will be identically located to those of ODBA/VeDBA: B8
Equating to zero the partial derivatives of ζ with respect to α and β gives B9and B10These expressions reduce to B11and B12The first condition is satisfied only if α = β, which also ensures the second is satisfied. When α ≠ β, the second condition holds if α + β = −3a_{x} which cannot be true given that neither a_{y} nor a_{z} are negative, which implies that α + β ≥ −2a_{x}. This leaves only α = β and hence a_{y} = a_{z} after restoring moduli. For a_{x} → 0, one then finds that (ODBA/VeDBA)^{2} can be no less than 2. This can be seen by letting γ = a_{y} = a_{z} and a_{x} = ε with : B13
However, as a_{x} → a_{y} = a_{z} one finds that (ODBA/VeDBA)^{2} can be no greater than 3. This can be seen by letting a_{x} = γ ± ε where, once again, and γ = a_{y} = a_{z}: B14
Hence, ODBA/VeDBA is maximized when the acceleration components on each axis have an identical magnitude. Because the minimum is already known to occur when two of the components are zero, one may conclude that ODBA is confined to the range B15
It is interesting to ask whether VeDBA can be recovered from historical records of ODBA. This is possible only in an approximate statistical sense, and best results would be obtained when data have been collected from animals whose orientation in space varies considerably or in situations where the orientation of the accelerometer itself is free to drift. In order to determine the rescaling factor, it is first necessary to calculate the mean exaggeration of ODBA relative to the true acceleration magnitude. Due to symmetry, it is sufficient to consider the solid angle Ω = π/2 corresponding to the octant 0 ≤ θ ≤ π/2, 0 ≤ ϕ ≤ π/2. Although the mean value of VeDBA is simply the vector length r, the mean value of ODBA is B16 B17 B18 B19 B20
Therefore, on average, ODBA exaggerates the true acceleration magnitude by 50% and so, in some circumstances, ODBA data records can be translated into estimates of VeDBA simply using B21A highly significant linear relationship between these two measures has already been experimentally observed, with the best fit corresponding to VeDBA ≈ 0.014 + 0.6418 ODBA for units of g [27]. The slope of this empirical relationship agrees with the theoretical value to within 4%. In the same work, the envelope of figure 2 exhibits a wedge distribution whose upper and lower slopes are approximately unity and , corresponding to the anticipated range in error of ODBA due to changes in orientation.
When evaluating r.m.s. accelerations, there is no explicit requirement to calculate vector lengths, but we stress that the Euclidean formulation, as adopted by VeDBA, is implicit in the present modelling. While VeDBA and r.m.s. acceleration are in agreement here, and one would generally expect to find an excellent correlation between them, their magnitudes are anticipated to differ whenever there is any spread in the distribution of the dynamic acceleration data. For complex or aperiodic acceleration profiles, we recommend the r.m.s. method over VeDBA in all cases, because the ambiguity implies that no single rescaling factor will suffice. However, for sinusoidal motion ψ = β sin(t) along a straight line, a simple rescaling is possible. The timeaveraging used by VeDBA follows that of ODBA, namely the L^{1}norm or arithmetic mean. However, that of the r.m.s. prescription follows the L^{2}norm or quadratic mean which, as found previously in (3.11), yields . If follows from the properties of the generalized means that will never be smaller than and its value is B22
For this simplified situation, one finds and hence B23
In applications where approximate alignment of the accelerometer's zaxis with gravity is attempted, small changes in orientation occurring at frequencies that are not rejected by the sampling window will cause ODBA to increase even in the absence of vibration. Because the magnitudes of a_{x} and a_{y} are then much smaller than gravity one finds that B24
Hence, the lack of rotational invariance in ODBA (as depicted in figure 10) then causes it to respond linearly to tilting of the accelerometer, tending to accentuate the sensitivity of ODBA to activity. Dependent on the moment of inertia there can be a cost associated with periodic adjustments in attitude, as might occur in animals exercising on a treadmill, and ODBA may be sensitive to it. For arbitrary orientations with respect to gravity, ODBA can both exaggerate and underestimate changes in the acceleration magnitude caused by rotation. Almost any acceleration measure is vulnerable to errors when the rotation of a transducer is not limited to low frequencies because it is then impossible to accurately separate the static and dynamic accelerations without information from a gyroscope.
When significant rotation exists, in the absence of a gyroscope, some measure of the variability of the instantaneous magnitude of the total acceleration vector a could offer an alternative proxy for biomechanical effort. For instance, one might calculate the standard deviation of for all the raw outputs of the accelerometer falling within some time interval. This could prove informative whenever a significant component of animal activity involves body rotation. Because accelerometers remain sensitive to gravity even within buoyant media, this may be especially useful in the context of aquatic animals as they are not constrained by gravity when swimming underwater. In circumstances where the dynamic acceleration is extremely small, activity can also be estimated by quantifying the rate of body reorientation dϕ/dt. This can be approximated as Δϕ/Δt using pairs of static acceleration vectors and separated by a fixed time interval Δt appropriate to the rotation rates of interest: B25and B26
 Received May 2, 2013.
 Accepted July 1, 2013.
© 2013 The Authors. Published by the Royal Society under the terms of the Creative Commons Attribution License http://creativecommons.org/licenses/by/3.0/, which permits unrestricted use, provided the original author and source are credited.