Optimizing heat shock protein expression induced by prostate cancer laser therapy through predictive computational models

Marissa Nichole Rylander; Yusheng Feng; Yongjie Zhang; Jon Bass; Roger Jason Stafford; Andrei Volgin; John D. Hazle; Kenneth R. Diller

doi:10.1117/1.2241310

1 July 2006 Optimizing heat shock protein expression induced by prostate cancer laser therapy through predictive computational models

Marissa Nichole Rylander, Yusheng Feng, Yongjie Zhang, Jon Bass, Roger Jason Stafford, Andrei Volgin, John D. Hazle, Kenneth R. Diller

Author Affiliations +

Journal of Biomedical Optics, Vol. 11, Issue 4, 041113 (July 2006). https://doi.org/10.1117/1.2241310

Abstract

Thermal therapy efficacy can be diminished due to heat shock protein (HSP) induction in regions of a tumor where temperatures are insufficient to coagulate proteins. HSP expression enhances tumor cell viability and imparts resistance to chemotherapy and radiation treatments, which are generally employed in conjunction with hyperthermia. Therefore, an understanding of the thermally induced HSP expression within the targeted tumor must be incorporated into the treatment plan to optimize the thermal dose delivery and permit prediction of the overall tissue response. A treatment planning computational model capable of predicting the temperature, HSP27 and HSP70 expression, and damage fraction distributions associated with laser heating in healthy prostate tissue and tumors is presented. Measured thermally induced HSP27 and HSP70 expression kinetics and injury data for normal and cancerous prostate cells and prostate tumors are employed to create the first HSP expression predictive model and formulate an Arrhenius damage model. The correlation coefficients between measured and model predicted temperature, HSP27, and HSP70 were 0.98, 0.99, and 0.99, respectively, confirming the accuracy of the model. Utilization of the treatment planning model in the design of prostate cancer thermal therapies can enable optimization of the treatment outcome by controlling HSP expression and injury.

1. Introduction

Prostate cancer is the second leading cause of cancer-related deaths in the United States.¹ Thermally driven treatment options include minimally invasive energy-based therapies such as thermal ablation, local hyperthermia, hyperthermia sensitization as an adjuvant to radiotherapy, chemotherapy or brachytherapy, and thermally mediated drug or gene deliveries. ^{2, 3} To ensure maximum therapy effectiveness, it is vital to provide a reasonable prediction of the outcome at the time of treatment. Knowledge of the temperature time history during and following therapy enables prediction of thermal necrosis for regions where damage is severe. However, the results of the thermal treatment are often difficult to predict in regions where the temperatures are insufficient to coagulate proteins. Applied thermal stress tends to elicit the offsetting effects of overexpressed molecular chaperones known as heat shock proteins (HSP) and hyperthermia mediated cell necrosis.^{4, 5}

Hyperthermia treatments (low temperature $\sim 42$ to $44 ° C$ ) are typically performed in a single application, but are generally employed as an adjuvant to brachytherapy, radiotherapy, chemotherapy, or gene therapy to increase therapeutic efficacy. However, sublethal temperatures experienced during hyperthermia can induce HSP expression, enhancing cell viability in surviving tumor cells following thermal therapy and imparting increased tumor resistance to subsequent adjuvant treatments thereby reducing the effectiveness of the overall therapy.⁶ Alternatively, HSP induction in the normal tissue surrounding the tumor volume may enhance survival at higher thermal doses, permitting delivery of greater energy to the target tumor region and increase resistance of healthy cells to subsequent chemotherapy and radiation treatments.

Thermal ablative therapy (i.e., high temperature $T > 55 ° C$ based tissue coagulative treatments) of the prostate using interstitial applicators typically requires several treatments to ablate the entire gland, which can be repeated immediately or upon recurrence.⁷ Effective destruction of the tumor requires complete elimination of the HSP protective effects within the targeted region by either extensive thermal injury, causing denaturation of proteins, or careful coordination of repeated ablative therapies or adjuvant therapies (in the case of hyperthermia) with periods of minimal expression in the HSP expression cycle. Thermal induction of HSP expression in the normal tissue during the initial treatment or through preconditioning can provide protection against repeated ablative episodes, reducing morbidity associated with these therapies. Therefore, knowledge of the thermal dose necessary to activate or deactivate HSP expression in the prostate can be critical in planning and implementing an effective thermal treatment for thermal ablation, hyperthermia alone, or as an adjuvant to radiotherapy or chemotherapy.

Molecular chaperons such as HSP assist in refolding and repair of denatured proteins and aid in the synthesis of new proteins in response to damage in both normal and stressed cells.^{8, 9, 10} The HSP family of proteins is classified according to their molecular weight. In this study, models are employed for prediction of the thermally induced kinetics of the two HSPs most associated with poor prognosis and progression in prostate cancer: HSP27 and HSP70 (molecular mass of 27 and 70 kilodaltons, respectively). Increased expression of HSP27 is linked to induction of apoptosis and modulating oxygen species in a manner that enhances tumor survival against chemotherapeutic agents.^{11, 12, 13} Elevated levels of HSP70 have been observed in breast and cervical cancers¹⁴ and have a proven role in cell proliferation and drug resistance.¹⁵ HSP70 and HSP27 overexpression has also been linked to the synergistic effect of hyperthermia with radiation or chemotherapy.¹⁶

The present study was designed to build on the current understanding of HSP27 and HSP70 expression kinetics and injury resulting from thermal stress in normal prostate cells (RWPE-1), prostate cancer cells (PC3),¹⁷ and PC3 tumors,¹⁸ and to then incorporate these dynamics into a finite element model to predict tissue response to laser therapy. Implementation of this model will facilitate development of dosimetry guidelines and optimized hyperthermia protocols by controlling HSP expression and tissue injury in both the targeted prostate tumor region and the healthy surrounding tissue. Utilization of this therapy planning model will permit design of hyperthermia therapies that will achieve maximum tumor destruction and preservation of healthy tissue.

2. Computational Methods

Optimal laser protocols may be defined in terms of the most advantageous energy deposition pattern producing a temperature profile, which induces the desired HSP expression and tissue damage distributions. The therapeutic protocol is described by the appropriate laser source parameters consisting of wavelength, optical fiber orientation in the tissue, number of laser probes, laser power, and pulse duration. This work focuses on protocol manipulation by varying laser power, pulse duration, and number of laser probes while keeping the other laser parameters constant. The defined energy dissipation parameters will ultimately be applied to the laser source to produce destruction of the tumor by minimizing HSP expression in the targeted region and permit preservation of the healthy tissue by inducing elevated HSP expression in these regions.

To determine optimal laser parameters, it is essential to understand the levels of HSP expression in various states and the desired expression level in the targeted normal and tumor tissue following laser therapy. Figure 1 illustrates the HSP expression levels in a normal prostate with an interior tumor for nonstress conditions, following traditional laser therapy, and the desired HSP expression for an optimal laser therapy. Intratumoral heating is the method of laser irradiation depicted in Fig. 1. The basal level of HSP expression is represented by a normalized expression level of $1 mg ∕ ml$ . Under nonstress conditions, the HSP expression in the normal tissue is represented by the basal level of expression, whereas the tumor tissue exhibits a higher level of expression with increasing levels corresponding to more progressive grades of prostate cancer.¹⁶ Traditional laser therapies induce significant thermal damage in the tumor nearest the probe; however, insufficient temperature elevation to coagulate proteins along the tumor border often results in minimal injury and elevated HSP expression. HSP expression in the tumor region can result in tumor recurrence due to enhanced cell viability and increased resistance to subsequent chemotherapy and radiation treatments. Due to insufficient temperature elevation throughout the tumor in traditional laser therapies, the HSP expression in the normal tissue is unaltered, preventing utilization of the protective and curative effects of HSP expression, which could decrease morbidity to subsequent thermal exposures or adjuvant therapies. In designing the optimal laser therapy, the objective should be to control the thermally induced HSP expression in the tumor and healthy tissue region to permit complete destruction of the tumor and preservation of the healthy surrounding tissue. The desired treatment outcome can be characterized by extensive thermal injury and elimination of HSP expression within the targeted tumor region, thereby preventing enhanced cell viability and resistance to chemotherapy and radiation associated with elevated HSP expression. The healthy tissue should exhibit elevated levels of HSP expression along the healthy tissue bordering the tumor to diminish morbidity associated with repeated thermal therapies or adjuvant treatments, thereby permitting improved patient recovery of tissue function. This work focuses on laser therapy protocol design using computational models based on measured injury and HSP27 and HSP70 expression data^{17, 18} to explore the tissue response to a variety of laser parameters.

Fig. 1

HSP expression (mg/ml) in normal and cancerous prostate tissue for nonstress conditions and following both traditional and optimal laser therapy.

2.1.

Model Organization

The specific constitutive models applied and their corresponding inputs and outputs are depicted in Fig. 2 . The source parameter variables are defined as $λ$ , the wavelength of the laser source; $P$ , the laser power; $\vec{o}$ , the optical fiber orientation; $m$ , the number of laser probes; and $t$ , the duration of thermal exposure. The functions $Φ (\vec{r}, t)$ , ${HSP}_{27, 70} (T (\vec{r}, t))$ , and $Ω (T (\vec{r}, t))$ are the spatiotemporal dependent distributions of fluence, HSP27 and HSP70 expression, and damage, respectively. The dependence is implicit based on temperature $T (\vec{r}, t)$ , where $\vec{r} = (x, y, z)$ is the position vector. The bidirectional arrows indicate that both forward and inverse solution techniques can be employed to obtain the associated inputs or outputs. The model may be exercised as a closed loop control process to adjust the source parameters for improved results.

Fig. 2

Major model components for developing a laser ablation system for treatment of tumors with maximum injury and HSP expression elimination.

To design the optimal laser therapy based on a desired tissue response, an adaptive, finite element model was developed to predict the temperature history, damage, and HSP27 and HSP70 distributions associated with defined laser therapies. The numerical methods of choice in the simulations involved utilization of $h$ - $p$ adaptive finite element methods in the spatial dimension and fully implicit Crank—Nicholson scheme in the temporal dimension. The $h$ - $p$ adaptive capability refers to the ability to adaptively refine the mesh by subdividing element size (decreasing $h$ ) or increasing the polynomial order (enlarging $p$ ). This permits exponential convergence by optimizing both mesh size and the polynomial order of the elements such that the numerical error is reduced to a specified precision. An element-by-element residual method was employed to estimate the error. The code automatically evaluates the solution gradient to determine whether to subdivide or increase polynomial order.

Various tissue/tumor geometries were generated by ${Hypermesh}^{®}$ (Altair Engineering, Ann Arbor, MI) and then implemented on the ${ProPHLEX}^{®}$ (Altair Engineering) platform, which is an $h$ - $p$ adaptive, finite element package. Three types of system geometries, quarter sphere, magnetic resonance imaging (MRI)-derived tumor, and prostate with an interior tumor, were employed in the simulation to illustrate various capabilities of the model.

Initial mesh, boundary conditions, and material properties were applied to all three system geometries in ${Hypermesh}^{®}$ and then imported into the ${ProPHLEX}^{®}$ software. This software allows determination of the temperature distribution during laser heating and study of the sensitivity of the thermal behavior to manipulation of individual source parameters. ${ProPHLEX}^{®}$ also permits specification of equations and boundary conditions; definition of the systems including irregular geometrical domains, optimization of the finite element mesh, running of the simulation, visualization of the temperature, damage, and HSP27 and HSP70 expression profiles throughout the tissue; and optimization of the parameter analysis to identify target irradiation parameter values.

2.2.

System Geometries

2.2.1.

Quarter sphere tumor

In previous experiments, PC3 cells in the backs of female NOD, CB17 PRKDC-SCID (severe combined immunodeficiency) mice formed hemispherical-shaped tumors as shown in the MR image in Fig. 3 . Naturally occurring tumors within the prostate tend to be lobular and can also be reliably modeled as oblate spheroids.

Fig. 3

MR image of mouse with prostate tumor.¹⁸

Due to the inherent symmetry of a hemisphere, a quarter sphere can be considered a valid initial approximation having a simple geometry. A hexahedron 3-D finite element mesh was defined by dividing the prostate tumor into distinct incremental volumes via the meshing operations in ${Hypermesh}^{®}$ .

2.2.2.

Magnetic-resonance-imaging-derived tumor

To illustrate use of the model for providing a prediction of the tissue response based on patient-specific data, MR images of the mouse prior to laser irradiation were employed to develop an accurate model of the tumor.

To build a hexahedral mesh for the tumor and the tissue from the MR images, the boundaries were detected using segmentation techniques. A semiautomatic segmentation method was adapted to find the boundaries of the tumor and the tissue on each imaging slice manually. Cubic spline and lofting methods were applied to obtain smooth boundary surfaces from which hexahedral meshes were extracted from the segmented MRI data.

An extended dual contouring isosurface extraction method^{19, 20, 21} generated hexahedral meshes for the tumor and the tissue conforming to the segmented boundaries. Boundary surfaces were further smoothed and improved using geometric flows.²² There were three steps in the quality improvement scheme. The first step involved denoising the surface mesh through vertex adjustment in the normal direction with volume preservation. Surface diffusion flow was employed to smooth the surface, and the discretized Laplacian-Beltrami operator was solved using geometric partial differential equations (GPDE). The second step was devoted to improving the aspect ratio of the surface mesh through vertex adjustment in the tangent direction with feature preservation. Surface features were preserved, since movement in the tangent plane does not change the surface shape.^{23, 24} The third step focused on improving the aspect ratio of the volumetric (hexahedral) mesh with vertex adjustment inside the volume. Generated hexahedral meshes were imported into ${ProPHLEX}^{®}$ to simulate the temperature, HSP70, HSP27 expression, and damage distributions in the tissue using the $h$ - $p$ finite element method.

2.2.3.

Prostate with an interior tumor

To simulate a more clinically relevant scenario where healthy tissue surrounds target tumor volumes, ${Hypermesh}^{®}$ was employed to create a prostate with an interior tumor. This system geometry was imported into ${ProPHLEX}^{®}$ for simulation and determination of the solution. Using a posteriori error estimation in ${ProPHLEX}^{®}$ , all three system geometries were adaptively refined by one-irregular meshing techniques with highest polynomial order up to six $(p = 6)$ to improve the accuracy of the solution to the preset tolerance. An example of the mesh for the computational geometry for the prostate (outer region) and interior tumor (inner region) is depicted in Fig. 4 .

Fig. 4

Meshed geometry representing the prostate (outer region) with an interior tumor (interior region).

2.3.

Thermal model

The mathematical representation of the temperature distribution in the tissue incorporates both the Penne’s bioheat transfer term^{25, 26} for the thermal effects of local blood perfusion and an expression for laser energy absorption.

Eq. 1

c_{t} ρ_{t} \frac{\partial T}{\partial t} = \nabla • [k (T) \nabla T] - ω_{b} (T) c_{b} (T - T_{a}) + Q (\vec{r}, t),

where

ρ_{t}

,

c_{t}

,

T_{a}

, and

c_{b}

are the density and specific heat of the tissue, arterial blood temperature, and specific heat of the blood, respectively. The temperature-dependent thermal conductivity of the tissue and blood perfusion rate within the tumor are denoted by

k

and

ω_{b}

, respectively. The nonlinear temperature dependence of thermal conductivity and blood perfusion are incorporated into the model to give a more precise measure of the temperature distribution. The mathematical formulation employed for the nonlinear effects of the temperature-dependent blood perfusion in the tumor²⁷ is shown in Eq. 2.

Eq. 2

ω_{b} = {\begin{matrix} 0.833 & T < 37.0 \\ 0.833 - {(T - 37.0)}^{4.8} ∕ 5.438 \times 10^{3}, & 37.0 \leq T \leq 42.0 \\ 0.416 & T > 42.0 \end{matrix}} (kg ∕ s ∕ m^{3}) .

The temperature dependence of thermal conductivity for tissue has not been characterized. The thermal variation properties of water are well known in the range of 20 to $100 ° C$ and are important since the thermal properties of tissue are dependent on the water content. The temperature-dependent thermal conductivity for water is shown in Eq. 3, where $ζ$ is a temperature-dependent dimensionless correction factor and $w_{c}$ is the water content of prostate tissue, 0.511.²⁸

k (T) = 4.19 (0.133 + 1.36 ζ w_{c}) * 10^{- 1}, (W ∕ mK)

Eq. 3

where ζ = 1 + 1.78 * 10^{- 3} (T - 20 ° C) .

The nonlinear temperature dependencies of density and specific heat were neglected because they are insignificant in comparison with the nonlinearity associated with blood perfusion and thermal conductivity. Sensitivity analysis showed that the temperature dependence of thermal conductivity and blood perfusion caused a 5% decrease in predicted tissue temperature for the laser irradiation protocols.

We hypothesized that repeated thermal treatments may cause dehydration of the tissue due to an elevated vapor pressure of water at higher temperatures, leading to greater surface evaporation. With evaporation, there would be a change in the thermal conductivity, density, and specific heat of the tissue, thereby affecting the model’s prediction of temperature. We considered these effects, but found them to be insignificant. As a result, the model neglects the minimal effects of dehydration.

The rate of laser energy absorbed per unit volume within the tissue is given by the following equation:

Eq. 4

Q (\vec{r}, t) = μ_{a} Φ (\vec{r}, t),

where

μ_{a}

(m^{- 1})

, and

Φ

(W ∕ m^{2})

are the irradiation absorption coefficient and fluence distribution, respectively.^{29, 30} Diffusion theory provides an accurate approximation to the radiative transport equation in this application, since the wavelength employed for irradiation is in the infrared region of the spectrum where the scattering coefficient of tissue is much larger than the absorption coefficient. Since intratumoral heating is simulated in all cases assuming a spherical tumor geometry, the diffusion approximation yields the following expression for the fluence distribution:

Eq. 5

Φ (\vec{r}, t) = \sum_{i = 1}^{m} ϕ_{i} (\vec{r}, t) = \sum_{i = 1}^{m} \frac{3 P (t) μ_{a} μ_{t r}}{4 π ∥ \vec{r} - {\vec{r}}_{i} ∥} \exp (- μ_{eff} ∥ \vec{r} - {\vec{r}}_{i} ∥),

where

ϕ_{i}

(W ∕ m^{2})

is the fluence associated with each laser probe,

P (t)

(W) is the laser power,

\vec{r}

is the position vector from the origin,

{\vec{r}}_{i}

(i = 1 \dots m)

is the position vector of each laser probe (where

m

is the number of laser probes), and

∥ ∙ ∥

denotes the Euclidean norm. The transport attenuation and effective irradiation coefficients, represented as

μ_{t r}

(m^{- 1})

and

μ_{e f f}

(m^{- 1})

, respectively, are determined by Eqs. 6, 7. The scattering coefficient and anisotropy factor are denoted as

μ_{s}

(m^{- 1})

and

g

, respectively. ^{29, 30} In Eq. 5, we assume that the distances between each laser applicator are greater than the maximum diameter of the influence region of heat applied by each applicator, so that the assumptions in Refs. 29, 30 are valid.

Eq. 6

μ_{t r} = [μ_{a} + μ_{s} (1 - g)],

Eq. 7

μ_{eff} = {(3 μ_{a} μ_{t r})}^{1 ∕ 2} .

The optical parameters

μ_{a}

,

μ_{s}

, and

g

employed in the simulation were measured for canine prostate at a wavelength of

810 nm

. Thermal properties including

ρ_{t}

and

c_{t}

of liver were utilized, since no measured properties exist for prostate. Parameters employed in the simulation are shown in Table 1 .

Table 1

Optical and thermal properties employed in the model.31

Parameter	Symbol	Value
Diode laser wavelength	$λ$	$810 nm$
Absorption coefficient ofcanine prostate	$μ_{a}$	$4.6 m^{- 1}$
Scattering coefficient ofcanine prostate	$μ_{s}$	$1474.4 m^{- 1}$
Anisotropy factor ofcanine prostate	$g$	0.9
Density of liver tissue	$ρ_{t}$	$1045 kg ∕ m^{3}$
Density of blood	$ρ_{b}$	$1058 kg ∕ m^{3}$
Specific heat of liver tissue	$c_{t}$	$3600 J ∕ kgK$
Specific heat of blood	$c_{b}$	$3840 J ∕ kgK$
Arterial blood temperature	$T_{a}$	$310 K$

The boundary conditions applied on the quarter sphere tumor were imposed on the surface $Γ_{QS}$ , consisting of the boundaries on all the surfaces associated with the quarter sphere such that $Γ_{B} \cup Γ_{M} \cup Γ_{S} = Γ_{Q S}$ . The defined boundary surfaces and their associated boundary conditions are shown in Fig. 5 , where $T_{s}$ is the temperature of the boundary tumor surface, $T_{\infty}$ is the temperature of the air surrounding the tumor, and $h$ is the convection coefficient. Insulating boundary conditions were imposed on surfaces $Γ_{B}$ and $Γ_{M}$ , since the tumor is symmetric about these boundaries. A convective boundary condition was imposed on the curved outer surface $Γ_{S}$ to account for air flow over the tumor.

Fig. 5

The defined boundary surfaces and their associated boundary conditions for the quarter sphere tumor geometry.

The MRI-derived tumor is composed of the surface $Γ_{MRI}$ consisting of all the boundaries such that $Γ_{D} \cup Γ_{T} = Γ_{M R I}$ . Boundary surfaces and their associated boundary conditions are shown in Fig. 6 . An insulating boundary condition was imposed on surface $Γ_{D}$ , since the tumor is symmetric about this boundary. A convective boundary condition was imposed on the curved top surface $Γ_{T}$ to account for air flow over the tumor.

Fig. 6

The defined boundary surfaces and their associated boundary conditions for the MRI-derived tumor geometry.

The prostate with interior tumor depicted in Fig. 1 was assumed to be exposed to air flow on all external surfaces except for the surface between the healthy tissue and tumor interface. Therefore, a convective boundary condition was imposed on the surface $Γ_{P}$ , which consists of all boundary surfaces noncontiguous with the tumor according to:

k \frac{\partial T (x, y, z, t)}{\partial n} = h (T_{s} - T_{\infty}) on Γ_{P} .

The expressions employed for calculating the average Nusselt number $({Nu}_{Dav})$ and the convection coefficient needed for quiescent air surrounding tissue surfaces with convective boundaries are shown in Eqs. 8, 9, 10.

Eq. 8

{Nu}_{Dav} = {0.60 + \frac{0.387 {Ra}_{D}^{1 ∕ 6}}{{[1 + {(0.59 ∕ \Pr)}^{9 ∕ 16}]}^{8 ∕ 27}}}^{2},

where Pr is the Prandtl number and

{Ra}_{D}

is the Rayleigh number as defined next:

Eq. 9

{Ra}_{D} = \frac{G β (T_{s} - T_{\infty}) D^{3}}{ν α},

where

G

is the local acceleration due to gravity (m/s),

β

is the volume coefficient of air expansion (1/K) evaluated at the mean value of the surface and air temperatures,

υ

is the kinematic viscosity

(m^{2} ∕ s)

, and

α

is thermal diffusivity

(m^{2} ∕ s)

. The values for

T_{s}

and

T_{\infty}

were stipulated as 310 and

295 K

, respectively. The parameters Pr,

υ

, and

α

were evaluated for air at

295 K

with values of 0.849,

1.55 \times 10^{- 5} m^{2} ∕ s

, and

2.18 \times 10^{- 5} m^{2} ∕ s

, respectively.³² The convection coefficient for quiescent air surrounding the tumor is determined by the following equation.

Eq. 10

h = \frac{{Nu}_{Dav} k}{D},

where

k

is the thermal conductivity of the tissue and

D

is the diameter of the tumor. The radius of the simulated tumor was

6 mm

, since typical PC3 tumors grown on mice attain a maximum spherical volume of

1 {cm}^{3}

prior to exhibiting necrotic cores.

To formulate the temperature profile in Eq. 1 by using adaptive finite element methods, the functional space of admissible functions $V$ was defined as $V = {v ∊ [S^{1} (ψ)]^{3}; v = 0 o n Γ}$ , where $S^{1} (ψ)$ is the Sobolev space of functions having distributional derivatives of order 1, which is $L^{2} (ψ)$ integrable, where $ψ$ is the computational domain of interest and $Γ$ represents the surface of interest ( $Γ_{QS}$ , $Γ_{MRI}$ , or $Γ_{P}$ ). Furthermore, a $θ$ family of approximations for integrating the solution forward in time was introduced. The weak formulation of Eq. 1 can be characterized by the following bilinear form: given $T^{0}$ , find $T^{n + 1}$ $(n = 0, 1, 2, 3, \dots)$ such that

11.

B (T^{n + 1}, v) = B (T^{n}, v) + F (v) \forall v ∊ V, where

B (T^{n + 1}, v) = \int_{D} ρ c (T^{n + 1} ∕ Δ t) v d x + \int_{D} k θ (\nabla T^{n + 1} • \nabla v) d x + \int_{D} c_{b} ω_{b} θ T^{n + 1} v d x

Eq. 12

B (T^{n}, v) = \int_{D} ρ c (T^{n} ∕ Δ t) v d x - \int_{D} k (1 - θ) (\nabla T^{n} • \nabla v) d x - \int_{D} c_{b} ω_{b} (1 - θ) T^{n} v d x,

Eq. 13

F (v) = \int_{D} (Q - c_{b} ω_{b} T_{a}) v d x .

Let

P_{s}

denote a partition of

ψ

with

N_{e}

elements

e_{i}

(i = 1, \dots N_{e})

. The partition is assumed to satisfy

\bar{ψ} = \underset{{\bar{e}}_{i} ∊ P_{s}}{\cup} e_{i}

, with

e_{i} \cap e_{j} = empty (when i \neq j)

. We associate with the partition the parameter

s

, where

s = \max e_{i} ∊ P_{s} (s_{e i})

and

s_{e i} = d i a (e_{i})

. In addition, the boundary of an element

e_{i}

will be denoted by

\partial e_{i}

. Then, a finite element space

V^{s}

can be defined as:

Eq. 14

V^{s} = {v ∊ S^{1} (ψ) : v = 0 on \partial ψ, and v ∣_{e i} = \hat{v} \circ F_{e i}^{- 1}, \hat{v} ∊ I_{h p} (\hat{e}) \forall e_{i} ∊ P_{s}},

where

V^{s} \subset V

,

F_{e i}

is the affine map from the master element

\hat{e}

to the element

e_{i}

in the partition, and

I_{h p}

is the space of polynomial functions defined on

\hat{e}

with greatest degree of

p_{e i}

.

2.4.

Cell Viability and Heat Shock Protein Expression Measurement Methods

The computational models for cellular HSP expression and damage were based on previously measured thermally induced cell injury and HSP27 and HSP70 expression for normal (RWPE-1 cells) and cancerous prostate cells (PC3 cells) heated with a water bath.¹⁷ Cells were cultured in flasks and upon reaching confluence were submerged in a water bath for predetermined temperatures ranging from 44 to $60 ° C$ and heating durations for 1 to $30 \min$ . The maximum experimental temperature caused complete cell death for the shortest heating duration. Following heating, the flasks were returned to a $37 ° C$ incubator for subsequent manifestation of damage and HSP elevation. HSP expression was measured with western blotting 16-h postheating (shown to be an effective evaluation period for measuring maximum HSP 70 expression³³). Cell viability was determined with propidium iodide staining and flow cytometry 72-h postheating to permit adequate time for the extent of injury to be measured as described in previous work.^{17, 34}

The accuracy of the HSP expression and damage models were further refined by measuring the thermally induced HSP expression and tissue damage following laser irradiation of prostate tumors as described in previous work.¹⁸ The PC3 prostate cancer cells were inoculated in the backs of SCID mice and grown to a tumor burden of approximately $1.0 cc$ , at which time tumors were laser irradiated $(λ = 810 nm)$ according to a variety of heating protocols (laser fluence of 3 to $5 W ∕ {cm}^{2}$ and pulse durations of 1 to $4 \min$ ). The tumors were sectioned and stained with fluorescently tagged antibodies for both HSP27 and HSP70 following 16-h postheating to enable comparison with cellular work and adequate time for maximum HSP expression induction. The thermally induced spatiotemporal HSP27 and HSP70 expression distributions were visualized and quantified with a 3-D laser scanning confocal microscope and Leica Lite image processing software, respectively. The measured HSP expression data for the irradiated tumors complied with the proposed model, but new expression kinetics parameters were determined. Tissue damage was characterized with hematoxylin and eosin (H and E) staining as described in previous work.¹⁸ The H and E staining values permitted correlation of the HSP expression maps with tissue injury, but did not yield any quantitative information for determining new injury parameters.

2.5.

Cell Damage Model

Measured cell viability data¹⁷ enabled determination of the constitutive parameter values for an Arrhenius damage model³⁵:

Eq. 15

Ω (τ) = \ln (C_{o} ∕ C_{τ}) = A \int_{0}^{τ} \exp {- [E_{a} ∕ R T (t)]} d t,

where damage

Ω

is defined as the logarithm of the ratio of the initial concentration of healthy cells

C_{o}

to the concentration of healthy cells remaining after thermal stimulation

C_{τ}

for a heating duration of

τ

(s)

. A

(1 ∕ s)

is a scaling factor,

E_{a}

(J/mol) is the injury process activation energy,

R

(J/molK) is the universal gas constant, and

T

(K) is the instantaneous absolute temperature of the cells during stress, which is a function of time,

t (s)

. The Arrhenius damage integral was fit to the cell and tissue injury data to characterize the response as a function of thermal stress temperature and duration. At each temperature, the threshold time

(τ)

was determined for

Ω = 1

, for which

C_{τ} = 1 ∕ e

of

C_{o}

. For isothermal stress conditions, when

Ω = 1

the damage equation simplifies to the logarithmic form,

Eq. 16

\ln (τ) = [(E_{a} ∕ R) (1 ∕ T)] - \ln (A) .

The thermal damage kinetic coefficients of $A$ and $E_{a}$ were determined from the intercept and slope, respectively, for the best-fit linear function of the experimental cell viability data as described in previous work. ¹⁷ The Arrhenius damage parameters employed in the model are shown in Table 2 . There existed a breakpoint at $54 ° C$ for both cell types where unique values of $E_{a}$ and $A$ were determined. The break point may occur as a consequence of different thermal injury mechanisms for temperatures above and below $54 ° C$ . In a similar study by Bhowmick and Bischof, AT-1 cells were described as having an injury process break point at $55 ° C$ , wherein it was also hypothesized to be a consequence of a change in the mechanism of thermal injury.³⁶ The damage parameters utilized in the computational model compare well with those determined for AT-1 Dunning rat prostate tumor cells, where cell viability was measured with calcein leakage or propidium iodide (PI) staining. The damage parameters for AT-1 cells for $27 ° C < T \leq 55 ° C$ are $E_{a} = 4.78 \times 10^{4} J ∕ mol$ and $A = 1.14 \times 10^{5} s^{- 1}$ (calcein leakage), and $E_{a} = 6.58 \times 10^{6} J ∕ mol$ and $A = 1.16 \times 10^{48} s^{- 1}$ (50% survival using PI staining). For $T > 55 ° C$ , the measured AT-1 cell damage parameters were $E_{a} = 8.78 \times 10^{4} J ∕ mol$ and $A = 2.4 \times 10^{11} s^{- 1}$ (calcein leakage), and $E_{a} = 9.22 \times 10^{4} J ∕ mol$ and $A = 6.11 \times 10^{12} s^{- 1}$ (50% survival using PI staining).³⁶

Table 2

Values of activation energy (Ea) and frequency factor (A) calculated by fitting the Arrhenius damage model to measured injury data for both cell types.17

	T≤54°C		T> 54°C
Cell type	$E_{a}$ (kcal/mole)	A $(s^{- 1})$	$E_{a}$ (kcal/mole)	A $(s^{- 1})$
PC3	$2.38 \times 10^{5}$	$1.80 \times 10^{36}$	$1.24 \times 10^{5}$	$7.00 \times 10^{17}$
RWPE-1	$2.49 \times 10^{5}$	$1.03 \times 10^{38}$	$5.88 \times 10^{4}$	$5.65 \times 10^{7}$

According to the Arrhenius formulation, increasing temperature and heating duration will cause the damage value to rise indefinitely with complete cell damage represented by infinity. In order to more meaningfully represent the damage, another parameter was employed for determining injury in the finite element simulations. The damage fraction $F_{D}$ of damaged to undamaged tissue is given by

Eq. 17

F_{D} = \frac{C_{o} - C_{f}}{C_{o}} = 1 - \exp^{- Ω} .

Native tissue is represented by $F_{D} = 0$ $(Ω = 0)$ and tissue with complete cell death is denoted as $F_{D} = 1$ (when $Ω \to \infty$ ).²⁸ A damage fraction value of 0.6 typically yielded sufficient denaturation and injury to eliminate HSP expression, as evidenced in the results in Sec. 3.

2.6.

Heat Shock Protein Expression Model

The measured HSP expression kinetics data for PC3 and RWPE-1 cells and prostate tumors enabled formulation of an empirical computational model for HSP expression prediction associated with a given thermal stimulus. The proposed model describes HSP expression as a function of only temperature and heating duration, which is adequately supported by our previously measured experimental data. ^{17, 18} Employing ${Maple}^{®}$ permitted a wide array of functions to be explored for determination of the most appropriate mathematical formulation to accurately fit the entire dataset for all measured temperatures and heating durations. HSP expression induced by a transient temperature field was found to be proportional to its concentration $H = H (T, t)$ and can be described in a general form:

Eq. 18

\frac{\partial H (t, T)}{\partial t} = f (t, T) ∙ H (t, T),

where

f (t, T)

is a general rate function. We select

f (t, T) = (α - β_{1} t^{γ - 1})

, which captures the characteristics of HSP expression induction, where

α

,

β

, and

γ

are parameters that are independent of time, but are dependent on temperature, with

γ > 1

.

H (t, T)

was found to be represented by the following function¹⁷:

Eq. 19

H (t, T) = A \exp (α t - β t^{γ}),

where

A

is a temperature-dependent constant. Since the basal value of

H (t, T) = 1

at

t = 0

due to normalization,

A = 1

for the measured dataset. All HSP expression parameters were determined by the least square approach. The HSP expression model was integrated into the finite element model to enable prediction of the thermally induced HSP response.

2.7.

Model Validation Methods

To verify the accuracy of the computational model for prediction, it is necessary to implement in vivo methods of validation. The precision of the model for predicting temperature and HSP expression was confirmed through magnetic resonance temperature imaging (MRTI) and immunofluorescence staining, respectively.

2.7.1.

Magnetic resonance temperature imaging model validation of model

Prostate tumors (PC3) grown on the backs of SCID mice were irradiated externally with an 810-nm laser and the temperature distribution was measured by magnetic resonance temperature imaging (MRTI) using the proton-resonance frequency-shift method.^{37, 38, 39} The MRTI measured temperature was then compared with the model predicted temperature following laser irradiation using identical laser parameters.

All imaging studies were performed in a 1.5-T MR scanner (Signa Echospeed) equipped with high-performance gradients ( $23 mT ∕ m$ maximum amplitude and $120 T ∕ m ∕ s$ maximum slew rate) and fast receiver hardware (bandwidth $\pm 500 MHz$ ). Mice were imaged with a 3-in.-diam receive-only spiral surface coil specially designed for small animal imaging (courtesy of Giaquinto, General Electric Global Research Center, Schenectady, New York). MRTI was performed by using a complex phase-difference technique with a fast, 2-D rf-spoiled gradient-recalled echo sequence ( $TR ∕ TE = 49.5 ms ∕ 20 ms$ , flip angle=30 deg, and $bandwidth = 9.62 kHz$ ). To achieve a $5 s$ per image scan time for the small view $(4 \times 2 {cm}^{2})$ required for these experiments, a partial-Fourier imaging method was used for the phase-encode direction ( $256 \times 48$ encoding matrix) with phase-preserving reconstruction and a reduced bandwidth to minimize gradient heating limitations. The acquired voxel size was approximately $0.16 \times 0.31 \times 3.0 {mm}^{3}$ . The in-plane resolution was interpolated to a $0.16 \times 0.16 {mm}^{2}$ pixel size by using zero padding before extrapolation of the temperature data. The change in temperature from baseline after $N$ images $(Δ T_{N})$ was extrapolated from the complex-valued MRTI data $(F)$ by using the temperature dependence of the proton resonance frequency shift⁴⁰ and an assumed temperature sensitivity $(χ)$ of $- 0.01 ppm ∕ ° C$ according to the following equation.

Eq. 20

Δ T_{N} = \sum_{i = 1}^{N} δ T_{i} = \frac{\sum_{i = 1}^{N} δ φ_{i}}{2 π γ B_{o} χ T E} = \frac{\sum_{i = 1}^{N} \arg (F_{i} - i F_{i})}{2 π γ B_{o} χ T E},

where

δ T_{i}

is the temperature difference measured between the

i

and

(i - 1)

images,

δ φ_{i}

is the change in phase measured between the

i

and

(i - 1)

images,

γ B_{0}

is the resonance frequency

(63.87 MHz)

, and TE is the sequence echo time.⁴¹

2.7.2.

Validation of heat shock protein expression model with in vivo data

The accuracy of the model prediction of HSP expression was verified by comparing measured HSP expression data (acquired through immunostaining) from laser irradiated mouse prostates tumors with the model simulated HSP expression employing the same laser parameters. The laser-induced HSP expression was quantified at three tissue depths within five tumors.¹⁸

3. Computational Results

3.1.

Model Verification

The accuracy of the model’s prediction was verified by comparing the MRTI measured and model predicted temperature following laser irradiation. PC3 tumors were externally irradiated with a laser power of $3 W$ for a pulse duration of $3 \min$ and the temperature data were acquired. The finite element model was run with identical laser parameters. Figure 7 depicts the MRTI measured and model predicted temperature at the surface of the tissue. The curve generated for the MRTI measured temperature is a mean curve representing thermal profiles for five mice. Standard deviations for each measured point were not shown for clarity but were in the range of 1.0 to $2.2 ° C$ . The precision of the fit between the MRT measured and model predicted temperature was determined by the correlation coefficient of 0.98 (value of 1 denotes a perfect fit). There is close correspondence between the measured and model predicted temperature, verifying the accuracy of the computational model.

Fig. 7

Comparison between MRTI measured and model predicted temperature.

The model was further validated by comparing the model predicted HSP27 and HSP70 expression with measured expression data from immunofluorescence and confocal microscopy analysis following irradiation of PC3 tumors in vivo. The same laser parameters were employed in irradiating tumors as discussed in the temperature validation. HSP27 and HSP70 expression at three tumor depths, $0$ (surface), 1, and $2 mm$ were measured and compared to the model predicted values as shown in Fig. 8 . The curves denoting the measured HSP27 and 70 expression represent data from five mice with a standard deviation of $\pm 0.2$ and $\pm 0.27 mg ∕ ml$ for HSP27 and HSP70, respectively. The correlation coefficient between measured and model predicted HSP27 and HSP70 are both 0.99, verifying the accuracy of the model for HSP expression prediction.

Fig. 8

Comparison between measured and model predicted HSP27 and HSP70 expression.

Before discussing laser treatment simulation results for a variety of system geometries and laser parameters, it is important to clearly define the criteria used for evaluating the success of each thermal therapy. In developing the most appropriate thermal therapy, the desired treatment outcome is characterized by maximum tumor destruction and minimum healthy tissue injury. The highest priority was to design a therapy that is characterized by complete tumor destruction to prevent recurrence of tumor growth as described in the design criterion 1 next. The second most important objective is minimizing damage to healthy surrounding tissue as described in criterion 2. HSP27 or HSP70 expression in the tumor region following therapy will compromise the effectiveness of the therapy, while elevation of HSP expression in the surrounding healthy tissue will reduce morbidity associated with subsequent thermal or adjuvant therapies. In order to achieve the ideal therapy, two criteria must be satisfied.

1. Tumor must experience significant cell injury (e.g., $F_{D} \geq 0.60$ ) and HSP27 and HSP70 expression must be diminished below their basal levels of 1 $({HSP}_{27, 70} \leq 1.0 mg ∕ ml)$ . The desired $F_{D}$ value was chosen based on the minimal amount of injury necessary to maintain HSP expression below its basal level.
2. Healthy tissue must receive minimal thermal damage (e.g., $F_{D} \leq 0.01$ ) and induction of increased levels of HSP27 and HSP70 expression $({HSP}_{27, 70} > 1.0 mg ∕ ml)$ .

3.2.

Quarter Sphere Tumor

Since PC3 cells form approximately hemispherical tumors, the initial simulations were run modeling the tumor as a quarter sphere for computational efficiency. For all simulations in this work, the laser source was positioned within the tumor and the laser wavelength employed was $810 nm$ . Although numerous simulations were investigated to look at a wide range of laser parameters, Fig. 9 illustrates the ability of the finite element model to predict the temperature, damage fraction, and HSP27 and HSP70 expression distributions in the tumor for specified laser parameters of laser power of $3 W$ and pulse duration of $1 \min$ . The normal body temperature is prescribed a value of $310 K$ . For all simulations, HSP expression was normalized to a basal level of expression denoted as 1. As a result, HSP expression higher than 1 represents an elevation in expression. Native tissue is represented by an $F_{D}$ value of 0, and tissue with complete cell death is denoted by an $F_{D}$ value of 1.0 in all simulations.

Fig. 9

Quarter sphere tumor model of radius $6.3 mm$ depicting (a) temperature, (b) HSP70, (c) HSP27, and (d) damage fraction distributions following laser irradiation intratumorally.

The quarter sphere simulation shows that significant temperature elevation and damage were induced in the center of the tumor volume, with the greatest injury occurring near the laser probe and diminishing with distance from the source. There are three observable zones of interest with regard to HSP27 and HSP70 expression in the tumor. The tissue region nearest the probe where HSP27 and HSP70 expression is below the basal level represents a zone where significant damage was induced, causing denaturation or inactivation of proteins and molecules necessary for HSP synthesis. With greater distance from the source, a distinct region of increased HSP27 and HSP70 expression is encountered, in which HSP expression is elevated above the basal level of 1. This is the thermal regime where temperatures caused significant induction of HSP27 and HSP70 expression and minimal thermal injury, providing enhanced tumor cell viability and resistance to subsequent therapy. To prevent tumor recurrence in this zone of concern, HSP expression must be eliminated. The outermost tumor region is the zone where temperature elevation was below the threshold to cause significant cell damage and limited HSP expression. The two outermost regions of the tumor are characterized by insufficient thermal injury as well as increased HSP expression in an extensive volume. As a result, this therapy would not satisfy the first and most important criterion applied to the tumor volume and would be deemed a highly unsuccessful therapy.

3.3.

Magnetic-Resonance-Imaging-Derived Tumor

To create a more realistic tumor model to mimic the current state of experimentation and improve the accuracy of the finite element prediction, the mesh built from the mouse MRI data was employed as the system geometry. Figure 10 illustrates the hexahedral mesh of the tumor and the tissue employed in the simulation.

Fig. 10

Hexahedral mesh of the tumor and the tissue. The tumor mesh consists of 296 vertices and 160 hexahedra. The tissue mesh is composed of 7522 vertices and 6147 hexahedra.

The tumor component of the model developed from the MRI data was imported into ${Hypermesh}^{®}$ , where boundary conditions were imposed. The tumor was then imported into the ${ProPHLEX}^{®}$ finite element software, where the temperature, HSP27 and HSP70 expression, and damage fraction distribution were predicted. A laser power of $5 W$ and heating duration of $2 \min$ were specified. To permit visualization of the heated region of interest, a $y$ -stack series through the tumor is shown in Fig. 11 .

Fig. 11

MRI-derived tumor model depicting (a) temperature, (b) HSP70, (c) HSP27, and (d) damage fraction distributions following laser irradiation intratumorally.

Significant temperature elevation and thermal injury were induced near the probe interface, denaturing all HSP. Throughout the remainder of the tumor volume, minimal injury was incurred, leading to a lack of HSP induction in the tumor periphery. HSP expression will have no impact on the outcome of this therapy, since temperatures associated with this therapy were not conducive to their induction ( ${HSP}_{27,70} < 1$ ). This therapy will yield a suboptimal outcome due to insufficient thermal injury throughout the entire tumor volume, thereby violating criterion 1.

3.4.

Prostate with Interior Tumor

Previous simulations have only considered irradiation of the tumor volume without inclusion of surrounding healthy tissue. A prostate model represented by a normal prostate (denoted by a coarsely meshed region) with an interior tumor (denoted by a finely meshed interior region) is employed for all subsequent simulations, as shown in Fig. 12 . This system geometry demonstrates a more realistic application of the predictive model, in which healthy tissue surrounds the targeted tumor volume. The laser power employed was $3.5 W$ with a pulse duration of $1 \min$ . Significant temperature elevation and damage were induced near the laser probe, but due to the large size of the tumor, only a portion of it is affected by the laser heating. Overall, minimal thermal injury is induced throughout the tumor volume. HSP27 and HSP70 are eliminated close to the probe due to thermal injury, but elevated levels of expression are induced in a portion of the periphery of the thermally treated zone. Due to minimal temperature elevation and insignificant injury in the outermost zone of the tumor, there was no induction of HSP expression. The current therapy would result in minimal tumor destruction and a high likelihood of tumor recurrence due to elevated levels of HSP expression, thereby violating both criteria.

Fig. 12

Prostate (coarsely meshed region) with interior tumor (finely meshed region) depicting (a) temperature, (b) HSP70, (c) HSP27, and (d) damage fraction distributions following laser irradiation intratumorally, employing a single source with a laser power of $3.5 W$ and pulse duration of $1 \min$ .

In order to more effectively satisfy the therapy design criteria, the laser power was increased to $6 W$ with a pulse duration of $2 \min$ , as shown in Fig. 13 . More substantial temperature elevation occurred, creating a larger zone of thermal injury. There still exists a region in the periphery of the tumor with high levels of HSP27 and HSP70 expression and minimal thermal injury. This therapy still unsuccessfully satisfies both criteria. A single laser probe cannot effectively coagulate the entire tumor volume to permit complete tumor destruction due to the limits of thermal diffusion.

Fig. 13

Prostate (coarsely meshed region) with interior tumor (finely meshed region) depicting (a) temperature, (b) HSP70, (c) HSP27, and (d) damage fraction distributions following laser irradiation intratumorally, employing a single source with a laser power of $6 W$ and pulse duration of $2 \min$ .

3.5.

Multiple Laser Probes

In the previous simulations, a single source was insufficient to effectively destroy a tumor of this magnitude. As a result, the use of several laser probes was investigated, similar to clinical treatments employed in interstitial laser thermotherapy. Four laser probes each with a power of $1.5 W$ and pulse duration of $1 \min$ were employed to permit substantial heating of the entire tumor. Figure 14 illustrates the temperature, HSP27 and HSP70 expression, and damage fraction distribution following laser irradiation. The entire tumor experienced significant temperature elevation. Extensive protein denaturation and cell damage occurred near the probe, resulting in diminished HSP expression in these regions. However, the thermal insult induced high levels of HSP expression in a majority of the tumor, violating the first criterion. This therapy would be highly unsuccessful due to the level of protection afforded by these proteins. Criterion 2 has also been violated, since no HSP expression is induced in the healthy tissue.

Fig. 14

Prostate with interior tumor (finely meshed region) depicting (a) temperature, (b) HSP70, (c) HSP27, and (d) damage fraction distributions following laser irradiation intratumorally, employing four laser probes each with a laser power of $1.5 W$ and pulse duration of $1 \min$ .

Increasing the laser power and pulse duration for each of the four sources to $4 W$ for $3 \min$ permits extensive thermal damage throughout the entire tumor, causing complete protein denaturation and elimination of HSP27 and HSP70 expression, thereby achieving criterion 1. Criterion 2 has also been satisfied, since HSP27 and HSP70 expression is elevated along the healthy tissue bordering the tumor, imparting enhanced protection and recovery to injured normal prostate tissue. This would be considered a highly successful therapy with minimal tumor recurrence and effective mitigation of healthy tissue destruction, as shown in Fig. 15 .

Fig. 15

Prostate with interior tumor (finely meshed region) depicting (a) damage fraction and (b) HSP70 distributions following laser irradiation intratumorally, employing four laser probes each with a laser power of $4 W$ and pulse duration of $3 \min$ .

4. Discussion

State of the art computational techniques and measured cellular and tissue data were coupled to facilitate creation of an accurate therapy planning model for prediction of the temperature, HSP70 and HSP27 expression, and damage fraction distribution in prostate tumors following laser irradiation. Utilization of this computational model permits the physician to explore various sets of laser parameters to observe the tissue response prior to performing the laser therapy.

Although other researchers have developed computational models to assist with treatment planning for prostate cancer thermal therapy,^{28, 31} the proposed model incorporates several unique aspects. This model utilizes a large amount of measured cellular and tissue data in the tissue of interest, which increases the accuracy of the prediction and prevents a purely theoretical analysis. The damage predictive model was based on measured cell viability profiles using propidium iodide staining in normal prostate (RWPE-1) and prostate cancer (PC3) cells following water bath heating for thermal stimulation temperatures of 44 to $60 ° C$ . Employing the same cell lines and hyperthermia induction techniques, HSP27 and HSP70 expression kinetics were determined using western blotting techniques.¹⁷ The HSP expression predictive model was further refined to reflect differences in in vitro and in vivo responses by characterizing the HSP27 and HSP70 distribution through immunostaining and confocal microscopy in laser irradiated prostate tumors. H and E staining was also performed on these tumors to gain a better understanding of the correlation between temperature and tissue injury associated with laser heating.¹⁸ A lower thermal threshold was observed for destruction of PC3 tumors in vivo compared to their in vitro counterparts under similar conditions due to the presence of the vascular network in vivo as observed by other researchers.^{36, 42} The HSP expression and injury models based on cellular data were valid for the PC3 tumors, but the HSP expression and injury parameters employed in the model were not identical.^{17, 18}

Several types of system geometries were employed to illustrate the utility of the model for predicting the tissue response associated with laser irradiation. Mice inoculated with PC3 cells formed hemispherical tumors, therefore initial simulations were conducted using a quarter sphere geometry for computational simplicity. Our laser heating experiments in vivo thus far have employed only external irradiation. As a result, the accuracy of the computational model was verified by comparing model predicted temperature and HSP27 and HSP70 expression with in vivo temperature and HSP expression measurement using MRTI and immunofluorescence, respectively, following external laser irradiation of prostate tumors. Intratumoral heating is more clinically relevant, therefore, following model validation this form of irradiation was employed for all subsequent simulations. MR images acquired of prostate tumors grown on the back of mice were employed to create mouse-specific tumor geometries through the use of the extended dual contouring isosurface extraction method and ${Hypermesh}^{®}$ software. Utilization of this method illustrates the ability of the computational model to predict the tissue response for patient-specific data. To mimic a more realistic scenario, a prostate model containing an interior tumor was also developed. Initially, a large tumor was simulated within the prostate and treated with a single laser source. It was evident that a single source was insufficient to effectively irradiate a tumor of this magnitude. We explored the use of multiple laser probes, as is performed in interstitial laser thermotherapy. Inclusion of four laser probes in the tumor with identical source characteristics and computational iteration of source parameters enabled development of an optimal therapy with maximum tumor destruction and preservation of healthy surrounding tissue.

Thus far, our laser heating experiments in vivo have employed only external irradiation. In the future, we will perform intratumoral heating and use the associated MRTI data to provide a more accurate verification of the precision of the finite element model. External irradiation was initially employed to permit measurement of HSP expression arising from only thermal stimulation. Utilization of intratumoral laser irradiation, though more clinically relevant, may induce HSP expression from the mechanical shear stresses associated with insertion of the laser probe and related inflammation response. Once the thermally induced HSP by the laser is fully characterized, all further experimentation and refined verification will involve intratumoral heating.

A wide array of laser parameters was simulated for the various system geometries, and for the majority of cases, significant zones of HSP expression were identified. These zones represent regions of the tumor with enhanced tumor cell viability and increased resistance to subsequent chemotherapy and radiation treatments. Ideally, significant protein denaturation must be induced to reduce HSP expression below its basal level in the entire tumor to prevent compromising the therapy outcome. Induction of HSP expression in the healthy tissue bordering the tumor would enhance recovery and reduce morbidity associated with repeated thermal therapies or adjuvant treatments. In practice, healthy tissue injury may be unavoidable in completely eradicating the tumor, but HSP induction in these injured regions may help mitigate injury due to repeated thermal therapies.

The unique contribution of the proposed computational model is the integration of the first predictive model for thermally induced HSP expression based on measured cellular and tissue HSP kinetics data into a finite element model capable of temperature and damage prediction. Currently, the model is founded on a much larger repository of measured cellular HSP expression and viability data, which was induced by water bath heating. Water bath heating has a much longer time constant (approximately $1 \min$ to reach the maximum temperature of $60 ° C$ ) than the laser source, bringing into question the accuracy of measurements for heating on the order of 1 to $2 \min$ due to nonisothermal heating.

The small size of the mouse prostate prevents its use in this study. The measured HSP expression tissue data has been acquired through correlation of MRT measured temperature with HSP expression through immunostaining of laser irradiated prostate tumors grown on the back of mice.¹⁸ More extensive tissue data must be acquired in naturally occurring prostate tumors on a larger scale through the use of canine models, coupled with the investigation of intratumoral laser heating.

The largest source of error associated with the computational model is due to a lack of available parameters values for prostate tissue at the desired temperature and wavelength investigated. Thermal properties for prostate were not available for a wide range of temperatures, so properties for liver were often employed. The available measured optical properties for prostate tumor and tissue were limited to a few wavelengths. Optical properties for native canine prostate were specified for normal human prostate tumor at a wavelength of $810 nm$ . To diminish error associated with input parameters, optical and thermal measurements must be performed for the wavelength and temperatures studied. Currently, the model does not incorporate the dynamic optical properties associated with denaturation of proteins during the laser heating process, which may lead to alterations in tissue absorption and scattering properties. The model does include the nonlinear temperature dependence of perfusion and thermal conductivity, which was found to decrease the predicted temperature by 5%. The damage parameters $E_{a}$ and $A$ were determined from cell viability measurements following water bath heating of PC3 prostate cancer (PC3) and normal prostate (RWPE-1) cells. Although H and E staining has been conducted in the laser irradiated prostate tumors, quantitative values cannot be derived for determination of damage parameters. A difference between in vitro and in vivo response to hyperthermia has been documented in AT-1 Dunning rat prostate tumor models, where prostate tumors experienced significant thermal damage compared to previous studies involving AT-1 cells. This discrepancy is thought to be due to the presence of fragile and poorly organized tumor vasculature in vivo, which is easily destroyed during thermal therapy leading to hypoxia, ischemia, and accelerated tumor death.^{36, 42} As a result, additional methods for measurement of apoptosis and necrosis in laser irradiated prostate tumors will be essential to improve the accuracy of the cell injury prediction component of the computational model.

5. Conclusion

An accurate finite element model is developed that enables prediction of the temperature, damage fraction, and HSP27 and HSP70 expression distributions in laser irradiated prostate tissue. This computational model is unique in that it incorporates the first HSP expression predictive model created for thermal therapy planning. An Arrhenius injury model is also included in the computational model, and enables the damage in the tumor to be predicted with precision based on a given temperature profile in the tumor.

The integration of thermal, damage, and HSP expression models into a single finite element model enables prediction of the prostate tumor and tissue response to a given laser therapy. Simulations using an idealized quarter sphere, MRI derived tumor, and prostate with interior tumor enables investigation of the temperature, damage fraction, and HSP expression distribution for a wide range of laser parameters typically employed in surgical procedures. Definite zones of HSP expression are identified for a large number of therapies. The existence of these regions of elevated HSP expression is of great concern, and minimization of HSP expression in the tumor should be considered in patient treatment design to prevent tumor recurrence. Utilization of this predictive model will enable a physician to investigate the HSP expression response to candidate treatments to better tailor a patient-specific therapy to achieve maximum tumor destruction and minimal healthy tissue injury.

Acknowledgments

This research was funded by the Abell-Hanger Foundation, the National Science Foundation Award numbers CTS-0332052 and CNS-0540033, and the Robert and Prudie Leibrock Professorship in Engineering at the University of Texas at Austin.

References

1.

R. F. Gittes, “Carcinoma of the prostate,” N. Engl. J. Med., 324 236 –245 (1991). 0028-4793 Google Scholar

2.

S. R. Denmeade and J. T. Issacs, “Development of prostate cancer treatment: The good news,” Prostate, 58 211 –224 (2004). 0270-4137 Google Scholar

3.

M. B. Garnick, “Prostate cancer: screening diagnosis, and management,” Ann. Intern Med., 118 804 –818 (1993). 0003-4819 Google Scholar

4.

S. Madersabacher, M. Grobl, G. Kramer, S. Dirnhofer, G. Steiner, and M. Marberger, “Regulation of heat shock protein 27 expression of prostatic cells in response to heat treatment,” Prostate, 37 174 –181 (1998). 0270-4137 Google Scholar

5.

N. B. Gibbons, R. W. G. Watson, R. Coffey, H. P. Brady, and J. M. Fritzpatrick, “Heat shock proteins inhibit induction of prostate cancer cell apoptosis,” Prostate, 45 58 –65 (2000). 0270-4137 Google Scholar

6.

J. Roigas, E. S. Wallen, S. A. Loening, and P. L. Moseley, “Effects of combined treatment of chemotherapeutics and hyperthermia on survival and the regulation of heat shock proteins in Dunning R3327 prostate carcinoma cells,” Prostate, 34 195 –202 (1998). 0270-4137 Google Scholar

7.

M. D. Sherar, M. R. Gertner, C. K. K. Yue, M. E. O’Malley, A. Toi, A. S. Gladman, S. R. H. Davidson, and J. Trachtenberg, “Interstitial microwave thermal therapy for prostate cancer: method of treatment and results of a phase I/II trial,” J. Urol. (Baltimore), 166 1707 –1714 (2001). https://doi.org/10.1016/S0022-5347(05)65658-3 0022-5347 Google Scholar

8.

J. G. Kiang and G. C. Tsokos, “Heat shock protein 70kDa: molecular biology, biochemistry, and physiology,” Pharmacology, 80 183 –201 (1998). 0031-7012 Google Scholar

9.

J. Martin, A. Horwich, and F. U. Hartl, “Prevention of protein denaturation under heat stress by the chaperonin hsp60,” Science, 258 995 –958 (1992). 0036-8075 Google Scholar

10.

H. Weich, J. Buchner, R. Zimmermann, and U. Jakob, “HSP90 chaperones protein folding in vitro,” Nature (London), 358 169 –170 (1992). https://doi.org/10.1038/358169a0 0028-0836 Google Scholar

11.

P. Mehlen, C. Kretz-Remy, X. Preville, and A. P. Arigo, “Human hsp27, drosophilla hsp27, and human

α β

-crystallin expression-mediated increase in gluthathione is essential for the protective activity of these proteins against TNF-

α

-induced cell death,” EMBO J., 15 2695 –2706 (1996). 0261-4189 Google Scholar

12.

S. Oesterreich, C. N. Weng, M. Qiu, S. G. Hilsenbeck, and S. A. W. Fuqua, “The small heat shock protein hsp27 is correlated with growth and drug resistance in human breast cancer cell lines,” Cancer Res., 53 4442 –4448 (1993). 0008-5472 Google Scholar

13.

E. H. Richards, E. Hickey, L. Weber, and J. R. Masters, “Effects of over-expression of the small heat shock protein HSP27 on the heat and drug sensitivities of human testis tumor cells,” Cancer Res., 56 2446 –2451 (1996). 0008-5472 Google Scholar

14.

C. Georgopoulous and W. J. Welch, “Role of the major heat shock proteins as molecular chaperones,” Annu. Rev. Cell Biol., 9 601 –634 (1993). 0743-4634 Google Scholar

15.

E. A. Craig, J. S. Weissman, and A. L. Horwich, “Heat shock proteins and molecular chaperones: mediators of protein conformation and turnover in the cell,” Cell, 78 365 –372 (1994). 0092-8674 Google Scholar

16.

P. A. Cornford, A. R. Dodson, K. F. Parsons, A. D. Desmond, A. Woolfenden, M. Fordham, J. P. Neoptolemos, Y. Ke, and C. S. Foster, “Heat shock protein expression independently predicts clinical outcome in prostate cancer,” Cancer Res., 60 7099 –7105 (2000). 0008-5472 Google Scholar

17.

M. N. Rylander, Y. Feng, K. Diller, and S. Aggarwal, “Thermally induced HSP27, 60, and 70 expression kinetics and cell viability in normal and cancerous prostate cells,” Google Scholar

18.

M. N. Rylander, R. J. Stafford, A. Eliott, A. Shetty, J. Hazle, and K. R. Diller, “Characterization of HSP27 and 70 expression in laser irradiated prostate tumors with nanoshell inclusion,” Google Scholar

19.

Y. Zhang, C. Bajaj, and B. S. Sohn, “Adaptive and quality 3D meshing from imaging data,” 286 –291 (2003). Google Scholar

20.

Y. Zhang, C. Bajaj, and B. S. Sohn, “3D finite element meshing from imaging data,” Comput. Methods Appl. Mech. Eng., 194 5083 –5106 (2005). 0045-7825 Google Scholar

21.

Y. Zhang and C. Bajaj, “Adaptive and quality quadrilateral/hexahedral meshing from volumetric data,” 365 –376 (2004). Google Scholar

22.

Y. Zhang and C. Bajaj, “Adaptive and quality quadrilateral/hexahedral meshing from volumetric data,” Comput. Methods Appl. Mech. Eng., (0045-7825) Google Scholar

23.

Y. Zhang, C. Bajaj, and G. Xu, “Surface smoothing and quality improvement of quadrilateral/hexahedral meshes with geometric flow,” 449 –468 (2005). Google Scholar

24.

G. Sapiro, “Geometric partial differential equations and image analysis,” 71 –134 (2001) Google Scholar

25.

H. H. Pennes, “Analysis of tissue and arterial blood temperatures in the resting forearm,” J. Appl. Physiol., 1 93 –122 (1948). 0021-8987 Google Scholar

26.

E. H. Wissler, “Pennes paper revisited,” J. Appl. Physiol., 85 35 –41 (1998). 8750-7587 Google Scholar

27.

J. Lang, B. Eridmann, and M. Seebass, “Impact of nonlinear heat transfer on temperature control in regional hyperthermia,” IEEE Trans. Biomed. Eng., 46 1129 –1138 (1999). https://doi.org/10.1109/10.784145 0018-9294 Google Scholar

28.

D. Zhu, Q. Luo, G. Zhu, and W. Liu, “Kinetic thermal response and damage in laser coagulation of tissue,” Lasers Surg. Med., 31 313 –321 (2002). https://doi.org/10.1002/lsm.10108 0196-8092 Google Scholar

29.

W. Star, “Diffusion theory of light transport,” Optical-Thermal Response of Laser-Irradiated Tissue, 166 –169 Plenum Press, New York (1995). Google Scholar

30.

B. Kim, S. L. Jacques, S. Rastegar, S. Thomsen, and M. Motamedi, “Nonlinear finite-element analysis of the role of dynamic changes in blood perfusion and optical properties in laser coagulation of tissue,” IEEE J. Sel. Top. Quantum Electron., 2 922 –933 (1996). https://doi.org/10.1109/2944.577317 1077-260X Google Scholar

31.

B. Anvari, S. Rastegar, and M. Motamedi, “Modeling of intraluminal heating of biological tissue: implications for treatment of benign prostatic hyperplasia,” IEEE Trans. Biomed. Eng., 41 854 –864 (1994). https://doi.org/10.1109/10.312093 0018-9294 Google Scholar

32.

F. Incropera and D. DeWitt, Fundamentals of heat and mass transfer, 482 –505 John Wiley and Sons, New York (1996). Google Scholar

33.

S. Wang, K. R. Diller, and S. Aggarwal, “Heat shock protein 70 expression kinetics,” J. Biomech. Eng., 125 794 –797 (2003). https://doi.org/10.1115/1.1632522 0148-0731 Google Scholar

34.

M. N. Rylander, S. Wang, S. Aggarwal, and K. R. Diller, “Correlation of HSP 70 expression and cell viability following thermal stimulation of bovine aortic endothelial cells,” J. Biomech. Eng., 127 751 –757 (2005). https://doi.org/10.1115/1.1993661 0148-0731 Google Scholar

35.

F. C. Henriques, “Studies of thermal injury, the predictability and the significance of thermally induced rate processes leading to irreversible epidermal injury,” Arch. Pathol., 43 489 –502 (1947). 0363-0153 Google Scholar

36.

S. Bhowmick and J. C. Bischof, “Supraphysiological thermal injury in Dunning AT-1 prostate tumor cells,” J. Biomech. Eng., 122 51 –59 (1998). https://doi.org/10.1115/1.429627 0148-0731 Google Scholar

37.

J. De Poorter, C. De Wagter, Y. De Deene, C. Thomsen, F. Ståhlberg, and E. Achten, “The proton resonance frequency shift method compared with molecular diffusion for quantitative measurement of two dimensional time dependent temperature distribution in phantom,” J. Magn. Reson., Ser. B, 103 234 –241 (1994). https://doi.org/10.1006/jmrb.1994.1035 1064-1866 Google Scholar

38.

Y. Ishihara, A. Calderon, H. Watanabe, K. Okamoto, Y. Suzuki, K. Kuroda, and Y. Suzuki, “A precise and fast temperature mapping using water proton chemical shift,” Magn. Reson. Med., 34 814 –823 (1995). 0740-3194 Google Scholar

39.

J. Olsrud, R. Wirestam, S. Brockstedt, A. Nilsson, K. Tranberg, F. Ståhlberg, and B. Persson, “MRI thermometry in phantoms by use of the proton resonance frequency shift method: application to interstitial laser thermotherapy,” Phys. Med. Biol., 43 2597 –2613 (1998). https://doi.org/10.1088/0031-9155/43/9/012 0031-9155 Google Scholar

40.

J. C. Hindman, “Proton resonance shift of water in the gas and liquid states,” J. Chem. Phys., 44 4582 –4592 (1966). https://doi.org/10.1063/1.1726676 0021-9606 Google Scholar

41.

L. R. Hirsch, R. J. Stafford, J. A. Bankson, “Nanoshell-mediated near-infrared thermal therapy of tumors under magnetic resonance guidance,” Proc. Natl. Acad. Sci. U.S.A., 100 13549 –13554 (2003). https://doi.org/10.1073/pnas.2232479100 0027-8424 Google Scholar

42.

S. Bhowmick, J. E. Coad, D. J. Swanlund, and J. C. Bischof, “In vitro thermal therapy of AT-1 Dunning prostate tumours,” Int. J. Hyperthermia, 20 73 –92 (2004). https://doi.org/10.1080/0265673031000111932 0265-6736 Google Scholar

Citation Download Citation

Marissa Nichole Rylander, Yusheng Feng, Yongjie Zhang, Jon Bass, Roger Jason Stafford, Andrei Volgin, John D. Hazle, and Kenneth R. Diller "Optimizing heat shock protein expression induced by prostate cancer laser therapy through predictive computational models," Journal of Biomedical Optics 11(4), 041113 (1 July 2006). https://doi.org/10.1117/1.2241310

Published: 1 July 2006

Access the abstract

JOURNAL ARTICLE
16 PAGES

DOWNLOAD PAPER SAVE TO MY LIBRARY