We propose a variational form of the BDF2 method as an alternative to the commonly used minimizing movement scheme for the time-discrete approximation of gradient flows in abstract metric spaces. Assuming uniform semi-convexity --- but no smoothness --- of the augmented energy functional, we prove well-posedness of the method and convergence of the discrete approximations to a curve of steepest descent. In a smooth Hilbertian setting, classical theory would predict a convergence order of two in time, we prove convergence order of one-half in the general metric setting and under our weak hypotheses. Further, we illustrate these results with numerical experiments for gradient flows on a compact Riemannian manifold, in a Hilbert space, and in the $L^2$-Wasserstein metric.