The two-dimensional transverse-field Ising model on the square lattice is non-integrable. As a consequence the total energy is the only conserved quantitiy under unitary time evolution, thus it is possible that the system thermalizes, but non-integrability is not equivalent to thermalization. For this reason it is necessary to check whether the system really thermalizes after quantum quenches. In order to do so we compute the time evolution of the magnetization and correlation functions after quenches of the coupling between the spins (interaction quenches) in the paramagentic phase or of the external transverse field (field quenches) in the ferromagnetic phase and compare the time averages of the observables to ensemble averages of the system in thermal equilibrium according to the canonical Gibbs ensemble at a temperature defined by the quench parameters. The conditional equation for the so-called effective temperature of the system after the quenches is given by the condition that the thermal excess energy in the system at the effective temperature has to be equal to the excess energy due to the quench. As the two-dimensional transverse-field Ising model cannot be solved exactly, we apply real-time Variational Monte Carlo to compute its time evolution after the quenches. This method allows us to compute both the time evolution of the observables as well as the time-dependent distributions of the eigenvalues of the operators describing the observables for large system sizes up to 16x16 and long time scales. We find that after the interaction quenches the system thermalizes, while after the field quenches there are significant deviations to the thermal system which increase with the quench strength and which become especially clear considering the distributions of the eigenvalues.