We propose a new method for calculating multipole matrix elements between wave eigenfunctions of the one-dimensional Schrödinger equation. The method is based on the transition to the auxiliary third- and fourth-order equations, to which an analogue of the Laplace transform is then applied. The resulting recursive procedure allows us to evaluate matrix elements starting with a number of eigenvalues that are assumed to be known and several basis matrix elements. As an example, we consider the multipole matrix elements between the wave functions of the harmonic and nonharmonic oscillators.