MD simulations were performed on PTEN mutant systems using the updated CHARMM program with the modified all-atom force field (version 36m) (Brooks et al., 2009 (link); Huang et al., 2017 (link); Klauda et al., 2010 (link)). Our computational studies closely followed the same protocol as in our previous works (Grudzien et al., 2022 (link); Haspel et al., 2021 (link); Jang et al., 2016a (link); Jang et al., 2019 (link); Jang et al., 2016b (link); Jang et al., 2021 (link); Jang et al., 2020 (link); Liao et al., 2020 (link); Liu et al., 2022b (link); Liu et al., 2022c (link); Maloney et al., 2021 (link); Maloney et al., 2022 (link); Weako et al., 2021 (link); Zhang et al., 2021a (link); Zhang et al., 2021b (link)). Prior to productions runs, a series of minimization and dynamics cycles were performed for the solvents including ions and lipids with a harmonically restrained protein backbone until the solvent reached 310 K. Next, preequilibrium simulations with dynamic cycles were performed while gradually releasing the harmonic restraints on the backbones of PTEN mutants. The particle mesh Ewald (PME) method was used to calculate the long-range electrostatic interaction, and the van der Waals (vdW) interactions using switching functions with the twin range cutoff at 12 Å and 14 Å were calculated for the short-range interaction between atoms. In the production runs, the Nosé-Hoover Langevin piston control algorithm was used to sustain the pressure at 1 atm, and the Langevin thermostat method was employed to maintain the constant temperature at 310 K. The SHAKE algorithm was applied to constrain the motion of bonds involving hydrogen atoms. Simulations were performed for eight mutant systems each with 1 μs, and additional simulations for the same systems were also performed to check reproducibility. The production runs were performed with the NAMD parallel-computing code (Phillips et al., 2005 (link)) on a Biowulf cluster at the National Institutes of Health (Bethesda, MD). The result analysis was performed in the CHARMM program (Brooks et al., 2009 (link)). To determine the most populated conformation, the ensemble clustering in Chimera (Pettersen et al., 2004 (link)) was implemented to obtain the conformational representatives. The weighted implementation of suboptimal path (WISP) (Van Wart et al., 2014 (link)) algorithm was used to identify the allosteric signal propagation pathways through the protein. To observe conformational changes in proteins, the normal mode analysis (NMA) and principal component analysis (PCA) were conducted by the ProDy program (Bakan et al., 2011 (link)). In the analysis, averages were taken afterward discarding the first 200 ns trajectories.