We propose and study a scheme combining the finite element method and machine learning techniques for the numerical approximations of coupled nonlinear forward-backward stochastic partial differential equations (FBSPDEs) with homogeneous Dirichlet boundary conditions. Precisely, we generalize the pioneering work of Dunst and Prohl [SIAM J. Sci. Comp., 38(2017), 2725--2755] by considering general nonlinear and nonlocal FBSPDEs with more inclusive coupling; self-contained proofs are provided and different numerical techniques for the resulting finite-dimensional equations are adopted. For such FBSPDEs, we first prove the existence and uniqueness of the strong solution as well as of the weak solution. Then the finite element method in the spatial domain leads to approximations of FBSPDEs by finite-dimensional forward-backward stochastic differential equations (FBSDEs) which are numerically computed by using some deep learning-based schemes. The convergence analysis is addressed for the spatial discretization of FBSPDEs, and the numerical examples, including both decoupled and coupled cases, indicate that our methods are quite efficient.