In this study, we used in situ measurements for the first time to analyze the applicability and effectiveness of evaluating groundwater storage (GWS) changes across China using Gravity Recovery and Climate Experiment (GRACE) satellite products and hydrological data derived from the WaterGap Global Hydrological Model (WGHM), Global Land Data Assimilation System (GLDAS) and eartH2Observe (E2O). The results show that the GWS derived from GRACE JPL Mascons products combined with GLDAS Noah V2.1 data most accurately reflect the overall distribution of GWS changes in China and the correlation coefficient between the in situ measurements reaches 0.538. The empirical orthogonal function decomposition for GWS indicates clear interannual variation and seasonal variation in China. The trends of China’s GWS changes showed a clear regional characteristic from 2003 to 2016. The GWS in the northeast, central-south, and western junction of Xinjiang-Qinghai-Tibet had increased significantly, and the North China Plain (NCP) had a severe decline. The correlation coefficient between the annual trends of precipitation and GWS was 0.57, and it reached 0.73 when four provinces (Beijing, Tianjin, Shanxi, Hebei) that are wholly or partially located in the NCP were excluded. The seasonal variability of GWS in China was obvious and the volatilities in Jiangxi, Hunan and Fujian provinces were the highest, reaching 6.39 cm, 6.33 cm and 5.20 cm, respectively. The empirical orthogonal function decomposition for GWS and precipitation over China indicated seasonal consistency with a correlation coefficient of 0.76. The awareness of areas with significant depletion and large seasonal fluctuation of GWS help adaptations to manage local GWS situation.