PostseismicGlobal Navigation Satellite System (GNSS) time seriesfollowed by megathrust earthquakes can be interpreted as a result of afterslip on the plate interface especially in its early phase. Afterslip is a stress release process accumulated by adjacent coseismic slip andcan be considered a recovery process for future events during earthquake cycles. Spatio-temporal evolution of afterslip often triggers subsequent earthquakes through stress perturbation. Therefore, it is important toquantitativelycapture the spatio-temporal evolution of afterslip and related postseismic crustal deformation and to predict their future evolution with a physics-based simulation. We developedanadjoint data assimilation method, which directly assimilates GNSS time series into a physics-based model to optimize the frictional parameters that control the slip behavior on the fault.The developed method was validated with synthetic data. Through the optimization of frictional parameters, the spatial distributions of afterslip can be roughly reproduced but not in detail if the observation noise is included. The optimization of frictional parameters provides not only the reproduction ofpostseismic displacements used for the assimilation but also the improvement in the prediction skill of the following time series. Then, we appliedthe developed method to the observed GNSS time series for the first 15 d following the 2003 Tokachi-oki earthquake. The frictional parameters in the afterslip regions were optimized to A-B ~ O(10 kPa), A ~ O(100 kPa), and L ~ O(10 mm). The large afterslip is inferred on the shallower side of the coseismic slip area. The optimized frictional parameters quantitatively predicted the postseismicGNSS time series for the following 15 d. These characteristics can be also detected if the simulation variables can be simultaneously optimized. The developed data assimilation method, which can be directly applied to GNSS time series following megathrust earthquakes, isan effective quantitative evaluation method for assessing risks of subsequent earthquakes and for monitoring the recovery process of megathrust earthquakes.