An implicit, multi-step transversal linearization (MTL) family of methods is proposed for accurate, efficient and numerically stable integration of non-linear oscillators of interest in structural dynamics. The presently developed method is a multi-step extension of the locally transversal linearization (LTL) method proposed earlier by Roy (2001). The MTL-based linearization is achieved through a non-unique replacement of the nonlinear part of the vector field by a conditionally linear interpolating expansion of known accuracy, whose coefficients contain the discretized state variables defined at a set of discretization points. Any available interpolating expansion (such as the ones based on polynomials, wavelets or distributed approximating functionals) may be used to achieve linearization and, indeed, an appropriate choice may be based on the kind of system trajectories being simulated. It is shown that the tangent spaces of the non-linear and linearized systems are transversal almost everywhere in the associated phase space, and, in particular, at the points of discretization. The unknown state vectors are thus determined by constructing a set of algebraic constraint equations, which ensure transversal intersections of the linearized and non-linear solution manifolds. Since an exact solution of the linearized dynamical system can be constructed, the formal accuracy of the MTL method as a function of the time step-size depends only on that of the interpolating expansion, used to replace the nonlinear terms. Moreover, if the original system is Hamiltonian, then so is its MTL-based linearized form. Finally, a limited numerical illustration of the method is provided for a few low dimensional nonlinear oscillators in their periodic and chaotic regimes.