For pulsar projects it is often necessary to predict the pulse phase in advance, for example when preparing for new observations. Interpolation of the pulse phase between existing measurements is also often required, for example when folding X-ray or gamma-ray observations according to the radio pulse phase. Until now these procedures have been done using various ad hoc methods. The purpose of this paper is to show how to interpolate or predict the pulse phase optimally using statistical models of the various noise processes and the phase measurement uncertainty. © 2012 The Authors Monthly Notices of the Royal Astronomical Society © 2012 RAS.