Modelling and control of non-linear tissue compression and heating using an LQG controller for automation in robotic surgery

Abstract

Robot-assisted surgery is being widely used as an effective approach to improve the performance of surgical procedures. Autonomous control of surgical robots is essential for tele-surgery with time delay and increased patient safety. In order to improve safety and reliability of the surgical procedure of tissue compression and heating, a control strategy for simultaneously automating the surgical task is presented in this paper. First, the electrosurgical procedure such as vessel closure that involves tissue compression and heating has been modelled with a multiple-input–multiple-output (MIMO) non-linear system for automation simultaneous. After linearizing the models, the linear-quadratic Gaussian (LQG) is used to control the tissue compression process and tissue heating process, and the particle swarm optimization (PSO) algorithm was used to choose the optimal weighting matrices for the LQG controllers according to the desired controlling accuracy. The LQG controllers with optimal weights were able to track both the tissue compression and temperature references in finite time horizon and with minimal error (tissue compression – the max absolute error was $9.6057 \times 10^{- 5}$ m and temperature – the max absolute error is 0.4561°C). Compared with a LQG controller with weighting matrices chosen by trial and error, a PSO-based optimized controller provided the least error and faster convergence. We have developed a control framework for simultaneously automating the surgical tasks of tissue compression and heating in robotic surgery, and modelled the automation of electrosurgical task using LQG controller with optimal weighting matrix obtained using a PSO algorithm.

Keywords

LQG MIMO PSO robot-assisted surgery tissue compression and heating

Introduction

Robot-assisted surgery is increasingly being used in various surgical procedures (Nezhat et al., 2005). The benefits of robot-assisted surgery include reduced patient trauma, shortened recovery time, increased dexterity, elimination of tremors, motion scaling and increased ergonomics that reduces fatigue for the operating surgeon (Feng et al., 2012). Though current robot-assisted surgical procedures are performed with a robot faithfully following the surgeon’s hand motion in a master–slave teleoperation system, semi or fully autonomous control of surgical tasks will allow surgeons to let the robot finish tasks with more precision or accuracy, such as applying forces and temperatures to tissues, and eventually improving the safety and reliability of surgical procedure (Feng et al., 2012; Kanno et al., 2013; Low and Phee, 2006). Previous work on automation in robot-assisted surgery has focused on tissue compression (Sankaranarayanan et al., 2011; Yu, 2008; Yu and Chizeck, 2008a, 2008b; Yu et al., 2007) or tissue retraction (Patil and Alterovitz, 2010). However, simultaneous compression and heating of tissue are common electrosurgical procedures, such as vessel closure, dissection and coagulation.

In order to automate a surgical procedure, two steps need to be followed. First, a system model needs to be derived to describe the process of interest mathematically (Patil and Alterovitz, 2010; Sankaranarayanan et al., 2011; Yu, 2008; Yu and Chizeck, 2008a, 2008b; Yu et al., 2007). For the models of surgical tasks, there have been some research results that mainly include tissue compression (Patil and Alterovitz, 2010; Sankaranarayanan et al., 2011; Yu, 2008), tissue impedance (Belmont et al., 2013; Dodde, 2011) and tissue temperature models (Mandel et al. 2013). These models are non-linear, the equations of which are obtained based on experimental data. However, there is a lack of a systematic modelling or comprehensive application for the overall surgical procedure. In this paper, these models and experiment data will be comprehensively used to design a controller for surgical automation. The second step is the design of control algorithms based on the system model (Patil and Alterovitz, 2010; Sankaranarayanan et al., 2011; Yu, 2008; Yu et al., 2008). Relevant control approaches studied in the references include robot force control, robot impedance control, automatic needle placement, automatic ligation and force control of knot tying (Yu, 2008). However, the system models of surgical tasks, especially the temperature model, are rarely seen in the control design in these references. Our preliminary work on designing a PID controller for tissue compression and heating was the first work of its kind in automation of a surgical task (Sinha et al., 2014). Although some optimal control algorithms have been used in automatic control design (Lewis et al., 2012), they have rarely been seen in surgical robot design to obtain a more effective, robust and accurate controller.

In this paper, we present the control strategy for simultaneously automating the surgical task of tissue compression and heating to improve safety and reliability of the surgical procedure. Figure 1 shows a general overview of the environment to which our autonomous control strategy was applied. A robotic gripper clamps down on a piece of tissue and applies an AC current in order to desiccate it. Our goal is to develop a controller that can drive the gripper to compress the tissue to a desired thickness and then heat it up to the desiccation temperature. However, the model of tissue compression and heating is a non-linear multiple-input–multiple-output (MIMO) system. In this work, feedback linearization was used to cancel the non-linearities and a linear-quadratic Gaussian (LQG) controller was used to track the desired tissue compression and temperature. LQG is the combination of a Kalman filter, such as a linear-quadratic estimator (LQE), with a linear-quadratic regulator (LQR) (Hespanha, 2014; Rowley and Batten, 2014). The controller based on LQG is one of the most fundamental optimal algorithms. However, selecting the optimum values of the weighting matrices for the LQG is the key procedure for the LQG controller development. The common approach to select the weighting matrices is via trial and error (Das et al., 2013; Hespanha, 2014; Mohd et al., 2012; Rowley and Batten, 2014), but this method is time consuming, cumbersome and results in a non-optimized performance. Therefore, optimal selections of weighting matrices have been used in many prior works in recent years (Das et al., 2013; Mobayen et al., 2011; Mohd et al., 2012). In this paper, a PSO method was used to find the best weighting matrices for making the surgical robot finish tasks with a desired accuracy (Meng et al., 2014; Mobayen et al., 2011; Mohd et al., 2012).

Figure 1.

Electrosurgery on a block of tissue.

System modelling

In order to automate the electrosurgical task, which is a non-linear MIMO system with both tissue compression and heat constraints, it is broken down into three components: a tissue compression model, tissue impedance model and tissue temperature model. According to the electrosurgical task in this paper and the related research results (Belmont et al., 2013; Berjano, 2006; Dodde, 2011; Kim, 2008; Sommer et al., 2013; Yu, 2008), the following assumptions were made: movement of the mass is restricted to one direction and only resistive heating is modelled. The following sections go through the process of creating the models of tissue compression, impedance and temperature, which are then applied to the controller development.

Tissue compression model

Soft tissues exhibit an exponential-like response when forces are applied to them (Yu, 2008). Therefore, the one-dimensional (1D) tissue dynamics can be modelled through a mass-spring-damper system with a non-linear spring, as shown in Figure 2.

Figure 2.

Mass-spring-damper model.

The non-linearity of the spring in this model can be described by the Blatz-Ko model, which was created by William L. Ko and P.J. Blatz for applications in elastic theory (Rosen et al., 2008; Yu, 2008). Therefore, the tissue dynamics model can now be written as:

u_{1} = m \overset{\cdot\cdot}{p} + d \overset{\cdot}{p} - \frac{s α}{β + 1} ((1 - \frac{p}{L_{0}}) e^{β ({(1 - \frac{p}{L_{0}})}^{2} - 1)} - \frac{1}{{(1 - \frac{p}{L_{0}})}^{2}} e^{β \frac{p}{L_{0} - p}})

(1)

where $u_{1}$ is the force applied to the tissue, $p$ is the position of the robotic gripper, $m$ is the mass of the tissue, $d$ is the viscosity of the tissue, $s$ is the contact surface area between one side of gripper and tissue, $α$ and $β$ are parameters related to the stiffness of the tissue, and $L_{0}$ is the initial thickness of the tissue.

For the implementation of the controller, let $y_{1} = x_{1}$ , $x_{1} = p$ , $x_{2} = \overset{\cdot}{p}$ ; then Equation (2) can be rewritten in the state-space form as:

\begin{matrix} \overset{\cdot}{x} = [\begin{matrix} {\overset{\cdot}{x}}_{1} \\ {\overset{\cdot}{x}}_{2} \end{matrix}] = [\begin{matrix} x_{2} \\ - \frac{d}{m} x_{2} + \frac{s α}{m (β + 1)} ((1 - \frac{x_{1}}{L_{0}}) e^{β ({(1 - \frac{x_{1}}{L_{0}})}^{2} - 1)} - \frac{1}{{(1 - \frac{x_{1}}{L_{0}})}^{2}} e^{β \frac{x_{1}}{L_{0} - x_{1}}}) \end{matrix}] + [\begin{matrix} 0 \\ \frac{1}{m} \end{matrix}] u_{1} \\ y_{1} = x_{1} \end{matrix}

(2)

Tissue impedance model

The impedance of the tissue is formulated based on the results of research on the bioimpedance of soft tissue under compression (Dodde, 2011). The tissue impedance can be represented by a Cole–Cole model as shown in Figure 3, where $R_{ext}$ is the resistance of the extra-cellular fluid, $R_{int}$ is the resistance of the intracellular fluid and $C_{m}$ is the bulk capacitance of the cellular membranes. It is important to note that the tissue impedance is related to geometrical dimensions and compression ratio of the tissue. However, the tissue manipulated by robotic gripper has almost the same geometrical dimensions. Therefore the $R_{ext}$ , $R_{int}$ and $C_{m}$ could hardly been changed by them. In other words, the $R_{ext}$ , $R_{int}$ and $C_{m}$ change non-linearly with respect to the compression ratio of the tissue in this paper.

Figure 3.

Soft-tissue impedance model.

The total impedance of the tissue can be written as:

Z = {(R_{ext}^{- 1} + {(R_{int} + \frac{1}{j ω C_{m}})}^{- 1})}^{- 1}

(3)

where $ω$ is the circular frequency of the electricity applied to the tissue.

In order to characterize the non-linearity of $R_{ext}$ , $R_{int}$ and $C_{m}$ , experimental data from Dodde (2011) is used. Due to the lack of any mathematical equations corresponding to each circuit component, the method of non-linear curve fitting is used to find the best-fit polynomial equations for the data presented in Dodde (2011).

Figure 4 shows the resistance or capacitance of each component with respect to compression ratio, and the maximum compression ratio in the plots above is limited to 0.8. This is because, at very high compression ratios, the tissue undergoes structural breakdown, which affect its electrical properties. Moreover, the amount of force needed to accomplish high compression ratio exceeds the capabilities of a surgical robotic grippers. Therefore, any data gathered for compression ratio above 0.8 is disregarded. The following polynomial Equations (4)–(6) are the best-fit curve for each component.

R_{ext} = - 171970 ε^{6} + 339800 ε^{5} - 243000 ε^{4} + 71810 ε^{3} - 7450 ε^{2} + 590 ε + 1230

(4)

R_{int} = 10203 ε^{3} - 6559 ε^{2} + 1980 ε + 2220

(5)

C_{m} = (42.40 ε^{3} - 21.06 ε^{2} - 2.04 ε + 8.9) \times 10^{- 12}

(6)

where $ε$ is the compression ratio of tissue and $ε = \frac{y_{1}}{L_{0}} \times 100 %$ .

Figure 4.

Best-fit curves for each component of Cole–Cole model: (a) best-fit curve for $R_{ext}$ ; (b) best-fit curve for $R_{int}$ ; (c) best-fit curve for $C_{m}$ .

Tissue temperature model

The tissue temperature model describes the voltage-to-temperature map. With the impedance $Z$ , calculated in Tissue impedance model section, the power dissipated in the tissue can now be written as:

P = Re {\frac{V^{2}}{Z}}

(7)

where $V$ is the AC voltage applied to the tissue. In this paper, we are assuming that the operating frequency is 200 kHz, which is a commonly used during electrosurgical procedures. In order to obtain the relative-accuracy temperature model, heat dissipation from tissue should be taken into account. It is given by:

q = h_{c} A_{t} | T - T_{r} |

(8)

where $q$ is the dissipated energy per unit time, $h_{c}$ is the heat transfer coefficient, $A_{t}$ is the surface area subject to convection, $T$ is the temperature of the tissue and $T_{r}$ is the ambient temperature.

Combining Equations (7) and (8), the temperature model of tissue can be written as:

\overset{\cdot}{T} = \frac{P - q}{mc} = \frac{V^{2} Re {\frac{1}{Z}} - h_{c} A_{t} | T - T_{r} |}{mc}

(9)

where $m$ is the mass of the tissue heated and $c$ is the specific heat capacity of the tissue.

Let $y_{2} = x_{3}$ , $x_{3} = T$ and $u_{2} = V^{2}$ , then Equation (9) can be rewritten in the state-space form as:

\begin{matrix} {\overset{\cdot}{x}}_{3} = \frac{- h_{c} A_{t} | x_{3} - T_{r} |}{mc} + \frac{Re {\frac{1}{Z}}}{mc} u_{2} \\ y_{2} = x_{3} \end{matrix}

(10)

With compression, impedance and temperature models accounted for, the overall system model has been obtained. It is obvious that the system is a non-linear MIMO system, where the inputs are force and voltage and the outputs are compressed displacement, impedance and temperature.

Controller design

Because precision or accuracy is very important to surgery, the key aspect of this work is to develop a controller that can successfully control surgical tool attached to the robotic end effectors to finish surgical tasks with minimal deviation and in a finite time horizon according to the compression and temperature trajectories. In order to realize the demand, a modular approach is taken for developing the controller, i.e. the tissue compression and tissue heat dissipation models are separately linearized and then LQG controllers are applied to each part.

Linearization of system models

Tissue compression model defined in Equation (2) is a second order non-linear system with a relative degree of 2. Thus, it has no zero dynamics and can be fully linearized. By utilizing the state feedback control law (Khalil and Grizzle, 2002), the linearized tissue compression model can be given in state-space form by:

\begin{matrix} \overset{\cdot}{x} = [\begin{matrix} 0 & 1 \\ 0 & 0 \end{matrix}] [\begin{matrix} x_{1} \\ x_{2} \end{matrix}] + [\begin{matrix} 0 \\ 1 \end{matrix}] v_{1} \\ y_{1} = x_{1} \end{matrix}

(11)

where $v_{1}$ is the error between the desired displacement and feedback value of tissue compression, it can be given by:

v_{1} = r_{1} - K_{1} [\begin{matrix} {\hat{x}}_{1} \\ {\hat{x}}_{2} \end{matrix}]

(12)

where $r_{1}$ is the desired displacement of tissue compression, $K_{1}$ is the state-feedback gain matrices, ${\hat{x}}_{1}$ and ${\hat{x}}_{2}$ are actual state values, $K_{1} {[\begin{matrix} {\hat{x}}_{1} & {\hat{x}}_{2} \end{matrix}]}^{T}$ is the feedback value of tissue compression.

Combining Equations (2), (11) and (12), the original input (force applied tissue) $u_{1}$ can be written as:

u_{1} = {mv}_{1} + d {\hat{x}}_{2} - \frac{s α}{β + 1} ((1 - \frac{{\hat{x}}_{1}}{L_{0}}) e^{β ({(1 - \frac{{\hat{x}}_{1}}{L_{0}})}^{2} - 1)} - \frac{1}{{(1 - \frac{{\hat{x}}_{1}}{L_{0}})}^{2}} e^{β \frac{{\hat{x}}_{1}}{L_{0} - {\hat{x}}_{1}}})

(13)

Similarly, the tissue temperature model defined in Equation (10) is a first-order non-linear system with a relative degree of 1. Therefore, it has no zero dynamics and is fully linearizable. The feedback control law is applied once again and the linearized tissue temperature model can be given in state-space form:

\begin{matrix} {\overset{\cdot}{x}}_{3} = v_{2} \\ y_{2} = x_{3} \end{matrix}

(14)

where $v_{2}$ is the error between the desired and feedback values of tissue temperature, given by:

v_{2} = r_{2} - K_{2} {\hat{x}}_{3}

(15)

where $r_{2}$ is the desired value of tissue temperature, $K_{2}$ is the state-feedback gain, ${\hat{x}}_{3}$ is actual state value and $K_{2} {\hat{x}}_{3}$ is the feedback value of tissue temperature.

Combining Equations (10), (14) and (15), the original input (voltage applied tissue) $u_{2}$ can be written as:

u_{2} = \frac{mc}{Re {\frac{1}{Z}}} v_{2} + \frac{h_{c} A_{t} | {\hat{x}}_{3} - T_{r} |}{Re {\frac{1}{Z}}}

(16)

Optimum controller design based on LQG

The system models have been linearized in Linearization of system models section. Now the LQG controlling strategy can be used for designing the controllers applied to tissue compression and tissue heating models separately. LQG regulator is formed by connecting the Kalman filter or estimator and the optimal state-feedback gain designed with LQR. The Kalman estimator can not only filter the disturbance and noise from inputs and measurements but can also estimate the system state values according to inputs and outputs. Therefore, controllers based on LQG are dynamic and rely on noisy output measurements to generate controlling commands. Combining the linearized system models and LQG design method (Hespanha, 2014), the control scheme based on LQG for the overall system model can be given, as shown in Figure 5, where $r_{1}$ is the input of the desired compression displacement trajectory, $r_{2}$ is the input of the desired temperature reference trajectory, $d_{1}$ and $d_{2}$ are disturbances from inputs, $K_{1}$ and $K_{2}$ are the optimal state-feedback gains, $n_{1}$ and $n_{2}$ are noises from measurements, ${\hat{y}}_{1}$ is the measured compression displacement (on which disturbance and noise have been filtered), ${\hat{y}}_{2}$ is the measured temperature (on which disturbance and noise have been filtered), and ${\hat{x}}_{1}$ , ${\hat{x}}_{2}$ and ${\hat{x}}_{3}$ are the estimated actual state values of the systems.

Figure 5.

Cartoon control scheme based on linear-quadratic Gaussian (LQG) for the overall system model.

Determination of optimal weight matrices ( Q and R ) for LQG based on PSO

With the control strategy in this paper, the most important part is to find the optimal state-feedback gains ( $K_{1}$ and $K_{2}$ ) for tissue compression and heating. However, before proceeding to obtain the optimal feedback gain, an LQR controller should first be implemented. An LQR is an optimal controller minimizes its cost function. Assuming a linear system is represented by Equation (17) with state-space form, the cost function can be given in Equation (18).

\begin{matrix} \overset{\cdot}{x} = Ax + Bu \\ y = Cx + Du \end{matrix}

(17)

J = \int_{0}^{\infty} (x^{T} Qx + u^{T} Ru) d t

(18)

where $Q$ is symmetric positive semi-definite matrix, $R$ is a positive definite matrix.

According to the LQR law (Hespanha, 2014), the optimal state feedback can be calculated:

u = - Kx

(19)

where $K$ is the state-feedback gain and it can be given as:

K = R^{- 1} B^{T} P

(20)

$P$ is obtained by solving the algebraic Riccati Equation (21) that is derived from Equations (17) and (18).

PA + A^{T} P - {PBR}^{- 1} B^{T} P = - Q

(21)

In order to choose the appropriate weight matrices ( $Q$ and $R$ ) used in Equation (21) to solve $P$ and eventually calculate the $K$ , the PSO algorithm is used in this work. The concept of the PSO is developed according to the social behaviour of a swarm of fishes, birds, bees and other animals. It is a robust population-based search algorithm inspired by the performance of a swarm (Mobayen et al., 2011). Each individual within the swarm is represented by a particle in search space. This particle has one assigned position $x_{i}$ (solution) and a vector $v_{i}$ , which determines the next movement of the particles, called the velocity vector. All particles simultaneously search for the global optimal position by updating the velocity and position of each particle by iteration. In the iteration process, each particle placed in the problem space can find a fitness value evaluated by an objective function to be optimized at its current location. The particles in a local neighbourhood share memories of their best passing positions. These positions are used to adjust the velocities and positions of all particles based on Equations (22) and (23) (Mobayen et al., 2011; Mohd et al., 2012).

v_{i} (t + 1) = w (t) v_{i} (t) + c_{1} ϕ_{1} (p_{i} - x_{i} (t)) + c_{2} ϕ_{2} (p_{g} - x_{i} (t))

(22)

x_{i} (t + 1) = x_{i} (t) + v_{i} (t + 1)

(23)

where $w (t)$ is the inertia weight, $c_{1}$ , $c_{2}$ are positive acceleration constants, $ϕ_{1}$ , $ϕ_{2}$ are random variables with uniform distribution between 0 and 1, $p_{i}$ is the best position that the particle ever passed, and $p_{g}$ is the neighbourhood position that has resulted in the best fitness so far.

$Q$ and $R$ are the entries that form a particle in the PSO algorithm. In order to simplify the problem and make the elements of matrices have actual meaning, the $Q$ and $R$ matrices are chosen to be diagonal matrices. Due to the accuracy is the most important performance for surgery, the maximum absolute error ( $MAE = max {| y_{d} (t) - y_{a} (t) |}$ ) between the desired and actual values of outputs is selected to be the fitness value calculation (objective function). $y_{d} (t)$ is the desired displacement of tissue compression or desired value of tissue temperature. $y_{a} (t)$ is the actual displacement of tissue compression or actual value of tissue temperature, which is determined by the controller or the optimal state-feedback gains ( $K_{1}$ and $K_{2}$ ). According to Equations (18) and (19), it is obvious that the optimal state-feedback gains ( $K_{1}$ and $K_{2}$ ) is related to $Q$ and $R$ matrices, and it is important that $Q > 0$ and $R \geq 0$ must be satisfied during searching the solutions. According to the above contents, the calculation process of the weight matrices based on the PSO can be described with a flow chart in Figure 6.

Figure 6.

Flow chart of the calculation process of $Q$ and $R$ based on the particle swarm optimization (PSO).

Results

In order to verify the performance of the LQG controller implemented in this paper, giving the desired compression reference trajectory in Figure 7(b), the temperature reference trajectory (thermo coagulation is usually expected at temperature of above 60–70 °C (Mandel et al., 2013)) in Figure 9(b), system parameters of in vitro pig model that are taken from the studies listed in Table 1, the desired maximum absolute error less than 0.0001 m for the tissue compression process and less than 1°C for the heating process, the solutions of $Q$ and $R$ are solved and the entire controlling process is simulated with MATLAB.

Figure 7.

Compression trajectory tracking: (a) compression trajectory tracking unfiltered while using the linear-quadratic Gaussian (LQG) controller with optimal weighting matrix obtained by a particle swarm optimization (PSO) algorithm; (b) compression trajectory tracking after filtering.

Table 1.

Parameters for models and the controllers.

Symbol [units]	Value
$d$ [N·s/m]	0.285 (Kerdok, 2006)
$m$ [kg]	0.001302 (Dodde, 2011)
$s$ [m²]	$7.15 \times 10^{- 5}$ (Rosen et al., 2008)
$α$ [Pa]	3289.4 (Rosen et al., 2008)
$β$	20.63 (Rosen et al., 2008)
$L_{0}$ [m]	0.0124 (Dodde, 2011)
$A_{t}$ [m²]	$1.43 \times 10^{- 4}$ (Rosen et al., 2008)
$c$ [J/(kg·°C)]	3600 (Mandel et al., 2013)
$h_{c}$ [W/(m²·°C)]	10 (Mandel et al., 2013)
$T_{r}$ [°C]	25 (Mandel et al., 2013)

Using the PSO algorithm, the optimal weighting matrices for the tissue compression and heating are obtained, as shown in Table 2. In order to compare the performance of the approaches of selecting weighting matrices, the weighting matrices have also been chosen through a cumbersome process of trial and error.

Table 2.

Linear-quadratic regulator (LQR) weighting matrices.

Approach	Tissue compression		Tissue heating
Approach	$Q_{1}$	$R_{1}$	$Q_{2}$	$R_{2}$
PSO	$diag [3.777 2.944 \times 10^{- 6}]$	$1.027 \times 10^{- 7}$	$0.3012$	$6.248 \times 10^{- 4}$
Trial and error 1	$diag [1.0 1.0]$	$1.0$	$1.0$	$1.0$
Trial and error 2	$diag [1.0 1.0]$	$1.0 \times 10^{- 5}$	$1.0$	$1.0 \times 10^{- 4}$
Trial and error 3	$diag [1.0 1.0 \times 10^{- 5}]$	$1.0$	$0.1$	$1.0$
Trial and error 4	$diag [1.0 1.0 \times 10^{- 5}]$	$1.0 \times 10^{- 5}$	$0.1$	$1.0 \times 10^{- 4}$
Trial and error 5	$diag [10.0 1.0 \times 10^{- 5}]$	$1.0 \times 10^{- 5}$	$0.5$	$1.0 \times 10^{- 4}$
Trial and error 6	$diag [5.0 1.0 \times 10^{- 5}]$	$1.0 \times 10^{- 5}$	$0.25$	$1.0 \times 10^{- 4}$

PSO, particle swarm optimization.

The performance of the controller in tracking the compression reference trajectory is shown in Figure 7(b). The actual trajectory converges to the desired tissue compression reference trajectory with negligible deviation (the max absolute error was $9.6057 \times 10^{- 5}$ m). It can also be clearly seen from the Figure 7(b) and Table 3 that the $Q_{1}$ and $R_{1}$ selecting through PSO can provide better controlling performance for the tissue compression process. Moreover, Figures 7(a) and (b) show the disturbance and noise can be filtered effectively. The force (control signal) applied to the tissue to achieve desired compression is shown in Figure 8.

Table 3.

Max absolute error between the desired and actual values.

Approach	Max absolute error between the desired and actual values
Approach	Displacement of tissue compression [m]	Tissue temperature [°C]
PSO	$9.6057 \times 10^{- 5}$	$0.4561$
Trial and error 1	$2.3972 \times 10^{- 3}$	$9.8894$
Trial and error 2	$1.9573 \times 10^{- 3}$	$4.4724$
Trial and error 3	$2.3931 \times 10^{- 3}$	$24.0026$
Trial and error 4	$2.8189 \times 10^{- 4}$	$1.4142$
Trial and error 5	$4.2095 \times 10^{- 4}$	$3.1628$
Trial and error 6	$2.3684 \times 10^{- 4}$	$1.0005$

PSO, particle swarm optimization.

Figure 8.

Force applied tissue while using the linear-quadratic Gaussian (LQG) controller with optimal weighting matrix obtained by a particle swarm optimization (PSO) algorithm.

Figure 9 shows the results of tracking the temperature trajectory. It can be seen that the controller was able to track the reference with almost no deviation (the max absolute error was 0.4561°C) and it also converges to the desired temperature in finite time. It can also be seen clearly from the Figure 9(b) and Table 3 that the $Q_{2}$ and $R_{2}$ selected using PSO method provided better controlling performance for the tissue heating process. Figure 10 shows the voltage (control signal) applied to the tissue to achieve the desired temperature.

Figure 9.

Temperature tracking: (a) temperature tracking unfiltered while using the linear-quadratic Gaussian (LQG) controller with optimal weighting matrix obtained by a particle swarm optimization (PSO) algorithm; (b) temperature tracking after filtering.

Figure 10.

Voltage applied tissue while using the linear-quadratic Gaussian (LQG) controller with optimal weighting matrix obtained by a particle swarm optimization (PSO) algorithm.

Conclusion

In this paper, a non-linear MIMO model representing the variation of tissue compression and conductance under applied force was used to automate the electrosurgical procedure simultaneously. The MIMO model with the two inputs (force and voltage) and the two outputs (compression and temperature) were separately linearized using state feedback linearization techniques. Two separate LQG controllers were then used to track the ramp and hold reference trajectories. During the process of designing the controller, the PSO algorithm was used to choose the optimal weighting matrices for LQG regulator to ensure the maximum absolute error was less than the desired value. The approach to select the weighting matrices solved the time-consuming and cumbersome problem of the traditional trial and error method. Moreover, the PSO method can provide better controller performance for the MIMO system. MATLAB simulation results also showed that the controllers successfully tracked the trajectories with minimal deviation and in a finite time horizon, and the disturbances and noises also can be filtered effectively. Therefore, the controller for the MIMO system developed in this paper can be used to drive a surgical robot autonomously and perform electrosurgical procedures such as coaptic vessel closures.

The next stage of this work is to model the thermal property of the tissue more accurately, which would require in vivo experimentation on porcine tissue. Moreover, in order to test the robustness of the system, the simulation needs to take other possible faults into account. In this work, the dynamics of the operating robot was not considered, which will be part of our future work. Eventually the controller will be implemented in a robotic surgical system like the RAVEN II open-source surgical robot (Hannaford et al., 2013).

Footnotes

Declaration of Conflicting Interests

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The authors disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work was supported by the RPI Office of Research Seed Grant.

References

Belmont

Dodde

Shih A

(2013) Impedance of tissue-mimicking phantom material under compression. Journal of Electrical Bioimpedance 4(1): 2–12.

Berjano

(2006) Theoretical modeling for radiofrequency ablation: state-of-the-art and challenges for the future. Biomedical Engineering Online 5(1): 24.

Das

Pan

Halder

. (2013) LQR based improved discrete PID controller design via optimum selection of weighting matrices using fractional order integral performance index. Applied Mathematical Modelling 37(6): 4253–4268.

Dodde

(2011) Bioimpedance of Soft Tissue Under Compression and Applications to Electrosurgery. Ann Arbor, MI: University of Michigan.

Feng

Pan

. (2012) Development of a medical robot system for minimally invasive surgery. The International Journal of Medical Robotics and Computer Assisted Surgery 8(1): 85–96.

Hannaford

Rosen

Friedman

. (2013) Raven-II: an open platform for surgical robotics research. IEEE Transactions on Biomedical Engineering 60(4): 954–959.

Hespanha

(2014) Lecture notes on LQR/LQG controller design. http://www.uz.zgora.pl/∼wpaszke/materialy/kss/lqrnotes.pdf (accessed November 2014).

Kanno

Haraguchi

Tadano

. (2013) Design of a 4-DOF forceps manipulator for robotic surgery. System Integration (SII). 2013 IEEE/SICE International Symposium on IEEE 2013, pp. 778–783.

Kerdok

(2006) Characterizing the Nonlinear Mechanical Response of Liver to Surgical Manipulation. Boston, MA: Harvard University.

10.

Khalil

Grizzle

(2002) Nonlinear Systems. Upper Saddle River, NJ: Prentice Hall.

11.

Kim

(2008) An ancient soft tissue characterization algorithm from in vivo indentation experiments for medical simulation. The International Journal of Medical Robotics and Computer Assisted Surgery 4(3): 277–285.

12.

Lewis

Vrabie

Syrmos

(2012) Optimal Control. Hoboken, NJ: John Wiley & Sons.

13.

Low

Phee

(2006) A review of master–slave robotic systems for surgery. International Journal of Humanoid Robotics 3(4): 547–567.

14.

Mandel

Malki

Adawi

. (2013) Hemorrhage control of liver injury by short electrical pulses. PloS One 8(1): e49852.

15.

Meng

Chen

(2014) Research on optimal control for the vehicle suspension based on the simulated annealing algorithm. Journal of Applied Mathematics.

16.

Mobayen

Rabiei

Moradi

. (2011) Linear quadratic optimal control system design using particle swarm optimization algorithm. International Journal of the Physical Sciences 6(30): 6958–6966.

17.

Mohd

NAS

Selamat

Alimin

(2012) Optimal controller design for a railway vehicle suspension system using particle swarm optimization. Jurnal Teknologi 54(1): 71–84.

18.

Nezhat

Saberi

Shahmohamady

. (2005) Robotic-assisted laparoscopy in gynecological surgery. Journal of the Society of Laparoendoscopic Surgeons 10(3): 317–320.

19.

Patil

Alterovitz

(2010) Toward automated tissue retraction in robot-assisted surgery. IEEE International Conference on Robotics and Automation (ICRA) 2010, pp. 2088–2094.

20.

Rosen

Brown

. (2008) Biomechanical properties of abdominal organs in vivo and postmortem under compression loads. Journal of Biomechanical Engineering 130: (021020–1)-(021020–17).

21.

Rowley

Batten

(2014) Dynamic and closed-loop control. http://cwrowley.princeton.edu/papers/RowleyBatten09.pdf (accessed November 2014).

22.

Sankaranarayanan

Arikatla

(2011) A simulation framework for tool tissue interactions in robotic surgery. Studies in Health Technology and Informatics 173: 440–444.

23.

Sinha

Sankaranarayanan

(2014) Modeling and control of tissue compression and temperature for automation in robot-assisted surgery. 36th Annual International Conference of the Engineering in Medicine and Biology Society (EMBC), pp. 366–370.

24.

Sommer

Eder

Kovacs

. (2013) Multiaxial mechanical properties and constitutive modeling of human adipose tissue: a basis for preoperative simulations in plastic and reconstructive surgery. Acta Biomaterialia 9(11): 9036–9048.

25.

(2008) Control Methods for Automated Surgery. Seattle, WA: University of Washington.

26.

Chizeck

(2008a) Equilibrium selection in automated surgery. 2nd IEEE RAS & EMBS International Conference on Biomedical Robotics and Biomechatronics, BioRob 2008, pp. 740–746.

27.

Chizeck

(2008b) Lower bounds on the optimal control of soft tissue grasping. 30th Annual International Conference of the Engineering in Medicine and Biology Society, EMBS 2008, pp. 1943–1947.

28.

Chizeck

Hannaford

(2007) Comparison of transient performance in the control of soft tissue grasping. IEEE/RSJ International Conference on Intelligent Robots and Systems, IROS 2007, pp. 1809–1814.