Detecting anomalous changes in remote sensing images is a challenging problem, where many approaches and techniques have been presented so far. We rely on the standard field of multivariate statistics of diagnostic measures, which are concerned about the characterization of distributions, detection of anomalies, extreme events, and changes. One useful tool to detect multivariate anomalies is the celebrated Cook's distance. Instead of assuming a linear relationship, we present a novel kernelized version of the Cook's distance to address anomalous change detection in remote sensing images. Due to the large computational burden involved in the direct kernelization, and the lack of out-of-sample formulas, we introduce and compare both random Fourier features and Nystrom implementations to approximate the solution. We study the kernel Cook's distance for anomalous change detection in a chronochrome scheme, where the anomalousness indicator comes from evaluating the statistical leverage of the residuals of regressors between time acquisitions. We illustrate the performance of all algorithms in a representative number of multispectral and very high resolution satellite images involving changes due to droughts, urbanization, wildfires, and floods. Very good results and computational efficiency confirm the validity of the approach.