Optimization of Photothermal Therapy Treatment Effect under Various Laser Irradiation Conditions

The photothermal effect refers to a phenomenon in which light energy is converted into heat energy, and in the medical field, therapeutics based on this phenomenon are used for anticancer treatment. A new treatment technique called photothermal therapy kills tumor tissue through a temperature increase and has the advantages of no bleeding and fast recovery. In this study, the results of photothermal therapy for squamous cell carcinoma in the skin layer were analyzed numerically for different laser profiles, intensities, and radii and various concentrations of gold nanoparticles (AuNPs). According to the heat-transfer theory, the temperature distribution in the tissue was calculated for the conditions under which photothermal therapy was performed, and the therapeutic effect was quantitatively confirmed through three apoptotic variables. In addition, the laser intensity and the volume fraction of AuNPs were optimized, and the results provide useful criteria for optimizing the treatment effects in photothermal therapy.


Introduction
The photothermal effect refers to the phenomenon in which light energy is converted into thermal energy [1,2]. This effect is utilized in various applications, e.g., for measuring thermal properties [3][4][5] and for laser hair removal, scar removal, and anticancer treatment in the medical field [6][7][8]. In photothermal therapy, tumor tissue is killed by increasing its temperature using a laser (light energy) in the presence of a photothermal conversion agent. One advantage of this method is that there is no bleeding, and the patient recovers faster compared with conventional treatment techniques based on incision [9,10].
Because photothermal therapy causes tumor death by increasing the temperature, it is necessary to understand the death mechanism of the biological tissue with respect to the temperature. Depending on the temperature, the mechanism of cell tissue death is either necrosis or apoptosis [11,12]. Necrosis occurs at ≥50 • C, and cell death occurs through the leakage of intracellular contents to the outside. In the case of a tumor, if it dies because of necrosis, the intracellular contents leak into the surrounding tissues, leading to a risk of cancer metastasis. Conversely, apoptosis refers to a form of self-destruction without affecting surrounding tissues and occurs between 43 and 50 • C. If the tumor tissue undergoes necrosis, as mentioned previously, there is a risk of cancer metastasis. Therefore, in photothermal therapy, it is important to kill the tissue through apoptosis by maintaining an appropriate temperature range.
In photothermal therapy, biological tissues absorb different levels of laser energy according to the wavelength bands [13][14][15]. For visible-light lasers, biological tissue absorbs large amounts of energy owing to its large light-absorption coefficient; therefore, it exhibits a large temperature increase. In addition, the temperature of the surrounding tissue excessively increases through conduction heat transfer, which causes the tissue temperature to increase beyond the temperature band corresponding to apoptosis to the temperature band where necrosis occurs. In contrast, while using an infrared laser, the light-absorption coefficient of the biological tissue is small; accordingly, the increase in the temperature of the tissue is small, owing to the minimal absorption of light energy. However, because this leads to a situation where the temperature range corresponding to apoptosis cannot be reached, photothermal agents are administered to the tumor tissue for increasing the temperature to the desired range. The photothermal agents are materials that increase the light-absorption coefficient at a specific wavelength and have the advantage of being able to increase the absorption rate for a desired laser wavelength. They are composed of various materials, and among them, gold nanoparticles (AuNPs) are widely used for reasons such as their harmlessness to the human body [16,17].
Considering these factors, many researchers investigate photothermal therapy. Fan et al. [18] confirmed the effect of photothermal therapy on subcutaneous cervical cancer in mice using indocyanine green (IR820) as a photothermal and fluorescence agent. To improve the photothermal effect and biocompatibility of IR820, the amphiphilic polymer Pluronic F-127 was used to form IR820@F-127 nanoparticles through self-assembly. The IR820@F-127 nanoparticles exhibited a high photothermal conversion efficiency of 35.2% and excellent photothermal stability as well as excellent biocompatibility, cell absorption ability, and bioimaging ability for tumor sites. The study was conducted in vitro and in vivo, and the therapeutic effect on cervical cancer was confirmed using a 793 nm laser. Chen et al. [19] developed plasmonic nanoparticle clusters (PNCs) to fabricate highly efficient photothermal agents. Au-based PNCs using FeOOH were used, and the resonance frequency and absorption efficiency were calculated by combining the effective medium approximation theory and full-wave electrodynamic simulations. In in vitro and in vivo experiments targeting human breast cancer cells (MCF-7), mouse breast cancer cells (4T1), and human cervical cancer cells (HeLa), a photothermal conversion efficiency of up to 84% was achieved under the optimized cluster growth condition. Li et al. [20] developed new nanoparticles to increase the efficiency of photothermal treatment and address the issue of cancer recurrence owing to the lack of photothermal nanomaterials with high photothermal conversion efficiencies. Nanoparticles (Ag@Ag 2 O/LDHs-U) in which Ag@Ag 2 O core-shell and ultrathin layered double hydroxide (LDH) structures were combined to improve the photothermal performance significantly in the NIR-II region were developed. The results of high-resolution transmission electron microscopy, extended X-ray absorption fine structure spectroscopy, and X-ray photoelectron spectroscopy confirmed that the ultrafine Ag@Ag 2 O core-shell nanoparticles were highly dispersed and firmly fixed on the ultrathin LDH nanosheet. The particles exhibited a photothermal conversion efficiency of 76.9% at a wavelength of 1064 nm, and in vitro and in vivo experiments confirmed that Ag@Ag 2 O/LDHs-U had excellent biocompatibility and phototherapy efficacy in the NIR-II region. Ren et al. [21] conducted a numerical analysis about irradiating a laser onto biological tissue in photothermal therapy using gold nanoparticles. Through the Monte Carlo method and Beer's law, the heat generation of biological tissue to which gold nanoparticles were injected was calculated, and the temperature distribution inside the tissue was calculated. Wang et al. [22] created a numerical analysis model for the tumor injected with GNPs and confirmed the temperature distribution of the tumor located on the skin surface and surrounding normal tissues through the Monte Carlo method through numerical analysis. Maksimova et al. [23] confirmed the effect of photothermal therapy by confirming the temperature distribution and the absorbed laser distribution when 808 nm and 810 nm lasers were irradiated to animals treated with a silica/gold nanoshell and malignant tumors through experiments and numerical analysis. This study confirmed the destruction of tumor cells according to the intensity and type of various lasers, and through this, it was confirmed that photothermal therapy can effectively destroy tumor cells.
Recent research on photothermal therapy has mainly been directed toward developing new photothermal agents to overcome the low photothermal conversion efficiency of existing photothermal agents. Although photothermal therapy utilizes the photothermal effect, which is based on a heat-transfer theory, related theoretical studies are insufficient.
In particular, lasers do not have a uniform light intensity in the radial direction owing to their optical characteristics and have various light intensity distributions. However, studies for simultaneously confirming the photothermal therapy effect according to the various profiles, intensities, and radii of the laser and concentration of photothermal agents are insufficient. Additionally, in photothermal therapy, the treatment effect was simply confirmed phenomenologically, or the degree of apoptosis in tumors, and the amount of thermal damage to surrounding normal tissues was not quantified. Therefore, in the present study, the amounts of laser heat absorbed and scattered in the medium when various lasers were applied to the skin layer containing squamous cell carcinoma were investigated through the Monte Carlo method, and the temperature distribution in the medium was calculated through a numerical analysis using the thermal-diffusion equation. In addition, using the apoptotic variable proposed by Kim et al. [24], the extent to which the temperature range corresponding to apoptosis in the tumor was maintained and the amount of thermal damage to the surrounding normal tissue were quantified to determine the conditions that produce the optimal treatment effect under various laser irradiation conditions. Figure 1 shows the temperature of the tumor center with respect to the treatment time when the f v of injected AuNPs was 10 −6 , and the laser intensity was 100 mW. Figure 1a,b presents the results for the Gaussian and top-hat profiles, respectively. When the laser radius ratio (ϕ r ) was 1, the radii of the tumor and laser were equal, and when ϕ r was 0.5, the radius of the laser was 50% of that of the tumor. As shown in the graph, the temperature at the tumor center varied with respect to ϕ r . In general, the temperature of the tumor center was higher when ϕ r was lower. This is because the tumor tissue absorbs more energy because of the increase in the energy per unit area as ϕ r decreases. In addition, the temperature was lower for the top-hat profile than for the Gaussian profile because the energy density of the center was higher in the former case. Therefore, the temperature distribution of the tumor and surrounding normal tissues was confirmed for different laser profiles, laser radii and intensities, and f v values of injected AuNPs. Furthermore, the maintenance degree of the apoptosis temperature range during the treatment period as well as the amount of thermal damage to surrounding normal tissues was quantitatively analyzed to suggest treatment conditions.

Apoptosis Retention Ratio
In this study, the apoptosis retention ratio (θ * A ) was used to determine the degree of apoptosis in tumor tissue quantitatively [24]. θ * A is obtained by calculating the ratio of the total tumor volume to the volume corresponding to the apoptosis temperature range of 43-50 • C for each total treatment time and then calculating the average of the values. Through this variable, the degree of maintenance at 43-50 • C (a temperature range known to induce apoptosis in the given treatment period) was quantitatively confirmed. In this study, the time used for the calculation of θ * A was selected as 600 s, which is the total treatment time. Figure 2 shows the results for θ * A with respect to the ϕ r for each f v of AuNPs when a laser with a Gaussian profile was used. The intensity of the laser (P l ) at which θ * A was maximized increased with ϕ r for all f v values. This is because the apoptosis temperature band can only be induced by applying a high P l because the energy absorbed per unit area decreases as ϕ r increases. In addition, the maximum value of θ * A decreased as ϕ r decreased. This is because as ϕ r decreased, the temperature increased excessively and exceeded the apoptosis temperature range owing to intensive heating in a narrow area of the tumor tissue. Additionally, the maximum value of θ * A increased as f v decreased. This is because as f v decreased, the depth to which the laser penetrated the medium increased, and heating occurred in a wide area in the depth direction, and excessive energy absorption did not occur because the absorption coefficient was relatively small. For the laser with a Gaussian profile, the apoptosis temperature band was maintained for the longest time under the conditions of ϕ r = 1.75, f v = 10 −6 , and P l = 210 mW.
Apoptosis retention ratio (θ * A ) for different radius ratios (ϕ r ) (Gaussian profile). Figure 3 shows the results for θ * A with respect to ϕ r for each f v of AuNPs when a laser with a top-hat profile was used. The overall trend is similar to that for the Gaussian profile, and the maximum value of θ * A was larger. This is because the tumor tissue was uniformly heated along the radial direction, as energy was uniformly absorbed along each laser radius. For the laser with a top-hat profile, the apoptosis temperature band was maintained for the longest time under the conditions of ϕ r = 1.25, f v = 10 −6 , and P l = 170 mW. In contrast to the case of the Gaussian profile, the optimal point for θ * A depended on the ϕ r . This is because the amount of energy absorbed per unit area decreased as the ϕ r increased, and after the optimal radius ratio was reached, the laser energy was insufficient for achieving the apoptosis temperature range.

Thermal Hazard Retention Value
When a NIR laser is applied to tumor tissue for treatment, the laser energy is absorbed by the tumor tissue owing to the increased light-absorption coefficient resulting from the AuNP injection. Subsequently, the absorbed heat is transferred to the surrounding normal tissue through conduction, increasing the temperature of the surrounding normal tissue.
To perform selective treatment of the tumor tissue, the thermal damage to the surrounding normal tissue must be minimized. Therefore, in this study, the amount of thermal damage to normal tissue near the tumor was quantitatively analyzed using the thermal hazard retention value (θ * H ) [24]. θ * H uses weights given to phenomena expressed in each temperature band in biological tissue as shown in Table 1. In addition, as in Equation (1), the ratio of the weighted sum of the volumes corresponding to each temperature band to the total volume of the surrounding normal tissue at each time step is calculated. After that, the calculated values as in Equation (2) are taken as the average of the total treatment time.
θ * H quantitatively identifies the average amount of thermal damage during the treatment as the ratio of the volume of surrounding normal tissues to the volume of thermal damage after weighting various phenomena occurring in biological tissues according to the temperature. Figure 4 shows θ * H with respect to ϕ r for each f v of AuNPs and each laser profile. In the case of the Gaussian profile, the thermal damage was maximized when ϕ r was 0.75. For ϕ r ≤ 0.75, the amount of thermal damage to normal tissue around the tumor tissue was relatively small because the laser energy was focused on the center of the tumor tissue. Therefore, θ * H was smaller than that for ϕ r = 0.75. For ϕ r ≥ 0.75, the heat transfer amount to the surrounding normal tissues was reduced because the applied laser energy per unit area decreased. In the case of the top-hat profile, similar trends were observed, and the thermal damage was maximized when ϕ r was 1. Additionally, θ * H decreased as f v decreased. This is because the medium absorbed a smaller amount of heat owing to the reduced light-absorption coefficient as f v decreased.
For both laser profiles, θ * H was minimized when ϕ r was 1.75, that is, when the radius of the laser significantly exceeded that of the tumor tissue. However, from a treatment perspective, θ * A should also be considered simultaneously. For Gaussian lasers, the thermal damage to the surrounding normal tissues is minimized when ϕ r is 1.75, and θ * A is maximized when ϕ r is 1.75, even when considering the degree of apoptosis in tumor tissues. However, for the top-hat laser, the thermal damage to the surrounding normal tissues is minimized when ϕ r is 1.75, but θ * A is maximized when ϕ r is 1.25. However, θ * H increases with the laser intensity, and θ * A has the intensity of the laser showing the maximum value for each ϕ r . Accordingly, because the intensities for both laser profiles must be considered simultaneously, the optimal treatment conditions must be identified by simultaneously considering θ * A , θ * H , and the laser intensity.

Effective Apoptosis Retention Ratio
Previously, θ * A and θ * H were used to confirm quantitatively the maintenance of the apoptosis temperature in tumor tissues and the amount of thermal damage to surrounding normal tissues. However, because apoptosis in the tumor tissue and thermal damage to the surrounding normal tissues occur simultaneously during actual treatment, the results were confirmed through the effective apoptosis retention ratio (θ * e f f ) [24]. θ * e f f was calculated as the ratio of θ * A to θ * H , and using this parameter, the optimal treatment conditions were confirmed for different laser conditions and f v values of injected AuNPs. Figure 5 shows θ * e f f with different ϕ r values for each f v of AuNPs and each laser profile. Simultaneously confirming the degree of apoptosis in tumor tissues and the thermal damage amount of surrounding normal tissues revealed that θ * e f f was maximized when ϕ r was 1 for both laser profiles. This indicates that θ * A was maximized when ϕ r was 1.75 using the Gaussian profile and 1.25 using the top-hat profile, respectively, under the two laser profile conditions, but high thermal damage to the surrounding normal tissues occurred owing to the high P l at the optimal point. For the Gaussian profile, the optimal treatment conditions were as follows: ϕ r = 1, f v = 10 −6 , and P l = 130 mW. For the top-hat profile, they were ϕ r = 1, f v = 10 −5 , and P l = 120 mW. From a heat-transfer viewpoint, setting ϕ r to 1 for both laser profiles was identified as the optimal treatment condition. However, if minimizing thermal damage is prioritized over the therapeutic effect from a clinical perspective, a condition in which the amount of thermal damage is smaller may be preferential even if less therapeutic effects are obtained. In cases where the amount of thermal damage is minimized, treatment may be impossible. Therefore, clinically, information on the reference point for the amount of thermal damage and required treatment time is necessary, and the results of this study are useful for establishing such reference points.

Monte Carlo Simulation and Heat-Transfer Model
In Monte Carlo simulation, laser particles absorbed and scattered in a medium, e.g., biological tissue, can be considered, and the behavior of the laser particles can be estimated stochastically using random numbers [27]. The Monte Carlo simulation was performed by calculating the moving angle and distance of a particle. First, the angle of movement of the particles can be determined using the deflection angle in Equation (3) and the azimuth angle in Equation (4), where ξ, θ, g, and ψ represent selected random numbers, deflection angles, anisotropy coefficients, and azimuth angles, respectively.
After the deflection angle and azimuth angle are calculated, they are converted into a vector representing the direction in which the particle moves in the Cartesian coordinate system through Equations (5)- (7), where µ x , µ y , and µ z represent the direction cosines of the axes.
If the calculated angle is perpendicular to the tissue surface, it can be converted into Equations (8)-(10) along each axis.
The moving distance of the particle is calculated using Equation (11) using the total attenuation coefficient of the medium and the selected random number. The total attenuation coefficient of the medium is the sum of the absorption coefficient and scattering coefficient of the medium, as indicated by Equation (12), where S, µ tot , µ a , and µ s represent the distance traveled by photons once, total attenuation coefficient, light-absorption coefficient, and light-scattering coefficient, respectively.
Finally, when the angle and distance of movement are calculated, the ratio of energy attenuation when the particle moves once is calculated using Equation (13). Here, W is the energy weight of the particle, and according to the attenuation ratio based on the absorption of the laser particle, the particle moves until the energy of the particle converges to zero. The light-absorption distribution of the laser in the medium can be calculated by repeatedly calculating the number of specified particles in this process.
The temperature distribution inside the medium was calculated using the thermaldiffusion equation, i.e., Equation (14), utilizing the light-absorption distribution calculated via the Monte Carlo simulation and the thermal properties of the medium. Here, k, ρ, and c v represent the thermal conductivity, density, and specific heat, respectively, and τ, F, and P l represent time, the fluence rate, and the laser intensity, respectively [28].

Optical Properties of AuNPs and Medium
Biological tissues exhibit a small light-absorption coefficient for lasers in the NIR region. This results in the poor absorption of heat from the laser which becomes an obstacle to reaching the temperature range required for tumor death. To compensate for this, in photothermal therapy, a photothermal agent that increases the light-absorption coefficient at a specific wavelength is injected into tumor tissue to enhance the laser heat absorption. In this study, AuNPs, which increase the light-absorption coefficient at a specific wavelength via localized surface plasmon resonance, were used as photothermal agents [29].
The optical properties of AuNPs can be calculated using Equation (15) [30], where f v , Q, and r eff represent the volume fraction of injected AuNPs, optical efficiency, and lightabsorption area, respectively. The optical efficiency of the AuNPs was calculated using the discrete dipole approximation [31,32]. In this study, rod-type AuNPs with an effective light-absorption area of 20 nm and an aspect ratio of 6.67 were used. After the optical properties of AuNPs are calculated, the optical properties of the biological tissue injected with AuNPs can be calculated using Equation (16) [33].

Numerical Investigation
First, verification of numerical analysis modeling was performed. The validation of the numerical analysis model was performed through two methods: comparison of results with previous studies and an experiment using a biomimetic phantom. Comparison with previous researchers was compared with the results of Crochet et al. [34]. This study assumed that a spherical tumor with a radius of 5 mm existed in the center of a tissue with a radius of 30 mm and a length of 30 mm as shown in Figure 6a. The radius of the irradiated laser was set to 5 mm and 15 mm, and the intensity was fixed to 1.7 W. The laser irradiation time was fixed at 800 s, and the results of comparison with previous studies are shown in Figure 6b. As shown in the figure, as a result of comparison with previous studies, it was confirmed that the results were consistent within a maximum error of 2%. The verification of numerical analysis modeling through experiments was confirmed in the author's previous study through a biomimetic phantom [24]. The temperature results were compared between the experiment and numerical analysis, and the average root-meansquare error was 0.1677; therefore, the simulation accurately reflected the actual condition. In this study, numerical analysis modeling was conducted for a condition under which squamous cell carcinoma occurred inside the skin layer, which was composed of the epidermis, papillary dermis, reticular dermis, and subcutaneous fat. It was assumed that squamous cell carcinoma with a radius of 5 mm and depth of 2 mm was located at 0.1 mm below the skin surface inside normal tissue with a radius of 15 mm and depth of 20 mm, as shown in Figure 7. A 1064 nm wavelength laser was used as the heat source, and the thermal and optical properties of each skin layer and tumor are presented in Table 2.  In the numerical analysis, the profile and radius of the laser and the f v of injected AuNPs were changed to identify the optimal treatment conditions for photothermal therapy under various conditions. The lasers had Gaussian and top-hat profiles, as shown in Figure 8. The typical laser profile has a Gaussian shape, and in the case of the top-hat profile, an additional device is installed in front of the laser beam to distribute the energy density uniformly according to the beam radius. As described previously, the top-hat laser applies a uniform energy according to the laser radius. However, because the amount of energy absorbed per unit area varies according to the beam radius, for a given energy level, an appropriate laser radius must be selected. In the case of the Gaussian laser, the energy absorbed by the tumor tissue varies according to the radial direction because of the nonuniformity of the energy density along the radial direction, as shown in Figure 9. If the diameter of the laser is significantly larger than that of the tumor tissue, uniform heating is induced in the entire tumor tissue because the profile of the laser is gradually flattened. However, because the surrounding normal tissue also absorbs laser heat, an appropriate laser radius must be selected to maintain the temperature in the range where apoptosis occurs. In this study, the laser beam diameter was selected based on 1/e 2 width. Accordingly, in this study, the degree of apoptosis temperature band maintenance in tumor tissue was confirmed by varying the radius ratio of the Gaussian and top-hat lasers to the tumor. In addition, the therapeutic effect of photothermal therapy was quantitatively confirmed by changing the f v of injected AuNPs and the intensity of the laser. The treatment time was fixed at 600 s, and the numerical-analysis parameters are presented in Table 3. The laser profile was set to the Gaussian and top-hat distributions, and the ratio of the laser radius to the tumor radius (ϕ r ) was varied from 0.5 to 1.75 at intervals of 0.25. The laser intensity was varied in the range of 0-500 mW at 10-mW intervals, and the f v of injected AuNPs was varied from 10 −3 to 10 −6 at intervals of 10 −1 . When AuNPs are injected into tumor tissue, the optical properties of the entire tissue change according to the f v . In this study, rod-type AuNPs with an effective radius of 20 nm and an aspect ratio of 6.67 were used, and the optical properties of the tumor tissue containing AuNPs for each f v are presented in Table 4.

Conclusions
A heat-transfer-based numerical study was performed on photothermal therapy using AuNPs on a skin layer in which squamous cell carcinoma occurred. The behavior of the laser was analyzed using a Monte Carlo simulation, and the temperature distribution in the medium was calculated using the thermal-diffusion equation. Using the calculated temperature distribution, the treatment conditions for photothermal therapy were optimized through θ * A that quantitatively reflects the degree of apoptosis in tumor tissue, θ * H that quantifies the thermal damage to surrounding normal tissues, and θ * e f f that simultaneously considers the above two situations.
A numerical analysis was performed by changing the profile, radius, and intensity of the laser and the f v of injected AuNPs. When a laser with a Gaussian profile was used, the optimal treatment effect was obtained when the ϕ r was 1; the f v of injected AuNPs was 10 −6 ; and laser intensity was 130 mW. For a top-hat laser, the optimal treatment effect was obtained when the ϕ r was 1; the f v of injected AuNPs was 10 −5 ; and laser intensity was 120 mW. Therefore, it can be used as a criterion for improving the treatment effects in photothermal therapy under various laser conditions. In the future, the feasibility of clinical photothermal therapy must be confirmed by verifying the numerical analysis results through in vitro and in vivo experiments. The results of this study can be used as criteria for clinically selecting the degree of thermal damage and treatment effect.