where 𝐸 >7Z denotes the lowest energy in the training set and 𝛥𝐸 defines the range of favorably weighted energies, which was set to 10 𝑘 [ 𝑇 for the fitting of two-body interactions and to 5 𝑘 [ 𝑇 for the fitting of three-body interactions upon careful experimentation. Both linear and nonlinear parameters were obtained by minimizing the (regularized) weighted sum of squared residuals (Eq. 11). The linear parameters (𝐶 V ) were obtained through singular value decomposition and nonlinear parameters were optimized using the simplex algorithm. Simulations Using the Many-body Potential. The three-body potentials were incorporated into LAMMPS 70 through MBX software, where forces and energy of NPs are calculated based on the developed potential. All simulations were carried out in Langevin dynamics at a temperature of 𝜀/𝑘 * , and a damping parameter of 0.0864 (𝑚𝜎 ! 𝜀 ⁄ ) 4/! was used so that the dynamics of the NPs matched those of polymer-grafted NPs in the CG MD simulations (Fig. ). Periodic boundary conditions were employed in all three directions of the simulation box. A velocity-Verlet algorithm with a time step of 0.02 (𝑚𝜎 ! 𝜀 ⁄ ) 4/! was used for integrating the equations of motion.