In this study, a supervised machine learning approach called Gaussian process regression (GPR) was applied to approximate optimal bi-impulse rendezvous maneuvers in the cis-lunar space. We demonstrate the use of the GPR approximation of the optimal bi-impulse transfer to patch points associated with various invariant manifolds in the cis-lunar space. The proposed method advances preliminary mission design operations by avoiding the computational costs associated with repeated solutions of the optimal bi-impulsive Lambert transfer because the learned map is computationally efficient. This approach promises to be useful for aiding in preliminary mission design. The use of invariant manifolds as part of the transfer trajectory design offers unique features for reducing propellant consumption while facilitating the solution of trajectory optimization problems. Long ballistic capture coasts are also very attractive for mission guidance, navigation, and control robustness. A multi-input single-output GPR model is presented to represent the fuel costs (in terms of the ΔV magnitude) associated with the class of orbital transfers of interest efficiently. The developed model is also proven to provide efficient approximations. The multi-resolution use of local GPRs over smaller sub-domains and their use for constructing a global GPR model are also demonstrated. One of the unique features of GPRs is that they provide an estimate of the quality of approximations in the form of covariance, which is proven to provide statistical consistency with the optimal trajectories generated through the approximation process. The numerical results demonstrate our basis for optimism for the utility of the proposed method.