Improving Pellet Quality in a Pharmaceutical Hot Melt Extrusion Process via PID Control and LOLIMOT-Based MPC

The aim of this paper is the development of a process control concept for a hot melt extrusion (HME) and pelletization process. The new concept should improve the particle size distribution of the pellets produced. Production of pellets containing an active pharmaceutical ingredient (API) can be achieved by means of HME, followed by a pelletization process step. The quality of pellets produced depends on the strand temperature at the pelletizer’s inlet and the pelletizer’s intake speed. This paper presents a strategy for the strand diameter and temperature control based on adjusting the cooling intensity on the cooling track between the HME and the pelletization step and altering the pelletizer’s intake speed. Two concepts are presented and compared to the open-loop operation of the system: the first one is model predictive control (MPC) in combination with a model based on the local linear model tree (LOLIMOT) algorithm, and the second one is PID control. The quality of the pellets produced was analyzed in terms of particle size distribution (PSD). By implementation of the two control concepts, strand temperature and diameter could be kept close to the desired set points. Consequently, the presented concepts yielded pellets with a narrower particle size distribution than the open-loop operation of the plant. The application of the presented control strategies can improve the quality of the pellets produced by an HME and pelletization system in terms of their particle size distribution.


Introduction
This paper presents a model-based control system for the strand temperature and strand diameter in a pharmaceutical hot melt extrusion (HME) process. HME processes are used in pharmaceutical manufacturing to produce advanced dosage forms, e.g., amorphous solid dispersions (ASDs) for improving solubility [1] and for producing oral dosage forms containing poorly soluble active pharmaceutical ingredients (APIs) [2] or for sustained-release tablets [3] and pellets [4]. HME processes do not require solvents and are reliable and reproducible [5]. After the extrusion step, the extruded strands can be cut into pellets of desired size, which can be filled into capsules or blended with excipients and pressed into tablets. There are several methods of producing pellets from the extrudates [6]. In this study, strand pelletization was used, which typically requires a cooling track between the extruder outlet and the pelletizer inlet to achieve suitable strand temperature. The cooling power and the inlet speed of the pelletizer influence the quality of pellets produced. As is the case with many other continuous pharmaceutical manufacturing processes [7][8][9][10][11][12], appropriate control strategies must be defined. Model-based approaches for control strategy development were successfully demonstrated in the pharmaceutical industry already [8,10,11,13,14]. For HME, only few publications are available which address control strategy development [5,15]. This study focuses on the cooling and pelletization step after the HME process itself. In [16], a study of a process control approach in a similar setup was reported. However, no details on the development of the feedback control were provided.
This study proposes two different control strategies that automatically adjust both of the above-mentioned parameters (i.e., strand temperature and strand diameter) to increase the pellet quality. One of them, i.e., the application of a datadriven modeling approach for controller synthesis in a pharmaceutical manufacturing setup, is new to the knowledge of the authors. The results achieved were compared to the process data and analytical results determined offline. The purpose of this manuscript is not to demonstrate that one or the other feedback control strategy performs better, but to show the achievable performance increase by implementation of feedback control compared to open-loop operation of the system. Quantitative measures are provided to compare the different techniques, and a discussion highlighting the pros and cons of the individual approaches is given. The application of model-based control design approaches to this specific application in pharmaceutical manufacturing is novel. In comparison to traditional approaches, which typically run the system in an open-loop fashion, i.e., keeping intake speed and air pressure constant, the implementation of feedback control significantly improves pellets quality in terms of their particle size distribution. A narrow particle size distribution is desired, since, on the one hand, homogeneously sized pellets provide a well-defined release profile of the final drug product, and on the other hand, the option of volumetric filling of pellets into capsules delivers more homogeneous capsule masses compared to pellets with wider particle size distribution. The demonstrated approaches allow for a systematic design of feedback control systems, which are able to narrow the particle size distribution of the pellets.
The paper is structured as follows: Section "Process Description" describes the plant setup and the equipment used. Section "Model Development" details the selected modeling approach, which is used in Sect. "Controller Design" for the controller design. Section "Results" discusses the results, and Sect. "Discussion and Outlook" concludes the paper.

Process Description
The pellet production line investigated consists of a twin screw extruder (Coperion ZSK18 twin screw extruder, https://www.coperion.com) equipped with two feeders (K-Tron KT20, https://coperion.com; Brabender MTS-Hyg, https://www.brabender-technologie.com), a strand pelletizer (Primo 60 E, http://www.maag.com) and a custom-made cooling track. A sketch of the manufacturing process is shown in Fig. 1, and a photo is provided in Fig. 2. As the model formulation, 40% (w/w) of vitamin B1 (thiamine mononitrate, pharma grade, DSM http://dsm.com) and 60% (w/w) of Eudragit E PO (pharma grade, Evonik https://evonik.de) were used as API and excipient, respectively. The total mass flow was set to 2 kg/h. One goal of the process is to manufacture homogeneously sized pellets. The pellets produced are cylindrically shaped. Their length (determined based on the ratio between the intake roll speed and the knife speed) is adjusted in the pelletizer. Given a constant mass flow at the extruder outlet, the strand diameter is a function of the pelletizer's intake roll speed. Under the assumption of constant strand density, the volume flow needs to be constant, and therefore, higher intake speeds cause lower strand diameters. Consequently, the strand diameter can be altered during runtime (by varying the intake speed of the pelletizer) and measured using a 3-axis laser measurement head ODOC 13TRIO and the corresponding processing unit USYS 200 (http://www.zumbach.com). For decreasing the strand temperature before it enters the pelletizer, a cooling track between the extruder's outlet and the pelletizer's inlet is required. The conveyor belt, whose speed is kept constant, transports the strand to the pelletizer. The quality of pellets produced greatly depends on the temperature of the strand at the pelletizer's inlet. Strand temperatures that are too high or too low can impair the cutting quality of the pelletizer (e.g., strands sticking to the intake rolls if the temperature is too high or breaking due to brittleness if the temperature is too low). To adjust the temperature of the strand at the inlet of the pelletizer, a cooling track consisting of a pressure regulator (FESTO VPPE-3-1-1, https://festo.com) and several air distribution nozzles was constructed. The pressure regulator sets the air pressure at the distribution nozzles, which affects the air mass flow and the heat transfer from the strand to the cooling air. The extruder die temperature may also influence the strand temperature, which is measured using an infrared pyrometer micro-epsilon CT-SF22-C3 with CF lens (https://www.microepsilon.de/) with a measurement spot of diameter 0.6 mm.
The equipment shown in Fig. 1 is connected to a SCADA system containing the software PLC (XAMControl, evon GmbH, https://www.evon-automation.com). The software PLC executed the control algorithm with a sampling interval of T s = 0.35 s. By using this sampling interval, the time required to compute the control algorithms remains below 50% of T s , i.e., more than 50% of the sampling interval are available as buffer time, before the next computation needs to be started.

Model Structure and LOLIMOT Algorithm
The inputs and outputs of the plant model are shown in Fig. 3. In the process investigated, the pelletizer's intake speed u 1 (in %; actual intake speed v in m/min is given by v ¼ 0:4521 m=min % u 1 À1:293 m=min; the relation between actuating signal u 1 and actual intake speed v was obtained experimentally) and air pressure u 2 (in % of 6 bar) are deemed actuating signals. Due to construction and process limitations, their values are restricted according to Eqs. (1) and (2). The excitation run shown in Sect. "Identification Data" allows lower values of u 2 . Due to a lack of process robustness at such low air pressure values, the lower limit of 10 % given in (2) was selected for the controller design presented in Sect. "Controller Design": Strand temperature y 1 (in°C) and diameter y 2 (in mm) are the controlled variables. At this point, it should be noted that the die temperature d 1 is typically set by the operator of the extruder. The selection of that temperature is not within the scope of the current paper, and it should not be altered by the feedback control strategies. However, in order to cover the typical range of die temperatures during operation, several different values are considered during the identification experiments; see Sect. "Identification Data". Therefore, in terms of controller synthesis, die temperature d 1 (in°C) is considered a measurable disturbance.
The suggested modeling approach relies on the so-called local linear model tree (LOLIMOT) algorithm [17,18]. The underlying model structure consists of multiple local linear 1 models (LLM). The output of each local model is computed via a static relation from the model's inputs. The local model outputs are weighed via the so-called validity function Φ and summed up to compute the global model output y. Figure 4 shows the model structure and the respective equations. The encircled part indicates the "core" model, and the shift operators denoted as "z −1 " and the feedback of model output ŷare used in order to model the system's dynamics. A shift operator "z −i " delays the input signal by i times the sampling interval, which is in our setup equal to i • 350 ms. Two models were developed: one for the strand temperature and the other for the strand diameter. For both of them, local model order n = 5 was chosen, as indicated by the maximum delay of output y. This value turned out to be a good compromise between model complexity and prediction accuracy. The number of local models, M, was restricted to be less than or equal to 10 in order to avoid overfitting. The local model parameters w i0 , …, w im and the validity function parameters c j1 , …, c jm and σ j1 , …, σ jm are determined by the LOLIMOT algorithm, which is responsible for partitioning the input space (inputs ûi in Fig. 4) into M local models.
The steps of the LOLIMOT algorithm are depicted in Fig. 5 for a system with 2 inputs: u1 and u2. Beginning with one local model 1-1, a split in each dimension is tested. The model parameters are identified based on the experimental data for both possible splits. Subsequently, the model error is computed. The split leading to the smallest error is chosen for further processing. The next step is to determine which of the local models, 2-1 or 2-2, contributes to a larger error. The one that does (2-1 in the example shown in Fig. 5) is divided further. More details of the LOLIMOT algorithm and further variants of the algorithm can be found in [18,19].

Identification Data
Experimental data for executing the LOLIMOT algorithm for system identification were generated by running specific experiments. The plant was excited by inputs u 1 , u 2 , and d 1 .
Responses y 1 and y 2 were captured (Fig. 6). The input signals were arranged for the "steady-state" periods to be available (see time interval between 120 minutes and 175 minutes) and the dynamic excitation of the system to be present (see randomly selected inputs between 0 minute and 120 minutes). For the "steady states", one input remained constant, while the other one was varied according to a staircase function. The duration of the staircase steps was selected such that a steady state will be reached at each level. Since disturbance input d 1 is not expected to undergo rapid variations, no dynamic excitation was performed. Obviously, the raw data 1 Strictly speaking, the model consists of locally-affine models. Fig. 1 Schematic of the HME process investigated captured have considerable noise. One reason for the fluctuations is the strand guidance, which allows for a minor radial movement that shifts the temperature measurement spot on the strand and affects the thickness measurement. In order to dampen the noise, a 50-element moving average filter was used. The filtered data used for the model identification (see the following subsection) and prediction by the identified model are presented in Fig. 6. Fig. 7 is providing a detailed view of one section of the identification data. Furthermore, a detailed view of a validation dataset is provided in Fig. 8. There, the measured validation data is compared to the LOLIMOT prediction and to the prediction obtained by the linear time varying (LTV) system which was used in the MPC setup; see Sect. "MPC Controller". Obviously, the LTV system output almost coincides with the LOLIMOT prediction.

Model Identification
LMNTool [20] was used to identify two dynamic models: one for predicting strand temperature y 1 and the other for predicting strand diameter y 2 . For both models, the inputs/ local model order and the maximum number of local models were chosen as described in Sect. "Model Structure and LOLIMOT Algorithm". Figure 6 and Fig.  7 present a comparison between the measured and predicted outputs, showing a sufficiently good agreement between the model and the real plant.

Controller Design
Two control concepts were designed and implemented at the real plant. The first one consists of two PID controllers, and the second one relies on a model predictive controller (MPC) based on the process model presented in Sect. "Model Development". These two concepts were chosen, because the PID approach is very common in industrial control applications, which is often used in comparative studies with respect to other control strategies. The LOLIMOT and MPC approach was chosen, because it allows for a systematic controller design, and it is applicable to processes, where no mechanistic model is available. Further, it offers the inherent capability of considering constraints (e.g., actuator saturation) and control of multi-input, multi-output systems. The following subsections outline these two concepts. The results are presented and discussed in Sect. "Results".

PID Controllers
Two independent PID controllers are designed for the strand temperature and thickness control. Since strand diameter y 2 is mainly affected by the pelletizer's intake speed u 1 , one feedback loop is designed for this actuator/control signal pairing. The other PID controller adjusts air pressure u 2 in order to control strand temperature y 1 . The control system is depicted in Fig. 9. The and the corresponding controller parameters are provided in Table 1.
An anti-windup strategy was implemented via integrator clamping [21]. The controller tuning was performed using simulation studies based on the process model. Next, a finetuning of the parameters was performed at the real plant. The final controller settings are summarized in Table 1. The controller gains are all negative since the plant gains from u 1 to y 2 and from u 2 to y 1 are both negative (the higher the pelletizer's intake speed, the lower the strand diameter; the higher the air pressure, the lower the strand temperature).

MPC Controller
A model predictive controller [22,23] with objective function was applied. Vectors r k ,ŷ k , and u k are composed of reference signals, controlled variables, and manipulated variables at time instant k, i.e.,  Vector u Δ, k denotes the difference between two consecutive actuator signals, i.e., The positive (semi-) definite matrices Q ≽ 0, R ≻ 0 and R Δ ≽ 0 are used for weighing the control error, the actuating effort and the actuating variable variation, respectively. Their values (see Table 2) are established in the simulation studies, followed by fine-tuning at the real plant. The weight of R is chosen much lower than of the other matrices, which allows the controller to fully exploit the admissible range of actuators.
The following constraints have to be taken into account when minimizing the objective function (4): Measureable disturbance d 1 is assumed to remain constant along the prediction horizon. Constraints (8) and (9) describe the system's dynamics. System matrices A k , B u,k , B d,k , and C k are computed using the LOLIMOT model, as described in [24]. Since these are not constant along the prediction horizon (they are computed separately for each prediction step k + i in (4)), a nonlinear optimization problem in the form of min u k ;…;u kþnc−1 J k ; subject to : 7 ð Þ; 8 ð Þ; 9 ð Þ ð10Þ has to be solved. The computation time of real-time execution of the optimization task monitored using the software PLC is around 30-40% of the sampling interval T s .
To account for the model-plant mismatch, the term r kþi −ŷ kþi in the objective function (4) is replaced by Fig. 9 Control concept using two separate PID controllers  (8) and (9)) of a validation dataset. The temperature of the validation dataset was corrected by an offset of 12°C. This was necessary due to a fact that the temperature sensor was slightly misaligned during the creation of the validation dataset. However, the dynamics of the temperature is still captured accurately, and a constant offset between measurement and model prediction is inherently compensated by the proposed control concept (see Eq. (11)) It is assumed that the difference between the model output y k and measurement y k remains constant along the prediction horizon.

Results
This section presents the results obtained using the proposed control concepts (i.e., PID and MPC) and open-loop experiments (i.e., no feedback control and constant u 1 and u 2 values). In order to compare those, process data (temperature y 1 and diameter y 2 ) and offline analysis of the pellets produced in terms of particle size distribution (PSD) were applied. For testing the feedback control loop, reference profiles for temperature and diameter were specified, and samples were taken, when the process was run at steady state of the specific set points; see Fig.  10 and Fig. 11. The open-loop operation was tested for one constant setting of the actuating signals. The plant inputs u 1 and u 2 were chosen to be around the center region of the feasible values and are provided in Sect. "Process Data".

Process Data
Reference profiles for the strand temperature and strand diameter were created, and tracking performance of the two control strategies was assessed. The profiles, the controlled variables, and the actuating signals are depicted in Fig. 10 and Fig. 11 for the MPC and PID control, respectively. The MPC tightly tracks the diameter profile, primarily by adjusting the pelletizer's intake speed. From 1100 s to 1200 s, a deviation in the diameter is noticeable. The controller reacts by increasing the pelletizer's intake speed until the diameter is close to the reference again. The temperature profile is not tracked with the same accuracy. There are three main reasons for the larger deviations: First, the weight of the temperature control in Q is much lower than that of the diameter control. Consequently, the controller prioritizes control actions that keep the diameter close to the reference. This can clearly be observed between 1000 s and 1200 s. Second, the selected temperature profile cannot be met at all times due to constraints of the plant (e.g., minimum/maximum cooling power). Obviously, between 800 s and 1000 s, the temperature does not reach the desired value of 75°C, although the cooling power is already reduced to its minimum value of 10%. Third, the temperature measurement in the process setup applied is very sensitive to the radial strand position. Minor deviations from the nominal strand positioning shift the small measurement spot away from the strand's center, causing a change in the temperature reading. Additionally, since the strand can freely move on the conveyor belt (perpendicularly to the transport direction), the airstream hits the strand under varying conditions, changing the temperature. This deviation between the model and the real plant can impair the control performance as well. Another parameter not  In addition to the closed-loop experiments discussed above, the process was operated in the open-loop configuration at the following set point: u 1 = 50%, u 2 = 30%; see Fig.  12. Pellet samples were taken from the closed-loop experiments (see gray-shaded numbering in Fig. 10 and Fig. 11) and the open-loop experiment. Their analysis is presented in the following subsection. From the process data, the root mean squared error between reference signal and actual diameter and temperature, respectively, was computed. N denotes the number of samples used for computation. The samples were selected according to the intervals used for taking offline samples; see grayshaded areas in Fig. 10 and Fig. 11. Table 3 summarizes the results.

Offline Analysis Data
The PSD measurements were performed using the Microtrac PartAn 3D dry particle image analyzer (http://microtrac.com). For PSD analysis, the area-equivalent diameter is applied. Area A is the average area in a sequence of images of a specific particle [25].
In Compared to the open-loop operation, the PSD slope is higher under both control concepts, indicating a narrower PSD using the feedback control. To quantitatively compare the pellets properties, the second central moment μ 2 [26] of each of the PSDs was computed and is provided in Table 3. Compared to the second central moment of the open-loop experiment, which is Table 3 RMSE for the proposed control concepts. Provided sample numbers correlate to the time intervals used for sampling; see Fig. 10 and Fig. 11. Further, the second central moment μ 2 of the measured PSD is provided (see Sect. "Offline Analysis Data")

Discussion and Outlook
Our results demonstrate the advantages of a closed-loop control of strand diameter and strand temperature over an openloop system operation. The PSD of pellets produced can be narrowed by applying any of the suggested control strategies. The two concepts have both pros and cons. While the "LOLIMOT+MPC" approach requires a plant model, tuning is very intuitive and can easily be accomplished by selecting the weights. Moreover, the control concept inherently takes interactions between actuating signals and controlled variables into account: e.g., the influence of a higher intake speed on the temperature (and not only on the diameter) is considered in the MPC strategy. Furthermore, rate constraints of the actuating signals can be implemented straightforwardly.
Measureable disturbance d 1 is also considered in the computation of actuating signals. While the PID concept does not require a plant model, tuning of the parameters is less intuitive.
Typically, tuning rules [21,27,28] can be used in order to select the appropriate parameters. However, the performance of the diameter control via PID control is better for most of the investigated steady-state points. In contrast, the temperature control is better using the "LOLIMOT+MPC" approach for most of the steady-state points. Taking into account the complete trial, including also the transients between the different set-points, the RMSE of the PID strategy is lower.
Our promising results suggest that the proposed methods can be applied on continuous pharmaceutical manufacturing lines, delivering a better pellet quality in terms of their PSD than a conventional operation with constant cooling track and intake speed settings. Future work will focus on improving the strand guidance and the IR-temperature sensor in order to achieve a more robust temperature measurement. Furthermore, the mass flow of the material could be considered in the modeling approach and the MPC control concept, allowing mass flow adjustments to the HME process during its runtime.
Funding Information Open access funding provided by Graz University of Technology. This work was partially funded by Austrian Research Promotion Agency FFG, program "Produktion der Zukunft," project number 858704.

Compliance with Ethical Standards
Conflict of Interest The authors declare that they have no conflict of interest.
Ethical Approval This paper does not contain any studies with human participants or animals performed by any of the authors.
Nomenclature A, Area (mm 2 ); D a , Area-equivalent diameter (mm); e, Control error, e = r − y; K P , Proportional gain PI controller (%/°C for PID 1, %/mm for PID 2); K I , Integral gain PI controller (%/(°C s) for PID 1, %/(mm s) for PID 2); K D , Derivative gain PI controller (% s/°C for PID 1, % s/mm for PID 2); μ 2 , Second central moment of PSD (mm 2 ); n c , Control horizon; n p , Prediction horizon; Q, Weighing matrix control error; RMSE d , Root mean square error of the diameter (mm); RMSE T , Root mean square error of the temperature (°C); R, Weighing matrix actuating effort; R Δ , Weighing matrix actuating variable variation; r 1 , Reference value strand temperature (°C); r 2 , Reference value strand diameter (mm); T s , Sampling interval (s); y 1 , Strand temperature, measured (°C); y 2 , Strand diameter, measured (mm);ŷ 1 , Strand temperature, predicted (°C);ŷ 2 , Strand diameter, predicted (mm); u 1 , Pelletizer intake speed (%); u 2 , Air pressure (in % of 6 bar); v, Pelletizer intake speed (m/min) Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.