Remote sensing images play a significant role in global land cover monitoring. However, due to the influence of cloud contamination, optical remote sensing images inevitably contain a large number of missing data, which severely limits their applicability. Existing cloud removal methods generally use only the effective information from a single band of temporally close known images, which is insufficient to predict accurately the changes between the known and cloudy images. In this paper, we proposed a spatial-spectral random forest (SSRF) method for thick cloud removal by gap filling, which uses spatially adjacent and multispectral information of known images simultaneously based on random forests. With its capability to fit nonlinear relations and adaptively assign variable contributions, SSRF can handle the potentially complex relationship between the known and cloudy images, thus, producing more accurate predictions. Based on the Landsat 8 OLI and Sentinel-2 MSI data in 13 regions, the effectiveness of SSRF was demonstrated through experiments on both simulated and real cloudy images. The results show that in urban areas with strong heterogeneity and agricultural areas with complex temporal changes, SSRF can yield satisfactory predictions both visually and quantitatively. Moreover, SSRF is more accurate than two popular benchmark methods. In addition, SSRF is less affected by cloud size and the time interval between the known and cloudy images, and can still produce reliable predictions when the used known images are also contaminated by clouds which reduce the amount of available neighborhood information. SSRF is also less affected by the omission error in cloud detection caused by thin clouds. SSRF is simple to implement and has great potential for widespread application.