113
research article
A statistical model describing combined irreversible electroporation and electroporation-induced blood-brain barrier disruption
Shirley Sharabi12, Bor Kos3, David Last1, David Guez1, Dianne Daniels12, Sagi Harnof24, Yael Mardor12, Damijan Miklavcic3
1	The Advanced Technology Center, Sheba Medical Center, Ramat-Gan, Israel
2	Sackler Faculty of Medicine, Tel-Aviv University, Tel-Aviv, Israel
3	University of Ljubljana, Faculty of Electrical Engineering, Ljubljana, Slovenia
4	Neurosurgery Department, Sheba Medical Center, Ramat-Gan, Israel
Radiol Oncol 2016; 50(1): 28-38.
Received 23 October 2015 Accepted 3 January 2016
Correspondence to: Dr. Yael Mardor, Ph.D., The Advanced Technology Center, Sheba Medical Center Tel-Hashomer, 52621 Israel. Phone: +972 3 530 2993; Fax: 972 3 530 3146; E-mail: yael.mardor@sheba.health.gov.il
Disclosure: We declare that Damijan Miclavcic holds patents in the area of electroporation and is consulting for various companies with financial interest in the area of electroporation use in medicine. The other authors declare no competing financial interests.
Background. Electroporation-based therapies such as electrochemotherapy (ECT) and irreversible electroporation (IRE) are emerging as promising tools for treatment of tumors. When applied to the brain, electroporation can also induce transient blood-brain-barrier (BBB) disruption in volumes extending beyond IRE, thus enabling efficient drug penetration. The main objective of this study was to develop a statistical model predicting cell death and BBB disruption induced by electroporation. This model can be used for individual treatment planning.
Material and methods. Cell death and BBB disruption models were developed based on the Peleg-Fermi model in combination with numerical models of the electric field. The model calculates the electric field thresholds for cell kill and BBB disruption and describes the dependence on the number of treatment pulses. The model was validated using in vivo experimental data consisting of rats brains MRIs post electroporation treatments. Results. Linear regression analysis confirmed that the model described the IRE and BBB disruption volumes as a function of treatment pulses number (r2 = 0.79; p < 0.008, r2 = 0.91; p < 0.001). The results presented a strong plateau effect as the pulse number increased. The ratio between complete cell death and no cell death thresholds was relatively narrow (between 0.88-0.91) even for small numbers of pulses and depended weakly on the number of pulses. For BBB disruption, the ratio increased with the number of pulses. BBB disruption radii were on average 67% ± 11% larger than IRE volumes.
Conclusions. The statistical model can be used to describe the dependence of treatment-effects on the number of pulses independent of the experimental setup.
Key words: electroporarion; blood brain barrier; Peleg-Fermi
Introduction
Electroporation (EP) is a physical phenomenon in which electric fields make cell membranes transiently permeable to ions and macromolecules which are otherwise deprived of or have limited
trans-membrane transport mechanisms.1-3 Electric pulses applied to the tissue induce an electric field which in turn induces a change in cell membrane potential. This change depends on various tissue related parameters such as tissue type and cell size as well as pulse parameters including pulse ampli-
Radiol Oncol 2016; 50(1): 113-120.
doi:28.1515/raon-2015-0015
Sharabi S et al. / A statistical model describing brain electroporation
29
tude, shape, duration, number of pulses, and pulse repetition frequency. As a function of the induced electrical field, electric pulses can either: reversibly permeabilize the cell membrane (reversible EP) or permeabilize the cell membrane in a manner that leads to cell death (irreversible EP).4 It was recently demonstrated that when applying EP to brain tissue it also induces reversible disruption of the blood-brain barrier (BBB).5-7
Both irreversible EP (IRE), and reversible EP combined with chemotherapy, also known as electrochemotherapy (ECT), are emerging as new treatment techniques for solid tumors.38-16 ECT uses EP to allow increased uptake of chemothera-peutic drugs into tumor cells12 and IRE is a method aimed at inducing tumor ablation without thermal damage.1718 Brain tumors are excellent candidates for local EP treatment. Glioblastoma multiforme (GBM) is the most frequent and most aggressive primary brain tumor with an average survival of 14 months from diagnosis. Existing treatments offer poor prognosis for GBM mainly due to tumor infiltration into the surrounding brain, high resistance to therapeutic apoptotic stimuli and poor BBB penetration of most therapeutic agents.1920 A combined approach, consisting of inducing significant/ rapid necrosis in the tumor mass and simultaneous delivery of high chemotherapy doses to the tumor and surrounding infiltrating zone is suggested as a treatment strategy. EP-induced tissue necrosis within the massive region of the tumor and surrounding BBB disruption, enabling efficient local delivery of systemically administered chemotherapy was recently demonstrated.721
Individual treatment planning is an important key for EP-based treatment success.22 Treatment parameters should be chosen in such a manner that will induce maximal damage to the tumor while sparing surrounding healthy tissue. This is usually done by numerical models. Several numerical models describing the electric filed distribution in the tissue have been introduced, and are applied for predicting treatment outcome and planning the electrodes placement to ensure full tumor coverage by electrical fields higher than the EP threshold.23-28 These models are usually based on experimental data. Treatment volumes calculated from MRI29 or histological data3031 are incorporated into computerized models together with the organ characteristics and the electrodes configuration. These calculations traditionally use deterministic models, i.e all the cells exposed to electrical fields higher than a specific threshold, known in the literature, will be irreversibly/reversibly electropo-
rated. Nevertheless, live tissues are more complex, especially malignant tissues which are inherently inhomogeneous, and therefore assuming a statistical effect of EP parameters maybe more appropriate.3233 For this reason we chose to apply a statistical model to describe reversible/irreversible effects in vivo.
The Peleg-Fermi model is the most widely used mathematical model for describing cell death as a consequence of IRE in medicine.32-34 Although several other models have been proposed35 the Peleg-Fermi model seems the most adequate since it includes dependency on the number of pulses as well as in electrical field. For this reason we decided to apply it on our experimental data and further extend it to irreversible and reversible EP effects in vivo.
The Peleg-Fermi statistical model was first introduced as a model describing the survival of bacteria after exposure to pulsed electrical fields.36 Later on it was suggested that this model can be adapted to describe the effects of IRE.32-34 Goldberg and Rubinsky32 extrapolated experimental data obtained using prostate cancer cells and demonstrated the feasibility of applying this model to describe the effects of IRE for up to 10 treatment pulses. Garcia et al. extended the model up to 90 pulses, by theoretical analysis that is yet to be confirmed with experimental data.34
Treatment parameters such as pulse shape, amplitude, frequency, duration and number of puls-es37 38 affect treatment outcome. Here, we chose to study and model the effect of number of pulses while other pulse parameters remain fixed.
A numerical model describing electric field distribution in the brain tissue based on the applied voltage, tissue and electrodes electrical properties and electrodes configuration was constructed. The calculated electrical field was then implemented in the statistical model that was estimating the effect of the number of pulses on the outcome- irreversible damage and BBB disruption.
The first goal of our study presented below was to extend the Peleg-Fermi model to describe a wider range of the number of treatment pulses in vivo and to validate the extended model using experimental data obtained from naive rats treated with EP in the brain.
The second goal was to adapt the statistical Peleg-Fermi model to describe the effects of pulse parameters on BBB disruption. BBB disruption is a vital key in treating brain tumors since it is important to disrupt a large enough volume surrounding the tumor mass to enable efficient drug penetration
Radiol Oncol 2016; 50(1): 28-38.
30
Sharabi S et al. / A statistical model describing brain electroporation	30
into the infiltrating zone. Once established, models describing both IRE and BBB disruption can be implemented to provide a complete treatment planning for brain tumors with EP.
Materials and methods
Animal experiments
The study was approved by and performed in accordance with the guidelines of The Animal Care and Use Committee of the Sheba Medical Center, which is approved by the Israeli authorities for animal experimentation.
We have recently presented the results of an animal experiment designed to study both IRE and BBB disruption using the same experimental setup.67 Here we describe in detail the aspects relevant to our statistical model which are based on that experimental data. Our unique electrode setup employs a single insulated intracranial needle electrode with an exposed tip placed in the target tissue and an external surface electrode pressed against the skin. The electric field produced by this electrode configuration is highest at the exposed tip of the intracranial electrode tissue interface and then decays with the square of the distance. Therefore, the electric fields surrounding the needle electrode tip induce nearly spherical IRE effects at the target tissue and gradually decrease further away to reversible EP effects which induce BBB disruption. Regions of interest (ROIs) plotted on MR images acquired post EP treatments with various pulse parameters were used for calculating the tissue damage and BBB disruption radii. We then studied the correlation between the experimental radii and the extended statistical model.
Animal model and procedure
The study was performed by treating 46 male Spring Dawly rats with 50 |js monopolar electric pulses at 1 Hz and 600 V, as previously described.7 The rats were divided into seven groups of 5-7 rats each, treated with varying number of pulses (N = 10, 45, 90, 180, 270, 450 and 540).
MR imaging
Rats were scanned 30 minutes post treatment and periodically thereafter up to 2 weeks post treatment, using a 1.5 T GE Optima MR system (Optima MR450w, General Electric, Milwaukee). The MR sequences included contrast-enhanced T1-
weighted MRI for depiction of BBB disruption and T2-weighted MRI for depiction of tissue response. Gradient echo (GE) MRI was acquired to assess possible procedure-related bleeding.
The damage radius induced by IRE (rd) (in mm) for each rat was calculated from the hyper-intense regions on T2-weighted MR images acquired two weeks post treatment. This time point was previously determined by histology as adequate to describe IRE.7 BBB disruption radius (rb), referring to the maximal radius of tissue in which the BBB was breached, was calculated from enhancing regions on contrast-enhanced T1-weighted MR images acquired 30 minutes post EP treatment.
In both cases the radii were calculated by delineating ROIs over the entire enhancing region in each slice (excluding the ventricles). The number of pixels in the ROIs was then counted and multiplied by the volume of a single pixel to receive the ROI volume. The slice thickness was 2 mm and in-plane pixel size was 0.3 X 0.3 mm. Radii of each slice was then extracted by calculating the biggest radius based on the Euclidean distance transform of the corresponding slice. The biggest radii computed over all slices were chosen as IRE radius and BBB disruption radius.
The radii rd and rb where then plotted as a function of the number of treatment pulses (N) to determine the dependence of the radii on the number of treatment pulses.
Numerical modeling
The mathematical models were based on a two-dimensional finite element model (assuming spherical symmetry of the produced IRE lesions and BBB disruption) (Figure 1) that was implemented in the COMSOL software package (Comsol Multiphysics, v.4.2a; Stockholm, Sweden) as previously described.734
The rat head and chest were modeled as a 30 x 15 mm ellipse (Figure 1C) with an initial conductivity of 0.258 S/m to match the conductivity used by Sel et al.37 The electric field was described by the Laplace equation for electric potential distribution in a volume conductor:
V-(a(E)V<p) = 0	[1]
where a is the electric conductivity of the tissue, E is the applied electric field and ty is the potential. The a(E) dependence of brain tissue was described by an smoothed Heaviside function using 500 V/ cm and 700 V/cm as reversible and irreversible
Radiol Oncol 2016; 50(1): 28-38.
Sharabi S et al. / A statistical model describing brain electroporation	55
thresholds.3739 These values where used traditionally for 90 pulses and are recalculated using the Peleg-Fermi model for the different number of pulses in this paper.
Dirichlet boundary condition was applied to the surface of the electrode:
<P = %
and to the ground <p = 0
[2]
[3]
where (£>0 is the applied potential on the intracranial electrode.
The boundaries where the analyzed domain was not in contact with an electrode were treated as electrically isolative and Neumann boundary condition was set to zero on the outer border of the model:
fr = 0
[4]
V• (/cVT) + wbcb(Ta -T) + Qmet + q"
Qmef = a|Vp|2
[6]
where k is the thermal conductivity of the tissue, T is the temperature, wb is the blood perfusion, cb is the heat capacity of the blood, Ta is the arterial temperature, q''' is the metabolic heat generation, p is the tissue density, cp is the heat capacity of the tissue and is the local voltage amplitude. Qmel accounts for Joule heating, where ^ is the electrical potential and a is the electrical conductivity of the
®
▼ 2.3312x10"" ▲ 38.835
where n denotes the normal to the boundary.
Thermal modeling
Control of the temperature during EP treatments is important in order to avoid damage to unwanted regions. The goal is to achieve complete coverage of the targeted tissue with sufficiently high electric field while ensuring that the temperature increase during the procedure does not generate thermal damage. The thermal effects of EP were determined from solution of the modified Pennes' bioheat equation (equation [5]) in the 2D numerical model with the inclusion of the Joule heating source term. A duty-cycle approach was used, in which a time dependent solver for the duration of the treatment was applied and the thermal dissipation was multiplied by the pulse length.
lp dt	[5]
®
0 10 x(mm)
FIGURE 1. Simulation results. (A) Electric field distribution in the numerical model. The shape of the field assumes a nearly spherical shape. (B) Temperature distribution after 540 pulses. (C) Model geometry including the location of the electrodes (red arrows)
tissue. The initial brain temperature was set to 37°C to match human temperature although anesthe-sized rat temperature is around 32°C. The parameter values utilized in the bioheat equation were taken from the literature40 and were used by others to follow/measure temperature increase and determine possible thermal damage due to EP treat-ments.4142 All parameters used in the simulations are summarized in Table 1. The thermal properties of the silver plating and copper were taken from the Comsol Multiphysics database.
Statistical modeling
The original Peleg-Fermi model computes the ratio (S) of surviving bacteria after EP. Here we extended this model to describe the effects of EP on brain tissue as following:
Radiol Oncol 2016; 50(1): 28-38.
32
Sharabi S et al. / A statistical model describing brain electroporation	32
TABLE 1. Material properties used for numerical model
Brain	a - basic conductivity	0.258[S/m]
	k - Thermal conductivity	0.0565[W/(m*K)]
	Cp - Heat capacity	3680 [J/(kg*K)]
	p - density	1039 [kg/mA3]
	Q'"- metabolic heat generation	10437 [W/mA3]
	T - temperature	37°C
Blood	Cp-heat capacity	3840 [J/(kg*K)]
	p density	1060 [kg/mA3]
	Wb-Perfusion rate	7.15E-3 [1/s]
copper	a - basic conductivity	5.998E7 [S/m]
	k - thermal conductivity	400 [W/(m*K)]
	Cp heat capacity	385 [J/(kg*K)]
	p - Density	8700 [kg/mA3]
Silver	a - basic conductivity	6.273E7 [W/mA3]
	k - thermal conductivity	429 [W/(m*K)]
	Cp heat capacity	234 [J/(kg*K)]
	p - Density	10500 [kg/mA3]
S(E,N) =
1 + exp
E-Ec(N) A(N)
[7]
First, the model was adapted to predict tissue damage (cell death) probability induced by EP. In the Peleg-Fermi model the probability for cells survival is given by:
1
where E is the electrical field, N is the number of pulses, Ec is the critical electric field in which 50% of the cells are killed and A is a kinetic constant which defines the slope of the curve.
The electric field calculated using the numerical model was exported to Matlab (R2011a, Mathworks, USA) and was implemented in the Peleg-Fermi model.
We have previously shown that the hyper-intense regions on T2-weighted MRI obtained 14 days post treatment were significantly correlated with rarified regions in histology, confirming that these regions represent damaged tissues.7 Based on this, rd was set as S(E,N) = 0, assuming over 99.99% of the cells were irreversibly electroporated.
An optimization based on A Nelder-Mead simplex algorithm43 with added constrains was applied to Equation [7] for each group treated with N pulses, calculating a map of S(E), until r(S = 0) matched r . The coefficients E and A for the differ-
dc
ent number of pulses were extracted and behavior equations were fitted to Ec(N) and A(N).
For each group treated with N pulses, Electric field distribution (E) was calculated using COMSOL Multiphisics and extracted to Matlab. The map E, along with the equation [7] allows to associate a map of S with any pair of the Fermi distribution (Ec,A). From the S map, the two isocontours of S = 0.9999 and S = 0.0001 are fitted to circles. An optimization based on A Nelder-Mead simplex algorithm on Ec and A as variables is used to find the (Ec,A) pair of parameters best corresponding to rd/rb.
The process of extracting r from S is nonlinear as it is based on fitting S = 1 iso-contour to a circle. Therefor the global dependency between Ec, A and rb/rd is noisy. This noisiness could potentially cause problems with computation of derivatives. Additionally, since S is monotonous there is no risk of the simplex finding a local minimum.
In order to assess whether the goodness of the Ecd(N) and Ad(N) (Ec and A for IRE) fits to the experimental data, r(S = 0) for different number of treatment pulses was calculated and compared to rd.
Although the Peleg-Fermi model was originally used to describe cell death, here we adapted it to describe BBB disruption as well and calculated the relevant coefficients. For this purpose the model was fitted to the radii calculated from contrast-enhanced T1-weighted MR Images. This time rb was set as BBB(E,N) = 1, meaning less then 0.001b% of the BBB was disrupted in radii larger than rb.
After determining Ecb and Ab (Ec and A for BBB disruption) for each N and the behavior equations Ec(N) and A(N), the goodness of the fit to the experimental data was evaluated by recalculating r(BBB = 1) and fitting it to rb.
The electrical field threshold for cell kill, i.e. IRE extent and for BBB disruption, i.e. reversible EP extent, for different N values were then extracted from the results of the model and compared with thresholds reported in the literature.
Results
MR images of 46 rats that were previously treated with EP as described above were included in the current analysis. Treatment parameters were 600 V, 50 |js pulses at 1 Hz with varying number of treatment pulses from 10 to 540 pulses. The extent of tissue damage and BBB disruption, i.e. rd- the irreversible damage radius and rb - the BBB disruption radius were calculated from the MR images acquired 30 minutes post treatment and 2 weeks post treatment as described in the Methods section.
Radiol Oncol 2016; 50(1): 28-38.
Sharabi S et al. / A statistical model describing brain electroporation
33
Dependence on the number of treatment pulses
The average radius of each treatment group as calculated from the MR images is presented in Table 2.
The dependence of r on N has been previously described by both logarithmic and power functions.44 Here, by fitting the mean rd of each treatment group to the number of electric pulses - N, we found the logarithmic function to provide a better fit to the data, resulting in the following dependence of r, on N:
rd(N) = 0.3267 • ln(0.8123 • N)
[8]
r = 0.84, p < 0.005.
Similarly, by fitting the mean rb of each treatment group to N, the dependence of rb on N was found to be:
\ «ft. w
®
^ I &f
w W<W
®
FIGURE 2. MRI example. (A) 3 slices of contrast-enhanced T1-weighted MR images of a rat treated with 45 electroporation pulses. The MR images was obtained 30 min post treatment. Each slice is 2mm thick. The enhancing region represents BBB disruption. (B) ROI (green) plotted in the MR image to mark the enhancing region.
r (A/) = 0.4213ln(1.535 A/)
[9]
r = 0.96, p < 0.0001.
rd and rb can be seen in Figure 3. The average ratio between rb(N)and rd(N) was found to be 1.67 ± 0.11 (s.e.m), confirming the coverage of significant volumes surrounding the IRE with BBB disruption. The small error suggests that the ratio between rd(N) and rb(N) is not affected by the number of applied pulses. The ratio between rd(N) and rb(N) plotted as a function of the number of treatment pulses supports this observation (Figure 3B). The coefficients of the empirical function for the BBB disruption are higher, because the BBB is disrupted by electric fields lower than those required for IRE ablation.
Irreversible damage model
The coefficients Ecd and Ad of equation [7] were calculated for each value of N as shown in Figure 4A-B. In order to find Ecd and Ad we used rd values obtained from equation [8] rather than using the average values obtained from the experiments, as this equation describes the dependence of rd on the number of treatment pulses based on the experimental data. Although Ec(N) is traditionally described with an exponential function we chose to describe it here using a power function as it fitted the data considerably better (r2 was considerably larger: 0.89 for the power function versus 0.5 for the exponential function), especially in the high N range. Still, when fitting the optimization results of Ec(N) of only the first 90 pulses to an exponential function, r2 increased to 0.83 (Figure 4C).
1 H 0.8 -0.6 -
o 0.6 -T T
2 0.4 -I T
0.2 -0
200 400 600 number of pulses
0
200 400 number of pulses
600
FIGURE 3. (A) radii of irreversible damage and BBB disruption calculated from the MRIs, as a function of the number of treatment pulses, and the logarithmic equations fits (B) ratio between rb(N) and rd(N) as a function of number of the number of treatment pulses.
1,500
Ec= 150475*N-0.185
r2 = 0.8906	lj000
500 0
An = 1119.9*e-o.oo2N r* = 0.9611
® 0 200 400 600 number of pulses
0 200 400 600 number of pulses
^ 1500 X
-ë 1000
— 500 o
Ec = 106228 e-000™ r2 = 0.8322
50
number of pulses
100
2.5
Ê 2 JL 1.5
? 1 oo
S" 0.5 0
®
R(S=0) =0.8743rb + 0.1984
Ia = 0.7883	$
05
1	15
rd (mm)
FIGURE 4. Dependence of Ecd (A) and Ad (B) on the number of treatment pulses. (C) Exponential dependence of Ecd on the number of treatment pulses with N limited to 90 pulses. (D) Correlation between radii obtained from experimental data and radii obtained from the statistical model for IRE. Error bars represent 95% confidence level.
Radiol Oncol 2016; 50(1): 28-38.
34
Sharabi S et al. / A statistical model describing brain electroporation	34
TABLE 2. Average radii of IRE and BBB disruption for each treatment group. Each group of 5-7 rats was treated with different number of pulses (10-540) at 600V, 50|js pulses at 1Hz
# of pulses 10	45	90 180	270	450 540
IRE rad'us 0.62 ± 0.15	1.35 ± 0.18	0.89 ± 0.20 1.42 ± 0.15	1.37 ± 0.16	1.92 ± 0.07 1.80 ± 0.21 (mm)
BBB disruption 1.25 ± 0.06	1.74 ± 0.04	1.84 ± 0.07 2.54 ± 0.15	2.19 ± 0.14	2.84 ± 0.04 2.69 ± 0.12 radius (mm)
6000 -, 5000 -4000 -(3000 -' 2000 -1000 -0
An = 48.072e-a004N r2 = 0.9742
200 400 Number of pulses
600
®
0
200 400 Number of pulses
600
FIGURE 5. Dependence of Ecb (A) and Ab (B) on the number of treatment pulses for BBB disruption. Error bars represent 95% confidence level.
=0) = 200326x-0.245 r2 = 0.9647
♦BBB=1 =0
FIGURE 6. Electrical field thresholds. (A) IRE thresholds. Dashed line represents published IRE thresholds for white matter for 80 50 js pulses at 4 Hz. (B) BBB disruption thresholds. Dashed line represent previously published threshold for 90 50 js pulses at 4 Hz.5 (C) Thresholds for E(S = 0) for the IRE and E(S = 1) for BBB disruption. (D) Ratio between E(S = 1) and E(S = 0) for IRE and E(BBB = 0) and E(BBB =1 ) for BBB disruption. Error bars are smaller than markers.
Linear regression analysis confirmed that r(S = 0), calculated from the Peleg-Fermi equation with Ecd(N) and Ad(N)described well the rd obtained from the experimental data: F(1,5) = 45, p < 0.008, r2= 0.79. The resulting regression equation was: rd=
0.19 + 0.87 x, (x = r(S = 0)).
BBB disruption model
The same optimization method that was used to calculate Ecd and Ad was applied to the BBB disruption data with rb set to BBB(E,N) = 1, meaning that for radii larger than rb the BBB was not breeched(less than 0.001% ). Ecb and Ab were calculated for each value of N as can be seen in Figure 5A-B.
As for the IRE models, the goodness of the fit to the experimental data was also evaluated. r(BBB = 1) and r(BBB = 0) were calculated from BBB(E,N) for each value of N using Ab(N) and Ecb(N).
Next we evaluated the correlation between rb(N) obtained from the experimental data, and r(BBB = 0) linear regression analysis confirmed that r(BBB = 1), calculated from the extended Peleg Fermi model with Ecb(N) and Ab(N) described well the behavior of rb obtained from experimental data: F(1,5) = 45, p < 0.001, r2= 0.91. The regression equation was: rb = 0.19 + 0.87 x (x = r(S = 0)).
Electric field thresholds
The electrical field thresholds for E(S = 0) and E(S = 1) were calculated from the model for cell death. E(BBB = 0) and E(BBB = 1) were calculated for BBB disruption. In the cell death model, E(S = 0) represents the threshold needed for over 99.99% cells death whereas electric field lower then E(S = 1) will cause cell death lower than 0.001%. In the BBB disruption model E(BBB = 0) represents the threshold needed for over 99.99% of the BBB to be breeched while electric field lower then E(BBB = 1) will not disrupt the BBB(BBB disruption lower than 0.001%). Thresholds are presented in Figure 6.
The ratio between S(E,N) = 0 and S(E,N) = 1 thresholds, representing the transition zone between over 99.99% cell death and no cell death (S(E,N) = 0 / S(E,N) = 1) thresholds was calculated. The ratio is relatively high (between 0.88 and 0.91) even for small numbers of pulses and depends only weakly on the number of treatment pulses. This means that the transition between 99.99% cell death threshold and no cell death threshold is narrow and gets even narrower for large numbers of
Radiol Oncol 2016; 50(1): 28-38.
Sharabi S et al. / A statistical model describing brain electroporation
35
treatment pulses. This is not the case with BBB disruption, where the ratio between the thresholds increases with the number of treatment pulses eventually converging to one (Figure 6). The average ratio between BBB disruption ratio and damage ratio is 1.67 ± 0.11 (s.e.m).
Thermal model
The initial temperature of the rats' brain in the simulation was set to 37°C to show that the treatment should not induce thermal damage in clinical use. The maximum temperature reached in the tissue after treatment at 600 V was 38.9°C (Figure 1B). This temperature was reached using 540 pulses. Although temperature was not measured during EP treatment, histological analysis of brains extracted 60 min post treatment revealed no signs indicative of thermal damage such as coagulation, or extensive hemorrhages.45 Connective tissue and blood vessels were preserved in the treated area suggesting damage induced only by IRE.46
Discussion
When treating tumors by EP, it is important to deliver the electric pulses so that the entire tumor volume will be treated to avoid recurrence. It is also vital to treat the infiltrating zone surrounding the tumor mass with high efficacy while preserving the healthy tissue. This is especially important in the case of brain tumors, where the infiltrative zone is relatively large47 and the preservation of healthy brain tissue is of critical importance.
The electrode configuration in this experiment, which consists of a single intracranial insulated electrode with an exposed tip combined with an external surface ground electrode, provides an electric field distribution that is strongest at the in-tracranial needle tip tissue interface and decreases with the square of the radius. This setup provides a well-controlled region of permanent damage induced by IRE, further surrounded by significant BBB disruption zone. This combined response offers the potential of this configuration for the treatment of brain tumors combining IRE and chemotherapy. This setup induces rapid tissue damage in the tumor mass surrounded by significant BBB disruption, thus potentially enabling efficient drug delivery of systemically administered drugs to the infiltration zone surrounding the tumor. Since GBM cells are highly resistant to therapeutic ap-optotic stimuli, however, they exhibit a paradoxi-
cal propensity for extensive cellular necrosis1920, IRE may be efficient for treating the tumor mass. Disrupting the BBB in the local vicinity of the tumor can also improve drug intake since peripheral administration of therapeutic agents is inefficient due to poor penetration of most drugs across the BBB. As clinical trials using IRE or ECT for the treatment of deep seated tumor are becoming more common8-16, a model that can predict treatment outcome and enable individual treatment planning is increasingly being recognized as a need.
The statistical model of EP-induced cell death used in this manuscript was originally suggested by Golberg et al.32 who validated the model using experimental data in vitro. Here, we present for the first time an experimentally validated statistical model for tissue IRE where the Peleg-Fermi model was extended to a wide range of number of treatment pulses r and to BBB disruption.
The results of this study demonstrate the feasibility of applying the Peleg-Fermi model for describing irreversible EP in the brain and for treatment planning. Furthermore, since our model is based on experiments with up to 540 pulses we were able to extend the model beyond the up to 90 traditional pulses used for IRE. This is important since protocols outside the traditional 100 pulses are being evaluated48 49 and a tool to evaluate protocols with higher number of pulses is needed.
The results of the model indicate, as expected, that with increasing number of treatment pulses it is possible to treat larger volumes of tissue and that the IRE threshold decreases with the number of pulses. This is however true up to a limited extent since both Ec(N) and the thresholds eventually plateau. This suggests that although increasing the number of pulses while lowering the treatment voltage may represent a safe way to avoid thermal damage while still achieving large enough treatment volumes, there is an upper limit for this effect. In addition, when using higher voltages, raising the number of pulses will eventually lead to increased damage induced by Joule heating but will not increase the damage induced by IRE. This plateau phenomenon does not only result from the logarithmic behavior of r(N),as can be seen in Equations [8], [9] (Figure 6) but also in Ecd(N) (Figure 4B) and Ecb(N) (Figure 5B).
Although the behavior of the equations describing Ec(N) was previously described as expo-nential32,36, which supports the claim that larger number of pulses increases lethality, we found that Ec(N) is better described by a power function. As power functions plateau faster than exponential
Radiol Oncol 2016; 50(1): 28-38.
36
Sharabi S et al. / A statistical model describing brain electroporation	36
functions it further supports the limited effect of increasing the number of treatment pulses. One explanation for the difference might be that previously the model was limited to 10-90 pulses3236, and therefore the plateau effect was not yet reached. In a paper published about evaluation of the Fermi equation as a model of dose-response curves on dose response the author described Ec as a Weibull function suggesting exponential function is just a simplification for limited range of pulses.50
Once we found that the Peleg-Fermi model can be used to describe IRE, we continued to further extend the model to describe BBB disruption induced by EP. For the model, we correlated radii calculated from contrast-enhanced Tl-weighted MRI with BBB(E,N) = 1 since even a relatively small BBB disruption can be visible. We found that the extended Peleg-Fermi model describes well not only the behavior of IRE radii but also that of BBB disruption induced by EP with high statistical significance. This also indicates that there are possibly similar underlying mechanisms at play, which cause the effects.
The combination of the two models can be used for efficient treatment planning for brain tumors where IRE is used for ablating the tumor mass while BBB is disturbed in the rims and infiltrating zone thus allowing efficient access of therapeutic agents. The rims in this setup are on average 1.67 + 0.11 mm wider than the damage, with no correlation to number of pulses. This suggests that the volume of BBB disruption is over 4 times larger than the volume of IRE.
Although both the IRE and the BBB disruption models were constructed separately, when planning a treatment protocol for brain, both should be used since BBB disruption with no irreversible damage only occurs in relatively low electric fields and when higher voltages or higher number of pulses are used, irreversible damage is difficult to avoid.
The electric field thresholds for IRE and reversible EP are mostly limited to the traditional treatment protocol, i.e. 90 pulses, but when using protocols that include a different number of pulses, different thresholds should be used.51 In this study we calculated the thresholds needed for the different number of pulses and fitted them to a power function. The thresholds we found for 90 pulses fit well within the thresholds previously reported in the literature for IRE in white matter52 and BBB dis-ruption5 as can be seen in Figure 6.
The ratio between the thresholds of BBB(E,N) = 0 and BBB(E,N) = 1 were found to increase with
the number of pulses suggesting that the window between 99.99% of BBB disruption and no BBB disruption narrows with the number of pulses. This suggests that while increasing the number of pulses will eventually not lead to bigger radius of BBB disruption, larger percentage of the BBB will be disrupted thus improving drug penetration to the tissue. This is not the case for IRE where the ratio between the thresholds for S(E,N) = 0 and S(E,N) = 1 seems to be nearly independent on the number of treatment pulses. This is somewhat surprising but could be explained by the fact that the ratio is relatively high to begin with (between 0.88-0.91) and that our dataset starts with 10 pulses, however the range of pulse numbers in this study covers the most commonly used IRE protocols in clinical practice. This is also consistent with previous publications saying there is a sharp delineation of IRE treated and healthy tissue53 and demonstrates that the sharp delineation is maintained even for high number of pulses.
The ratio between rb and rd was found to be nearly independent on the number of pulses. This is further supported by the relatively constant ratio between the thresholds that where calculated from the model. Thus, during treatment planning it might be sufficient to calculate one radius. It also indicates that BBB disruption may be used as a safety limit for irreversible EP. Though when using other electrode configurations, caution is needed. If thermal damage occurs, typically at high voltages or high number of pulses, it may influence the ratio between cell death and BBB disruption.
The thermal model showed only a mild increase in brain temperature. The maximal temperature at the end of 540 pulses reached 38.9°C. Since 42°C is often considered the thermal damage threshold if sustained for long durations42, it is safe to assume that the tissue damage found in our experiments was induced solely by EP and not by thermal effects. This was also confirmed by histology7 showing no signs of thermal damage, although a temperature assessment in real time is advisable.
Despite our understanding that this model may be used by physicians and researchers for the selection of treatment protocols, a model that also incorporates dependence on additional treatment parameters such as frequencies and pulse durations should be developed.28 Such all-inclusive model would enable physicians to choose the safest and most efficient protocol on a per-patient bases. Another point to bear in mind is that although the electrode configuration suggested in this paper produces very low Joule heating, using other elec-
Radiol Oncol 2016; 50(1): 28-38.
Sharabi S et al. / A statistical model describing brain electroporation
37
trode configurations with high number of pulses might induce thermal damage in addition to IRE.54 Although this study indicates that the combination of IRE and BBB disruption may be applied for the treatment of brain tumors, experimental validation using animals bearing intracranial tumors is yet to be done.
In conclusion, the results of our study indicate that it is possible to apply high voltage electric pulses in a manner that induces localized focused irreversible damage in the brain surrounded by a larger volume of BBB disruption while using a single minimally invasive intracranial electrode. We used existing statistical models of cell kill by electric pulses that were based on theoretical cases and validated them using in vivo experimental data and extended the knowledge of EP thresholds beyond the traditional 90 pulses protocol used in IRE. We further extended the model to describe BBB disruption induced by EP. These models can assist physicians and researchers in selecting optimal treatment protocols allowing them to achieve the desired outcome in treating brain tumors. Although validation of the model in tumors is yet to be done, the results confirm that treatments outside the most commonly used protocols can achieve expected outcome.
Acknowledgements
This work was performed in partial fulfillment of the requirements for a Ph.D. degree of Shirley Sharabi, Sackler Faculty of medicine, Tel Aviv University, Israel and was supported by COST TD1104 STSM number 010315-057446. This joint paper is a result of networking efforts within COST Action TD1104 (www.electroporation.net). The work was supported by Slovenian Research Agency under various grants. The work was performed in the scope of LEA EBAM. This research was supported by the Israel Science Foundation.
References
1.	Weaver JC. Electroporation of biological membranes from multicellular to nano scales. IEEE Trans Dielectr Electr Insul 2013; 10: 754-68.
2.	Kotnik T, Kramar P, Pucihar G, Miklavcic D, Tarek M. Cell membrane elec-troporation-Part 1: The phenomenon. IEEE Elect Insul Mag 2012; 28: 14-23.
3.	Yarmush ML, Golberg A, Sersa G, Kotnik T, Miklavcic D. Electroporation-based technologies for medicine: principles, applications, and challenges. Annu Rev Biomed Eng 2014; 16: 295-320.
4.	Davalos RV, Mir LM, Rubinsky B. Tissue ablation with irreversible electroporation. Ann Biomed Eng 2005; 33: 223-31.
5.	Garcia PA, Rossmeisl JH, Jr., Robertson JL, Olson JD, Johnson AJ, Ellis TL, et al. 7.0-T magnetic resonance imaging characterization of acute blood-brain-barrier disruption achieved with intracranial irreversible electroporation.
PloS One 2012; 7: e50482.
6.	Hjouj M, Last D, Guez D, Daniels D, Sharabi S, Lavee J, et al. MRI study on reversible and irreversible electroporation induced blood brain barrier disruption. PloS One 2012; 7: e42817.
7.	Sharabi S, Last D, Guez D, Daniels D, Hjouj MI, Salomon S, et al. Dynamic effects of point source electroporation on the rat brain tissue. Bioelectrochemistry 2014; 99: 30-9.
8.	Edhemovic I, Gadzijev EM, Brecelj E, Miklavcic D, Kos B, Zupanic A, et al. Electrochemotherapy: a new technological approach in treatment of metastases in the liver. Technol Cancer Res Treat 2011; 10: 475-85.
9.	Kwon D, McFarland K, Velanovich V, Martin RC, 2nd. Borderline and locally advanced pancreatic adenocarcinoma margin accentuation with intraoperative irreversible electroporation. Surgery 2014; 156: 910-20.
10.	Linnert M, Iversen HK, Gehl J. Multiple brain metastases - current management and perspectives for treatment with electrochemotherapy. Radiol Oncol 2012; 46: 271-8.
11.	Mevio N, Bertino G, Occhini A, Scelsi D, Tagliabue M, Mura F, et al. Electrochemotherapy for the treatment of recurrent head and neck cancers: preliminary results. Tumori 2012; 98: 308-13.
12.	Miklavcic D, Mali B, Kos B, Heller R, Sersa G. Electrochemotherapy: from the drawing board into medical practice. Biomed Eng Online 2014; 13: 29.
13.	Pech M, Janitzky A, Wendler JJ, Strang C, Blaschke S, Dudeck O, et al. Irreversible electroporation of renal cell carcinoma: a first-in-man phase I clinical study. Cardiovasc Intervent Radiol 2011; 34: 132-8.
14.	Philips P, Hays D, Martin RC. Irreversible electroporation ablation (IRE) of unresectable soft tissue tumors: learning curve evaluation in the first 150 patients treated. PloS One 2013; 8: e76260.
15.	Scheffer HJ, Nielsen K, de Jong MC, van Tilborg AA, Vieveen JM, Bouwman AR, et al. Irreversible electroporation for nonthermal tumor ablation in the clinical setting: A systematic review of safety and efficacy. J Vasc Interv Radiol 2014; 25: 997-1011.
16.	Um SJ, Choi YJ, Shin HJ, Son CH, Park YS, Roh MS, et al. Phase I study of autologous dendritic cell tumor vaccine in patients with non-small cell lung cancer. Lung Cancer 2010; 70: 188-94.
17.	Jiang C, Davalos RV, Bischof JC. A review of basic to clinical studies of irreversible electroporation therapy. IEEE Trans Biomed Eng 2015; 62: 4-20.
18.	Rossmeisl JH, Jr., Garcia PA, Pancotto TE, Robertson JL, Henao-Guerrero N, Neal RE 2nd, et al. Safety and feasibility of the NanoKnife system for irreversible electroporation ablative treatment of canine spontaneous intracranial gliomas. J Neurosurg 2015; 123: 1008-25.
19.	Brat DJ, Van Meir EG. Vaso-occlusive and prothrombotic mechanisms associated with tumor hypoxia, necrosis, and accelerated growth in glioblastoma. Lab Invest 2004; 84: 397-405.
20.	Raza SM, Lang FF, Aggarwal BB, Fuller GN, Wildrick DM, Sawaya R. Necrosis and glioblastoma: a friend or a foe? A review and a hypothesis. Neurosurgery 2002; 51: 2-12; discussion 12-3.
21.	Agerholm-Larsen B, Iversen HK, Ibsen P, Moller JM, Mahmood F, Jensen KS, et al. Preclinical validation of electrochemotherapy as an effective treatment for brain tumors. Cancer Res 2011; 71: 3753-62.
22.	Miklavcic D, Snoj M, Zupanic A, Kos B, Cemazar M, Kropivnik M, et al. Towards treatment planning and treatment of deep-seated solid tumors by electrochemotherapy. Biomed Eng Online 2010; 9: 10.
23.	Neal RE, 2nd, Garcia PA, Kavnoudias H, Rosenfeldt F, Mclean CA, Earl V, et al. In vivo irreversible electroporation kidney ablation: experimentally correlated numerical models. IEEE Trans Biomed Eng 2015; 62: 561-9.
24.	Pavliha D, Kos B, Marcan M, Zupanic A, Sersa G, Miklavcic D. Planning of electroporation-based treatments using Web-based treatment-planning software. J Membr Biol 2013; 246: 833-42.
25.	Pavliha D, Music MM, Sersa G, Miklavcic D. Electroporation-based treatment planning for deep-seated tumors based on automatic liver segmentation of MRI images. PloS One 2013; 8: e69068.
26.	Zupanic A, Kos B, Miklavcic D. Treatment planning of electroporation-based medical interventions: electrochemotherapy, gene electrotransfer and irreversible electroporation. Phys Med Bio 2012; 57: 5425-40.
Radiol Oncol 2016; 50(1): 28-38.
38
Sharabi S et al. / A statistical model describing brain electroporation
27.	Groselj A, Kos B, Cemazar M, Urbancic J, Kragelj G, Bosnjak M, et al. Coupling treatment planning with navigation system: a new technological approach in treatment of head and neck tumors by electrochemotherapy. Biomed
Eng Online 2015; 14 Suppl 3: S2.
28.	Miklavcic D, Davalos RV. Electrochemotherapy (ECT) and irreversible electroporation (IRE) -advanced techniques for treating deep-seated tumors based on electroporation. Biomed Eng Online 2015; 14 Suppl 3: 11.
29.	Kranjc M, Markelc B, Bajd F, Čemažar M, Serša I, Blagus T, et al. In situ monitoring of electric field distribution in mouse tumor during electroporation. Radiology 2015; 274: 115-23.
30.	Qin Z, Jiang J, Long G, Lindgren B, Bischof JC. Irreversible electroporation: an in vivo study with dorsal skin fold chamber. Ann Biomed Eng 2013; 41: 619-29.
31.	Miklavcic D, Semrov D, Mekid H, Mir LM. A validated model of in vivo electric field distribution in tissues for electrochemotherapy and for DNA electrotransfer for gene therapy. Biochim Biophys Acta 2000; 1523: 73-83.
32.	Golberg A, Rubinsky B. A statistical model for multidimensional irreversible electroporation cell death in tissue. Biomed Eng Online 2010; 9: 13.
33.	Dermol J, Miklavčič D. Predicting electroporation of cells in an inho-mogeneous electric field based on mathematical modeling and experimental CHO-cell permeabilization to propidium iodide determination. Bioelectrochemistry 2014; 100: 52-61.
34.	Garcia PA, Davalos RV, Miklavcic D. A numerical investigation of the electric and thermal cell kill distributions in electroporation-based therapies in tissue. PloS One 2014; 9: e103083.
35.	Dermol J, Miklavcic D. Mathematical models describing Chinese hamster ovary cell death due to electroporation in vitro. J Membr Biol 2015; 248: 865-81.
36.	Peleg M. A model of microbial survival after exposure to pulsed electric fields. J Sci Food Agric 1995; 67: 93-9.
37.	Sel D, Lebar AM, Miklavcic D. Feasibility of employing model-based optimization of pulse amplitude and electrode distance for effective tumor electropermeabilization. IEEE Trans Biomed Eng 2007; 54: 773-81.
38.	Miklavcic D, Towhidi L. Numerical study of the electroporation pulse shape effect on molecular uptake of biological cells. Radiol Oncol 2010; 44: 34-41.
39.	Corovic S, Lackovic I, Sustaric P, Sustar T, Rodic T, Miklavcic D. Modeling of electric field distribution in tissues during electroporation. Biomed Eng Online 2013; 12: 16.
40.	Elwassif MM, Kong Q, Vazquez M, Bikson M. Bio-heat transfer model of deep brain stimulation-induced temperature changes. J Neural Eng 2006; 3: 306-15.
41.	Garcia PA, Rossmeisl JH, Jr., Neal RE, 2nd, Ellis TL, Olson JD, Henao-Guerrero N, et al. Intracranial nonthermal irreversible electroporation: in vivo analysis. J Membr Biol 2010; 236: 127-36.
42.	Garcia PA, Rossmeisl JH, Jr., Neal RE, 2nd, Ellis TL, Davalos RV. A parametric study delineating irreversible electroporation from thermal damage based on a minimally invasive intracranial procedure. Biomed Eng Online 2011; 10: 34.
43.	Lagarias JC, Reeds JA, Wright MH, Wright PE. Convergence properties of the Nelder--Mead simplex method in low dimensions. SIAM J Optim 1998; 9: 112-47.
44.	Pucihar G, Krmelj J, Rebersek M, Napotnik TB, Miklavcic D. Equivalent pulse parameters for electroporation. IEEE Trans Biomed Eng 2011; 58: 3279-88.
45.	Sherar M, Moriarty J, Kolios M, Chen JC, Peters RD, Ang LC, et al. Comparison of thermal damage calculated using magnetic resonance ther-mometry, with magnetic resonance imaging post-treatment and histology, after interstitial microwave thermal therapy of rabbit brain. Phys Med Bio 2000; 45: 3563-76.
46.	Maor E, Ivorra A, Leor J, Rubinsky B. The effect of irreversible electropora-tion on blood vessels. Technol Cancer Res Treat 2007; 6: 307-12.
47.	Nieto-Sampedro M, Valle-Argos B, Gomez-Nicola D, Fernandez-Mayoralas A, Nieto-Diaz M. Inhibitors of glioma growth that reveal the tumour to the immune system. Clin Med Onco 2011; 5: 265-314.
48.	Faroja M, Ahmed M, Appelbaum L, Ben-David E, Moussa M, Sosna J, et al. Irreversible electroporation ablation: Is all the damage nonthermal? Radiology 2013; 266: 462-70.
49.	Olweny EO, Kapur P, Tan YK, Park SK, Adibi M, Cadeddu JA. Irreversible electroporation: evaluation of nonthermal and thermal ablative capabilities in the porcine kidney. Urology 2013; 81: 679-84.
50.	Peleg M. Evaluation of the Fermi equation as a model of dose-response curves. Appl Microbiol Biotechnol 1996; 46: 303-6.
51.	Pucihar G, Krmelj J, Rebersek M, Napotnik T, Miklavcic D. Equivalent pulse parameters for electroporation. IEEE Trans Biomed Eng 2011; 58: 3279-88.
52.	Garcia PA, Neal RE, Rossmeisl JH, Davalos RV. Non-thermal irreversible electroporation for deep intracranial disorders. Conf Proc IEEE Eng Med Biol Soc 2010; 2010: 2743-6.
53.	Ellis TL, Garcia PA, Rossmeisl JH, Jr., Henao-Guerrero N, Robertson J, Davalos RV. Nonthermal irreversible electroporation for intracranial surgical applications. Laboratory investigation. J Neurol 2011; 114: 681-8.
54.	Kos B, Voigt P, Miklavcic D, Moche M. Careful treatment planning enables safe ablation of liver tumors adjacent to major blood vessels by percutaneous irreversible electroporation (IRE). Radiol Oncol 2015; 49: 234-41.
Radiol Oncol 2G16; 50(1): 28-38.