As commonly understood, the noise spectroscopy problem-characterizing the statistical properties of a noise process affecting a quantum system by measuring its response-is mathematically ill-posed, in the sense that no unique noise process leads to a set of responses unless extra assumptions are taken into account. Ad-hoc solutions assume an implicit structure, which is often never determined. Thus, it is unclear when the method will succeed or whether one should trust the solution obtained. Here, we propose to treat the problem from the point of view of statistical estimation theory. We develop a Bayesian solution to the problem which allows one to easily incorporate assumptions which render the problem solvable. We compare several numerical techniques for noise spectroscopy and find the Bayesian approach to be superior in many respects. initial states, to apply suitable control sequences, and to measure a set of convenient observables. The main difficulty is that noise correlations influence the dynamics of the quantum system in a highly nonlinear way. Thus, inferring these correlations in detail from the response of the quantum system is generally an ill-posed problem, unless a priori information on the noise is assumed. Even when standard assumptions, such as Gaussian noise or a dephasing coupling are satisfied, the problem remains nonlinear and inverting it carries along a set of non-trivial complications that in turn constrain the type of noise that can be characterized. For example, in [8-12] a control-induced frequency comb approach is used in order to overcome the nonlinear character of the problem but it comes at the cost of being only effective when the noise correlations are smooth functions in frequency space.We propose that many of these problems can be ameliorated, or at least properly quantified, using a statistically principled approach. Within the statistical phrasing of the problem we provide a Bayesian solution [17], complete with a numerical implementation. We show the problem can be solved analytically in the limit of large of amounts of experimental data. At the other extreme-the small-data limit-a numerically stable Monte Carlo algorithm [18] approximates the full Bayesian solution. Our two approaches provide a robust solution to the software side of the noise spectroscopy problem, and can be used to improve state-of-the-art spectroscopy protocols [8][9][10][11][12]. These two regimes are schematically depicted in figure 1. Though the physical model we consider is simplified to illustrate the novel aspects of this work, we note that our approach is very versatile. We discuss this later.We summarize the performance of our approach for the large-and small-data limits in section 6, with complete details including all source code in the supplementary material. In the large data limit, our approach can give almost an order of magnitude improvement in performance over state-of-the-art estimation strategies. By contrast, in the small-data regime, we can even achieve two orders of magnitude improv...