Classification and Diagnosis of Bioprocess Cell Growth Productions

Jun 22, 2019 - (7−10) For example, principal component analysis (PCA) and partial least-squares ... (13) Related research work and variants were als...
0 downloads 0 Views 3MB Size
Article pubs.acs.org/IECR

Cite This: Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Classification and Diagnosis of Bioprocess Cell Growth Productions Using Early-Stage Data Yuan Jin,† S. Joe Qin,*,† Qiang Huang,‡ Victor Saucedo,¶ Zheng Li,¶ Angela Meier,§ Siddhartha Kundu,∥ Bri Lehr,⊥ and Salim Charaniya□

Downloaded via GUILFORD COLG on July 31, 2019 at 22:52:02 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.



Mork Family Department of Chemical Engineering and Materials Science, University of Southern California, Los Angeles, California 90089, United States ‡ Department of Industrial and Systems Engineering, University of Southern California, Los Angeles, California 90089 United States ¶ Process Development Engineering, Genentech, South San Francisco, California 94080, United States § Late Stage Cell Culture, Genentech, South San Francisco, California 94080, United States ∥ MSAT, Genentech, Vacaville, California 95688, United States ⊥ MSAT, Genentech, Oceanside, California 92056, United States □ Global MSAT, Genentech, South San Francisco, California 94080, United States ABSTRACT: Many industrial pharmaceutical manufacturing processes are composed of multiple-step batch operations. However, uncontrolled variations often occur during operations that affect the cell growth performance. Because of the complexity of biological processes, one leading challenge in process operation is the identification of potential causes of undesirable process variabilities. In this paper, we propose a classification and diagnosis strategy to analyze cell culture manufacturing variability in bioprocesses with the objective of unveiling hidden factors affecting process yield and performance. The proposed strategy includes two parts: (i) a clustering method performed in the principal component residuals that effectively separates the low lactate batches into different clusters and (ii) a fault diagnosis method based on regularized LDA contribution analysis for exploring the leading contributors to each of the low performance classes. The proposed strategy is applied to industrial production data collected over 8 years from a batch process. The effectiveness of the proposed approach is demonstrated on this data set, allowing the classification and diagnosis of the sources of the low performance situations.



INTRODUCTION Multiple-step batch processes are widely used in biopharmaceutical manufacturing in the utilization of mammalian cells for the production of recombinant protein therapeutics. Because of the lack of process observability and a means to control the disturbances, process variability is often experienced in these complex bioprocesses.1,2 Data of enormous volume and variety are generated and collected from different cell culture process stages, but they are not fully utilized in practice. Therefore, mining the historical data is pivotal in uncovering the reasons that cause the process variability and diagnosing them in early stages for improved process performance.1,3 Low titers in CHO monoclonal antibody (mAb) production are often correlated with reduced cell growth and with high lactate production.1 Several studies have employed lactate metabolism models to study this phenomenon.4−6 For example, a steady state topological metabolic model has been proposed,5 indicating the mechanism of lactate consumption and accumulation that explain the cause of high lactate. Nevertheless, these first principles based models are built with lactate metabolic pathway simulation. It is difficult to apply them to large-scale manufacturing processes because the © 2019 American Chemical Society

concentrations of additional metabolites are needed but are not routinely available in commercial manufacturing. Batch process monitoring is widely applied in the chemical, semiconductor, and pharmaceutical industries and is useful for ensuring product quality and monitoring process operation anomalies. In a vast body of literature, data-driven modeling is a popular approach to exploring knowledge from bioprocess operation data. In particular, process data analytics (PDA) has been successfully applied to process monitoring and fault diagnosis for bioprocesses.7−10 For example, principal component analysis (PCA) and partial least-squares (PLS) algorithms have been used to model variable correlations using historical data. Specifically, multiway PCA (mPCA) and multiway PLS (mPLS) were used to deal with batch process data. The key idea behind mPCA and mPLS is to unfold the Special Issue: Dominique Bonvin Festschrift Received: Revised: Accepted: Published: 13469

March 6, 2019 May 30, 2019 June 22, 2019 June 22, 2019 DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Article

Industrial & Engineering Chemistry Research

discriminant analysis (LDA) for further diagnosis. This strategy includes the following two steps:

three-way batch data tensor into a two-way matrix and apply PCA or PLS on the unfolded data. Nomikos et al. applied mPCA and mPLS to extract the information in the batch trajectory data and developed simple monitoring charts for process fault detection.11,12 Lee et al. used multiway kernel PCA to capture the nonlinear characteristics within normal batch processes and detect faults.13 Related research work and variants were also reported.14,15 To take dynamics into account, batch dynamic PCA (BDPCA) and batch dynamic PLS (BDPLS) were developed.16 Another approach to process monitoring is to use clustering methods to group different process runs into subgroups according to their similarity on the basis of a distance measure.1,17 For example, a clustering technique was used to isolate the variable contributors that were faulty by calculating the Euclidean distance between pairs of contributors.17 Singhal et al. developed a modified k-means methodology on the basis of PCA to cluster multivariate time-series data using similarity factors.18 Husson et al. combined PCA, hierarchical clustering, and partitional clustering to enrich the description of data and obtain better visualization.19 Classification based methods can be applied when the fault data with class information is available.14,20−22 Support vector data description (SVDD) is proposed20 and applied to monitor the batch process inspired by support vector machines (SVM).21,22 This method has no Gaussian assumption of the process data and is effective for nonlinear data modeling. Furthermore, Yao et al. proposed functional SVDD (FSVDD) based batch process monitoring, which incorporated functional data analysis (FDA).23 Le et al. analyzed batch process data from a large-scale pharmaceutical manufacturing process.3 In their study, both PLS and SVM were employed to predict the final process outcome. With known fault types such as bias and variance increases, Wang and Zhao24 proposed a nested-loop Fisher discriminant analysis algorithm combined with a measure of relative changes to discriminate among them. With the advent of big data analytics and the successful use of machine learning in other industries, one motivation of this paper is to expand the volume and variety of data for the diagnosis and classification of critical biopharmaceutical manufacturing anomalies. Through a collaboration between university and industrial practitioners, we studied the classification and diagnosis of a type of industrial bioprocess on the basis of 8 years of off-line measured production quality related data. Although there are plenty of process measurement data, these off-line measurement data are critical, sparse, and expensive to measure. Therefore, they are valuable for discovering the contributing features of the uncontrolled variations that lead to poor performance. In this paper, we propose an integrated data analytics strategy to analyze data from the last step of the bioprocess manufacturing chain. On the basis of data collected from 241 batches over 8 years, we first classified the batches into high performance (low lactate) and low performance (high lactate) ones. Then, we used this classification information to select relevant features that best separate these two classes using PCA. Because the batch duration is over 11 days and the overall batch performance is largely determined in the first 6 days into the batch, we extracted features in the first 144 h (6 days), which we used to predict the final batch lactate levels. In this feature space, we performed clustering of the low performing batches and proposed the use of regularized linear

• the use of hierarchical clustering on the useful features found by PCA to cluster the low performing batches into different groups • the use of regularized LDA to visualize and diagnose the major contributors to each of the high lactate batch clusters The organization of the paper is as follows. Following the introduction, we present the description of the large-scale cell culture process and the problem to be studied. Then, we introduce the proposed strategy of the data analytics approaches for bioprocess variability source analysis. Case studies on data collected over 8 years from an industrial manufacturing bioprocess are then presented to demonstrate the effectiveness of the proposed analytics strategy. Conclusions and discussions are given in the last section.



BIOPROCESS AND PROBLEM DESCRIPTION Cell culture is a process by which cells are grown under controlled conditions. The layout of a large-scale cell culture process developed by Genentech is well discussed in Li et al.,25 which shows the large-scale cell culture processing steps, consisting of the cell bank, seed train, inoculum train, and final production run. The production run is the last step, N, with a manufacturing scale of 12 000 L per run. Every step has its key process performance indicators (KPI), such as cell growth rate, final viability, and titer, along with appropriate product quality attributes. The success of process scale up is usually measured by these KPIs and product quality attributes to meet predefined criteria.25 In this study, data from 8 years of historical batch runs of a large-scale cell culture process from Genentech were collected. Product quality of bulk drug substance of all selected batches met release specifications. Process data were obtained from the inoculum train, which was cultured for approximately 4 days with one measurement per day. The production run, which is the focus of this study, was cultured for approximately 11 days with two measurements per day. For each stage, the data include 11 off-line measurements, 11 initial batch conditions, and 7 media preparation measurements (e.g., media preparation osmolality and media preparation duration). These variables are highly correlated and mainly consist of physiological parameters, such as viable cell density (VCD), viability, and packed cell volume (PCV), and physio-chemical environment variables, such as buffer pH and osmotic pressure. In this long processing time, the pH value is adjusted in realtime by the control system. Midcourse corrections can be made around the 72 h mark. At this point, extra glucose and nutrition can be added to the batch on the basis of standard operating procedures. Any corrections attempted after 144 h will have virtually no effect on the batch. Further, the production process must run to finish even though high lactate or poor performance is anticipated at an early stage. As an important performance metric, the lactate concentration at the end of the production cell culture process was measured and scaled to [0, 1]. The distribution of the scaled final lactate concentration is shown in Figure 1. It is observed that the distribution of the lactate concentration is bimodal, with a larger mode for the low lactate runs and a smaller mode for the high lactate runs. On the basis of the process knowledge, if the normalized final lactate value of a run is 13470

DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Article

Industrial & Engineering Chemistry Research

paper, data from 8 years of manufacturing needed to be classified into various types of low performing situations, and their contributing features needed to be diagnosed. Therefore, we designed the analytics method to achieve the following tasks: (i) selection of pertinent features for differentiation between low and high lactate batches; (ii) application of hierarchical clustering (HC) to the selected features that are sensitive to faults, essentially modifying the dissimilarity measures for HC and categorizing the high lactate batches; (iii) application of LDA for visual confirmation of the clustering of all high lactate batches; and (iv) development of regularized LDA based contribution analysis between normal and high lactate clusters to diagnose the contributing features. The flowchart of the proposed strategy is introduced at the end of this section. PCA Based Feature Selection and Monitoring. PCA is a popular and basic method of unsupervised multivariate analytics. Given the scaled process matrix X ∈ n × m consisting of m variables and n samples, PCA aims to maximize the variance of the projected scores. With a number of principal components, PCA decomposes the data into a principal portion and a residual portion as follows:

Figure 1. Scaled final lactate concentration distribution.

larger than 0.4, this run is considered a high lactate run, which is classified as poor performance. Otherwise, it is considered a low lactate run. It is further observed that the lactate near the end of the batch run fluctuates, especially for high lactate runs. To refine the classification, the lactate concentration at the end of the batch and the average of the final three lactate measurements were used for classification. The results are shown in Figure 2. If the final measurement of a run was higher than 0.3, and the average of the final three measurements was greater than 0.25, this run was recognized as a high lactate run (red curves). Similarly, if the final measurement of a run was smaller than 0.2, and the average of the final three measurements was less than 0.18, it was recognized as a low lactate run (blue curves). The rest of the runs were in-between trace runs (green curves). One goal of this study was to detect and diagnose the potential causes for high lactate runs in the first 144 h, so that effective corrections could be implemented in future operations to control the production cell culture performance. Work on monitoring on the basis of multiphase partition of batch processes has been reported,26 but the scope included analysis of the entire batch. The work in this paper focuses on the early stage (first 144 h) of the batch process to diagnose early symptoms of batches with various unknown causes that lead to low-value products.

̃ ̃T X = TPT + TP

(1)

where T and P are the principal score matrix and loading matrix, respectively. T̃ and P̃ are the residual score matrix and loading matrix, respectively. The two portions of scores are orthogonal and so are the two portions of loadings, which leads to P̃ P̃ T = I − PPT. For a sample, xi, its PCA model estimate, x̂i, and residual, x̃i, are x̂ i = PPTx i ̃ ̃ Tx i = (I − PPT)x i x̃ i = PP

(2)

PCA based fault detection assumes that the statistical properties of the principal component (PC) scores and residuals are within a control limit for all normal operations. The principal components and residuals are monitored with different statistics;29,30 that is, Hotelling’s T2 is used to monitor the PC scores, and SPE is used to monitor the residuals. The T2 index and SPE index for the ith sample are defined as



PROPOSED ANALYTICS FOR CLASSIFICATION AND DIAGNOSIS In this section, we briefly discuss the proposed methods for the batch classification and diagnosis analytics of historical biopharmaceutical manufacturing operation data. Unlike the existing batch process monitoring work in the literature,27,28 which is based on a PCA model of the normal class and in which the T2 and squared prediction error (SPE) indices are used for the detection and diagnosis of anomalies, in this

Ti 2 = x i TPΛ−1PTx i SPEi = x i T(I − PPT)x i

(3)

where Λ represents the covariance of the principal scores, T. The control limits with (1 − α) confidence are31,32

Figure 2. Batch lactate profiles of the classified high and low lactate runs (left) and unclassified intermediate lactate profiles (right). 13471

DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Article

Industrial & Engineering Chemistry Research

Ti 2 ≤ χl ; α 2 SPEi ≤ gχh; α 2

First, it is an unsupervised approach; it treats all included variables with equal importance without the ability to select critical ones. If irrelevant variables are included in HC, they can distort the clustering results significantly. Second, it is scale dependent. In this paper, we use PCA to select relevant features on the basis of their ability to separate low lactate batches from high lactate ones and then apply HC to categorize the high lactate batches into appropriate clusters. From eq 5, if we use the PCA residuals as the features for further clustering, we see the connection between residuals based WSS and SPE for the normal batches as follows:

(4)

where coefficient g and degrees of freedom h are calculated by m

g=

∑l + 1 λi 2 m ∑l + 1 λi

m

and h =

2

(∑l + 1 λi) m ∑l + 1 λi 2

.

It is commonly accepted that the principal components are more likely to be able to separate normal data from abnormal situations because the principal component scores represent the largest variability in the data. For the process data studied in this paper, however, we will show the opposite; that is, the variability in the principal components is least useful for separating normal batches from abnormal ones; it is the residuals that are effective in the discrimination between them. Thus, PCA residuals will be used in the subsequent hierarchical clustering, which will be discussed next. Hierarchical Clustering. The goal of data clustering is to group similar data points into the same group, so that data points in the same group have similar properties. Hierarchical clustering (HC) is an unsupervised learning approach that does not require class information nor a prespecified number of clusters.33 The hierarchy of clusters is represented as a dendrogram with the height of the dendrogram corresponding to the dissimilarity of the clusters.33,34 In our proposed strategy, the dissimilarity is defined as the Euclidean distance with complete linkage because it tends to yield a more balanced dendrogram. The total within-cluster sum of square (WSS) for cluster Ck is defined as follows: WSS =

∑ x ∈ Ck



WSS0 =

x i ∈ C0

x̃ i

2

=



SPEi (6)

x i ∈ C0

This relation helps us see the connection between SPE for the normal batches and the corresponding WSS measure. LDA and Regularized LDA Based Contributions. Although both LDA and PCA are commonly used linear transformation techniques, the major difference is that PCA is an unsupervised learning algorithm that does not need data labels. In contrast, LDA is a supervised learning algorithm that seeks to maximize the separation of different known classes. Let Xk represent the set of samples in class k; the scatter matrix for class k is Sk =



(x − x̅k )(x − x̅k )T (7)

x ∈ Ck

where x̅k is the sample mean of the kth cluster. The within-class scatter matrix is

∥x − μk ∥2

K

(5)

Sw =

where μk is the mean vector of all data samples in cluster Ck. Although hierarchical clustering seems to be a well-defined approach, it can lead to meaningless results if used carelessly.

∑ Sk

(8)

k=1

and the between-class scatter matrix is 13472

DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Article

Industrial & Engineering Chemistry Research K

Sb =

∑ nk(x̅k − x̅)(x̅k − x̅)T k=1

(9)

where x̅ is the overall mean of the data, and K is the total number of classes. The direction along which to best separate the different classes can be found by maximizing the ratio of the between-class scatteredness versus the within-class scatteredness: J(ϕ) =

ϕ T S bϕ ϕTSw ϕ

(10)

The optimal direction is obtained by solving the generalized eigenvalue problem: S bϕ = λSw ϕ

(11)

When there exists strong collinearity in the within-class data, Sw will be ill-conditioned; thus, the optimal direction found by LDA is very sensitive to noise in the data. In order to avoid this issue, regularized LDA is proposed by adding a regularization term to Sw:35 Sw = γ Sw + (1 − γ )I

γ ∈ [0, 1]

Figure 3. Flowchart of the proposed classification and diagnosis strategy.



(12)

CASE STUDY: CLASSIFICATION AND DIAGNOSIS OF THE BIOPROCESS DATA In this section, we present the results of applying the proposed strategy to large-scale manufacturing bioprocess data. We start from data preprocessing and then discuss how the two objectives in our strategy are achieved for the industrial process data classification and diagnosis. The advantages of the proposed approach are shown by comparison with other alternatives. Data Preprocessing. The batch process data is a threedimensional array that is organized into a cube, X(N × M × K), as shown in Figure 4. After preprocessing for all batches collected over 8 years, the final data set contains 241 batch runs, where each of them has 10 off-line measured process variables. Each batch is sampled twice a day over the batch length of 11 days. Although 23 samples are collected for each batch, we focus on the analysis of the first 13 samples, which spread over 144 h. The variabilities during this period largely determine the final batch performance, because subsequent changes to the process after this stage do not alter the outcome. To analyze the batch-wise data, we unfold X into a matrix (N × (M × K)), following Nomikos et al.12,37 The new matrix has dimensions of 241 × 130, with each column representing a physical variable at a certain time. Because of the low sampling frequency, process dynamics is not a consideration in this paper. PCA for Feature Selection. We apply PCA on the first 144 h of the data for all the low lactate batches shown in Figure 2 and check against the high lactate batches to find the best features to distinguish the two classes. Figure 5 shows the cumulative proportion of variance (CPV) versus the number of PCs for the normal data. As typical for unfolded batch data, CPV does not show a sharp change in the rate of increments. Therefore, the CPV method is inconclusive for this work. To examine which part of the PCA model is more effective to separate the high lactate batches from the low lactate ones, we show in Figure 6 the T2 and SPE charts for 70, 80, 85, and 90% CPV captured. The corresponding 95% T2 control limits

where the hyperparameter γ is selected on the basis of crossvalidation. Note that if γ = 1, this is the standard LDA; if γ = 0, then eq 10 reduces to the objective of PCA. The regularized LDA method can be used for fault diagnosis to characterize the variable contributions to a fault cluster.36 With a class of normal process data, X0, and a class of fault data, Xk, for this biclass LDA problem, we only need to find one optimal direction to discriminate fault cluster Xk from normal cluster X0. This direction is defined as the fault direction for Xk, and the weights of variables projected onto this direction are defined as the variable contributions to fault class k. This approach extends to the regularized LDA based contribution analysis as well, which is shown in Algorithm 1. A comparison of LDA based contribution analysis and regularized LDA based contribution analysis is discussed in the next section. Proposed Flowchart for Classification and Diagnosis. The flowchart of the proposed strategy to analyze the variability for the bioprocess is shown in Figure 3. There are two major objectives in this strategy: (i) clustering of the low performance batches into appropriate groups with similar characteristics based on carefully selected features and (ii) diagnosis of the contributing features of the high lactate clusters. We analyzed early stage batch data focusing on the first 144 h of the bioprocess to discover sensitive process features that best separated low lactate batches from high lactate ones. We first unfolded the three-way data into a two-way matrix and then constructed a PCA model for the normal batches, which were the low lactate batches. Then we selected features that best characterized the separation of normal and high lactate batches. After that, we perform hierarchical clustering on the selected features (i.e., the residuals in this case) to discover groups of the uncontrolled process variations. The results of how different clusters are separated are visualized with LDA directions. Further, regularized LDA based contribution analysis is applied to diagnose the contributing features for each high lactate cluster. 13473

DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Article

Industrial & Engineering Chemistry Research

Figure 4. Three-way process data unfolding.

substance product qualities of all these batches met registered product specifications. From Figure 6, we observe that in the principal component subspace, the high lactate batches and the low lactate batches show little difference in terms of T2 values, whereas using SPE in the residual subspace, the high lactate batches and the low lactate batches are easily separated. Furthermore, CPV of 90% gives the best separation between the low and high lactate batches, which corresponds to 38 PCs. Therefore, in the subsequent analysis of the high lactate batches, we chose to use the residual components after excluding the first 38 PCs. Because there are overall 130 dimensions in the unfolded data, 92 dimensions are left in the residuals for subsequent fault classification and diagnosis. To further illustrate how the PCs are incapable of separating the high and low lactate batches, we show in Figure 7 the scores of the leading two PCs for the high and low lactate batches. It is clear from this figure that the leading PCs are not informative at all in separating them. Hierarchical Clustering in the PCA Residual Space. We use hierarchical clustering on the PCA residuals after removing 38 PCs to study if there are multiple groups within the high lactate batches. To calculate WSS, we adopt the complete linkage,34 which uses all members in a cluster to calculate the dissimilarity with another cluster. Figure 8 shows the dendrogram obtained by hierarchical clustering for the high lactate batches. From Figure 8, we chose a dissimilarity

Figure 5. CPV explained by the PCA model.

are 22.36, 33.92, 42.56, and 53.38, respectively, and the corresponding 95% SPE control limits are 53.47, 34.86, 25.83, and 17.87, respectively. These limits are shown as the horizontal red lines in the figure. The first 160 data points in the figure are for the low lactate batches, and the rest are for the high lactate batches. It should be noted that those high lactate batches exceeding the control limits do not mean that they are outside the regulatory requirements. The drug

Figure 6. T2 and SPE charts for 70, 80, 85, and 90% CPV captured. The residual space of the PCA model with 90% CPV best separates the low and high lactate batches. 13474

DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Article

Industrial & Engineering Chemistry Research

The two large clusters, clusters 1 and 2, have multiple batches that were manufactured consecutively in time, which are shown in parentheses. This information could be useful for further root cause diagnosis. Visualization of the Clusters in the LDA Space. The five clusters of high lactate batches relative to the low lactate cluster are visualized using an overall LDA of all clusters. Figure 9 shows how these batches distribute in the five leading LDA directions, ϕ1 through ϕ5, labeled as LD1 through LD5, in separation of the six clusters. It appears that clusters 3 and 4 are far away from the low lactate cluster, and they dominate the LD1 and LD2 directions, respectively, as shown in the LD1−LD2 plot in Figure 9. This dominance makes the leading discriminant directions incapable of separating other clusters. Therefore, the leading LDA directions are not necessarily the best directions to visualize all clusters. In fact, the subsequent LD3, LD4, and LD5 directions are more effective in visualizing the separation of the high lactate clusters relative to the low lactate cluster, as shown in Figure 10. In this figure, arrows are used to connect the center of the low lactate cluster to those of the high lactate clusters. The separations of the low lactate cluster from all high lactate clusters are clearly visualized except for that of cluster 3, which dominates LD1. This evidence shows that the clusters are separable by LDA, and further diagnosis of the contributing features to the five clusters is feasible. Diagnosis Based on Pairwise LDA Directions. To examine how the high lactate clusters departed from the low lactate cluster, pairwise LDA is performed between the low lactate cluster and each of the five high lactate clusters to obtain the directions (i.e., ϕ1) along which the low lactate cluster departed. These directions help us diagnose the contributing variables of the high lactate clusters. Figure 11 shows the histograms of the LDA scores along the first LDA direction between each of the five high lactate clusters and the low lactate cluster. It is shown in the figure that clusters 1 and 2 have many batches, but clusters 3 to 5 have a few in each of them. The unbalanced distribution of data in different clusters can lead to ill-conditioning for regular LDA, so we apply regularized LDA in this work. We used 3-fold cross-validation with the average test accuracy as the performance measure. We found that the best value of γ is 0.6 because its average test accuracy is the highest. Figure 11 also shows that clusters 3 to

Figure 7. Scores of the low lactate and high lactate batches for the leading two PCs.

limit that leads to five clusters for all high lactate batches, which are labeled in the figure. This choice seems to give a good gap between a finer clustering and a coarser one. The batch IDs for the five high lactate clusters are summarized as follows. cluster 1: 479 521 876 878 (884 886 890 891 895 898 899 902 904) cluster 2: 480

536

566 596

(686 688

690

693

698

708 710

714

715

721

723

724

727 729) (756 760

762

766

768) 786 807

889) 897

905 913

(872 874) (888 (959 960) 981

1000 cluster 3: 523 725 cluster 4: 524 597 cluster 5: 804 879 892 906

Figure 8. Hierarchical clustering of the high lactate batches in the PCA residual subspace. 13475

DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Article

Industrial & Engineering Chemistry Research

Figure 9. Visualization of the low lactate cluster and five high lactate clusters by an overall LDA.

Figure 10. Visualization of the low lactate cluster and five high lactate clusters by LD3, LD4, and LD5.

5 are far away from the low lactate cluster and easily separable. Cluster 2 is marginally separated from the low lactate cluster. Figure 12 shows the leading discriminant direction between each of the five high lactate cluster and the low lactate cluster, with the bar charts in the left column depicting the loading contributions using the standard LDA, and the bar charts in the right column depicting the leading contributions using regularized LDA. To visualize the difference more clearly, we show in Figure 13 the contributions of each variable averaged through the batch time. It is observed that whereas the standard LDA spreads the contributions to almost all variables, the regularized LDA clearly identifies four variables with high

Figure 11. Histograms of the LDA scores along the first LDA direction between each of the five high lactate clusters and the low lactate cluster.

contributions: pO2, pCO2, ammonium, and sodium. The other six variables have almost no contributions to the fault situations. For each of the contributing variables, the peak contributions vary from cluster to cluster, showing that these clusters indeed represent different types of faults. 13476

DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Article

Industrial & Engineering Chemistry Research

Figure 12. Diagnosis of the five high lactate clusters by pair-wise LDA: contributions using standard LDA (left column) and contributions using regularized LDA (right column). Rows 1 through 5 correspond to clusters 1 to 5, respectively.

Because all five high lactate clusters have high contributions from the four variables, it is necessary to examine how different the directions are. If the regularized LDA directions for the high lactate clusters are nearly collinear or of the same direction, then they are not easily separable. On the other hand, if the angles between them are large, they are not collinear and are easily separable. Therefore, we calculate the angles between the leading directions ϕ1 of the five regularized LDA directions, which are shown in Table 1. It is observed that the smallest angle between them is greater than 55°, whereas a few are close to 90°. Therefore, it is evident that these high

lactate clusters are separable from each other because of different contributions from the four variables. In conclusion, by applying the proposed strategy, we are able to clearly separate low lactate batches from the clusters of high lactate batches. Although the five high lactate clusters are all related to unusual variations in pO2, pCO2, ammonium, and sodium, the regularized LDA directions for high lactate clusters are very different, as shown by the angles between them. For all five high lactate clusters, the regularized LDA outperforms standard LDA in focusing on contributing variables. 13477

DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Article

Industrial & Engineering Chemistry Research

Figure 13. Diagnosis of the five high lactate clusters by pairwise LDA: averaged contributions using standard LDA (left column) and averaged contributions using regularized LDA (right column). Rows 1 through 5 correspond to clusters 1 to 5, respectively.

shown that the PCA residual space is sensitive to the high and low lactate batches in the process, whereas the principal components are not useful for the discrimination. The contribution of this work is a novel application and combination of data analytics tools to a complex cell-growth process with limited measurements. Process operation knowledge is incorporated in selecting the early stage of the process data for discrimination and diagnosis. One of the findings in the work is that the normally active variations in the principal components are of little use in the discrimination of the various high and low lactate batches. By applying hierarchical clustering on the PCA residual subspace, various clusters of faults are clearly identified and visualized. Another finding is that the leading LDA directions can be distorted more by a cluster far away from the rest of the data, and the subsequent LDA directions offer a more balanced visualization of all clusters. A third issue in this application is the uneven number

Table 1. Angles between the Five Leading Directions of Regularized LDA for Each of the Five High Lactate Clusters vs the Low Lactate Cluster angles (°)



cluster cluster cluster cluster

1 2 3 4

cluster 2

cluster 3

cluster 4

cluster 5

63.1

55.4 56.5

80.8 62.8 78.9

87.9 73.8 74.8 83.7

CONCLUSIONS AND DISCUSSIONS The proposed classification and diagnosis strategy is successfully applied to an industrial biopharmaceutical batch process with data collected over 8 years. Data from the early stage of the batch process are shown as being able to distinguish low lactate batches from high lactate batches. It is 13478

DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Article

Industrial & Engineering Chemistry Research

(4) Zagari, F.; Jordan, M.; Stettler, M.; Broly, H.; Wurm, F. M. Lactate metabolism shift in CHO cell culture: the role of mitochondrial oxidative activity. New Biotechnol. 2013, 30, 238−245. (5) Mulukutla, B. C.; Yongky, A.; Grimm, S.; Daoutidis, P.; Hu, W.S. Multiplicity of steady states in glycolysis and shift of metabolic state in cultured mammalian cells. PLoS One 2015, 10, No. e0121561. (6) Larson, T. M.; Gawlitzek, M.; Evans, H.; Albers, U.; Cacia, J. Chemometric evaluation of on-line high-pressure liquid chromatography in mammalian cell cultures: Analysis of amino acids and glucose. Biotechnol. Bioeng. 2002, 77, 553−563. (7) Alcala, C. F.; Qin, S. J. Analysis and generalization of fault diagnosis methods for process monitoring. J. Process Control 2011, 21, 322−330. (8) Qin, S. J. Process data analytics in the era of big data. AIChE J. 2014, 60, 3092−3100. (9) Garcı ́a-Muñoz, S.; Kourti, T.; MacGregor, J. F.; Mateos, A. G.; Murphy, G. Troubleshooting of an industrial batch process using multivariate methods. Ind. Eng. Chem. Res. 2003, 42, 3592−3601. (10) Kirdar, A. O.; Green, K. D.; Rathore, A. S. Application of multivariate data analysis for identification and successful resolution of a root cause for a bioprocessing application. Biotechnology progress 2008, 24, 720−726. (11) Nomikos, P. Detection and diagnosis of abnormal batch operations based on multi-way principal component analysis World Batch Forum, Toronto, May 1996. ISA Trans. 1996, 35, 259−266. (12) Nomikos, P.; MacGregor, J. F. Multi-way partial least squares in monitoring batch processes. Chemom. Intell. Lab. Syst. 1995, 30, 97− 108. (13) Lee, J.-M.; Yoo, C.; Lee, I.-B. Fault detection of batch processes using multiway kernel principal component analysis. Comput. Chem. Eng. 2004, 28, 1837−1847. (14) Yu, J.; Qin, S. J. Multiway Gaussian mixture model based multiphase batch process monitoring. Ind. Eng. Chem. Res. 2009, 48, 8585−8594. (15) Wurl, R. C.; Albin, S. L.; Shiffer, I. J. Multivariate monitoring of batch process startup. Quality and Reliability Engineering International 2001, 17, 269−278. (16) Chen, J.; Liu, K.-C. On-line batch process monitoring using dynamic PCA and dynamic PLS models. Chem. Eng. Sci. 2002, 57, 63−75. (17) Cedeño, M. V.; Rodrı ́guez Aguilar, L. P.; Sánchez, M. C. Bioprocess statistical control: Identification stage based on hierarchical clustering. Process Biochem. 2016, 51, 1919−1929. (18) Singhal, A.; Seborg, D. E. Clustering multivariate time-series data. J. Chemom. 2005, 19, 427−438. (19) Husson, F.; Josse, J.; Pages, J. Principal component methods− hierarchical clustering−partitional clustering: why would we need to choose for visualizing data?; Technical Report of the Applied Mathematics Department; Agrocampus, 2010. (20) Tax, D. M.; Duin, R. P. Support vector data description. Machine learning 2004, 54, 45−66. (21) Ge, Z.; Gao, F.; Song, Z. Batch process monitoring based on support vector data description method. J. Process Control 2011, 21, 949−959. (22) Ge, Z.; Song, Z. Bagging support vector data description model for batch process monitoring. J. Process Control 2013, 23, 1090−1096. (23) Yao, M.; Wang, H.; Xu, W. Batch process monitoring based on functional data analysis and support vector data description. J. Process Control 2014, 24, 1085−1097. (24) Wang, Y.; Zhao, C. Probabilistic fault diagnosis method based on the combination of nest-loop fisher discriminant analysis and analysis of relative changes. Control Engineering Practice 2017, 68, 32− 45. (25) Li, F.; Vijayasankaran, N.; Shen, A.; Kiss, R.; Amanullah, A. Cell culture processes for monoclonal antibody production. MAbs 2010, 2, 466−479. (26) Zhao, C. Quality-relevant fault diagnosis with concurrent phase partition and analysis of relative changes for multiphase batch processes. AIChE J. 2014, 60, 2048−2062.

of batches in the various clusters and collinearity in the data. To deal with this situation, regularized LDA is adopted, and contributions based on the regularized LDA directions are shown to be effective for clear diagnosis. These findings provide a clear characterization and understanding of the high lactate batches. Future work is planned to understand whether the high lactate batches are attributable to uncontrolled variations in the batch initial conditions and the preceding steps of the manufacturing chain. The analytics performed on the multiyear industrial data in this paper reinforced the importance of process knowledge in applying machine learning and data analytics to industrial problems. A critical issue for process data analytics is the design of an appropriate strategy based on the best understanding of machine learning tools and process knowledge. The recent overview paper by Qin and Chiang38 highlights the importance of aligning the appropriate machine learning tools with process knowledge to come up with effective process analytics solutions. This work can serve as one example in showing that standard analytics have to be tailored to provide effective solutions to process analytics problems. For instance, although hierarchical clustering appears to be a versatile and elegant clustering approach, it would not be effective without being fed appropriately selected variables or features; irrelevant features included in the scheme tend to distort the conclusion. As another instance, although the principal components contain the most variation, it is the residuals that are more informative for distinguishing normal batches from high lactate batches in this project. Lastly, the LDA based analysis is supposed to best separate classes with reduced dimensions, but the leading discriminant dimensions can be dominated by a class far away from the rest, making them least useful in visualizing the separations of the other classes.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

S. Joe Qin: 0000-0001-7631-2535 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work is supported in part by Genentech, a member of the Roche Group. We appreciate the detailed comments and suggestions from the reviewers that helped improve the presentation of the results of the paper.



REFERENCES

(1) Charaniya, S.; Le, H.; Rangwala, H.; Mills, K.; Johnson, K.; Karypis, G.; Hu, W.-S. Mining manufacturing data for discovery of high productivity process characteristics. J. Biotechnol. 2010, 147, 186−197. (2) Tôrres, A. R.; de Oliveira, A. D. P.; Grangeiro Júnior, S.; Fragoso, W. D. Multivariate statistical process control in annual pharmaceutical product review. J. Process Control 2018, 69, 97−102. (3) Le, H.; Kabbur, S.; Pollastrini, L.; Sun, Z.; Mills, K.; Johnson, K.; Karypis, G.; Hu, W.-S. Multivariate analysis of cell culture bioprocess dataâĂ Ť lactate consumption as process indicator. J. Biotechnol. 2012, 162, 210−223. 13479

DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480

Article

Industrial & Engineering Chemistry Research (27) Nomikos, P.; MacGregor, J. Multivariate SPC charts for monitoring batch processes. Technometrics 1995, 37, 41−59. (28) Ge, Z.; Gao, F.; Song, Z. Batch process monitoring based on support vector data description method. J. Process Control 2011, 21, 949−959. (29) Joe Qin, S. Statistical process monitoring: basics and beyond. J. Chemom. 2003, 17, 480−502. (30) Qin, S. J. Survey on data-driven industrial process monitoring and diagnosis. Annual reviews in control 2012, 36, 220−234. (31) Box, G. Some Theorems on Quadratic Forms Applied in the Study of Analysis of Variance Problems, I. Effect of Inequality of Variance in the One-way Classification. Ann. Math. Stat. 1954, 25, 290−302. (32) MacGregor, J. F.; Kourti, M. Statistical process control of multivariate processes. Control Engineering Practice 1995, 3, 403−414. (33) Johnson, S. C. Hierarchical clustering schemes. Psychometrika 1967, 32, 241−254. (34) James, G.; Witten, D.; Hastie, T.; Tibshirani, R. An Introduction to Statistical Learning with Applications in R; Springer, 2013. (35) Guo, Y.; Hastie, T.; Tibshirani, R. Regularized linear discriminant analysis and its application in microarrays. Biostatistics 2007, 8, 86−100. (36) He, Q. P.; Qin, S. J.; Wang, J. A new fault diagnosis method using fault directions in Fisher discriminant analysis. AIChE J. 2005, 51, 555−571. (37) Nomikos, P.; MacGregor, J. F. Monitoring batch processes using multiway principal component analysis. AIChE J. 1994, 40, 1361−1375. (38) Qin, S. J.; Chiang, L. H. Advances and opportunities in machine learning for process data analytics. Comput. Chem. Eng. 2019, 126, 465−473.

13480

DOI: 10.1021/acs.iecr.9b01175 Ind. Eng. Chem. Res. 2019, 58, 13469−13480