8694
Ind. Eng. Chem. Res. 2010, 49, 8694–8704
Multiblock-Based Qualitative and Quantitative Spectral Calibration Analysis Chunhui Zhao and Furong Gao* Department of Chemical and Biomolecular Engineering, The Hong Kong UniVersity of Science and Technology, Clear Water Bay, Kowloon, Hong Kong SAR
In this article, an improved spectral calibration and statistical analysis approach is proposed. Having realized the multiplicity of underlying spectral characteristics and their different effects on the quality of interpretation and prediction, the major task lies in how to qualify and quantify the descriptor-quality relationship more meaningfully over different wavelength subspaces. The underlying spectral information can be explored more comprehensively by a spectral subspace separation and multiblock modeling strategy. Unlike the specific purpose of quality predictions over different spectral subspaces, systematic information in both descriptor and quality spaces is decomposed into different parts under the control of different between-set relationships. Closely related variations and irrelevant ones are discriminated and evaluated separately, where, in particular, the orthogonal variations, as important systematic information, are quantified in specific model parameters, even though they are not predictive. In this way, the approach well tracks the wavelength-varying underlying characteristics, takes full advantage of underlying variations in both X and Y spaces, and allows for a more meaningful interpretation of the between-set relationship. This article theoretically and experimentally illustrates the efficiency of the proposed method in qualitative and quantitative spectral calibration analysis. 1. Introduction During the past decade or so, the use of spectroscopic information1-19 has received much attention and begun to emerge as an important technique that is being heavily encouraged and practiced for different purposes. Predicting a dependent variable from spectral measurements is a frequently encountered problem in chemometrics. To analyze the composition of a chemical mixture, a calibration model is often constructed using spectra of the mixture together with the reference concentrations of constituents to form a quantitative predictive relationship.10-19 Common multivariate calibration methods20-24 used for spectral analysis include principal component regression (PCR) and partial least-squares (PLS) regression. Such methods are based on the fact that variable collinearity is typical in spectral data, in which some wavelength variables can be represented as a linear combination of some latent variables (LVs). Spectral data often comprise more variables (wavelengths) than observations (samples), so to deal with this problem of high data dimensionality and redundancy, PCR and PLS reduce the number of spectral variables by means of LV extraction, so that the original spectral data space is shrunk to a feature subspace of smaller dimension. Those underlying features are then used for regression modeling. Moreover, although calibration models are typically developed by including all available wavelengths, both theoretical and experimental evidence exists to demonstrate that it is possible to enhance prediction performance through the implementation of variable selection,8,9,12,19,25,26 also termed wavelength selection. The assumption is that there might be parts of the wavelengths that contain little information about the analyte properties. When these wavelengths are included in the regression model, predictive performance on unseen test data might be poor. Therefore, variable selection in multivariate analysis is a very relevant step, because the removal of uninformative variables will produce better predictions and simpler models. In the past few years, several techniques devoted to feature selection in PLS models applied to spectral data have been presented, such as interactive variable selection (IVS),25 * To whom correspondence should be addressed. Tel.: +85223587136. Fax: +852-23580054. E-mail:
[email protected].
uninformative variable elimination (UVE),8 and iterative predictor weighting (IPW).26 The selection of variables for multivariate calibration can be considered as an optimization problem. The selection result directly condenses the calibration model structure by reducing the dimensionality of the regression coefficient matrix. However, with the elimination of some wavelengths that are judged to be uninformative for quality prediction, the model might suffer from potential information loss problem compared with the original spectral space. Actually, these removed wavelength variables might also cover some underlying information about quality description, though not very much. Sometimes, they might even be an indispensable part of the complete description of the quality properties. Compared with direct variable elimination, orthogonal signal correlation (OSC)27-30 was designed to collect underlying components scanned over all wavelengths that account for as much as possible of the variation in the descriptor matrix (X) while remaining orthogonal to the property matrix (Y). Consequently, OSC filtering indirectly predigests the calibration model by reducing the required number of latent variables, thereby making the interpretation of the model easier. Without curtailing any wavelengths from the predictors, this approach avoids the problem of missing information. Considering that the spectra of mixtures are actually often linear combinations (with coefficients corresponding to the proportions) of the spectra of their constituent species,10,31-34 it would be useful if the component spectra could be recovered from the spectrum of a mixture. The goal is the assessment and resolution of all possible conformations. None of the abovementioned methods can properly identify the unknown species in a mixture, as well as their concentrations. Recently, the techniques of evolving factor analysis (EFA)35,36 and independent component analysis (ICA)37 have been developed to allow for the resolution of concentration profiles and pure spectra of the different species from spectral measurement but using different calculation procedures. EFA can extract the equilibrium concentrations of the species present in a solution, as well as their spectra, when there exists a physically meaningful ordering (“evolution variable”) of rows or columns in the data matrix.
10.1021/ie100892y 2010 American Chemical Society Published on Web 08/11/2010
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
A brief description of the EFA calculation algorithm has been presented by Tomisic and Simeon.35 Moreover, details can be found in the book by Malinowski.36 ICA37 is not generally sufficient to consider only up to second-order information (correlations, covariances) when describing the measured variables. It attempts to recover statistically independent signal sources, instead of merely uncorrelated ones, given only observations that are assumed to be linear mixtures of the original signals. To achieve this purpose, ICA makes full use of higher-order statistics and applies criteria related to information theory and entropy. In probability theory, independence is a high-order statistic that guarantees uncorrelatedness so that it is a much stronger condition than uncorrelatedness. It is designed to extract the constituent species in the form of independent components (ICs), as well as determining their effects on the observed mixture spectra, which is called the blind source signal separation process. Previous works10,31,32,34 have applied independent component analysis (ICA) to spectral data and successfully demonstrated their effectiveness in recovering the components of interest from the spectra of mixtures and estimating their concentrations. The independent component regression (ICR) method was first proposed by Chen and Wang31 for application to near-infrared (NIR) spectral data analysis. The authors pointed out that, by comparing the spectra of separated ICs with a library of the spectra of pure substances, it was possible to identify unknown constituent species existing in the mixture. Ideally, if the separated ICs exactly matched the pure substances constituting the mixture, then the mixing matrix would agree well with the concentrations of the substances existing in the mixture. Considering that they did not match very well in practice, regression analysis could be readily conducted between the estimated mixing matrix and the real concentration measurements based on the simple leastsquares algorithm. Shao et al.32 further reported that, in addition to the completely equivalent quantitative prediction performance compared with PCR, ICs could provide more chemical explanations than principal components (PCs), which were found to be strongly correlated to the NIR spectra of the source components in the spectra. Zhao et al.34 further applied ICA for the analysis of interference-subject spectral analysis and developed a robust calibration modeling strategy. Its general principle was that interference-induced spectral variation can be regarded as the linear combination of some underlying interference sources and can be split from the common part over all spectral measurement cases. Then, the regression model created in the interference-filtered common subspace could be readily transferred over different measurement cases for quality prediction. Although many investigators have reported that ICs extracted from the spectra of mixtures are able to capture the essential structure and are chemically meaningful, some underlying problems have not yet been brought into focus. Commonly, evolving along the wavelength direction, the spectra of mixtures show significant fluctuations and dynamics that actually are alternately dominated by the spectral profiles of different constituent species. This reveals different underlying chemical characteristics and different influential relationships on quality properties. Here, it is called “multiplicity” of spectral characteristics. Without losing generality, over different subspaces, different underlying chemical characteristics and different effects on the final products are suggested. Located in different spectral regions, the influences on ICA feature extraction and regression performance are also different. For example, certain source species can be accurately extracted in some spectral subspaces
8695
by ICA, whereas in other spectral subspaces, they cannot be estimated reliably enough or even cannot be determined by ICA because they do not dominate the current mixture spectra. Clearly, a unified ICA feature extraction model over the fullrange spectral space is not sufficient and desirable enough to recover the real constituent features and track the varying chemical characteristics along the wavelength direction. A larger bias could be found when different underlying spectral information is mixed together. The resulting mixing matrix, which is directly related to the estimation of the constituent concentrations, can result in high prediction error when the recovered ICs poorly match the real pure substances existing in the chemical mixture. Considering the multiplicity of spectral characteristics, the idea of spectral subspace separation was first proposed by Zhao et al.38 Over different subspaces, different ICA decomposition relationships were determined, revealing different underlying chemical characteristics along the wavelength direction. The multiblock ICR model was thus constructed and was demonstrated to be able to give superior quality prediction results compared with the unified ICR model. However, the multiblock ICR model focused only on quality prediction and did not consider the interpretation of multivariate calibration models and the statistical analysis of the betweenset (X-Y) relationship, which, however, can be as important as predictive ability. In the present work, a multiblock version of the joint LV (MBJLV) modeling strategy is developed to establish qualitative and quantitative analysis for spectral calibration. The proposed algorithm highlights the facts that, along the wavelength direction, each spectral subspace can explain only one part of the quality property and all subspaces work together to describe the overall quality variations. This phenomenon is called “subspace-specific local and global effects” here. Once the argument is recognized, it is not advisible to strive for excessive prediction precision localized in a specific separate subspace. Therefore, our emphasis is on how to find the different local contribution powers of various subspaces and their collective efforts to express the final qualities. To achieve this purpose, our analysis is implemented in two steps. First, considering the specific effect of ICA on spectroscopy, it is used to identify multiple spectral subspaces from the original wavelength space. In different subspaces, which cover different spectral wavelengths and enclose different underlying chemical characteristics, both constituent spectra and their mixing relationships are determined, revealing different dominant source spectral profiles. In the second step, the mulitblock joint LV (MBJLV) modeling strategy is designed on the basis of subspace separation. In both the descriptor and quality spaces, the underlying systematic variations are decomposed into different parts under the supervision of each other. This yields superior interpretability because it is possible to “zoom in” on separate blocks and evaluate their specific relationship with Y in present of the other Xb blocks. This provides (a) estimates of the pure-constituent profiles over different spectral wavelength regions and also models (b) the closely related Xb-Y variations, (c) the Y-orthogonal variation in each Xb block, and (d) the Xb-orthogonal variation in Y. The covarying part portrays the prediction-useful descriptor and quality variations within each specific block. Particularly, the decomposition of orthogonal variations in specific model structures evaluates the information that does not contribute to the mutual prediction (Xb T Y), but might still be of importance for interpretation of the between-set relationship. This article is organized as follows: In the next section, the multiblock quantitative and qualitative analysis strategy is
8696
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
described in detail based on the separation of the spectral subspace. The underlying principle is explained, and knowledge of some related modeling methods is also introduced. The application to a real spectral case demonstrates the effectiveness of the proposed method. The results are discussed, with an emphasis on the suitability of the proposed method and its advantages for spectral calibration analysis and understanding. Finally, conclusions are drawn in the last section. 2. Proposed MBJLV Calibration Analysis Strategy As analyzed before, for spectral data, the wavelength variables are many in number and vary over different wavelength regions, and the underlying chemical characteristics are different for different dominating source spectra. A unified ICA model focusing on the full range of wavelengths cannot probe the varying underlying information comprehensively enough. A larger bias could be found when different underlying spectral information is mixed together over wavelengths. The resulting mixing matrix, which is directly related to the constituent concentrations, can result in high prediction error and worse calibration analysis performance when the recovered ICs provide poor matches for the real pure substances existing in the chemical mixture. Therefore, one spectral data space should be further divided into several different subspaces if the inherent spectral information changes. Without losing generality, over different subspaces, different underlying chemical characteristics and different relationships with the qualities are suggested. The important and meaningful issue is how to check their details and quantify them in specific model parameters. In the proposed MBJLV modeling strategy, blocking of spectral variables and partitioning of underlying variations can be combined to check the spectral space more comprehensively, reveal more spectral information, and thus obtain more meaningful statistical analysis results. The modeling details are described in the following subsections. 2.1. Spectral Subspace Separation. Considering the multiplicity of underlying spectral characteristics, it is necessary to divide a spectral wavelength region into several meaningful subspaces and analyze both their different underlying chemical information and their different effects on the final qualities over multiple spectral regions. Different subspaces can be identified by evaluating the resulting IC extraction performance and model representability. Here, a complete automatic spectral subspace separation algorithm is designed as described in Appendix A. Its goal is to find the segments along the wavelength direction that can be well approximated by one ICA model with sufficient representability. It performs a greedy search based on ICA model performance evaluation. A division performance index, the mean squared error (MSE), is set up, which is related to the unexplained squared error of each wavelength region after ICA model fitting. It is calculated as MSE )
1 NJ
N
J
∑ ∑ (x
i,j
- xˆi,j)2
i)1 j)1
where xi,j is the measured spectral variable (j ) 1, 2, ..., J) in each sample (i ) 1, 2, ..., N) and xˆi,j is its reconstructed value by ICA. During the separation process, every possible wavelength interval will be tried, and the best spectral separation points will be determined in terms of reconstruction power. Any possible spectral subdivision is evaluated in terms of the percentage reduction of MSE, which is used to check whether the subdivision can improve the representation performance of the ICA model. The subdivision with the greatest reduction of
Figure 1. Spectral subspace identification scheme.
the MSE is chosen, and the procedure proceeds recursively with the resulting segments (corresponding submatrices of spectral data), as shown in Figure 1. Upon application of the designed spectral subspace separation strategy (described in Appendix A), multiple specific regions are obtained, enclosing different spectral wavelengths. This approach provides different analysis scenes for ICA-based source signal extraction and improves the ICA model representability. The advantage of the subspace separation is that, in addition to a comprehensive view for the whole spectral region, one also obtains local analysis viewpoints for different spectral wavelength regions. It is easy to determine how they will act under the influence of each other. It is assumed that both their respective underlying characteristics and their contributions to quality properties tend to be hidden by each other to a certain extent when the whole range of spectral wavelengths is not divided. Focusing on separate spectral subspaces provides improved interpretation and potential for capturing the key factors in different subspaces for quality interpretation and improvement. In the separation procedure, the retained ICs are mainly used to aid in evaluating whether the separation could improve the ICA model representability with the same model order. Therefore, for the entire wavelength region, the number of retained ICs is determined by cross-validation so that a small reconstruction error (MSE) is obtained. This approach is implemented and evaluated by considering both training and testing data sets. Then, the same number of ICs is retained for modeling the resulting subspaces during the whole subspace separation procedure. 2.2. Multiblock-Based Calibration Analysis. Over different subspaces, different IC decompositions lead to different effects on qualities, providing different analysis platforms for quality description. One analysis approach is to treat each predictor block independently. However, this approach overlooks the correlations among different blocks and does not explore their joint effects. The multiblock40-44 modeling idea is expected to be able to treat them wisely, which can provide superior interpretability because it is possible to obtain local information (block scores) and global information (super scores) simultaneously from the data. Here, based on the separation of different spectral regions (i.e., blocking of spectral wavelength variables), the multiblock
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
calibration analysis strategy is formulated so that information from multiple subspaces can be combined sensibly, in addition to providing detailed comprehension of local spectral regions. Generally, it is assumed that product qualities depend not only on the local level of each spectral subspace but also on the overall effects of all subspaces within the entire wavelength region, in which different wavelength ranges might contribute to different parts of quality variations to different extents. Therefore, calibration analysis is interested in determining which parts of the descriptor variations are explaining which parts of the quality variations within each wavelength subspace. Another important issue is whether unique information useful for quality description exists in the separate blocks, that is, whether one subspace is sufficient or several or even all blocks are needed. If multiple subspaces are needed, calibration analysis determines how much of the quality information can be well explained by each block. These measures can be advantageous when several blocks are available to see. To achieve the above analysis purposes, over different spectral subspaces, both X and Y spaces are decomposed into two different systematic parts under the control of different between-set relationships, covering both closely related and irrelevant ones. First, the multiblock partial least-squares canonical correlation analysis (MBPLS-CCA) method is designed as described in Appendix B to extract the closely related part, whereas a postprocessing CCA46,47 is implemented on the PLS LVs and thus directly maximizes the X-Y correlations. Then, the second part is the mutual irrelevant systematic information obtained by orthogonal signal correlation (OSC) in each subspace, which is also of analytical interest to the comprehension of the X-Y relationship. The above two parts of systematic variations explain different underlying pieces of information in each subspace that function differently in relating X and Y. Let {X1,X2, ..., Xb, ..., XB} represent B (Jb × N) matrices of spectral measurements with the same N samples scanned over different spectral regions (where the subscript b is used to represent one specific block and Jb is the number of wavelength variables in each subspace), which all point to the same concentration matrix (Y) of constituents in a mixture. From a mathematical point of view, by performing ICA decomposition,39 their underlying source spectral profiles, located in different spectral subspaces, can be determined as ˆ 1 ) S1A1 X ˆ 2 ) S2A2 X l ˆ b ) SbAb X l ˆ B ) SBAB X
(1)
where different source spectra [Sb(Jb × Rb)] and multiple mixing relationships [Ab(Rb × N)] are decomposed (Rb is the number of retained ICs, which can be different for different spectral subspaces). Ab reveals the different contributions of the source spectral profiles to the spectra of mixtures in different subspaces/ blocks. Her,e it should be noted that the ICA model order (Rb), which is chosen by similar guidelines with respect to the smallest reconstruction error, might be different from that used in spectral subspace separation. All of the mixing matrices obtained from different spectral subspaces can thus be organized as multiple blocks, A ) [A1T A2T ... AbT ... ABT], that reveal different underlying information for quality description and meanwhile can correlate with each other over subspaces. The prepared ICA mixing matrices (A )
T
T
T
8697
T
[A1 A2 ... Ab ... AB ]) and the concentration measurements (Y) are input into the multiblock-based modeling algorithm described in Appendix B. Typically, A and Y are mean-centered and appropriately scaled, thus giving each equal weight. Over different subspaces, both closely related and irrelevant systematic parts are decomposed and distinguished in each descriptor block and the quality space under the supervision of each other, where the related statistics are calculated as Model of Ab
Tbc ) AbWbc Eb ) Ab - Abc ) Ab - Tbc(TbcTTbc)-1TbcTAb ) Ab - TbcPbcT Tbo ) AbWbo Ebf ) Eb - Abo ) Eb - Tbo(TboTTbo)-1TboTEb ) Eb - TboPboT
(2) Model of Y Ubc ) YCbc Fb ) Y - Ybc ) Y - Ubc(UbcTUbc)-1UbcTY ) Y - UbcQbcT Ubo ) YCbo Fbf ) Fb - Ybo ) Fb - Ubo(UboTUbo)-1UboTFb ) Fb - UboQboT
(3) where Wbc and Wbo are PLS-CCA weights and OSC weights, respectively, in Ab space, with Pbc and Pbo as the corresponding loadings, and Cbc and Cbo are PLS-CCA weights and OSC weights, respectively, in Y space, with the corresponding loadings Qbc and Qbc. Ebf and Fbf are the final residuals. Both PLS-CCA and OSC components can be calculated directly from the original data space (Ab and Y). Moreover, the correlated and orthogonal variations are orthogonal to each other. In summary, the MBJLV modeling result is formulated as Model of Ab Ab ) Abc + Abo + Ebf ) TbcPbcT + TboPboTP + Ebf
(4) Model of Y Y ) Ybc + Ybo + Fbf ) UbcQbcT + UboQboT + Fbf
(5) where the first parts (TbcPbcT and UbcQbcT) reveal the systematic pieces of information that are really related to each other and can be used to predict each other, revealing the between-set similarity. The second part reveals the orthogonal variations, TboPboT and UboQboT, or the systematic pieces of information that cannot be interpreted by each other; they are used to evaluate the between-set dissimilarity. Therefore, both betweenset similarity and between-set dissimilarity are represented and quantified in specific model parameters that play different roles in calibration analysis to evaluate the between-set (X-Y) relationship. Moreover, it is worth noting that Y is decomposed differently over different subspaces, because of the different effects of the subspaces on quality interpretation. From the above decomposition results, various model statistics can then be calculated to quantitively evaluate the different variations underlying both different spectral subspaces and quality space for improved model interpretation and spectral understanding. This aspect is further clarified in the next section. 3. Simulations and Discussions In this section, the performance of calibration modeling and quality prediction using the proposed method is illustrated through two case studies. Quantitative and qualitative analyses
8698
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
are performed to reveal the benefits resulting from subspace separation with respect to both quality prediction performance and statistical interpretations. 3.1. Experimental Data Set. Case Study 1. The data set used for case study 1 consisted of spectra from 80 samples of corn with wavelengths in the range 1100-2498 nm at 2-nm intervals (700 channels) scanned on an m5 NIR spectrometer, composing spectral measurement data X(700 × 80). The corresponding concentration values were provided in the response matrix Y(80 × 4), referring to the four constituents moisture, oil, protein, and starch. Fifty samples were used for model training, and the other 30 were used as testing data. The corn spectral data are available at the eigenvector Research homepage: http://www.eigenvector.com/data/Corn/. Case Study 2. The spectral data for case study 2, obtained from the literature,48 were recorded for ternary mixtures of water, ethanol, and 2-propanol in a 1-cm cuvette in the wavelength range of 580-1091 nm. The short-wave NIR spectra of 22 samples were obtained at 70 °C, and the temperature of the samples was controlled (∼0.2 °C variation). In addition to the spectral measurement matrix X(512 × 22), the quality matrix Y(22 × 3) describes different contents of ethanol, water, and isopropanol. Fifteen samples were selected for model training, and the other seven were used as testing data. 3.2. Simulation Methodology and Results. First, for both case studies, mixture spectra are shown in Figure 2, taking examples from the training samples. Along the wavelength direction, the spectral profiles show great fluctuations. As analyzed before, different basic chemical components dominate over different wavelength regions. By performing the proposed subspace separation, different spectral blocks are obtained, covering different wavelength ranges. As shown in Figure 2, for case study 1, five spectral subspaces are obtained, with the ranges 1100-1364, 1366-1468, 1470-1894, 1896-2392, and 2394-2498 nm; for case study 2, three spectral subspaces are obtained, with the ranges 1-149, 150-302, and 303-512 nm. Over different subspaces, the source spectra and their mixing coefficients are extracted. For both cases, by cross-validation, the proper number of ICs is retained in each subspace to recover the important source spectra, as reported in Table 1. Moreover, the reconstructed spectral variations were also calculated, indicating the model fitting competency in each subspace. Generally, for both cases, most of the underlying spectral variations can be modeled by using only a few ICs. Taking as an example the first two source components, the spectral profiles are illustrated in Figure 3 over different subspaces. For comparison, single ICA (SICA) decomposition results over the entire wavelength region are also displayed for both cases in Figure 4. It is clear that the recovered source component profiles based on spectral subspace separation are quite different from those decomposed by SICA. These spectral subspaces can correlate with and influence each other. Multiple mixing matrices are then collected and arranged into a predictor matrix, A ) [A1T A2T A3T A4T A5T] in case study 1 and A ) [A1T A2T A3T] in case study 2, and then related to the real concentration measurements (Y). First, by the standard ICR (SICR), multiblock PLS (MBPLS) modeling, and proposed MBJLV algorithms, quality predictions are performed for both case studies, as reported in Table 2. Here, the SICR method is designed as follows: First, a unified ICA decomposition is performed to decompose the entire wavelength range, and then single PLSCCA is used to fit the regression relationship. In this way, the comparison will focus on revealing the effects of spectral subspace separation. For the MBPLS and MBJLV methods, the
Figure 2. Mixture spectra trajectories and the subspace separation results for (a) corn data and (b) ternary mixture. Table 1. ICA Decomposition Results over Different Spectral Subspaces for (a) Corn Data and (b) Ternary Mixture (a) Corn Data spectral subspaces statistics
I
II
III
IV
V
model order ˆ b (%) ICA reconstruction R2X
6 88.58
5 98.80
6 99.52
5 99.86
4 99.91
(b) Ternary Mixture spectral subspaces statistics
I
II
III
model order ˆ b (%) ICA reconstruction R2X
4 99.99
5 99.98
5 97.16
comparison can focus on the modeling technique, that is, the conventional PLS approach and the current joint LV modeling approach. Two indices, the mean square error of calibration (MSEC) and the mean square error of prediction (MSEP), which are calculated based on training samples and testing samples, respectively, are used to evaluate the prediction performance. Generally, the results indicate that the method using the separation of the spectral subspace can better decompose the underlying source spectra so that the resulting ICA mixing matrices are more quality-informative, revealing better fitness and generalization adaptability for quality estimation. Com-
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
8699
Figure 3. Spectral profiles of the first two ICs over different spectral subspaces for (a) corn data and (b) ternary mixture.
paratively, with no subspace separation, the direct application of single ICA decomposition to the full-wavelength-range spectra might damage the prediction accuracy. This demonstrates that, because of the varying underlying spectral characteristics, a unified ICA decomposition model cannot recover source components well enough over the full range of wavelengths. Therefore, the undesirable extraction result, as the descriptors, will directly influence the subsequent regression modeling performance. Such a problem might be more apparent when the differences in the underlying chemical characteristics over different spectral subspaces are larger. This result can be understood from the fact that the more different the spectral subspaces, the larger their spectral variations, making the fitting of a unified ICA model to the full wavelength range at the same time less representative and more difficult. Generally, it is expected that spectral variations can be better fitted by dividing the original spectral space into different subspaces according to the changes in the underlying spectral characteristics and then modeling the different subspaces using multiple specific ICA models. The separated subspaces provide different statistical analysis platforms, which reveal different impacts of different substances on the spectra of mixtures. In particular, if the spectral library of pure substances is known, one simply needs to compare the uncovered ICs with the available spectra, and then it is easy to identify the unknown constituent species existing in the mixture and determine which chemical components dominate in each wavelength subspace. Moreover, compared with the proposed method, the MBPLS algorithm shows comparable quality prediction performance in both case studies, which agrees well with the report in previous works31,32 that the benefit of ICR lies in the improvement of chemical
explanations, rather than improvement of quality prediction performance) as the ICs were found to be strongly correlated to the NIR spectra of the source components in the spectra. In the following discussion, using the proposed method, quantitative and qualitative calibration analysis can be implemented to evaluate different roles of each spectral subspace. A summary of the modeling and analysis results is given in Table 3. In each specific spectral subspace, the systematic variations are partitioned into two parts, the X-Y covarying part and the X-Y orthogonal (or, more precisely, unrelated) part. Although only the X-Y covariations participate in the interpretation of the counterpart, the systematic orthogonal variations represented by the OSC structure are still of interest for the evaluation of between-set relationships. Generally, the descriptor variations are more complex, and therefore, more CCA and OSC components are needed in X space than in Y space. Four model statistical criteria are taken into account here in both X and Y spaces; they are calculated based on the modeling results shown in eqs 4 and 5: modeled correlation variations as indicated by Abc and Ybc, namely, R2Abc and R2Ybc, and modeled irrelevant variations as indicated by Abo and Ybo, namely, R2Abo and R2Ybo. Moreover, variations predicted by Tbc and Ubc are both ˆ b. Here, it should be noted that R2Y ˆb ˆ b and R2A calculated, R2Y indicates the prediction competency of each subspace, which is different from R2Ybc because R2Ybc denotes the variations modeled by the correlated part (Ybc) separated in Y space. These model statistics quantify the local effects of different subspaces with their correlations taken into consideration, for example, what roles various wavelength regions play concerning quality interpretation and prediction, how many spectral variations provide systematic information and how many take part in
8700
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010 Table 2. Quality Prediction Comparisons for (a) Corn Data and (b) Ternary Mixture method
MSEC
MSEP
(a) Corn Data SICR algorithm
0.2851 0.3739 0.2108 0.3711 0.0749 0.1442 0.0369 0.1038 0.0385 0.1004 0.0263 0.1279
MBPLS algorithm
MBJLV algorithm
1.0801 0.8699 0.5587 0.9413 0.3147 0.6241 0.3177 0.7130 0.3027 0.8086 0.2966 0.3318
(b) Ternary Mixture SICR algorithm
0.0664 0.0036 0.0457 0.0023 0.0006 0.0012 0.0049 0.0013 0.0020
MBPLS algorithm MBJLV algorithm
0.4557 0.0732 0.2673 0.0700 0.0343 0.1575 0.0863 0.0373 0.2380
Table 3. Summary of Quantitative Calibration Analysis Results for (a) Corn Data and (b) Ternary Mixture (a) Corn Data spectral subspaces modeling variations X model
Y model
number of CCA components OSC 2 R Abc (%) R2Abo (%) ˆ b (%) R 2Y number of CCA components OSC R2Ybc (%) R2Ybo (%) ˆ b (%) R 2A
I
II
III
IV
V
1 1 71.78 0.56 34.87 1 1 47.09 40.78 62.71
2 2 75.59 0.32 14.50 2 1 86.58 0.01 59.73
1 2 47.79 1.79 31.12 1 1 46.80 40.78 41.77
3 2 91.07 2.51 4.62 1 1 37.01 41.48 13.16
2 1 84.42 3.10 17.13 2 1 87.66 0.01 45.64
(b) Ternary Mixture spectral subspaces
Figure 4. Spectral profiles of the first two ICs using the single ICA method for (a) corn data and (b) ternary mixture.
quality prediction, how many quality variations can be explained by the current subspace, and which subspace is more important for quality interpretation and prediction. Based on the quantitative evaluation of these model statistics, we can further realize the underlying meaning of subspace separation and variation partition. Generally, over different subspaces, different types of variations dominate. For example, for case study 1, in X space, as indicated by R2Abc, the correlation variations take a great amount compared with the Y orthogonal variations as indicated by R2Abo. Over different subspaces, there are different predictive abilities with respect to quality variations. For ˆ b along the wavelength example, in case study 1, comparing R2Y direction, the first subspace is assumed to be relatively most quality-predictive, whereas the fourth subspace is most uncritical for quality prediction. Moreover, in case study 2, the three subspaces generally show similar quality prediction performances, whereas in case study 1, the contributions to quality prediction fluctuate greatly over subspaces. It should be noted that larger correlated variations (R2Abc and R2Ybc) do not indicate that they can predict greater amounts of variation in their ˆb counterparts, which can be clearly seen in comparison with R2A 2ˆ and R Yb for the two case studies. Generally, predicted variations
modeling variations X model
Y model
number of components R2Abc (%) R2Abo (%) ˆ b (%) R2 Y number of components R2Ybc (%) R2Ybo (%) ˆ b (%) R2 A
CCA OSC
CCA OSC
I
II
III
1 2 28.99 51.36 52.91 1 1 54.03 45.97 28.15
2 3 64.39 35.00 45.28 1 1 45.74 52.94 25.04
2 2 81.53 16.89 42.07 1 1 53.11 4.13 68.77
in each subspace are less than correlated variations in Y space ˆ b over multiple subspaces, we obtain (R2Ybc). Summing all R2Y values much greater than 100%, which demonstrates that there are overlapping parts with respect to prediction ability over various subspaces. In addition, aside from the correlated and orthogonal variations, as indicated by their sums in X and Y spaces, respectively, R2Abc + R2Abo < 100 and R2Ybc + R2Ybo < 100, there are still some variations that are not modeled and are thus left in the residuals. They might also cover some systematic information to some extent, in addition to measurement errors, but they are of no use for between-set relationship analysis. As mentioned before, each subspace might quantify only one fraction of quality variations, in which only those quality-
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
8701
Figure 5. Subspacewise regression coefficients with respect to different quality indices for (a) corn data and (b) ternary mixture.
responsible systematic variations are used to model quality variations. Moreover, the overall spectral properties should be comprehensively stacked and characterized for quality interpretation and prediction, revealing the subspacewise cumulative effects along the wavelength direction. To illustrate this, the final quality prediction is obtained by considering all subspaces ˆ values are 79.86% for case study 1 simultaneously. The R2Y ˆb and 99.71% for case study 2, which are both larger than R2Y over different subspaces. This well demonstrates the cumulative effects of multiple subspaces on quality prediction. Moreover, over different subspaces, as shown in Table 1, there are different descriptor variables in X space, and correspondingly, there will be 26 regressors in case study 1 and 14 regressors in case study 2. The corresponding regression coefficients with respect to all quality indices (four quality variables in case study 1 and three variables in case study 2) are shown in Figure 5, revealing different weight distributions over different regressors and different subspaces. Although it is not necessary that the size of the regression coefficient be consistent with the significance of the associated regressor, these coefficients might still reveal some information more or less from another viewpoint. For example, they might indicate that, in different subspaces, the extracted ICs, as regressors, might be different and also have different effects on qualities, negative and positive. From the plot shown in Figure 5a, it is interesting to find that, for quality
variables 1 and 2, the magnitudes of the regression coefficients are generally similar and the only difference is their influencing direction, positive or negative. A similar circumstance can be also found in the situation for quality variables 3 and 4 in case study 1. In conclusion, the above case studies have illustrated the theoretical analyses mentioned in section 2 concerning different chemical characteristics over wavelength regions and their different roles in quality prediction and interpretation. From the results, spectral subspace separation might have two beneficial aspects: First, it improves the quality prediction performance compared with that obtained by SICR. Second, the separated multiple subspaces provide a more detailed statistical analysis platform, from which one can expect to further check the local specific characteristics of different subspaces and obtain a more comprehensive understanding before regression modeling. Based on the variation separation, a quantitative statistical analysis performed in both regressor and response spaces over different subspaces has revealed their different roles in quality interpretation. The results of this study constitute a step forward toward spectral data analysis and calibration modeling for practical applications, revealing desirable improvement in model representation and quality interpretation. Also, this approach provides the basis for further work and improvement. Future investigations might profitably cover the following issues:
8702
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
(1) With respect to subspace separation, how can the separation algorithm be further improved to obtain more proper model representation? Instead of the searching algorithm described in subsection 2.1, there are many other optimization approaches, such as evolutionary programming, that can be used. It would be meaningful to check whether the subspace separation results can be further improved by alternative methods. (2) The subspace separation results provide a potential statistical analysis platform. The question of whether, instead of the current multiblock analysis strategy, these spectral subspaces can be made better for use in model design and quality interpretation is also a meaningful issue. (3) It has to be pointed out that ICA is not the only choice for the recovery of pure spectra of constituent species and their concentrations from the spectra of chemical mixtures. Other approaches, such as EFA as mentioned in the Introduction, can also be used as an alternative, which, however, is not addressed here. Further, considering that nonlinearity is not uncommon in spectral data, possibly a nonlinear version is required so as to better decompose the underlying nonlinear structure. 4. Conclusions In this article, a qualitative and quantitative calibration analysis method for spectral data is presented. Different spectral subspaces are separated from the original wavelength space based on changes in the underlying chemical characteristics. This provides a multiblock calibration analysis platform from which both local information and global information can be readily obtained. In each subspace, as indicated by their different roles in quality interpretation and prediction, the underlying systematic variations are decomposed into different parts and quantitatively evaluated. In this way, this method can provide more interesting and reliable model representation, better identify the underlying spectral information, and obtain enhanced statistical interpretation. Simulation examples demonstrate the effectiveness of the proposed method. It should be generally applicable to a broad range of spectral applications, aiding in more comprehensive calibration analysis and understanding. Acknowledgment This work was supported by China’s National 973 program (2009CB320603). Appendix A Spectral Subspace Separation Algorithm The input of the division algorithm involves: (a) the whole spectral wavelength space, (b) the minimum length of the separated wavelength region (minlen), (c) the improvement threshold of MSE to accept a subdivision, and (d) the maximum number of separated spectral subspaces (maxnum).
The recursive division procedure is implemented as follows: (1) Set the initial number of separated spectral regions as numreg ) 1. (2) Input the prepared spectral data set, perform the modified ICA algorithm39 on the data, and obtain the initial ICA model. Calculate the resulting MSE0 index value to evaluate the model representation performance.
(3) For each wavelength interval (j) within the current spectral data set, if the subdivision in j generates two segments with lengths higher than the predefined minimum region length (minlen), perform ICA on the two segments individually; calculate the resulting new MSE values for both segments, denoted as MSE1j_1 and MSE1j_2. (4) Find the sampling wavelength interval (j•) at which the average of MSE1j_1 and MSE1j_2 is lowest. Comparing both MSE1j_1 and MSE1j_2 with MSE0, if the improvement for either segment, MSE0-MSE1j_1 and MSE0-MSE1j_2, does not reach the predefined improvement threshold, then stop this branch. (5) Otherwise, accept the subdivision and update the number of separated wavelength regions numreg: ) numreg + 1. If the current numreg is no longer less than the predefined maxnum, stop the iteration procedure; output the division result. (6) Otherwise, recursively repeat steps 2-5 for either of the two resulting segments, each now employed as the new spectral input data space in step 2. The output is the spectral subspace separation result.
Appendix B
MBJLV Modeling Procedure Input Data. X ) [X1 X2 ... XB] and Y Step 1. Choose a starting matrix u and iterate the following equations until convergence of tT: (1.1) Perform the regular PLS-CCA algorithm45 between Xb and u with block weights wbc ) Rb,pls · wb,cca (where Rb,pls is a Jb × Ab,pls-dimensional PLS weight matrix, wb,cca is an Ab,plsdimensional CCA weight vector, Jb is the number of measurement variables in block data Xb, and Ab,pls is the number of retained block PLS LVs) and block scores tbc ) Xbwbc and T ) [t1c t2c ... tBc]. (1.2) Perform the regular PLS-CCA45 between T and u, with super weights wTc ) RT,pls · wT,cca (where RT,pls is a B × AT,plsdimensional PLS weights matrix, wT,cca is an AT,pls-dimensional weight vector, and AT,pls is the number of retained super PLS LVs) and super scores tTc ) T · wTc. (1.3) Obtain the new Y super score u
q)
YTtTc tTcTtTc
where u ) Yq, and then scale u to be of length 1. Step 2. Deflate the residuals
pbc )
XbTtbc tbcTtbc
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
Eb ) Xb - tbcpbc
T
E ) [E1 E2 ... EB] F ) Y - tTcqT From here, one can go to step 1 to implement the above procedures for the new LVs, where Xb and Y are both replaced by their corresponding residual matrices Eb and F. In this way, block canonical LVs, Tbc, are prepared. Step 3. In each block, perform the regular PLS-CCA algorithm45 on Y and Tbc, and obtain Y canonical scores, Ubc, Y weights Cbc and loadings Qbc, as well as the residual Fb ) Y - UbcQbcT. Step 4. In each block, perform Fearn’s OSC algorithm29 on Xb and Ubc to extract the Ubc-orthogonal components in Xb, Tbo, as well as the weights Wbo. Then, calculate the loadings PboT ) (TboTTbo)-1TboTEb, so as to model the Ubc-orthogonal variations from the PLS-CCA residual Eb. Also, perform Fearn’s OSC algorithm29 on Y and Tbc to extract the Tbc-orthogonal components in Y, Ubo, and the weights Cbo. Then, calculate the loadings QboT ) (UboTUbo)-1UboTFb so as to model the Tbc-orthogonal variations from the PLS-CCA residual Fb. Output Results. From PLS-CCA, one obtains block LVs, Tbc and Ubc; block weights,Wbc and Cbc; and block loadings, Pbc and Qbc. From OSC, one obtains block LVs, Tbo and Ubo; block weights,Wbo and Cbo; and block loadings, Pbo and Qbo. Here, some important points should be noted: (a) During the algorithm, to calculate each PLS-CCA block score and super score, multiple PLS scores are prepared as shown in step 1. This results from the consideration that PLS LVs prepare possible quality-related descriptor variations for the following CCA postprocessing and then only those closely quality-related systematic information will be extracted from them by CCA. (b) Because the super LVs (T) summarize the information contained in all blocks, deflation of Xb using T can make the information of the separate blocks get mixed up. This leads to block LVs that describe not only the information from their specific block but also information from other blocks, resulting in interpretation problems. Therefore, here, we uses block LVs to deflate Xb, which makes them orthogonal and prepares blockspecific local information useful for quality description. Moreover, these block LVs are extracted by considering the blockwise correlations, which give information for the specific block in relation to Y in the presence of the other blocks and make interpretation much easier. Literature Cited ¨ hman, J. Orthogonal signal (1) Wold, S.; Antti, H.; Lindgren, F.; O correction of near-infrared spectra. Chemom. Intell. Lab. Syst. 1998, 44, 175. (2) Bijlsma, S.; Louwerse, D. J.; Smilde, A. K. Rapid estimation of rate constants of batch processes using on-line SW-NIR. AIChE J. 1998, 44, 2713. (3) Westerhuis, J. A.; Gurden, S. P.; Smilde, A. K. Spectroscopic monitoring of batch reactions for on-line fault detection and diagnosis. Anal. Chem. 2000, 72, 5322. (4) Gurden, S. P.; Westerhuis, J. A.; Smilde, A. K. Monitoring of batch processes using spectroscopy. AIChE J. 2002, 48, 2283. (5) Othman, N. S.; Fevotte, G.; Peycelon, D.; Egraz, J. B.; Suau, J. M. Control of polymer molecular weight using near infrared spectroscopy. AIChE J. 2004, 50, 654.
8703
(6) Gabrielsson, J.; Jonsson, H.; Trygg, J.; Airiau, C.; Schmidt, B.; Escott, R. Combining process and spectroscopic data to improve batch modeling. AIChE J. 2006, 52, 3164. (7) Reis, M. M.; Arau´jo, P. H. H.; Sayer, C.; Giudici, R. Spectroscopic on-line monitoring of reactions in dispersed medium: Chemometric challenges. Anal. Chim. Acta 2007, 595, 257. (8) Ye, S. F.; Wang, D.; Min, S. G. Successive projections algorithm combined with uninformative variable elimination for spectra variable selection. Chemom. Intell. Lab. Syst. 2008, 91, 194. (9) Xu, H.; Liu, Z.; Cai, W.; Shao, X. A wavelength selection method based on randomization test for near-infrared spectral analysis. Chemom. Intell. Lab. Syst. 2009, 97, 189–193. (10) Zhao, C.; Gao, F.; Wang, F. Phase-based joint modeling and spectroscopy analysis for batch processes monitoring. Ind. Eng. Chem. Res. 2010, 49, 669–681. (11) Sjo¨blom, J.; Svensson, O.; Josefson, M.; Kullberg, H.; Wold, S. An evaluation of orthogonal signal correction applied to calibration transfer of near infrared spectra. Chemom. Intell. Lab. Syst. 1998, 44, 229. (12) Gusnanto, A.; Pawitan, Y.; Huang, J.; Lane, B. Variable selection in random calibration of near-infrared instruments: Ridge regression and partial least squares regression settings. J. Chemom. 2003, 17, 174. (13) Trygg, J. Prediction and spectral profile estimation in multivariate calibration. J. Chemom. 2004, 18, 166. (14) Andrew, A.; Fearn, T. Transfer by orthogonal projection: Making near-infrared calibrations robust to between-instrument variation. Chemom. Intell. Lab. Syst. 2004, 72, 51. (15) Chen, T.; Morris, J.; Martin, E. Gaussian process regression for multivariate spectroscopic calibration. Chemom. Intell. Lab. Syst. 2007, 87, 59–67. (16) Preys, S.; Roger, J. M.; Boulet, J. C. Robust calibration using orthogonal projection and experimental design. Application to the correction of the light scattering effect on turbid NIR spectra. Chemom. Intell. Lab. Syst. 2008, 91, 28. (17) Benoudjit, N.; Melgani, F.; Bouzgou, H. Multiple regression systems for spectrophotometric data analysis. Chemom. Intell. Lab. Syst. 2009, 95, 144. (18) Alciaturi, C. E.; Quevedo, G. J. Bayesian regularization: Application to calibration in NIR spectroscopy. J. Chemom. 2009, 23, 562. (19) Chen, T.; Martin, E. Bayesian linear regression and variable selection for spectroscopic calibration. Anal. Chim. Acta 2009, 631, 13. (20) Geladi, P.; Kowalski, B. R. Partial least-squares regressionsA tutorial. Anal. Chim. Acta 1986, 185, 1. (21) Brereton, R. G. Introduction to multivariate calibration in analytical chemistry. Analyst. 2000, 125, 2125. (22) Kleinbaum, D. G.; Kleinbaum, D. G. Applied Regression Analysis and Other MultiVariable Methods, 4th ed; Thomson Brooks/Cole: Belmont, CA, 2008. (23) Kutner, M. H.; Nachtsheim, C.; Neter, J. Applied Linear Regression Models, 4th ed.; McGraw-Hill/Irwin: Boston, 2004. (24) Ergon, R. Reduced PCR/PLSR models by subspace projections. Chemom. Intell. Lab. Syst. 2006, 81, 68. (25) Lindgren, F.; Geladi, P.; Rannar, S.; Wold, S. Interactive variable selection (IVS) for PLS. Part I: Theory and Algorithm. J. Chemom. 1994, 8, 349. (26) Forina, M.; Casolino, C.; Pizarro Millan, C. Iterative predictor weighting (IPW) PLS: A technique for the elimination of useless predictors in regression problems. J. Chemom. 1999, 13, 165. ¨ hman, J. Orthogonal signal (27) Wold, S.; Antti, H.; Lindgren, F.; O correction of near-infrared spectra. Chemom. Intell. Lab. Syst. 1998, 44, 175. (28) Sjo¨blom, J.; Svensson, O.; Josefson, M.; Kullberg, H.; Wold, S. An evaluation of orthogonal signal correction applied to calibration transfer of near infrared spectra. Chemom. Intell. Lab. Syst. 1998, 44, 229. (29) Fearn, T. On orthogonal signal correction. Chemom. Intell. Lab. Syst. 2000, 50, 47. (30) Westerhuis, J. A.; De Jong, S.; Smilde, A. K. Direct orthogonal signal correction. Chemom. Intell. Lab. Syst. 2001, 56, 13. (31) Chen, J.; Wang, X. Z. A new approach to near-infrared spectra data analysis using independent component analysis. J. Chem. Inf. Comput. Sci. 2001, 41, 992. (32) Shao, X. G.; Wang, W.; Hou, Z. Y.; Cai, W. S. A new regression method based on independent component analysis. Talanta. 2006, 69, 676. (33) Navea, S.; Tauler, R.; Juan, A. D. Monitoring and modeling of protein processes using mass spectrometry, circular dichroism, and multivariate curve resolution methods. Anal. Chem. 2006, 78, 4768. (34) Zhao, C. H.; Gao, F. R.; Yao, Y.; Wang, F. A robust calibration modeling strategy for analysis of interference-subject spectra data. AIChE J. 2010, 56, 196.
8704
Ind. Eng. Chem. Res., Vol. 49, No. 18, 2010
(35) Tomisic, V.; Simeon, V. Raman spectra of aqueous solutions of strong electrolytes: Evolving-factor analysis. Phys. Chem. Chem. Phys. 1999, 1, 299. (36) Malinowski, E. R. Factor Analysis in Chemistry, 2nd ed.; Wiley: New York, 1992. (37) Hyvarinen, A.; Oja, E. Independent component analysis: Algorithms and applications. Neural Networks 2000, 13, 411. (38) Zhao, C. H.; Gao, F. R.; Wang, F. L. Spectra Data Analysis and Calibration Modeling Method Using Spectra Subspace Separation and Multiblock Independent Component Regression Strategy. AICHE J., published online Jun 21, 2010, http://dx.doi.org/10.1002/aic.12333. (39) Lee, J.; Qin, S. J.; Lee, I. Fault detection and diagnosis based on modified independent component analysis. AIChE J. 2006, 52, 3501. (40) MacGregor, J. F.; Jaeckle, C.; Kiparissides, C.; Koutoudi, M. Process monitoring and diagnosis by multiblock PLS methods. AIChE J. 1994, 40, 826. (41) Kourti, T.; Nomikos, P.; MacGregor, J. F. Analysis, monitoring and fault-diagnosis of batch processes using multiblock and multiway PLS. J. Process Control 1995, 5, 277. (42) Westerhuis, J. A.; Kourti, T.; MacGregor, J. F. Analysis of multiblock and hierarchical PCA and PLS models. J. Chemom. 1998, 12, 301.
(43) Qin, S. J.; Valle, S.; Piovoso, M. J. On unifying multiblock analysis with application to decentralized process monitoring. J. Chemom. 2001, 15, 715. (44) Lopes, J. A.; Menezes, J. C.; Westerhuis, J. A.; Smilde, A. K. Multiblock PLS analysis of an industrial pharmaceutical process. Biotechnol. Bioeng. 2002, 80, 419. (45) Yu, H. L.; MacGregor, J. F. Post processing methods (PLS-CCA): Simple alternatives to preprocessing methods (OSC-PLS). Chemom. Intell. Lab. Syst. 2004, 73, 199. (46) Anderson, T. W. An Introduction to MultiVariate Statistical Analysis, 2nd ed.; Wiley: New York, 1984. (47) Burnham, A. J.; Viveros, R.; MacGregor, J. F. Frameworks for latent variable multivariate regression. J. Chemom. 1996, 10, 31. (48) Wulfert, F.; Kok, W. T.; Smilde, A. K. Influence of temperature on vibrational spectra and consequences for the predictive ability of multivariate models. Anal. Chem. 1998, 70, 1761.
ReceiVed for reView April 15, 2010 ReVised manuscript receiVed June 13, 2010 Accepted July 24, 2010 IE100892Y