A statistical model describing combined irreversible electroporation and electroporation-induced blood-brain barrier disruption

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 amplitude, 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.3,8–16 ECT uses EP to allow increased uptake of chemotherapeutic drugs into tumor cells12 and IRE is a method aimed at inducing tumor ablation without thermal damage.17,18 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.19,20 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.7,21

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 data30,31 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 electroporated. 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.32,33 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 pulses37,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 naïve 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 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.

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 μ s 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. r_d- the irreversible damage radius and r_b – 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.

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.

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μ s pulses at 1Hz

# of pulses	10	45	90	180	270	450	540
IRE radius (mm)	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
BBB disruption radius (mm)	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

The dependence of r on N has been previously described by both logarithmic and power functions.44 Here, by fitting the mean r_d 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_d on N: [8]rd(N)=0.3267⋅ln(0.8123⋅N)${r_d}\left( N \right) = 0.3267 \cdot \ln \left( {0.8123 \cdot N} \right)$r2 = 0.84, p < 0.005.

Similarly, by fitting the mean r_b of each treatment group to N, the dependence of r_b on N was found to be: [9]rb(N)=0.4213⋅ln(1.535⋅N)${r_b}\left( N \right) = 0.4213 \cdot \ln \left( {1.535 \cdot N} \right)$r2 = 0.96, p < 0.0001.

r_d and r_b can be seen in Figure 3. The average ratio between r_b(N)and r_d(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 r_d(N) and r_b(N) is not affected by the number of applied pulses. The ratio between r_d(N) and r_b(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.

Figure 3. (A)

Irreversible damage model

The coefficients E_cd and A_d of equation [7] were calculated for each value of N as shown in Figure 4A–B. In order to find E_cd and A_d we used r_d values obtained from equation [8] rather than using the average values obtained from the experiments, as this equation describes the dependence of r_d on the number of treatment pulses based on the experimental data. Although E_c(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 E_c(N) of only the first 90 pulses to an exponential function, r2 increased to 0.83 (Figure 4C).

Figure 4.

Figure 2.

Linear regression analysis confirmed that r(S = 0), calculated from the Peleg-Fermi equation with E_cd(N) and A_d(N)described well the r_d obtained from the experimental data: F(1,5) = 45, p < 0.008, r2= 0.79. The resulting regression equation was: r_d = 0.19 + 0.87 x, (x = r(S = 0)).

BBB disruption model

The same optimization method that was used to calculate E_cd and A_d was applied to the BBB disruption data with r_b set to BBB(E,N) = 1, meaning that for radii larger than r_b the BBB was not breeched(less than 0.001% ). E_cb and A_b were calculated for each value of N as can be seen in Figure 5A–B.

Figure 5.

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 A_b(N) and E_cb(N).

Next we evaluated the correlation between r_b(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 E_cb(N) and A_b(N) described well the behavior of r_b obtained from experimental data: F(1,5) = 45, p < 0.001, r2= 0.91. The regression equation was: r_b = 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.

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 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).

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 intracranial 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 apoptotic stimuli, however, they exhibit a paradoxical propensity for extensive cellular necrosis19,20, 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 E_c(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 E_cd(N) (Figure 4B) and E_cb(N) (Figure 5B).

Although the behavior of the equations describing Ec(N) was previously described as exponential32,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 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 pulses32,36, 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 T1-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 disruption5 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 r_b and r_d 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 electrode 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.