The T4 lysozyme system with toluene bound (as described in 2.3.1) was minimized in GROMACS 5.1 via steepest-descent, followed by 1 ns of NVT simulation and then 5 ns of NPT simulation at 1 atm and 300 K for equilibration. The leapfrog integrator was used with a 2 fs timestep and the bonds involving hydrogen constrained with LINCS. The system was then simulated for a total of 806 ns under the same NPT conditions, saving configurations to a trajectory file every 30ps. tICA was performed on the pairwise-distances of the toluene heavy-atoms and the alpha carbons of the binding site of the trajectory with a lagtime of 0.6ns to generate order parameters for MSM construction. Of the 210 initial dimensions, 22 dimensions were retainedenough to account for 95% of the kinetic variance in the data, and were scaled by the kinetic map scheme. These large number of initial dimensions were used to help ensure all relevant binding modes were separated. An initial MSM was estimated from the order parameters computed from the trajectory using PyEMMA, using 1000 microstates generated from k-means clustering and a lagtime of 6 ns. The MSM was coarse-grained into four macrostates using Perron-Cluster Cluster Analysis ++(PCCA+); full details are available in scripts deposited in the SI. Two random trajectory frames from each of the four macrostates were then used as the starting point for new simulations to further sample each identified binding mode and potentially generate additional transitions. These eight additional simulations were each run for 60 ns and combined with the longer run above to re-estimate a MSM, following the same sequence of steps, with these additional simulations added to better explore transitions out of each macrostate. The total amount of aggregated simulation time used for the final MSM was 1.286 µs spread across nine trajectories. Addi-tional simulation details can be found in the SI.