This paper presents a geostatistical methodology which accounts for spatially varying population size in the processing of cancer mortality data. The approach proceeds in two steps: (1) spatial patterns are first described and modeled using population-weighted semivariogram estimators, (2) spatial components corresponding to nested structures identified on semivariograms are then estimated and mapped using a variant of factorial kriging. The main benefit over traditional spatial smoothers is that the pattern of spatial variability (i.e. direction-dependent variability, range of correlation, presence of nested scales of variability) is directly incorporated into the computation of weights assigned to surrounding observations. Moreover, besides filtering the noise in the data the procedure allows the decomposition of the structured component into several spatial components (i.e. local versus regional variability) on the basis of semivariogram models. A simulation study demonstrates that maps of spatial components are closer to the underlying risk maps in terms of prediction errors and provide a better visualization of regional patterns than the original maps of mortality rates or the maps smoothed using weighted linear averages. The proposed approach also attenuates the underestimation of the magnitude of the correlation between various cancer rates resulting from noise attached to the data. This methodology has great potential to explore scale-dependent correlation between risks of developing cancers and to detect clusters at various spatial scales, which should lead to a more accurate representation of geographic variation in cancer risk, and ultimately to a better understanding of causative relationships.