We investigate the dynamics of kinetically constrained models of glass formers by analysing the statistics of trajectories of the dynamics, or histories, using large deviation function methods. We show that, in general, these models exhibit a first-order dynamical transition between active and inactive dynamical phases. We argue that the dynamical heterogeneities displayed by these systems are a manifestation of dynamical first-order phase coexistence. In particular, we calculate dynamical large deviation functions, both analytically and numerically, for the Fredrickson-Andersen model, the East model, and constrained lattice gas models. We also show how large deviation functions can be obtained from a Landau-like theory for dynamical fluctuations. We discuss possibilities for similar dynamical phase-coexistence behaviour in other systems with heterogeneous dynamics.