In the precise point positioning/real-time kinematic (PPP-RTK) technique, high-precision ionospheric delay correction information is an important prerequisite for rapid PPP convergence. The commonly used ionospheric modeling approaches in the PPP-RTKs only take the trend term of the ionospheric total electron content (TEC) variations into account. As a result, the residual ionospheric delay still affects the positioning solutions. In this study, we propose a two-step regional ionospheric modeling approach that involves combining a polynomial fitting model (PFM) and a Kriging interpolation (KI) model. In the first step, a polynomial fitting method is used to model the trend term of the ionospheric TEC variations. In the second step, a KI method is used to compensate for the residual term of the ionospheric TEC variations. Datasets collected from continuously operating reference stations (CORSs) in Hunan Province, China, are used to validate the PFM/KI method by comparing with a single PFM method and a combined PFM and inverse distance weighting interpolation (IDWI) method. The experimental results show that the two-step PFM/KI modeled ionospheric delay achieves an average root mean square (RMS) error of 1.8 cm, which is improved by about 48% and 23% when compared with the PFM and PFM/IDWI methods, respectively. Regarding the positioning performance, the PPP-RTK with the PFM/KI method takes an average of 1.8 min or 4.0 min to converge to a positioning accuracy of 1.3 cm or 2.5 cm in the horizontal and vertical directions, respectively. The convergence times are decreased by about 18% and 14% in the horizontal direction and 9% and 5% in the vertical direction over the PFM and the PFM/IDWI methods, respectively.