Frequency-response analysis is a crucial tool in the study of dynamical systems that manifest in many engineering problems. Evaluating it computationally can be expensive, especially if many frequency samples or multiple model parameters are to be considered. Interpolation-based model order reduction (MOR) methods offer effective strategies to reduce the computational cost of a frequency-response analysis by computing cheap-to-evaluate surrogate models. While effective, interpolatory MOR methods present several challenges which have, thus far, not been fully resolved. These include, for instance, (a) the automatic selection of expansion points in the parameter domain of interest and (b) the adaptive choice of the dimension of the surrogate model for a given tolerance in the approximation error. Existing solutions to the above issues require adept strategies to estimate the approximation error, as the solution of the full-order model is typically not available. Although estimating the approximation error is computationally inexpensive compared to computing the exact error, its repeated evaluation at many samples can cause a heavy computational burden. An adaptive error estimation strategy can offer a substantial gain in such a situation. To be practically applicable, a reduction algorithm should also require as few tuning parameters as possible. In this contribution, a framework is presented for the automatic and adaptive computation of reduced-order models using interpolation-based MOR methods. The individual parts of the algorithm, i.e. basis computation, error estimation, and parameter sampling, are modular and can be exchanged. This offers a thorough and flexible comparison of the individual approaches. Their versatility and efficacy are assessed and guidelines are identified towards establishing a robust black-box reduction approach. Various parametrized dynamical systems modeling problems from applications such as vibro-acoustics and poro-elasticity are used as benchmark examples.