Introduction

Irreversible electroporation (IRE) is a new minimally invasive surgical technique to ablate undesirable tissue (1). The technique is easy to apply, can be monitored and controlled, is not affected by local blood flow, and does not require the use of adjuvant drugs. The procedure involves placing electrodes into or around the targeted area to deliver a series of short and intense electric pulses that induce irreversible structural changes in the cell membranes of the targeted tissue through a phenomenon known as irreversible electroporation (2, 3), causing cell death.

Electroporation, which results in an increase in the permeabilization of the cell membrane, can be initiated by exposing cells or tissues to electric pulses (4, 5). As a function of the resulting transmembrane potential (the electric potential difference across the plasma membrane) the electroporation pulse can either: have no effect on the cell membrane, reversibly open the cell membrane after which cells can survive, or irreversibly open the cell membrane, which leads to cell death. This potential is normally on the order of 70mV. If the potential drop across the membrane exceeds approximately 1V, structural rearrangement of the lipid bilayer occurs, creating permanent aqueous pathways or pores for ions and macromolecules to pass through, i.e., irreversible electroporation (6).

This induced potential is dependent on a variety of conditions such as tissue type and cell size. There are also a number of pulse parameters that affect the resulting transmembrane potential and, thus, the degree of electroporation, such as the pulse shape, duration, number, and repetition rate. Readers interested in electroporation at the single cell level are referred to (7, 8, 9, 10). Weaver et al. illustrate how an applied electric field in tissue affects the transmembrane potential of its cells and the equivalent circuit of tissue undergoing electroporation (11). The microscopic electric field distribution is difficult to determine analytically due to tissue anisotropy, the locations of large blood vessels, electrode geometry and orientations, and the distribution of electrical impedance within the tissue. The issue is further complicated because the bioelectrical properties will change due to minor localized heating of the tissue from the pulses and due to the cells becoming permeable as a result of treatment (12, 13). However, for a specific tissue type, the macroscopic electric field distribution to which the tissue is exposed is the primary parameter (14). Therefore, knowledge of the electric field distribution is critical in optimizing pulse parameters for in vivo applications.

Reversible electroporation of tissue has been used in conjunction with membrane impermeant anticancer drugs for almost twenty years as a form of cancer therapy known as electrochemotherapy (ECT) (15, 16, 17). Once the non-permeable substance is injected into the body, electrodes placed into or around the targeted tissue generate reversibly permeabilizing electric fields throughout the targeted tissue to facilitate the entry of the drug into the cells of the targeted area (18) to enhance the effectiveness of the anticancer drug. The cell membrane subsequently reseals itself with the drug trapped inside the cell. Typical parameters for electrochemotherapy are 4-8 100 μs pulses at a frequency of 1Hz with a voltage-to-distance ratio of 1300V/cm (19). ECT benefits from the ease of application and the area treated is comparable to both high temperature treatment therapies and non-selective chemical therapies. In addition, because the electric field is not strongly affected by local blood flow, control over the extent of the affected tissue is possible, unlike in many other diffuse modes such as thermal and non-selective chemical therapies. ECT is beneficial due to its selectivity; however, a disadvantage is that by its nature, it requires the combination of chemical agents with an electric field to have the desired effect.

Researchers in electrochemotherapy have demonstrated the importance of the electric field distribution to induce reversible electroporation to enhance the effectiveness of drug delivery (20). Miklavcic et al. (1998) recognized the importance of the electric field distribution when designing electroporation protocols and showed that minimally invasive needle electrodes produce an inherently inhomogeneous field that must be taken into consideration when planning a treatment. By comparing their models with experimental data, they have found that for a specific tissue type and a specified set of voltage parameters, there is a strong correlation between the electric field and the tissue?s permeabilization state (21). They estimated the electric field threshold for reversible electroporation as 362 +/- 21 V/cm and the threshold for irreversible electroporation as 637 +/- 43 V/cm for rat liver tissue using eight 100 μs pulses at a frequency of 1 Hz. Such studies are valuable for determining the electric field thresholds to induce reversible and irreversible electroporation for a specific tissue and set of voltage parameters. The culmination of the research by Miklavcic et al. illustrates the effectiveness of modeling to select appropriate pulse parameters and electrode geometries to successfully treat a targeted area for ECT (21).

Recently, Davalos et al. developed models to study the potential thermal effects associated with reversible electroporation. They analytically investigated the effects of various parameters on the resultant temperature distribution resulting from reversible electroporation pulses (22). They have shown what parameters are optimal to ensure negligible heating of the tissue. For example, electrodes with smaller diameters result in higher and less uniform distribution of electric fields. Other researchers showed that the temperature increase and the thermal damage is typically insignificant unless the procedure is conducted over the skin (23, 24). When designing protocols for transdermal drug delivery via electroporation, they showed that it is important to recognize that a large amount of the voltage drop occurs over the stratum corneum (25, 26). When the electrodes are placed over the skin, significant heating can occur due to the high electrical resistance of the skin, unless active cooling is applied (23). Furthermore, improperly prepared skin can drastically reduce the actual electric field delivered within the tissue when electrodes are applied at the body surface.

Similarly, in order to ablate only the targeted tissue with irreversible electroporation, there are specific parameters that must be chosen. Through modeling and simulations, one can accurately predict the necessary electrode voltage as well as optimize the electrode placements, pulse patterns, timing, et cetera, for each unique treatment to effectively treat the targeted area. Previous models have shown that IRE can ablate substantial volumes of tissue without inducing any relevant thermal effects (1, 2, 3). The research showed that the most critical parameters affecting the extents of ablation are impedance distribution, pulse characteristics, and electrode configuration (size, spacing, and number) (1).

The primary goal of this paper is to provide interested physicians and researchers with the tools and understanding necessary to design appropriate protocols to use of IRE in cancer therapy. To this end, we first present the fundamental theory that determines how electric fields will result from a chosen electrode configuration, pulse characteristics, and the electrical properties of the tissue. We then proceed to describe the other important theoretical paradigm for modeling IRE, heat diffusion via the Pennes bioheat equation, and we detail the connection between them, namely Joule heating.

To provide an intuitive understanding of how to design an effective IRE procedure, we then proceed to discuss in detail two important cases that come up frequently in IRE treatments: coaxial disk electrodes and parallel needle electrodes. To expound upon these fundamental 2-electrode cases across a range of relevant electrode sizes and inter-electrode gaps, we provide the reader with a functional set of figures into which they can plug their scenario specific parameters. We then discuss the use of cell constants to predict power consumption and measure bulk conductivity, explain how to predict and control the size of the treated region, and provide the tools necessary to estimate the extents, if any, of thermal damage at various points within the treated tissue. We also discuss various complexities of IRE such as heterogeneous tissues and the changes in tissue conductivity that result from electroporation.

The synthesis of these various aspects will provide physicians with a framework to design IRE protocols for their particular application. If the physician knows the electrode configuration and geometry as well as the electrical impedance distribution of the tissue, pulse parameters can be designed to ensure ablation with irreversible electroporation throughout the entire targeted region.

Materials and Methods

In order to design protocols for an IRE procedure, the electric field distribution must be determined, which is dependent on the procedure?s specific electrode geometry and tissue impedance distribution. By predicting the electric field distribution for a specific scenario, the physician can optimize the electrode geometry to ablate the entire targeted region while minimally affecting the surrounding tissue. Furthermore, to verify that a specific protocol does not induce thermal effects, the temperature distribution can be calculated from the electric field distribution, the electric pulse parameters and tissue properties. The methodology below illustrates how to obtain this information and presents the results through two basic two-dimensional models: coaxial disk electrodes and parallel needle electrodes.

Estimating the Electric Field Distribution Due to Electroporation

The electrical field distribution associated with an electroporation pulse is determined by solving the Laplace equation:



where φ is the electrical potential and σ is the electrical conductivity of the tissue. The solution of the Laplace equation yields the distribution of the potential within the tissue. The electric field is then the gradient of the potential, i.e., E = ∇φ.

The typical formula to approximate the induced transmembrane potential (V_m) resulting from an applied electric field, which determines whether a cell will undergo electroporation, is:



where λ is the shape factor of the cell (1.5 for spherical cells), r is the radius of the cell, E_a is the applied electric field, θ is the angle between electric field and the vector from the cell center to any point on its surface, f_s is approximately equal to the frequency where the beta dielectric dispersion occurs (below which the cell membrane charge is in step with the electric field), and f is the frequency of the assumed sinusoidal E_a (27). Since even a DC pulse can be decomposed into a superposition of various sine waves at different frequencies, this equation enables the transmembrane potential resulting from a short pulse to be computed in an approximate sense.

For most cases, the transient terms can be neglected because the electroporation pulse (100 μs-50 ms) is much larger than the membrane charging time (about 1 μs) (28). The equation above reveals that if the pulses are ultrashort (i.e., less than several microseconds), the problem is no longer static and the charging time of the membrane must be considered. Under such circumstances, it is hypothesized that the lipid bilayer membrane is not electroporated but that intracellular organelles are affected by the pulsed electric field (29).

Estimating the Temperature Distribution Due to Electroporation

Solving the Laplace equation enables one to calculate the associated Joule heating (p), the heat generation rate per unit volume from an electric field:



The heating of the tissue resulting from each procedure was calculated by adding the Joule heating source term to the Pennes Bioheat transfer equation (30). The Pennes Bioheat equation, which accounts for metabolism and blood flow, is the most commonly used equation to solve heat transfer problems in the body. Our modified Pennes Bioheat equation has the following form:



where k is the thermal conductivity of the tissue, T is the temperature, w_b is the blood perfusion, c_b is the heat capacity of the blood, T_a is the arterial temperature, q''' is the metabolic heat generation, ρ is the tissue density, and c_p is the heat capacity of the tissue.

For software that cannot handle the blood perfusion term in Equation [4], a variable can be substituted to circumvent this limitation. Let:



Therefore, the transformation is:



After canceling out like-terms and dividing through by the exponential, the equation reduces to:



The equation above resembles the general heat conduction equation, which is readily solvable with most finite element packages, with the term as the heat source The same variable substitution is applied to the boundary conditions to solve the problem. The solution, U, is then substituted back into [5] to determine the temperature distribution in the analyzed domain.

Case Studies

Two fundamental models that are representative of common conditions associated with irreversible electroporation are used to illustrate the critical parameters in designing treatments. These models are depicted in Figure 1. The first configuration involves two coaxial disk electrodes and the second contains two parallel needle electrodes. These are common electrode geometries and many more geometries are of course possible, yet through these two basic examples, we illustrate the trends that occur as a function of critical parameters such as electrode size and shape.

The models used to generate the results in Figures 3, 4, and 6 are two-dimensional, as depicted in Figure 2. For the disk electrodes, as long as conductivity is uniform, no loss of generality is incurred by using a cylindrically symmetric model. However, the 2D simplification for the needle electrodes implies infinitely long electrodes. As will be discussed in the following sections, this causes the cell constant, which can be used to predict power consumption or to estimate bulk tissue conductivity, to be overestimated, the treatment volume to be underestimated, and has minimal impact on the predicted maximal temperature rise due to the pulse. These errors become substantial when the ratio of electrode length (L) to electrode gap (d) is small, yet it will be shown that the 2D model is still a useful approximation for the purpose of IRE treatment planning.

Boundary Conditions

The models are fully defined and readily solvable using a numerical method once an appropriate set of boundary conditions and the properties of the tissue are defined. Boundary conditions most often include surfaces where electric potential is specified, as in the case of a source or sink electrode, or surfaces that are electrically insulating, as on the free surfaces of the tissue, for example. For each electrode configuration, the surface of one electrode is assumed to have a prescribed voltage, and the other electrode is set to ground. Specifically, the boundary condition of the tissue that is in contact with a charged electrode is defined as:



where V_o is the applied voltage and the electrical boundary condition at the interface of a grounded electrode and the tissue is defined as:



The boundaries where the analyzed domain is not in contact with an electrode are treated as electrically insulating:



Other possibilities for electrode boundary conditions include those for passive electrodes, where a film conductance is specified, or those where the electrode interface impedance is included in an approximate manner, as in the complete electrode model (31).

The computation of φ and, hence, the electric field can then be used to estimate the variation of temperature, where the tissue is initially at physiological temperature (37 °C). Several thermal boundary conditions can be explored to study their heat exchange between the electrodes and the tissue (1, 2, 22, 24); however, in our study, the boundary of the analyzed domain and the surfaces of the electrodes are taken to be adiabatic to produce an upper limit to the calculated temperature distribution in the tissue:



Numerical Modeling

In this study, the computations were performed with a commercial finite element package (FEMLab, Comsol AS, Stockholm, Sweden), whereas data analysis was done in Matlab (The MathWorks, Inc., Natick, MA). The finite element method generates approximate solutions that consist of the value of electric potential on every node within a mesh that fills the geometric confines of the tissue in question. Analytical solution of electrical conduction on a non-trivial geometry usually requires the use of techniques of numerical analysis, such as the finite difference method or the finite element method. Even where analytical solutions exist, it is often difficult to extract accurate predictions due to the frequent occurrence of infinite integrals that are difficult to evaluate (even numerically). Furthermore, analytical solutions typically fail when the tissue conductivity is heterogeneous. In our study, the analyzed domain extends far enough (typically at least ten times the electrode gap for the case of needle electrodes and five times for disk electrodes) from the area of interest (i.e., near the electrodes) that the electrically and thermally insulating boundary at the edges of the domain does not significantly influence the results in the treatment zone.

Results and Discussion

Example of Electric Field Distribution for Each Electrode Configuration

Figure 1 depicts the electric field that was computed for the two cases with a scalar and uniform electrical conductivity across the tissue, and Figure 2 displays six additional representative cases in 2D. The first case (Fig. 1a) was computed for a pair of coaxial disk electrodes with an electrode-electrode gap (d) of 4 mm and a disk diameter (D) of 10 mm, and where a pulse of 400 V was applied across the electrodes. The second case (Fig. 1b) was computed for a pair of 1-mm diameter (D) parallel needle electrodes with a center-to-center distance (d) of 10 mm, and a pulse of 1000 V was applied between the two needles. In each case, the voltage-to-distance ratio is 1000 V/cm, but it is important to notice that the actual electric field that develops in the area between the electrodes varies strongly with position, especially in the case of needle electrodes. In fact, the regions with the largest electric fields (1600 V/cm) occur as the pair of thin red annular rings in the first case and the thin red covering over the needle electrodes in the second case. These high fields are due to the sharp electrode borders in the first case and the electrode geometry in the second case. It is also apparent from the cases reported in Figure 2 that electrode diameter greatly affects the region of tissue that will experience any particular electric field strength (larger electrodes mean larger treatment volumes if all other parameters are fixed).


Figure 1: Two commonly employed electrode designs ? (a) coaxial disk electrodes and (b) parallel needle electrodes. Full 3D solution is depicted for each case, where the aspect ratio (d/D = interelectrode gap/electrode diameter) is 0.4 for (a) and 10 for (b), and where the applied voltage to distance ratio (V/d) is 1000 V/cm.

The results presented in Figures 1 and 2 are for a specific set of conditions. Figures 3, 4, and 6 provide vital statistics and trends for the possible variations on these two basic designs.

Cell Constants, Power Consumption, and Bulk Conductivity

The cell constant (K) is numerically equal to the product of the bulk conductivity (in S/m) of the tissue (σ) and the magnitude of electrical impedance (Z, in Ohms) between the electrodes, and has units of inverse length. Figure 3 presents the cell constants (K) computed for the two electrode configurations over a range of d/D values typically encountered. To obtain these data, the applied current was calculated by integrating the normal current density over the surface of one electrode. Most importantly, K provides: (i) an estimate of σ when Z has been measured and (ii) an estimate of the power consumption when σ is known ahead of time (P = σV²/K).

For example, consider the two sets of cases presented in Figure 2. For the disk electrode case where d/D = 0.4 (Fig. 2c), KD ∼ 0.37 (K = 37 m^-1) from Figure 3a. Similarly, for the needle electrode case where d/D = 10 (Fig. 2e), KL ∼ 0.96 (K = 96 m^-1 for L = 10 mm). If we assume σ ∼ 0.2 S/m, then the power during the pulse will be around 0.86 kW (P ∼ 0.2·400²/37) for the disk electrode case in Figure 2c and about 2.1 kW (P ∼ 0.2·1000²/96) for the needle electrode case in Figure 2e. For the needle electrode configuration, however, this 2D simplification has resulted in a lower than actual rate of power consumption during the pulse being predicted. The full 3D model predicts that the pulse will in fact consume 2.8 kW, so that K = 72 m^-1. This discrepancy can be understood by considering that the impedance in the circuit falls when current can flow not only within the 10 mm thick section of tissue that bounds the electrodes, but also above and below it; since adding impedances in parallel will reduce the overall impedance (Z), the 3D model will always have higher power consumption and a lower cell constant than would be predicted from the 2D model. This effect becomes more pronounced as L/d is reduced.


Figure 2: Cross-sectional electric field distribution for some representative cases. Graphs show simplified 2D models where the same surfaces of constant electric field from Figure 1 are represented with contours. Electrode diameters are 2.5, 5, and 10 mm for (a)-(c) and 0.5, 1, and 2 mm for (d)-(f), respectively. All other parameters are identical with those listed for the corresponding graphs of Figure 1, so that (c) and (e) are for the same two cases as shown in Figure 1. For the disk electrodes in (a)-(c), cylindrical symmetry ensures the accuracy of the 2D model; however, the 2D simplification for the needle electrodes in (d)-(f), where the electrodes are perpendicular to the plane of the model, implies infinitely long electrodes.

Table I examines this discrepancy between the cell constant predicted for the case of 1 mm diameter parallel needle electrodes by the 2D model and that for the full 3D model. From these computed data, it can indeed be seen that as L increases for constant d, the 2D model becomes more accurate. Conversely, the 2D model becomes less accurate as d is increased across the three columns that correspond with the d/D points marked as circles in Figure 3b. By adjusting the data in Figure 3b according to the trends in Table I, it is possible to generate reasonably accurate predictions of the cell constant with just the 2D model.




Figure 3: Cell constant (?) as a function of aspect ratio (d/D) plotted for (a) coaxial disk electrodes and (b) parallel needle electrodes. Cell constant is non-dimensionalized in the first case by multiplying with the diameter (D) of the disk electrodes and in the second case by multiplying with the active electrode length normal to the page (L). Circles in (b) correspond with the three columns in Table I that report discrepancies between the 2D assumption (infinitely long parallel needle electrodes) and the fully 3D case.

Continued

Ecommerce Solutions
Online Payment Solutions