6104
Ind. Eng. Chem. Res. 2003, 42, 6104-6111
PROCESS DESIGN AND CONTROL Modeling and Control of a Microthermal Cycler for DNA Polymerase Chain Reaction Su Whan Sung* Corporate R&D, LG Chem Ltd./Research Park, 104-1, Moonji-dong, Yuseong-gu, Daejeon 305-380, Korea
Sunwon Park Department of Chemical & Biomolecular Engineering and Center for Ultramicrochemical Process Systems, Korea Advanced Institute of Science and Technology, 373-1 Guseong-dong, Yuseong-gu, Daejeon 305-701, Korea
Dae Sung Yoon, Youngsun Lee, and Geunbae Lim Biochip Project Team and MEMS Laboratory, Samsung Advanced Institute of Technology, San 14 Nongseo-Ri, Kiheung-Eup, Yongin-City, Kyunggi-Do 449-712, Korea
We have developed a microthermal cycler composed of a microreactor integrated with a platinum sensor/heater and automation equipment/software such as a data acquisition system, a control program, and a graphical user interface. From the control point of view, we have analyzed the dynamic characteristics of the fabricated microreactor and found some interesting dynamic features. On the basis of the analysis, we suggest an appropriate model structure and estimate the model parameters using the prediction error identification method. Requirements for highperformance operation are discussed, and a nonlinear control strategy is proposed to linearize the nonlinear dynamics of the microthermal cycler. We determine the parameters of the nonlinear controller using the optimal tuning method. The developed microthermal cycler shows excellent control performances and performs a successful DNA polymerase chain reaction. 1. Introduction The miniaturized processes integrated with microreactors, microsensors, and microactuators using semiconductor fabrication technologies can realize extremely high-efficiency and high-throughput operation, as well as surprising reductions of reagents, cost, reaction time, power, etc. The microthermal cycler is one of the most typical miniaturized devices. In particular, it is useful for the DNA polymerase chain reaction (PCR), which requires rapid temperature control to shorten the total running time and precise temperature control for high efficiency. Many types of microthermal cyclers have been proposed. Belgrader et al.1 developed a compact, battery-powered fluorometric thermal cycler that consisted of two reaction modules for multiplex real-time PCR. Northrup et al.2 designed a portable thermal cycle system including silicon-based reaction chambers with integrated heaters for efficient temperature control and optical windows for real-time fluorescence monitoring. Also, a silicon-based thermal microsystem has been fabricated, and its precise temperature control and rapid heating and cooling have been demonstrated.3 In the latter work, a gain-scheduling algorithm was used for * To whom all correspondence should be addressed. Email:
[email protected]. Tel.: +82-42-866-5785. Fax: +82-42-863-7466.
the proportional-integral (PI) controller to incorporate the nonlinearity of the thermal cycler. In this research, we mainly focus on the systematic modeling and control of the microthermal cycler in the belief that a system-level analysis and optimization will contribute to maximization of the performance of the microthermal cycler and provide some insights into its optimal operation. We fabricated a silicon-based microreactor integrated with a platinum sensor/heater and hardware/software for data acquisition, control, and power supply. We then analyzed dynamic thermal characteristics of the microthermal cycler, especially from the control and modeling points of view. We propose an appropriate model structure on the basis of the dynamics analysis and estimate the model parameters using the prediction error identification method. Requirements for the high-performance operation of this device are discussed, and a nonlinear control strategy linearizing the nonlinear dynamics of the microthermal cycler is proposed. We use the optimal tuning method to obtain the adjustable parameters of the controller. The developed microthermal cycler shows excellent control performances and performs a successful DNA polymerase chain reaction. 2. Construction of Microthermal Cycler We manufactured the microthermal cycler as shown in Figure 1a. It is composed of six parts: a silicon-based
10.1021/ie030031d CCC: $25.00 © 2003 American Chemical Society Published on Web 10/31/2003
Ind. Eng. Chem. Res., Vol. 42, No. 24, 2003 6105
Figure 1. Microthermal cycler: (a) overall scheme, (b) top view of the microreactor.
6106 Ind. Eng. Chem. Res., Vol. 42, No. 24, 2003
microreactor, a cooling fan, an amplification circuit, a data acquisition system, an external power supply, software for the automatic control algorithm, and a graphical user interface. The control signal of the microthermal cycler flows as follows: The amplification circuit amplifies the voltage of the platinum sensor (equivalently, the temperature of the microreactor) and transfers the amplified voltage to the analog input of the data acquisition system as shown in Figure 1a. Then, the automatic control algorithm adjusts the analog output to control the temperature as rapidly and precisely as possible, and the graphical user interface graphically displays the temperature and analog output on screen. Subsequently, the external power supply powers the platinum heater in proportion to the analog output of the data acquisition system. Then, the voltage of the platinum sensor (equivalently, the temperature of the microreactor) changes and passes through the amplification circuit again. The whole procedure is repeated for each sampling time. The amplification circuit plays an important role in increasing the resolution in reading the temperature. The platinum heater connected to the external power supply is to heat the microreactor, while the cooling fan is for rapid cooling of the microreactor. We use a nonlinear proportional-integral (PI) controller as the control algorithm. The graphical user interface is a useful tool for users who are unfamiliar with the system to operate the microthermal cycler as they want. The graphical user interface includes various functions for the user’s convenience such as easy real-time scheduling of the desired temperature profile, manual setting of the controller tuning parameters, real-time plotting of the temperature, and so on. The silicon-based microreactor was integrated with a thin-film platinum sensor and heater as shown in Figure 1b. The microchip was fabricated through several steps as shown in Figure 2. At step 4, the silicon was etched to a depth of 100 µm with tetramethylammoniumhydroxide (TMAH) for the microreactor and the microchannel. The thermal oxide film at step 5 served as an electrical insulation layer. At step 7, titanium and platinum film were deposited by dc off-axis magnetron sputtering. To guarantee high resolution in the data acquisition, the following devices were integrated to form the data acquisition system: the AD7715 analog-digital converter of Analog Devices, Inc., for the 16-bit analog input; the DAC8043 digital-analog converter of Analog Devices, Inc., for the 12-bit analog output; and the AT89C2051 microcontroller of ATMEL, Inc., for communication with the personal computer. 3. Modeling of Microthermal Cycler In this section, we develop a model to mathematically represent the dynamics between the reactor temperature and the voltage applied to the heater. After emphasizing several interesting dynamic characteristics of the microthermal cycler, we determine the model structure on the basis of the dynamic characteristics. Then, we estimate the model parameters using the prediction error identification method and demonstrate the model performances. 3.1. Determining Model Structure. Before estimating the model parameters, we should construct the
Figure 2. Steps in the fabrication of the microreactor.
model structure. The following five items are considered for the best model structure selection: (1) For electrical heating systems such as microthermal cyclers, the voltage has frequently been chosen as the input of the model. However, this choice is not recommended because the temperature is proportional to the electrical power [i.e., the square of the voltage (v2) divided by the resistance of the heater (R)] rather than the voltage. If the voltage were as the model input, the square nonlinearity between the model output and the model input could not be avoided. Thus, the square of the voltage should be chosen as the model input to avoid the nonlinearity. (2) A step test is one of the simplest techniques for detecting the dynamic characteristics of processes.4 Figure 3 shows the step response of the microthermal cycler. It should be noted that there is a sudden jump from the initial temperature (67.9 °C) to around 80 °C at the instance of the input change. This phenomenon is due to the position of the sensor in the middle of the heater. Because the sensor is very close to the heater,
Ind. Eng. Chem. Res., Vol. 42, No. 24, 2003 6107
mance of the feedback control system.5 Thus, this observation of a nearly zero time delay is very favorable for the design of a high-performance feedback controller. (4) Another dynamic characteristic of the step response is that the jump is not a vertical line and the transition from the initial jump to the next dynamic response of the temperature is smooth rather than clearly separated. This means that the order of the dynamic system is at least 2. (5) Finally, it should be noted that there is always a small input nonlinearity for various reasons such as the resistance variation of the heater during heating, nonlinear characteristics of the data acquisition system and the power supply, etc. Putting all the pieces together, we should consider the following model requirements to choose an appropriate model structure for the microthermal cycler: (1) The input and output of the dynamic model should be the square of the voltage and the microreactor temperature, respectively. (2) The dynamic model should contain a small negative zero. (3) The time delay can be chosen as zero. (4) The order of the dynamic microthermal cycler is at least 2. (5) Input nonlinearity should be included. Especially if the model does not satisfy the first and second requirements, we cannot achieve acceptable model performances by increasing the model parameters. We suggest the following nonlinear model structure that satisfies the above model requirements
Figure 3. Step response of the microthermal cycler: (a) step input, (b) temperature.
the rapid heat transfer from the heater to the sensor is initially dominant, whereas the slow heat transfer from the center to the outside of the heater becomes dominant after a while, resulting in the sudden jump. This is strong evidence that the transfer function from the input (i.e., the square of the voltage) to the output (i.e., the temperature) includes a small negative “zero” (here, zero is defined as the root that makes the transfer function equal to 0). For a simple justification of the abovestatement, consider the following second-order transfer function for example
y(s) )
5s + 1 u(s) 4s + 10s + 1 2
(1)
where u(s) and y(s) are the Laplace transforms of the input and output, respectively. The transfer function has a negative zero of -0.2. Note that the process output y(s) includes the time derivative of the input [i.e., 5su(s), or equivalently, 5 du(t)/dt]. The time derivative [du(t)/ dt] of the input jumps from zero to infinity at the instance of the step input change and again becomes zero after the instance, which results in the response of Figure 3. However, the sudden jump cannot be expected if the derivative term of the input is not included. (3) We also recognize from the step response of Figure 3 that there is nearly no time delay between the input and the output. The existence of time delays seriously deteriorates the maximum achievable control perfor-
u(t) ) v(t)2
(2)
q(t) ) u(t) + p1u2(t) + ‚‚‚ + pm-1um(t)
(3)
dx(t) ) Ax(t) + Bq(t) dt
(4)
y(t) ) Cx(t) + e(t)
(5)
where u(t) and y(t) denote the model input (the square of the voltage) and the model output (the temperature), respectively. e(t) is a white measurement noise. x(t) is the n-dimensional state. System matrices A-C have the following forms
[
0 1 0 A) l 0 0
0 0 1 l 0 0
0 0 0 l 0 0
‚‚‚ ‚‚‚ ‚‚‚ l ‚‚‚ ‚‚‚
0 0 0 l 0 1
-an -an-1 -an-2 l -a2 -a1
]
(6)
B ) [bn bn-1 bn-2 ‚‚‚ b2 b1 ]T
(7)
C ) [0 0 0 ‚‚‚ 0 1 ]
(8)
Here, eq 3 is the nonlinear static function, and the system of eqs 4 and 5 is the linear dynamic system. Thus, the model from eqs 3-5 is a Hammerstein-type model.6 The chosen dynamic model of eqs 2-8 satisfies the model requirements as follows: Equation 2 is to satisfy the first requirement. If we did not fix bi, i ) 1, 2, ..., n - 1 at zeros, the second requirement would be satisfied.
6108 Ind. Eng. Chem. Res., Vol. 42, No. 24, 2003 Table 1. Estimated Model Parameters and Modeling Error model parameters model structure
a1
a2
b1
b2
p1
p2
ISEa
n ) 2, b1 ) 0, u(t) ) m)1 n ) 2, u(t) ) v(t), m ) 1 n ) 2, u(t) ) v2(t), m ) 1 n ) 2, u(t) ) v2(t), m ) 2 n ) 2, u(t) ) v2(t), m ) 3
177.4 5.972 2.246 2.657 2.604
1.028 0.195 0.281 0.272
56.70 9.090 9.781 9.790
84.98 27.92 1.666 2.236 2.194
0.0129 0.0142
0.0010
1232.4 1464.5 122.7 18.43 17.83
v2(t),
a
ISE ) integral of the square error.
Because eq 4 does not include time delay, the third requirement is satisfied. The process order of the dynamic system composed of eqs 4 and 5 is n, so we should have n g 2 to satisfy the fourth requirement. The fifth requirement is incorporated by the nonlinear form of eq 3. 3.2. Estimating Model Parameters. Above, we itemized the model requirements for the microthermal cycler and established the model structure that satisfies the requirements. In this section, we use the prediction error identification method7 to estimate the model parameters for the given model structure of eqs 2-8. This approach estimates the model parameters by minimizing the cost function as follows
{
min V(P ˆ ,A ˆ ,B ˆ) ) P ˆ ,A ˆ ,B ˆ
0.5 N
[y(ti) - yˆ (ti)]2 ∑ N i)1
}
(9) Figure 4. Performances of the proposed model.
subject to
u(t) ) v2(t)
(10)
qˆ (t) ) u(t) + pˆ 1u2(t) + ‚‚‚ + pˆ m-1um(t)
(11)
dxˆ (t) )A ˆ xˆ (t) + B ˆ qˆ (t) dt
(12)
yˆ (t) ) Cxˆ (t)
(13)
t0 ) 0 < t1 < ‚‚‚ < tN-1 < tN
(14)
where y(t) and yˆ (t) denote the measured process output and the predicted model output, respectively. To solve the optimization problem, we use the LevenbergMarquardt optimization method because its convergence rate is high and it is robust. All equations needed to calculate the derivatives of the cost function in the Levenberg-Marquardt method can be easily derived from eqs 9-13, as is done in Sung et al.7 3.3. Model Performances. We activated the microthermal cycler using a roughly tuned PI controller to generate test data. Table 1 shows the model performances of three linear model types: the simple linear second-order model with v2(t) as the model input [n ) 2, b1 ) 0, u(t) ) v2(t), m ) 1], the linear second-order model that has a zero and the voltage as the model input [n ) 2, u(t) ) v(t), m ) 1], and the linear second-order model that has a zero and v2(t) as the model input [n ) 2, u(t) ) v2(t), m ) 1]. These estimation results are consistent with the expected results from the model requirements in subsection 3.1. Both of the second-order models with the second-order nonlinear polynomial function [n ) 2, u(t) ) v2(t), m ) 2] and the third-order nonlinear polynomial function [n ) 2, u(t) ) v2(t), m ) 3] provide satisfying identification results. We identified the third-order model also. Its performance was found to be almost the same as that of the second-order model.
We finally choose the simplest model [n ) 2, u(t) ) v2(t), m ) 2] because its performance is close to the best. The performance of this model is shown in Figure 4. 4. Control Strategy In this section, we will establish an automatic control strategy using the nonlinear proportional-integral (PI) controller shown in Figure 5 on the basis of the identified nonlinear model. It satisfies the following three control requirements for a high-performance microthermal cycler. 4.1. Linearization. It should be noted that, because the PI controller is linear, the corresponding process also should be linear to achieve the full performance of the PI controller. To satisfy this requirement, the proposed control strategy linearizes the nonlinear dynamics of the microthermal cycler using the identified nonlinear model shown in Figure 5a, where the dotted-line box is equivalent to the linear process of eqs 4 and 5 if the model is perfect. That is, the linear PI controller controls the linear process as shown in Figure 5b, which makes the tuning of the PI controller easier and also guarantees high performances as in the linear case. 4.2. Tuning of the PI Controller. For the tuning of the proposed controller in Figure 5a, we do not have to consider the nonlinear parts of eqs 10 and 11 because Figure 5a is essentially equivalent to Figure 5b. That is, we can tune the PI controller for the linear model of eqs 12 and 13. To tune the PI controller, we estimate the tuning parameters by minimizing the following integral of the time-weighted square value of the error (ITSE) for the unit step set point change8 N
min{V(kp,ki) ) kp,ki
ti[yd(ti) - yˆ (ti)]2∆t} ∑ i)1
(15)
Ind. Eng. Chem. Res., Vol. 42, No. 24, 2003 6109
Figure 5. Linearization of the nonlinear dynamics of the microthermal cycler: (a) proposed nonlinear control strategy, (b) equivalent control system.
subject to N
[ys(ti) - yˆ (ti)]∆t ∑ i)1
qˆ (t) ) kp[ys(t) - yˆ (t)] + ki
ys(t) ) 1
(16)
and
eqs 12 and 13 where kp, ki, and ∆t represent the proportional gain, integral gain, and sampling time, respectively. ys(t) and yˆ (t) are the set point and the predicted model output, respectively. Equation 16 is a mathematical representation of the PI controller. yd(t) denotes the user-specified desired trajectory, which has the following form
yd(s) )
(
)
1 ys(s) τ s + 1.414τs + 1 2
(17)
As we decrease the time constant τ, the closed-loop response becomes faster but more sensitive to the modeling error and shows a larger overshoot. To solve the optimization problem, we use the LevenbergMarquardt optimization method. For details, refer to Sung et al.8 4.3. Integral Windup. The maximum achievable control performance for the microthermal cycler is very high because the cycler has almost negligible time delay and no unstable zeroes.9 On the other hand, the voltage range is confined between 0 and a maximum voltage (which cannot be sufficiently large because of resolution and cost problems). Therefore, we cannot avoid “integral windup” due to actuator saturation for such a highperformance operation. In this research, we stop the integral action of the PI controller to prevent the integral windup phenomenon when the control output is saturated. 4.4. Control Results. We choose τ ) 0.2 s as the time constant of the desired trajectory and the sampling time of ∆t ) 0.055 s. The estimated optimal tuning results for the desired trajectory are kp ) 0.2107 and ki ) 0.2107/0.2242. The control performance is excellent, as shown in Figure 6. It shows remarkable heating and cooling rates of approximately 36 and -22 °C/s, respectively, about 15 times faster than commercial PCR machines. The overshoot is less than 0.8 °C, and the steady-state error is less than (0.1 °C. Figure 7 demonstrates the results of using only a linear PI controller without linearizing the nonlinear dynamics
Figure 6. Control performances of the proposed control strategy.
of the microthermal cycler. The linear PI controller is tuned for the operating region of 55 °C. As a result, it shows acceptable performance for the corresponding region, but it shows poor control performances for the other operating regions of 72 and 95 °C because the linear controller cannot remove the nonlinear dynamics. 5. Experimental Results of DNA PCR Until now, we experimentally demonstrated that the temperature control performance of the developed microthermal cycler is superior to that of the conventional machine. In this section, we experimentally confirm that
6110 Ind. Eng. Chem. Res., Vol. 42, No. 24, 2003
Figure 8. Quantitative analyses of PCR products: (a) conventional bulk PCR, (b) microthermal cycler PCR.
Figure 7. Control performances of the linear PI controller.
the proposed microthermal cycler can successfully perform the DNA polymerase chain reaction. Conventional bulk PCR was performed using a commercial thermal cycler, and silicon-glass microreactor PCR was performed using our microthermal cycler. The PCR mixture [PCR Core System II (Promega Corp., Madison, WI)] contained a 10 × PCR buffer [10 mM Tris-HCl (pH 9.0), 50 mM KCl and Triton(r)X-100], 1.5 mM MgCl2, 200 µM dNTP, 1 µM upstream and downstream control primer, 1 ng/50 µL positive control plasmid DNA (supplied by kit), and 1.25 units/µL Taq DNA polymerase. In the case of the microthermal cycler, bovine serum albumin (BSA) was added to the PCR mixture to impede surface adsorption of the enzyme, primer, and template onto the wall of the microreactor. The concentration of BSA was adjusted to be 100 µg/mL, which is the best concentration found among our PCR experiments. The PCR buffer mixture was injected into the microreactor through the inlet using a syringe, and the inlet and outlet passages were sealed using acrylic plastic blocks and a VersaChem epoxy (ITW Performance Polymers, FL). After a preheating treatment for perfect denaturation of DNA, 30 thermal cycles were performed for PCR. The microreactor was heated to 95 °C for 2 min and then cycled for 30 cycles: 30 s at 95 °C, 30 s at 55 °C, and 1 min at 72 °C. It takes about 62 min to complete 30 cycles. A final extension was performed at 72 °C for 7 min.
The conventional bulk PCR was done in a GeneAmp(r) tube (Perkin-Elmer) with a DNA thermal cycler (Eppendorf, Inc). The bulk PCR using the Eppendorff mastercycler was conducted with the same thermal schedule as used in the microthermal cycler. However, it takes over 90 min to complete 30 cycles because the temperature ramp rate of the conventional PCR machine is much lower than that of the microthermal cycler. After the PCRs, the samples of both cases were extracted, and quantitative analyses were conducted using an Agilent 2100 bioanalyzer. Figure 8 presents the capillary electrophoresis data of the PCR products. It confirms that DNA is successfully amplified by the developed microthermal cycler without undesirable PCR bands. 6. Conclusions We developed a microthermal cycler composed of a silicon-based microreactor and other equipment/software required for full automation. The dynamic characteristics of the microthermal cycler were analyzed, and on the basis of this analysis, the model structure was determined. We identified the model parameters using the prediction error identification method and demonstrated the excellent model performance. A nonlinear control strategy was proposed to linearize the nonlinear dynamics of the microthermal cycler, which makes the tuning of the PI controller easier and guarantees high performances as in the linear case. We demonstrated high-performance operation with the antiwindup nonlinear PI controller tuned by the optimal tuning method. Finally, we confirmed that the proposed microthermal cycler could successfully perform the DNA polymerase chain reaction within a much shorter reaction time compared to the conventional bulk PCR machine. Acknowledgment This work was supported by the BK21 Project and Center for Ultramicrochemical Process Systems sponsored by KOSEF
Ind. Eng. Chem. Res., Vol. 42, No. 24, 2003 6111
Literature Cited (1) Belgrader, P.; Young, S.; Yuan, B.; Primeau, M.; Christel, L. A.; Pourahmadi, F.; Northrup, M. A. A Battery-Powered Notebook Thermal Cycler for Rapid Multiplex Real-Time PCR Analysis. Anal. Chem. 2001, 73, 286. (2) Northrup, M. A.; Benett, B.; Hadley, D.; Landre, P.; Lehew, S.; Richards, J.; Stratton, P. A Miniature Analytical Instrument for Nucleic Acids Based on Micromachined Silicon Reaction Chambers. Anal. Chem. 1998, 70, 918. (3) Lao, A. I. K.; Lee, T. M. H.; Hsing, I.; Ip, N. Y. Precise Temperature Control of Microfluidic Chamber for Gas and Liquid Phase Reactions, Sens. Actuators 2000, 84, 11. (4) Seborg, D. E.; Edgar, T. F.; Mellichamp, D. A. Process Dynamics and Control; John Wiley & Sons: New York, 1989 (5) Sung, S. W.; Lee, I. Limitations and Countermeasures of PID Controllers. Ind. Eng. Chem. Res. 1996, 35, 2596.
(6) Haber, R.; Unbehauen, H. Structure Identification of Nonlinear Dynamic SystemssA Survey on Input/Output Approches. Automatica 1990, 26, 651. (7) Sung, S. W.; Lee, I. Prediction Error Method for ContinuousTime Processes with Time Delay. Ind. Eng. Chem. Res. 2001, 40, 5743. (8) Sung, S. W.; Lee, T.; Park, S. Optimal PID Controller Tuning Method for Single-Input-Single-Output (SISO) Processes. AIChE J. 2002, 48, 1358. (9) Morari, M.; Zafiriou, E. Robust Process Control; PrenticeHall: Englewood Cliffs, NJ, 1989.
Received for review January 13, 2003 Revised manuscript received June 23, 2003 Accepted September 26, 2003 IE030031D
