A search is performed for dark matter (DM) produced in association with a single top quark or a pair of top quarks using the data collected with the CMS detector at the LHC from proton-proton collisions at a center-of-mass energy of 13 TeV, corresponding to 138 fb−1 of integrated luminosity. An excess of events with a large imbalance of transverse momentum is searched for across 0, 1 and 2 lepton final states. Novel multivariate techniques are used to take advantage of the differences in kinematic properties between the two DM production mechanisms. No significant deviations with respect to the standard model predictions are observed. The results are interpreted considering a simplified model in which the mediator is either a scalar or pseudoscalar particle and couples to top quarks and to DM fermions. Axion-like particles that are coupled to top quarks and DM fermions are also considered. Expected exclusion limits of 410 and 380 GeV for scalar and pseudoscalar mediator masses, respectively, are set at the 95% confidence level. A DM particle mass of 1 GeV is assumed, with mediator couplings to fermions and DM particles set to unity. A small signal-like excess is observed in data, with the largest local significance observed to be 1.9 standard deviations for the 150 GeV pseudoscalar mediator hypothesis. Because of this excess, mediator masses are only excluded below 310 (320) GeV for the scalar (pseudoscalar) mediator. The results are also translated into model-independent 95% confidence level upper limits on the visible cross section of DM production in association with top quarks, ranging from 1 pb to 0.02 pb.