The progression of breast cancer involves cancer-cell invasions of extracellular matrices. To investigate the progression, 3D cell cultures are widely used along with different types of matrices. Currently, the matrices are often characterized using parallel-plate rheometry for matrix viscoelasticity, or liquid-like viscous and stiffness-related elastic characteristics. The characterization reveals averaged information and sample-to-sample variation, yet, it neglects internal heterogeneity within matrices, experienced by cancer cells in 3D culture. Techniques using optical tweezers and magnetic microrheometry have measured heterogeneity in viscoelasticity in 3D culture. However, there is a lack of probabilistic heterogeneity quantification and cell-size-relevant, microscale-viscoelasticity measurements at breast-tumor tissue stiffness up to ≃10 kPa in Young's modulus. Here, we have advanced methods, for the purpose, which use a magnetic microrheometer that applies forces on magnetic spheres within matrices, and detects the spheres displacements. We present probabilistic heterogeneity quantification using microscale-viscoelasticity measurements in 3D culture matrices at breast-tumor-relevant stiffness levels. Bayesian multilevel modeling was employed to distinguish heterogeneity in viscoelasticity from the effects of experimental design and measurement errors. We report about the heterogeneity of breast-tumor-relevant agarose, GrowDex, GrowDex-collagen and fibrin matrices. The degree of heterogeneity differs for stiffness, and phase angle (i.e. ratio between viscous and elastic characteristics). Concerning stiffness, agarose and GrowDex show the lowest and highest heterogeneity, respectively. Concerning phase angle, fibrin and GrowDex-collagen present the lowest and the highest heterogeneity, respectively. While this heterogeneity information involves softer matrices, probed by ≃30 μm magnetic spheres, we employ larger ≃100 μm spheres to increase magnetic forces and acquire a sufficient displacement signal-to-noise ratio in stiffer matrices. Thus, we show pointwise microscale viscoelasticity measurements within agarose matrices up to Young's moduli of 10 kPa. These results establish methods that combine magnetic microrheometry and Bayesian multilevel modeling for enhanced heterogeneity analysis within 3D culture matrices.