Mathematical models promise new insights into the mechanisms underlying the emergence of collective behaviour in fish. Here, we establish a mathematical model to examine collective behaviour of zebrafish, a popular animal species in preclinical research. The model accounts for social and hydrodynamic interactions between individuals, along with the burst-and-coast swimming style of zebrafish. Each fish is described as a system of coupled stochastic differential equations, which govern the time evolution of their speed and turn rate. Model parameters are calibrated using experimental observations of zebrafish pairs swimming in a shallow water tank. The model successfully captures the main features of the collective response of the animals, by predicting their preference to swim in-line, with one fish leading and the other trailing. During in-line swimming, the animals share the same orientation and keep a distance from each other, owing to hydrodynamic repulsion. Hydrodynamic interaction is also responsible for an increase in the speed of the pair swimming in-line. By linearizing the equations of motion, we demonstrate local stability of in-line swimming to small perturbations for a wide range of model parameters. Mathematically backed results presented herein support the application of dynamical systems theory to unveil the inner workings of fish collective behaviour.