We show here results of four-component calculations of nuclear magnetic resonance σ for atoms with 10 ≤ Z ≤ 86 and their ions, within the polarization propagator formalism at its random phase level of approach, and the first estimation of quantum electrodynamic (QED) effects and Breit interactions of those atomic systems by using two theoretical effective models. We also show QED corrections to σ(X) in simple diatomic HX and X2 (X = Br, I, At) molecules. We found that the Z dependence of QED corrections in bound-state many-electron systems is proportional to Z5, which is higher than its dependence in H-like systems. The analysis of relativistic ee (or paramagneticlike) and pp (or diamagneticlike) terms of σ exposes two different patterns: the pp contribution arises from virtual electron-positron pair creation/annihilation and the ee contribution is mainly given by 1s → ns and 2s → ns excitations. The QED effects on shieldings have a negative sign, and their magnitude is larger than 1% of the relativistic effects for high-Z atoms such as Hg and Rn, and up to 0.6% of its total four-component value for neutral Rn. Furthermore, percentual contributions of QED effects to the total shielding are larger for ionized than for neutral atoms. In a molecule, the contribution of QED effects to σ(X) is determined by its highest-Z atoms, being up to -0.6% of its total σ value for astatine compounds. It is found that QED effects grow faster than relativistic effects with Z.