The 2004 Mw 9.1?9.3 Sumatra-Andaman earthquake is one of the largest earthquakes of the modern instrumental era. Despite considerable efforts to analyze this event, the different available observations have proven difficult to reconcile in a single finite-fault slip model. In particular, the critical near-field geodetic records contain variable and significant postseismic signal (between 2 weeks’ and 2 months’ worth), while the satellite altimetry records of the associated tsunami are affected by various sources of uncertainties (e.g., source rupture velocity and mesoscale oceanic currents). In this study, we investigate the quasi-static slip distribution of the Sumatra-Andaman earthquake by carefully accounting for the different sources of uncertainties in the joint inversion of available geodetic and tsunami data. To this end, we use nondiagonal covariance matrices reflecting both observational and modeling uncertainties in a fully Bayesian inversion framework. Modeling errors can be particularly large for great earthquakes. Here we consider a layered spherical Earth for the static displacement field, nonhydrostatic equations for the tsunami, and a 3-D megathrust interface geometry to alleviate some of the potential epistemic uncertainties. The Bayesian framework then enables us to derive families of possible models compatible with the unevenly distributed and sometimes ambiguous measurements. We infer two regions of high fault slip at 3°N?4°N and 7°N?8°N with amplitudes that likely reach values as large as 40 m and possibly larger. These values are a factor of 2 larger than typically found in previous studies?potentially an outcome of commonly assumed forms of regularization. Finally, we find that fault rupture very likely involved shallow slip. Within the resolution provided by the existing data, we cannot rule out the possibility that fault rupture reached the trench.