This article presents a tool to quantify uncertainties in magnitude–depth () estimates for earthquakes associated with macroseismic intensity data. The tool is an open‐source code written in Python and is named quantifying uncertainties in earthquakes’ magnitude and depth (QUake‐MD). In QUake‐MD, uncertainties are propagated from the individual intensity data point (IDP) to the final magnitude (M), depth (H), epicentral intensity () solution. It also accounts for epistemic uncertainties associated with the use of different intensity prediction equations (IPEs). For each IPE, QUake‐MD performs a sequential least‐square inversion process to estimate the central M, H value. QUake‐MD then explores the uncertainties around this central M, H solution by constructing a probability density function (PDF) constrained to be consistent with the range of plausible epicentral intensity , a plausible depth range, and IDP uncertainties. The resulting PDFs of all IPEs provided to QUake‐MD are then stacked to obtain a final PDF of possible M, H, solutions representative of both data quality and IPE epistemic uncertainties. This tool is geared toward end users who would like to grasp a more complete understanding of the uncertainties associated with historical earthquake parameters beyond the classical standard deviation values proposed today in parametric earthquake catalogs. We apply QUake‐MD to two events of the SisFrance macroseismic database to illustrate the challenges involved in building realistic spaces of M, H, solutions reflecting the quality of the data and the epistemic uncertainties in IPEs.