Even though the effects of benthic fauna on aquatic biogeochemistry have been long recognized, few studies have addressed the combined effects of animal bioturbation and metabolism on ecosystem–level carbon and nutrient dynamics. Here we merge a model of benthic fauna (BMM) into a physical–biogeochemical ecosystem model (BALTSEM) to study the long-term and large-scale effects of benthic fauna on nutrient and carbon cycling in the Baltic Sea. We include both the direct effects of faunal growth and metabolism and the indirect effects of its bioturbating activities on biogeochemical fluxes of and transformations between organic and inorganic forms of carbon (C), nitrogen (N), phosphorus (P) and oxygen (O). Analyses of simulation results from the Baltic Proper and Gulf of Riga indicate that benthic fauna makes up a small portion of seafloor active organic stocks (on average 1 %–4 % in 2000–2020) but contributes considerably to benthic–pelagic fluxes of inorganic C (23 %–31 %), N (42 %–51 %) and P (25 %–34 %) through its metabolism. Results also suggest that the relative contribution of fauna to the mineralization of sediment organic matter increases with increasing nutrient loads. Further, through enhanced sediment oxygenation, bioturbation decreases benthic denitrification and increases P retention, the latter having far-reaching consequences throughout the ecosystem. Reduced benthic–pelagic P fluxes lead to a reduction in N fixation and primary production, lower organic matter sedimentation fluxes, and thereby generally lower benthic stocks and fluxes of C, N and P. This chain of effects through the ecosystem overrides the local effects of faunal respiration, excretion and bioturbation. Due to large uncertainties related to the parameterization of benthic processes, we consider this modelling study a first step towards disentangling the complex ecosystem-scale effects of benthic fauna on biogeochemical cycling.