  We study the cosmology of Galileon modified gravity models in the linear
perturbation regime. We derive the fully covariant and gauge invariant
perturbed field equations using two different methods, which give consistent
results, and solve them using a modified version of the {\tt CAMB} code. We
find that, in addition to modifying the background expansion history and
therefore shifting the positions of the acoustic peaks in the cosmic microwave
background (CMB) power spectrum, the Galileon field can cluster strongly from
early times, and causes the Weyl gravitational potential to grow, rather than
decay, at late times. This leaves clear signatures in the low-$l$ CMB power
spectrum through the modified integrated Sachs-Wolfe effect, strongly enhances
the linear growth of matter density perturbations and makes distinctive
predictions for other cosmological signals such as weak lensing and the power
spectrum of density fluctuations. The quasi-static approximation is shown to
work quite well from small to the near-horizon scales. We demonstrate that
Galileon models display a rich phenomenology due to the large parameter space
and the sensitive dependence of the model predictions on the Galileon
parameters. Our results show that some Galileon models are already ruled out by
present data and that future higher significance galaxy clustering, ISW and
lensing measurements will place strong constraints on Galileon gravity.
