In this paper the combined integral method (CIM) is applied to non-classical two-phase Stefan problems with delayed onset of phase change. This can occur if the phase change is caused by a heat-flux or Robin boundary condition. The method requires choosing an approximating function, typically a polynomial, but it is not clear what should be used as the exponent in the highest order term. Previous studies have determined exponents either from exact solutions or from expansions valid over short time scales; neither approach is satisfactory and can be very inaccurate for larger times. We combined the heat balance and refined integral methods to determine this exponent as part of the solution process, allowing the exponent to be time-dependent. From comparing the approximate solutions with numerical and exact analytical solutions whenever possible, we show that the CIM greatly improves the accuracy on standard heat balance integral methods, without detracting from its simplicity. (C) 2014 Elsevier B.V. All rights reserved.