Importance sampling can be highly efficient if a good importance sampling density is constructed. Although the parametric sampling densities centered on the design points are often good choices, the determination of the design points can be a difficult and inefficient task itself, especially for problems with multiple design points, or highly nonlinear limit state functions. This paper introduces a nonparametric importance sampling method based on the Markov chain simulation and maximum-entropy density estimation (MEDE). In the proposed method, Markov chain simulation is utilized to generate samples that distribute asymptotically to the optimal importance sampling density. A nonparametric estimation of the optimal importance sampling density is then obtained using the MEDE technique. The conventional MEDE method is difficult for multi-dimensional problems as it needs to solve a set of simultaneous nonlinear integral equations. This paper developed a new MEDE technique for multivariate dataset. The method starts with using histogram to approximate a density. The multi-dimensional histogram is converted into a series of one-dimensional conditional PDFs in each dimension and the density is reconstructed by means of orthogonal expansion. Thus, the solution of MEDE is converted to a set of coefficients of the Legendre polynomials. The new importance sampling method is illustrated and compared with the classical kernel-based importance sampling using a number of numerical and structural examples.