<p>This study introduces a novel method for calculating the thermal conductivity of graphene using a Monte Carlo approach to evaluate anisotropic three-phonon interactions. The phonon dispersion relation is derived using a force constant model that incorporates up to fifth-order nearest neighbor interactions, while the phonon density of states (DOS) is computed via a generalized Gilat–Raubenheimer method. A quantitative relationship for the scaling exponent of the specific heat capacity at low temperatures is established, emphasizing the unique two-dimensional characteristics of graphene. To address anisotropic effects, the Monte Carlo approach efficiently identifies three-phonon combinations that adhere to the conservation laws of energy and momentum. The findings highlight the pivotal role of anisotropic phonon interactions in graphene’s thermal conductivity. The thermal conductivity values obtained through the iterative method exhibit strong agreement with previous three-phonon calculations, thereby validating the model. Nevertheless, discrepancies with experimental data suggest that incorporating higher-order phonon processes, such as four-phonon scattering, may further improve predictive accuracy.</p>