We model the spatial and seasonal variability of three key water quality variables (water temperature and concentration of nitrates and phosphates) for entire river networks in a large area in northern Spain. Models were developed with the Random Forest technique, using 12 (water temperature and nitrate concentration) and 15 (phosphate concentration) predictor variables as descriptors of several environmental attributes (climate, topography, land-uses, hydrology and anthropogenic pressures). The effect of the different predictors on the response variables was assessed with partial dependence plots and partial correlation analysis. Results indicated that land-uses were important predictors in defining the spatial and seasonal patterns of these three variables. Water temperature was positively related with air temperature and the upstream drainage area, whereas increases in forest cover decreased water temperature. Nitrate concentration was mainly related to the area covered by agricultural land-uses, increasing in winter, probably because of catchment run-off processes. On the other hand, phosphate concentration was highly related to the area covered by urban land-uses in the upstream catchment and to the proximity of the closest upstream effluent. Phosphate concentration increased notably during the low flow period (summer), probably due to the reduction of the dilution capacity. These results provide a large-scale continuous picture of water quality, which could help identify the main sources of change in water quality and assist in the prioritization of river reaches for restoration projects.