The Two-Thirds Power Law Derived from a Higher-Derivative Action

Abstract

The two-thirds power law is a link between angular speed $\omega$ and curvature $\kappa$ observed in voluntary human movements: $\omega$ is proportional to ${\kappa }^{2/3}$. Squared jerk is known to be a Lagrangian leading to the latter law. However, it leads to unbounded movements and is therefore incompatible with quasi-periodic dynamics, such as the movement of the tip of a pen drawing ellipses. To solve this drawback, we give a class of higher-derivative Lagrangians that allow for both quasi-periodic and unbounded movements, and at the same time lead to the two-thirds power law. The current study extends this framework and investigates a wider class of Lagrangians admitting generalised conservation laws.

1. Introduction

The Newtonian approach to mechanics relies on the fundamental equation $\overrightarrow{F}=m\overrightarrow{a}$ in Cartesian coordinates, where the total external force $\overrightarrow{F}$ acting on a point-like body of inertial mass m gives it an acceleration $\overrightarrow{a}=\ddot{\overrightarrow{x}}$ (with the dot on top denoting the time (t) derivative), which in turn yields the velocity $\overrightarrow{v}\left(t\right)=\dot{\overrightarrow{x}}$ and position $\overrightarrow{x}\left(t\right)$ of the body after successive integrations and upon imposing initial conditions. By invoking a principle of least action, Joseph-Louis Lagrange replaced Newton’s equations by the search for a single function L that is called the Lagrangian. In its original formulation, the Lagrangian L depends on the configuration (position) space variables, $q^{\alpha }\left(t\right)$, and their first-order time derivative ${\dot{q}}^{\alpha }\left(t\right)$, where $\alpha \in \{1,\dots ,n\}$ and the integer n is the dimension of the configuration space, i.e., the number of degrees of freedom of the system once all the holonomic constraints have been taken into account. Typically, $L=E_{\mathrm{k}}-E_{\mathrm{p}}$, the difference between the kinetic and potential energy functions of the system. Lagrange’s postulate is that the dynamical system, during its temporal evolution between two (initial and final) instants $t_{\mathrm{i}}$ and $t_{\mathrm{f}}$, minimises the action $S\left[q^{\alpha }\right]={\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}L(q^{\alpha },{\dot{q}}^{\alpha })dt$. The motion is then given by Hamilton’s variational principle $\delta S=0$ subjected to the boundary conditions $\delta q^{\alpha }\left(t_{\mathrm{i}}\right)=0=\delta q^{\alpha }\left(t_{\mathrm{f}}\right)$. Hamilton’s principle has then been generalised by Mikhail Ostrogradsky [1] to encompass Lagrangian functions depending on higher than one time derivatives of the dynamical variables $q^{\alpha }\left(t\right)$, therefore involving adapted boundary conditions. Note that by the term “Lagrangian”, one understands the function appearing under the integral in the expression of the action functional $S\left[q^{\alpha }\right]$, be it of the original form $L(q^{\alpha },{\dot{q}}^{\alpha })$ or in a higher-derivative form $L(q^{\alpha },{\dot{q}}^{\alpha },{\ddot{q}}^{\alpha },\dots )$. We say that an action is “higher-derivative” if it is given by the integral, along a given path in configuration space, of a higher-derivative Lagrangian, $L(q^{\alpha },{\dot{q}}^{\alpha },{\ddot{q}}^{\alpha },\dots )$.

In the realm of motor control, which resides at the intersection of biomechanics, neuroscience, and mathematics, there has been a paradigm shift akin to the one from Newtonian to Lagrangian mechanics. Instead of merely describing motion, researchers have begun to ask why organisms move the way they do. Indeed, certain actions or behaviours are repeated in a consistent manner. The concept of stereotyped movements was observed centuries ago, among others by Charles Scott Sherrington in 1906 [2]. The answer may lie in the minimisation of a cost function, a mathematical representation of the effort, energy, or some other metric that a biological system tries to optimise. The interested reader may find reviews about optimal control theory in motor control in Refs. [3,4]. By quantifying the cost linked with diverse trajectories or control strategies, these functions illuminate pathways to optimal movement strategies.

Traditionally applied to describe the dynamics of non-living physical systems, the recent introduction of Lagrangians as cost functions offers new insights accompanied by new challenges. One frequently cited cost function in motor control relies on what is called the jerk, $\overrightarrow{j}=\dot{\overrightarrow{a}}$, a vectorial quantity related to movement smoothness [5]. Research suggests that when humans make smooth, unperturbed movements, the trajectory they follow tends to minimise the integrated squared jerk, i.e., the movement minimises the following higher-derivative action [6]

$$S_J\left[\overrightarrow{x}\left(t\right)\right]={\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{L}_{J}dt={\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{\parallel \overrightarrow{j}\parallel }^{2}dt.$$

(1)

An apparent paradox that is going to be solved in the current study is that, from a mechanical point if view, the higher-derivative Lagrangian, $L_J={\parallel \overrightarrow{j}\parallel }^2$, does not lead to bounded trajectories. Indeed, the solutions $\overrightarrow{x}\left(t\right)$ of the Euler–Lagrange equations ${\overrightarrow{x}}^{\left(4\right)}=\overrightarrow{0}$ that minimise the action $S_J$ based on the Lagrangian $L_J$ are polynomials of degree three in the time parameter t; the polynomials are not periodic solutions. The jerk has however been considered in the study of bounded trajectories, which is problematic for describing voluntary movements such as the drawing of ellipses on a sheet of paper. In the present paper, we propose another framework in which higher-derivative actions produce bounded trajectories satisfying the two-thirds law (see below), and in which the minimisation of jerk appears to be a dynamical consequence of the variational principle. We reiterate that, by “higher-derivative”, we mean a Lagrangian $L({\overrightarrow{x}}^{\left(0\right)},{\overrightarrow{x}}^{\left(1\right)},\dots ,{\overrightarrow{x}}^{\left(N\right)})$ depending on the first $N>1$ derivatives of the dynamical variables $\overrightarrow{x}\left(t\right):={\overrightarrow{x}}^{\left(0\right)}\left(t\right)$, i.e., depending on ${\overrightarrow{x}}^{\left(0\right)}\left(t\right)$, ${\overrightarrow{x}}^{\left(1\right)}\left(t\right):=\dot{\overrightarrow{x}}\left(t\right)$, ${\overrightarrow{x}}^{\left(2\right)}\left(t\right):=\ddot{\overrightarrow{x}}\left(t\right)$, and so on. Three-dimensional vectors are assumed throughout. The study of such higher-derivative classical systems with finitely many degrees of freedom is achieved here by resorting to Ostrogradsky’s approach [1], which pioneered the study of higher-derivative actions.

Although cost functions can be central in elucidating optimal movement trajectories, cost functions are not the only framework in motor control research. Another approach is to identify conserved quantities or invariants. Invariants, as the term implies, are quantities whose value stays constant during dynamical evolution. A particularly illuminating example of this field is the empirical two-thirds power law, hereafter referred to as “2/3-PL”. It has been found in Ref. [7] and can be written as

$$v=C{\kappa }^{-1/3},$$

(2)

with velocity $\overrightarrow{v}=\dot{\overrightarrow{x}}$, speed $v=\parallel \overrightarrow{v}\parallel$ and $C\in {\mathbb{R}}_0^{+}$ a constant. One speaks of “two-thirds” because this law’s original formulation involved the angular speed $\omega =v\kappa$ where $\kappa =1/R$ is the curvature (the inverse of the radius R of the osculating circle) of the trajectory at a given time leading to $\omega =C{\kappa }^{2/3}$, i.e., the $2/3$ exponent. The initial observation was that the speed of a drawing or writing movement is related to the curvature of the drawing. This law has since been observed in a wide range of planar movements [8], and especially in elliptic trajectories. Note that a more general link of the type may be observed in an even wider class of trajectories [9], though here we do not consider the generalisation,

$$v=C{\kappa }^{\beta }.$$

(3)

There has been quite a debating on the 2/3-PL, with claims that the law (3) is mostly an artefact due to fitting, while others (including the authors of this paper) argue that this law is indeed a behavioural consequence of a fundamental law in human; see [8] for a detailed review on this topic.

Results from various studies challenge a purely kinematic interpretation of the 2/3-PL and highlight the role of the central nervous system in motion planning, which leads to Equation (2), using the motor imagery paradigms for covert movements [10,11].

The present investigation extends the understanding of cost functions through the lens of higher-derivative Lagrangians. The study proposes that a class of higher-derivative actions broader than Equation (1) leads to the 2/3-PL, providing new insights into cost functions in human motion. This is developed in Section 3 after general considerations about the 2/3-PL in Section 2.

Complementing the Lagrangian perspective, the Hamiltonian formalism provides a phase-space representation of the dynamics of a system as well as a general way to compute invariants through action-angle variables [12]. The proposed higher-derivative Lagrangians give rise to corresponding Hamiltonian functions through Ostrogradsky’s procedure [1], after which the invariants are computed in Section 4. A contextualisation of the results obtained in the framework of motor control is then discussed in Section 5.

2. Two-Thirds Law: Kinematical Considerations

In the context of three-dimensional motion in Euclidean space, the curvature and torsion of a given trajectory are obtained through standard formulae (see, e.g., [13,14])

$$\kappa ={\displaystyle \frac{\parallel \overrightarrow{v}\times \overrightarrow{a}\parallel }{v^3}},\tau ={\displaystyle \frac{\overrightarrow{j}·(\overrightarrow{v}\times \overrightarrow{a})}{{\parallel \overrightarrow{v}\times \overrightarrow{a}\parallel }^2}},$$

(4)

where $\overrightarrow{a}=\ddot{\overrightarrow{x}}$ and where the symbol × stands for the vector product in three-dimensional space. From the first Equation (4) giving the curvature of the trajectory in terms of the velocity and acceleration, one can write

$$v={\parallel \overrightarrow{v}\times \overrightarrow{a}\parallel }^{1/3}{\kappa }^{-1/3}.$$

(5)

Hence, the 2/3-PL is valid if and only if the norm ${\ell }_2$ of the vector ${\overrightarrow{\ell }}_2=\overrightarrow{v}\times \overrightarrow{a}$ is constant. A sufficient although not necessary condition is ${\dot{\overrightarrow{\ell }}}_2=\overrightarrow{0}$. Using ${\dot{\overrightarrow{\ell }}}_2=\overrightarrow{a}\times \overrightarrow{a}+\overrightarrow{v}\times \overrightarrow{j}$ along with $\overrightarrow{a}\times \overrightarrow{a}\equiv 0$, one can find that for Equation (2) to hold true, it suffices to have

$$\overrightarrow{j}=\gamma \overrightarrow{v}.$$

(6)

The coefficient function $\gamma$ may explicitly depend on time and on the various derivatives ${\overrightarrow{x}}^{\left(i\right)}$. Since $\overrightarrow{j}={\overrightarrow{x}}^{\left(3\right)}$ and $\overrightarrow{v}={\overrightarrow{x}}^{\left(1\right)}$, Equation (6) can be rewritten as ${\overrightarrow{x}}^{\left(3\right)}-\gamma ({\overrightarrow{x}}^{\left(i\right)},t){\overrightarrow{x}}^{\left(1\right)}=\overrightarrow{0}$. That is, a trajectory satisfying this last relation also satisfies Equation (2). The simplest choice is that of a constant function $\gamma$ leading to elliptic trajectories. In Ref. [15], different choices of the form $\gamma =\gamma \left(t\right)$ are explored leading to trajectories that all comply with the 2/3-PL.

Let us note that a straightforward way to obtain dynamics leading to a conserved norm ${\ell }_2$ of the vector ${\overrightarrow{\ell }}_2$ is proposed in Ref. [16]. The Lagrangian considered is actually proportional to ${\ell }_2=\parallel \overrightarrow{v}\times \overrightarrow{a}\parallel$ itself. Since the Lagrangian does not explicitly depend on the evolution parameter, via Noether’s theorem, there exists the conserved quantity $B:=\overrightarrow{v}·{\overrightarrow{p}}_0+\overrightarrow{a}·{\overrightarrow{p}}_1-L$, where, in this case, ${\overrightarrow{p}}_0={\displaystyle \frac{\partial L}{\partial \overrightarrow{v}}}-{\displaystyle \frac{d}{dt}}{\displaystyle \frac{\partial L}{\partial \overrightarrow{a}}}$ and ${\overrightarrow{p}}_1={\displaystyle \frac{\partial L}{\partial \overrightarrow{a}}}$. In the case when $L={\ell }_2$, one has $B=2L=2{\ell }_2$, ensuring a conservation of the norm of the vector ${\overrightarrow{\ell }}_2$, which we recall is a sufficient and necessary condition for the 2/3-PL (see Equation (5)). However, the dynamics associated with the Lagrangian $L={\ell }_2$ is ${\overrightarrow{x}}^{\left(3\right)}=0$ leading to unbounded trajectories in phase-space, as is the case for the jerk Lagrangian $L_J$ discussed in Section 1.

Condition (6) implies that the motion is planar, so that the torsion $\tau$ defined in Equation (4) vanishes. Hence, non-planar trajectories should not be related to Equation (2). Although it is not the main topic of the current paper, it has been proposed in Ref. [17] that the law

$$v=C{\kappa }^{-1/3}{\tau }^{-1/6}$$

(7)

should hold in non-planar, three-dimensional motion, with C being a positive constant. From Equation (4), it can be deduced that $C={\parallel \overrightarrow{j}·(\overrightarrow{v}\times \overrightarrow{a})\parallel }^{1/6}$. In other words, the law (7) is valid if $\overrightarrow{j}·(\overrightarrow{v}\times \overrightarrow{a})$ is constant, or

$$\dot{\overrightarrow{j}}=f\overrightarrow{v}.$$

(8)

Since $\dot{\overrightarrow{j}}={\overrightarrow{x}}^{\left(4\right)}$ and $\overrightarrow{v}={\overrightarrow{x}}^{\left(1\right)}$, Equation (8) can be rewritten as

$${\overrightarrow{x}}^{\left(4\right)}-f({\overrightarrow{x}}^{\left(i\right)},t){\overrightarrow{x}}^{\left(1\right)}=\overrightarrow{0}.$$

(9)

Hence, any trajectory such that given by Equation (9) holds satisfies the law (7).

3. Dynamical Principle—Lagrangian Formalism

3.1. The Model

We propose that the actions associated with the higher-derivative Lagrangians are

$$L={\displaystyle \frac{\lambda }{2}}{\parallel {\overrightarrow{x}}^{\left(N\right)}\parallel }^2-{\displaystyle \frac{1}{2}}U\left({\parallel {\overrightarrow{x}}^{(N-1)}\parallel }^2\right),$$

(10)

where $N\ge 1$ and $\lambda \in {\mathbb{R}}_0^{+}$ are relevant cost function candidates that may lead to trajectories satisfying (2) for any “potential” function $U\left(z\right)$ where the variable z denotes the squared norm ${\parallel {\overrightarrow{x}}^{(N-1)}\parallel }^2$. Recall that the Lagrangian (10) is higher-derivative as soon as $N>1$.

From the variational principle based on the action functional $S\left[\overrightarrow{x}\left(t\right)\right]={\int }_{t_{{\prime }^{\mathrm{i}}}}^{t_{\mathrm{f}}}Ldt$, the equations of motion read as follows, for a generic higher-derivative Lagrangian [1,18]:

$$\overrightarrow{0}={\displaystyle \frac{\delta L}{\delta \overrightarrow{x}}}\equiv \sum _{j=0}^N{\left(-{\displaystyle \frac{d}{dt}}\right)}^j{\displaystyle \frac{\partial L}{\partial {\overrightarrow{x}}^{\left(j\right)}}}.$$

(11)

Equations (11) have to be satisfied along with the vanishing of the boundary terms defining the momenta ${\overrightarrow{p}}_i$ [1,18]:

$$\begin{array}{ccc}{\left.\sum _{i=0}^{N-1}\delta {\overrightarrow{x}}^{\left(i\right)}{\overrightarrow{p}}_i\right|}_{t_1}^{t_2}\hfill & =& 0,\hfill \end{array}$$

(12)

$$\begin{array}{ccc}\hfill {\overrightarrow{p}}_i& :=& {\displaystyle \frac{\delta L}{\delta {\overrightarrow{x}}^{(i+1)}}}\equiv \sum _{j=0}^{N-i-1}{\left(-{\displaystyle \frac{d}{dt}}\right)}^j{\displaystyle \frac{\partial L}{\partial {\overrightarrow{x}}^{(i+j+1)}}},i\in \{0,\dots ,N-1\}.\hfill \end{array}$$

(13)

So far, no assumption has been made to derive Equations (11) and (12). These equations are mere consequences of Hamilton’s variational principle, $\delta S=0$, for a generic higher-derivative action. One chooses to cancel the boundary terms (12) by imposing the following conditions at the boundaries of the integration domain:

$$\delta {\overrightarrow{x}}^{\left(j\right)}\left(t_{\mathrm{f}}\right)=0=\delta {\overrightarrow{x}}^{\left(j\right)}\left(t_{\mathrm{i}}\right),\forall j\in \{0,\dots ,N-1\}.$$

Since the Lagrangians under discussion do not explicitly depend on time, let us choose the initial date $t_{\mathrm{i}}=0$ for convenience. More specifically, the equations of motion computed from the Lagrangian (10) are

$$\lambda {\displaystyle \frac{d^N}{d{t}^N}}{\overrightarrow{x}}^{\left(N\right)}=-{\displaystyle \frac{d^{N-1}}{d{t}^{N-1}}}\left[U^{\prime }({\parallel {\overrightarrow{x}}^{(N-1)}\parallel }^2){\overrightarrow{x}}^{(N-1)}\right],$$

(14)

where the prime denotes the derivative over the squared norm ${\parallel {\overrightarrow{x}}^{(N-1)}\parallel }^2$.

Integrating Equation (14) ($N-1$) times leads to

$$\lambda {\overrightarrow{x}}^{(N+1)}=-U^{\prime }({\overrightarrow{x}}^{(N-1)2}){\overrightarrow{x}}^{(N-1)}+\sum _{j=0}^{N-2}{\displaystyle \frac{{\overrightarrow{b}}_j}{j!}}{t}^j,$$

(15)

with ($N-1$) constant vectors $\{{\overrightarrow{b}}_j\}$, $j=0,\dots ,N-2$, that can be fixed by initial conditions. As one can see from the last terms on the right-hand side of Equation (15), the Lagrangian field equation (14) leads to general solutions with a polynomial dependence on the evolution parameter, therefore signalling an instability that is a landmark of higher-derivative models.

One can find that in order to avoid any instabilities—also called run-away solutions in the context of field theory—one has to impose the initial conditions

$${\overrightarrow{b}}_j=\overrightarrow{0},\forall j\in \{0,\dots ,N-2\},$$

(16)

as soon as then, the field Equation (14) leads to the differential equation

$$\lambda {\displaystyle \frac{d}{dt}}{\overrightarrow{x}}^{\left(N\right)}+U^{\prime }({\overrightarrow{x}}^{(N-1)2}){\overrightarrow{x}}^{(N-1)}=0$$

(17)

that provides stable solutions and may imply the 2/3-PL in the cases discussed below. Note that Equation (16) only fixes ($N-1$) conditions among the $2N$ initial conditions needed to ensure a unique solution.

An equivalent way of presenting the initial conditions (16) is by first defining the ($N-1$) vectors

$${\overrightarrow{A}}_j^N:=\lambda {\displaystyle \frac{d^{j+1}}{d{t}^{j+1}}}{\overrightarrow{x}}^{\left(N\right)}+{\displaystyle \frac{d^j}{d{t}^j}}[U^{\prime }({\parallel {\overrightarrow{x}}^{(N-1)}\parallel }^2){\overrightarrow{x}}^{(N-1)}],j\in \{0,\dots ,N-2\},$$

(18)

and then setting the initial conditions

$${\left.{\overrightarrow{A}}_j^N\right|}_{t=0}\equiv {\overrightarrow{b}}_j=\overrightarrow{0},\forall j\in \{0,\dots ,N-2\}.$$

(19)

The point with the presentation (19) of the initial conditions (16) is that, for the considered Lagrangian (10), the set of vectors $\{{\overrightarrow{A}}_j^N\}$, $j\in \{0,\dots ,N-2\}$, is in one-to-one correspondence with the set of momenta $\{{\overrightarrow{p}}_j\}$, $j\in \{0,\dots ,N-2\}$; see Equation (13). Indeed, one readily sees from the Lagrangian (10) and the definition (13) of the momenta that

$${\overrightarrow{A}}_j^N={(-1)}^{j-1}{\overrightarrow{p}}_{N-2-j},j=0,\dots ,N-2.$$

(20)

Therefore, the choice of initial conditions (16), equivalently the conditions (19), is in turn equivalent to setting the following initial conditions on the first ($N-1$) momenta:

$${\overrightarrow{p}}_j(t=0)=\overrightarrow{0},j=0,1,\dots ,N-2.$$

(21)

To summarise, with the Lagrangian (10) and the above initial conditions (19), the Lagrangian field equations (14) lead to the differential equation

$$\lambda {\displaystyle \frac{d}{dt}}{\overrightarrow{x}}^{\left(N\right)}+U^{\prime }({\parallel {\overrightarrow{x}}^{(N-1)}\parallel }^2){\overrightarrow{x}}^{(N-1)}=0.$$

(22)

From Equation (22), one observes that the vector

$${\overrightarrow{\ell }}_N=\lambda {\overrightarrow{x}}^{(N-1)}\times {\overrightarrow{x}}^{\left(N\right)}$$

(23)

remains constant during the dynamical evolution, i.e., ${\dot{\overrightarrow{\ell }}}_N=\overrightarrow{0}$. Conversely, from the variational equations (14) and their consequence (15), one sees that imposing the condition ${\dot{\overrightarrow{\ell }}}_N=\overrightarrow{0}$ is equivalent to imposing the initial conditions (19). Therefore, a compact way of imposing the initial conditions (19) is by imposing ${\dot{\overrightarrow{\ell }}}_N=\overrightarrow{0}$.

Now, let us write the squared norm of the vector ${\overrightarrow{\ell }}_N$ (23) in a way that simplifies the expressions for the Lagrangian and the Hamiltonian. To this end, at any instant t, one decomposes the vector ${\overrightarrow{x}}^{(N-1)}$ into its norm $x^{(N-1)}$ times a unit vector ${\overrightarrow{1}}_t^{(N-1)}$. In the formula, one has ${\overrightarrow{x}}^{(N-1)}=x^{(N-1)}{\overrightarrow{1}}_t^{(N-1)}$. This allows one to write

$${\overrightarrow{x}}^{\left(N\right)}\equiv {\dot{\overrightarrow{x}}}^{(N-1)}={\dot{x}}^{(N-1)}{\overrightarrow{1}}_t^{(N-1)}+k_{(N-1)}{x}^{(N-1)}{\overrightarrow{1}}_n^{(N-1)},$$

(24)

where we use ${\displaystyle \frac{d}{dt}}{\overrightarrow{1}}_t^{(N-1)}=k_{(N-1)}{\overrightarrow{1}}_n^{(N-1)}$ in terms of some function $k_{(N-1)}$; we do not need to specify and of a unit vector ${\overrightarrow{1}}_n^{(N-1)}$ orthogonal to ${\overrightarrow{1}}_t^{(N-1)}$, as one straightforwardly verifies by taking the derivative of the scalar product ${\overrightarrow{1}}_t^{(N-1)}·{\overrightarrow{1}}_t^{(N-1)}=1$. From the decompositions (24) of the vectors ${\overrightarrow{x}}^{(N-1)}$ and ${\overrightarrow{x}}^{\left(N\right)}$, one readily finds ${\ell }_N^2={\lambda }^2{k}_{(N-1)}^2{\left(x^{(N-1)}\right)}^4$. In turn, this allows us to rewrite the Lagrangian (10) in the form

$$L={\displaystyle \frac{\lambda }{2}}{\left({\dot{x}}^{(N-1)}\right)}^2+{\displaystyle \frac{{\ell }_N^2}{2\lambda {\left(x^{(N-1)}\right)}^2}}-{\displaystyle \frac{1}{2}}U\left({\left(x^{(N-1)}\right)}^2\right)\equiv L(x^{(N-1)},{\dot{x}}^{(N-1)}).$$

(25)

Equation (25) is the higher-derivative generalisation of the reduced Lagrangian for a dynamical system with spherical symmetry.

The question of finding a Lagrangian whose equations of motion are Equations (9) arises naturally but, at this stage, it remains an open problem since we have not found any Lagrangian leading to Equation (9). In what follows, we choose not to investigate the case of non-planar trajectories further.

3.2. The Case $N=2$

Equations (6) and (22) surely coincide if $N=2$. In this case, the Lagrangian (10) reduces to the Flash–Handzel Lagrangian [19], i.e.,

$$L={\displaystyle \frac{\lambda }{2}}{\parallel \overrightarrow{a}\parallel }^2-{\displaystyle \frac{1}{2}}U({\parallel \overrightarrow{v}\parallel }^2).$$

(26)

If one imposes the initial condition ${\overrightarrow{b}}_0:=\lambda \overrightarrow{j}\left(0\right)+\overrightarrow{v}\left(0\right)U^{\prime }\left(v^2\left(0\right)\right)=\overrightarrow{0}$, the vector ${\overrightarrow{\ell }}_2=\lambda \overrightarrow{v}\times \overrightarrow{a}$ is conserved since Equation (22) now reads

$$\lambda \overrightarrow{j}=-U^{\prime }\left(v^2\right)\overrightarrow{v},$$

(27)

where we use the notation $V^2={\parallel \overrightarrow{V}\parallel }^2$ for the squared norm of vectors. The trajectories generated by any choice of U therefore satisfies 2/3-PL. Note that the velocity at initial time should not vanish, since in that case, the conserved vector ${\overrightarrow{\ell }}_2$ would vanish at all ulterior times leading to a vanishing constant of proportionality between the norm of the velocity and the curvature $\kappa$ to the power $p=-1/3$ (compare Equation (2) with Equation (5)), thereby enforcing a straight line as the only possible movement. For a movement with non-zero curvature $\kappa$, it is therefore crucial that the initial condition ${\overrightarrow{b}}_0=\overrightarrow{0}$ combines the two vectors $\overrightarrow{j}\left(0\right)$ and $\overrightarrow{v}\left(0\right)$ in such a way as to cancel their weighted sum, $\lambda \overrightarrow{j}\left(0\right)+U^{\prime }\left(v^2\left(0\right)\right)\overrightarrow{v}\left(0\right)$, without cancelling them separately.

3.3. The Case of a Linear Potential Function

Keeping N arbitrary but setting $U\left(z\right)=\lambda {\omega }^{2}z$ leads to Pais–Uhlenbeck oscillators [18] since Equation (10) now reads

$$L_{PU}={\displaystyle \frac{\lambda }{2}}\left({\overrightarrow{x}}^{\left(N\right)2}-{\omega }^2{\overrightarrow{x}}^{(N-1)2}\right).$$

(28)

After ($2N-2$) integrations, the equations of motions (14) reduce to

$$\ddot{\overrightarrow{x}}=-{\omega }^2\overrightarrow{x}+\sum _{j=0}^{2N-3}{\displaystyle \frac{{\overrightarrow{b}}_j}{j!}}{t}^j$$

(29)

whose general solution is

$$\overrightarrow{x}\left(t\right)={\overrightarrow{K}}_{1}cos\left(\omega t\right)+{\overrightarrow{K}}_{2}sin\left(\omega t\right)+\sum _{j=0}^{2N-3}{\overrightarrow{a}}_j{\displaystyle \frac{t^j}{j!}},$$

(30)

with ${\overrightarrow{a}}_{2N-3}={\overrightarrow{b}}_{2N-3}/{\omega }^2$, ${\overrightarrow{a}}_{2N-4}={\overrightarrow{b}}_{2N-4}/{\omega }^2$, the other ${\overrightarrow{a}}_j$ with $j\le 2N-5$ being recursively given by

$${\omega }^2{\overrightarrow{a}}_j+{\overrightarrow{a}}_{j+2}={\overrightarrow{b}}_j,j\in \{0,1,\dots ,2N-5\}.$$

(31)

Imposing the initial conditions ${\overrightarrow{b}}_j=\overrightarrow{0}$ for all j certainly implies that ${\overrightarrow{a}}_j=\overrightarrow{0}$ for all j and the dynamics reduce to elliptic trajectories $\overrightarrow{x}\left(t\right)={\overrightarrow{K}}_{1}cos\left(\omega t\right)+{\overrightarrow{K}}_{2}sin\left(\omega t\right)$ without any run-away modes. Indeed, with these initial conditions, the equations of motion reduce to $\ddot{\overrightarrow{x}}+{\omega }^2\overrightarrow{x}=0$; hence, ${\overrightarrow{\ell }}_2={\omega }^2{\overrightarrow{\ell }}_1$ and the angular momentum ${\overrightarrow{\ell }}_1$ is conserved, as in standard Newtonian mechanics with central forces.

3.4. The Case of a Vanishing Potential Function

In the case when $U\left(z\right)=0$, one recovers the mean squared (ms) derivative cost functions of [20]:

$$L_{\mathrm{ms}}={\displaystyle \frac{\lambda }{2}}{\parallel {\overrightarrow{x}}^{\left(N\right)}\parallel }^2,$$

(32)

whose equations of motion are ${\overrightarrow{x}}^{\left(2N\right)}=\overrightarrow{0}$ leading to

$$\overrightarrow{x}\left(t\right)=\sum _{j=0}^{2N-1}{\displaystyle \frac{{\overrightarrow{b}}_j}{j!}}{t}^j.$$

(33)

Such trajectories are unbounded and do not satisfy Equation (6). Still, Lagrangians of this type have been shown to successfully model pointing tasks since the seminal paper [21] using $N=3$.

4. Hamiltonian Formalism

4.1. Ostrogradsky’s Approach

Following Ref. [1], the position degrees of freedom are defined as ${\overrightarrow{q}}^j={\overrightarrow{x}}^{\left(j\right)}$ and the momenta ${\overrightarrow{p}}_j$ are defined by Equation (13) for $j=0,\dots ,N-1$. Specifically, for the Lagrangian (10), one has

$$\begin{array}{ccc}\hfill {\overrightarrow{p}}_j& =& {(-1)}^{N-j-1}\lambda {\overrightarrow{x}}^{(2N-j-1)}+{(-1)}^{N-j-1}{\displaystyle \frac{d^{N-j-2}}{d{t}^{N-j-2}}}\left(U^{\prime }({\overrightarrow{x}}^{(N-1)2}){\overrightarrow{x}}^{(N-1)}\right)\hfill \end{array}$$

(34)

$$\begin{array}{ccc}& =& {(-1)}^{N-j-1}{\displaystyle \frac{d^{N-j-2}}{d{t}^{N-j-2}}}{\overrightarrow{A}}_0^N\equiv {(-1)}^{N-j-1}{\overrightarrow{A}}_{N-j-2}^N,\hfill \end{array}$$

(35)

and the corresponding Hamiltonian reads

$$H=\sum _{j=0}^{N-2}{\overrightarrow{p}}_j·{\overrightarrow{q}}^{j+1}+{\displaystyle \frac{{\overrightarrow{p}}_{N-1}^2}{2\lambda }}+{\displaystyle \frac{1}{2}}U\left({({\overrightarrow{q}}^{N-1})}^2\right).$$

(36)

Hamilton’s equations ${\dot{\overrightarrow{q}}}^j={\displaystyle \frac{\partial H}{\partial {\overrightarrow{p}}_j}}$, ${\dot{\overrightarrow{p}}}_j=-{\displaystyle \frac{\partial H}{\partial {\overrightarrow{q}}^j}}$ lead to

$$\begin{array}{ccc}\hfill {\dot{\overrightarrow{q}}}^0& =& {\overrightarrow{q}}^1,\hfill \\ & ⋮& \end{array}$$

(37)

$$\begin{array}{ccc}\hfill {\dot{\overrightarrow{q}}}^{N-2}& =& {\overrightarrow{q}}^{N-1},\hfill \end{array}$$

(38)

$$\begin{array}{ccc}\hfill {\dot{\overrightarrow{q}}}^{N-1}& =& {\displaystyle \frac{{\overrightarrow{p}}_{N-1}}{\lambda }},\hfill \end{array}$$

(39)

$$\begin{array}{ccc}\hfill {\dot{\overrightarrow{p}}}_0& =& \overrightarrow{0},\hfill \end{array}$$

(40)

$$\begin{array}{ccc}\hfill {\dot{\overrightarrow{p}}}_1& =& -{\overrightarrow{p}}_0,\hfill \\ & ⋮& \end{array}$$

(41)

$$\begin{array}{ccc}\hfill {\dot{\overrightarrow{p}}}_{N-2}& =& -{\overrightarrow{p}}_{N-3},\hfill \end{array}$$

(42)

$$\begin{array}{ccc}\hfill {\dot{\overrightarrow{p}}}_{N-1}& =& -{\overrightarrow{p}}_{N-2}-U^{\prime }\left({({\overrightarrow{q}}^{N-1})}^2\right){\overrightarrow{q}}^{N-1}.\hfill \end{array}$$

(43)

On the conserved phase-space surface given by

$${\overrightarrow{p}}_0=\cdots ={\overrightarrow{p}}_{N-2}=\overrightarrow{0},$$

(44)

and by using Equation (39) along with Equation (24), one obtains the reduced Hamiltonian,

$$\tilde{H}(q^{N-1},p_{N-1})={\displaystyle \frac{p_{N-1}^2}{2\lambda }}+{\displaystyle \frac{{\ell }_N^2}{2\lambda {\left(q^{N-1}\right)}^2}}+{\displaystyle \frac{1}{2}}U\left({\left(q^{N-1}\right)}^2\right),$$

(45)

with $q^{N-1}=\parallel {\overrightarrow{q}}^{N-1}\parallel$ and ${\ell }_N$ is the norm of the constant vector ${\overrightarrow{\ell }}_N$ given in Equation (23). Let us recall that the conditions (44) are equivalent to the conditions (19).

The trajectory in the plane $(q^{N-1},p_{N-1})$ is such that $\tilde{H}=\tilde{\mathcal{E}}$ is constant—we use the notation $\tilde{\mathcal{E}}$ although the Hamiltonian does not a priori possess the dimension of an energy. For example, in the Pais–Uhlenbeck case discussed in Section 3.3, the trajectory in the plane $(q^{N-1},p_{N-1})$ is a closed loop whose equation is given by

$$p_{N-1}^2+{\displaystyle \frac{{\ell }_N^2}{{\left(q^{N-1}\right)}^2}}+{\lambda }^2{\omega }^2{\left(q^{N-1}\right)}^2=2\lambda \tilde{\mathcal{E}},q^{N-1}>0.$$

(46)

A stability analysis of trajectories can be carried out by applying the method of Ref. [22] (Chapter 7) to the Hamiltonian Equations (37)–(43). Equations (37)–(43) can be rewritten in a matrix form: $\dot{\zeta }=B\zeta$, where $\zeta$ is a vector containing the $2N$ coordinates in phase-space and B is a $2N$×$2N$ matrix. The eigenvalues of B are 0, with multiplicity $(2N-2)$, and $\pm \sqrt{U^{\prime }\left(0\right)/\lambda }$. The zero modes are global translation modes. If $U^{\prime }\left(0\right)<0$ as in the harmonic oscillator case, the non-zero eigenvalues are complex conjugated and all trajectories are bounded, global translation excepted. If $U^{\prime }\left(0\right)>0$, there necessarily exist unbounded trajectories even if the global translation mode is set to zero. The existence of unbounded trajectories preserving the 2/3-PL is a prediction that could be experimentally studied and is out of the scope of this paper.

4.2. Action Variables

Let us now restrict the discussion to potentials U for which there exist some values of $\tilde{\mathcal{E}}$ such that the equation $\tilde{\mathcal{E}}={\displaystyle \frac{{\ell }_N^2}{2\lambda {\left(q^{N-1}\right)}^2}}+{\displaystyle \frac{1}{2}}U[{(q^{N-1})}^2]$ has two finite, distinct solutions for $q^{N-1}>0$ leading to bounded trajectories that will appear as closed loops $\mathrm{\Gamma }$ in the plane $(q^{N-1},p_{N-1})$. The two distinct solutions are denoted $q_{\mathrm{i}}^{N-1}$ and $q_{\mathrm{f}}^{N-1}$, where one can choose $q_{\mathrm{i}}^{N-1}<q_{\mathrm{f}}^{N-1}$ without loss of generality. These correspond to the two turning points of the radial motion in the three-dimensional space with coordinates $(q_x^{N-1},q_y^{N-1},q_z^{N-1})$. The two turning points $q_{\mathrm{i}}^{N-1}$ and $q_{\mathrm{f}}^{N-1}$ are reached at the successive dates $t_{\mathrm{i}}$ and $t_{\mathrm{f}}$, respectively.

It is known that the action variable $I_N={\displaystyle \frac{1}{2\pi }}{\oint }_{\mathrm{\Gamma }}{p}_{N-1}d{q}^{N-1}$ is a constant of motion [12]. $I_N$ can be expressed in various ways, where in all cases one chooses to parametrise the closed curve $\mathrm{\Gamma }$ by the evolution parameter t. First, by using Equation (39),

$$I_N={\displaystyle \frac{1}{\pi }}{\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{p}_{N-1}{\dot{q}}^{N-1}dt={\displaystyle \frac{1}{\pi \lambda }}{\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{p}_{N-1}^{2}dt={\displaystyle \frac{\lambda }{\pi }}{\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{\left({\dot{q}}^{N-1}\right)}^{2}dt,$$

(47)

where ($t_{\mathrm{f}}-t_{\mathrm{i}}$) is half the period of the radial movement. By starting from the equivalent expression $I_N=-{\displaystyle \frac{1}{2\pi }}{\oint }_{\mathrm{\Gamma }}{q}^{N-1}d{p}_{N-1}$ and using the solution (17), and equivalently Equation (39) with Equation (43), one obtains

$$I_N={\displaystyle \frac{{\lambda }^2}{\pi }}{\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{\displaystyle \frac{{\left({\ddot{q}}^{N-1}\right)}^2}{U^{\prime }\left(q^{(N-1)2}\right)}}dt.$$

(48)

For potential functions U such that $U^{\prime }\left(0\right)\ne 0$, let us define

$$I_N^a={\displaystyle \frac{{\lambda }^2}{\pi U^{\prime }\left(0\right)}}{\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{\left({\ddot{q}}^{N-1}\right)}^{2}dt$$

(49)

for convenience.

In the Pais–Uhlenbeck case, $I_N\equiv I_N^a={\displaystyle \frac{\lambda }{\pi {\omega }^2}}{\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{\left({\ddot{q}}^{N-1}\right)}^{2}dt$ and $\tilde{H}=\omega (2I_N+{\ell }_N)$ [22] (Chapter 6); see also [23] for a discussion of action variables in the generic Pais–Uhlenbeck oscillator model. It can be expected that $I_N\approx I_N^a$ from the Taylor expansion $U\left(z\right)\approx U\left(0\right)+U^{\prime }\left(0\right)z+\dots$. It is worth pointing out that the mean squared derivative cost functions used in Ref. [20] are equal to $I_N^a$.

4.3. Application to the Pais–Uhlenbeck $N=2$ Oscillator

We apply the formalism exposed so far to the simplest higher-derivative case with $N=2$ and a linear potential, giving a Pais–Uhlenbeck Lagrangian $L_{PU}(\overrightarrow{x},\dot{\overrightarrow{x}},\ddot{\overrightarrow{x}})={\displaystyle \frac{\lambda }{2}}\left(a^2-{\omega }^2{v}^2\right)$. The Ostrogradsky momenta $({\overrightarrow{p}}_0,{\overrightarrow{p}}_1)$ read ${\overrightarrow{p}}_0=-\lambda (\overrightarrow{j}+{\omega }^2\overrightarrow{v})$ and ${\overrightarrow{p}}_1=\lambda \overrightarrow{a}$. With the notation $({\overrightarrow{q}}^0,{\overrightarrow{q}}^1)=(\overrightarrow{x},\overrightarrow{v})$, the expression for ${\overrightarrow{p}}_1$ allows us to express the acceleration as $\overrightarrow{a}={\dot{\overrightarrow{q}}}^1={\overrightarrow{p}}_1/\lambda$.

The Ostrogradsky Hamiltonian reads

$$\begin{array}{cc}\hfill H({\overrightarrow{q}}^0,{\overrightarrow{q}}^1,{\overrightarrow{p}}_0,{\overrightarrow{p}}_1)& ={\overrightarrow{p}}_0·{\overrightarrow{q}}^1+{\displaystyle \frac{1}{\lambda }}{\overrightarrow{p}}_1·{\overrightarrow{p}}_1-\frac{\lambda }{2}\left(\frac{p_1^2}{{\lambda }^2}-{\omega }^2{\left(q^1\right)}^2\right)\hfill \\ \hfill & ={\overrightarrow{p}}_0·{\overrightarrow{q}}^1+\frac{1}{2\lambda }{p}_1^2+\frac{\lambda {\omega }^2}{2}{\left(q^1\right)}^2.\hfill \end{array}$$

(50)

Hamilton’s equations ${\dot{\overrightarrow{q}}}^j={\displaystyle \frac{\partial H}{\partial {\overrightarrow{p}}_j}}$, ${\dot{\overrightarrow{p}}}_j=-{\displaystyle \frac{\partial H}{\partial {\overrightarrow{q}}^j}}$ for $j\in \{0,1\}$ yield

$${\dot{\overrightarrow{q}}}^0={\overrightarrow{q}}^1,{\dot{\overrightarrow{q}}}^1={\displaystyle \frac{{\overrightarrow{p}}_1}{\lambda }},{\dot{\overrightarrow{p}}}_0=0,{\dot{\overrightarrow{p}}}_1=-\lambda {\omega }^2{\overrightarrow{q}}^1.$$

(51)

Altogether, Hamilton’s equations (51) lead to the field equations ${\overrightarrow{x}}^{\left(3\right)}+{\omega }^2{\overrightarrow{x}}^{\left(1\right)}=0$, with the solution

$$\overrightarrow{x}\left(t\right)={\overrightarrow{K}}_{1}cos\left(\omega t\right)+{\overrightarrow{K}}_{2}sin\left(\omega t\right)+{\overrightarrow{a}}_0,$$

(52)

with ${\overrightarrow{K}}_1$, ${\overrightarrow{K}}_2$, ${\overrightarrow{a}}_0\in {\mathbb{R}}^3$. On the conserved surface ${\overrightarrow{p}}_0=\overrightarrow{0}$ in phase-space, the reduced Hamiltonian reads

$$\tilde{H}(q^1,p_1)={\displaystyle \frac{p_1^2}{2\lambda }}+{\displaystyle \frac{{\ell }_2^2}{2\lambda {\left(q^1\right)}^2}}+{\displaystyle \frac{\lambda {\omega }^2}{2}}{\left(q^1\right)}^2=\tilde{\mathcal{E}},$$

(53)

similarly to the more general case Hamiltonian (45). Finally, $I_2=I_2^a={\displaystyle \frac{\lambda }{\pi {\omega }^2}}{\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{q}^{\left(3\right)2}dt$ as a constant of motion and $\tilde{H}=\omega (2I_2+{\ell }_2)$. Solution (52) is graphically illustrated in Figure 1A, while the phase-space trajectory computed from Hamiltonian (53) is illustrated in Figure 1B.

5. Conclusions

5.1. Variational Principle for 2/3-PL

The application of Lagrangian and Hamiltonian formalisms in motor control demands a multidisciplinary approach that respects both an established mathematical formalism and the intricacies of human physiology. This paper proposed a broader class of higher-derivative Lagrangians that, upon defining appropriate initial conditions, lead to trajectories complying with the 2/3-PL, thus providing new insights into cost functions critical to human motion. These Lagrangians are given in Equation (10). A salient issue for the observation of the 2/3-PL is the necessity of setting accurate initial conditions: If a Lagrangian involving up to the $N\mathrm{th}$ time derivative ${\overrightarrow{x}}^{\left(N\right)}$ is used, it is necessary, in order to have a bounded motion, to set the Ostrogradsky momenta ${\overrightarrow{p}}_j\left(0\right)=\overrightarrow{0}$, with $0\le j\le N-2$. For the class of Lagrangians (10), the momenta are given by Equation (34). If these ($N-1$) initial conditions are unaligned with the natural capabilities of the human motor system, the considered Lagrangian is not qualified to model voluntary human movement. This consideration leads us to the conclusion that a minimal N is the most natural choice. Therefore, we consider that $N=2$ actions of the form

$$S={\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}\left[{\displaystyle \frac{\lambda }{2}}{\overrightarrow{a}}^2-{\displaystyle \frac{1}{2}}U\left({\overrightarrow{v}}^2\right)\right]dt$$

(54)

are favoured:

• The actions (54) naturally lead to the 2/3-PL, provided one individual is able to fix the initial condition ${\overrightarrow{b}}_0:=\lambda \overrightarrow{j}\left(0\right)+U^{\prime }({\overrightarrow{v}}^2\left(0\right))\overrightarrow{v}\left(0\right)=\overrightarrow{0}$. For a motion with vanishing initial speed, one just needs to impose $\overrightarrow{j}\left(0\right)=\overrightarrow{0}$, irrespective of the choice of potential function U.
• The actions (54) may lead to quite a large variety of trajectories satisfying the 2/3-PL according to the choice made for U. In the case of harmonic potential, elliptic trajectories are recovered, which are the best-known case in which this law appears.
• The action variable $I_2={\displaystyle \frac{\lambda }{2\pi }}{\oint }_{\mathrm{\Gamma }}adv\sim {\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{\parallel \overrightarrow{j}\parallel }^{2}dt$ makes the mean squared jerk explicitly appear, and it is known that minimising this function (maximising smoothness) is an experimentally observed principle in motor control [19]. Mechanics imposes that $I_2$ is constant but not necessarily minimal. However, provided $\lambda$ is a mass scale, $I_2$ has the dimension of the mechanical power. Minimising $I_2$ during motion is therefore a way to minimise power expenditure. Figure 2 gives a schematic representation of $I_2$ in this case.

Note that the 2/3-PL and minimal jerk, which are closely linked within our framework, are not the only proposed models to describe voluntary motion. Let us consider Lagrangian (54) in which the position degree of freedom is actually an angular degree of freedom describing a circular motion in a plane, i.e., a circular motion restricted to a single constant axis of rotation. In this case, the movement tends to minimise $I_2\sim {\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{\dot{\alpha }}^{2}dt$, with $\alpha$ the angular acceleration. By virtue of Euler equation with a single rotation axis, $\dot{\alpha }\sim \dot{M}$ with M the external torque applied on the considered body and $I_2\sim {\int }_{t_{\mathrm{i}}}^{t_{\mathrm{f}}}{\dot{M}}^{2}dt$. This is a simplified version of the minimal torque-change model [24] suggesting that voluntary human motion minimises the rate of change of total external torque. This last comment suggests that, beyond the 2/3-PL, the integration of Lagrangian and Hamiltonian formulations of mechanics may fit within a general motor control model of voluntary motion encompassing several kinds of cost functions.

5.2. First Principles Shaping Voluntary Motion

Beyond the confines of the 2/3-PL, humans exhibit various stereotyped behaviours, such as bell-shaped velocity profiles [25] and adherence to Fitt’s law [26]. In the former case, the kinematics of horizontal reaching movements invariably shows a peak velocity in the middle of the trajectory, with a slight asymmetry depending on whether the movement is upward or downward [27]. These phenomena, which describe kinematic features of horizontal reaching movements and the compromise between speed and accuracy, are underpinned by the central nervous system’s capacity to maintain consistency in the face of environmental variability. Our study here suggests that theoretical constructs like action variables may mirror such neural processes, serving as constants/constraints in the complex generation of human motion. The central nervous system, equipped with its intricate network of neurons and sensory receptors, plays a pivotal role in orchestrating movements. The brain, acting as the command centre, receives sensory feedback, processes it, and generates precise motor commands. It is conceivable that the brain, through its complex computations and feedback loops, sets dynamical principles and initial conditions that align with the physiological capabilities of the human body, ultimately giving rise to observed movement patterns, including the 2/3-PL.

We expect that assuming a higher-derivative action principle as the dynamical principle shaping voluntary human movement fits naturally within an optimal feedback control model; see, e.g., [3,28]. Although the detailed matching between both frameworks is the subject of a separate paper, the correspondence can be summarised as follows. A higher-derivative action can be seen as the optimal feedback control law, coming into play at the level of movement planification/task selection: knowing the initial and final points to reach, an action principle allows us to plan the optimal movement and its trajectory in phase-space. This leads to an estimation of the initial conditions, with the initial acceleration as a tunable parameter as a consequence of the higher-derivative nature of the action. It is indeed crucial that the acceleration, i.e., muscular force, can be tuned by one individual and is not just a consequence of the external forces as in standard Newtonian mechanics. Then, after movement is initiated and sensory feedback is collected, the phase-space trajectory can be used as optimal state estimator: The actual trajectory can be compared to the planned one, and corrections can be made by the individual. The invariance of action variables can be seen as a further control policy in the case of rhythmic movements. Also, identifying the neural correlates for these mechanical concepts is an open challenge. While certain brain regions responsible for sensory processing are understood well enough, those involved in integrating multimodal sensory information to estimate physical forces like gravity remain less defined [29].

The framework considered here can be experimentally studied. A first type of study is to compute trajectories generated by a particular form of potential term, then to ask participants to perform those trajectories and check that they actually follow the 2/3-PL. A second type of study is to identify movements that lead to self-intersecting trajectories in a (position, speed)-plane, hence incompatible with standard dynamics. Then, the underlying higher-derivative Lagrangian may be found by solving the inverse problem. As an example, such movements have been found in participants performing infinity-shaped rhythmic motion of the forearm in parabolic flight [30]. A last direction is to make one individual evolve in a time-changing environment, for example, in variable gravity. If an action principle shapes voluntary movement, the time-changing environment a priori causes some parameters of the Lagrangian to be variable. The corresponding modifications of trajectories can be predicted via the equations of motion, and the question of whether an individual actually perform these changes or not can be experimentally studied. We have found in Ref. [30] the first hints that such a study is actually the case in variable gravity and consider to pursue this research program in forthcoming investigations.