Numerical Study on Death of Squamous Cell Carcinoma Based on Various Shapes of Gold Nanoparticles Using Photothermal Therapy

Due to increased exposure to ultraviolet radiation caused by increased outdoor activities, the incidence of skin cancer is increasing. Incision is the most typical method for treating skin cancer, and various treatments that can minimize the risks of incision surgery are being investigated. Among them, photothermal therapy is garnering attention because it does not cause bleeding and affords rapid recovery. In photothermal therapy, tumor death is induced via temperature increase. In this study, a numerical study based on heat transfer theory was conducted to investigate the death of squamous cell carcinoma located in the skin layer based on various shapes of gold nanoparticles (AuNPs) used in photothermal therapy. The quantitative correlation between the conditions of various AuNPs and the laser intensity that yields the optimal photothermal treatment effect was derived using the effective apoptosis ratio. It was confirmed that optimal conditions exist for maximizing apoptosis within a tumor tissue and minimizing the thermal damage to surrounding normal tissues when using AuNPs under various conditions. Furthermore, it is envisioned that research result will be utilized as a standard for photothermal treatment in the future.


Introduction
Recently, ultraviolet exposure has increased due to increased outdoor activities as a result of global economic development. Accordingly, the incidence of skin cancer is increasing annually [1,2]. Skin cancer can be primarily classified into squamous cell carcinoma, basal cell carcinoma, and malignant melanoma, and the treatment of these skin cancers is performed using various methods such as chemotherapy, cryotherapy, and incision [3][4][5]. These therapy methods, however, pose various negative effects [6][7][8]. In particular, surgery via incision causes bleeding and possible secondary infection [9,10].
An alternative treatment (i.e., photothermal therapy) is garnering attention as a method that can alleviate these adverse effects [11,12]. Photothermal therapy is a treatment method that uses the photothermal effect, a phenomenon in which light energy is converted into thermal energy when light energy is irradiated onto a medium, to induce the death of tumor tissues by increasing the temperature of the target tumor tissue [13,14]. This treatment method affords quick recovery and minimal risk of subsequent infection [15,16]. Photothermal therapy primarily supplies heat to biological tissues via laser, thereby affording the easy control of the heating range and intensity [17,18]. When a visible ray laser is used (among lasers of various wavelengths), a significant amount of light absorption occurs not only in the target tumor tissue, but also in normal tissues. Hence, laser in the near-infrared region, which involves a low light-absorption coefficient in biological tissues, is applied in photothermal therapy [19]. However, because the tumor tissue has a low light absorption rate with respect to laser in the near-infrared region, a light absorption enhancer is injected into the tumor tissue to increase the light absorption rate of the tumor tissue to perform treatment. This causes a temperature increase in only the tumor tissue, enabling selective treatment [20]. Among the various light absorption enhancers, gold nanoparticles (AuNPs) are widely used because they are harmless to the human body and afford high surface workability [21][22][23]. AuNPs reach the tumor in various ways such as direct administration and intravenous injection, and are injected into the tumor through a process such as endocytosis and used for treatment [24][25][26].
Biological tissues including tumor tissues cause various types of death depending on temperature [27,28]. Apoptosis occurs between 43 • C and 50 • C, and is a form of selfdeath that does not affect the surroundings. Necrosis, in contrast, occurs at temperatures of 50 • C or above; it poses the risk of cancer cell metastasis and recurrence as it affects adjacent tissues when the tissue dies. Accordingly, the temperature range corresponding to apoptosis should be maintained to prevent necrosis. In photothermal therapy, treatment is performed by controlling the appropriate intensity of the heat source, injection amount of AuNPs, and type of AuNPs to minimize thermal damage to surrounding normal tissues while maintaining a temperature band corresponding to apoptosis [29].
Based on these factors, various studies regarding photothermal therapy are being conducted. Nam et al. [30] conducted an experimental study pertaining to photothermal therapy combined with chemotherapy. Polydopamine-coated spike-shaped AuNPs with excellent photothermal stability and optical efficiency were developed to enable photothermal therapy; consequently, successful treatment for 85% of CT26 colon carcinoma cases was confirmed. In addition, the therapeutic efficacy of TC-1 submucosa-lung metastasis was confirmed. Mackey et al. [31] conducted a study regarding photothermal therapy using rod-type AuNPs. The study was conducted both theoretically and in vitro, and information pertaining to the shapes of three sizes of gold nanorods that afforded the optimal therapeutic effect was obtained. Finally, it was demonstrated that gold nanorods with a length of 28 nm and a diameter of 8 nm were the most effective light absorption enhancers.
Broek et al. [32] conducted an experimental study regarding photothermal therapy using branched AuNPs. In this study, branched AuNPs were biofunctionalized with nanobodies of heavy chain-only antibodies and bound to HER2 antigen expressed in both breast and ovarian cancer cells. As a result of the treatment, it was confirmed that the tumor died when 690 nm continuous wave (CW) laser was irradiated with an intensity of 38 W/cm 2 for 5 min, and that the tumor was killed only under an optical density of 4 or higher when the corresponding AuNP was used. Xi et al. [33] developed new photothermal agents (PTAs) to improve the low photothermal conversion efficiency (PCE) of conventional PTAs. PCE is facilitated by the -CF3 moiety included in the meso-position of the BODIPY scaffold in the PTAs because the free rotation of -CF3 provides a pathway for efficient non-radiative decay. In addition, it was confirmed that the barrier-free rotation of -CF3 continued even when tfm-BDP was encapsulated with polymer nanoparticles. Excellent therapeutic effects were confirmed in both in vitro and in vivo experiments, and mouse experiments confirmed that tfm-BDPNP was efficiently accumulated in the tumor area and completely resected tumor tissue when near-infrared lasers of 0.3 W/cm 2 , 808 nm were used.
In summary, studies regarding photothermal therapy under various conditions have been conducted both experimentally and theoretically. However, treatment trends based on different types of AuNPs and laser irradiation conditions could not be determined, and studies suggesting optimal treatment conditions quantitatively when considering various AuNP concentrations and other conditions at the same time are insufficient. In addition, although photothermal therapy induces cell death through temperature increase due to the photothermal effect, studies pertaining to heat transfer are insufficient. Therefore, in this study, a numerical study of photothermal therapy based on heat transfer theory was conducted for actual skin layers including squamous cell carcinoma. The temperature distribution in tumor and normal tissues based on various AuNPs shapes, laser intensities, and volume fractions of injected AuNPs was obtained. In addition, the optimal treatment effect was proposed by quantitatively confirming conditions that maximize the apoptosis of tumor tissues while minimizing thermal damage to surrounding normal tissues, which is the main purpose of photothermal therapy based on the apoptotic variable proposed by Kim et al. [34].

Discrete Dipole Approximation (DDA) Method
In this study, the DDA method [35,36] was used to calculate the optical efficiency of various nanoparticles. In this method, the scattering and absorption characteristics are calculated after assuming that dipoles with polarized shapes are formed at regular intervals based on information regarding a specific shape. Compared with the Mie theory [37,38], which can be applied only to existing spheres or ellipses, the DDA method allows a wide range of shapes that can be calculated and applied. In addition, compared to the finite difference domain method, the DDA method has the advantage of being very fast in computation speed and requiring less memory for computation [39].
To perform calculations using the DDA method, the polarization vector P must be determined. P depends on the interaction between the dipole and local electric field E; it can be calculated using Equation (1), where α is the polarizability, and r is the position vector. In addition, the local electric field E can be calculated using Equation (2).
where k is the wavenumber of the radiation. The electric field at position i due to the dipole at position j is included in the second term on the right side of Equation (2). Meanwhile, A is the interaction matrix between the dipoles (Equation (4), where r ij is r i − r j , and A ij is the interaction matrix under the condition that i = j. If i and j are the same, then the interaction matrix can be simplified as α −1 i . Subsequently, Equation (2) can be transformed into a three-dimensional complex linear equation; hence, the P of each dipole can be calculated.
In the case of the absorption, attenuation, and scattering, cross-sections can be calculated as shown in Equations (5)-(7) using the calculated P, where * denotes complex conjugation.
Finally, the absorption, attenuation, and scattering efficiencies can be calculated using Equation (8), where r e f f is the effective radius of the particle, and V is the volume of the particle.
r e f f = 3V 4π

Heat Transfer Model and Optical Properties
In this study, the Pennes bioheat equation [40], which is widely used in the field of bio-heat transfer, was used for the thermal analysis of biological tissues. In this equation, it is assumed that the heat generated by blood and metabolism is uniformly generated in the biological tissue, as expressed in Equation (10).
where ρ is the density; c p is the specific heat; T is the temperature; and k m is the thermal conductivity of the medium. q b and q met are the heat generated by blood flow and metabolism, respectively. In this study, the skin surface was irradiated with a Gaussian profile laser, and the temperature was increased via the photothermal effect. Therefore, the heat generated by the laser must be simultaneously considered [41]. In addition, the heat generated by metabolism and the blood flow is insignificant compared with the heat generated by the laser, and Equation (10) is written as shown in Equation (11) to confirm the study results under steady-state conditions. − k m ∇ 2 T = q l (11) where µ abs is the light absorption coefficient of the medium; P l is the laser intensity; and r l is the laser radius. µ tot represents the total light attenuation coefficient of the medium and is the sum of the absorption coefficient (µ abs ) and attenuated scattering coefficient (µ sca ) of the medium. The optical properties of the medium should be calculated separately for normal and tumor tissues. For the former, because AuNPs are not injected into them, only their optical properties are considered. However, for the tumor tissues, because various shapes of AuNPs are injected inside them, the corresponding optical properties must be considered simultaneously.
µ abs = µ abs,n + µ abs,m , µ sca = µ sca,n + µ sca,m Equations (13)- (14) are equations that represent the optical properties of AuNPs [42], where f v is the volume fraction of AuNPs in the tumor; g is the anisotropy coefficient; and Q abs and Q sca are the absorption and scattering efficiencies of the particles, respectively. Finally, the optical properties of the medium containing AuNPs were calculated as the sum of the optical properties of the AuNPs and the medium, as shown in Equation (15). The optical property of the medium is not a calculated value, but an intrinsic property of the material. The optical efficiency of AuNPs was calculated using the DDA method described above, and it was assumed that the injected AuNPs were uniformly distributed in the medium [43].

Apoptotic Variables
In this study, three variables proposed by Kim et al. [34] were used to quantitatively confirm the therapeutic effect of photothermal treatment.
First, the apoptosis ratio (θ A ), which represents the quantity of tumor tissues that undergo apoptosis at a temperature range of 43 • C-50 • C, was calculated as the ratio of the volume of the tumor tissue to the tissue volume entering the temperature band corresponding to apoptosis. For example, if all regions of the tumor are in the temperature band of apoptosis, then θ A is 1.
Next, the thermal hazard value (θ H ), which quantitatively represents the amount of thermal damage to normal tissues, is the ratio between the normal tissue volume and weighted sum that affects each phenomenon based on the temperature of the biological tissue, as shown in Equation (17). The phenomena that occur in each temperature range in the biological tissue are shown in Table 1 [44,45]. The minimum value of θ H was 1, indicating that all regions of the normal tissue corresponded to 37 • C to 43 • C. As shown in Equation (17), the result for θ H varies depending on the range of normal tissues for calculating thermal damage. Accordingly, in this study, the range of the normal tissue surrounding the tumor tissue was selected as 50% of the length of the tumor tissue.
Finally, the effective apoptosis ratio (θ e f f ), which quantitatively represents the maximal occurrence of apoptosis in tumor tissues and the minimal occurrence of thermal damage to surrounding normal tissues simultaneously (i.e., the ultimate purpose of photothermal therapy), is the ratio between θ A and θ H , as shown in Equation (18).
Because θ e f f is the relative value between the ratio of the volume corresponding to the temperature band where apoptosis occurs in the tumor tissue and the amount of thermal damage to the surrounding normal tissues, more effective treatment is induced as θ e f f increases. Finally, the optimal conditions for photothermal therapy can be determined using θ e f f .

Validation of Numerical Model
To verify the numerical analysis model to be used in this study, the studies of Soni et al. [46] and Ren et al. [47] were used. The study assumed that a cylindrical tumor was located from the surface in cylindrical normal tissues as shown in Figure 1, and assumed that AuNPs in the tumor were uniformly distributed. Normal tissue has a radius of 20 mm and a depth of 10 mm, and a tumor tissue has a radius of 10 mm and a depth of 5 mm. The radius of the irradiated laser was selected to be 10 mm, which is the same as the radius of the tumor, and the intensity of the laser was 0.5 W/cm 2 . In addition, convection with air exists at the skin surface, and the initial temperature of all media was assumed to be 37 • C. The thermal and optical properties of normal and tumor tissues are summarized in Table 2.   Figure 2 is a graph of the validation results with the numerical analysis model in this study. The temperature along the radial direction when the depth was 0 mm and 5 mm in the central part of the tumor was confirmed, and the error was within about 1.4%. Through this, it was confirmed that the numerical analysis model used in this study was valid.

Numerical Investigation
In this study, a numerical analysis of photothermal therapy was conducted on a skin structure composed of four layers that included squamous cell carcinoma, as shown in Figure 3. The radius and depth of the entire normal tissue were 10 and 6 mm, respectively It was assumed that tumor tissues with a radius of 5 mm and a depth of 2 mm were presen

Numerical Investigation
In this study, a numerical analysis of photothermal therapy was conducted on a skin structure composed of four layers that included squamous cell carcinoma, as shown in Figure 3. The radius and depth of the entire normal tissue were 10 and 6 mm, respectively. It was assumed that tumor tissues with a radius of 5 mm and a depth of 2 mm were present on the skin surface. In addition, the irradiated laser exhibited a Gaussian distribution and had the same radius of 5 mm as the tumor tissue. The thickness of each skin layer and the thermal properties of all the media are summarized in Table 3. For the optical properties of the tumor and skin layer, the studies of Meglinski et al. [19] and Salomatina et al. [48] were referenced since values at various wavelengths were required.  Using the numerical analysis model, the thermal behavior of tumor tissues and the surrounding normal tissues for various laser intensities, volume fractions of injected AuNPs, and AuNP types were confirmed, and the conditions are shown in Table 4. The laser intensity was set from 0 to 1.2 W with 0.002 W intervals. Meanwhile, six types of AuNPs with diameters ranging from 10 to 50 nm were selected at intervals of 5 nm, and the volume fraction of injected AuNPs in the tumor was classified into four stages from 10 −3 to 10 −6 . Size of AuNP (r e f f ) 10 to 50 nm 9 Interval: 5 nm Volume fraction of AuNP ( f v ) 10 −3 to 10 −6 4 Interval: 10 −1 As described above, six AuNPs were injected into the tumor, as shown in Figure 4. In the rod-type AuNPs, the aspect ratio was fixed at 6.67; in the shell-type AuNPs, the difference between the outer diameter and the inner diameter was fixed at 4 nm; and in the prism-type AuNPs, the aspect ratio was fixed at 1.67. Finally, for various shapes of AuNPs, their absorption and scattering efficiencies were calculated in the wavelength range of 500 to 1500 nm using the DDA method. Finally, a numerical analysis of photothermal treatment was performed on the skin layer containing squamous cell carcinoma under the various conditions described above. For all cases, the temperature distribution of tumor tissues and the surrounding normal tissues was obtained, and the optimal treatment conditions were derived to maximize the temperature band corresponding to the apoptosis of tumor tissues while minimizing thermal damage to the surrounding normal tissues.

Derivation of Optical Properties for Various AuNP Types
The optical efficiencies of AuNPs with various shapes were calculated using the DDA method. Figure 5 shows the absorption efficiency and scattering efficiency of the rod-type AuNPs. Calculations were performed for a wavelength range of 500 to 1500 nm, and it was confirmed that a laser wavelength that results in the optimal efficiency exists for each r e f f . Analysis was performed on AuNPs of six different shapes, and the wavelength (λ max ) values with the maximum absorption efficiency at various shapes and r e f f were recorded. The absorption and scattering efficiencies at the corresponding wavelengths are summarized in Table 5. In the case of r e f f , after the size of the particles is determined, it is calculated through Equation (9), and the size of the particles having a value as close as possible to the r e f f range presented in this study was set. As shown in Table 5, the absorption efficiency of the rod type was generally greater than that of the other types of AuNPs. In addition, the absorbance efficiency increased with r e f f increased generally. However, for the rod and shell types, the r e f f with the optimal absorbance efficiency existed depending on r e f f . Finally, the temperature distribution inside the tumor and normal tissues when the volume fraction of AuNPs was changed to achieve the optimal optical efficiency conditions based on the shape of various AuNPs was obtained.  Figure 6 shows the temperature distribution of the normal and tumor tissues when the laser intensity was 0.4 W and the volume fraction of injected AuNPs was 10 −3 and 10 −6 , separately. The injected AuNPs were of the rod type-their aspect ratio and r e f f were 6.67 and 50 nm, respectively. When the volume fraction of AuNPs in the tumor was high ( f v = 10 −3 ), the light absorption coefficient of the tumor tissue including AuNPs was high. Therefore, a significant amount of laser energy was absorbed from the medium. Accordingly, the temperature of the tumor region increased significantly, as shown in Figure 6a. In addition, the temperature of the surrounding normal tissue increased because of heat conduction from the tumor tissue. In contrast, when the volume fraction of AuNPs in the tumor was low ( f v = 10 −6 ), the light absorption coefficient of the tumor tissue was low; therefore, even when a laser of the same intensity was used, the increase in temperature was insignificant, as shown in Figure 6b. The results confirmed that the temperature distribution of the tumor and surrounding normal tissues differed depending on the volume fraction of AuNPs in the tumor. Therefore, in this study, by obtaining the temperature distribution in the medium based on various shapes and volume fractions of AuNPs, the degree of thermal damage to the surrounding normal tissue as well as the corresponding temperature range that caused apoptosis were quantitatively confirmed.

Apoptosis Ratio
In photothermal therapy, laser is irradiated to the affected area to kill the tumor tissue by increasing the temperature via the photothermal effect. From a temperature perspective, the temperature band corresponding to apoptosis must be maintained to prevent necrosis; therefore, the extent to which the apoptosis temperature band is maintained must be determined quantitatively by verifying the temperature distribution in the tumor tissue. Figure 7 shows a graph of the apoptosis ratio (θ A ) as functions of P l and f v for the rod-type AuNPs with various r e f f . As shown in the graph, as r e f f increased, the intensity of the laser with the maximum θ A increased. This is because, when calculating the absorption coefficient, as Q a does not increase at the same rate as r e f f , the absorption coefficient decreases, and a higher laser intensity is required to achieve the corresponding temperature range in which apoptosis occurs. In addition, it was confirmed that as f v decreased, the intensity of the laser with the maximum θ A increased. This is because as f v decreased, the light absorption coefficient of the medium decreased; as such, the amount of heat absorbed by the medium decreased. Hence, the amount of heat generated by the laser must be increased to achieve the temperature band in which apoptosis occurs. Based on Figure 5, the condition in which the apoptosis ratio becomes 1 for all cases can be obtained.

Apoptosis Ratio
In photothermal therapy, laser is irradiated to the affected area to kill the tumor tissue by increasing the temperature via the photothermal effect. From a temperature perspective, the temperature band corresponding to apoptosis must be maintained to

Thermal Hazard Value
One of the most significant goals of photothermal therapy is to maintain a temperature range of 43 • C-50 • C in the tumor, which is known to cause apoptosis. Although direct heating using laser is not performed to the surrounding normal tissues, temperature increase is inevitable due to conduction heat transfer from the tumor tissue. Hence, the death of surrounding normal tissues due to temperature increase must be minimized by verifying the temperature distribution of not only the tumor tissue, but also the surrounding normal tissues. Accordingly, in this study, the amount of thermal damage to the surrounding normal tissues was quantitatively confirmed based on the thermal hazard value (θ H ) (Equation (17)). Figure 8 shows θ H as functions of P l and f v for sphere-type AuNPs with various r e f f . As shown in the graphs, as f v increased, θ H at the same P l was high. This is because when f v increased, the light absorption coefficient of the medium increased and the amount of heat absorbed by the tumor tissue increased; therefore, the temperature of the surrounding normal tissue increased higher, and the amount of thermal damage increased. Meanwhile, as r e f f increased, θ H decreased. This is because as r e f f increased, the light absorption coefficient of the medium decreased, as described in Section 3.3, which implies that the amount of heat absorbed by the medium decreased at the same laser intensity, and the temperature increased slowly.

Effective Apoptosis Ratio
Previously, the apoptosis ratio (θ A ) and thermal hazard value (θ H ) were used to define the conditions for maximizing the temperature range where apoptosis occurs in the tumor tissue and the relative ratio of thermal damage amount to the surrounding normal tissue. To confirm the distribution of the effective apoptosis ratio (θ e f f ), which combines the two variables above, calculations were performed based on the shape and size of the AuNPs as well as the volume fraction of the injected AuNPs. Figure 9 shows θ e f f as functions of P l and f v for shell-type AuNPs with various r e f f . As shown, the P l value that results in the maximum θ e f f existed for each f v , and that the P l with the maximum θ e f f decreased as f v increased. This is because, as described above, when f v increased, the light absorption coefficient of the medium increased and the amount of heat absorbed from the medium increased; therefore, the intensity of the laser required to maintain the temperature range where apoptosis occurs was decreased.
In addition, as r e f f increased, the intensity of the laser at which θ e f f was maximized increased. However, for the shell-type AuNPs, as shown in Figure 9, the intensity of the laser that afforded the optimal treatment effect was not constant because r e f f , which maximized the absorption efficiency for various r e f f , existed. In addition, it was confirmed that an excessive increase in the laser intensity decreased θ e f f due to the increase in θ H . The results confirmed that certain values of laser intensity and volume fraction of injected AuNPs yielded the optimal therapeutic effect for the different shapes and sizes of AuNPs. Figure 8 shows the derivation results of θ e f f for each r e f f for six types of AuNPs investigated in this study. Based on the analysis, the volume fraction of injected AuNPs that resulted in the optimal therapeutic effect for all AuNP shapes was derived based on f v = 10 −6 . As shown in Figure 10, the optimal values of r e f f and P l afforded the optimal therapeutic effect for various AuNP shapes. For the rod type, the absorption coefficient decreased as r e f f increased; consequently, P l increased, and the optimal θ e f f was obtained. For the shell type, it was confirmed that the P l value that afforded the optimal therapeutic effect was minimized as the absorption coefficient was maximized when r e f f was 25 nm.   For the pyramid and cube types, as r e f f increased, the absorption coefficient exhibited similar values; however, when r e f f was 40 nm, the absorption coefficient decreased and the optimal point of P l decreased. For the sphere and prism types, the calculated absorption coefficients were similar as r e f f increased; as such, the optimum points of P l based on r e f f were similar. Table 6 summarizes the f v and P l values that afforded the maximum θ e f f based on the shape and size of each AuNP.

Conclusions
In this study, a numerical study based on heat transfer theory was performed to investigate photothermal therapy using AuNPs on an actual skin layer containing squamous cell carcinoma. The temperature distribution inside the biological tissue was calculated using the Pennes bioheat equation, and the optical properties of the AuNPs were calculated using the discrete dipole approximation method.
The optimal treatment conditions for photothermal therapy were confirmed through the apoptosis ratio, which quantitatively confirms the volume ratio corresponding to the temperature band of apoptosis in tumors, the thermal hazard value that confirms the thermal damage of surrounding normal tissues, and the effective apoptosis ratio that considers the above two variables at the same time.
Finally, the AuNP shape, r e f f , f v , λ max , and P l that afforded the optimal therapeutic effect were obtained; hence, they can be utilized as the optimal treatment conditions for performing photothermal treatment in the future. In addition, future studies will not simply assume a cylindrical tumor, but rather conduct studies on tumors of various shapes.  Institutional Review Board Statement: Not applicable.

Informed Consent Statement: Not applicable.
Data Availability Statement: Data sharing is not applicable to this article.

Conflicts of Interest:
The authors declare that they have no known competing financial interests or personal relationships that could have influenced the work reported herein.

Abbreviations
The