Raw Material Variability and Its Impact on the Online Adaptive Control of Cohesive Powder Blend Homogeneity Using NIR Spectroscopy

: It is significant to analyze the blend homogeneity of cohesive powders during pharmaceutical manufacturing in order to provide the exact content of the active pharmaceutical ingredient (API) for each individual dose unit. In this paper, an online monitoring platform using an MEMS near infrared (NIR) sensor was designed to control the bin blending process of cohesive powders. The state of blend homogeneity was detected by an adaptive algorithm, which was calibration free. The online control procedures and algorithm’s parameters were fine-tuned through six pilot experiments and were successfully transferred to the industrial production. The reliability of homogeneity detection results was validated by 16 commercial scale experiments using 16 kinds of natural product powders (NPPs), respectively. Furthermore, 19 physical quality attributes of all NPPs and the excipient were fully characterized. The blending end time was used to denote the degree of difficulty of blending. The empirical relationships between variability of NPPs and the blending end time were captured by latent variable modeling. The critical material attributes (CMAs) affecting the blending process were identified as the particle shape and flowability descriptors of cohesive powders. Under the framework of quality by design (QbD) and process analytical technology (PAT), the online NIR spectroscopy together with the powder characterization facilitated a deeper understanding of the mixing process.


>
The near infrared (NIR) spectroscopy combined with modified adaptive algorithm was used to monitor the blending process of natural product powders and excipient systems. > The effectiveness of the online monitoring platform was justified through 16 batches of commercial blending operations. > An empirical relationship between raw material variability and detected blending end time were established. > The particle irregularity and powder flowability were identified as critical material attributes affecting the cohesive powder blending process.

Physical Characterization of Powder Materials
Nineteen physical parameters were used to comprehensively characterize the properties of NPPs and dextrin. Whenever possible, the pharmacopoeia test methods were preferred. Before measurement, all powders were sifted through a 24-mesh size sieve to remove any clumps. Each parameter value was the mean of three replicates.
Bulk density (Da): The bulk density was measured according to the method described in Section 2.9.34 of European Pharmacopoeia 9.0 (E.P. 9.0). In simple terms, 100 g (m) of powder was filled into a 250 ml graduated cylinder with a funnel, and the resulting volume (Va) was recorded. Bulk density was computed using Equation (1): Tapped density (Dc): The tapped volume (Vc) was obtained after 2500 strokes using a tap density tester (HY-100; Dandong Hengyu Instrument Ltd., Dandong, China). This parameter was computed using Equation (2): Based on the test results of bulk density and tapped density, inter-particle porosity (Ie), Carr index (IC) and Hausner ratio (IH) were calculated using Equations (3)-(5), respectively. The angle of repose (α) was measured according to the official method illustrated in Section 2.9.36 under E.P. 9.0. The sample was poured through a funnel onto a platform until a stable powder cone was formed. Then the height (h) and the diameter (r) of the powder cone was measured. Finally, the angle of repose was calculated using Equation (6): Flow time (t") was measured according to the official method illustrated in Section 2.9.16 under E.P. 9.0. A 100 g sample was added into a dry funnel whose bottom opening had been blocked. The times needed for the entire sample to flow out of the funnel were recorded after unblock the bottom.
Moisture content (MC) was measured using the Sartorius MA35 instrument (Sartorius AG, Goettingen, Germany). About 2.00 g of the test powder was tiled on a sample tray and heated at 105 • C. The percentage of moisture content was recorded until a constant weight was obtained.
The percentage of particles measuring <20 µm (pf ) was determined using a laser particle size analyzer (BT 2001; Dandong Bettersize Instrument Ltd., Dandong, China) according to the official method in Section 2.9.31 of E.P. 9.0. The particle size distribution of the powder was measured using air as a dispersion medium, and the percentage of the powder with a particle size less than 20 µm was calculated.
The homogeneity index (Iθ) test method is the same to the particle size analysis by laser light diffraction. Due to the small particle size of the NPP and dextrin, the four sizes used were 0.150, 0.075, 0.050, and 0.020 mm, and the percentage of powder in each particle size range was calculated. Iθ was calculated according to the Equation (7): where F m is the percentage of particles in the majority range, F m−1 is the percentage of particles in the range immediately below the majority range, F m+1 is the percentage of particles in the range immediately above the majority range, n is the order number of the fraction studied under a series with respect to the majority fraction, d m is the mean diameter of the particles in the majority fraction, d m−1 is the mean diameter of the particles in the fraction of the range immediately below the majority range, and d m+1 is the mean diameter of the particles in the fraction of the range immediately above the majority range. Hygroscopicity (H) was determined by pharmacopoeia method descripted in Section 5.11 under E.P. 9.0. Under a constant temperature, the ammonium sulfate saturated solution placed at the bottom of the dryer provided a constant humidity environment. Hygroscopicity could be determined by calculating the percentage of increasing weight of powders in the mentioned dryer after 24 h.
According to the particle size distribution curves tested by the laser diffraction instrument, D 10 , D 50 , and D 90 were the particle sizes at the 10th, 50th, and 90th percentiles, respectively. The width of powder size distribution (Span) was calculated by Equation (8): For the particle shape, the powder was pasted on a metal disc and placed on an ion sputtering apparatus (SC7620, Quorum, East Sussex, UK). Then, the surface morphology of the powder was observed using a scanning electron microscope (Quanta250, FEI, Hillsboro, OR, USA). Photographic images of about 40 particles for each sample were analyzed by Image Pro plus 6.3 software (Media Cybernetics, Rockville, MD, USA) to measure various particle shape descriptors, i.e., aspect ratio, roundness, irregularity, and radius ratio. The aspect ratio means a ratio between major axis and minor axis of ellipse equivalent to a particle; a perfect circle has an aspect ratio of 1. The radius ratio represents ratio between maximum radius and minimum radius, and the max or min radius is a maximum or minimum distance between each particle's centroid position. The irregularity is a ratio between the perimeter and the Feret diameter (the longest length of a particle); with an ideal ratio of π Processes 2019, 7, 568 6 of 23 meaning a perfect circle. The roundness can be calculated by Equation (9). A circular particle will have a roundness equal to 1, and the other shape will have a roundness greater than 1.

Near Infrared Spectroscopy
The Antaris Target Online Near Infrared Sensor (Thermo Nicolet, Madison, WI, USA) equipped with a semiconductor tunable light source (covering 5500-7500 cm −1 ) was used in the blending process (see Figure 1). For on-line detection, the sensor was firmly mounted onto a flush top of the blend bin (Wuhan Hengda Chang Machinery Equipment Co., Ltd., Wuhan, China) which was modified with a sapphire window. The rotating speed of bins for both the pilot scale and the commercial scale was set constant at 14 r/min. The position of the bin was continuously monitored. When the NIR sensor was rotated between 170 • and 190 • , the NIR sensor was triggered to collect the diffuse reflection spectra through the sapphire window. The sample was scanned 8 times to obtain an average spectrum during each rotation. A wireless access point enables successful communications between the computer and the sensor via the 802.11 wireless Ethernet standard. The RESULT software (Thermo Nicolet, Madison, WI, USA) used workflow controls to receive the spectral data and to synchronize data collection with the bin position. The parameters of the NIR were set as follows: Log (1/R) was absorbance data format, and the spectral scanning rang was within 7100-5500 cm −1 . The pilot experiment with a 75 L bin was carried out, and the resolution and magnification of NIR were, respectively, set to 8 cm −1 and 125. In the actual production process, a 3000 L bin was used. The magnification and resolution were, respectively, set to 625 and 4 cm −1 . The sequence NIR data files collected from time-based blending experiments was achieved in the specified folder on the hard disk. The real-time monitoring of the blending process as well as the detection of the blending endpoint was achieved by the adaptive modeling algorithm.

Near Infrared Spectroscopy
The Antaris Target Online Near Infrared Sensor (Thermo Nicolet, Madison, WI, USA) equipped with a semiconductor tunable light source (covering 5500-7500 cm −1 ) was used in the blending process (see Figure 1). For on-line detection, the sensor was firmly mounted onto a flush top of the blend bin (Wuhan Hengda Chang Machinery Equipment Co., Ltd., Wuhan, China) which was modified with a sapphire window. The rotating speed of bins for both the pilot scale and the commercial scale was set constant at 14 r/min. The position of the bin was continuously monitored. When the NIR sensor was rotated between 170° and 190°, the NIR sensor was triggered to collect the diffuse reflection spectra through the sapphire window. The sample was scanned 8 times to obtain an average spectrum during each rotation. A wireless access point enables successful communications between the computer and the sensor via the 802.11 wireless Ethernet standard. The RESULT software (Thermo Nicolet, Madison, WI, USA) used workflow controls to receive the spectral data and to synchronize data collection with the bin position. The parameters of the NIR were set as follows: Log (1/R) was absorbance data format, and the spectral scanning rang was within 7100-5500 cm −1 . The pilot experiment with a 75 L bin was carried out, and the resolution and magnification of NIR were, respectively, set to 8 cm −1 and 125. In the actual production process, a 3000 L bin was used. The magnification and resolution were, respectively, set to 625 and 4 cm −1 . The sequence NIR data files collected from time-based blending experiments was achieved in the specified folder on the hard disk. The real-time monitoring of the blending process as well as the detection of the blending endpoint was achieved by the adaptive modeling algorithm.

Principles of the Adaptive Modeling Algorithm
In this paper, a modified adaptive algorithm including three steps, i.e., dynamic data pretreatment, process variation identification, and reliable homogeneity judgement, was proposed. The detailed procedures of the adaptive algorithm are as follows.
Step 1: The NIR sensor is used to collect the spectral data along with the mixing process and the data are subjected to real-time dynamic pretreatment.

Principles of the Adaptive Modeling Algorithm
In this paper, a modified adaptive algorithm including three steps, i.e., dynamic data pretreatment, process variation identification, and reliable homogeneity judgement, was proposed. The detailed procedures of the adaptive algorithm are as follows.
Step 1: The NIR sensor is used to collect the spectral data along with the mixing process and the data are subjected to real-time dynamic pretreatment. S1-1: At time point h, all the raw near infrared data are combined into a matrix G (size m × n, where m is the number of samples, and n is the number of spectral variables). S1-2: The chemometric method, e.g., multiplicative scatter correction (MSC) or standard normal variate transformation (SNV), is used to preprocess the spectral matrix G to obtain the matrix X (size m × n). S1-3: At the time point h + 1 (or h + z, where z is an integer), the calculation process is repeated from step S1-1.
Step 2: Two windows are established to simultaneously move on the preprocessed NIR data. The data in one window (modeling window) are used to build the principal component analysis (PCA) model and to set the control limit. The data in the other window (prediction window) are monitored and the number of samples larger than critical value is counted. S2-1: On the matrix X, two parallel sample matrices A (size w 1 × n) and B (size w 2 × n) are set, where w 1 is the number of (m − w 1 − w 2 + 1) to (m − w 2 ) samples in matrix X, w 2 is the number of (m − w 2 + 1) to m samples in matrix X, and m≥ (w 1 + w 2 ). The B matrix is used to build the PCA model: In Equation (10), the matrix T (size w 2 × k) is the score matrix, the matrix P (n × k) is the loading matrix, the matrix E (w 2 × n) is the residual matrix, and k is number of principal components.
S2-2: Based on the PCA model, the critical value is set according to Equation (11): In Equation (11), F (k,w 2 −k,α) is the F-distribution critical value at degrees of freedom k and w 2 − k, and α is the significance level.
S2-3: The sample data in the matrix A are projected in the principal component space in the direction of the load P, and the score matrix T new (size w 1 × k) of the matrix A is calculated: S2-4: The Mahalanobis distance d i (i = 1, 2, 3, . . . , w 1 ) for each sample in matrix A is calculated using Equation (13): where the vector t i (size 1 × k) represents the score vector of the sample i in the matrix T, mu is the mean vector of the matrix B, and S is the covariance matrix of the matrix B. Comparing the values of di with D crit , the number of atypical samples in matrix A is counted.
Step 3: Along with the proceeding of the mixing process, the number of atypical samples generated at each time point is recorded and displayed in the graphical user interface (GUI). In order to enhance the reliability of end-point detection, a lag time (LT) could be set. When the number of atypical samples is kept at zero in the LT interval, the materials in the mixer are considered to reach the blend homogeneous state.
Several advantages can be found in this modified adaptive algorithm. First, at the current monitoring window, all previous NIR spectra are preprocessed to reflect the system dynamics in real time. Second, the modelling window is situated behind the prediction window, since the former window is considered to contain more variation information than the latter window along with the proceeding of mixing process. Third, the reliability of mixing endpoint detection is improved by setting the LT reasonably.

Reference Analysis
To validate the reliability of the blending homogeneity detected by the online NIR monitoring system, a stratified sampling plan was made. Immediately after the blending endpoint was found, the blender was stopped by the programmable logic controller (PLC). Samples were taken by a thief probe at 11 locations of the industrial bin. As shown in Figure 2, sampling locations 1, 2 and 3 stand for the upper powder bed, sampling locations 4, 5 and 6 stand for the bin corner, sampling locations 7, 8 and 9 represent the middle layer of the powder bed, sampling site 10 represents the lower surface of the powder bed and the tank blind angle, and sampling site 11 denotes the sampling site on the axis of the powder bed. Figure 2 shows 11 sampling locations. The contents of APIs were detected by an Agilent 1260 High Performance Liquid Chromatograph (Agilent, Santa Clara, CA, USA,) with a UV detection scheme. The chromatographic conditions referred to used a 250 mm × 4.6 mm column (Agilent Extend C-18, Agilent, Santa Clara, CA, USA,), and the column temperature was maintained at 25 • C. Other information, such as the mobile phase and the detection wavelength, are listed in Table  S1 (Supplementary Materials). All HPLC methods were validated in accordance with the Chinese Pharmacopoeia (ChP 2015).  The stratified sampling plan for the industrial production. The yellow frame represents the surface of the powder bed, and the blue stars represents stratified sampling sites.

Software
The homemade KRT software running on the MATLAB 2009a platform (MathWorks, Natick, MA, USA) could realize all steps of adaptive algorithms and can be operated at the manufacturing site. The partial least square (PLS) models were developed by SIMCA 13.0 software (Umetrics, Umea, Sweden). An analysis of variance (ANOVA) was performed on SAS 8.0 software (SAS, Carey, NC, USA).

Physical Properties of Powder Materials
The measured 19 physical quality attributes for 22 NPPs and dextrin (material No. 23) are shown in Table 2. All tested material properties were stored in the iTCM database [42]. The SEM images of all materials are displayed in Figures 3 and S1 (Supplementary Materials). This information reflected different aspects of powder properties, such as particle shape, particle size, density, flowability and stability.
For the particle size, the D10, D50 and D90 values of NPPs powders were relatively small compared to that of frequently used pharmaceutical excipients in the report [42,43]. Except for dextrin (Span value equal to 0.88), the Span values of NPPs were large, indicating that these materials had a wide Figure 2. The stratified sampling plan for the industrial production. The yellow frame represents the surface of the powder bed, and the blue stars represents stratified sampling sites.

Software
The homemade KRT software running on the MATLAB 2009a platform (MathWorks, Natick, MA, USA) could realize all steps of adaptive algorithms and can be operated at the manufacturing site. The partial least square (PLS) models were developed by SIMCA 13.0 software (Umetrics, Umea, Sweden). An analysis of variance (ANOVA) was performed on SAS 8.0 software (SAS, Carey, NC, USA).

Physical Properties of Powder Materials
The measured 19 physical quality attributes for 22 NPPs and dextrin (material No. 23) are shown in Table 2. All tested material properties were stored in the iTCM database [42]. The SEM images of all materials are displayed in Figure 3 and Figure S1 (Supplementary Materials). This information reflected different aspects of powder properties, such as particle shape, particle size, density, flowability and stability. indicating that the fine particle content (FPC) in these materials were high. The D50 values of these materials were within a range of 6.11-48.54 μm. Except for Scrophulariae Radix (No. 14) and Chuanxiong Rhizoma (No. 15), the D50 values of all powders were less than 30 μm. They can be classified into Geldart's group C powders that are cohesive. The range of D90 was within 12.56-127.80 μm, and only three materials' D90 values were larger than 100 μm. This indicates that most powders in the experiment had cohesive properties.  Table 1, and No. 12 is the dextrin image).  Table 1, and No. 12 is the dextrin image).
For the particle size, the D 10 , D 50 and D 90 values of NPPs powders were relatively small compared to that of frequently used pharmaceutical excipients in the report [42,43]. Except for dextrin (Span value equal to 0.88), the Span values of NPPs were large, indicating that these materials had a wide particle size distribution. The values of pf were large and the values of Iθ were extremely small, indicating that the fine particle content (FPC) in these materials were high. The D 50 values of these materials were within a range of 6.11-48.54 µm. Except for Scrophulariae Radix (No. 14) and Chuanxiong Rhizoma (No. 15), the D 50 values of all powders were less than 30 µm. They can be classified into Geldart's group C powders that are cohesive. The range of D 90 was within 12.56-127.80 µm, and only three materials' D 90 values were larger than 100 µm. This indicates that most powders in the experiment had cohesive properties.
With respect to density properties, the fluctuations of Da and Dc values were relatively large. Additionally, there was a strong linear relationship between Da and Dc (R 2 = 0.87). The values of Ie varied greatly from 0.31 to 3.15, which meant that different voids may be formed during powder accumulation. The values of IH had a relatively wide range (1.23-1.91), and the IC values (18.69-47.75) changed obviously. Both IH and IC are qualitative indicators of the cohesiveness of granular materials. In accordance with Geldart et al. [44], particles with IH values larger than 1.4 were classified into group C powders. From the perspective of IH, it was found that 18 batches of powders in the experiment were cohesive powders. According to the E.P. 9.0, when the IH was greater than 1.35 or the IC value was greater than 26, the flowability of powder was considered poor. Based on the experimental data and pharmacopoeia standards, the flow properties of most materials in the experiment were generally poor. As revealed by Capece et al. [45,46], the inter-particle cohesiveness correlated well to the powder flowability. Poor flowability indicated greater inter-particle cohesiveness. However, the IH values of processed Fructus Xanthii (No. 2) and Mume Fructus (No. 21) were less than 1.25, and the corresponding IC values were the smallest and the second smallest, respectively. This shows that the two powders behave more like group A powders, which are mildly cohesive. The same phenomenon was also observed by Berthiaux et al.-some powders like talc, whose cohesion was smaller than expected, cannot be classified easily with the comparison between cohesion and the Carr index [47]. The reason may be that the Geldart's classification was only based on particle size and density, and the powder flowability was influenced by a variety of material properties including particle shape and surface compositions.
With respect to the particle shape factors, it can be seen form Table 2 that the roundness values were in a range of 1.02-1.37. The roundness values of most powders were close to 1, indicating that the particles approached a circle, while the rest of the powders were not round. The values of the aspect ratio were within the range of 1.15-1.73, and radius ratios fluctuated from 1.27 to 2.79, indicating that some particles were close to circular and some had certain elongation. The descriptors of irregularity are ranged from 2.71 to 3.10, showing that the particles' morphology was different to a certain extent. It is especially noteworthy that the irregularity of processed Fructus Xanthii (3.10) and Mume Fructus (3.00) were closest to π, indicating that they were closer to the circle. This may be a reason for their good flowability. Additionally, it was clear from the SEM image ( Figure 3 and Figure S1 in Supplementary Materials) that the particles were of different sizes, and this could be confirmed by particle size distribution of these powders. Though these particles still exhibited inhomogeneous shapes, there were almost no acicular or sharp-edged structures. In particular, primary particles of dextrin (No. 12 in Figure 3) could be agglomerated to larger unit, which may have altered their size and shape.
For the stability factor, the MC values were in a range from 2.29 to 8.89. By contrast, the range of hygroscopicity (H) was extremely large from 6.05 to 39.57. The H values of all NPP materials were greater than 15%, indicating they were very hygroscopic according to the E.P. 9.0 standards in Section 5.11. The reason for the high hygroscopicity of the NPPs may be that they contain a lot of hydrophilic compounds, such as glycosides, organic acids, carbohydrates, and amino acids [49].

Principle Component Analysis of Physical Properties
In order to discuss the relationship among physical properties of various materials, a PCA model was established, and a detailed explanation can be seen below.
First principal component (PC1): Descriptors of the particles size and density. The first PC explained 35.0% and predicted 7.82% of the variability in the properties of the raw materials and excipient. Figure 4 shows that particle size related parameters such as D 10 , D 50 , D 90 and Iθ are located on the right side of the loading scatter plot. The pf is in the opposite direction. A powder with a high FPC will cause D 10 , D 50 , D 90 and Iθ to decrease. The density related descriptors such as Da and Dc had positive loading values to the first PC, revealing that powder with high density may have reduced the FPC. Figure 4 illustrates that most powder parameters not only changed with the PC1 but also with the second principal component (PC2). Nevertheless, the Iθ, pf, H and MC are rather unique parameters that combined a strong loading for PC1 with a neutral loading for PC2. ratio and the roundness exhibited moderately positive PC1 loadings and high positive PC2 loadings. As expected, the irregularity descriptor was located at the opposite direction of the loading scatter plot, i.e., the bottom left corner. Overall, the location of irregularity, other shape and flow-related parameters on the loading scatter plot demonstrated common pharmaceutical knowledge: Low irregularity meant that the particles were approximately spherical, and more spherical particles always flowed better.  PC2: Descriptors of the particles shape and flow. The second PC explained 27.0% and predicted 28.3% of the variability in the properties of the raw materials and excipient. The shape-related descriptors are displayed on the top right side in Figure 4, indicating that the radius ratio, the aspect ratio and the roundness exhibited moderately positive PC1 loadings and high positive PC2 loadings. As expected, the irregularity descriptor was located at the opposite direction of the loading scatter plot, i.e., the bottom left corner. Overall, the location of irregularity, other shape and flow-related parameters on the loading scatter plot demonstrated common pharmaceutical knowledge: Low irregularity meant that the particles were approximately spherical, and more spherical particles always flowed better.
The raw materials can be clustered in the score scatter plot if they possess similar properties. For example, Angelica sinensis (No. 20) and Dipsaci Radix (No. 21) in Figure 5 are located in the completely opposite direction, because the Ie values of the two powders are the maximum and the minimum, respectively. It can be seen from Figure 5 that the six batches of powders (No. 1-6) used in the pilot experiment are evenly distributed among materials used in the industrial experiment, indicating that the adaptive algorithm optimized by these six batches of powders has the potential of being transferrable to other NPP-dextrin binary systems.
Processes 2019, 7, x FOR PEER REVIEW 13 of 23 the adaptive algorithm optimized by these six batches of powders has the potential of being transferrable to other NPP-dextrin binary systems.  Table 1. No. 23 is dextrin. The red triangles represent the materials used in the pilot experiment, and green dots represent the materials used in the industrial production). Figure 6a shows the original NIR spectra for the dextrin and six kinds of NPPs used in the pilot experiment. It can be seen that the near infrared absorption spectra of the natural plant powders were different from the dextrin, indicating that the NIR sensor has the potential ability to distinguish the NPP and the excipient. In the pilot experiment, 840 NIR spectra for each sample blend were obtained during the mixing process, which lasted for one hour. Figure 6b shows changes of NIR spectra in mixing process of Lonicerae Japonicae Flos extract powder-dextrin. The zoom of a specific spectral region (7500~7100 cm −1 ) revealed that the NIR absorption spectrum was not smooth and appeared to fluctuate because of the large noise. Consequently, this regional band was cutoff in the following analysis. Then, the NIR spectra were preprocessed by the SNV method to eliminate the influence of the particle scattering and optical path difference, which was conductive to the optimization of the parameters of the adaptive algorithm.  Table 1. No. 23 is dextrin. The red triangles represent the materials used in the pilot experiment, and green dots represent the materials used in the industrial production). Figure 6a shows the original NIR spectra for the dextrin and six kinds of NPPs used in the pilot experiment. It can be seen that the near infrared absorption spectra of the natural plant powders were different from the dextrin, indicating that the NIR sensor has the potential ability to distinguish the NPP and the excipient. In the pilot experiment, 840 NIR spectra for each sample blend were obtained during the mixing process, which lasted for one hour. Figure 6b shows changes of NIR spectra in mixing process of Lonicerae Japonicae Flos extract powder-dextrin. The zoom of a specific spectral region (7500~7100 cm −1 ) revealed that the NIR absorption spectrum was not smooth and appeared to fluctuate because of the large noise. Consequently, this regional band was cutoff in the following analysis. Then, the NIR spectra were preprocessed by the SNV method to eliminate the influence of the particle scattering and optical path difference, which was conductive to the optimization of the parameters of the adaptive algorithm. The NIR data from the mixing process of Lonicerae Japonicae Flos extract powder and dextrin.

Parameters Optimization of the Adaptive Algorithm
Before industrial applications of the NIR spectroscopy, the parameters involved in the adaptive algorithm were optimized in order to make the algorithm sensitive toward process normal variations and robust enough toward signal disturbances. The significance level α, which is the F test level in

Parameters Optimization of the Adaptive Algorithm
Before industrial applications of the NIR spectroscopy, the parameters involved in the adaptive algorithm were optimized in order to make the algorithm sensitive toward process normal variations and robust enough toward signal disturbances. The significance level α, which is the F test level in the D crit control limit calculation formula, was set to 0.05 according to prevailing statistical considerations. Next, the number of principal components k and the window size w were fine-tuned. Six pilot blending experiments using, respectively, NPP materials No. 1-6 in Table 1 were carried out to investigate impacts of the window size and the principal component number on the performance of the adaptive algorithm. Each experiment was run for one hour, and a total of 840 NIR spectra were obtained. At given values of w and k, the total number of atypical samples during the blending process were recorded. By changing the values of w and k, a response surface for the total number of atypical samples can be plotted, as shown in Figure 7. All the six binary powder systems exerted a similar trend, in that the total number of atypical samples grew rapidly with a reduction of the number of principal components and an increase of the window sizes.
Six pilot blending experiments using, respectively, NPP materials No. 1-6 in Table 1 were carried out to investigate impacts of the window size and the principal component number on the performance of the adaptive algorithm. Each experiment was run for one hour, and a total of 840 NIR spectra were obtained. At given values of w and k, the total number of atypical samples during the blending process were recorded. By changing the values of w and k, a response surface for the total number of atypical samples can be plotted, as shown in Figure 7. All the six binary powder systems exerted a similar trend, in that the total number of atypical samples grew rapidly with a reduction of the number of principal components and an increase of the window sizes.
Generally speaking, the modeling window should be wide enough to provide the best possible representation of the changes in all stages of the mixing process [50]. The NIR datasets in the modeling window were then concentrated by the principal component analysis into the first few principal components, which could achieve a high signal to noise ratio. Besides, the combination of w and k with few atypical samples was favored. For the six experiments, the cumulative contribution rate of the two first principal components was more than 99.9%, which means that two first principal components could explain most of variations of the data. Consequently, the optimal number of principal component k of the adaptive algorithm was set to 2. As for the window size (w), it could be found that when the window size was 1-5, as the numbers of atypical samples were too small to reflect the mixing dynamics. When the window sizes exceed 15, the total number of occurrence samples dramatically increases. Considering the above criteria for parameters optimization, the window size (w) was set to 10. The optimal settings were combinations with few total atypical samples and a wider windows size.

R eal Time Detection of the Blending Homogeneity during Industrial Production
With the optimized parameters w being set to 10 and k being set to 2, the adaptive algorithm was utilized to monitor the industrial blending process of natural plant granules in real time. Since the Generally speaking, the modeling window should be wide enough to provide the best possible representation of the changes in all stages of the mixing process [50]. The NIR datasets in the modeling window were then concentrated by the principal component analysis into the first few principal components, which could achieve a high signal to noise ratio. Besides, the combination of w and k with few atypical samples was favored. For the six experiments, the cumulative contribution rate of the two first principal components was more than 99.9%, which means that two first principal components could explain most of variations of the data. Consequently, the optimal number of principal component k of the adaptive algorithm was set to 2. As for the window size (w), it could be found that when the window size was 1-5, as the numbers of atypical samples were too small to reflect the mixing dynamics. When the window sizes exceed 15, the total number of occurrence samples dramatically increases. Considering the above criteria for parameters optimization, the window size (w) was set to 10. The optimal settings were combinations with few total atypical samples and a wider windows size.

Real Time Detection of the Blending Homogeneity during Industrial Production
With the optimized parameters w being set to 10 and k being set to 2, the adaptive algorithm was utilized to monitor the industrial blending process of natural plant granules in real time. Since the mixing bin was 3000 L, the lag time was set to 300 s. If the number of atypical samples continued to be 0 during this lag time, a high assurance for accurate detection of the homogeneity would be provided.
Sixteen kinds of NPPs, i.e., No. 7 to the No. 22 in Table 1, were blended with dextrin in the GMP-certified room. Figure 8 shows an example of the monitoring results from the mixing process of Mume Fructus extract and dextrin. At the beginning of the mixing stage, the atypical samples were rapidly raised then promptly reduced to 0, indicating that the NPP and dextrin were mixed briskly and reached the homogeneous state. At this stage, the powder may have been in the action of strong convection and may have overcome the hierarchical trend caused by material property differences. With the mixing process progressing, the number of occurrences started to fluctuate within a small range, and the mixed powder at this stage may be in a state of mixing and segregation due to the random diffusion of the powder particles. Finally, in the late stage of the mixing process, the number of overrun samples continued to be 0 within the lag time, indicating that the NPP and dextrin reached the homogeneous state. The stop time for the blender was 33.47 min, and the corresponding endpoint time was 28.47 min. The on-line monitoring results of the other fifteen materials can be seen from Figure S2 in Supplementary Materials. The corresponding endpoint times are listed in Table 3. It can be seen that the mixing time of the Ecliptae herba extract powder and dextrin was the longest (36.07 min), while that of Cibotii Rhizoma extract powder was the least (2.60 min).
Processes 2019, 7, x FOR PEER REVIEW 16 of 23 mixing bin was 3000 L, the lag time was set to 300 s. If the number of atypical samples continued to be 0 during this lag time, a high assurance for accurate detection of the homogeneity would be provided. Sixteen kinds of NPPs, i.e., No. 7 to the No. 22 in Table 1, were blended with dextrin in the GMPcertified room. Figure 8 shows an example of the monitoring results from the mixing process of Mume Fructus extract and dextrin. At the beginning of the mixing stage, the atypical samples were rapidly raised then promptly reduced to 0, indicating that the NPP and dextrin were mixed briskly and reached the homogeneous state. At this stage, the powder may have been in the action of strong convection and may have overcome the hierarchical trend caused by material property differences. With the mixing process progressing, the number of occurrences started to fluctuate within a small range, and the mixed powder at this stage may be in a state of mixing and segregation due to the random diffusion of the powder particles. Finally, in the late stage of the mixing process, the number of overrun samples continued to be 0 within the lag time, indicating that the NPP and dextrin reached the homogeneous state. The stop time for the blender was 33.47 min, and the corresponding endpoint time was 28.47 min. The on-line monitoring results of the other fifteen materials can be seen from Figure S2 in Supplementary Materials. The corresponding endpoint times are listed in Table 3. It can be seen that the mixing time of the Ecliptae herba extract powder and dextrin was the longest (36.07 min), while that of Cibotii Rhizoma extract powder was the least (2.60 min).   In order to verify the accuracy of the detected homogeneity, samples were taken at 11 locations inside the mixing bin by thief sampling immediately after the blender was stopped, according to the predefined sampling protocol. Each location was sampled three times, and the API content of 33 samples were determined by HPLC. An analysis of variance (ANOVA) was used to determine whether the mean contents of 11 sampling sites were different. The mean square between group (MS b ) and the mean square within group (MS w ) were calculated to represent the variation between sample means and variation within the samples, respectively. For the one-way ANOVA, the ratio of the between-group variability to the within-group variability followed a probability distribution known as an F-distribution. If the probability (p) was low enough, the null hypothesis that all the group means are equal was rejected using the common significance level of 0.05. The ANOVA results are provided in Table 3. Except for material No. 8 and material No. 16, the p values were all larger than 0.05, revealing that there was no difference among the mean content and that the binary powder system reached a homogeneous state. As for material No. 8 and material No. 16, it could be seen that the MS w values were close to zero and were therefore sensitive to small changes in the MS b values. The relative standard deviation (RSD) values of the mean content of 11 sampling sites were 1.7% and 0.67% for material No. 8 and material No. 16, respectively, indicating that the between-group variability was small (<3.0%) and the uniformity was acceptable.

Effect of Powder Properties on Mixing Process
The endpoint time is a measure of the difficulty degree to which the powder system reached a homogeneous state. The longer the mixing time is, the more difficult it is mixed. Through correlation analysis, it could be seen that the relationship between mixing time and proportion of added dextrin was extremely poor (R 2 < 0.01). Therefore, the factors affecting the mixing time were further explored from the perspective of the physical properties of powdered materials. Because cohesive powders can hardly flow freely, the flow time (t") was not used as an independent variable to develop a model. The independent variables matrix was composed of 16 samples and numeric difference of 18 physical parameters (Da, H, MC, IC, Ie, α, IH, Iθ, pf, Dc, aspect ratio, radius ratio, roundness, irregularity, Span, D 10 , D 50 , and D 90 ) between the natural plant powder and the dextrin. Taking Lonicerae Japonicae Flos extract powder as an example, the numeric difference of bulk density (Da) could be obtained by calculating the difference of Da between NPP (Da = 0.446) and dextrin (Da = 0.439), and then by computing the absolute value (0.007). For ease of understanding, the independent variables matrix of physical properties difference is shown in Table S2 (Supplementary Materials). The dependent variables matrix was constructed from the endpoint time of 16 blends detected by NIR with the adaptive algorithm. Subsequently, a PLS model was established to correlate the independent variables and the response variable. From Table 4, it can be seen that the first three latent factors can explain 83.9% of the Y-variation (30.4%, 43.4%m and 10.1%). Each PLS factor can explain 33.4%, 10.6% and 20.1% of the X-variation. A PLS biplot ( Figure 9) was drawn to jointly visualize all samples and all variables. A significant correlation of variables can be suggested when variables are projected between the outer circle (r = 1) and the middle circle (r = 0.75). There was a strong relationship among flowability variables α, IH and IC, among the aspect ratio, the radius ratio and roundness and among variables D 50 , D 90 , pf and Iθ. Moreover, it can be seen from Figure 9 that α, IH, and IC are located near the mixing time, indicating that the differences of flowability properties were positively correlated with the mixing time. The greater the differences of these properties were, the longer the mixing time was needed. The critical raw material attributes (CMAs) influencing the mixing time were investigated by the variable importance in the projection (VIP) index. As shown in Figure 10, irregularity, α, H, and IC have VIP values larger than 1. The error bar of the hygroscopicity passed through zero, reflecting that the influence of H on the mixing process was uncertain. Therefore, irregularity, αand IC were identified as CMAs. A PLS biplot ( Figure 9) was drawn to jointly visualize all samples and all variables. A significant correlation of variables can be suggested when variables are projected between the outer circle (r = 1) and the middle circle (r = 0.75). There was a strong relationship among flowability variables α, IH and IC, among the aspect ratio, the radius ratio and roundness and among variables D50, D90, pf and Iθ. Moreover, it can be seen from Figure 9 that α, IH, and IC are located near the mixing time, indicating that the differences of flowability properties were positively correlated with the mixing time. The greater the differences of these properties were, the longer the mixing time was needed. The critical raw material attributes (CMAs) influencing the mixing time were investigated by the variable importance in the projection (VIP) index. As shown in Figure 10, irregularity, α, H, and IC have VIP values larger than 1. The error bar of the hygroscopicity passed through zero, reflecting that the influence of H on the mixing process was uncertain. Therefore, irregularity, α and IC were identified as CMAs.   Both α and IC were related to the flowability of powders; the quite different flow performance of dextrin and NPP means that the mixing process involved in the movement of high and low cohesive particles. In the blending process, high cohesive particles have more contact with other particles than those of low cohesion [51], causing a low particle velocity. This can seriously affect the movement between powders, leading to low particle diffusion and resulting in a longer mixing time. Among 16 binary powders, natural plant granules of Ecliptae Herba (No. 8) and Mume Fructus (No. 21) have the longest mixing times. It could be seen from Table 2 that the flowability indexes of the two materials were better than that of dextrin. This meant that the cohesion of dextrin may be stronger than that of two materials, and dextrin was hard to be de-agglomerated, making it difficult for dextrin to evenly distribute in binary powders. As for Radix Paeoniae Rubra (No. 13), Scrophulariae Radix (No. 14) and Chuanxiong Rhizoma (No. 15), the differences in the flowability indexes were small-hence, the mixing times were smaller than that of the samples No. 8 and 21. However, due to the biggest difference in the irregularity for the sample No. 13, the mixing time was consequently longer than that of the samples No. 14 and 15.
In the two batches binary mixtures for Alismatis Rhizoma (No. 11) and Rhei Radix Et Rhizoma (No. 18), there was no significant difference in flowability performance. However, the irregularity difference of the two systems was quite different in the Table S2 (Supplementary Materials), which Both α and IC were related to the flowability of powders; the quite different flow performance of dextrin and NPP means that the mixing process involved in the movement of high and low cohesive particles. In the blending process, high cohesive particles have more contact with other particles than those of low cohesion [51], causing a low particle velocity. This can seriously affect the movement between powders, leading to low particle diffusion and resulting in a longer mixing time. Among 16 binary powders, natural plant granules of Ecliptae Herba (No. 8) and Mume Fructus (No. 21) have the longest mixing times. It could be seen from Table 2 that the flowability indexes of the two materials were better than that of dextrin. This meant that the cohesion of dextrin may be stronger than that of two materials, and dextrin was hard to be de-agglomerated, making it difficult for dextrin to evenly distribute in binary powders. As for Radix Paeoniae Rubra (No. 13), Scrophulariae Radix (No. 14) and Chuanxiong Rhizoma (No. 15), the differences in the flowability indexes were small-hence, the mixing times were smaller than that of the samples No. 8 and 21. However, due to the biggest difference in the irregularity for the sample No. 13, the mixing time was consequently longer than that of the samples No. 14 and 15.
In the two batches binary mixtures for Alismatis Rhizoma (No. 11) and Rhei Radix Et Rhizoma (No. 18), there was no significant difference in flowability performance. However, the irregularity difference of the two systems was quite different in the Table S2 (Supplementary Materials), which resulted in the mixing time of sample No. 11 (9.77 min) being longer than that of Sample No. 18 (4.77 min). This may be due to fact that the shape of the particles is related to the flowability [52,53]. An irregular shape leads to particle mutual occlusion, hindering particle movement. This is not conductive to the rearrangement and the diffusion of particles.
Additionally, it can be seen from Table S2 in Supplementary Materials that the No. 20 sample had the largest Da difference. This indicates that higher bulk density ratios may lead to segregation [23], which will affect the mixing process, resulting in prolonged mixing time.
As for Taraxaci Herba (No. 16) and Cirsii Herba (No. 22), the contribution plots for two binary mixtures are shown in Figure 11. In Figure 11a,b, the dominating bars exhibit the H parameter deviates most from the model average, which indicates that the hygroscopicity of the powder can affect mixing time. Hygroscopic powders can increase cohesion and reduce flowability by adsorbing moisture in the air [54], leading to longer mixing time. These results reflected that the mixing process of cohesive powder systems was complex, and a distinction should be made between the primary and secondary influencing factors for a particular cohesive powder mixture.
Processes 2019, 7, x FOR PEER REVIEW 20 of 23 resulted in the mixing time of sample No. 11 (9.77 min) being longer than that of Sample No. 18 (4.77 min). This may be due to fact that the shape of the particles is related to the flowability [52,53]. An irregular shape leads to particle mutual occlusion, hindering particle movement. This is not conductive to the rearrangement and the diffusion of particles. Additionally, it can be seen from Table S2 in Supplementary Materials that the No. 20 sample had the largest Da difference. This indicates that higher bulk density ratios may lead to segregation [23], which will affect the mixing process, resulting in prolonged mixing time.
As for Taraxaci Herba (No. 16) and Cirsii Herba (No. 22), the contribution plots for two binary mixtures are shown in Figure 11. In Figure 11a,b, the dominating bars exhibit the H parameter deviates most from the model average, which indicates that the hygroscopicity of the powder can affect mixing time. Hygroscopic powders can increase cohesion and reduce flowability by adsorbing moisture in the air [54], leading to longer mixing time. These results reflected that the mixing process of cohesive powder systems was complex, and a distinction should be made between the primary and secondary influencing factors for a particular cohesive powder mixture.

Conclusions
In this paper, the blending process in the preparation of diverse natural product powders with typically cohesive characteristics was taken as the research object. The MEMS NIR sensor combined with the modified adaptive modeling algorithm was successfully used to detect the blending homogeneity of cohesive powders during their industrial production. Compared with the previous fixed mixing time during the production, the online MEMS NIR spectroscopy technique can significantly reduce mixing time, improve mixing efficiency, and save production costs.
The innovation of this work is that it investigates the critical material attributes affecting the mixing time by developing an empirical relationship. It was found that the particle shape and powder flowability descriptors play a significant role in mixing cohesive powders. Thereby, decreasing the difference of flowability or particle shape between powders would be beneficial in improving blending processability. In future studies, more materials could be used to verify the accuracy of the empirical model, and numerical simulation studies could be conducted to understand the mixing mechanism of cohesive powders. Additionally, this paper provided a reference for other types of miniaturized NIR spectrometer to be applied in the pharmaceutical blending process.

Conclusions
In this paper, the blending process in the preparation of diverse natural product powders with typically cohesive characteristics was taken as the research object. The MEMS NIR sensor combined with the modified adaptive modeling algorithm was successfully used to detect the blending homogeneity of cohesive powders during their industrial production. Compared with the previous fixed mixing time during the production, the online MEMS NIR spectroscopy technique can significantly reduce mixing time, improve mixing efficiency, and save production costs.
The innovation of this work is that it investigates the critical material attributes affecting the mixing time by developing an empirical relationship. It was found that the particle shape and powder flowability descriptors play a significant role in mixing cohesive powders. Thereby, decreasing the difference of flowability or particle shape between powders would be beneficial in improving blending processability. In future studies, more materials could be used to verify the accuracy of the empirical model, and numerical simulation studies could be conducted to understand the mixing mechanism of cohesive powders. Additionally, this paper provided a reference for other types of miniaturized NIR spectrometer to be applied in the pharmaceutical blending process.