Solute transport and mixing in heterogeneous porous media are important to many processes of practical applications. Most of the previous studies focused on solute transport in flow of Newtonian fluids, whereas there are many processes in which the phenomenon takes place in flow of a non-Newtonian fluid. In this paper, we develop a computational approach to evaluate and upscale dispersion of a solute in flow of a shear-thinning (ST) fluid in a heterogeneous porous medium. Our results indicate that the dispersivity is a non-monotonic function of the Péclet number and the shear rate, and this behavior is accentuated by the heterogeneity of the pore space and spatial correlations between the local permeabilities. As a result, solute transport in ST fluids deviates significantly from the same phenomenon in Newtonian fluids. Moreover, the shear-dependence of the dispersivity strongly influences the fate of solute transport in porous media at large length scales, including larger effluent concentration at the breakthrough point, which also occurs much faster than Newtonian fluids. To provide further evidence for the numerical findings, we compare dispersion in flow of a power-law fluid in a single tube with the same in a bundle of such tubes. Our results emphasize the shortcomings of the current theories of dispersion to account for the role of fluid rheology in solute mixing and spreading.