Abstract: Soil organic carbon (SOC) pool is one of the important pathways for “carbon neutral”. However, one of the difficulties is how to accurately estimate the amount and its change of soil carbon stock during a period of time within a pre-defined region. This depends on a practical scheme of spatial data analysis, as well as high quality datasets from sampling and chemical analysis. The projects of multi-purpose regional geochemical survey (MPRGS) and land quality geochemical survey (LQGS) carried out in recent decades in China provided densely sampled soil data. However, the traditional statistical methods are not appropriate for those data, because they are not randomly sampled and the autocorrelation exists in them. Based on the geostatistical theory, this study established a methodology for calculating soil carbon storage and its error, detecting the regions of soil source and/or sink, and estimating the amount of soil organic carbon change during a tertian period with an example based on the MPRGS and LQGS datasets for arable land in Shandong province. The results indicated that the block-kriging method is appropriate for regional soil organic carbon stock estimation in that it does not require the random sampling data, and it fully considered the autocorrelation of SOC density (SOCD); The regions of significant SOCD change, including increase or decrease, can be detected with the punctual kriging method and the equation provided by the authors. As for the example of Shandong province, the 0-20cm SOC storage in arable land is 172.85±0.38Tg; The arable land in the typical regions, including irrigated, rainfed, and paddy, totally acted as a carbon “sink” during 2003-2020; The average annual SOCD increase is 0.0878±0.0010、0.0442±0.0016 and 0.0484±0.0054 kg/(m² a) respectively.