Abstract: In this research, ground LiDAR point cloud was used as a data source; hilly terrain under the northern subtropical secondary forest was set as a research object. Six commonly used zonal interpolation methods were applied to create DEMs, global error and the relationship between the errors and levels of the terrain factor, ground interpolation point density and surface vegetation cover were analyzed. Then the uncertainty map of interpolation error was created by the random forest method. The research results demonstrated that: The predicted values on the interpolation surface were generally underestimated and the best output resolution was 2 m; the natural neighbor interpolation method showed the highest accuracy and the best visual quality, but the tension spline method had the poorest accuracy; the global error increased corresponding to the increment of slope but decreased by the increment of ground interpolation point density; the terrain interpolation errors were relatively high under the young and middle age natural secondary forest but received the highest error under mature forest. The global errors under shrubs were not high but the variations of errors were high. Meanwhile, there existed a relative high correlation between LiDAR-derived forest parameters and the errors of terrain interpolation, while no significant relationship was found between NDVI and the global errors; this again proved that LiDAR-derived forest variables could better reflect the terrain interpolation accuracy in the NDVI saturation regions.