Research on hypersonic vehicles has become increasingly important worldwide in recent years. However, accurately simulating the dynamics of the nonequilibrium high-temperature reactions that are in the hypersonic flow around the vehicles presents a significant challenge as a large number of states and transitions are accessible even for the smallest atom-diatom reaction systems. It is quite difficult, sometimes even impossible, to exhaustively investigate all relevant combinations or determine high-dimensional analytical representations for the state-to-state reaction probabilities. In this study, we used Gaussian process regression (GPR) to fit a model based on only 807 QCT data for training. The confidence interval of the GPR prediction and the Kullback-Leibler (KL) divergence were used to help minimize the sampling amount of data for fitting the converged GPR model. The model aims to predict the state-to-state integral cross section (ICS) of the O + O2 → 3O dissociation reaction under random initial conditions (Et, v, j). In total, it took almost a month to obtain this converged GPR model, but it took only a few seconds to predict the ICS value for any initial condition. For 330 initial conditions not included in the training set, the mean-square error (MSE) between the QCT-calculated ICSs and the GPR-predicted ones is only 0.08 Å2 and the R2 is 0.9986, indicating that the GPR model can replace the direct expensive QCT calculation with high accuracy. Finally, we calculated the equilibrium dissociation rate coefficients based on the StS ICS values predicted by the GPR model, and the results were in good agreement with available experimental and theoretical results. Thus, this study provides an effective and accurate approach to the extensive direct state-to-state reaction dynamic calculations.