High temperature superconducting (HTS) rotatory flux pump, or so called HTS dynamo, can output none-zero time-averaged DC voltage and charge the rest of the circuit if a closed loop has been formed. This type of flux pump is often employed together with HTS coils, where the HTS coils can potentially work in persistent current mode, and act like electromagnets with considerable magnetic field, having wide range of applications in industry. The output characteristics of HTS rotary flux pumps have been extensively explored through experiments and finite element method (FEM) simulations, yet the work on constructing statistical models as an alternative approach to capture key characteristics has not been studied and published. A 2D FEM program has been used to model the HTS rotatory flux pumps and evaluate the effects of different factors upon the output voltage through parameter sweeping and analysis of variance. Typical design considerations, including operation frequency, air gap, HTS tape width and remanent flux density have been investigated, in particular the bilateral effect of HTS tape width has been explained by looking at the averaged integration of the electric field over the tape. Based on the data obtained from various simulations, regression analysis has been conducted through a collection of machine learning methods and demonstrated that the output voltage of a rotary flux pump can be obtained promptly with satisfactory accuracy via Gaussian process regression, aiming to provide a novel approach for future research and powerful design tool for industrial applications using HTS rotary flux pump devices.