## Abstract

This paper proposes an improved reflectance reconstruction method by adaptively selecting training samples. Modified Principal Component Analysis estimation was proposed by orthogonal regression considering the system noise; deriving the optimum number of training samples by BP-Adaboost neural network; and grouping the representative samples together by hierarchical cluster analysis from a large database of samples. Finally, the training samples were selected by colorimetric subspace tracking. Experimental results indicated that the proposed method significantly outperforms the traditional methods in terms of both spectral and colorimetric accuracy, and our reflectance modeling is a reasonable and convenient tool to generate adaptive training sets.

© 2017 Optical Society of America

## 1. Introduction

One of the ultimate goals of multispectral imaging is to achieve high spectral accuracy when representing the original appearance under various illumination conditions, which is widely used in art archiving, restoration, high fidelity reproduction, etc. Several multispectral reconstruction techniques, such as Wiener estimation [1,2], pseudo-inverse method, Principal Component Analysis [3–5], etc., have been extensively studied to accurately estimate the spectral reflectance of object surfaces. Recently, a great deal of work has been conducted to investigate the impact of training sets on the accuracy of spectral reconstruction, which has been mostly concerned with representative samples. Hardeberg proposed a method to select the most significant target patches from Munsell samples by minimizing the condition number, thus selecting samples that are different from other training samples already selected [6]. Cheung proposed a MAXMINC algorithm to select different samples in CIELAB space from a large set of samples for color characterization [7]. Mohammadi used agglomerative cluster analysis to select 14 targets from a data set, pointing out that the spectral properties are more important than the number of samples [8]. Shen proposed a sequential method for target selection according to spectral similarities by spectral distance [9,10]. Kruschwitz selected ten color block mirrors (narrow notch mirror) from 10 “Block Dye” designs as specialized color targets for reflectance reconstruction, which reduced the number of color patches required [11]. Eckhard used a recursive top-down algorithm to reject clusters of training samples that do not enhance the performance of a linear least-square regression process [12]. However, most research efforts have mainly concentrated on spectral accuracy, whereas color accuracy has been ignored. Therefore, it is necessary to develop a new method for adaptively selecting training data sets based on the characteristics of multispectral response.

In this study, we propose a method to reconstruct spectral reflectance based on modified principle component analysis by orthogonal regression considering the system noise, followed by investigating whether the accuracy of reflectance reconstruction could be improved if the training sets were adaptively selected, representative and correlated in spectral colorimetric subspace. Finally, the ideal number of training sets was predicted using BP-Adaboost neural network. The novelty of the proposed method is mainly the manner of training sample selection. The performances of both the proposed and traditional methods are compared in the case of a wide-band multispectral imaging system, which uses a modified trichromatic camera combined with two absorption filters.

## 2. Spectral reconstruction model

The spectral reflectance can be reconstructed using prior knowledge by measured color patches and digital signals. Suppose that the spectral reflectance of color patches has the dimension *n ×* 1, where *n* is the number of discrete wavelengths uniformly sampled over the visible wavelength from 400 *nm* to 700 *nm*. Based on principal component analysis [6], the spectral reflectance can be estimated from *p* scalar eigenvalues and basis eigenvectors calculated from measured spectra [5,6]. When *F* is the first *p* eigenvectors [*e _{1}*,

*e*,

_{2}*e*,…

_{3}*e*] and

_{p}*a*is the scalar eigenvalues [

*a*,

_{1}*a*,…

_{2}*a*], then the spectral reflectance can be expressed as a linear combination of eigenvectors and scalar vector, as given by Eq. (1)

_{p}*c*represents the multispectral imaging response, and the signal can be rewritten as a scalar product in matrix notationwhere

*M*denotes the transform matrix of spectral and multi-channel responses of the imaging system. considering the disturbance of imaging system and nonlinearity of the principal component basis eigenvector [5]. The matrix

*MF*can be rewritten as$MF=M{F}_{0}+\Delta MF$, where$M{F}_{0}$ and$\Delta MF$are the undisturbed matrix and disturbed matrix, respectively. The vector

*c*can be rewritten as$c={c}_{0}+\epsilon $, where ${c}_{0}$is undisturbed response vector,

*ε*is the noise. The above equation can be rewritten as Let$A=[-c,MF]$be the augmented matrix and$B=[\epsilon ,-\Delta MF]$ be the noise matrix [13]. The solution to the above homogeneous Eq. (4) can be formulated as seeking a solution vector

*a*, that is, minimizing the Frobenius norm of the denoted noise matrix

*B*.

*H*denotes the conjugate transpose, ${b}_{ij}$ is element of matrix

*B*, The optimal solution vector can be solved as

*A*. the notion$\u2020$ is Moor-Penrose pseudo-inverse. If the error disturbance is unrelated and has the same variance, $E\left\{{\left(\Delta MF\right)}^{T}\Delta MF\right\}={\sigma}^{2}I$, the estimated spectral reflectance vector can be given by Eq. (7)

## 3. Optimized selection of training samples

Previous studies have shown that the key to multispectral imaging is to find a proper transformation matrix from the camera response to spectral reflectance, which depends on the imaging system, calibration target samples, and matrix coefficients. In the traditional spectral reconstruction method, the transformation matrix was usually calculated from ensemble samples, random samples or common calibration targets, e.g., ColorChecker DC, ColorChecker Rendition Chart, Esser Test Chart TE221, etc. However, most methods focused on the widely representative training sets and ignored the chromaticity characteristic. Thus, a new method for training sample selection is proposed, assuring representation across the entire sample set and improving the correlation with the color target.

#### 3.1 Training set grouping based on hierarchical clustering

Hierarchical clustering groups samples by creating a cluster tree such that the samples in the same cluster are similar, whereas objects in different clusters are distinct. First, all of the training samples are merged according to their similarities using Chebyshev’s distance because of its better performance for distinguishing camera response. Second, similar samples are grouped into sub-clusters using Ward’s distance linkage method between all pairs of items. The mathematical expression for Ward’s distance between the two clusters is given as

*r*and

*s*,

*n*and

_{r}*n*are the number of elements in clusters

_{s}*r*and

*s*. The clustering process can be implemented by the Statistics toolbox of Matlab 2012a. Finally, the representative training sample of each cluster was selected by colorimetric subspace tracking analysis.

#### 3.2 Training sample selection based on colorimetric subspace tracking

The goal of this part is to select training samples that have similar chromaticity with target patches by subspace tracking. Let${P}_{n}(\lambda )$ be the spectral power distribution of illuminant, ${r}_{n}(\lambda )$be the spectral reflectance, ${V}_{n}(\lambda )$be the color matching functions, *n* be the sampling interval, the tristimulus values ${t}_{c}$ can be expressed as

*r*is a

*n*× 1 reflectance matrix, $\Theta $is a

*n*× 3 transformation matrix with only 3 nonzero singular values because$\Theta $is a column full rank matrix, i.e., $rank(\Theta )=3$, it can be calculated by singular value decomposition

*V*is a$3\times 3$ real matrix, ${U}_{1}=\left[{u}_{1},{u}_{2},{u}_{3}\right]$ and${U}_{2}=\left[{u}_{4},{u}_{5}\cdots {u}_{m}\right]$are the singular vector matrices corresponding to the first three non-singular values and

*m-3*zero-singular values, respectively. Then, the projection matrix${P}_{\Theta}$ of $\Theta $can be defined as follows:

The most similar targets can be selected as representative samples using spearman rank correlation coefficient in colorimetric subspace. The Spearman’s rank correlation coefficient defined as

*r*and

_{s}*r*are coordinate-wise rank vectors of fundamental stimulus${r}_{s}^{*}$and${r}_{t}^{*}$respectively, i.e., ${r}_{s}=({r}_{s1},{r}_{s2},\cdots ,{r}_{sn}),$${r}_{t}=({r}_{t1},{r}_{t2},\cdots ,{r}_{tn})$,${r}_{sj}$and${r}_{ti}$is tiedrank computed by matlab function

_{t}*tiedrank*.

*n*is the number of discrete wavelengths (usually

*n*= 31),

*r*and

_{s}*r*are constant (supposing

_{t}*n*= 31, that $\overline{{r}_{s}}=\overline{{r}_{t}}=16$). The training sample set can be expressed asCombining Eqs. (7) and (15), the reconstructed spectral can be expressed as

## 4. The number of training samples prediction

The prediction of training number for spectral reconstruction is frequently limited because it is difficult to confirm the best number for different target. The five-layer BP-Adaboost backpropagation neural network was proposed, the core of this approach is to establish predicting model on determining the weight of weak predictors, according to the predicted result by weighted test samples, and then weak predictor sequence used as the strong predictor, the signals is processed from the input layer to output layer in order to perform an optimization mapping by minimize the *rms* error between the desired output and actual multi-channel response.

With regard to the BP-Adaboost neural network, the authors believe the prediction accuracy is closely related to the characteristic of multi-channel signal and correlation with the training sets. The BP-Adaboost artificial neural network was built and trained by the optimization Toolbox of Matlab 2012a with the aim of predicting the number of training samples. Figure 1 shows that the network structure is 9-6-6-1-1, nine neurons in the input layer, two hidden layers with six neurons in each layer owing to its superiority in fitting the non-linearity of complex mathematical issues, and one output neuron. The available input data consist of nine parameters computed in the domain of 6-channel camera response, i.e., the Euclidean minimum distance, maximum distance, mean distance, centroid distance, ward distance, angle values between target and all training samples, and average values, norms, mean absolute deviation of the multispectral response, where the output denotes the best training numbers. The data were normalized before training and testing, and all the networks were trained 100 times with the maximum number of training epochs set to 100. The learning rate and performance goal were set to 0.1 and 0.00003, respectively. The TANSIG function was used as a transfer function, whereas the TRAINLM function was employed as a training function, LEARNGDM. LEARNGDM and MSE were used as the learning and performance functions. After the training process, the optimized number of training samples could be implemented by these trained networks.

## 5. Experiment & Results

To verify the proposed approach, comparative experiments were implemented. In the multispectral imaging system, we used a trichromatic camera with the built-in IR filter removed (Sinarback 75H color filter array, Sinar, Swiss) with 16-bit digitization, two wideband filters (GG475 and BG39, product of U-optic company, China), two Elinchrom CANLITE1000 sources with soft box placed 45° on either side of the sample plane, 6-channel multispectral images were acquired by the imaging system, and a gray-scale card was used for calibrating the lighting spatial non-uniformity. 1641 mineral pigment patches were used as color targets, 844 patches were randomly selected for training and another 802 for testing. The reflectance was measured by an X-Rite Spectroeye spectrophotometer. The transform between camera response and reflectance was calculated from the training set and evaluated on the test set. The hierarchical clustering was carried out in the domain of entire training samples’ camera response,and 6-channel camera signals for each channel were corrected using the gain-offset-gamma(GOG) model, and then converted to tristimulus values. The estimated reflectance was evaluated using spectral as well as colorimetric metrics. Two different metrics [4], goodness fitting coefficient (GFC) and root-mean-square (*rms*), were used as spectral metrics, with CIE DE76 and CIE DE2000 (color difference) under D50 CIE standard illuminant as the colorimetric metric. The GFC ranges from 0 to 1, with 1 corresponding to the perfect estimation. The *rms*, CIE DE76 and CIE DE2000 are positive values from 0 and higher, with 0 corresponding to the perfect estimation. These metrics are given by the equations:

First, the reflectance of 562 patches were reconstructed using 844 training samples by the proposed method and nine parameters, where the number of training samples corresponding to the best reconstructed spectral were used as training parameters of BP-Adaboost neural network and the remaining 240 samples were used as test samples. The performance of the proposed method was compared with Mohammadi’s method [8], Hardeberg’s PCA estimation with 20 training samples [6], respectively, and Wiener estimation. Tables 1 and 2 give the performance of spectral and colorimetric accuracy. Obviously, our method outperformed the other methods, which clearly demonstrates that the proposed method approximately produces the minimum spectral *rms* error and color difference through the optimized training samples. Other methods exhibit a large discrepancy without considering the over-fitting effect caused by an excessive number of training samples, such as the Wiener estimation and the PCA algorithm. Figure 2 shows the boxplot distribution of the color difference and spectral *rms* error compared with other estimation methods, of which most of the CIEDE2000 and spectral *rms* are less than those of other traditional methods.

Figures 3-6 demonstrate the distribution of training samples in principal component coordinate space, multi-dimensional space, xyY space and a-b chromaticity space, respectively. The six subfigures were random samples (No. 90, 132, 143, 188, 202, 227) of all 240 reconstructed samples, corresponding to the number of training samples (7, 5, 12, 3, 12, 11) predicted by the neural network. The figure shows that the training samples selected by Mohammadi’s method are relatively concentrated in a-b space, located in green and blue areas, with a small similarity between the selected samples and target. In multidimensional scale space and principal component space, the selected training samples are more concentrically distributed and too many samples are correlated, which could cause over-fitting in spectral reconstruction and reduce the reconstruction accuracy, leading to a smaller mean distance between samples and stronger similarity. In addition, the samples selected present strong linear correlation in principal component coordinate space and multi-dimensional space. Our method achieves improved estimation accuracy by improving the representative and correlation in the response and reflectance matrices.

The accuracies of reconstructed spectral samples using different numbers training samples were investigated. Table 3 lists the spectral good fitting coefficient and CIEDE2000 color difference errors. Obviously, our method outperforms other methods with a fixed number. However, it is not safe to conclude that more training samples are needed to obtain more accurate results. Excessive training numbers would cause over-fitting and thus reduce accuracy. This finding is consistent with experimental data. A histogram exhibiting the distribution of ideal training sample numbers is depicted in Fig. 7(a), which clearly illustrates that the accuracy improved when the number of training samples was increased from 3 to 13.

The reflectance of the entire collection of 240 test samples was shown in Fig. 7(b). The reflectance was sampled at 10 *nm* intervals in the range of 400 *nm*-700 *nm*. The reconstructed reflectance of the proposed method, Wiener estimation, and Hardeberg PCA are shown in Fig. 8, where the black line is the measured spectral reflectance. The reflectance reconstructed by the proposed method was found to be more accurate than those of the traditional method.

## 6. Conclusion

This paper proposed a method for spectral reconstruction by optimizing the training samples. The procedure was divided into three main stages. In the first stage, the ideal number of training samples was predicted by BP-Adaboost neural network. In the second stage, the training samples were grouped into sub-clusters. In the third stage, the training samples were selected by colorimetric subspace tracking. The feasibility of the proposed technique was examined in experiments using multispectral imaging with mineral pigment patches. The results show that in terms of spectral and colorimetric metrics, the proposed technique is superior to Mohammadi’s method and other traditional methods, such as Wiener estimation and PCA estimation. However, there are still problems to be solved in the future. Our method is more computationally expensive than traditional methods in training procedure, which is the most common shortcoming for nearly all adaptive methods. The proposed method has potential applications in textile, printing, and other industries.

## Funding

National Natural Science Foundation of China (NSFC) (61405104, 61505149); Key Basic Research Project of QFNU (Xkj201321); Chinese National Undergraduate Innovation Project (201610446070).

## Acknowledgments

Thanks are due to X. X. Wan for valuable discussion and to J. X. Liang for assistance with the experiments.

## References and links

**1. **J. H. Yoo, D. C. Kim, H. G. Ha, and Y. H. Ha, “Adaptive spectral reflectance reconstruction method based on wiener estimation using a similar training set,” J. Imaging Sci. Technol. **60**(2), 020503 (2016). [CrossRef]

**2. **Y. Murakami, M. Yamaguchi, and N. Ohyama, “Piecewise Wiener estimation for reconstruction of spectral reflectance image by multipoint spectral measurements,” Appl. Opt. **48**(11), 2188–2202 (2009). [CrossRef] [PubMed]

**3. **X. Zhang and H. Xu, “Reconstructing spectral reflectance by dividing spectral space and extending the principal components in principal component analysis,” J. Opt. Soc. Am. A **25**(2), 371–378 (2008). [CrossRef] [PubMed]

**4. **R. Shrestha, A. Mansouri, and J. Y. Hardeberg, “Multispectral imaging using a stereo camera: concept, design and assessment,” Eurasip. J. Adv. Sig. Pr **2011**(1), 57 (2011). [CrossRef]

**5. **Z. Liu, X. X. Wan, X. G. Huang, Q. Liu, and C. Li, “The study on spectral reflectance reconstruction based on wideband multi-spectral acquisition system,” Guangpuxue Yu Guangpu Fenxi **33**(4), 1076–1081 (2013). [PubMed]

**6. **J. Y. Hardeberg, “Acquisition and reproduction of color images: colorimetric and multispectral approaches,” Ph.D dissertation (Ecole Nationale Superieure des Telecommunications, 1999).

**7. **V. Cheung and S. Westland, “Methods for optimal color selection,” J. Imaging Sci. Technol. **50**(5), 481–488 (2006). [CrossRef]

**8. **M. Mohammadi and M. Nezamabadi, “A Prototype Calibration Target for Spectral Imaging”, in 10th Congress of the International Colour Association (2005), pp. 387–390.

**9. **H. L. Shen, P. Q. Cai, S. J. Shao, and J. H. Xin, “Reflectance reconstruction for multispectral imaging by adaptive Wiener estimation,” Opt. Express **15**(23), 15545–15554 (2007). [CrossRef] [PubMed]

**10. **H. L. Shen, H. G. Zhang, J. H. Xin, and S. J. Shao, “Optimal selection of representative colors for spectral reflectance reconstruction in a multispectral imaging system,” Appl. Opt. **47**(13), 2494–2502 (2008). [CrossRef] [PubMed]

**11. **J. D. T. Kruschwitz, “Specialized Color Targets for Spectral Reflectance Reconstruction of Magnified Images,” Ph.D dissertation (Rochester Institute of Technology, 2015).

**12. **T. Eckhard, E. M. Valero, J. Hernández-Andrés, and M. Schnitzlein, “Adaptive global training set selection for spectral estimation of printed inks using reflectance modeling,” Appl. Opt. **53**(4), 709–719 (2014). [CrossRef] [PubMed]

**13. **M. D. Rahman and Y. Kai-Bor, “Total least squares approach for frequency estimation using linear prediction,” IEEE. T. Acoust. Spee. Signal. Proces **35**(10), 1440–1454 (1987). [CrossRef]