  Large longitudinal studies provide lots of valuable information, especially
in medical applications. A problem which must be taken care of in order to
utilize their full potential is that of correlation between intra-subject
measurements taken at different times. For data in Euclidean space this can be
done with hierarchical models, that is, models that consider intra-subject and
between-subject variability in two different stages. Nevertheless, data from
medical studies often takes values in nonlinear manifolds. Here, as a first
step, geodesic hierarchical models have been developed that generalize the
linear ansatz by assuming that time-induced intra-subject variations occur
along a generalized straight line in the manifold. However, this is often not
the case (e.g., periodic motion or processes with saturation). We propose a
hierarchical model for manifold-valued data that extends this to include trends
along higher-order curves, namely B\'ezier splines in the manifold. To this
end, we present a principled way of comparing shape trends in terms of a
functional-based Riemannian metric. Remarkably, this metric allows efficient,
yet simple computations by virtue of a variational time discretization
requiring only the solution of regression problems. We validate our model on
longitudinal data from the osteoarthritis initiative, including classification
of disease progression.
