Oceanic eddy kinetic energy (EKE) is a key quantity for measuring the intensity of mesoscale eddies and for parameterizing eddy effects in ocean climate models. Three decades of satellite altimetry observations allow a global assessment of sea surface information. However, the subsurface EKE with spatial filter has not been systematically studied due to the sparseness of subsurface observational data. The subsurface EKE can be inferred both theoretically and numerically from sea surface observations but is limited by the issue of decreasing correlation with sea surface variables as depth increases. In this work, inspired by the Taylor-series expansion of subsurface EKE, a multiple-input neural network approach is proposed to reconstruct the subsurface monthly mean EKE from sea surface variables and subsurface climatological variables (e.g., horizontal filtered velocity gradients). Four neural networks are trained on a high-resolution global ocean reanalysis dataset, namely, surface-input fully connected neural network model (FCNN), surface-input Residual neural network model (ResNet), multiple-input fully connected neural network model (MI-FCNN), and multiple-input residual neural network model (MI-ResNet). The proposed MI-FCNN and MI-ResNet models integrate the surface input variables and the vertical profiles of subsurface variables. The MI-ResNet model outperforms the FCNN, ResNet, and MI-FCNN models, and traditional physics-based models in both regional and global reconstruction of subsurface EKE in the upper 2000 m. In addition, the MI-ResNet model performs well for both regional and global observational data based on transfer learning. These findings reveal the potential of the MI-ResNet model for efficient and accurate reconstruction of subsurface oceanic variables.