In addition to studies of sea level change and mantle rheology, reliable Glacial Isostatic Adjustment (GIA) models are necessary as a background model to correct the widely used Gravity Recovery and Climate Experiment (GRACE) monthly gravity solutions to determine subsecular, nonviscous variations. Based on spherical harmonic analyses, we developed a method using degree-dependent weighting to assimilate the Global Positioning System (GPS) derived crustal uplift rates into GIA model predictions, in which the good global pattern of GIA model predictions and better local resolution of GPS solutions are both retained. Some systematic errors in global GPS uplift rates were also corrected during the spherical harmonic analyses. Further, we used the refined GIA uplift rates to infer the GIA-induced rates of Stokes coefficients (complete to degree/order 120) relying on the accurate relationship between GIA vertical surface deformation and gravitational potential changes. The results show notable improvements relative to GIA model outputs, and may serve as a GIA-correction model for GRACE time-variable gravity data. more interested in the rate of crust uplift . u GIA ≡ ∂u GIA (θ, λ, t)/∂t and rate of geopotential change∂t rather than u GIA and δV GIA themselves. Reliable . u GIA can be used to infer mantle viscosity, and δ .V GIA are necessary as a background model to correct the widely used Gravity Recovery and Climate Experiment (GRACE) monthly gravity solutions for these secular variations [9][10][11][12].Currently, the GIA-induced uplift rates can be obtained either from a certain GIA model [13][14][15] or from Global Positioning System (GPS) measurements [16]. However, the agreements between model-based and GPS-based results are usually limited due to a number of factors (some are discussed in the following sections). It is believed that the GIA models are good at describing large-scale and especially global GIA patterns while GPS data are better at local details although the GPS rates can also be affected by unrelated local effects [17,18].Arguing that the GIA-induced crust uplifts are dominated by the Earth's viscosity mode M0 [19], Wahr et al. [20,21] proposed an approximate theory to infer changes in Stokes coefficients from the measured GIA uplift rates. Later, Purcell et al. [22] refined the theory relying on viscous load Love numbers (also see Section 2.3 for some details). While most studies [23][24][25] of the GIA-induced geopotential change are on the basis of the approximate theory of Wahr et al., Purcell et al. [22] and Jia et al. [26] pointed out that this method will lead to a~15% relative error, which would weaken GRACE's monitoring of global changes. However, both of the two methods ignore possible systematic errors in measured uplift rates as discussed in this study.On the other hand, Argus et al. [13] and Peltier et al. [14,15] had assimilated GPS uplift rates into their GIA numerical simulations, leading to a notable improvement to their GIA models. However, only a few tens of GPS sites were used...