This paper considers the problem of identifying bilinear systems from input-output measurements. A novel approach to this problem is proposed, one based upon the so-called Expectation Maximisation algorithm, wherein maximum likelihood estimates are generated iteratively without the need for a gradient-based search algorithm. This simple method is shown to perform well in simulation and it allows a multivariable bilinear system to be estimated directly in state-space form without the need for explicit parameterisation.