This research examines the behavior of a class of lattice Boltzmann (LB) models designed to simulate immiscible multiphase flows. There is some debate in the scientific literature as to whether or not the "color gradient" models, also known as the Rothman-Keller (RK) models, are able to simulate flow with density variation. In this paper, we show that it is possible, by modifying the original equilibrium distribution functions, to capture the discontinuity present in the analytical momentum profile of the two-layered Couette flow with variable density ratios. Investigations carried out earlier were not able to simulate such a flow correctly. Now, with the proposed approach, the new scheme is compatible with the analytical solution, and it is also possible to simulate the two-layered Couette flow with simultaneous density ratios of 0(1000) and viscosity ratios of 0(100). To test the model in a more complex flow situation, i.e. with non-zero surface tension and a curved interface, an unsteady simulation of an oscillating bubble with variable density ratio is undertaken. The numerical frequency of the bubble is compared to that of the analytical frequency. It is demonstrated that the proposed modification greatly increases the accuracy of the model compared to the original model, i.e. the error can be up to one order of magnitude lower with the proposed enhanced equilibrium distribution functions. The authors believe that this improvement can be made to other RK models as well, which will allow the range of validity of these models to be extended. This is, in fact, what the authors found for the method proposed in this article. (C) 2013 Elsevier Ltd. All rights reserved.