Modeling the boiling phase change is particularly challenging for Lagrangian particle methods due to a high density ratio and dramatic volume expansion. In this study, the incompressible moving particle semi-implicit (MPS) method and the weakly compressible smoothed particle hydrodynamics (SPH) method are coupled to develop an incompressible–compressible particle method for modeling a multiphase flow with boiling. The coupling strategies developed by Lind et al. (JCP, 2016) are adopted. A high speed of sound must be employed in SPH for incompressible bubbles in a heavy liquid under gravity, which can cause severe spurious pressure fluctuations. In this situation, a technique for detecting bubble connectivity is developed to average the pressure inside each bubble for stable coupling. This pressure averaging strategy is physically consistent due to the high density ratio. Because the pressure coupling strategy tends to break the force balance between the surface tension and the pressure jump across the interface, the surface tension model must be based only on the liquid particles near the interface. The boiling mass transfer is modeled in a straightforward manner by injecting gas particles towards the gas phase from the liquid interface particles. The gas volume expansion is naturally considered by the coupling method, as the gas particle distribution can be automatically regulated in the SPH calculation. Rising bubble simulations with different surface tension coefficients and density ratios verify the proposed method for a bubble flow. The Stefan and sucking problems as well as a horizontal film boiling flow are simulated to verify the boiling model and to demonstrate its straightforward capability to handle gas volume expansion.