An important issue in the numerical simulation of geothermal reservoirs is the problem of scales. Data are collected at a scale usually smaller than the one used to discretise the sedimentological units in the numerical model. For instance, thermal conductivities sampled from field scale cores have measurement support in the order of centimeters to meters, whereas numerical models for heat flow require conductivities representative of scales ranging between tens to hundreds of meters. We present a study aimed at demonstrating the upscaling of thermal conductivities. Based on the spatial characteristics of a large sample data set of thermal conductivities of permo-carboniferous sedimentary rocks, 10 different realizations of the system are randomly generated at a fine scale of resolution and are then upscaled to four different resolutions using diverse averaging procedures (based on arithmetic, geometric, or harmonic averaging) as well as renormalization. Results show that upscaling based on harmonic averaging of local values is superior in reproducing the original values while renormalization gives the poorest results. Generally it is demonstrated that the specific kind of upscaling has only a small impact on the resulting temperature distribution. Due to the diffusive character of heat conduction all results tend towards the arithmetic mean value associated with the data.