We prsent a multiscale method for computing the effective behavior of a class of stiff and highly oscillatory ordinary differential equations. The oscillations may be in resonance with one another and thereby generate hidden slow dynamics. The proposed method relies on correctly tracking a set of slow variables whose dynamics is closed up to $\epsilon$ perturbation, and is sufficient to approximate any variable and functional that are slow under the dynamics of the ODE. This set of variables is detected numerically as a preprocessing step in the numerical methods. Error and complexity estimates are obtained. The advantages of the method is demonstrated with a few examples, including a commonly studied problem of Fermi, Pasta, and Ulam.