## 1 Introduction  
Hyperspectral images (HSIs) with high spatial resolution play an increasingly important role in feature extraction and terrain classification of the underlying materials [1; 2; 3]. However, hyperspectral sensors have some limitations regarding spatial resolution, such as hardware and technology. Furthermore, the spatial and spectral resolutions are closely linked; i.e., better spectral resolution automatically leads to coarser spatial resolution and vice versa. Therefore, it is very difficult to obtain high-spatial-resolution HSIs directly. As an effective and economical approach to image enhancement or super-resolution, fusion of hyperspectral and multispectral images (MSIs) has attracted extensive attention in remote sensing, and substantial research resultshave emerged [4; 5]. Nevertheless, there are still many problems to be solved [1; 4; 6], such as spatial resolution, spectral distortion and anti-noise performance. Pansharpening, originally designed to solve the fusion of MSIs and panchromatic images, can be adapted to handle the fusion of HSIs and MSIs using the classical information processing methods [5], such as the component substitution method [7; 8] and multi-resolution analysis [9; 10]. For example, the component substitution method [8] constructed the spectral band mapping between the MSI and clustered HSI, and then replaced the spatial component of the HSI with the corresponding part of the MSI for the reconstructed image. The multi-resolution-analysis-based pansharpening methods [7] linearly combined multispectral band images to synthesize a high-spatial-resolution image for each HSI band. However, these methods do not fully exploit the intrinsic relationship between the HSI and MSI, which often results in serious spectral distortion and performance degradation. In addition, with the development of artificial intelligence, there are many fusion algorithms based on deep neural networks [11; 12; 13]. For example, [PERSON] et al. [11] proposed a fusion method based on a deep convolutional neural network (CNN) with two branches for HSI and MSI separately, by exploiting spectral correlation. Multi-branch back propagation neural networks were employed by [PERSON] et al. [13] to propose a cluster-based fusion method to address the non-linear spectral mapping for each cluster. Generally, supervised learning needs a training database based on images or pixels. Data fusion aims at integrating the observed low-spatial-resolution HSI data and high-spatial-resolution MSI data to obtain the reconstructed super-resolution image with high spatial and spectral resolutions. Conversely, the observed MSIs and MSIs can be regarded as spectral and spatial degradation of the high-spatial-resolution HSI, which is called the observation model. Data fusion is often an inverse problem, made so by fusing the observed HSI and MSI data to enhance the spatial resolution of reconstructed images. Furthermore, the inverse problem of data fusion is ill-posed. To address the difficulties, Bayesian-based methods dedicated to the fusion of HSIs and MSIs have been designed, which adopted the priors (e.g., Guassian prior) to design the algorithms in the principal component subspace or wavelet domain [14]. Another popular approach is hyperspectral super-resolution (HySure) [15], in which the fusion problem was reformulated by introducing the vector-total-variation regularization into convex data fitting terms, and was solved by the split augmented Lagrangian shrinkage algorithm. In addition, the research by [PERSON] et al. [16] developed a Bayesian-based sparse representation (BSR) method, in which a proper posterior distribution was refined by the decomposition of the image on a set of dictionaries. However, this kind of fusion method typically involves a high computational load, so some measures were taken to reduce the large computational complexity of Bayesian-based fusion approaches; for example, a Sylvester equation-based fast fusion [17].