The availability of reliable and long time series of runoff data is fundamental for most of the hydrological analyses and for the assessment and the management of water resources even in condition of global climatic change. However, hydrologic data sets are often characterized by a short duration and also suffer from missing data values, mainly due to malfunctioning of gauging stations for a specific period. In order to overcome this problem and obtain long and continuous runoff time series, different models and methods have been previously developed and proposed. While some models, used to extent the streamflow record, are conceptual, empirical, regressive models based on the rainfall input, other models are based on the derivation of runoff maps at different time scale; these maps allow the runoff estimation in gauged basins characterized by the absence of data or in ungauged basins. Aim of this paper is the derivation of a map relative to the mean annual runoff at regional scale using a stochastic approach derived from the kriging interpolator. This approach can be assimilated to a kriging system, which considers explicitly the areal nature of runoff variable by imposing the constraint of the water balance; it allows to derive gridded annual runoff maps with finer and finer resolution. The methodology has been applied to 23 main hydrographic basins of Sicily, Italy using the mean annual runoff dataset provided by Osservatorio delle Acque. All these basins have been previously grouped in three homogeneous zones, as suggested by previous studies. A cross-validation procedure has been performed in order to validate the procedure for each homogeneous zone.