In the expression above, the volume of the system shows negligible fluctuations during the simulation, hence it is considered constant (the average volume of the simulation box was used to calculate rates). We use a linear fit of the transient period at the start of the simulation (until the number of clusters reach a maximum) and then we calculate J from the slope. This corresponds to the rate of nucleation after it is normalized by V. Then, we model the decay in the number of clusters by fitting an exponential function. The growth rate is calculated by the exponent normalized by V. We use NumPy for the fitting procedure. At last, we consider fitting errors as the root mean square deviation of the linear fit for the nucleation rate, and the one standard deviation error on the exponent for the growth rate. Principal Component Analysis. We calculate different sets of data which characterize our analysis of the early stages of MOF self-assembly. A dimensionality reduction can be performed while retaining meaningful information by projecting data on principal components that possess most of the variation of the dataset . In this effort, we consider the average values of quantities corresponding to the largest cluster along the trajectory of each simulation. We standardize data to have a mean of zero and standard deviation equal to the unity. This is done to combine data that have different units and magnitude . Then, we identify eigenvectors and eigenvalues of the covariance matrix . At last, we project data on the eigenvectors (principal components) with the highest eigenvalues.