We investigate the level of galaxy assembly bias in the Sloan Digital Sky Survey (SDSS) main galaxy redshift survey using ELUCID, a state-of-the-art constrained simulation that accurately reconstructed the initial density perturbations within the SDSS volume. On top of the ELUCID haloes, we develop an extended HOD model that includes the assembly bias of central and satellite galaxies, parameterized as $Q_\mathrm{cen}$ and $Q_\mathrm{sat}$, respectively, to predict a suite of one- and two-point observables. In particular, our fiducial constraint employs the probability distribution of the galaxy number counts $N_8^g$ and the projected correlation functions of quintiles of galaxies selected by $N_8^g$. We perform extensive tests of the efficacy of our method by fitting the same observables to mock data using both constrained and non-constrained simulations. We discover that in many cases the level of cosmic variance between the two simulations can produce biased constraints that lead to an erroneous detection of galaxy assembly bias if the non-constrained simulation is used. When applying our method to the SDSS data, the ELUCID reconstruction effectively removes an otherwise strong degeneracy between cosmic variance and galaxy assembly bias in SDSS, enabling us to derive an accurate and stringent constraint on the latter. Our fiducial ELUCID constraint, for galaxies above a stellar mass threshold $M_* = 10^{10.2} , h^{-2} , M_\odot$, is $Q_\mathrm{cen}{=}{-}0.09\pm{0.05}$ and $Q_\mathrm{sat}{=}0.09\pm{0.10}$, indicating no evidence for a significant galaxy assembly bias in the local Universe probed by SDSS. Finally, our method provides a promising path to the robust modelling of the galaxy-halo connection within future spectroscopic surveys like DESI and PFS.