We apply a recently developed time-dependent one-zone leptohadronic model to study the emission of the blazar Mrk 421. Both processes involving proton-photon interactions, i.e. photopair (Bethe-Heitler) and photopion, have been modeled in great detail using the results of Monte Carlo simulations, like the SOPHIA event generator, in a self-consistent scheme that couples energy losses and secondary injection. We find that TeV gamma-rays can be attributed to synchrotron radiation either from relativistic protons or, alternatively, from secondary leptons produced via photohadronic processes. We also study the variability patterns that each scenario predicts and we find that while the former is more energetically favored, it is the latter that produces, in a more natural way, the usual quadratic behavior between X-rays and TeV gamma-rays. We also use the obtained SEDs to calculate in detail the expected neutron and neutrino fluxes that each model predicts.