We propose an accurate multi-step scheme on time-space grids for solving backward stochastic differential equations (BSDEs). For high accuracy, we preprocess the target equation before discretization, changing it from a stochastic problem to a deterministic problem. Then, we apply a multi-step method in time and the Gauss-Hermite quadrature approximation in space to construct the multi-step scheme. Error estimates are rigorously proved for the semi-discrete version of the proposed scheme for BSDEs with certain type of simplified generator functions. For efficiency, we improve our scheme by applying a new kind of random process, called the Gauss-Hermite process, by which the time-space domain needed for solving BSDEs at a particular point of interest can be significantly compared to the domain used by the multi-step scheme. Error estimates are rigorously derived to verify that the improved scheme can keep the essentially the same accuracy as that of the multi-step scheme. Furthermore, in order to break the curse of dimensionality, we apply hierarchical bases and associated sparse grids to effect the spatial interpolation. Combining these techniques, high dimensional BSDEs can be solved efficiently with acceptable accuracy.