Score–moment combined linear discrimination analysis (SMC-LDA) as an improved discrimination method

Abstract

A new discrimination method called the score–moment combined linear discrimination analysis (SMC-LDA) has been developed and its performance has been evaluated using three practical spectroscopic datasets. The key concept of SMC-LDA was to use not only the score from principal component analysis (PCA), but also the moment of the spectrum, as inputs for LDA to improve discrimination. Along with conventional score, moment is used in spectroscopic fields as an effective alternative for spectral feature representation. Three different approaches were considered. Initially, the score generated from PCA was projected onto a two-dimensional feature space by maximizing Fisher's criterion function (conventional PCA-LDA). Next, the same procedure was performed using only moment. Finally, both score and moment were utilized simultaneously for LDA. To evaluate discrimination performances, three different spectroscopic datasets were employed: (1) infrared (IR) spectra of normal and malignant stomach tissue, (2) near-infrared (NIR) spectra of diesel and light gas oil (LGO) and (3) Raman spectra of Chinese and Korean ginseng. For each case, the best discrimination results were achieved when both score and moment were used for LDA (SMC-LDA). Since the spectral representation character of moment was different from that of score, inclusion of both score and moment for LDA provided more diversified and descriptive information.

---

Introduction

Principal component analysis (PCA) combined with linear discriminant analysis (LDA) is a well-accepted chemometric discrimination method.^1–6 The score from PCA effectively represents spectral variation using principal components (factors or eigenvectors) with far fewer variables. It is frequently utilized for the analysis of diverse vibrational spectroscopic data. Additionally, LDA is a popular method for maximizing discrimination among clusters or groups.

Even though PCA-LDA has been widely used, improved discrimination is always needed, especially for differentiation of minute spectral features. With this in mind, a new algorithm, called score–moment combined linear discrimination analysis (SMC-LDA), was developed. This method utilizes the moment of spectrum as an alternative descriptor of spectral features along with PCA score for LDA. Moment is frequently adopted to describe the texture of an image in image processing fields, such as image retrieval systems.^7–11 It effectively describes the outline and overall shape of the object using a few scalar values (moments). Also, it is known to yield a relatively stable description under situations involving a minor shift or rotation of the image. The concept of moment has been extended to vibrational spectroscopic fields in this study and utilized as an alternative and effective descriptor of spectral features.

To demonstrate the discrimination ability of the proposed method, three spectroscopic datasets were used: (1) infrared (IR) spectra of normal and malignant stomach tissue, (2) near-infrared (NIR) spectra of diesel and light gas oil (LGO) and (3) Raman spectra of Chinese and Korean ginseng. First, moments and principal component (PC) scores of spectra were used separately for LDA. Then, both moments and scores were utilized for LDA (SMC-LDA). Discrimination performance was improved for all three datasets when using SMC-LDA.

Score–moment combined linear discrimination analysis (SMC-LDA)

The main concept of the proposed method was to combine the score and moment as descriptors (representation method) of spectral features and inputs for linear discriminant analysis (LDA). PCA involves the representation of a spectrum as scores that are projected from the principal component (PC) and that describe the largest spectral variations in a given dataset. Since it is widely known, PCA will not be discussed in detail here. However, moment is not frequently used in the field of vibrational spectroscopy, so it will be addressed at length.

Moment is used to describe image texture in the fields of image retrieval and image processing. In this study, moment is proposed as an alternative method of spectral representation that could be comparable to score from PCA. Detailed mathematical descriptions of moment follow. Initially, raw spectra (or preprocessed data such as baseline corrected), x(λ), are Fourier transformed to achieve X(jω) by use of eqn (1):

where and ω is angular frequency. Maclaurin's expansion is applied to the exponential part in eqn (1) to yield:where m_n is the nth moment of x(λ). The term m_n can be derived as:

The inverse Fourier transform of eqn (2) is then performed to yield eqn (4):

Raw spectra x(λ) can now be expressed as moments (m₁, m₂, … , m_n) that represent the spectral features as a few scalar values. The meanings of individual moments are difficult to translate for spectroscopic data. However, in the image processing field, it is known that the first, second, third and fourth moments represent average, standard deviation, skew and slope of a shape, respectively. Higher order moments are harder to conceptualize. Nevertheless, moment can be used as an effective “spectral feature descriptor” in spectroscopic fields.

A noticeable feature of moment is that it is calculated in the Fourier-transformed domain where the original spectrum is decomposed into individual sinusoidal components (angular frequency in eqn (2)). In this domain, minor differences in sinusoidal components derived from the original spectrum can be observed. This is the most distinguishable portion compared to the PCA score.

LDA is used to maximize the discrimination between classes. In the case of using moment, each moment is projected onto two-dimensional (2D) feature space by maximizing Fisher's criterion function. The Fisher's criterion function J(d) is given by:

B and W are distances between classes and within class, respectively (d^T: transpose of d). B and W are defined as:

where μ is the mean of m_ij. In these equations, N_c is the number of classes (C₁, C₂,…, C_Nc) in the dataset. Also, N_i is the number of samples in a given class and m_ij is the moment corresponding to the ith sample and the jth dimension. The vector d, that maximizes the criterion function J, is a series of eigenvectors corresponding to eigenvalues of the matrix W⁻¹B. d₁ is the eigenvector corresponding to the maximum eigenvalue of matrix W⁻¹B and forms a one dimensional feature axis. To describe more spectral variations, d₂ is also utilized and then forms an additional one dimensional feature axis. Finally, matrix D (feature space) is defined as:

D = [d₁d₂] (8)

Moment is projected onto the feature space (formed by d₁ and d₂) and then the projected value y on the feature space is described as:

y_ij = D^Tm_ij (9)

Consequently, two dimensional plot of projected value y is generated viaeqn (9) and this projection plot is utilized for discrimination of classes throughout this study.

The same procedure [eqn (5)–(9)] was performed for score as well as score–moment combinations. From this point, three abbreviations will be used: S-LDA, M-LDA and SMC-LDA. S-LDA, that is PCA-LDA, uses score for LDA and M-LDA utilizes only moment for LDA. SMC-LDA corresponds to the use of both score and moment for LDA.

Experimental

Spectral collection and data processing

IR spectra were collected using IR microscopy equipped with an MCT detector (SensIR Technologies, Danbury, CT). The aperture size was 25 µm. Stomach tissues of four different subjects containing both normal and malignant parts were acquired from Guro Medical Center of Korea University, Seoul, Korea. Tissue samples were freeze dried, cut (10 µm in thickness) and transferred onto a slide glass. IR spectra were acquired by accurately focusing on different areas of the sample based on pathological diagnostic information and with the help of a visual image. Each IR spectrum corresponded to 32 scans with a resolution of 4 cm⁻¹. A total of 87 spectra and 134 spectra were collected for normal and malignant tissues, respectively.

Forty-seven diesel and light gas oil (LGO) samples were acquired from a refinery in Korea.¹² NIR spectra were collected over the 1100 to 1700 nm region using a dispersive NIR spectrometer (Foss NIRSystems, Silver Spring, MD) equipped with a tungsten halogen lamp, PbS detector, and a fiber optic interactance/reflectance probe. The resolution of the collected spectra was 10 nm with 2 nm intervals between data. NIR spectra were collected by positioning the fiber optic probe into a sealed bottle that contained the sample. Each NIR sample spectrum consisted of 16 co-added scans.

Forty-two Korean and Chinese ginseng samples were acquired from the National Agriculture Production Inspection Office (NAPIO) in Seoul, Korea.¹³ All samples were powdered using a cyclone mill fitted with a 1 mm screen. The particle size of the powder was maintained below 20 mesh. FT-Raman spectra were collected on a Bruker (Karlsruhe, Germany) IFS55 spectrometer equipped with a FRA 106 Raman module, CaF₂ beam splitter and a cryogenically cooled Ge detector. The excitation source was a diode-pumped Nd:YAG laser operated at 1064 nm and a power of 300 mW. Each spectrum was collected with a resolution of a 4 cm⁻¹ in a 180° back-scattering arrangement. Descriptions of these three spectroscopic datasets are presented in Table 1.

Table 1 Descriptions of three vibrational spectroscopic datasets

Spectroscopy used	Sample	Number of spectra
IR	Stomach (normal)	87
	Stomach (malignant)	134
NIR	Diesel	47
	Light gas oil (LGO)	47
Raman	Chinese ginseng	42
	Korean ginseng	42

---

All algorithms were developed using Matlab Version 6.5 (The MathWorks Inc., MA, USA). For PCA, the commercial Matlab program package was used. All the pre-processing methods such as baseline correction or normalization were also accomplished using Matlab.

Results and discussion

Spectral features

Fig. 1(a) shows averaged IR spectra (1800–900 cm⁻¹) of normal and malignant stomach tissues from four different patients. The spectra were baseline-corrected at 1800 and 900 cm⁻¹ and normalized by the area under the 1800–900 cm⁻¹ range. By using normalization, unwanted spectral variations from different sample thicknesses were reduced. As shown, overall spectral features between normal and malignant tissues were generally similar to each other. However, moderately different spectral features were observed in the 1350 to 900 cm⁻¹ range. No noticeable spectral differences were observed in the 1800–1350 cm⁻¹ range where amide I (centered at 1650 cm⁻¹) and amide II (centered at 1550 cm⁻¹) are located. It is beyond the scope of this paper to describe the spectral differences based on medical and pathological origins. However, it is apparent that there are small spectral differences between normal and malignant tissues. Accurate discrimination of these differences is important for biomedical diagnosis based on spectroscopy.

Fig. 1 Three different types of vibrational spectra used in this study: IR spectra of normal and malignant stomach tissues (a), NIR spectra of LGO (light gas oil) and diesel (b), and Raman spectra of Korean and Chinese ginseng samples (c).

Fig. 1(b) represents selected NIR spectra of LGO and diesel. Diesel is a blended product of three to five components. LGO is a major component constituting around 80% of diesel. Therefore, NIR spectral features of LGO and diesel are similar. The bands centered at 1210 nm and 1400 nm corresponded to the second overtone and combination bands of CH vibration. As shown in Fig. 1(b), only small spectral differences were observed at 1210 nm.

Fig. 1(c) shows selected Raman spectra of Korean and Chinese ginseng samples. Even with the use of a low-energy near-infrared laser, the baselines of the spectra varied slightly due to weak fluorescence. Spectral differences between Korean and Chinese samples were observed at 940, 870 and 500 cm⁻¹. Also, the peaks in the 1500–1000 cm⁻¹ range differed slightly from each other. Before performing the spectral processing, the weak fluorescence background that was superimposed on the Raman spectra was effectively removed using Wavelet Transformation (WT).^14–16 A detailed investigation of the chemical and botanical origins of the differences is beyond the scope of this paper. Overall, the spectral features of each case are similar to each other.

Spectral representation by principal component scores and moment

To demonstrate the spectral description character of a moment, an extra dataset was used that had a quantitative property for each sample. The samples in the three datasets used in this study did not have quantitative values. It was easier to visualize the spectral description characteristics between score and moment using a dataset with quantitative properties. NIR spectra of gasoline samples (225 samples collected over a period of 12 months, research octane number (RON) range: 85.6–94.0) were examined. The main purpose of collecting the samples over a year was to include compositional variations based on four seasonal grade changes. Two-dimensional (2D) sample plots were constructed by two moments (moments 3 and 4) and two principal component (PC) scores (scores 1 and 2). These results are shown in Fig. 2. On each plot, all gasoline samples are represented by moments and PC scores, respectively. Two samples located at the edges of the cluster and one sample in the middle of the cluster were selected (3 total samples designated as A, B and C). The two samples (A and C) located on both edges of the cluster represented considerable spectral and RON differences. The corresponding three spectra for each case are also shown in Fig. 2. Spectral variation was successfully represented with the use of moments. RON variations and spectral variations of the three samples were distinct since they were located far apart in the cluster. The aromatic bands centered around 4600 cm⁻¹ were clearly increasing with respect to the RON. Aromatic compounds usually have high RON values. Similar results were observed for the PC scores, but the trend was different from that of the moment. Overall, the results in Fig. 2 show that the moment successfully represented the spectral variation.

Fig. 2 Two-dimensional sample plots constructed by two moments (moments 3 and 4) as well as two principal component (PC) scores (scores 1 and 2) for gasoline samples. Three samples were selected and the corresponding spectra are shown on the right plots.

Fig. 3 shows two-dimensional (2D) sample plots constructed by two scores and two moments for three cases of IR, NIR and Raman data. By using these plots, the difference in spectral description between score and moment was compared. The scores and moments were arbitrarily selected to show the different characteristics of spectral representation between score and moment in two dimensional space. For each 2D plot for score, 5 samples that showed large differences were chosen. The same samples were used in the 2D plot for moment. In the moment 2D plot for IR stomach data (open circles: normal tissue, filled circles: malignant tissue), the selected 5 samples were well-described and widely spread. The moment distribution of these samples was different from that of score, as seen in Fig. 2. Similar results were observed for the NIR (open circles: diesel, filled circles: LGO) and Raman data (open circles: Korean ginseng, filled circles: Chinese ginseng). Overall, results support the moment as an alternative, effective and comparable spectral description method along with score from PCA. It is similar to photographing an object at a different angle by a different photographer. Therefore, when score and moment are combined for discrimination, more descriptive and diverse information would be available. This will eventually lead to improvements in discrimination performance.

Fig. 3 Two-dimensional sample plots constructed by two scores and two moments for the three cases of IR, NIR and Raman data. Open and filled circles in IR stomach data correspond to normal tissue and malignant tissue, respectively. Open and filled circles in NIR data represent diesel and LGO, respectively. Open and filled circles in Raman data represent Korean and Chinese ginseng, respectively.

Comparison performance of S-LDA, M-LDA and SMC-LDA

To compare the discrimination ability of the three methods, an objective quantitative parameter for judging the degree of discrimination was needed. When two classes (A and B) are encountered in discrimination as in this study, it is possible to describe the degree of discrimination using the distance between two mean points of each class as well as sample distances from the mean point of each class. Mean points of the two classes were calculated and designated as C_A and C_B. For all the samples, Mahalanobis distances from C_A (MD_A) and C_B (MD_B) were determined. Then MD_B ÷ MD_A was calculated. When the MD_B ÷ MD_A value was less than one, a sample could be assigned to Class B. When the MD_B ÷ MD_A value was over one, it was assigned to Class A. Finally, discrimination accuracy was obtained by calculating the percentage of samples with correct assignment.

The cross validation method was essentially used to minimize possible over-fitting especially in the course of LDA. For this purpose, each data set was divided into 5 segments (80% in calibration and 20% in prediction) and cross-validated. Then, the corresponding discrimination accuracies from the calibration and prediction sets were evaluated.

Initially, the optimal number of scores for S-LDA and moments for M-LDA needed to be determined. To determine the optimal number of scores, the first ten scores were generated and continuously added one by one from the first score. By evaluating the pattern of increasing discrimination accuracy as a function of the number of scores, the optimal number of scores was determined when the best discrimination was achieved. The same procedure was used to determine the optimal number of moments. Fig. 4 shows discrimination accuracy plotted as functions of the number of factors (a) and moment (b) used to discriminate between normal and malignant tissues (IR data). Additionally, total percent variance (TPV) is also displayed in the plot (a). TPV is an indicator of how much variation is accounted for by factors (principal components, PCs). Factors represent the variation in the spectral data set, while eigenvalues are the relative weights of each individual PC. By summing eigenvalues and representing as a percentage, it can be estimated how much variance is described by the PCs. For S-LDA, four factors were chosen as an optimum since no significant improvement in discrimination accuracy as well as TPV after four factors. For M-LDA, the similar trend of improving discrimination accuracy was observed by adding additional moments and the use of six moments provided the optimal discrimination accuracy without over-fitting. The same procedures were accomplished to determine the numbers of scores and moment for other two data sets (NIR and Raman data sets).

Fig. 4 Discrimination accuracy plotted as functions of the number of factors (a) and moment (b) used to discriminate between normal and malignant tissues (IR data). Additionally, total percent variance (TPV) is also displayed in the plot (a).

Even though the spectral description characters between score and moment are different, possibly there is an overlap in feature description or some interactions between them. Therefore, the different approach has been accomplished to determine the optimal number of scores and moments for SMC-LDA. Initially, the number of scores is fixed as determined in Fig. 4 and then individual moment is added as combinatorial fashion. For example, 15 combinations of moment (₆C₂) are possible when 2 moments are used out of total 6 moments. Twenty combinations of moment (₆C₃) are possible when 3 moments are used out of total 6 moments. The maximum number of moments used corresponded to the optimal number of moments that was determined in M-LDA for each case. At each combination of moments with fixed number of scores, the corresponding discrimination accuracy was evaluated and the best combination was determined.

Table 2 shows the overall results of cross-validated discrimination accuracy for S-LDA, M-LDA and SMC-LDA. The numbers in parentheses from S-LDA and M-LDA results correspond to the optimal number of scores and moments, respectively. The numbers in parentheses from SMC-LDA results correspond to the optimal number of moments with the use of fixed number of scores from S-LDA. For IR data, the same number of moments as used in M-LDA was used for SMC-LDA. In the case of SMC-LDA for NIR data, the use of only one moment (1st moment) fairly improved the discrimination accuracy. In the case of SMC-LDA for Raman data, the use of three additional moments (1st, 6th and 7th moments) along with three scores led to 100% discrimination accuracy. Although it is difficult to justify how these moments selectively help to improve the discrimination, it is clear, however, based on the result, that the use additional moments provides complementary spectral information to score from PCA.

Table 2 Overall results of discrimination accuracy for S-LDA, M-LDA and SMC-LDA. The numbers in parentheses correspond to the optimal number of scores and moments, respectively

		IR data (1800–900 cm^−1)	NIR data (1100–1600 nm)	Raman data (1800–200 cm^−1)
S-LDA	Calibration	96.8 (4)	94.4 (3)	99.7 (3)
	Prediction	96.8	92.8	97.6
M-LDA	Calibration	96.4 (6)	94.4 (6)	98.8 (7)
	Prediction	96.3	93.7	98.9
SMC-LDA	Calibration	98.4 (6)	97.3 (1)	100.0 (3)
	Prediction	97.3	96.9	100.0

---

The discrimination accuracies from S-LDA and M-LDA are similar each other, or slightly improved when M-LDA is used based on the results from both NIR and Raman data. Therefore, moment can be a successful spectral representation tool compatible with conventional score. When scores and moments were combined for LDA (SMC-LDA), the best discrimination accuracies were achieved for all three data sets. As shown in Fig. 2 and 3, spectral description characteristics of score and moment differed from each other. The discrimination accuracy was improved by the input of complementary information into the LDA. The corresponding 2D projection plots generated by using the eqn (9) for IR (discrimination between malignant and normal stomach tissues), NIR (discrimination between LGO and diesel) and Raman (discrimination between Chinese and Korean ginseng) data are shown in Fig. 5, 6 and 7, respectively.

Fig. 5 Two-dimensional projection plots using S-LDA, M-LDA and SMC-LDA for the discrimination between normal and malignant stomach tissue using IR spectroscopy. Open and filled circles correspond to normal and malignant tissues, respectively.

Fig. 6 Two-dimensional projection plots using S-LDA, M-LDA and SMC-LDA for the discrimination between LGO and diesel using NIR spectroscopy. Open and filled circles correspond to diesel and LGO, respectively.

Fig. 7 Two-dimensional projection plots using S-LDA, M-LDA and SMC-LDA for the discrimination between Chinese and Korean Ginseng using Raman spectroscopy. Open and filled circles correspond to Chinese and Korean ginseng, respectively.

Conclusion

A new discrimination method, SMC-LDA, which can improve prediction accuracy over PCA-LDA was presented. Moment was uniquely utilized as an alternative spectral representation tool and its performance was demonstrated. Three different vibrational spectroscopic datasets acquired from diverse fields served to evaluate the performance and applicability of SMC-LDA. Nonetheless, additional studies are required to assure the improved performance of SMC-LDA. SMC-LDA performance will be investigated further using systematically designed datasets, such as signal-to-noise ratio change, baseline variations and bandwidth dependency, that are encountered frequently in general spectroscopic data. Also, SMC-LDA will be applied to more diverse spectral datasets, such as natural products, where minute spectral discrimination is critically necessary under highly complex chemical environments.

Acknowledgements

This work was supported by a grant (#20050401034632) from Biogreen 21 Program, Rural Development Administration, Republic of Korea.

References

1. H. Martens and T. M. Naes, Multivariate Calibration, John Wiley and Sons, New York, 1989.
2. K. R. Beebe, R. J. Pell and M. B. Seasholtz, Chemometrics: A Practical Guide, John Wiley and Sons, New York, 1998.
3. W. J. Krzanowski, J. Chemom., 1992, 6, 97–102.
4. Y. Mallet, D. Coomans and O. Devel, Chemom. Intell. Lab. Syst., 1996, 35, 157–173.
5. A. Baldovin, W. Wen and D. L. Massart, Chemom. Intell. Lab. Syst., 1997, 38, 25–37.
6. Z. Ramadan, M. Mulholland and D. B. Hibbert, Chemom. Intell. Lab. Syst., 1998, 40, 165–174.
7. M. Gibaldi and D. Perrier, Noncompartmental analysis based on statistical moment theory, in Pharmacokinetics, Marcel Dekker Inc., New York, 1982, pp. 409–416.
8. L. Gupta and M. D. Srinath, Pattern Recognit., 1987, 20(3), 267–272.
9. J. Tajima, Pattern Recognit., 1981, 14, 211–217.
10. D. W. Thomas and B. R. Wilkins, Pattern Recognit., 1972, 4(4), 379–389.
11. M. F. Zakaria, L. J. Vroomen, P. J. A. Zsombor, J. M. Murray and H. M. Kessel, Pattern Recognit., 1987, 20(6), 639–643.
12. H. Chung, H. J. Choi and M. S. Ku, Bull. Korean Chem. Soc., 1999, 20(9), 1021–1025.
13. Y. A. Woo, H. J. Kim and H. Chung, Analyst, 1999, 124, 1223–1226.
14. C. K. Chui, Introduction to Wavelets, Academic Press, San Diego, CA, 1991.
15. B. Walczak and D. L. Massart, Chemom. Intell. Lab. Syst., 1997, 36, 81–94.
16. X. Shao, A. K. Leung and F. Chau, Acc. Chem. Res., 2003, 36, 276–283.

---