S U M M A R YThis paper applies Bayesian inversion, with evaluation of data errors and model parametrization, to produce the most-probable shear-wave velocity profile together with quantitative uncertainty estimates from microtremor array dispersion data. Generally, the most important property for characterizing earthquake site response is the shear-wave velocity (V S ) profile. The microtremor array method determines phase velocity dispersion of Rayleigh surface waves from multi-instrument recordings of urban noise. Inversion of dispersion curves for V S structure is a non-unique and non-linear problem such that meaningful evaluation of confidence intervals is required. Quantitative uncertainty estimation requires not only a non-linear inversion approach that samples models proportional to their probability, but also rigorous estimation of the data error statistics and an appropriate model parametrization. This paper applies a Bayesian formulation that represents the solution of the inverse problem in terms of the posterior probability density (PPD) of the geophysical model parameters. Markov-chain Monte Carlo methods are used with an efficient implementation of Metropolis-Hastings sampling to provide an unbiased sample from the PPD to compute parameter uncertainties and inter-relationships. Nonparametric estimation of a data error covariance matrix from residual analysis is applied with rigorous a posteriori statistical tests to validate the covariance estimate and the assumption of a Gaussian error distribution. The most appropriate model parametrization is determined using the Bayesian information criterion, which provides the simplest model consistent with the resolving power of the data. Parametrizations considered vary in the number of layers, and include layers with uniform, linear and power-law gradients. Parameter uncertainties are found to be underestimated when data error correlations are neglected and when compressional-wave velocity and/or density (nuisance) parameters are fixed in the inversion.Bayesian inversion of microtremor array data is applied at two sites in British Columbia, the area of highest seismic risk in Canada, to study the ability to recover an accurate V S profile in relatively deep and shallow geological settings on the Fraser River delta in Greater Vancouver and in Victoria, respectively. A well-resolved V S profile to at least 110 m depth is determined at the Fraser River delta site for a power-law gradient parametrization. At the Victoria site, a layer with low V S and a weak linear gradient is indicated to 15-18 m depth, above much higher velocity material. Invasive V S measurements from seismic cone penetration testing and/or surface-to-downhole methods are used to assess the reliability of the Bayesian microtremor inversion results, with excellent agreement obtained at both sites.