Most measurements of mass in Astronomy that use kinematics of stars or gas rely on assumptions of equilibrium that are often hard to verify. Instead, we develop a novel idea that uses the clustering in action space, as a probe of underlying gravitational potential: the correct potential should maximize small-scale clustering in the action space. We provide a first-principle derivation of likelihood using the two-point correlation function in action space, and test it against simulations of stellar streams. We then apply this method to the 2nd data release of Gaia, and use it to measure the radial force fraction f h and logarithmic slope α of dark matter halo profile. We investigate stars within 9-11 kpc and 11.5-15 kpc from Galactic centre, and find ( f h , α) = (0.391 ± 0.009, 1.835 ± 0.092) and (0.351 ± 0.012, 1.687 ± 0.079), respectively. We also confirm that the set of parameters that maximize the likelihood function do correspond to the most clustering in the action space. The best-fit circular velocity curve for Milky Way potential is consistent with past measurements (although it is ∼ 5-10% lower than previous methods that use masers or globular clusters).Our work provides a clear demonstration of the full statistical power that lies in the full phase space information, relieving the need for ad hoc assumptions such as virial equilibrium, circular motion, or steam-finding algorithms.