Multispectral (MS) remote sensed imagery, reflecting the radiance from different land covers by more spectral bands, has the performance of accurately mapping the land surface composition. However, these multispectral sensors usually have lower spatial resolution, which limits their applications in mapping the complex land surface morphological structure. High spatial resolution remote sensed imagery, obtained from commercial satellite sensors, has the potential to give more accurate descriptions of urban surface, and has been used extensively in the fields of urban planning, urban building extraction and decision supporting   . Therefore, there is a desire to integrate the high spatial and high spectral information from these two kinds of imageries to give the most complete and accurate description of the study scene  .
Image fusion or pan-sharpening method is a technique producing images with high spatial and spectral resolution simultaneously, by injecting the spatial detail information in higher resolution panchromatic (PAN) image into the MS channels  .
Pan-sharpening means to use a panchromatic image to sharpen the multispectral images. There are several steps in a pan-sharpen algorithm. Firstly, registration between the PAN and MS images is made to get the spatial aligned images, which is a pivotal process to attain effective fusion results. Secondly, spatial information is extracted from the high resolution using a certain algorithm such as wavelet transform, intensity-hue-saturation (IHS) transform and principal component transform. Thirdly, the extracted spatial information is injected into the MS images to sharpen the spatial resolution meanwhile preserve the spectral information contained in the MS images. Finally, assessment will be made to evaluate the effectiveness of the pan-sharpen results. Another key point in this process lies in the mechanism of extraction and injecting of spatial information, which has become the hotspot issue in the applications of remotely sensed imageries.
Lots of pan-sharpening methods have been proposed in the past twenty years. These algorithms can be categorized into four categories: projection substitution methods, numerical methods, multi-resolution analysis based methods and hybrid methods.
IHS transformation and PCA transformation are two representative methods in the projection substitution fusion methods. In IHS transform, MS images are converted from Red-Green-Blue (RGB) color space into the intensity-hue-saturation color space, and then the intensity image, which mainly contains the low resolution spatial detail information, is substituted by the histogram matched PAN image. Fusion results are attained by inverse transforming from IHS to RGB color space     . This algorithm provides an effective and fast implementation for sharpening the MS images. However, it is reported that there is significant spectral distortion in the results which may be induced by adding inappropriate spatial information. A fast IHS algorithm with spectral adjustment for IKONOS imagery fusion is proposed by Tu et al.  to avoid the spectral distortion as far as possible.
PCA-based fusion was commonly used due to the uncorrelated property among the principal components after the PCA transform. The first principal component, which was considered as containing enough spatial information due to the largest variance compared with the remains principal components, was replaced by the histogram matched PAN image    . The new first principal component and all the other principal components preserving the spectral information were converted back to get the fused MS images with higher spatial resolution compared with the original MS images. However, “a higher variance of the first PC does not necessarily mean it has higher correlation with the PAN image”  . Therefore, several modified PCA-based fusion algorithms have been proposed recently to improve the effectiveness of algorithm    .
In the numerical fusion algorithms, PAN image is assumed as the linear combination of the original high resolution MS bands, such that the combination coefficients will be estimated using the degraded low resolution MS bands  . Brovey method  , color normalized, and P+XS  are the algorithms belong to this set. Disadvantages of such algorithms lie in the assumption of linear combination, which is inappropriate in reality and will lead to incorrect fusion results. Recently, Garzelli et al.  suggested an optimal algorithm, which based on the minimum mean-square-error (MMSE) sense, to sharpen the MS images. In this algorithm, the fused high resolution MS images are assumed as the weighted combination of low resolution MS images and PAN image, in which the weight coefficients can be estimated using the Least-Square (LS) algorithm. Another model in this literature is called Band-Dependent Spatial-Detail model, which used the assumption that spatial information can be induced from the difference between PAN image and the sum of LRMS images.
Lots of attentions have been paid on the multi-resolution analysis based methods. Idea behind such methods is that the missing spatial information in MS images can be inferred from the high frequencies, which is the foundation of ARSIS concept    . ARSIS comes from the French acronym for “Amélioration de la Résolution Spatiale par Injection de Structures”(Improving Spatial Resolution by Structure Injection)  . Multi-resolution analysis methods such as Wavelet analysis       , Pyramid decomposition, Contourlet analysis    -  and Shearlet analysis  are used to induce a scale-by-scale description of the information content of the PAN and MS images  . Among these algorithms, the key points are how to extract the spatial information as far as possible, and how to define a fusion rule to integrate the spatial information and the spectral information. Although different kinds of rules have been tested, a thorough investigation of this kind of algorithm is necessary to assess their performances.
Due to the limitation among different kinds of fusion algorithms, hybrid algorithms such as IHS and Wavelet, PCA and Wavelet, IHS and Contourlet, are used to give better fusion results. Intensity image or the first principal component will be extracted using the corresponding transform, and then wavelet decomposition will be used on the intensity image and PAN image simultaneously. The wavelet coefficient corresponding to the approximant part of the intensity image will be replaced by PAN image’s approximant wavelet coefficients. The fused MS image will be induced by inverse wavelet transform. Usually, better fusion effectives will be obtained using the hybrid algorithms.
As pointed by Tu   , the key point in IHS-based fusion algorithm lies in the extraction of spatial information, which can be deduced from the difference between PAN and Intensity images. Therefore, the new intensity image can be seen as the linear combination of PAN and the original Intensity image. The idea behind this algorithm motivates us to use the so called minimum mean-square error (MMSE) method to induce the new intensity image after an optimization calculation. Therefore, the proposed novel hybrid fused method is based on the IHS transform and the MMSE optimal algorithm.
Outline of this paper is as follows. A brief introduction is given in the first section. Then, the proposed hybrid fusion algorithm is introduced in section 2. Numerical experiments and results are shown in section 3. Section 4 gives the discussion of the experimental results and conclusions are made in Section 5.
2. The Hybrid Pansharpen Method
Based on the fact that the fused high resolution MS images contain the spatial information coming from low resolution MS images and the panchromatic image, the proposed hybrid pansharpen method utilized the optimal component coefficients of the MS images and the panchromatic image to get the optimum fusion result. The flowchart of the hybrid pan-sharpen method is shown in Figure 1. There are two key steps in the hybrid pansharpen method: IHS transformation and optimization calculation. IHS transformation is used to get the intensity image, which contains the spatial information of the MS images. Optimization calculation is used to get the final intensity image by calculating the optimal component coefficients.
2.1. IHS-Based Fusion Method
IHS transform is extensively used to convert the MS images from RGB color space into the IHS color space. The Intensity image contains most of the spatial information of the scene, while hue image and saturation image reflect the spectral information of the same land cover. Compared with the PAN image, Intensity image has lower spatial resolution, which makes the MS images shortage of spatial information. Therefore, usually, Intensity image is replaced by the histogram matched PAN image to increase the spatial structure of the MS images. Standard IHS-based fusion algorithm is introduced briefly as the following four steps.
Firstly, band combination of MS images is used to form the RGB components, and then, the low spatial resolution RGB images are resized by upsampling to match the size of the high spatial resolution PAN image  .
Figure 1. The logic flow of Hybrid IHS Pan-sharpen method.
Secondly, IHS transform is made to convert the images from RGB color space into IHS color space using Equation (1).
where and are the variables in the computation. Hue and saturation components in the IHS space are given as
Thirdly, the Intensity image I is replaced by the histogram matched PAN image. Finally, inverse IHS transform is used to get the fused MS images using Equation (3).
where and are the fused MS images.
Tu  introduced a computationally efficient method by rewrite the previous two equations, and the new formulation is given as
It was found that the spectral distortion mainly due to the change of saturation value, i.e., “the saturation value is expanded and stretched , when the PAN value is less than its corresponding I value; the saturation value is compressed when PAN value is larger than the I value”  . To avoid the change of saturation value among different land surface materials, we suggest using the optimization version of Intensity image as the replacement of the original Intensity image.
2.2. The Logic Flow of the Hybrid IHS Pan-Sharpen Method
To give a concise description of the hybrid pan-sharpen algorithm, some symbols are used to refer to the images and variables. Let , which has size of , denote the ith band of N MS images. Matrix P is the PAN image which has the size of , where r is the ratio of spatial resolution between the PAN image and the MS image. For example, r equals four for QuickBird sensor. and P also is used to denote the lexicographically ordered vector which have the size of and , respectively.
According to the ratio r of spatial resolution between PAN image and MS image, low resolution MS images are upsampled to get new MS images which have same size with that of PAN image. Then, IHS transform is used to convert the new MS images from RGB color space into IHS color space to get three component images: intensity image I, hue image H, saturation image S.
New PAN image P1 can be deduced using histogram match by Equation (6).
where are mean value of PAN image and Intensity image respectively, and are standard deviation of PAN image and intensity image, respectively.
The new intensity image , which will be estimated using optimization algorithm, can be written as
where and are coefficients to be defined.
This formulation is similar to the single spatial-detail (SSD) model given by Garzelli et al. for image fusion  . But, there are two differences between the two models after an in-depth investigation. Firstly, SSD model is used to describe the relationships among estimated HRMS image, LRMS image and PAN image, i.e., the ith band of LRMS can be expressed as
where is the parameter to be estimated. While, our model is used to depict the relationship among new intensity image, original intensity image and histogram matched PAN image. Secondly, in SSD model, parts of spatial information in PAN image are added into the low resolution MS images to get the high resolution MS images, which is necessary to enhance the spatial structure in the high resolution MS images. Whereas, in our model, the new intensity image is estimated as the linear combination of intensity image I and the histogram matched PAN image, which is better than the situation in which I is added or replaced totally by PAN image, due to the fact that the spatial information in intensity image will be lost.
To estimate the parameters and , we employ the least-square criteria, i.e., to minimize the following object function:
where denote the square of 2-norm of a vector.
Least-square solution can be deduced by calculating the partial derivative, and the solution can be expressed as
where denote the transpose matrix of matrix M, and denotes the inverse matrix of matrix M.
To give a better fusion result, the optimized calculation can be implemented in the sliding window which has the size of , i.e., the parameters should be estimated in each non-overlapped sliding window. The proposed algorithm include the following three steps:
Step 1: Let be the initialized matrix of new intensity image; Iteration Times, and Tolerance ;
Step 2: Calculate the parameters and using Equation (10); Estimate intensity image using Equation (7);
Step 3: if , output as the estimated intensity image; otherwise, go back to Step 2.
3. Numerical Experiments
3.1. Experiment Data
The QuickBird images are downloaded from http://www.digitalglobe.com. DigitalGlobe company provide commercial satellite QuickBird images, which contain one 0.6 m spatial resolution panchromatic image (450 - 900 nm) and four 2.4 m MS images: blue band (450 - 520 nm), green band (520 - 600 nm), red band (630 - 690 nm) and near infrared band (760 - 900 nm). A subset images which has the size of are cut from the original QuikBird images and are used as the experiment images. The MS images have been resampled to the same pixel size of PAN image. The experiment images are shown in Figure 2.
The second remote sensed images used in this paper is Landsat ETM+ images,
Figure 2. QuickBird experiment data. Left: PAN image; right: MS image, R: band 3, G: band 2, B: band 1.
which is downlodad at https://www.usgs.gov/. Panchromatic image of ETM+ sensor has the spatial resolution of 15 meter, while the multi-spectral bands have the spatial resolution of 30 meter. So, it is necessary to merge the abund spectral information in the mulit-spectral images into the panchromatic image to get the high resolution multi-spectral images. The images are shown in Figure 3.
3.2. Assessment Index
To give an objective assessment, correlation coefficients are used to assess the spectral distortion between the fused MS images and the up-sampled MS images, due to the shortage of original high resolution MS images. Correlation coefficient is defined as
Correlation coefficient measure the similarity degree of the same spectral band between fused image and original image. Its value should be as close to 1 as possible.
Another index is ERGAS (Erreur Relative Globale Adimensionnelle de Synthese)   or relative dimensional global error, which is defined as
where is the resolution of PAN image, is the resolution of MS image, is the mean radiance of each spectral band, RMSE is the root mean square error calculated using
Figure 3. Landsat ETM+ images: Panchromatic image (left), multi-spectral images (Right): R: band 5, G: band 4, B: band 3.
where NP is the total number of pixels in the original and fused image, and is the radiance value of pixel j in the ith band of the original image and the fused image, respectively. ERGAS is used to assess the spectral quality in the fused image, and the lower the value of the ERGAS, the higher the spectral quality of the merged image  .
4. Results and Discussion
4.1. Results of Fusing QuickBird Images
The outputs of applying different fusion methods to QuickBird images are shown in Figure 4. Firstly, it can be found by visual interpretation that there are more spatial detail information in the fused results compared with the original multi-spectral images. Spatial resolution of the fused results are much higher than original MS images. Most of the detail spatial structure in PAN image has been merged into the fused results. Some of the small spatial structure details which cannot be discerned from the original MS image (Figure 2), can be identified in the fused results. The results of fast IHS and wavelet fusion method have more sharp edge and texture than the results of Brovey fusion method, PCA fusion method and the proposed hybrid IHS fusion method, which can be verified by the correlation coefficients between the fused MS images and the PAN image.
It also can be found that the results from Brovey fusion, PCA fusion and fast IHS fusion, are severely disturbed by spectral distortion, which can be testified using the correlation coefficients between the fused MS images and the up-sampled MS images in Table 1. By preserving more spatial structure information in the fused images, wavelet fusion and the proposed hybrid IHS fusion generated more better results compared with the other fusion methods.
It can be found that the fast IHS results display higher spectral distortion compared with the results derived from hybrid IHS. The reason of spectral distortion of fast IHS has been investigated by Tu   . In the proposed hybrid IHS algorithm, the difference between PAN image and intensity image has been optimized selected to decrease the change of saturation image, which is critical to preserve the spectral information contained in the original MS images. Therefore, there is similar spectral characteristic between the hybrid IHS results and the original MS images.
In addition to the visual inspection, the performance of these two methods is further quantitatively analyzed using the assessment indexes. Firstly, the correlation coefficients verified that results from hybrid IHS have higher similarity to the original MS images compared with the results from fast IHS. Little spectral distortion emerged in the results of the proposed method, which can be seen by visual investigation. There is major difference in RMSE and ERGAS between the results derived from different methods. Results from HIHS have smaller RMSE and ERGAS than that of results from fast IHS, which demonstrate that hybrid IHS’s results have higher quality. Correlation coefficient to PAN of the hybrid IHS is less than that of fast IHS’s result, which demonstrates that there is short-
Brovey Fusion Wavelet Fusion PCA Fusion
Fast IHS Fusion Hybrid IHS Fusion
Figure 4. Image fusion results of QuickBird images using different methods.
Table 1. Values of different indexes to evaluate the quality of the fused QuickBird images.
age of spatial information in the results of hybrid IHS compared with that of fast IHS.
4.2. Results of Fusing Landsat ETM+ Images
In this subsection, the proposed hybrid IHS method together with other fusion methods are used to fuse the MS images and panchromatic image taken from Landsat ETM+ sensor. The fused Landsat ETM+ images are shown in Figure 5. It can be found by visual interpretation that results from wavelet fusion and hybrid IHS fusion induced better fusion images compared with the results from the
Brovey Fusion Wavelet Fusion PCA Fusion
Fast IHS Fusion Hybrid IHS Fusion
Figure 5. Landsat ETM+ fusion results using different fusion methods.
Table 2. Values of different indexes to evaluate the quality of the fused Landsat ETM+ images.
other methods. Spatial information in the fused images had been increased in some degree. However, significant spectral distortion emerged in the fused results of Brovey method, PCA method and fast IHS fusion method.
To give a throughout investigation of the proposed hybrid IHS fusion method, different indexes are used to assess the performance of these methods (Table 2). Correlation coefficients, which is used to assess the correlation relationship between the fused high resolution MS images and the low resolution MS images, demonstrated the degree of spectral similarity between two images. It was found that results from wavelet fusion and the proposed hybrid IHS fusion, which had slight spectral distortion, outperformed the other fusion methods. The other three indexes such as RMSE, correlation ceofficients to panchromatic image and ERGAS have demonstrated that the proposed hybrid IHS outperformed the other fusion methods.
In this paper, we give a hybrid of IHS and Minimum Mean-Square-Error for fusing low resolution multi-spectral and Panchromatic images from same scene. IHS is one of the commonly used fusion algorithms to merge the spatial information in PAN image and spectral information in LRMS images. However, spectral distortion phenomenon in IHS method seriously deteriorates the quality of the fused images. Reason of spectral distortion is due to the process of adding the difference between PAN image and intensity image directly into the original RGB images. To avoid or mitigate the influence of pixels that has bigger value compared with the ordinary pixels in the difference image, MMSE model is utilized to estimate the new intensity image from PAN and intensity images.
QuickBird PAN image and LRMS images are fused to evaluate the performance of our proposed algorithm. FIHS are used as reference to analyze the results from HIHS method. The comparison confirms that results from HIHS preserve most of spectral information with little spectral distortion, while results from FIHS have significant spectral distortion which has worse fusion quality. Therefore, the proposed hybrid method outperforms the commonly used FIHS method by providing higher quality fusion results.