Accurate wavefront integration based on gradient fields is crucial for various indirect measurement techniques, such as Shack-Hartmann sensing, shearography, and the fringe reflection technique. In this paper, a higher-order iterative compensation algorithm is proposed to enhance the reconstruction accuracy for the finite-difference-based least-squares integration (FLI) method. In this method, higher-order gradient fields are reconstructed and the calculated residual gradient fields compensate the truncation error with the traditional FLI by iterations. A comparison of different FLI methods, including traditional FLI, iterative FLI, higher-order FLI, and the proposed FLI method, is conducted. The result shows that the reconstructed wavefront with the proposed method is more accurate than those with other FLI methods. In addition, the impact of the gradient measurement noise is also discussed.