In this work, a fourth order Cartesian grid finite difference method with multigrid acceleration is introduced for solving parabolic interface problems. In comparison to other existing methods in this field, the proposed Augmented and Matched Interface and Boundary (AMIB) method includes novel developments to boost the order of convergence from two to four in both space and time. In temporal discretization, the acceleration in convergence is realized by applying the Richardson extrapolation to the Crank–Nicolson scheme. The bottleneck of the order enhancement in space in the previous parabolic study was due to the multigrid solver. To overcome this difficulty, a recently developed fourth order multigrid method for elliptic interface problems has been reformulated for the present parabolic problem. Moreover, this multigrid method has been further generalized so that it can handle any boundary condition, including Dirichlet, Neumann, Robin, and their mixed combinations. Numerical experiments indicate that the proposed AMIB multigrid method achieves a fourth order of accuracy in accommodating complex interfaces and variable coefficients, while maintaining unconditional stability. Moreover, the multigrid solver significantly improves the computational efficiency.