## 2 Materials and Methods  
### Monte Carlo Method and Heat Transfer Model  
When a laser is irradiated on a tissue, the laser is absorbed in the vicinity of the irradiation, and scattering occurs simultaneously by the medium. To analyze these behaviors, the Monte Carlo method is used in the field of bioheat transfer [26]. This technique probabilistically calculates the behavior of photons irradiated through a number selected at random, considering not only photons absorbed by the medium but also photons scattered at the same time.  
Variables computed through random numbers are the distance and angle that the photon travels. A photon's moving distance per 1 step is determined as in Equation (1) through the attenuation coefficient, which is the sum of absorption coefficient and scattering coefficient, where \(S\) is the photon's moving distance per 1 step, \(\xi\) is the random number generated between 0 and 1, and \(\mu_{t},\mu_{a},\ \mu_{s}\) are the attenuation coefficient, absorption coefficient, and scattering coefficient of the medium, respectively.  
\[S=\frac{-\ln(\xi)}{\mu_{t}}\ (\mu_{t}=\mu_{a}+\mu_{s}) \tag{1}\]  
When a photon travels 1 step, the energy of the photon decreases due to absorption in the medium. The energy reduction of a photon is determined by the ratio of the absorption coefficient and the attenuation coefficient of the medium as shown in Equation (2). Here, \(W\) is the energy weight of the photon.  
\[\Delta W=W\frac{\mu_{a}}{\mu_{t}} \tag{2}\]  
In the case of the angle in which photon travels, the deflection angle and the azimuthal angle are calculated as Equations (3) and (4), respectively, by the anisotropy coefficient, where \(g\) is the anisotropy coefficient, \(\theta\) represents deflection angle, and \(\psi\) is the azimuthal angle.  
\[\cos\theta=\begin{cases}\frac{1}{2g}\bigg{\{}1+g^{2}-\left[\frac{1-g^{2}}{1-g+2g \xi}\right]^{2}\bigg{\}}&if\ g>0\\ 2\xi-1&if\ g=0\end{cases}\ (0\ <\theta<\pi) \tag{3}\]  
\[\psi=2\pi\xi(0\ <\ \xi\ <\ 1) \tag{4}\]  
Once deflection angle and azimuthal angle are determined, angles the photon travels in the Cartesian coordinate system can be determined as in Equations (5)-(7), where, \(\mu_{x}\mu_{y}\mu_{z}\) are the directional cosines for each axis.  
\[\mu_{x}^{\prime}=\frac{\sin\theta}{\sqrt{1-\mu_{z}^{2}}}(\mu_{x}\mu_{z}\cos \psi-\mu_{y}\sin\psi)+\mu_{x}\cos\  
When the final path of photon movement is determined, the absorption distribution of heat in the medium is determined according to the intensity of the irradiated laser per unit area and the absorption coefficient of the medium. This allows the temperature distribution of the medium to be calculated. The temperature distribution is calculated using the thermal diffusion equation as shown in Equation (11), where \(q\) is the volumetric heat source, \(k\) is the thermal conductivity of the medium, \(\rho\) is the density, and \(c_{\eta}\) is the heat capacity.  
\[\frac{\partial T}{\partial t}=\frac{q+\
abla\cdot(k\
abla T)}{\rho c_{\eta}} \tag{11}\]  
The equation is calculated through the explicit finite element method as shown in Equation (11) below, where \(P\) is the intensity of the laser, \(F\) is Fluence rate, and \(dx,\;dy,\;dz\) are infinitesimal length for each axis [27].  
\[\begin{split}\Delta T=\frac{\Delta t}{\rho c_{\eta}}\Big{(}& P\mu_{a}Fdxdydz+(T_{x^{-}}-T)\frac{2 kk_{x^{-}}}{k+k_{x^{-}}}\frac{dydz}{dx}+(T_{x^{+}} -T)\frac{2 kk_{x^{+}}}{k+k_{x^{+}}}\frac{dydz}{dx}\\ &+(T_{y^{-}}-T)\frac{2 kk_{y^{-}}}{k+k_{y^{-}}}\frac{dxdz}{dy}+(T_ {y^{+}}-T)\frac{2 kk_{y^{+}}}{k+k_{y^{+}}}\frac{dydz}{dy}\\ &+(T_{z^{-}}-T)\frac{2 kk_{z^{-}}}{k+k_{z^{-}}}\frac{dxdy}{dz}+(T _{z^{+}}-T)\frac{2 kk_{x^{+}}}{k+k_{z^{+}}}\frac{dxdy}{dz}\Big{)}\end{split} \tag{12}\]