Principal wavelengths in the formation of spectral images of natural scenes

Abstract. Considering the high degree of correlation in the visible spectrum, the principal wavelengths from spectral measurements of radiance recorded in spectral images were selected using a method based on principal components analysis (PCA). It seems to be that this is the first time that, instead of using spectra, data is taken directly from the “slices” of spectral images; the method has the advantage of preserving the structure of the original data in the reduced data set. A “true” dimensionality of five wavelengths resulted for all the analyzed images. The averages of the selected wavelengths for 10 spectral images produced good results for a human observer. These results were possible using only four wavelengths. Though PCA by itself is not able to include the impact of specific sensors on the selection of basis functions, results suggest that the variable selection method used in this work (which is not just PCA) yielded objective information of the structure of the physical stimuli (i.e., the spectral structures) that have been shaping the visual systems of animals and insects since many years ago.


Introduction
Since the beginnings of color science, two basic tendencies have been present with its development. The first has to do with the idea that color vision is almost exclusive to humans, and the second implies that only a reduced part of the visible spectrum is required to achieve the perfect perception of colors. Both tendencies have been supported by relevant findings and theories. This work is concerned with the second tendency.
Color science was founded on human perception, even in the definition of light, but especially with Maxwell's determination of the color-matching functions of the normal human eye. 1,2 Fortunately, by the end of the 18th century, it was recognized that light is part of a wider spectrum and that "plant and animal colors do not exist simply for human pleasure." 3,4 The idea that combining a reduced number of basic colors is enough to generate any other color seems to have been introduced by Thomas Young. He stated that "it is almost impossible to conceive each sensitive point of the retina to contain an infinite number of particles, each capable of vibrating in perfect unison with every possible undulation, it becomes necessary to suppose the number limited, for instance, to three principal colors, red, yellow, and blue, etc." 5 The same idea was also suggested by Maxwell in his pioneering studies on color. Maxwell concluded that all the colors of the spectrum could be compounded from red, green, and blue, which he called the primaries. 1 However, since long ago, studies on color have shown that there are no three spectral "primary" colors from which all other colors can be compounded. 6 Apparently Maxwell was the first to consider the problem of determining the laws of the composition of colors, so that the number of "standard" colors was reduced to a minimum. The determination of the minimum number of components required for color vision has been assessed in different ways and from different perspectives. One approach has been supported by physiological studies on animal vision in general, considering humans as a particular case. 3,7 On the other hand, principal components analysis (PCA) has been successfully used for determining the level of redundancy or the degree of importance of visual information related to color vision. For instance, Judd et al. 8 applied PCA to study the spectral distribution of daylight from experimental determinations. They found that only two components were required to reconstitute the measured curves with a good approximation. Cohen 9 made use of PCA to obtain the best linear model that fit a subset of 150 spectral reflectances of Munsell color samples. He found that the first three characteristic vectors accounted for 99.18% of the total variance of the original data. Maloney 10 used PCA to show that within the visible spectrum, three components are sufficient to represent most of the variation in spectra of "natural formations." But he concluded that at least seven components are required to optimally reconstruct the original data. Parkkinen and Jaaskelainen 11 showed that color spectra can be accurately reconstructed using a few of what they called "principal spectra." Wandell 12 proposed a method for the analysis and synthesis of color images taking advantage of the high degree of correlation across the visible wavelengths. He suggested the use of PCA to model spectral reflectance and illumination. Parkkinen et al. 13 analyzed 1257 reflectance spectra in the visible region of the Munsell color chips. They found that eight characteristic spectra are needed to achieve a good representation for all spectra. Oxtoby and Foster 14 did psychophysical discrimination experiments using images of Mondrian-like patterns of Munsell surfaces and their spectral approximations produced by PCA. Observers required at least five basis functions for discrimination performance to be a chance occurrence.
It is worth noting that, except for Oxtoby and Foster, all the authors mentioned above worked with spectra, while the data analyzed in this work was taken directly from the spectral images.
As mentioned above, several research studies have found a high degree of correlation in the visible spectrum. Taking this as starting point, this work is not simply an attempt to assess the optimum dimensionality for the formation of spectral images of natural scenes, but also to determine specifically which of the several combinations of the available wavelengths are the most important for image formation. The used statistical methods do not consider specific observers. However, due to the importance of human vision, several results are given as could be perceived by a human observer.
An important point, already mentioned above, is that in this work, instead of using spectral data, data was taken directly from the "slices" of spectral images. To our knowledge, this is the first time an analysis of this kind has been conducted.

Spectral Images and the Color Signal
The spectral images used in this work were acquired as described by Brelstaff et al., 15,16 and made available by the same authors. Each spectral image consists of a sequence of 31 chromatically narrow-band-filtered 8-bit, 256 × 256 pixel images. These images were acquired sequentially using a set of 31 optical interference filters in the range from 400 to 700 nm with 10-nm spacings. The images include common natural objects like plants, trees, flowers, grasses of different colors, the sky, shapes and textures.
The images produced by the camera represent spectral measurements of radiance, that is, measurements of the power of light reflected by the surfaces of the objects in the scenes at specific wavelengths. In this work, the color signal corresponds to the radiance measured with the camera.
Using the notation suggested by Wandell, 12 color information is acquired at wavelengths λ n , for n ¼ 1, N. The color signal is obtained by multiplying the spectral power distribution of the illuminant at point x, E x ðλ n Þ times the surface spectral reflectance at x, S x ðλ n Þ. Therefore, the color signal is Considering that the response of the photoreceptor of a spectral camera is Rðλ n Þ, the response of the camera can be modeled as But as shown by Brelstaff et al., 16 the spectral response of the spectral camera used to acquire the images is practically constant from 400 to 700 nm. Given that, the spectral response of the camera can be considered as simply the color signal given in Eq. (1), but multiplied by a scalar. This scalar multiplication has no effect in the reconstruction of the RGB images, due to the normalization that is part of such a process.

Determination of the Principal Wavelengths
All the authors mentioned in the introduction found a high degree of redundancy across the visible wavelengths. This suggests that dimensionality reduction in spectral images is possible. However, a major deficiency of PCA for dimensionality reduction of multivariate data is that, while the dimension of the data space may be reduced from, say, p to q, all the p original variables are, in general, needed to define the q new variables. Thus, dimensionality reduction along with the identification of the most important, or principal, variables is highly desirable, especially if that reduced number of variables conveys the main features of the whole sample.

Principal Components Analysis
PCA is frequently used to uncover patterns that may be exhibited in a given set of data. This technique is also useful in determining if any given information is unnecessary or redundant. When unnecessary or redundant variables are identified, dimensionality reduction is possible. For instance, if p variables x 1 , x 2 . . . , x p , are observed on n individuals (n > p), their corresponding values are grouped in an (nxp) data matrix X.
The sample covariance matrix is given by ðn − 1Þ −1 X 0 X, and the sample correlation matrix can be written in the same form if the columns of X have previously been standardized. Principal component analysis generates p new variables y 1 , y 2 . . . , y p , which are uncorrelated and ordered in such a way that the first few explain most of the variation (i.e., useful information) contained in all the original variables. All this is achieved without disturbing the overall features of the sample. This statistical method implies the calculation of the eigenvalues and eigenvectors of the sample covariance or correlation matrix. The first k eigenvalues indicate the amount of variance of the original data that is explained when only k < p variables are considered. The corresponding eigenvectors are used to obtain the new variables.

Determination of the Minimum Number of Components to be Retained
In a given multivariate study, the dimension of the original data is given by the number of variables that are considered at the outset. When using PCA, the first step for selecting the most important variables is to determine the minimum number of principal components to be retained. This number serves as a guide for dimensionality reduction and for determining the number of important variables (wavelengths). Several criteria have been suggested to determine the number of principal components to be retained so that the original data is appropriately represented. Among the most used are: retaining components with eigenvalues greater than 1.0 (λ s > 1.0), the Scree plot, the broken stick model, and, keeping components with eigenvalues totaling to a fixed amount of the total variance. 17 For the sake of objectivity, in this work, the stopping method proposed by Eastment and Krzanowski was applied. 18 This method is based on a cross-validation procedure which assumes that a suitable choice is a number of components for which adequate prediction of the original data matrix is possible. The method is as follows: assuming that the observations of p variables (x 1 , x 2 . . . , x p ) on n individuals (n > p) are the columns of a (nxp) data matrix X, the singular value decomposition (SVD) of X is defined by and S ¼ diagðs 1 ; s 2 ; : : : ; s p Þ. If u ij is the (i, j)'th element of U, then the (i, j)'th element x ij of X is given by where v ij is the (i, j)'th element of V 0 and s i is the corresponding principal (or eigen-) values. When PCA is used to reduce the dimensionality of a problem, determining how many components are needed to reproduce the data with the desired accuracy is required. Thus, retaining the first M components can be considered as modeling the elements of X by where ε ij is an error or residual term. The aim of the method is to determine the optimal value of M using cross-validation. In the classic cross-validation procedure, the data matrix is divided into subgroups. Each subgroup is deleted from the data in turn, and the parameters of the predictor are estimated from the remainder of the data (for the model in turn). Finally, the deleted values are predicted for the model. A suitable object function relating actual and predicted values is used and the model that optimizes such a function is selected. In our case, for each possible choice of M, a predicted valuex M ij of x ij is calculated (i ¼ 1, . . . , n; j ¼ 1, . . . , p), and the discrepancy between actual and predicted values is calculated with To avoid using the same data for its own prediction, the corresponding row and column are deleted and the prediction is calculated using the following equation: whereû it ,ŝ t ands t ,v tj are calculated using the Singular Value Decomposition when the i'th row and the j'th column of X are eliminated, respectively. At this point, it is worth noting that Eastment and Krzanowski, 19 due to hardware limitations, used an updating algorithm for estimating the singular value decomposition. In this work, the singular value decomposition was actually calculated for predicting each x ij . That is, 1040 singular value decompositions of a 130 × 8 matrix were calculated.
To select the optimum value of M, the following statistic is calculated: where D M represents the degrees of freedom required to fit the M'th component and D R represents the remaining degrees of freedom after fitting the M'th component. Taking into account the number of parameters to be estimated and the constraints on the eigenvectors at each step, D M ¼ n þ p − 2M. Given that there are np − p degrees of freedom at the outset (considering the columns of X being mean-centered), D R can be easily calculated at each stage by successive subtraction. Eastment and Krzanowski 19 suggested that the optimum value for M is that at which W M becomes greater than unity or when the value of W M stabilizes. W M represents the increase in predictive information supplied by the M'th component divided by the average information in each of the remaining components.

Selection of Wavelengths
The variable selection method proposed by Krzanowsky 19 was used in this work because its optimality criterion implies the best subset of variables which reproduce as closely as possible the general features of the complete original data. Such optimality criterion implies direct comparisons between the multidimensional locations of individual points of the principal components of the subset configuration and the corresponding locations of points of the principal components configuration produced by the complete original data. This type of comparisons is called Procrustes analysis. 20 In our context, wavelengths correspond to variables and the respective reflectances recorded at each pixel of the images correspond to responses to those variables. Thus, redundant wavelengths are those which do not convey additional visual (color) information to that provided by the principal wavelengths. The purpose of this work is not only to determine how many but also which are the most important wavelengths in the image formation process. Thus, we are considering light as a compound stimulus such that several of its components are redundant for image formation.

Implementation of the Variable Selection Method
For determining the principal wavelengths in the spectral images, representative 64 × 64 pixels subimages were carefully selected from each "slice" of the images (Fig. 1). This means that, for each spectral image, 31 samples (each of 4096 pixels) were used. This is equivalent to considering the application of each of the 31 variables (wavelengths) to 4096 cases. The subimages containing the greatest amount of colors from the imaged scene were selected. For building the covariance matrix required for PCA, the pixel values of the subimage for each wavelength were entered as a column of the matrix following a lexicographic order in the subimage. The result was a 4096 × 31 covariance matrix for each spectral image. The same matrix was used for determining the minimum number of components to be retained and for the selection of variables (wavelengths). This approach is different from those used in previous publications where spectra were used in the analysis.

Image Reconstruction for a Human Observer
Once the principal wavelengths for each spectral image were determined, RGB images were reconstructed using only the required information for the reduced number of wavelengths. This was done to provide a "visual" representation of the results, and to investigate the suitability of the resulting images for human vision.
For generating the RGB images, the spectral images were structured as 3-D matrices whose dimensions were R × C × L. The first two were spatial dimensions (R ¼ 256 pixels, C ¼ 256 pixels, or R ¼ 64 and C ¼ 64 for subimages), while the third was the spectral dimension (L ¼ 31 wavelengths, or a much smaller number of selected wavelengths). When the matrix containing the CIE 1931 color-matching functions 21 was used to generate the RGB images, only the rows that corresponded to the selected wavelengths were used together with the respective "slices" of the spectral images. This way, only the information related to the selected wavelengths was introduced for each reconstruction. In Fig. 1 is shown a grayscale version of an RGB image that was reconstructed using the 31 "slices" of the Table 1 The first five eigenvalues of 10 spectral images, and the corresponding percentages of variance that are explained by principal component. spectral image and the corresponding information of the CIE 1931 color-matching functions for the respective 31 wavelengths.

Results
The procedure for determining the principal wavelengths described above was applied to 10 spectral images of natural scenes. These 10 images were in turn selected from a set of 29 spectral images. As will be shown, the obtained results were consistent for all the 10 spectral images. However, for simplicity, some discussions will refer to the results of two representative images used as examples.
When the stopping method proposed by Eastment and Krzanowski 18 was used, a "true" dimensionality of five variables resulted for the 10 analyzed images. Table 1 shows the first five eigenvalues of the covariance matrix of the 10 spectral images used in this work. The first five eignenvalues correspond to the image shown in Fig. 1. As can be seen, the first three principal components of each image explain around 95% of the total variance, and the first five explain around 98% of the total variance. This supports the results of the stopping method which, in all cases, suggested the retention of five principal components.
As mentioned above, the statistical method that was used for selecting the most important wavelengths has the advantage of preserving the structure of the original data. To show that the method actually worked as expected, in Fig. 2 are plotted the first versus the second and third principal components of the covariance matrix of the spectral image represented in Fig. 1. For Fig. 2(a), the principal components were calculated using the entire 4096 × 31 covariance matrix, while for Fig 2(b), the principal components were calculated using a reduced 4096 × 5 covariance matrix after applying the variable selection method. As can be seen in these figures, the structure of the data is well preserved after removing 26 wavelengths of the original data set.
At this point it is convenient to note that all images were constructed using the same exact procedure. The reference images [Figs. 3(a) and 4(a)] were reconstructed using the information associated with the 31 images, one per wavelength, and the corresponding 31 rows of the CIE 1931 color-matching functions, from 400 to 700 nm.
The results of applying the variable selection method are given in Table 2. As can be seen, the selected wavelengths for each spectral image are more or less evenly distributed along the visible spectrum. That is, no significant concentration of selected wavelengths was obtained in any range of the spectrum.
However, it is worth noting the constant appearance of wavelengths close to the long wave limit of the visible spectrum. In Figs. 3 and 4 below, the reconstructions of two of the 10 analyzed spectral images are shown.
For Figs. 3 and 4, Figs. 3(a) and 4(a) are the RGB reconstructions that resulted when all the 31 wavelengths were used, and are considered to be the reference images for a human observer. Figures 3(b) and 4(b) are the RGB reconstructions obtained  when only the five wavelengths selected by the variable selection method were used. For instance, for reconstructing Fig. 3(b), only the data related to 550, 580, 600, 680, and 690 nm was used. As can be seen in Figs. 3(b) and 4(b), the wavelengths chosen with the selection method produced poor "visual results" (poor RGB images). Considering that for all the analyzed images the selected wavelengths were rather evenly distributed in the visible spectrum, only the averages of the selected wavelengths (from the shortest to the longest for the 10 images) were used. The averages were calculated following the order shown in Table 2. In the first step, the magnitudes of the shortest wavelengths for each of the 10 images were added and the result was divided by 10, then the next 10 wavelengths were added and the result was divided by 10. The same procedure was applied to the remaining wavelengths. The results were: 470, 505, 557, 637, and 672 nm. Using these averages as a reference, the images that appear in Figs. 3(c) and 4(c) were produced using the wavelengths 470, 500, 560, 640, and 670 nm, for which the required information for image reconstruction was available. As can be seen, the resulting images represent a considerable improvement from the point of view of a human observer, though they are still easily discriminated from the reference images. These results suggested that few wavelengths evenly distributed in the visible spectrum could produce acceptable RGB images. Therefore, the images shown in Figs. 3(d) and 4(d) were generated using the data corresponding to 450, 500, 550, 600, and 650 nm. In all cases, a noticeable improvement was obtained, at least for a human observer (i.e., in the RGB images).
In several applications of PCA, dimensionality reduction is considered appropriate even when only a little more than 80% of the data variability is explained by the retained principal components. The percentages given in Table 1 suggest that less than five wavelengths could be enough to produce better results. Consequently, the images of Figs. 3(e) and 4(e) were obtained using only four wavelengths. The best improvement was obtained in Fig 3(e), where only the information relative to 450, 500, 550, and 600 nm was used. This result was even better than the results obtained with five wavelengths. Several tests were completed to generate RGB images with only three wavelengths, but no improvement, at least in the RGB images, was obtained. The best results with this approach are shown in Figs. 3(f) and 4(f). The wavelengths used for Fig. 3(f) were 450, 550, and 650 nm. The wavelengths used for Fig. 4(f) were 450, 600, and 660 nm. Those wavelengths were chosen because they are close to the wavelengths that produce the maximum individual responses of the human cones (L, M, and S).
A probable alternative to the PCA variable selection method could imply the analysis of the energies of the "slices" that constitute each spectral image. To investigate this, the energies of the 31 "slices" of the spectral image represented in Fig. 1 were calculated and plotted in Fig. 5. As can be seen, the wavelengths selected with the variable selection method for that spectral image (550, 580, 600, 680, and 690 nm) correspond to "slices" with the highest energies of their neighboring wavelengths.   5 Total energy of each of the 31 "slices" that constituted the spectral image represented in Fig. 1.
However, it is clear that it would be difficult to decide which five to choose just from the information of the energies of the slices. This shows the importance of the role of the variable selection method.

Discussion
The method for dimensionality reduction suggested that the "true" dimensionality of all the analyzed spectral images was five.
That is, five wavelengths should be enough to reconstruct any of the 10 spectral images with an acceptable quality. Taking into account that only 10 images were sampled and that the intensity of the "slices" of the spectral images were not manipulated to improve the results, it can be said that the images shown in Figs. 3(d), 3(e), 4(d), and 4(e) support this result. This result coincides with that of Oxtoby and Foster, 14 who did not analyze spectra, but did psychophysical discrimination experiments, presenting images to human observers that were reconstructed in the same way it was done here.
On the other hand, assuming that the recorded reflectances in the spectral images are perceived by many species of animals, it is worth mentioning that vertebrates have five classes of visual photopigments. 22 Furthermore, several groups, including birds and many lizards, retain all five classes. Consequently, considering that the spectral sensitivity of a photopigment is determined by its peak, λ max , 23 it can be said that these findings also support our dimensionality reduction results. Moreover, it has been consistently confirmed that spectral information can be coded with, at most, five types of receptors. [24][25][26] In addition, models indicate that animal color vision, involving five or fewer broadly tuned receptors, is well matched to most natural spectra. 27 With regard to the consistent appearance of a selected wavelength close to the 700 nm limit of the visible spectrum, it is worth noting that, according to Endler, 28 the reflectance spectra of green leaves shows one peak at about 550 nm and another one at about 700 nm. The fact that all the analyzed images in this work contain a considerable proportion of green areas explains the appearance of those wavelengths. This also partly explains why the images reconstructed with the five averaged wavelengths appear rather greenish [Figs. 3(c) and 4(c)]. Interestingly, for many species, green is one of the main colors associated with their food, which has always been a strong factor for adaptation.
According to Lythgoe and Partridge, 24 the long-wavelength limit for vision is likely to be set by the absorption of porphyropsin with a λ max of about 630 nm, which allows useful sensitivity to light longer than 740 nm.
An important aspect of this study was that the intensity of the slices of the spectral images was used as recorded by the spectral camera. Better results can probably be obtained, even with less than five wavelengths, by manipulating the proportions with which each slice contributes to the reconstruction of the RGB images.
Finally, though PCA by itself is not able to include the impact of specific sensors on the selection of basic functions, 29 our results suggest that the variable selection method that was used in this work (which is not just PCA) yielded objective information of the structure of the physical stimuli (i.e., the spectral structures) that have been shaping the visual systems of animals and insects since many years ago.

Conclusions
The variable selection method used here achieved the identification of redundant wavelengths and the selection of subsets of wavelengths without disturbing the overall features of the original data.
Considering that the response of the spectral camera used to acquire the images has a constant response all along the visible spectrum, our results suggest that five wavelengths, evenly distributed in the visible spectrum, should be enough to achieve color vision by different living beings.
Fine tuning in the perception of color might imply mechanisms that modulate the degree of participation of individual wavelengths depending on their mutual interaction.