Aiming at the problems of bispectral analysis when applied to machinery fault diagnosis, a machinery fault feature extraction method based on sparseness-controlled non-negative tensor factorization (SNTF) is proposed. First, a non-negative tensor factorization(NTF) algorithm is improved by imposing sparseness constraints on it. Secondly, the bispectral images of mechanical signals are obtained and stacked to form a third-order tensor. Thirdly, the improved algorithm is used to extract features, which are represented by a series of basis images from this tensor. Finally, coefficients indicating these basis images' weights in constituting original bispectral images are calculated for fault classification. Experiments on fault diagnosis of gearboxes show that the extracted features can not only reveal some nonlinear characteristics of the system, but also have intuitive meanings with regard to fault characteristic frequencies. These features provide great convenience for the interpretation of the relationships between machinery faults and corresponding bispectra.