Gaussian Process (GP) model interpolation is used extensively in geostatistics. We investigated the effectiveness of using GP model interpolation to generate maps of cortical activity as measured by Near Infrared Spectroscopy (NIRS). GP model interpolation also produces a variability map, which indicates the reliability of the interpolated data. For NIRS, cortical hemodynamic activity is spatially sampled. When generating cortical activity maps, the data must be interpolated. Popular NIRS imaging software HomER uses Photon Migration Imaging (PMI) and Diffuse Optical Imaging (DOI) techniques based on models of light behaviour to generate activity maps. Very few non-parametric methods of NIRS imaging exist and none of them indicate the reliability of the interpolated data. Our GP model interpolation algorithm and HomER produced activity maps based on data generated from typical functional NIRS responses. Image results in HomER were taken as the bench mark as the images produced are commonly considered to be representative of the true underlying hemodynamic spatial response. The output from the GP approach was then compared to these on a qualitative basis. The GP model interpolation appears to produce less structured image maps of hemodynamic activity compared to those produced by HomER, however a broadly similar spatial response is compelling evidence of the utility of GP models for such applications. The additional generation of a variability map which is produced by the GP method may have some utility for functional NIRS as such information is not explicitly available from standard approaches. GP model interpolation can produce spatial activity maps from coarsely sampled NIRS data sets without any knowledge of the system being modelled. While the images produced do not appear to have the same feature resolution as photonic model-based methods the technique is worthy of further investigation due to its relative simplicity and, most intriguingly, its generation of ancillary information in the form of the variability map. This additional data may have some utility in NIRS optode design or perhaps it may have application as additional input for response classification purposes. This GP technique may also be of use where model information is inadequate for DOI techniques.