Component Detection Weighted Index of Analogy: Similarity

Feb 17, 2005 - Multivariate data analysis in pharmaceutics: A tutorial review. Tarja Rajalahti , Olav M. Kvalheim. International Journal of Pharmaceut...
2 downloads 8 Views 909KB Size
Anal. Chem. 2005, 77, 1607-1621

Component Detection Weighted Index of Analogy: Similarity Recognition on Liquid Chromatographic Mass Spectral Data for the Characterization of Route/Process Specific Impurities in Pharmaceutical Tablets Simeone Zomer,† Richard G. Brereton,*,† Jean-Claude Wolff,‡ Christian Y. Airiau,§ and Caroline Smallwood‡

Centre for Chemometrics, School of Chemistry, University of Bristol, Cantock’s Close, Bristol BS8 1TS, U.K., GlaxoSmithKline, Gunnels Wood Road, Stevenage, Hertfordshire SG1 2NY, U.K., and GlaxoSmithKline, Old Powder Mills, Tonbridge, Kent TN11 9AN, U.K.

Detection and identification of impurities in pharmaceuticals is an essential task for determining the possible infringement of a patent. This article reports a multivariate analysis method to distinguish between tablets of the same substance on the basis of their origin, by characterizing route/process specific impurities via diagnostic ion chromatograms, using liquid chromatography/mass spectrometry (LC/MS). The approach is based on the formulation of a novel index that quantifies the similarity between LC/ MS samples, named the component detection weighted index of analogy. The index estimates similarity by fully exploiting the two-dimensional nature of the data, where the relative contribution of chromatograms relates to their quality and noise level. Results show that well-defined clusters are formed according to the origin of tablets; a series of ions are identified as characterizing each class and can be used to predict the origin of unknown tablet samples. The method presented is designed for analysis of larger data sets and can be suitable for exploratory analysis where any a priori knowledge on the data is scarce or absent, hence requiring the acquisition of chromatograms in a broad m/z range. Identification of impurities in pharmaceuticals is traditionally a vital task in the drug development process, but it can also be essential for the protection of intellectual property of marketed products. Typically, for a given pharmaceutical there would be a limited number of drug substance manufacturers, but a larger number of drug product manufacturers and distributors. Classifying the distributors with respect to the origin of the drug substance used and then determining the synthetic chemistry process is of great importance to the pharmaceutical industry to protect their assets and intellectual property. This can be achieved because in principle every synthetic chemistry route/process gives * To whom correspondence should be addressed. E-mail: r.g.brereton@ bris.ac.uk. † University of Bristol. ‡ GlaxoSmithKline, Stevenage. § GlaxoSmithKline, Tonbridge. 10.1021/ac048504t CCC: $30.25 Published on Web 02/17/2005

© 2005 American Chemical Society

rise to a unique impurity fingerprint where markers can be detected, although the task is complex in nature by the fact that these generally appear at a very low level, i.e., ppm or ppb of mass fraction levels relative to the drug substance. Because of its high sensitivity, a preferred analytical technique often used in such a context is liquid chromatography/mass spectrometry (LC/MS). A wide range of applications for the identification of impurities or degradants are reported in the literature1-4 with fewer examples that are also targeted on the determination of such route/process specific impurities (RSI).5,6 However, generally the interpretation of the LC/MS data is difficult because the instrument produces a vast amount of data in the form of several hundred chromatograms that also contain a variable amount of noise. Moreover, the analysis is not trivial since such impurities may reflect changes that are featured in only a very minor portion of the data, and any prior information to drive the analyst is often scarce or absent. As consequence, the development and application of nonstandard techniques of data analysis rooted in chemometrics7-9 can help the analysis. This paper describes a chemometric procedure for investigating the LC/MS fingerprint/signature of pharmaceutical tablets from the same formulated drug product, to characterize and differentiate their origin on the basis of RSI. The procedure relies on the formulation of a novel index of similarity for samples analyzed by LC/MS, named the component detection weighted index of analogy (CWIA). The approach proposed is holistic in the sense that it takes into account all the information available without explicitly extracting the spectral features of RSI but identifies the (1) Wu, Y. H. Biomed. Chromatogr. 2000, 14, 384-396. (2) Niessen, W. M. A. Chimia 1999, 53, 478-483. (3) Qin, X. Z.; Ip, D. P.; Chang, K. H. C.; Dradransky, P. M.; Brooks, M. A.; Sakuma T. J. Pharm. Biomed. 1994, 12, 221-233. (4) Volk, K. J.; Hill, S. E.; Kerns, E. H.; Lee, M. S. J. Chromatogr., B 1997, 696, 99-115. (5) Sheldon, E. M.; Downar, J. B. J. Pharm. Biomed. 2000, 23, 561-572. (6) Rudaz, S.; Souverain, S.; Schelling, C.; Deleers, M.; Klomp, A.; Norris, A.; Vu, T. L.; Ariano, B.; Veuthey, J. L. Anal. Chim. Acta 2003, 492, 271-282. (7) Brereton, R. G. Chemometrics Data Analysis for the Laboratory and Chemical Plant; Wiley: Chichester, U.K., 2003. (8) Lavine, B. K.; Workman, J. Anal. Chem. 2002, 74, 2763-2769. (9) Bessant, C.; Brereton, R. G.; Dunkerley, S. Analyst 1999, 124, 1733-1744.

Analytical Chemistry, Vol. 77, No. 6, March 15, 2005 1607

Figure 1. Data analysis flowchart.

chromatograms of the ions (i.e., certain m/z values) that are the most characteristic and useful for discrimination. Focusing on an optimal subset of ions is significant because it can help in speeding up any subsequent analysis of the same kind while also providing clues for a chemical interpretation of the phenomenon (e.g., which synthetic chemistry route was used or which were the reaction conditions of the process). The flowchart of the procedure is illustrated in Figure 1. The first step consists of the application of curve resolution techniques (CRTs) to remove the features due to the most intense peaks from the data: major chemicals relate to the tablet formulation (e.g., the drug substance itself) that is common to all samples, whose presence may potentially hide RSI. However, because they occur in relatively high concentrations, their traces can be accurately modeled and removed. In the second step, the component detection algorithm (CODA)10 is applied to estimate the quality/noise content of all the chromatograms recorded. Chromatographic mass spectral data often consist of several hundred chromatograms, some of which associate with noise or background rather than genuine peaks. These chromatograms worsen the quality of data and CODA operates by selecting an optimal subset where better quality information lies. Because of its simplicity and effectiveness, the use of CODA is widespread11-14 but suffers from the fact that it requires the user to define a subjective quality threshold to accept/reject chromatograms, and most importantly, when comparing several samples, each one may account for a different number and type of chromatogram retained, hence complicating (10) Windig, W.; Phalp, J. M.; Payne, A. W. Anal. Chem. 1996, 68, 3602-3606. (11) Barbarin, N.; Henion J. D.; Wu Y. H. J. Chromatogr., A 2002, 970, 141154. (12) Williams, A. Curr. Top. Med. Chem. 2002, 2, 99-107. (13) Windig, W.; Smith, W. F.; Nichols, W. F. Anal. Chim. Acta 2001, 446, 467476. (14) Windig, W.; Marchincin, T. F.; Meyer, G. N. Appl. Spectrosc. 2003, 57, 1575-1584.

1608

Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

any comparison that fully exploits the two-dimensional nature of LC/MS data. The final step of the flowchart is similarity recognition by means of CWIA: this index overcomes the limitations of CODA because it establishes a measure of similarity between pairs of samples by taking into account all chromatograms that contribute according to their quality. During the similarity recognition phase, CWIA is initially applied to all pairs of samples in the data set to produce a similarity matrix that can be used as input to a clustering algorithm, hence allowing the identification of the major groups. Subsequently, the use of information regarding the origin of tablets allows the identification of those ions that characterize each group upon which a reduced model for predicting the origin of unknown tablet samples can be constructed. As a result, well-defined clusters emerge according to the different sources investigated, and a series of ions, which possibly relate to RSI, are identified as characterizing each of the tablet sources. The proposed method for similarity recognition based on CWIA is designed for the handling of entire batches of samples (about 20-30 in this study) and can be suitable for exploratory analysis where a priori knowledge of the data is scarce or absent, hence requiring the acquisition of chromatograms in a broad m/z range. EXPERIMENTAL SECTION Materials. Drug product samples from different sources/ manufacturers, which may or may not use the same synthetic process or formulation (i.e., excipients may vary between the groups of samples), were fingerprinted by LC/MS. A proprietary drug of molecular formula C9H7Cl2N5 (nominal mass 255, [M + H]+ ) 256) was used as a model compound. The sample set (Table 1) consisted of drug product (i.e., tablets) from known sources G, I, and C, and an unknown sample, Un, that originates from one of the known sources. To simulate the real acquisition of new tablets that typically occurs at different stages over time, samples were analyzed at three different time points over a period of 3-4 weeks: the first batch consists of 27 samples (9 G, 8 I, 9 C, and 1 Un), the second of 26 (14 G, 5 I, 6 C, and 1 Un), and the third 25 (14 G, 4 I, 6 C, and 1 Un). The first batch also contains eight replicates to verify instrumental reproducibility and tabletto-tablet variation, and some samples are reanalyzed in each batch to verify the stability of the response over time. Sample Preparation. The equivalent tablet mass containing 5 mg of active pharmaceutical ingredient (i.e., drug substance) was dissolved in 1 mL of 75:25 volume fraction 0.05 M aqueous ammonium acetate and acetonitrile adjusted to pH 4 with acetic acid. The resulting solution was centrifuged and the supernatant used for analysis by LC/MS. Liquid Chromatography. An Agilent 1100 chromatography system (Agilent Technologies, Stockport, Cheshire, U.K.) was used. The column was a Hypurity Elite C18, 5 µm, 3 × 150 mm (Phenomenex, Macclesfield, Cheshire, U.K.). The flow rate was 0.6 mL/min, and the column temperature was set to 40 °C. Eluent A consisted of 0.05 M ammonium acetate in water, and eluent B consisted of acetonitrile. Elution started at 10% volume fraction of B, held for 2 min, and then linearly increased to 80% B over 23 min. The volume fraction of B was held at 80% for 5 min before re-equilibration at the initial conditions for 8 min.

Table 1. Data Set Analyzeda sample no, and origin G1 G2 G3 G4 G5 G6 G7 G8 G9 G 10 G 11 G 12 G 13 G 14 I1 I2 I3 I4 I 5 (replicate a) I 5 (replicate b) I 5 (replicate c) I 5 (replicate d) I 5 (replicate e) C1 C 2 (replicate a) C 2 (replicate b) C 2 (replicate c) C 2 (replicate d) C 2 (replicate e) C3 C4 C5 C6 Un

batch 1

batch 2

batch 3

x x x x x x x x x

x x x x x x x x x x x x x x x x x x x

x x x x x x x x x x x x x x

x x

x x

x x x x x x x x x x x x x x x x x x

x x x x x x

x x x x

x x x x x

a Batch 1 contains replicates to examine reproducibility. Many samples are common to all batches to study the stability of the instrumental response over time, while some others are unique of batch 2 and batch 3 in order to provide a test set for the model.

Mass Spectrometry. LC/MS experiments were carried out on a Waters/Micromass Q-ToF quadrupole orthogonal acceleration time-of-flight mass spectrometer (Waters UK Ltd., Wythenshawe, U.K.), equipped with the LockSpray dual electrospray ion source.15 The instrument was operated in positive ion electrospray (+ve ESI) mode, with both the reference and the sample sprays at 3 kV. The nitrogen desolvation and nebulizer gas flow rates were set to 400 and 90 L/h, respectively. The source temperature was set to 100 °C and desolvation temperature to 300 °C. The cone voltage was 25 V. The instrument was calibrated over 501200 Da mass range using [Glu1]-fibrinopeptide B (250 µmol in 0.1% acetic acid in water/methanol 1:1 volume fraction, Sigma Chemical Co., St. Louis, MO) in MS/MS mode. To match the calibration conditions, the collision gas (argon) was left on for all MS experiments. Spectra were acquired in the range 100-1000 Da over 30 min (a total of 1284 acquisitions/sample) in centroid mode and exported from Masslynx v4.0 (Waters UK Ltd.) raw data files into ASCII files using “Databridge” provided with Masslynx software. Finally, ASCII files were imported for processing into Matlab v6.0 (Mathworks Inc.). (15) Wolff, J.-C.; Eckers, C.; Sage, A. B.; Giles, K.; Bateman, R. Anal. Chem. 2001, 73, 2605-2612.

THEORY Removal of Main Components. A pharmaceutical tablet is a mixture that contains two primary types of ingredients: the drug substance that is the active pharmaceutical ingredient and the excipients (inert ingredients) that make up the composition of the formulated product and are critical to the proper delivery of the drug actives.16 In the LC/MS data acquired from a tablet sample, the traces of possible route/process specific impurities can be hidden by these major ingredients, particularly when coelution occurs. However, because of their relative abundance, their trace can be accurately modeled and removed. This goal is pursued with the use of curve resolution techniques:17-19 these are a set of chemometric methods that rely on a certain mathematical decomposition to deconvolve the two-way LC/MS data from its unresolved multicomponent mixture into factors of single species (their spectral and concentration profiles). The methods are based on the assumption that in each point in time along the profile the total spectrum recorded can be modeled as a linear combination of spectra of the coeluting species, each of which contributes in proportion to its concentration. The relationship can be written as

X ) C‚S

(1)

where X is the LC/MS data matrix with columns corresponding to chromatograms recorded for different m/z values and rows corresponding to mass spectra recorded at subsequent points in time, C is the concentration matrix that consists of the columns of the elution profiles c of all the species detected, and finally S is the spectral matrix that consists of the rows of the pure spectral profiles s of these species. The basic principle of CRTs is to seek for pure vectors c and s characterizing each compound in the mixture, that give the best fit to X by minimizing a quantity such as

E ) |X - C‚S|

(2)

where | | indicates the norm. To achieve accurate results in CRTs, it is critical to estimate the correct number of chemical components in the system initially, namely, the rank of matrix X, which can be a difficult task particularly when in the presence of high noise and correlation (i.e., similar elution and spectral profiles associate to the components).20 Rank estimation is often performed by principal component analysis (PCA),7,21 which decomposes the data matrix X according to the relationship

X ) T‚P

(3)

where T and P are the scores and loadings matrices and are composed of vectors (along column for T and along rows for P) (16) Rubbinstein M. H. In Pharmaceutics: The Science of Dosage Form Design; Aulton, M. E., Ed.; Churchill Livingstone: Edinburgh, Scotland, 1998. (17) Jiang, J. H.; Liang, Y.; Ozaki, Y. Chem. Intell. Lab. Syst. 2004, 71, 1-12. (18) Shen, H. L.; Grung, B.; Kvalheim, O. M.; Eide, I. Anal. Chim. Acta 2001, 446, 311-326. (19) Liang, Y. Z.; Kvalheim, O. M. Anal. Chim. Acta 1994, 292, 5-15. (20) Shen, H. L.; Wang, J. H.; Liang, Y. Z.; Pettersson, K.; Josefson, M.; Gottfries, J.; Lee F. Chem. Intell. Lab. Syst. 1997, 37, 261-269.

Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

1609

Figure 2. (a) Original TIC of sample C2. (b) TIC after the application of CRTs: a number of minor peaks can now be determined. (c) Ratio between the first and second eigenvalue of the residual data matrix X as a function of the number of iterations. As the trace due to the main component (i.e., active pharmaceutical ingredient) decreases, the ratio λ1/λ2 decreases because while λ1 associates to this main component, λ2 refers to other minor ones.

that are mutually orthogonal. The new components fit a model that maximizes the variance in X and are ordered in term of size, given by the sum of squares of the columns of T, which are known as eigenvalues λi. By examining the relative size of the eigenvalues it is possible to estimate the number of significant components that correspond to the rank of matrix X. In this application, a simplified version of the curve resolution technique known as alternating least squares (ALS)22 was used: at the beginning the original algorithm requires an estimate of the number of chemical components in the system described in X and a first guess of either the concentration or the spectral profiles of each component. Having X available from experiments and, for example, a first guess on the concentrations C ˆ , the spectral profile can be estimated by the pseudoinverse on eq 1:

Sˆ ) (C ˆ T‚C ˆ )-1‚C ˆ T‚X

(4)

and in turn, the concentration matrix C can be updated to

C ˆ ) X‚Sˆ T‚(Sˆ ‚Sˆ T)-1

(5)

where Sˆ and C ˆ are the estimates and the symbol T denotes transposition. The ALS algorithm iterates between eqs 4 and 5 until the reconstructed matrix from the estimates minimizes E in eq 2. Typically a series of constraints derived from a priori knowledge on the system such as nonnegativity on concentrations and spectra and unimodality of concentration profiles are also introduced.23 In the case under investigation, the task is simpler because X is largely dominated by the presence of a single compound (i.e., (21) Xu, C. J.; Liang, Y. Z.; Li, Y.; Du Y. P. Analyst 2003, 128, 75-81. (22) Karjalainen, E. J. Chem. Intell. Lab. Syst. 1989, 7, 31-38. (23) Gemperline P. J.; Cash E. Anal. Chem. 2003, 75, 4236-4243.

1610

Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

the active pharmaceutical ingredient) that we want to model and remove and that appears as a very large peak in the total ion current chromatogram (TIC) profile around elution time of 9 min (Figure 2a). Hence, given an initial rank equal to 1, a guess of the concentration profile is taken as the profile of the most abundant ion. Application of eqs 4 and 5 leads to some Sˆ and C ˆ estimates and a matrix of residuals XR that is calculated as

XR ) X - C ˆ ‚Sˆ

(6)

where XR should contain information on the remaining chemical components. Nonnegativity is imposed on XR by iteratively multiplying the product C ˆ ‚Sˆ by a constant 0.9 until such condition is satisfied. As consequence, some of the information held in XR may still refer to some residual interference due to the main component. The procedure therefore iterates over eqs 4-6 on XR until the matrix shows a significant variation in its rank, namely, when more than a significant chemical component emerges (Figure 2b). This can be performed by applying PCA to XR and determining the relative increase of the size of PC2 compared to PC1 by considering the eigenvalue ratio λ1/λ2 until this converges (Figure 2c). For those samples where a relatively abundant peak also related to presence of the excipient occurs, the initial rank is set equal to 2, and the procedure iterates analogously until the ratio λ2/λ3 converges. Denoising. The removal of the trace due to the main components results in structured data in X corresponding to a variety of minor compounds. However, X still contains noise both of instrumental and chemical origin in addition to genuine information, whose presence should be minimized. Since the inception of LC/MS, a number of methods have been proposed to address the problem in the context of spectrum extraction for compound identification and quantification. These methods are

Figure 3. Reduced TIC profile of the sample in Figure 1 after application of CODA using a MCQ threshold equal to 0.8. The profile is now cleaned from background and spikes with a number of peaks that clearly emerge. However, information related to minor analytes may also have been “thrown away” together with those chromatograms of poorer quality.

based on a variety of approaches such as simultaneous maximization of spectral peaks,24 backfolding,25 and model peak methods where spectra for individual components are extracted on the basis of an estimated shape similarity in different chromatograms.26,27 However, simpler approaches can be used when the analysis focuses on looking for characteristic m/z values featured in a certain class of samples rather than determining spectral components of the mixture, or more generally, the analysis involves evaluating the similarity of within/between sample classes, for example, by simply visualizing the TIC profile. Such methods are more holistic in the sense that they apply to entire portions of X (e.g., entire chromatograms or adjacent spectra) without explicit focus on chemical components. For example, sequential pair covariance (SPC)28 operates by multiplying intensities of adjacent unscaled un-normalized mass spectra: large intensities at each point will arise only if adjacent spectra have common features, while noise that is uncorrelated between neighbor spectra will be suppressed by the operation. CODA10 calculates for each chromatogram a quality index (MCQ) that is higher for those displaying genuine peaks, consequently those below a chosen threshold for MCQ are discarded. The morphological score (MPS)18 allows selection of low noise m/z values by measuring the oscillation of a profile around its mean. Finally, improved versions of SPC have been reported, as higher-order SPC29 that considers correlation between more than two consecutive scans and windowed mass selection30 that better deals with the influence due to background. (24) Biller, J. E.; Biemann, K. Anal. Lett. 1974, 7, 515-528. (25) Pool, W. G.; deLeeuw, J. W.; vandeGraaf, B. J. Mass Spectrom. 1997, 32, 438-443. (26) Dromey, R. G.; Stefik, M. J.; Rindfleisch, T. C.; Duffield, A. M. Anal. Chem. 1976, 48, 1368-1375. (27) Stein, S. E. J. Am. Soc. Mass Spectrom. 1999, 10, 770-781. (28) Muddiman, D. C.; Rockwood, A. L.; Gao, Q.; Severs, J. C.; Udseth, H. R.; Smith, R. D.; Proctor, A. Anal. Chem. 1995, 67, 4371-4375. (29) Muddiman, D. C.; Huang, B. M.; Anderson, G. A.; Rockwood, A.; Hofstadler, S. A.; WeirLipton, M. S.; Proctor, A.; Wu, Q. Y.; Smith, R. D. J. Chromatogr., A 1997, 771, 1-7. (30) Fleming, C. M.; Kowalski, B. R.; Apffel, A.; Hancock, W. S. J. Chromatogr., A 1999, 849, 71-85.

In this study, denoising relies on CODA because recognition in the data does not explicitly focus on spectral component detection. MPS is not exhaustive because the background is not automatically detected,31,32 while it was shown that SPC and related methods have drawbacks because signal suppression also occurs for ions that are relatively small, whose detection is essential in this application, but simultaneously those with initial large intensities are enhanced.30 In this paper, CODA was also suggested as the preferred method when data have consistent background and peaks are smooth and continuous. In CODA, the quality index of each chromatogram is calculated by

MCQj )

1

r-w+1

xr - w

i)1

∑ a(π) a(w,s) ij

ij

(7)

where a(π)ij is the length scaled data at time point i and chromatogram j, a(w,s)ij is the smoothed mean-subtracted data, w is the window size for smoothing, and r is the number of rows in the matrix, with MCQ taking values in the range from 0 to 1. The underlying principle is that chromatograms of good quality (higher MCQ) will display genuine peaks without spikes and background, and this type of profile a(π)j will be similar to the corresponding smoothed, mean subtracted version a(w,s)j. This is because a significant number of spikes, typically due to instrumental noise, will result in dissimilarity compared to the smoothed version of the profile, while high background due to chemical noise will result in dissimilarity compared to the meancentered version; hence, in both cases MCQ will tend to decrease. CODA is therefore applied by selecting a subportion of chromatograms displaying MCQ over a preselected threshold and discarding the remaining ones from further analysis. Figure 3 illustrates (31) Shen, H.; Carter, J. F.; Brereton, R. G.; Eckers, C. Analyst 2003, 128, 287292. (32) Zomer, S.; Brereton, R. G.; Carter, J. F.; Eckers, C. Analyst 2004, 129, 175181.

Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

1611

the TIC obtained from applying CODA with a threshold of 0.8 on the same sample of Figure 2. A number of peaks can now be clearly recognized. Component Detection Weighted Index of Analogy. 1. Similarity Recognition. A first method to compare samples after the application of CODA is to consider the reduced TIC profiles determined using an equal MCQ threshold for all samples. For each pair of profiles x and y, it is possible to calculate a similarity index q of the form I

∑x ‚y i

q)

i

i)1

x∑ ∑ I

xi2

i)1

with

0eqe1

(8)

I

yi2

i)1

This is a simple approach for similarity recognition that does not require the tuning of additional parameters with the advantage that it can be easily extended to larger data sets: q can be estimated for each pair of samples, hence producing a similarity matrix (a square matrix with a size matching the size of the data set where each entry corresponds to the similarity measure for each pair of samples) that can be used as input to classification and clustering algorithms.7 Because in this application profiles are largely dominated by the peak of the active ingredient (Figure 2a) that characterizes the entire data set, alignment can be performed by taking the position of this peak as reference. However, such an approach for determining similarity is limited because it relies on TIC profiles where information is lost along the mass direction and it does not take full advantage of the twodimensional nature of the data. This is because if an equal MCQ threshold in CODA is used for all samples, different types and numbers of chromatograms may be retained in each case, making difficult any comparison on a two-dimensional basis. This can be serious when there is instrumental instability, where noise changes over the time of analysis, or when minor analytes appear in certain samples, mixing their traces with noisy chromatograms. In the latter case, the MCQ will increase because peaks will now appear close to the noise level. Moreover, another disadvantage of considering the reduced TICs is that this does not provide information about which features in the data (i.e., certain m/z values) are responsible for the clustering of samples. Other methods have also been proposed13,14 that are able to identify characteristic chromatograms: here, after the application of CODA, the search is performed on the union of the chromatograms selected, namely, those appearing in at least one sample. Then, investigation of additional features of the signal such as intensity, peak width, and peak area ratio between different profiles allows the identification of the most characteristic ions. These methods are effective and provide a valuable tool when comparing a small number of samples, but unlike the analysis on reduced TICs, they are not designed for the analysis of larger data sets. Finally, in both the methods discussed above, a general drawback in the use of CODA as a primary step is that the selection of the optimal MCQ threshold remains subjective: in the original paper10 a value equal to 0.85 was suggested, but elsewhere30 authors pointed out that lower values are more 1612

Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

appropriate when the intensities of analyte ions are low and get closer to the noise level and that ultimately the choice should depend on the consistency of background and noise level. A drawback in this application is that the imposition of such threshold may result in discarding chromatograms related to impurities close to noise level, which can still be important for characterization. Here we propose a novel approach for similarity recognition, CWIA, that goes beyond the limitations illustrated above as follows: (a) In contrast to both the approaches, CWIA uses all the information available. CODA is still applied as primary step but without the imposition of a rejection threshold; hence, all the chromatograms recorded contribute to determine similarity. (b) In contrast to the analysis of the reduced TICs, CWIA extracts from the data the features that are responsible for differentiating groups of samples. (c) In contrast to the methods in refs 13 and 14, CWIA is designed for the analysis of larger data sets. The proposed procedure develops along three steps: the first is the detection of the clusters in the data set; the second is the detection of the ions that are responsible for the generation of such clusters; the third is the use of the optimal subset of ions for predicting the origin of the tablet samples. 2. Clustering. The formulation of CWIA is as follows: J

∑N ‚q j

CWIAX,Y )

j

j)1

with

J

∑N

0 e CWIAX,Y e 1

(9)

j

j)1

where qj is the similarity index of eq 8 calculated for the m/z value j between samples X and Y, J is the total number of ion chromatograms recorded, and Nj is the weighting factor of the m/z value j that is calculated as follows:

Nj )

exp(τI‚min (MCQxj,MCQyj)) - 1 exp(τI) - 1

(10)

In practice, Nj penalizes a m/z value when either one of related chromatograms is of poor quality through the “min” operator applied to the MCQ values, while it favors in an exponential manner m/z values whose chromatograms are of high quality. The parameter τI is tuneable and defines the weights over the MCQ range (Figure 4): the higher this is, the more the relative weight ascribed to higher MCQs. Positive values for τI should be preferred, so that m/z values related to high-quality chromatograms are overweighted while low-quality ones remain underweighted (a method for the fine-tuning of τI is discussed below). The CWIA values can then be used to construct a similarity matrix Ψ, where each entry corresponds to the CWIA for each pair of samples. Ψ can be used as input to clustering algorithms, from which samples will emerge grouped together in clusters according to their similarity. 3. Characterization. As formulated in eq 9, CWIA does not provide explicit insight as to which features in the data are responsible for the clustering of samples. This is important as

Figure 4. Weights assumed by each m/z value as function of the quality parameter MCQ (eq 10) and coverage (eq 14) of the related chromatogram for different values of τI andτII, respectively. τI is set equal to +3 (dashed line), and τII is set equal to -5 (dotted line).

often the analyst wants to know what are the features to look at in a certain group of samples, i.e., the ion chromatograms that discriminate between classes or that are characteristic for each one of them. CWIA can be used for this purpose if a single m/z value each time is given as input to the algorithm and information on class membership (i.e., origin of tablet samples) is available. For a single m/z value j eq 9 becomes

CWIAxj,yj ) Nj‚qj

(11)

Then, using the information on class membership, the similarity matrix Ψj derived from CWIA for a m/z value j can be compared with a target matrix Ω. For example, if a batch is composed of 10 samples, 3 of class membership C1, 5 of membership C2, and 2 of membership C3, Ω assumes the form

[ ]

1 1 1 0 0 Ω) 0 0 0 0 0

1 1 1 0 0 0 0 0 0 0

1 1 1 0 0 0 0 0 0 0

0 0 0 1 1 1 1 1 0 0

0 0 0 1 1 1 1 1 0 0

0 0 0 1 1 1 1 1 0 0

0 0 0 1 1 1 1 1 0 0

0 0 0 1 1 1 1 1 0 0

0 0 0 0 0 0 0 0 1 1

0 0 0 0 0 0 0 0 1 1

(12)

The Ω matrix is the target in the sense that, in a perfect case, the similarity matrix Ψj will perfectly resemble Ω. This is because CWIA takes values in the range [0,1], and while 1 represents a perfect match (pair of samples having same membership), 0 represents a perfect mismatch (pair of samples having different membership). If one wants to model one class, a sub-block of Ω can be considered in alternative, for examples for group C1 the first three rows. Following this scheme, the new index that

estimates the relevance of a certain m/z value j in discriminating classes can be formulated as follows: Ω )‚CFj dj ) (1 - EΨj

(13)

Ω is the error between the two matrices Ψj and where the term EΨj Ω, calculated as the root-mean-squared of the residuals; therefore, the discriminatory power dj increases linearly with the inverse of the error, and CFj is a coverage factor that takes into account the occurrence of the m/z value j in the data set, namely, for how many samples the related chromatogram generates a signal above a threshold. Here thresholding is due to the curve resolution method applied in the initial step, which has also the effect of bringing down to zero signals of those chromatograms arising from m/z values that are at a very low level and that originate from pure noise. The coverage factor is introduced because the most useful m/z values are those that minimize the error but that simultaneously appear in a significant proportion of the data set, hence providing diagnostic information for the majority of the samples. In such circumstances, the calculation of CWIA can be extended to a larger fraction of the data set, and comparison between the similarity matrix Ψj and the target matrix Ω is maximized. CFj is defined as

CFj )

exp(τII‚(hj/n)) - 1 exp(τII) - 1

(14)

where hj is the number of hits for the m/z value j (number of times the related chromatogram generates a signal above the threshold) in the samples of the data set that consists of n samples; hence the ratio hj/n represents the actual coverage, and τII is a tuneable parameter. Equation 14 operates on the coverage in manner similar to eq 10 on the MCQs: an exponential weight is applied where τII determines its distribution as in Figure 4. In Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

1613

contrast to τI, here negative exponential weights should be preferred so that higher discriminatory power d is assigned to those m/z values that give partial coverage and can still be important for diagnosis, while still underweighting those giving minimum coverage. Alternatively to the use of CF, a presence/ absence criterion could be established, namely, a weight that takes into account when certain m/z values frequently appear in one class but not in another; however, this can be dangerous particularly when the drug product manufacturing process does not follow a strict standard: conditions may change over time and different impurities may arise in the same batch, which means that a certain m/z value can appear in one sample and not in a second one, even though both relate to the same batch/class. Also, this can be a problem when dealing with relatively small data sets where certain classes may be underrepresented and consequently results can be unstable. 4. Modeling and Validation. Because CWIA provides results in the form of a similarity matrix Ψ, a method to investigate the structure in the data is hierarchical agglomerative cluster analysis.7 Here objects successively merge into clusters of bigger size until all converge into a single one, and the whole process that can be graphically presented in the form of a dendrogram. Other clustering techniques such as Fuzzy C-Means (FCM)33 can also be used: in this method, the number of clusters is preselected, e.g., equal to the number of classes; samples are initially allocated randomly and redistribute at each iteration in that cluster whose center is the closest. If the features that describe the samples relate well to class membership (i.e., tablet origin), in each batch three classes should be well distinguished. However, better results are expected with a better subset of features. In fact, features used to cluster objects are often not all equally informative: some may be redundant or irrelevant for the clustering, and often many candidate features are included since the relevant features are unknown a priori.34,35 Generally, the identification of relevant features can save future measurement time and costs or can be used for chemical interpretation of the phenomenon, which in this case would enable relating m/z values to markers that in turn may be related to one of the manufacturing processes. The selection of the optimal subset of ions occurs by tuning a third parameter τIII that corresponds to a threshold for the discriminatory power d (eq 13), by choosing only those for which d is above τIII. Once the optimal subset of ions has been selected, three questions should be posed. (a) Does the optimal subset allow a better description of data? Better clustering should be expected on samples from which the ions were extracted (batch 1), e.g.. a dendrogram where the three classes separate earlier. Moreover, on single samples, the reduced TICs should display peaks characteristic of each group. (b) Is the optimal subset predictive? For samples from other batches (2 and 3), it should help predicting the origin of tablets, hence providing a form of validation. (c) Does the optimal subset relate to impurities? Findings of CWIA could be compared to other sources of information when available, e.g., comparison with the spectral profiles of known impurities. (33) Bezdek J. C. Pattern Recognition with Fuzzy Objective Function Algorithms; Plenum Press: New York, 1981. (34) Milligan, G. W. Psychometrika 1980, 45, 325-342. (35) Questier, F.; Walczak, B.; Massart, D. L.; Boucon, C.; de Jong, S. Anal. Chim. Acta 2002, 466, 311-324.

1614 Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

Table 2. Statistics of Fuzzy C-Means Clustering for Different Values of τIa τI

no. of errors

distance between cluster centers

-10 -9 -8 -7 -6 -5 -4 -3 -2 -1 1 2 3 4 5 6 7 8 9 10

8.52 8.27 8.11 8.11 8.11 8.00 8.00 8.00 8.00 2.03 2.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00

1.92 × 10-4 2.08 × 10-4 2.22 × 10-4 2.45 × 10-4 2.82 × 10-4 3.33 × 10-4 4.24 × 10-4 6.13 × 10-4 1.17 × 10-3 1.25 2.08 2.31 2.48 2.59 2.65 2.67 2.66 2.62 2.58 2.54

a The number of errors corresponds to samples misallocated using a model with three clusters (i.e., the three classes C, I, and G) and are expressed in fractions because results are averaged over 100 repetitions of the algorithm.

RESULTS AND DISCUSSION Clustering of Samples. The application of CWIA in the first phase of clustering requires the tuning of τI, which is the parameter responsible for distributing the weight among m/z values depending on the quality of the related chromatograms (eq 10), hence determining their relevance in the measure of correlation between pair of samples. Rationally, τI should be a positive exponential so that high-quality chromatograms have high weights while low-quality ones remain of low significance and weight but still contribute. Empirically, a method for obtaining the optimal τI can be that of performing clustering with varying values of this parameter in a broad range. Table 2 illustrates the results of FCM by preselecting three clusters (G, I, and C) using the CWIA matrix obtained from samples in batch 1 and using values for τI in the range -10 to +10. The results are averaged on 100 repetitions of the algorithm, because the outcome of FCM is partly influenced by the initial random allocation of samples into three clusters. Batch 1 consists of 9 samples from group G, 8 samples from group I, and 9 samples from group C. For negative values of τI it can be observed that the performance of FCM is poor. Until τI ) -2, samples G and I are indistinguishable: all 8 samples from group I fall into G cluster, plus occasionally one from group C, from which the average error results between 8 and 9. In the range from τI ) -1 to τI ) +2, some error still occurs but is limited to one or two samples, while τI ) +3 is the value at which errors are no longer observed. Another useful diagnostic is the distance between the centers of the clusters: this increases after τI ) +3 to reach a maximum at τI ) +6. Therefore, a sensible choice for the optimal τI might be a value in the range +3 to +6. In this application, we selected the lower limit τI ) +3 because a lower value would distribute weights more evenly on m/z values whose chromatograms are of medium quality (MCQ between 0.5

Figure 5. (a) Dendrogram generated using the q index (eq 8) between reduced TICs after the application of CODA with threshold 0.8. (b) Dendrogram generated using the CWIA index with τI ) +3. In the latter case, all replicates cluster at earliest stages and groups emerge according to the origin of tablets. The method of “single linkage” was used as aggregation criteria between clusters.

and 0.7) where impurities at lower concentration may occur. Figure 5 compares the dendrograms on the samples that compose batch 1 using the similarity matrix generated from the reduced TICs method of eq 8 and the similarity matrix from CWIA of eq 9 with τI ) +3. By comparing replicates of the series I5 and C2, it can be seen that the first method generates rather unstable results as I5a, I5b, and C2a do not cluster together with other replicates. These samples are in fact expected to cluster together at very early stages because they are identical, and differences should relate only to instrumental and tablet-to-tablet variations within the same lot. There also appear to be some clusters with intercalation between samples of different origin: C4 and C5 group together with class I, which in turn is split in two subclusters, and G9 clusters far from all samples of its true class. The use of CWIA improves performance in both aspects. Replicates cluster

together, and the three classes relating to the origin of tablets can be readily distinguished. Some variability can still be observed within a class; for example, samples I1 and I2 cluster later than the other I samples as well as C4 and C5 compared to samples from class C. Sample Un appears to belong to class G as it associates well within this cluster, in the center. Independent information other than from LC/MS analysis of the samples suggested that sample Un originated from manufacturer G; thus, the methodology described here assigned sample Un to the correct class. In principle, better results are expected because CWIA better exploits the two-dimensional information in the data matrix, as it considers all chromatograms and information is not lost along the mass direction. However, at this stage CWIA does not provide information about which are the features in the data responsible for such clustering. Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

1615

Table 3. Number of Ions Selected as Optimal Subset and Minimum Coverage (min h/n%) Occurring in Such Optimal Sub-selection for Different Values of τII and τIIIa 0.40

τIII

0.45

0.50

0.55

0.60

0.65

0.70

0.75

τII

no ions

min h/n%

no ions

min h/n%

no ions

min h/n%

no ions

min h/n%

no ions

min h/n%

no ions

min h/n%

no ions

min h/n%

no ions

min h/n%

-10 -9 -8 -7 -6 -5 -4 -3 -2 -1 1 2 3 4 5 6 7 8 9 10

659 650 635 614 596 549 522 481 445 391 312 287 271 260 249 246 243 239 237 233

14.8 18.5 18.5 22.2 22.2 29.6 37.0 44.4 48.1 53.3 77.8 81.5 81.5 92.6 92.6 92.6 92.6 95.6 96.3 96.3

523 509 481 458 437 412 388 353 301 261 217 207 204 200 198 195 191 190 188 188

17.4 22.2 22.2 29.6 29.6 40.7 44.4 48.1 53.3 74.1 77.8 81.5 92.6 92.6 92.6 96.3 96.3 96.3 96.3 96.3

170 165 159 152 140 132 122 112 103 98 83 81 80 75 73 71 71 71 71 71

22.2 29.6 44.4 44.4 44.4 44.4 48.1 53.3 74.1 74.1 81.5 92.6 92.6 92.6 96.3 96.3 96.3 96.3 96.3 96.3

67 67 65 64 60 59 55 54 48 47 38 36 34 33 33 33 33 32 31 31

44.4 44.4 44.4 48.1 53.3 53.3 53.3 74.1 77.8 77.8 81.5 92.6 92.6 96.3 96.3 96.3 96.3 96.3 99.2 99.2

30 30 30 30 30 28 27 27 26 23 21 20 20 20 20 18 18 18 18 18

53.3 53.3 53.3 53.3 53.3 53.3 77.8 77.8 77.8 81.5 92.6 96.3 96.3 96.3 96.3 100 100 100 100 100

16 16 16 16 16 16 15 15 15 15 13 13 13 11 10 10 10 10 10 10

77.8 77.8 77.8 77.8 77.8 77.8 81.5 81.5 81.5 81.5 96.3 96.3 96.3 100 100 100 100 100 100 100

8 8 8 8 8 8 7 7 7 6 5 4 4 4 4 4 4 4 4 4

81.5 81.5 81.5 81.5 81.5 81.5 81.5 81.5 81.5 96.3 96.3 100 100 100 100 100 100 100 100 100

4 4 4 4 4 4 3 3 3 3 3 3 3 3 3 3 3 3 3 3

81.5 81.5 81.5 81.5 81.5 81.5 100 100 100 100 100 100 100 100 100 100 100 100 100 100

a A decrease of τ determines the selection of more ions, while a decrease of τ increases the minimum coverage occurring in the subselection. III II The optimal for τII and τIII can be defined by choosing the preferred size and minimum coverage. For example, by imposing no ions < 30 and Min h/n% > 50, 28 ions are selected with a minimum coverage equal to 53.3% (boldface type), which in turn correspond to τII ) -5 and τIII ) 0.6.

Table 4. Distribution of the Ions Selected as Optimal Subset between the Three Classes G, I, and C, Using Samples from Batch 1, with CWIA Parameter Values τI ) +3,τII ) -5, and τIII ) 0.6a G

I

C

rank

m/z

dj

EΨjΩ

hj/n%

m/z

dj

EΨjΩ

hj/n%

m/z

dj

EΨjΩ

hj/n%

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

429 431 299 434 433 301 342 428 432 430 344 314 407 382 452 405 279 300 316 409

0.785 0.771 0.753 0.751 0.746 0.740 0.665 0.651 0.651 0.609 0.595 0.572 0.563 0.549 0.546 0.540 0.535 0.533 0.533 0.531

0.215 0.229 0.247 0.241 0.253 0.260 0.335 0.349 0.349 0.391 0.405 0.427 0.420 0.451 0.451 0.453 0.465 0.467 0.463 0.421

100 100 100 81.5 96.3 100 100 100 100 100 100 96.3 66.7 100 88.9 77.8 100 100 85.2 48.1

330 270 272 332 271 273 434 433 346 344 640 318 315 639 314 316 343 525 392 286

0.675 0.661 0.659 0.639 0.636 0.617 0.613 0.613 0.611 0.607 0.605 0.597 0.589 0.583 0.577 0.575 0.574 0.571 0.570 0.568

0.325 0.339 0.341 0.354 0.357 0.383 0.380 0.386 0.389 0.393 0.354 0.394 0.400 0.416 0.422 0.421 0.426 0.428 0.429 0.432

100 100 100 81.5 81.5 96.3 81.5 96.3 100 100 53.3 77.8 74.1 95.6 96.3 85.2 100 96.3 96.3 100

402 299 404 301 415 398 342 330 300 302 344 328 399 470 429 403 427 517 445 469

0.731 0.720 0.700 0.694 0.672 0.655 0.648 0.626 0.619 0.616 0.612 0.612 0.608 0.607 0.597 0.595 0.592 0.589 0.583 0.583

0.268 0.280 0.298 0.306 0.327 0.336 0.352 0.374 0.381 0.384 0.388 0.388 0.392 0.392 0.403 0.394 0.405 0.411 0.417 0.416

96.3 100 92.6 100 99.2 77.8 100 100 100 100 100 100 100 96.3 100 74.1 88.9 100 100 96.3

a According to eqs 13 and 14, d inversely relates to the error E j ΨjΩ and directly relates in an exponential manner to the coverage hj/n. Those ions for which the discriminant power dj is above the threshold τIII ) 0.6, are selected as optimal subset (boldface type).

Detection of Diagnostic Ions. The most useful ions for discriminating between the classes of interest can be detected by providing the algorithm with information on the origin of tablets. Once this information is used, the relevant m/z values can be identified by comparison with the target similarity matrix Ω (eq 11, via eq 12). This phase of characterization requires the tuning of the remaining two parameters: τII that distributes the weight among m/z values depending on the coverage of the related 1616 Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

chromatograms in the data set (eq 14); τIII that corresponds to the threshold for the discriminatory index (eq 13) to accept/reject m/z values in the optimal subset of features. Table 3 presents results for τII in the range from -10 to +10 and τIII in the range from +0.40 to +0.75, reporting the number of m/z values selected as top ranked features, and the minimum coverage percentage found in the m/z values that compose this subselection. It can be noted that m/z values with lower coverage fall progressively into the optimal subset as τII decreases, while a decrease of τIII causes

Figure 6. Resemblance between the target similarity matrix Ω and the CWIA similarity matrix Ψ for different m/z values. Ψ402 and Ψ433 indicate these ions to be diagnostic, Ψ344 to be limitedly diagnostic, Ψ259 to be irrelevant, and Ψ448 to be potentially relevant but not diagnostic because of its limited coverage.

a selection of more m/z values because the demand in terms of discriminatory power is lessened. The optimal tuning for both τII and τIII can be determined simultaneously in a simple manner by the analyst: for example, if we are interested in looking at less than 30 chromatograms, all of which characterize at least half of the data set, we implicitly select τII ) -5 and τIII ) 0.6 that in Table 3 correspond to 28 chromatograms and a minimum coverage of 53.3%. Table 4 shows how the ions of the selected chromatograms distribute between the classes, ranked according to their discriminatory power dj (eq 13): 10 characterize class G, 11 characterize class I, and 14 characterize class C, with 7 repetitions, which makes 10 + 11 + 14 - 7 ) 28 in total. Notice that those ions minimizing the error can be lower ranked because of their limited coverage. An example is m/z ) 409 for class G that produces a relatively low error of 0.421 but is relegated from 14th position in the ranking to the 20th because it covers less than a half of the data set (48.1%). However, the value τII ) -5 (Figure 4) guarantees substantial redistribution only when the coverage is very limited. For example m/z ) 434 with coverage of 81.5%, loses only one

position in the ranking for classes G and I, while m/z ) 398 with coverage of 77.8% maintains its original position for class C. A major strength of CWIA is that it allows visualizing results in matrix form, hence providing an insight of the relevance of a certain ion for characterizing entire data sets rather than single samples, upon a comparison with the target similarity matrix Ω shown in Figure 6. Here, each pixel corresponds to the CWIA value between pairs of samples: the brighter, the higher the similarity. In the target matrix Ω pixels are either white or black as the “target” corresponds to a maximum within groups similarity (white ) 1) and minimum between groups similarity (black ) 0), respectively. Corresponding characteristic peaks can then be observed in detail on the profiles of single samples (Figure 7). These figures show some representative examples: m/z ) 402 is highly diagnostic for group C because Ψ402 reproduces the feature of Ω in the bottom-right corner. The respective chromatograms on single samples show that, at m/z ) 402, a genuinely intense peak occurs for class C with no peaks for classes G and I. m/z ) 433 is highly diagnostic for class G and partly diagnostic for class I as Ψ433 reproduces Ω well in the top-left corner but Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

1617

Figure 7. Ion chromatograms of samples from each class. m/z ) 402 is clearly diagnostic for class C, m/z ) 433 is relevant for class G and to a lesser extent for I, m/z ) 344 is partly diagnostic for I, while m/z ) 259 does not provide any information for discriminating between the three classes.

only partly in the center. Translated into chromatograms, m/z ) 433 is present as a well-defined peak for all the samples in class G, with another peak for some of the samples in class I but at different elution time, therefore uncorrelated with the former one. m/z ) 344 shows a different case: most of the features of Ω are represented in Ψ344, but the difference between classes is less evident. Particularly, classes G and C show some correlation as unwanted features, in the top-right and bottom-left corners. By examining the chromatograms, it can be seen that in all cases the classes display some peaks along the profile with one at elution 1618 Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

time around 20 min, common to all samples, another around 14 min that is common to G and C and two different ones for I. Overall m/z ) 344 is diagnostic for class I while being substantially unable to distinguish between classes G and C. Finally, m/z ) 259 is an example of a nondiagnostic ion as its map Ψ259 does not reproduce Ω with the corresponding m/z value that are inconclusive, while m/z ) 448 is a case where some correlation appears but overall it is not indicative because the coverage is small. In fact m/z ) 448 generates chromatograms for only few of the samples in the data set, hence being scarcely diagnostic.

Figure 8. (a) Comparison between the original TICs and the reduced TICs using the chromatograms of the optimal subset of ions selected with CWIA, on three samples from batch 1 representative of each class. A series of peaks relating to the classes can be determined. (b) Dendrograms using CWIA with the optimal subset of ions: the clusters emerge earlier and more clearly than in Figure 5b where all chromatograms were used.

Prediction of Unknown Samples and Consistency. A reduced model including only the optimal subset of selected ions is expected to enhance and highlight differences between the three classes G, I, and C. With respect to the samples that compose batch 1, it is more appropriate to indicate this property in terms of improved description rather than prediction, because batch 1 acts as training set on which information on class membership is provided to guide selection and determine the optimal subset of ions. Figure 8a on the right shows the reduced TIC profiles on three samples from batch 1 for each of the three classes C, I, and G, by using the chromatograms of the top-ranked 28 ions as in Table 4 with CWIA parameters set to τI ) +3, τI )

-5, and τIII ) 0.6. In comparison to the original TICs (on the left), a number of characteristic peaks can now be clearly determined in relation to each class. The dendrogram in Figure 8b demonstrates that such improvement applies to the entire data set. Clusters associate according to the origin of tablets earlier rather than what was observed in Figure 5, where all chromatograms where taken into account as input for CWIA. The use of this optimal subset of ions results in a model that shows also prediction ability. Those samples in batches 2 and 3 for which information on class membership was not provided for the selection of the optimal subset (i.e., samples G10, G11, G12, G13, G14, I3, and C6) group according to the correct origin Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

1619

Figure 9. Clustering of batches 2 and 3 using CWIA with the chromatograms of the optimal subset of ions selected from batch 1. Samples that are unique to batches 2 and 3 (underlined) cluster according to their class.

(Figure 9). This provides a form of model validation and further indicates that the method can be used to predict the origin of unknowns. Finally, the consistency of the findings using CWIA can be verified if independent information on the spectral features of impurities is available from other sources. In this study, the active pharmaceutical ingredient analyzed contains two chlorine atoms and therefore related impurities will show the distinctive isotope pattern for chlorine. Figure 10 shows the mass spectrum with its distinctive isotope pattern (dichloro compound) obtained for the compound eluting at 20.9 min, which corresponds to the peak observed in the chromatogram for m/z ) 402 (Figure 7). Hence m/z ) 402, which is highly characteristic for class C is also associated with m/z ) 403 and 404. It should be noted that the way Table 4 presents results generated by the CWIA method 1620 Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

Figure 10. Mass spectrum obtained from the chromatographic peak eluting at 20.9 min identified in the ion chromatogram corresponding to m/z ) 402 in Figure 7 for a sample of class C.

can therefore be redundant to some extent, because the m/z values identified do not correspond to any specific compound, since in the mass spectra each compound will result in several

significant ions. Interestingly, not all the ions of the isotope cluster related to this compound appear in Table 4. It is possible that this is because some are not diagnostic for one class alone but also arise from different compounds at different elution times that characterize other classes. For example, m/z ) 405 (which is part of the isotope pattern) is not only diagnostic for class C, as it appears in Table 4 in the list of the characteristic m/z values for class G. Similarly, the diagnostic m/z values 428, 429, 430, 431, 432, 433, and 434 for class G are from an isotope cluster from a single tetrachloro compound. However, this type of interpretation requires knowledge from independent sources that in many cases is not available and therefore is not implemented in the CWIA method. CONCLUSIONS In certain applications, pattern recognition on LC/MS data is a daunting task because the information of interest (route/process specific impurities) is frequently confounded with information that is not relevant for the target study (e.g., the main components) and with a variable quantity of noise both of chemical and instrumental origin. The situation is further complicated by the fact that, over time, some irreproducibility may occur due to variation in instrumental conditions and actual changes in tablet composition. Chemometrics provides a variety of tools to assist the analysis: curve resolution helps in removing the irrelevant information arising from the dominant compounds in the formulation, and the component detection algorithm is effective for denoising in the sense that it isolates the remaining chemical information from noise. However, for recognition of entire data sets, there are limitations due to the fact that CODA acts on individual samples, where a different part of the data (both the number and type of chromatograms) can be selected each time, and the most

appropriate threshold to extract this fraction can change due to variations of the system, for example, when minor impurities at noise level only appear in a portion of the samples. Nevertheless, determination of characteristic ions useful for discrimination can successfully be obtained on a few samples by considering for example the union of the selected m/z values.13 The CWIA method proposed in this paper overcomes these limitations for the analysis of larger data sets. By removing the need of a threshold when using CODA, CWIA fully exploits the two-dimensional information in LC/MS data, taking into account for similarity recognition between samples all chromatograms on the basis of their quality (eqs 9 and 10). The major groups can be visualized using clustering techniques, and the features responsible for such associations can be characterized with respect to the entire data set in the form of a matrix (Figure 6). Finally, CWIA allows the identification of an optimal subset of ions that can highlight differences between the classes in the form or reduced TICs (Figure 8a) and that can be used for prediction of unknowns (i.e., origin of tablets). The procedure is only semiautomated because it requires setting of three parameters, but it may be suitable for other applications involving exploratory studies of larger data sets when prior information on the data is scarce or absent, hence requiring the acquisition of chromatograms in a broad m/z range. ACKNOWLEDGMENT S.Z. gratefully acknowledges GlaxoSmithKline and EPRSC for financial support.

Received for review December 23, 2004.

October

10,

2004.

Accepted

AC048504T

Analytical Chemistry, Vol. 77, No. 6, March 15, 2005

1621