Frontera is made possible by National Science Foundation award OAC-1818253.
This research was carried out using resources provided by the Open Science
Grid, which is supported by the National Science Foundation and the U.S.
Department of Energy Office of Science. Additional work used ABACUS2.0, which
is part of the eScience center at Southern Denmark University. Simulations
were also performed on the SuperMUC cluster at the LRZ in Garching, on the
LOEWE cluster in CSC in Frankfurt, on the HazelHen cluster at the HLRS in
Stuttgart, and on the Pi2.0 and Siyuan Mark-I at Shanghai Jiao Tong
University. The computer resources of the Finnish IT Center for Science (CSC)
and the Finnish Computing Competence Infrastructure (FCCI) project are
acknowledged. This research was enabled in part by support provided by Compute
Ontario (http://computeontario.ca), Calcul Quebec (http://www.calculquebec.ca)
and Compute Canada (http://www.computecanada.ca).

The EHTC has received generous donations of FPGA chips from Xilinx Inc., under
the Xilinx University Program. The EHTC has benefited from technology shared
under open-source license by the Collaboration for Astronomy Signal Processing
and Electronics Research (CASPER). The EHT project is grateful to T4Science
and Microsemi for their assistance with Hydrogen Masers. This research has
made use of NASA’s Astrophysics Data System. We gratefully acknowledge the
support provided by the extended staff of the ALMA, both from the inception of
the ALMA Phasing Project through the observational campaigns of 2017 and 2018.
We would like to thank A. Deller and W. Brisken for EHT-specific support with
the use of DiFX. We thank Martin Shepherd for the addition of extra features
in the Difmap software that were used for the CLEAN imaging results presented
in this paper. We acknowledge the significance that Maunakea, where the SMA
and JCMT EHT stations are located, has for the indigenous Hawaiian people.

##  Appendix A Estimating Static Image Posteriors with Themis

###  A.1 The Themis Image Model

Themis is a general sampling-based parameter estimation framework developed
for comparing parameterized models with the VLBI data produced by the EHT
(Broderick et al., 2020a). Implemented within Themis is a wide variety of
image models, ranging from phenomenological geometric models (e.g., Gaussians,
rings, etc.) to physically motivated ray-traced radiative transfer models
(e.g., semi-analytic accretion flow models). Of relevance here is a set of
splined raster models, described in detail in Broderick et al. (2020b), and
applied to polarized reconstructions in M87∗ Paper VII. In these, the image
reconstruction process is replaced with a model reconstruction problem in
which the intensity at control points on a rectilinear raster grid are varied,
with the final image produced at arbitrary locations via a cubic spline.

The image model presented in Broderick et al. (2020b) has a fixed field of
view in the two cardinal directions. As described in M87∗ Paper VII, this
restriction has been subsequently relaxed, with a more general set of models,
adaptive splined rasters, that permit a rotation of the raster and a resizing
of the field of view among its principle axes. In this way, not only can the
brightness at each control point vary, but the locations of the control points
themselves can evolve. This results in a very flexible image model even when
the dimension of the raster is small. The splined raster and adaptive splined
raster models may be combined with any other Themis models, and typically
large-scale features will be absorbed by a large-scale Gaussian component. The
dimension of the underlying image model is, then, the sum of the number of
control points, $N_{x}\times N_{y}$, two fields of view and the orientation of
the raster, and the number of parameters associated with a potential geometric
addition.

The Themis model has been fit to a variety of EHT data types. The bulk of the
Themis-based analyses presented here employed the light-curve-normalized, LMT-
calibrated complex visibilities directly. The residual complex station gains
are reconstructed during each model evaluation and marginalized over via the
Laplace approximation (see Section 6.8 of Broderick et al., 2020a). Log-normal
priors are imposed on the station gain amplitudes with uncertainties that
depend on the a priori calibration estimates (see Section 2.2). Network
calibrated sites (ALMA, APEX, JCMT, SMA) are assumed to have a log-normal
1$\sigma$ uncertainty of 1%; following the LMT gain amplitude correction it is
assumed to have a residual log-normal 1$\sigma$ uncertainty of 20%; and the
remaining sites (PV, SMT, SPT) are assigned a log-normal 1$\sigma$ uncertainty
of 10%. In this way Themis fully explores the potential station gain space
during its exploration of the underlying image space. Some early analyses made
use of visibility amplitudes and closure phases. For these the station gain
amplitudes were reconstructed in a similar fashion as described above, though
with a larger log-normal prior on LMT gain of 100% and 20% on all other
stations.

###  A.2 Scattering and Variability Mitigation

Application of Themis-based imaging methods to Sgr A∗ is complicated by the
presence of the confounding factors described in Section 3. Unless otherwise
stated, diffractive scattering is mitigated via the scattering model from
Johnson et al. (2018) directly to the model visibilities, i.e., instead of
deblurring the data, the model is blurred. Refractive scattering and
variability are mitigated via a scheme similar to that described in Section
3.1.2. Unlike the description there, for the Themis analyses, the additional
noise terms are treated as parameterized contributions to the visibility
uncertainties as described in detail in Georgiev et al. (2021).

Explicitly, the Themis analyses assume that the visibility uncertainties are
described by

$\sigma_{j}^{2}(\sigma_{\rm ref},f,a,b,c,u_{0})=\\\ \sigma_{j,\rm
th}^{2}+\sigma_{\rm ref}^{2}+f^{2}|V|^{2}+\sigma_{\rm var}^{2}(a,b,c,u_{0}),$
(A1)

where $\sigma_{j,\rm th}$ is the intrinsic thermal uncertainty, $|V|$ is the
measured visibility amplitude, and $\sigma_{\rm var}$ is given in Equation 2.
Within this prescription, the $\sigma_{\rm ref}$ is conceptually identical to
the Const approach to refractive scattering mitigation described in Section
3.1.2, with the exception that it is not fixed. Similarly, the fractional
contribution, $f|V|$ is conceptually identical to the systematic noise term,
again with the sole distinction being to treat $f$ as a free parameter.

###  A.3 Likelihood and Parameter Priors

For analyses that fit the visibility amplitudes and closure phases, the Themis
log-likelihood is given by the appropriate combination of log-likelihoods from
Sections 6.3 and 6.8 from Broderick et al. (2020a). For analyses that fit the
full complex visiblities with the noise modeling, the log-likelihood is

$\mathcal{L}({\bf p},{\bf q})=\\\ -\sum_{\rm
data}\bigg{\\{}\frac{|\hat{V}_{j}-\hat{V}({\bf u}_{j};{\bf
p})|^{2}}{2\sigma_{j}^{2}({\bf
q})}+\frac{1}{2}\ln\left[2\pi\sigma_{j}^{2}({\bf q})\right]\bigg{\\}},$ (A2)

where $\\{\hat{V}_{j}\\}$ is the light-curve-normalized visibility data, ${\bf
p}$ are the image-model parameters and ${\bf q}=\\{\sigma_{\rm
ref},f,a,b,c,u_{0}\\}$ are the noise-model parameters. The novel terms on the
right are simply the $\sigma$-dependent normalization; while normally
constant, in this instance they induce a penalty term that constrains the
${\bf q}$. When descattered, the diffractive scattering kernel is applied to
the model, as opposed to deblurring the data, which maintains a consistent
treatment of the “noise” model among the descattered and on-sky analyses.

Priors on each of the image model parameters are listed in Table 6. In
summary, uninformative logarithmic priors are assumed on the brightness at the
control points, the raster size and orientation, and a potential shift from
the correlation phase center. Similarly uninformative priors are imposed on
the refractive scattering mitigation and systematic uncertainty terms of the
noise model. In contrast, informed priors are imposed on the remaining
parameters of the noise model, set by the pre-modeling considerations from
Georgiev et al. (2021) and applied to Sgr A∗ in Paper IV. The relevant
interquartile ranges are data-set specific, and are listed for the Themis
analyses reported here in Table 2 for Sgr A∗ and the various synthetic data
sets in Section 5.

###  A.4 Reconstructing the Posterior with Themis

Themis provides a number of sampling schemes with which to explore the
likelihood surface. The bulk of the Themis-based analyses presented here
employ a set of Markov Chain Monte Carlo (MCMC) samplers to explore the large-
dimensional parameter space, $\\{{\bf p},{\bf q}\\}$. To efficiently traverse
the prior, identify multiple modes, and rapidly locate the global maximum,
Themis makes use of parallel tempering, in which multiple, tempered copies of
the likelihood, $\mathcal{L}_{\beta}=\beta\mathcal{L}$ are simultaneously
explored; higher “temperature” copies ($\beta\rightarrow 0$) become the
reference distribution which we take to be uniform in the prior bounds, while
the bottom “temperature” ($\beta=1$) is the posterior to be sampled. These are
coupled via the deterministic even–odd swap tempering scheme (Syed et al.,
2019), which efficiently moves information through the multiple, tempered
levels. Periodically, the $\beta$ are optimized to efficiently move
information from the prior to the posterior following the scheme from (Syed et
al., 2019). Each individual tempering level is sampled via the Hamiltonian
Monte Carlo sampling kernel from the Stan package (Carpenter et al., 2017).

Convergence is assessed via chain statistics, such as, integrated
autocorrelation time, approximate split-$\hat{R}$ tests, parameter rank
distributions, and visual inspection. Reconstruction quality is determined by
inspection of the best-fit sample’s residuals, and “reasonableness” of the
associated image and gain reconstructions. Pathological solutions are
uncommon, and are typically indicative of poor sampler initialization and/or
errors in the input. The number of tempering levels is chosen to exceed the
global communication barrier $\Lambda$ (Syed et al., 2019), typically near 20,
and is run for $5\times 10^{4}$ to $10^{5}$ MCMC steps per tempering level.

##  Appendix B Refractive Noise Levels Anticipated for EHT Data

(a) On-sky images

(b) Refractive noise levels

Figure 29: Dependence of refractive substructures on the power-law index
$\alpha$ and inner scale of the turbulence $r_{\mathrm{in}}$. (a) Example
realizations of the on-sky (i.e., scattered) images of the crescent model in
Section 5. Each row shows a realization for a different power-law index
$\alpha$ (1.3, 1.4 and 1.5), while each column shows a different $r_{in}$
(600, 800 and 1200 km). (b) The corresponding levels of the flux-normalized
refractive noise anticipated for $(u,v)$-coverage of April 7 Sgr A∗ data. The
corresponding normalized visibility amplitudes are shown in Figure 3.

To assess the effects of refractive scattering on our EHT measurements, we use
the stochastic optics module implemented in eht-imaging library (Johnson,
2016; Chael et al., 2016). The stochastic optics allows the simulation of
scattering realizations for a given intrinsic image and generation of the
corresponding synthetic data sets. Here, we focus only on the refractive noise
that would affect our imaging, i.e., the substructures on the scattered image.
The other effects — centroid shift due to the position wander or total flux
modulation— do not affect the image reconstruction of data calibrated by self-
calibration and obtained from a single on-sky realization of scattering
(Johnson et al., 2018).

The refractive noise level at each $(u,v)$-coordinate is determined by the
standard deviation of the flux-normalized visibilities for all realizations,
where the centroid of each realization is shifted to the phase center to
mitigate the effects of position wander. We adopt a dipole model for the
magnetic field wandering in the scattering screen (Johnson et al., 2018). Note
that the refractive noise does not depend strongly on the choice of the model
for field wander (Psaltis et al., 2018) for the source size anticipated from
EHT data and scattering parameters constrained by Johnson et al. (2018).

We first compare the refractive effects for different combinations of the
power-law spectral index $\alpha$ of the density fluctuation of the ionized
plasma (see Section 3.1.1) and the inner scale of turbulence $r_{in}$. Our
refractive noise prescriptions J18model1 and J18model2 are derived with
$\alpha=1.38$ and $r_{in}=800\,{\rm km}$, which are the best-fit values in
Johnson et al. (2018). These measurements have the uncertainties of
$1.3<\alpha<1.5$ and $600\,{\rm km}<r_{in}<1500\,{\rm km}$ (Johnson et al.,
2018). Figure 29 shows an example gallery of scattered images and
corresponding noise levels in the visibility domain for each combination of
$\alpha$ and $r_{in}$ in these ranges. The level of refractive substructure as
well as the kernel size of diffractive angular broadening increase with
$\alpha$ and $r_{in}$ (Psaltis et al., 2018). Based on the results shown in
Figure 29, we adopt a scaling factor of up to 2 to cover the cases of highest
refractive noise level in the possible ranges of these two parameters.

Figure 30: The levels of the flux-normalized refractive noises as a function
of $(u,v)$-distance for different intrinsic structures. The refractive noise
level is computed on $(u,v)$-coverage of April 7 Sgr A∗ data. The
corresponding normalized visibility amplitudes are shown in Figure 3.

Refractive noise is expected to also be dependent on the intrinsic source
structure (e.g. Johnson & Narayan, 2016). To assess this dependence, we
estimate the refractive noise levels based on a circular Gaussian and the
time-averaged images of the seven geometric models introduced in Section 5.1,
for $\alpha=1.38$ and $r_{in}=800\,{\rm km}$. As described in Section 5.1,
these geometric models provide a broad consistency with the visibility
amplitudes profile of Sgr A∗ data. The size of the circular Gaussian is
adjusted to have the FWHM of 45 $\mathrm{\mu as}$, consistent with the
effective FWHMs of other seven geometric models estimated from their second
moments. In Figure 30, we show the refractive noise levels for the different
intrinsic source models. The refractive noise level is only a few percents of
the total flux density regardless of the intrinsic source structures. The
standard deviations across different geometric models are typically
$\lesssim$0.4% of the total flux density.

J18model2 was derived by taking the mean of the refractive noise levels of
seven geometric models. J18model1, on the other hand, is based on the
refractive noise level of the the circular Gaussian model. However, the
circular Gaussian model has the lowest refractive noise level in the most of
$(u,v)$-coverage of Sgr A∗ data, indicating that the refractive noise level
from the circular Gaussian model, if not scaled properly, may underestimate
that of Sgr A∗ data. We derived a scaling factor to make the median of the
$\chi^{2}$ of all seven geometric models unity to minimize the effects of its
systematically low refractive noise level.

##  Appendix C First Sgr A∗ Images from Blind Imaging

Independently performed comparison between synthetic and real observational
data can help to identify which image features are likely intrinsic and which
are most likely spurious (see this approach in Bouman (2017) for synthetic
data and M87∗ Paper IV for observed EHT data of M87∗).

To initiate the exploration of possible Sgr A∗ images while minimizing
influence from collective bias, in our first stage of analysis we
reconstructed images of Sgr A∗ in five independent teams. The teams were blind
to each others’ work and prohibited from discussing their imaging results.
This procedure was similar to that done in M87∗ Paper IV for first images of
M87∗, with two notable differences. First, since the participants in all teams
were already aware of problematic data from previous EHT imaging work, the
team activities were not blind to previously-identified problems in the data
that could be affecting results. Second, as teams approached the date to
submit results it was made known to the organizers that teams would not be
comfortable submitting only a single image or movie, and the instructions were
changed to submit three representative images or movies per team.

No restrictions were imposed on the data pre-processing or imaging procedures
used by each team. Teams 1 and 2 focused primarily on RML methods. Teams 3 and
4 primarily on CLEAN methods. Team 5 used a Bayesian method. All teams used an
early-release engineering data set. While the reduction procedures and
metadata corresponded to a preliminary version of the calibration pipeline,
the data products was deemed mature enough for the blind imaging consistency
test. The April 7 data set was selected for the first comparison, as it had
the best $(u,v)$-coverage and a very stable light curve over the full
observation (Wielgus et al., 2022).

Figure 31: The first EHT images of Sgr A∗, blindly reconstructed by five
independent imaging teams using an early, engineering release of Stokes I data
from the April 7 observations. Images from Teams 1 and 2 used RML methods (no
restoring beam); images from Teams 3 and 4 used CLEAN (restored with a
circular 20 $\mathrm{\mu as}$ beam, shown in the lower right); images from
Team 5 used Themis. Each Team shows three representative images from their
reconstructions. Some images with the label “Movie” are time-averaged images
of movie reconstructions. Many images show ring-like morphology with the
diameter of $\sim 50$ $\mathrm{\mu as}$, while some reconstructions show non
ring-like morphology. The significant differences in brightness temperature
between images are caused primarily by different assumptions regarding the
total compact flux density (see Table 8) and also restoring beams applied only
to CLEAN and Themis images. Table 8: Image Properties and Data Consistency
Metrics for the First Sgr A* Images

|  | Team 1 | Team 2 | Team 3 | Team 4 | Team 5
---|---|---|---|---|---|---
Image Properties |
Method | Image ID | RML | RML | CLEAN | CLEAN | THEMIS
$F_{\rm cpct}$ (Jy) | 1 | 2.37 | 2.07 | 1.90 | 1.64 | 1.59
| 2 | 2.38 | 2.24 | 1.29 | 1.18 | 1.61
| 3 | 2.38 | 1.99 | 1.25 | 0.69 | 1.73
Engineering Data (10-sec avg., stokes I, 0% sys. error) |
$\chi^{2}_{\text{CP}}$ | 1 | 2.27 (1.85) | 2.35 (1.75) | 4.52 (4.20) | 2.21 | 7.78
| 2 | 2.30 (1.81) | 2.08 | 2.51 | 3.31 | 6.79
| 3 | 2.81 (1.80) | 2.08 | 2.30 | 13.80 | 13.58
$\chi^{2}_{\text{log\,CA}}$ | 1 | 3.44 (1.77) | 4.67 (2.19) | 2.76 (2.98) | 2.06 | 5.90
| 2 | 2.44 (2.07) | 2.00 | 2.49 | 2.84 | 4.88
| 3 | 4.73 (2.08) | 6.44 | 2.93 | 1.91 | 4.95
Science Release (60-sec avg., Stokes $I$, 0% sys. error) |
$\chi^{2}_{\text{CP}}$ | 1 | 5.55 (2.53) | 6.09 (1.96) | 25.05 (21.77) | 4.58 | 54.11
| 2 | 6.08 (2.32) | 3.89 | 7.27 | 12.11 | 45.09
| 3 | 8.78 (2.31) | 4.11 | 5.73 | 108.88 | 109.10
$\chi^{2}_{\text{log\,CA}}$ | 1 | 7.12 (1.52) | 10.75 (1.69) | 8.31 (5.92) | 2.79 | 24.22
| 2 | 4.02 (1.46) | 2.81 | 5.43 | 9.66 | 19.94
| 3 | 6.47 (1.79) | 16.95 | 8.21 | 7.45 | 10.01
Science Release (60-sec avg., Stokes $I$, 10% sys. error) |
$\chi^{2}_{\text{CP}}$ | 1 | 2.42 (1.33) | 2.62 (1.18) | 7.80 (6.92) | 2.32 | 10.24
| 2 | 2.41 (1.25) | 2.09 | 2.85 | 5.93 | 8.02
| 3 | 4.24 (1.17) | 2.07 | 2.56 | 13.59 | 16.69
$\chi^{2}_{\text{log\,CA}}$ | 1 | 3.36 (0.74) | 4.73 (0.89) | 3.98 (2.86) | 1.32 | 10.70
| 2 | 1.89 (0.76) | 1.26 | 2.49 | 3.88 | 7.55
| 3 | 3.57 (0.88) | 6.88 | 3.72 | 3.59 | 4.66

† $F_{\rm cpct}$, $\chi^{2}_{\text{CP}}$, and $\chi^{2}_{\text{log\,CA}}$ mean
the compact total flux density of the VLBI scale images, $\chi^{2}$ to closure
phases and log closure amplitudes, respectively. For movies, the values in the
parentheses show $\chi^{2}$ to the original movies, while the values outside
of them show $\chi^{2}$ to the time-averaged images.

Figure 31 shows images of the three representative submissions made by each
team for Sgr A∗. In cases where dynamic movie reconstructions of Sgr A∗ were
submitted, we show the time-averaged image and indicate that the original
submission was a movie. Although most submissions contain flux with a $\sim$50
$\mathrm{\mu as}$ separation, not all submissions contain a clear ring
feature. Additionally, of the submissions that do contain a ring feature, the
azimuthal flux distribution varies drastically. For example, Team 1 submitted
three movies: one containing a ring with a position angle of $\sim
30^{\circ}$, another containing a ring with position angle of
$\sim-70^{\circ}$, and the third containing no ring at all.

This initial imaging stage indicated that a $\sim$50 $\mathrm{\mu as}$ feature
was likely to exist in the image, but there is significant uncertainty in the
presence of a ring structure, as well as the flux distribution around a
possible ring structure. In Table 8, we show $\chi^{2}$ values to closure
phases ($\chi^{2}_{\rm CP}$) and log closure amplitudes ($\chi^{2}_{\rm\log
CA}$) for both engineering data used for imaging and science release data
described in Section 2, as well as the compact total flux densities ($F_{\rm
cpct}$). There are no significant differences seen in $\chi^{2}$ values of
images regardless of the presence of a ring structure and between different
azimuthal intensity distributions of the ring structure. These ambiguous
results motivated a series of systematic tests to identify if source evolution
(similar to that expected to be present in Sgr A∗) could explain the multiple
solutions recovered.

##  Appendix D Imaging Pipelines

###  D.1 DIFMAP Pipeline

For this second stage of the imaging process we developed a scripted version
of CLEAN in DIFMAP (version higher than 2.5k), together with a python wrapper
used over the script, for carrying out an extensive parameter search. The
initial manual CLEAN analysis (see Appendix C) together with the pre-imaging
considerations (see Section 3) provided the set and range of parameters to be
explored on both the synthetic data sets (described in Section 5) and the
actual April 6 and 7 Sgr A∗ data, processed with both fringe fitting pipelines
HOPS and CASA, with high and low bands combined.

Prior to CLEAN imaging, the EHT data was pre-processed using the pre-imaging
pipeline discussed previously. This includes intra-site normalization to the
Sgr A∗ light curve, coherent averaging of the visibilities using 10 and 60
seconds, the inclusion of an additional systematic error by a factor of 0, 2,
and 5%, and the scattering and noise mitigation prescriptions with the common
parameters listed in Section 6.1 (see also, Table 3).

The first step in the DIFMAP imaging process is the gridding of the
visibilities in the $(u,v)$-plane, that was done by using two different
approaches: uniform and natural weighting. The next step is to down-weight the
data on all baselines to phased-ALMA in the self-calibration process. As found
for the case of M87 (M87∗ Paper IV), this prevents the ALMA baselines to
dominate the phase and amplitude self-calibration process due to their
significantly higher S/N. We have tested scaling factors for ALMA baselines
weights during self-calibration of 0.1 and 0.5.

Atmospheric fluctuations at the short wavelength of 1.3 mm severely limit the
S/N and the capability to reliably measure the visibility phases. We therefore
rely on the closure phase measurements and the reconstruction of station
phases using the Cornwell Wilkinson hybrid phase self-calibration approach.
However the self-calibration with CLEAN model solely relies on the visibility
phases and amplitudes, which forces an initial self-calibration of the phases
(visibility amplitudes are much better constrained; Paper II), a process which
is anchored on the more robustly measured closure phases.

To mitigate any possible bias in our choice for the initial model to self-
calibrate the phases, and therefore to minimize the chances of reaching a
local minimum in our manifold of image models, we have explored different
initial models for the first phase self-calibration. The “fiducial” initial
model has been chosen, based on the reduced-$\chi^{2}$ of closure phases with
the first reconstructed CLEAN model after the phase self-calibration with the
initial model among three different types: a Gaussian with a FWHM of 15
$\mu$as (i.e., unresolved symmetric model), a uniform disk with a size ranging
between 56 and 84 $\mu$as (in steps of 4 $\mu$as), and a uniform ring with
sizes ranging between 36 and 68 $\mu$as (also in steps of 4 $\mu$as, no
width).

A distinctive feature of CLEAN imaging is the use of masks (also known as
cleaning windows) to define the area of the image to be cleaned. For this
purpose we used a set of centered circular-disk-shaped CLEAN windows (i.e.,
with no hole in the middle to avoid any biasing).

The size of the cleaning windows should match the expected emission extend of
the imaged source. To avoid introducing any prior bias we have surveyed for a
relatively large range of mask sizes, from 80 to 110 $\mu$as, that covers our
prior Sgr A∗ image size constrains (Paper II). We note that due to the limited
$(u,v)$-coverage larger masks ($\geq$ 100 $\mu$as) may pick up some emission
from the main side-lobe of the interferometer, which may be particularly
problematic for the point source model given its diffuse emission, more
difficult to recover with the CLEAN algorithm.

Contrary to the case of M87 (M87∗ Paper IV), we do not expect a large extended
missing flux in Sgr A∗ (Paper II), hence our CLEAN stopping criterion is based
on setting a minimum threshold for the relative decrease in the rms of the
residual image over the image noise estimated from the visibilities.

As a result, the total number of surveyed parameters has been 1680 for the on-
sky data set, and 6720 for descattered data. Note that the descattered results
are based on a 4 times larger number of parameters due to the scattering
mitigation (i.e., refractive noise floor).

###  D.2 eht-imaging Pipeline

The eht-imaging software library implements the RML imaging technique
described in Section 4.2. As a RML method, a successful image reconstruction
is achieved when the proper hyperparameters balance the contribution of the
selected data products and regularizers to the minimization of Section 4.4.1.
Thus, we designed a framework around our imaging pipeline to explore an
extensive range of imaging and data pre-calibration parameters.

The pipeline pre-process the input data following the procedure described in
Section 6.1 as a first step. Only low- and high-band light-curve normalized
data sets were used for static imaging to minimize the effects of the source
variability. We opted for an integration time of 60 s and explored all the
non-closing systematic error, scattering and intra-day variability noise
budgets tabulated in Table 4. The algorithm is then initialized using a
symmetric flat disk image contained in a square-shaped grid of 80$\times$80
pixels, which corresponds to a FOV of 150 $\mathrm{\mu as}$. The disk diameter
is surveyed as a parameter with values of 70, 80, and 90 $\mathrm{\mu as}$.
Note that, among the different imaging regularizers used in the pipeline, the
maximum entropy regularizer (MEM) rewards similarity to a prior image, which
in this case is assumed to be equal to the initialization image.

Full closure quantities and visibility amplitudes are used as data products.
Contrary to M87∗, the relative weighting of the visibility amplitudes adopted
range from a less to an equally important contribution as closure quantities,
since complex station gains were derived from calibrator sources beforehand
(see Section 5.1.3 in Paper II). Several imaging regularizers and their
relative weighting are employed to constraint on the image properties, such as
MEM, total variation (TV), and total squared variation (TSV), which enforce
similarity to a prior image and/or smoothness over neighboring pixels,
respectively.

Convergence to an optimal image which minimizes Section 4.4.1 is then carried
out in an iterative process in which the reconstructed image is blurred using
a Gaussian kernel with the FWHM of the array nominal resolution ($\sim$24
$\mathrm{\mu as}$), and then used as initialization for the next one to
prevent the algorithm of being caught in local minima. In contrast to the M87∗
pipeline, no self-calibration was performed to the data, since the $\chi^{2}$
statistics obtained for the output images of just one cycle of this iterative
procedure were already good enough.

###  D.3 SMILI Pipeline

The SMILI imaging pipeline was designed to reconstruct images with RML imaging
techniques utilizing three imaging regularizers: weighted-$\ell_{1}$
($\ell_{1}^{w}$), TV, and TSV, similar to M87∗ Paper IV. Before the imaging
process, both low- and high-band EHT data are pre-processed as described in
Section 6.1. First, data are normalized by the time-dependent intra-site flux
density, and coherently time-averaged at an integration time of 120 sec. The
intra-scan fluctuations in LMT baselines are corrected by self-calibrating the
shortest VLBI baseline between SMT and LMT to a circular Gaussian with the
total flux density of 1 Jy and the FWHM size of 60 $\mathrm{\mu as}$. Finally,
the pipeline applies the parameterized mitigation schemes for the non-closing
errors, scattering effects and intra-day variations as described in Section
6.1.

The pipeline reconstructs an image from pre-processed data at both low and
high bands jointly with the field of view of $150\ \mu\text{as}$ discretized
by 2 $\mathrm{\mu as}$ pixels. Throughout the imaging, the pipeline adopts a
prior image of a circular Gaussian with the FWHM size of 140, 160 or 180
$\mathrm{\mu as}$ for $\ell_{1}$ prior. The imaging process consists of in
total nine imaging cycles, including a single self-calibration of complex
visibilities. Each imaging cycle performs 5000 iterations of the L-BFGS-B
algorithm used for the gradient-descent optimization in the SMILI’s image
solver.

In the first eight cycles of imaging, the images are reconstructed from
visibility amplitudes, closure phases and log closure amplitudes. To account
for residual errors in the amplitude calibration, the fractional 5 % and 30 %
uncertainties are added in quadrature to the errors of visibility amplitudes
on non-LMT and LMT baselines, respectively. The first cycle starts with a
circular disk image with a $70\ \mu\text{as}$ diameter with the total flux of
1 Jy, while the later cycles subsequently adopt the image from the previous
cycle as the initial image after re-centering its center of the mass and
blurring with a $5\ \mu\text{as}$ circular Gaussian.

After the eight cycles of imaging, complex visibilities are self-calibrated
with its output image. Since the budgets of scattering and intra-day
variations are not fractional to visibility amplitudes, these two budgets will
not be properly scaled with gains solved with a self-calibration. Therefore,
to reflect these two budgets accurately on self-calibrated data, the SMILI
pipeline updates these two budgets added to complex visibilities and closure
quantities after self-calibration. Finally, the final image is reconstructed
by a single imaging cycle using visibility amplitudes and closure quantities.

###  D.4 Themis Raster Dimension Survey

The sole hyperparameters in the Themis analyses not self-consistently explored
during the construction of the posterior are the raster dimensions, $N_{x}$
and $N_{y}$. While reasonable initial guesses may be made based on the
diffraction limit, in practice we employ data-driven method: comparing the
Bayesian evidence (or appropriate proxies) of different dimensions (Broderick
et al., 2020b).

Due to the computational expense of performing Themis posterior
reconstructions, we limit the number of raster surveys performed to three: for
the Sgr A∗ April 6 data, the Sgr A∗ April 7 data, and the GRMHD synthetic data
set presented in Section 6.4.2. Of these, we present only that on the Sgr A∗
April 7 data here for brevity; surveys on the other data sets were conducted
similarly. Furthermore, while we did experiment with $N_{x}\neq N_{y}$, given
the apparent symmetry of the Sgr A∗ image, and not wishing to introduce any
potential biases away from symmetry in the model specification, we restrict
ourselves here to square rasters, i.e., $N_{x}=N_{y}$.

Figure 32: Mean images from the Themis posteriors for $(N_{x},N_{y})=(4,4)$,
$(5,5)$, $(6,6)$, $(7,7)$ and $(8,8)$ from left to right. In each the
difference in the logarithm of the Bayesian evidence, $\Delta\log Z$, relative
to the $(7,7)$ model is listed in the lower left.

Five individual analyses were performed on the April 7, combined high- and
low-band complex visibility data151515A similar study was performed using only
the low-band data with similar results. as described in Appendix A, differing
in the raster dimensions: $(N_{x},N_{y})=(4,4)$, $(5,5)$, $(6,6)$, $(7,7)$ and
$(8,8)$. For each, after convergence the Bayesian evidence, $Z$, was computed
via thermodynamic integration across tempering levels (Lartillot & Philippe,
2006). Mean images from each analysis and the relative Bayesian evidence are
shown in Figure 32.

For raster resolutions that are too small, the fit quality is poor. For raster
resolutions that are too large, the added model complexity is not justified by
the fit improvement. Importantly, by selecting on $\log Z$, we avoid potential
complications associated with the modifications of the likelihood due to noise
modeling, correlations between the high- and low-band gains, and non-
Gaussianity of the posterior.

As roughly anticipated by the diffraction limit, the preferred raster size is
$(7,7)$, which we adopt henceforth. This is accompanied by a convergence of
the image structure for raster dimensions of $(5,5)$ and larger. Smaller and
larger raster dimensions are overwhelmingly disfavored. Nevertheless, by
$(4,4)$ the locations of the knots are identified, though the model
misspecification prevents a faithful recovery of their relative brightness.

An identical procedure applied to the Sgr A∗ April 6 data finds that a smaller
$(6,6)$ raster is preferred, ostensibly due to the smaller data volume.
Similarly, when applied to the GRMHD synthetic data set in Section 6.4.2, a
raster dimension of $(6,6)$ is again preferred, possibly due to the
competition between the simplicity of the ring structure.

##  Appendix E Top Set Selection

###  E.1 The Effect of Relaxation in $\rho_{NX}$ Criteria

As discussed in Section 6.2, we relaxed the normalized cross-correlation
criteria for each image by multiplying the value of $\rho_{NX}$ obtained for
$\alpha=24$ $\mathrm{\mu as}$ beam by a relaxation factor of 0.95. This
relaxation factor allows us to recover a large enough number of Top Set
parameters for evolving synthetic data with Sgr A∗ $(u,v)$-coverage. This
contrasts with the Top Set selection for M87 (M87∗ Paper IV), where no
relaxation factor was considered. The necessity of the relaxation factor
originates from the fact that temporal variations of the source structure
cause a loss in the image fidelity of the reconstructions.

(a) Descattered reconstructions

(b) On-sky reconstructions

Figure 33: The worst (i.e., lowest) $\rho_{\rm NX}$ images of synthetic
geometric models among those reconstructed from the Top Sets of imaging
parameters (a) with and (b) without scattering mitigation.

To ensure that this relaxed criteria is still able to distinguish different
morphologies, we show the Top Set images with the worst (i.e. lowest)
$\rho_{NX}$ for each geometric model in Figure 33. Each descattered or on-sky
reconstruction in Figure 33 shows the image that has the worst $\rho_{NX}$
among the Top Set images for April 7. These images still retain resemblance to
the ground truth morphology, demonstrating that Top Set selection with a
relaxation factor of 0.95 is still capable of reconstructing the basic
structures of synthetic images with varied morphologies.

###  E.2 $\chi^{2}$ Distribution of Sgr A∗ Images

Figure 34: $\chi^{2}$ distributions of the images from each pipeline and epoch
which passed the $\rho_{\rm NX}$ criteria. $\chi^{2}$ of closure phases and
log closure amplitudes are showed in blue and orange respectively. We assumed
$1\ \rm{\%}$ systematic noise and representative variability noise for
scattered images and we also added refractive noise when calculating
$\chi^{2}$ of deblurred images.

Although the Top Set selection described in Section 6.4.1 is sorely based on
the image fidelity of the reconstructions from the synthetic data sets, Top
Set images of Sgr A∗ provide reasonable fits to EHT measurements. In Figure
34, we show $\chi^{2}$ distributions of all the Sgr A∗ images passing the
$\rho_{NX}$ criterion using synthetic data reconstructions. $\chi^{2}$ values
were computed for closure phases and log closure amplitudes formed from data
after averaged at 1 min, based on the same definition with M87∗ Paper IV. To
account for non-closing effects, refractive scattering and time-variability,
we include 1 % of the fractional errors and the representative budgets for
scattering and temporal variability (see Section 3.1 and Section 3.2),
respectively – in particular, the J18model1 refractive noise model and a
variability model with parameters $a=0.02$, $u_{0}=2$, $b=2.5$, and $c=2$ are
used. As shown in Figure 34, most of the images have $\chi^{2}$ less than
unity, and all the images provide $\chi^{2}<2$, demonstrating that all Top Set
images provide reasonable fits to the Sgr A∗ data within the anticipated
deviations of time-variable on-sky images from its intrinsic time-averaged
structure. Therefore, we adopt all the parameter sets satisfying the
$\rho_{NX}$ criteria as the Top Sets without any cut by $\chi^{2}$.

###  E.3 Top Sets of Imaging Parameters for April 6

Table 9: Parameters in the DIFMAP Pipeline Top Set on April 6

April 6 (8400 Param. Combinations; 365 in Top Set)
---
Systematic | 0 | 0.02 | 0.05
error | 15.6% | 17.3% | 67.1%
$\bf Ref\ Type$ | No | Const | 2$\times$Const | J18 | 2$\times$J18
| 18.1% | 23.0% | 14.2% | 24.1% | 20.5%
$\bf a_{\bf\rm psd}$ | No | 0.015 | 0.02 | 0.025
| 21.9% | 26.8% | 25.8% | 25.5%
$\bf b_{\bf\rm psd}$ | No | 1 | 3 | 5
| 21.9% | 17.3% | 24.4% | 36.4%
$|u|_{0}$ | No | 2
| 21.9% | 78.1%
Time average | 10 | 60
(sec) | 49.9% | 50.1%
ALMA weight | 0.1 | 0.5
| 39.5% | 60.5%
UV weight | 0 | 2
| 72.6% | 27.4%
Mask Diameter | 80 | 85 | 90 | 95 | 100 | 105 | 110
($\mathrm{\mu as}$) | 10.7% | 19.2% | 36.2% | 17.0% | 6.0% | 7.9% | 3.0%

Note.  In each row, the upper line with bold text shows the surveyed parameter
value corresponding to the parameter of left column, while the lower line
shows the number fraction of each value in Top Set. The total number of
surveyed parameter combinations and Top Set are shown in the first row.

Table 10:  Parameters in the eht-imaging Pipeline Top Set on April 6

April 6 (112320 Param. Combinations; 1415 in Top Set)
---
Systematic | 0 | 0.02 | 0.05
error | 19.0% | 39.7% | 41.3%
$\bf Ref\ Type$ | No | Const | 2$\times$Const | J18 | 2$\times$J18
| 7.8% | 29.2% | 26.4% | 17.0% | 19.6%
$\bf a_{\bf\rm psd}$ | No | 0.015 | 0.02 | 0.025
| 13.0% | 34.2% | 29.4% | 23.4%
$\bf b_{\bf\rm psd}$ | No | 1 | 2 | 3 | 5
| 13.0% | 16.3% | 25.4% | 22.1% | 23.1%
$|u|_{0}$ | No | 2
| 13.0% | 87.0%
TV | $\bf 0$ | 0.01 | 0.1 | 1
| 12.4% | 27.8% | 51.8% | 8.0%
TSV | $\bf 0$ | 0.01 | 0.1 | 1
| 29.0% | 34.8% | 35.8% | 0.4%
Prior size | 70 | 80 | 90
($\mathrm{\mu as}$) | 15.3% | 27.1% | 57.6%
MEM | $\bf 0$ | 0.01 | 0.1 | 1
| 3.4% | 15.3% | 77.9% | 3.5%
Amplitude | 0 | 0.1 | 1
weight | 0% | 4.2% | 95.8%

Note.  Same as Table 9

Table 11: Parameters in the SMILI Pipeline Top Set on April 6

April 6 (54000 Param. Combinations; 292 in Top Set)
---
Systematic | 0 | 0.02 | 0.05
error | 34.2% | 31.8% | 33.9%
$\bf Ref\ Type$ | No | Const | 2$\times$Const | J18 | 2$\times$J18
| 7.9% | 14.0% | 18.8% | 24.3% | 34.9%
$\bf a_{\bf\rm psd}$ | No | 0.015 | 0.02 | 0.025
| 0.3% | 23.3% | 46.2% | 30.1%
$\bf b_{\bf\rm psd}$ | No | 1 | 2 | 3 | 5
| 0.3% | 28.8% | 37.0% | 28.1% | 5.8%
$|u|_{0}$ | No | 1 | 2
| 0.3% | 44.5% | 55.1%
TV | $\bf 10^{2}$ | $\bf 10^{3}$ | $\bf 10^{4}$ | $\bf 10^{5}$
| 7.2% | 92.5% | 0.3% | 0.0%
TSV | $\bf 10^{2}$ | $\bf 10^{3}$ | $\bf 10^{4}$ | $\bf 10^{5}$
| 41.4% | 58.6% | 0.0% | 0.0%
Prior size | 140 | 160 | 180
($\mathrm{\mu as}$) | 31.5% | 41.1% | 27.4%
$\bf\ell_{1}$ | 0.1 | 1 | $\bf 10$
| 99.7% | 0.3% | 0.0%

Note.  Same as Table 9

In Section 6.2.2, we show the results of the Top Set selection for April 7
data. Here, we show the summary of Top Set parameters of each pipeline for
April 6 data in Table 9, Table 10 and Table 11, which are selected based on
the method identical to that of April 7 data. Although the number of Top Set
parameters for April 6 exceeds 100, this number is significantly smaller than
the Top Set size identified for April 7 data. In particular, the number of Top
Set parameters is reduced to only $22\ \%$ (DIFMAP), $25\ \%$ (eht-imaging)
and $11\ \%$ (SMILI) of the April 7 Top Set size, as shown in Table 3, Table 4
and Table 5. This significant reduction is most likely due to poorer
$(u,v)$-coverage on April 6.

##  Appendix F Classification of Ring Images

As described in Section 5.2 and Section 7.2, we categorize GRMHD and Sgr A∗
images into clusters of images that share similar morphologies. Clustering of
images is performed separately for on-sky and descattered reconstructions from
each pipeline.

Prior to clustering, the images are aligned with an iterative method described
in below. We first derive the averaged image of all images. Each image is then
aligned to maximize the cross-correlation with the averaged image. After
aligning all images, the averaged image is recomputed with the aligned ones,
and used to align each image again. We repeat this procedure to align each
image for a total of three times, providing the convergence of the alignment.
After the above relative alignments between images, all images are centered
using the average image; the images are uniformly shifted with the same amount
of positional offset to maximize the cross-correlation between the average
image and the time-averaged image of Simple Disk model.

The clustering of images has two major steps: the identification of ring
images and non-ring images, and the clustering of the ring images by the peak
position angles. The ranges of the peak position angles are described in
Section 5.2 for the GRMHD model and Section 7.2 for the Sgr A∗
reconstructions. The separation of ring and non-ring images is based on two
morphological criteria as described in below.

First, we identify non-ring images by the degree of the central depression
seen in the image. We measure a typical brightness of the central region by
taking the median of the intensity within a radius of 10 $\mathrm{\mu as}$,
defined by

$I_{\rm
c}\equiv\operatorname*{median}_{r<10,\,\theta\in[-\pi,\pi)}I(r,\theta),$ (F1)

where $I(r,\theta)$ denotes the intensity at the radius of $r$ in $\mathrm{\mu
as}$ and the position angle of $\theta$ in radians. To measure the degree of
the central depression, we compare the measured central brightness $I_{\rm c}$
with the maximum of the azimuthally-averaged intensity at the outer area
within a radius of 10 to 40 $\mathrm{\mu as}$, given by

$I_{\rm
o}\equiv\max_{r\in[10,40]}\left(\frac{1}{2\pi}\int^{\pi}_{-\pi}I(r,\theta)\,d\theta\right).$
(F2)

We define the degree of the central depression by $f_{\rm cd}\equiv
1-I_{c}/I_{o}$, and classify as non-ring images those with $f_{\rm cd}<0.2$
and $f_{\rm cd}<0.15$ for descattered and on-sky reconstructions,
respectively. The slightly lower threshold for on-sky images takes account of
the angular broadening effects due to scattering, causing a systematic
decrease in the central depression of a ring emission. This criterion can
effectively distinguish point-like images from ring images.

The second criterion identifies non-ring images by the smoothness of the
azimuthal intensity distributions. To extract the azimuthal profile of the
ridge intensity for the outer emission, we take the maximum intensity of the
outer area within the radius of 10-40 $\mathrm{\mu as}$ for each position
angle, given by

$I_{p}(\theta)\equiv\max_{r\in[10,40]}I(r,\theta).$ (F3)

For each position angle $\theta$, we evaluate the difference between the ridge
brightness $I_{p}(\theta)$ to the central brightness $I_{c}$ normalized by the
maximum intensity of the outer area defined by

$f_{p}(\theta)\equiv\frac{I_{p}(\theta)-I_{c}}{\max_{\theta\in[-\pi,\pi)}I_{p}(\theta)}.$
(F4)

$f_{p}(\theta)$ allows to assess whether the azimuthal distribution of the
outer area has a dark gap area comparable or lower than the central
depression. We identify a continuous dark area with $f_{p}(\theta)<0.2$ and
$f_{p}(\theta)<0.15$ over a range of position angles broader than $70^{\circ}$
as a gap in the azimuthal intensity distribution of descattered and on-sky
reconstructions, respectively. If an image has more than two distinct gaps
separated by $>70^{\circ}$, we classify it as a non-ring image with two or
more distinct blobs not smoothly connected with each other. We also classify
an image as a non-ring one with an incomplete ring if it does not have a
continuous bright area without any gaps over the range of position angles less
than $180^{\circ}$ (i.e. not completing more than a half circle). The above
criteria effectively exclude multiple blob or highly corrupted ring-like
images.

The morphological criteria discussed above provides a classification of ring
and non-ring images that is broadly consistent with human perception. However,
as noted in Section 6.4.2, the classification of images that are borderline
between ring and non-ring classification are sensitive to the exact criteria
used. Therefore ring definitions that make use of slightly different criteria
can lead to classification that still largely aligns with human perception but
varies in the ring classification percentages quoted in this paper. In this
work, motivated by method interpretability, we chose to make use of the
simplest classification criteria that still largely aligned with human
perception.

##  Appendix G Validation of Static Imaging Results with Full-track Dynamic
Imaging

Table 12: Parameters in the SMILI dynamic Pipeline Top Set on April 7

April 7 (7776 Param. Combinations; 345 in Top Set)
---
Systematic | 0 | 0.02 | 0.05
error | 32.8% | 35.4% | 31.9%
$\bf Ref\ Type$ | No | J18model1
| 18.8% | 81.2%
TV | $\bf 10^{2}$ | $\bf 10^{3}$ | $\bf 10^{4}$ | $\bf 10^{5}$
| 9.0% | 22.0% | 69.0% | 0%
TSV | $\bf 10^{2}$ | $\bf 10^{3}$ | $\bf 10^{4}$ | $\bf 10^{5}$
| 32.2% | 44.6% | 23.2% | 0%
$\bf\ell_{1}$ | 0.1 | 1 | $\bf 10$
| 11.9% | 88.1% | 0%
Prior size | 140 | 160 | 180
($\mathrm{\mu as}$) | 29.6% | 36.5% | 33.9%
$\bf R_{t}$ | $\bf 10^{4}$ | $\bf 10^{5}$ | $\bf 10^{6}$
| 39.7% | 24.1% | 36.2%
$\bf R_{i}$ | $\bf 10^{4}$ | $\bf 10^{5}$ | $\bf 10^{6}$
| 14.2% | 40.9% | 44.9%

Note.  In each row, the upper line with bold text shows the surveyed parameter
value corresponding to the parameter of left column, while the lower line
shows the number fraction of each value in Top Set. The total number of
surveyed parameter combinations and Top Set are shown in the first row.

Figure 35:  The distribution of Sgr A∗ Top Set images on April 7 reconstructed
with the SMILI dynamic imaging pipeline. We show the distribution of images
for each cluster in the same convention with Figure 14 with three horizontal
panels separated by horizontal lines. The top panel shows individual images
randomly sampled from different clusters. The middle and bottom panels
visualize the distributions of reconstructed descattered and on-sky Sgr A∗
images for each cluster, respectively. In each panel, from top to bottom, we
show the average of each cluster, the distributions of the radial profiles,
and the distributions of azimuthal intensity profiles.

As described in Section 3.2, the temporal variation of Sgr A∗ is anticipated
to cause significant deviations of visibilities from those of the time-
averaged morphology, of which levels can be well characterized by a broken
power-law model. As described in Section 7.5.3, the variability noise model
allows us to enhance the overall fidelity of synthetic data reconstructions,
and enable more sets of the imaging parameter combinations to be selected as
Top Sets. However, the prescription for the temporal variability in Section
3.2 include some simplifications of its characteristic properties: for
instance, the circular symmetry assumed for the levels of the variability
amplitudes in Fourier space and no inclusion of the correlated variations
considered in data metrics as covariant components. Here, to assess the
dependence of the mitigation scheme for temporal variability, we show the
time-averaged images identified by full-track dynamic imaging not relying on
the variability noise model. We emphasize that the goal of this section is not
to characterize the dynamic evolution of Sgr A∗ on short timescales, which is
the primary focus of Section 9.

As shown in Figure 35, we find the primary three ring morphologies identified
in Section 7 with similar azimuthal variations in the most of Top Set images,
while non-ring structures are also reconstructed in a small fraction of Top
Sets. The results are broadly consistent with those from the imaging survey
presented in Section 7, and strongly indicate that our results described in
Section 7 are resilient to the methods to recover the time-averaged
morphology. We briefly describe the imaging process in Section G.1 and the
results in Section G.2.

###  G.1 SMILI Dynamic Imaging Pipeline and Top Set Selections

Similar to the RML and CLEAN parameter surveys described in Section 6.2, we
conducted a large imaging survey with RML dynamic imaging methods (see Section
4.4) implemented in a scripted pipeline using SMILI. The survey was performed
on all seven geometric models and Sgr A∗ data.

After completing the common pre-imaging process of data (Section 6.1), the
pipeline reconstructs the time-averaged images on $\sim 8000$ sets of the
imaging parameter combinations across a broad parameter space, as outlined in
Table 12. With the exactly same criteria described in Section 6.2, Top Sets of
the imaging parameters are then selected based on the fidelity of the
synthetic data reconstructions.

The SMILI dynamic imaging pipeline shares the same procedures with the SMILI
static imaging pipeline described in Section D.3, except for a major
difference. Instead of the variability noise model used in Section 6.2, we
allow temporal variations on time scales of typical scan intervals — for each
set of parameters, the pipeline reconstructs a movie with the frame interval
of 1 h, and then time-averages the movie to obtain the resultant
reconstruction of the time-averaged morphology. We utilize two temporal
regularizers, denoted by $R_{t}$ and $R_{i}$, enforcing the continuity of the
frame-to-frame intensity variations and the continuity between each frame and
time-averaged intensity distributions, respectively, based on the Euclidian
distance between images (see Johnson et al., 2017, for details). We note that
here we only explore J18model1 for the scattering mitigation, given
consistency among the different scattering mitigation schemes (Section 7.5.2).

###  G.2 Results

As shown in Figure 35, Top Set images from SMILI dynamic imaging pipeline can
be categorized into four clusters shown in Section 7.2. Similar to the RML,
CLEAN and Themis static imaging pipelines, the ring morphology with the
diameter of $\sim$50 $\mu$as was found as the dominant feature, while the non-
ring images were found in a small ($\sim 5$%) fraction of Top Sets. The
differences in the image appearance between on-sky and descattered
reconstructions are also consistent with Section 7.5.2; the on-sky
reconstructions are more blurry than the descattered ones. We note that
similar to Top Sets presented in Section 5 and Section 7, here Top Set images
from SMILI dynamic imaging pipeline do not constitute a likelihood and
therefore the fractions should not represent our degree of certainty.

##  Appendix H Reconstructions of a Best-bet GRMHD Model in Paper V

In Section 6, we utilize synthetic data based on a GRMHD model selected from a
library of time-dependent GRMHD models presented in Paper V (see Section 5.2)
to assess the performance of our imaging procedures on a physically motivated
evolving source. While this GRMHD model is broadly consistent with our
criteria in Section 5 based on 1.3 mm EHT data and light curves, it is not
identified as a “best-bet" model in Paper V; “best-bet" models satisfy
heterogeneous constraints derived from 1.3 mm EHT data, 86 GHz VLBI
observations with the GMVA, 2.2 ${\rm\mu m}$ flux density and x-ray
luminosity. Paper V identifies three GRMHD models that lay in a “best-bet”
region of strongly magnetized (MAD) models at low inclination with prograde
spin.

Figure 36:  Image and visibility characteristics of the best-bet GRMHD models
in Paper V (MAD, $a_{*}=0.5,i=30^{\circ},R_{\rm high}=160$). Left two panels
are the snapshot and averaged images of the ground-truth movie. Right panel
shows the simulated visibility amplitudes (red) and real Sgr A∗ measurements
(black) as a function of projected baseline length.  Figure 37:  The
distribution of reconstructed images with a best-bet GRMHD model identified in
paper V. The distribution of images from each pipeline for each cluster is
shown with the same convention as Figure 12.

Here we present example reconstructions of synthetic data derived from the
“best-bet" model shown in Figure 16 of Paper V, with positive spin of
$a_{*}=0.5$ and electron temperature of $R_{\rm high}=160$ viewed at
$i=30^{\circ}$. In Figure 36, we show a snapshot and time-averaged images of
this GRMHD model, as well as the amplitude versus $(u,v)$-distance
corresponding to synthetic data generated with April 7 coverage (generated in
the same manner with Section 5). As shown in Figure 36, visibility amplitudes
are broadly consistent with EHT Sgr A∗ data. Figure 37 shows descattered and
on-sky images, which are reconstructed with DIFMAP, eht-imaging, and SMILI
pipelines using their Top Set parameters and then further categorized by the
same clustering method presented in Appendix F.

As shown in Figure 37, the distributions of the Top Set images share key
results from those of the GRMHD reconstructions described in Section 6.4.2
(see also Figure 12). The vast majority of the reconstructions identify a ring
morphology with a diameter of $\sim$50 $\mu$as, consistent with the ground-
truth model. However, a small fraction of images have a non-ring morphology.
Furthermore, ring reconstructions have multiple azimuthal intensity modes. In
particular, ring images that have a peak PA of $\sim-154^{\circ}$ consistent
with the ground-truth image do not appear as the most popular modes in eht-
imaging and SMILI reconstructions. The broad consistency with the results in
Section 6.4.2 suggest that our main results (e.g., that for the RML and CLEAN
pipelines our methods produce a small fraction of non-ring modes for an
underlying ring model, and that the most popular mode reconstructed is not
always that true mode) likely generalize across GRMHD models that are in a
broad agreement with 1.3 mm EHT data.

##  Appendix I Ring fitting parameters

In Table 13 and Table 14 we list the diameter $d$, width $w$, position angle
$\eta$, asymmetry $A$, and fractional central brightness $f_{c}$ measured from
Top Set Sgr A∗ images for each identified cluster (see Section 8)
corresponding to the descattered and on-sky images, respectively.

Table 13: Mean and standard deviation of ring parameters, diameter $d$, width $w$, position angle $\eta$, asymmetry $A$, and fractional central brightness $f_{c}$ measured from Top Set or posterior descattered Sgr A∗ images for each cluster. |  | $d\ \rm{(\mu as)}$ | $w\ \rm{(\mu as)}$ | $\eta\ \rm{(^{\circ})}$ | $A$ | $f_{c}$
---|---|---|---|---|---|---
DIFMAP |  |  |  |  |  |
April 6 Ring | REx | $46\pm 4.1$ | $33\pm 3.5$ | $-100.6\pm 73.4$ | $0.15\pm 0.08$ | $0.47\pm 0.14$
| VIDA | $51\pm 3.1$ | $33\pm 3.1$ | $-108.3\pm 79.2$ | $0.26\pm 0.12$ | $0.46\pm 0.11$
April 7 Ring 1 |  | $50\pm 1.9$ | $31\pm 2.7$ | $-92.0\pm 42.0$ | $0.06\pm 0.04$ | $0.40\pm 0.08$
|  | $51\pm 1.3$ | $33\pm 1.4$ | $-82.7\pm 40.5$ | $0.07\pm 0.07$ | $0.37\pm 0.08$
April 7 Ring 2 |  | $49\pm 3.2$ | $32\pm 2.9$ | $-22.5\pm 54.4$ | $0.08\pm 0.06$ | $0.45\pm 0.12$
|  | $50\pm 1.7$ | $32\pm 1.4$ | $-27.0\pm 47.4$ | $0.12\pm 0.11$ | $0.42\pm 0.11$
April 7 Ring 3 |  | $49\pm 1.8$ | $32\pm 2.9$ | $20.5\pm 47.2$ | $0.08\pm 0.05$ | $0.44\pm 0.09$
|  | $51\pm 1.6$ | $33\pm 1.8$ | $5.3\pm 38.9$ | $0.12\pm 0.10$ | $0.41\pm 0.10$
eht-imaging |  |  |  |  |  |
April 6 Ring | REx | $56\pm 4.5$ | $24\pm 2.4$ | $32.3\pm 86.4$ | $0.17\pm 0.05$ | $0.16\pm 0.14$
| VIDA | $59\pm 11.3$ | $30\pm 10.4$ | $70.0\pm 88.3$ | $0.24\pm 0.14$ | $0.23\pm 0.19$
April 7 Ring 1 |  | $55\pm 1.9$ | $26\pm 2.3$ | $-140.6\pm 62.9$ | $0.11\pm 0.04$ | $0.18\pm 0.08$
|  | $56\pm 5.4$ | $27\pm 2.5$ | $-156.0\pm 60.6$ | $0.15\pm 0.13$ | $0.18\pm 0.09$
April 7 Ring 2 |  | $53\pm 1.8$ | $27\pm 2.7$ | $-60.1\pm 35.8$ | $0.12\pm 0.05$ | $0.23\pm 0.10$
|  | $53\pm 3.6$ | $27\pm 4.0$ | $-71.7\pm 33.1$ | $0.15\pm 0.11$ | $0.22\pm 0.10$
April 7 Ring 3 |  | $55\pm 1.5$ | $27\pm 3.0$ | $170.4\pm 101.4$ | $0.12\pm 0.07$ | $0.31\pm 0.17$
|  | $57\pm 6.4$ | $26\pm 4.5$ | $179.7\pm 67.1$ | $0.18\pm 0.21$ | $0.27\pm 0.11$
SMILI |  |  |  |  |  |
April 6 Ring | REx | $57\pm 3.4$ | $24\pm 1.9$ | $-27.7\pm 43.8$ | $0.23\pm 0.09$ | $0.03\pm 0.10$
| VIDA | $46\pm 12.0$ | $50\pm 16.6$ | $-71.6\pm 123.8$ | $0.18\pm 0.12$ | $0.59\pm 0.33$
April 7 Ring 1 |  | $52\pm 5.0$ | $26\pm 2.0$ | $151.9\pm 75.8$ | $0.12\pm 0.04$ | $0.22\pm 0.10$
|  | $52\pm 3.8$ | $27\pm 3.9$ | $175.9\pm 77.6$ | $0.10\pm 0.08$ | $0.22\pm 0.09$
April 7 Ring 2 |  | $53\pm 4.0$ | $25\pm 2.9$ | $-39.3\pm 51.2$ | $0.13\pm 0.06$ | $0.23\pm 0.11$
|  | $52\pm 6.2$ | $29\pm 8.4$ | $-59.2\pm 64.9$ | $0.15\pm 0.12$ | $0.27\pm 0.18$
April 7 Ring 3 |  | $51\pm 4.0$ | $26\pm 1.9$ | $109.1\pm 55.2$ | $0.13\pm 0.04$ | $0.19\pm 0.09$
|  | $51\pm 3.0$ | $27\pm 3.8$ | $116.2\pm 65.0$ | $0.10\pm 0.08$ | $0.19\pm 0.08$
Themis |  |  |  |  |  |
April 6 Ring | REx | $51\pm 3.9$ | $25\pm 1.2$ | $-128.6\pm 10.0$ | $0.20\pm 0.04$ | $0.27\pm 0.10$
| VIDA | $54\pm 0.9$ | $24\pm 0.9$ | $-121.0\pm 21.7$ | $0.27\pm 0.06$ | $0.34\pm 0.08$
April 7 Ring 1 |  | $53\pm 0.5$ | $23\pm 1.4$ | $-37.5\pm 11.3$ | $0.14\pm 0.01$ | $0.09\pm 0.07$
|  | $56\pm 1.4$ | $27\pm 0.9$ | $-37.6\pm 7.6$ | $0.30\pm 0.05$ | $0.19\pm 0.07$
April 7 Ring 2 |  | $53\pm 0.7$ | $22\pm 0.5$ | $-12.5\pm 8.3$ | $0.21\pm 0.02$ | $0.05\pm 0.05$
|  | $56\pm 1.2$ | $27\pm 0.7$ | $-20.6\pm 6.1$ | $0.36\pm 0.04$ | $0.15\pm 0.03$
April 7 Ring 3 |  | - | - | - | - | -
|  | - | - | - | - | -
Table 14: Mean and standard deviation of ring parameters, diameter $d$, width $w$, position angle $\eta$, asymmetry $A$, and fractional central brightness $f_{c}$ measured from Top Set and posterior on-sky Sgr A∗ images for each cluster. |  | $d\ \rm{(\mu as)}$ | $w\ \rm{(\mu as)}$ | $\eta\ \rm{(^{\circ})}$ | $A$ | $f_{c}$
---|---|---|---|---|---|---
DIFMAP |  |  |  |  |  |
April 6 Ring | REx | $46\pm 3.0$ | $34\pm 4.3$ | $97.9\pm 88.9$ | $0.11\pm 0.05$ | $0.50\pm 0.15$
| VIDA | $47\pm 1.9$ | $39\pm 3.4$ | $94.8\pm 87.0$ | $0.14\pm 0.07$ | $0.56\pm 0.10$
April 7 Ring 1 |  | $47\pm 3.0$ | $38\pm 2.7$ | $-111.6\pm 39.7$ | $0.06\pm 0.03$ | $0.59\pm 0.07$
|  | $49\pm 2.2$ | $37\pm 1.9$ | $-99.9\pm 40.2$ | $0.09\pm 0.06$ | $0.58\pm 0.08$
April 7 Ring 2 |  | $48\pm 2.4$ | $38\pm 2.6$ | $31.3\pm 81.2$ | $0.06\pm 0.04$ | $0.58\pm 0.07$
|  | $49\pm 2.0$ | $37\pm 2.1$ | $16.0\pm 75.7$ | $0.10\pm 0.06$ | $0.56\pm 0.07$
April 7 Ring 3 |  | $48\pm 2.7$ | $37\pm 2.9$ | $40.4\pm 63.5$ | $0.05\pm 0.03$ | $0.57\pm 0.08$
|  | $50\pm 2.1$ | $38\pm 1.9$ | $27.6\pm 60.9$ | $0.08\pm 0.06$ | $0.55\pm 0.08$
eht-imaging |  |  |  |  |  |
April 6 Ring | REx | $49\pm 3.9$ | $28\pm 3.6$ | $13.8\pm 123.8$ | $0.15\pm 0.05$ | $0.33\pm 0.15$
| VIDA | $50\pm 5.6$ | $41\pm 6.7$ | $110.8\pm 133.2$ | $0.20\pm 0.10$ | $0.55\pm 0.13$
April 7 Ring 1 |  | $51\pm 2.4$ | $32\pm 2.7$ | $-165.9\pm 45.0$ | $0.08\pm 0.03$ | $0.42\pm 0.10$
|  | $52\pm 4.5$ | $35\pm 3.5$ | $-170.7\pm 41.1$ | $0.11\pm 0.10$ | $0.41\pm 0.08$
April 7 Ring 2 |  | $50\pm 2.3$ | $32\pm 2.7$ | $-30.4\pm 68.1$ | $0.07\pm 0.03$ | $0.44\pm 0.07$
|  | $49\pm 2.6$ | $35\pm 4.0$ | $-71.6\pm 53.5$ | $0.06\pm 0.08$ | $0.42\pm 0.08$
April 7 Ring 3 |  | $49\pm 2.4$ | $33\pm 3.3$ | $161.6\pm 67.3$ | $0.10\pm 0.02$ | $0.50\pm 0.08$
|  | $53\pm 6.6$ | $35\pm 7.5$ | $173.3\pm 70.8$ | $0.19\pm 0.14$ | $0.52\pm 0.10$
SMILI |  |  |  |  |  |
April 6 Ring | REx | $43\pm 0.4$ | $28\pm 3.1$ | $163.1\pm 10.1$ | $0.10\pm 0.04$ | $0.57\pm 0.01$
| VIDA | $39\pm 3.9$ | $46\pm 8.8$ | $-98.6\pm 93.3$ | $0.17\pm 0.14$ | $0.67\pm 0.16$
April 7 Ring 1 |  | $43\pm 4.9$ | $33\pm 1.7$ | $127.7\pm 29.9$ | $0.09\pm 0.08$ | $0.51\pm 0.07$
|  | $45\pm 1.8$ | $34\pm 3.1$ | $136.8\pm 36.8$ | $0.07\pm 0.06$ | $0.48\pm 0.07$
April 7 Ring 2 |  | $47\pm 4.5$ | $33\pm 2.8$ | $58.4\pm 57.5$ | $0.10\pm 0.05$ | $0.48\pm 0.11$
|  | $46\pm 5.2$ | $37\pm 7.4$ | $59.8\pm 81.7$ | $0.08\pm 0.08$ | $0.50\pm 0.13$
April 7 Ring 3 |  | $48\pm 4.3$ | $32\pm 2.9$ | $89.6\pm 34.9$ | $0.11\pm 0.04$ | $0.45\pm 0.14$
|  | $47\pm 3.5$ | $35\pm 4.3$ | $84.8\pm 48.2$ | $0.11\pm 0.08$ | $0.46\pm 0.11$
Themis |  |  |  |  |  |
April 6 Ring | REx | $46\pm 2.1$ | $30\pm 1.6$ | $-139.4\pm 52.9$ | $0.15\pm 0.05$ | $0.34\pm 0.11$
| VIDA | $47\pm 3.1$ | $33\pm 2.0$ | $-127.4\pm 53.0$ | $0.17\pm 0.08$ | $0.42\pm 0.07$
April 7 Ring 1 |  | $51\pm 0.8$ | $28\pm 1.3$ | $-37.8\pm 53.6$ | $0.15\pm 0.01$ | $0.21\pm 0.10$
|  | $55\pm 1.6$ | $32\pm 2.4$ | $-41.4\pm 52.6$ | $0.25\pm 0.05$ | $0.30\pm 0.07$
April 7 Ring 2 |  | $50\pm 1.5$ | $26\pm 1.1$ | $-16.0\pm 8.9$ | $0.21\pm 0.03$ | $0.20\pm 0.08$
|  | $56\pm 2.6$ | $32\pm 0.9$ | $-25.2\pm 7.3$ | $0.35\pm 0.03$ | $0.29\pm 0.06$
April 7 Ring 3 |  | - | - | - | - | -
|  | - | - | - | - | -

##  Appendix J Dynamic Imaging and Snapshot Model Fitting tests

###  J.1 Selection of Time Windows with the Best $(u,v)$-coverage

Figure 38: Normalized metric computations for every scan of the April 7, 2017
EHT coverage of Sgr A∗. 0:00 GMST marks the day change from April 7 to April
8. Though the metrics have different considerations, all highlight the region
(labelled “II”) from ${\sim}$1:30 GMST to ${\sim}$3:10 GMST as a candidate
region for dynamic imaging.

Farah et al. (2021) demonstrates that the changing $(u,v)$-coverage created by
the Earth’s rotation during the aperture synthesis process leads to regions of
time that produce dynamic reconstructions of varying quality. The quality of a
short-timescale reconstruction is partially determined by the snapshot
$(u,v)$-coverage geometry, which introduces certain artifacts during the
imaging process.

The scale and severity of these artifacts can be predicted by quantitatively
scoring the $(u,v)$-coverage as a function of time, which can be done in a
number of ways. Some metrics examine how much of the Fourier plane is covered
by an interferometer (e.g., Palumbo et al. 2019), while others look at gaps
created by the sparse coverage (e.g., Wielgus et al. 2020). In addition to
these metrics, Farah et al. (2021) derives a novel metric that probes both the
anisotropy and radial homogeneity of the coverage.

By applying these metrics to the April 6 and 7 EHT $(u,v)$-coverage on Sgr A∗,
we can assess the scan-by-scan performance and identify regions of time which
are likely to produce the best reconstructions, independent of the underlying
source structure. The result of such an analysis for April 7 is shown in
Figure 38, and two candidate regions are highlighted. The metrics predict
dynamic imaging reconstructions will have the highest quality in the region
from 1.5-3.2 GMST (Region II); the reconstructions will produce substantially
worse results in the region from 19.4-21 GMST (Region I). We validate this
prediction by testing on high S/N data in Farah et al. (2021) and show that
Region II indeed allows for significantly better recovery of the source
variability than Region I. Therefore, based only on the EHT’s
$(u,v)$-coverage, we focus on dynamic imaging/modeling Region II throughout
Section 9.

###  J.2 Starwarps Temporal Regularizer Normalization

Figure 39:  April 7 EHT visibility amplitudes over the selected window for
dynamic imaging and modeling (black points). The horizontal lines show the
expected standard deviation of the visibility amplitudes in StarWarps
reconstructions for different values of $\beta^{-1}_{Q}$ in units of
$(\text{Jy/pixel})^{2}$: $5\times 10^{-4}$ (magenta), $5\times 10^{-6}$
(green), and $5\times 10^{-8}$ (cyan).

Temporal regularization in StarWarps is controlled by a parameter
$\beta_{Q}^{-1}$. This parameter corresponds to the variance of the
conditional distribution of pixel intensities for a given snapshot holding the
previous snapshot fixed:
$p(I_{k}|I_{k-1})=\mathcal{N}(I_{k-1},\beta^{-1}_{Q}\mathbbm{1})$. The units
of $\beta_{Q}^{-1}$ are $(\text{Jy/pixel})^{2}$. In the main text, we consider
values $\beta^{-1}_{Q}\in\\{5\times 10^{-4},5\times 10^{-6},5\times
10^{-8}\\}$. Larger values of $\beta^{-1}_{Q}$ correspond to less temporal
regularization, as the conditional distribution $p(I_{k}|I_{k-1})$ becomes
wider.

We can also interpret the values $\beta^{-1}_{Q}$ in visibility space. The
Fourier transform of the image $I_{k}$ is given by a $N_{\rm pix}\times N_{\rm
pix}$ matrix $\mathbf{F}$:

$V_{k}=\mathbf{F}I_{k}.$ (J1)

In our convention, the pixel values in $I_{k}$ have units Jy/pixel, so the
entries of $\mathbf{F}$ are pure phase terms without a $1/\sqrt{N}$
normalization (so that, for instance, the zero-baseline visibility in Jy is
just the sum of the pixel intensities). As a result
$\mathbf{F}\mathbf{F}^{\dagger}=N_{\rm pix}\mathbbm{1}$. Because Equation J1
is a linear transformation, $p(V_{k}|I_{k-1})$ is also a normal distribution,
with a mean $V_{k-1}$ and a covariance:

$\mathbf{\Sigma}=\mathbf{F}\left[\beta^{-1}_{Q}\mathbbm{1}\right]\mathbf{F}^{\dagger}=\beta^{-1}_{Q}N_{\rm
pix}\mathbbm{1}.$ (J2)

Thus, $\sigma_{\rm vis}\equiv\sqrt{\beta^{-1}_{Q}N_{\rm pix}}$ is the standard
deviation of a snapshot visibility measurement in StarWarps, holding the
previous frame fixed.

In Figure 39, we compare this quantity to the measured EHT visibility
amplitudes in the selected dynamic imaging window on April 11. The StarWarps
movie reconstructions in the main text have $N_{\rm pix}=40\times 40=1600$.
GRMHD simulations (see Paper IV and Georgiev et al. 2021) and the light curve
of Sgr A∗ (see Wielgus et al., 2022) suggest that the variations on minute
timescales should have a zero-baseline standard deviation of $\sigma_{\rm
vis}\sim 10\,{\rm mJy}$. Thus, reconstructions with $\beta^{-1}_{Q}\sim
10^{-7}({\rm Jy/pixel})^{2}$ are expected to give variability that is
consistent with what is measured in Sgr A∗.

Larger values of $\beta^{-1}_{Q}$ correspond to lower temporal regularization
and allow for larger variations in the visibility amplitudes. For instance,
our reconstructions in Section 9 with $\beta^{-1}_{Q}=5\times 10^{-6}$ permit
somewhat more variability than is seen in simulations and observations of Sgr
A∗. Nevertheless, we have also found that allowing excess variability helps to
trace evolution in tests on synthetic data from GRMHD simulations.

###  J.3 Testing $(u,v)$-coverage Effects

As discussed in Section 9.2, the geometry of the $(u,v)$-coverage can have an
effect on the recovered image structure, especially in cases where the
coverage is extremely sparse. To study the effects of $(u,v)$-coverage on
dynamic fits to Sgr A∗ data, we perform a number of tests on synthetic data
sets and study the effect of different $(u,v)$ baselines on fits to the real
data.

####  J.3.1 Recovering the Position Angle of a Static Crescent

As most of our analysis of the dynamic structure of Sgr A∗ revolves around
tracking the position angle (PA) of brightness around the ring, it is
important to assess our ability to recover the PA accurately in realistic
synthetic data. To that end, we constructed synthetic EHT data sets from 4
static crescent models with peak brightness points rotated at 60∘ increments
around the ring. The brightness ratio of each crescent model was chosen to
roughly match the 1.5:1 ratio recovered from geometric model fitting to Sgr A∗
data. Figure 40 shows the imaging and geometric modeling results obtained by
fitting to these synthetic data sets in the selected 1.7 h region. Note that,
for both approaches, the true PA is recovered as the primary mode for most of
the crescents. The imaging methods contain temporal regularization, which
likely makes it easier to recover a static underlying structure; however, the
geometric modeling results do not assume any temporal regularization.

####  J.3.2 Uniform Ring Synthetic Data

The interplay between the source size and sidelobes in the dirty beam pattern
from sparse coverage can cause imaging artifacts that appear in the form of
bright “knots" around ring sources. Computing the dirty image of an underlying
uniform ring source reveals the location of these knots when using calibrated
visibilities. To assess the impact of these knot artifacts on our analysis we
performed imaging and geometric model fitting on data generated from a uniform
ring (with no brightness changes in azimuth) with diameter of 49 $\mu$as
(refer to the uniform ring in Figure 7). As can be seen in Figure 40, both the
imaging and modeling results indicate an image structure with a preferred PA –
$\sim 0^{\circ}$ or $\sim 90^{\circ}$ for imaging and $\sim 100^{\circ}$ for
modeling. However, the associated asymmetry of the recovered uniform ring is
very low, a brightness ratio of less than 1.1:1 compared to 1.5:1 for Sgr A∗.
Thus, in combination with the results of Section J.3.1, we conclude that
although the $(u,v)$-coverage will bias the PA in the limit of low image
asymmetry, for the level of image asymmetry recovered in Sgr A∗ this bias
should have a small effect.

####  J.3.3 Baseline Test

In order to evaluate the contribution of each baseline to the recovered
evolution in Sgr A∗, we compared results obtained on data sets modified to
remove a particular baseline. In particular, we compared the PA posteriors
obtained using geometric modeling on 11 different data sets – 10 data set each
with a single baselines removed and one complete data set. As can be seen in
Figure 41, we find that most baselines do not heavily affect the trends we see
on April 7. However, there are two baselines that appear to have a significant
effect on the results: Chile–Hawaiʻi and LMT–Hawaiʻi. Without the
Chile–Hawaiʻi baseline we are not able to discriminate between the northwest
and southeast PA; without the LMT–Hawaiʻi baseline we do not recover as
significant of a PA shift. Upon inspecting the $(u,v)$-coverage of these
baselines it becomes apparent that these two baselines probe the northwest to
southeast orientations that we are interested in, and thus without them we are
unable to properly discriminate between these two PA orientations. It is also
worth noting that removal of the Chile-SPT baseline appears to “clean up" the
modeling results, suggesting that small-scale features probed by this baseline
may not be properly captured in our geometric model fits.

Figure 40: Position angle (PA) recovered from differently oriented static
crescents synthetic data sets and a uniform ring synthetic data set, using
both dynamic imaging (green) and geometric modeling (blue) techniques. The
crescents’ ground truth PA is shown as a vertical red line. Imaging uses a
prior image $\mu$ of a uniform ring convolved with a 25 $\mu$as beam. Figure
41: Position angle (PA) recovered using DPI geometric modeling after removal
of a particular baseline from Sgr A∗ data on April 7. The flagged baseline
appears in the title of each panel. The Chile-Hawaiʻi and Hawaiʻi-LMT
baselines are highlighted and their location in the $(u,v)$ plane is shown.
These two baselines probe the east–west and southeast–northwest orientations.

###  J.4 Testing Scattering Mitigation Strategies

In producing dynamic reconstructions and model fits of the Sgr A∗ data, the
choice of scattering mitigation strategy is a potential source of uncertainty.
We have explored the sensitivity of our dynamic imaging and snapshot model
fitting results to the same five scattering mitigation strategies we consider
in the static image surveys in Section 6 and Section 7.5.2. Namely, we produce
reconstructions and model fits to the unmodified data (i.e., on-sky with no
descattering), as well as with visibilities deblurred by the Sgr A∗
diffractive scattering kernel and with the thermal noise error bars inflated
by four models of the refractive noise: the Const model, the J18model1 model,
and then double the additional error tolerance from each of these models
(2$\times$Const, 2$\times$J18model1). Based on the analysis done in Section
3.1, this selection is conservative and spans our uncertainty in Sgr A∗’s
refractive noise.

Figure 42: The position angle (PA) recovered using dynamic imaging and
geometric modeling techniques under different scattering mitigation
assumptions. These range from no scattering mitigation whatsoever (yellow), to
different amounts of refractive noise added to the deblurred data: a constant
noise floor (blue), the refractive noise model J18model1 (red), and these two
scaled by a factor of two (green and magenta, respectively). Imaging results
were obtained using an initialization/prior image of a uniform ring blurred
with a Gaussian kernel with FWHM of 25 $\mu$as. Modeling results were obtained
from the DPI pipeline using the $m$-ring 2 geometric model. The gray band at
roughly 2.6 GMST indicates the region where the LMT has dropped out and data
coverage is poorer.

Figure 42 presents comparisons of the StarWarps reconstructions and snapshot
model fitting results on the April 6 and 7 Sgr A∗ data with all five
scattering mitigation strategies. We find that the general trends in the ring
position angle we discuss in Section 9.4 are not significantly changed by any
of the five scattering mitigation strategies we explore for geometric
modeling, although the position angle posteriors are significantly broader
(sometimes spanning a full 360∘) when using the larger refractive noise
budgets of 2$\times$Const and 2$\times$J18model1. In contrast, for imaging, on
April 7, we observe a transition from positive to negative PA to be the most
commonly recovered trend with all of the on-sky, Const and J18model1
scattering mitigation strategies when using a ring prior. However, when we add
a very large amount of refractive noise tolerance to the error bars in the
2$\times$Const and 2$\times$J18model1 models, the PA becomes more stable over
the observation window. This is due to the interplay of temporal
regularization with an increased flexibility in fitting the data with a static
model due to the expanded noise budget. In this figure, imaging with StarWarps
makes use of the ring$\textasteriskcentered$25$\mu$as ring
prior/initialization and modeling with DPI uses a second order $m$-ring
($m=2$) model with a parameterized central Gaussian floor.

###  J.5 Testing the Effects of Different Imaging Priors and Model
Specifications

As overviewed in Section 9.4.1, the recovered PA of the azimuthal brightness
distribution in the ring-like morphology of Sgr A∗ is sensitive to the
modeling choices made in both imaging and geometric modeling. In this section
we go into further detail on some of the effects seen.

####  J.5.1 Imaging

##### Temporal regularization

The level of continuity enforced between recovered frames is controlled by
temporal regularization. In particular, StarWarps encourages frames to be
similar by probabilistically modeling each frame $I_{k}$ as being a sample
from a normal distribution with mean $I_{k-1}$ and covariance
$\beta_{Q}^{-1}\mathbbm{1}$ (see Equation 7). Thus, decreasing the multiplier
$\beta_{Q}^{-1}$ will increase recovered continuity between frames. This is
seen, as expected, in the recovered movies of Sgr A∗ visualized in Figure 28;
recovered movies with low temporal regularization experience fast and drastic
variability on large scale features, while movies with high temporal
regularization experience slow yet steady variability on large-scale features
and absorb the remaining variability in the data with small scale
fluctuations. Due to the extreme sparsity of data on each snapshot, although
these movies contain significantly different levels of recovered variability,
they all fit the data in terms of $\chi^{2}$ fairly well. The only substantial
difference between data fits can be seen when inspecting the SMT-LMT-Hawaiʻi
closure triangle. It is worth noting that the positive to negative flip on the
April 7’s closure triangle is not reproduced by the recovered video with high
temporal regularization. Nonetheless, since all movies still match all
remaining baselines indistinguishably well it is difficult to form any solid
conclusions on the type of variability in Sgr A∗ based on this one closure
triangle.

#####  StarWarps spatial prior images

To explore the sensitivity of results to the ring features encouraged during
imaging, we introduce 5 different images used as both the initialization and
mean prior image $\mu$ in StarWarps. Figure 43 shows the mean prior images
explored for imaging the original data and deblurred data. Note that recovered
flux is constrained to evolve within the regions that have flux in the prior
image, therefore the puffier rings and the tapered disk with more extended
flux are less constraining during imaging. As expected, movies recovered with
the less constraining puffier prior images result in recovered movies with a
less clear underlying ring structures. Nonetheless, in all these cases (with a
ring init/prior), the same general trend in PA is recovered in one of the
modes, even when the central indent is very weak. When a disk prior is used
with no central indent whatsoever the same PA trend is not recovered; instead,
the PA trend appears to be reversed in sign (as discussed in Figure 27).
Figure 44 shows more detailed comparisons of StarWarps movies reconstructed
with ring and disk mean prior images on April 6 and 7, including data fits to
representative closure phases.

Figure 43: Position angle (PA) reconstructed using StarWarps under different
init/prior assumptions for both April 6 (blue) and 7 (green). These include
uniform rings with increasing Gaussian blurring and a uniform disk with no
central brightness depression. Figure 44: Comparing data fit of descattered
dynamic imaging versus static imaging on April 7. A representative descattered
image from the eht-imaging static imaging pipeline is shown. In the selected
closure plots below the measured data (60 s avg. without a variability noise
budget) is shown overlaid with the corresponding closure phases.

####  J.5.2 Geometric Modeling

For snapshot geometric modeling, we have two competing effects. One is that we
require a geometric model that can adequately explain the on-sky image.
However, given the sparseness of the $(u,v)$-coverage for each snapshot, the
risk of over-fitting the data is considerable and potentially leads to
artificially un-informative posteriors. However, under-fitting the data can
lead to large biases in the recovered parameters, and artificially narrow
posteriors. To find the preferred model we use relative measures. That is, we
don’t compare the absolute fit quality using a metric like the $\chi^{2}$
statistic, but rather how well a model does compared to the others considered.
For this purpose we use the Bayesian evidence (also referred to as evidence),

$Z({\rm data}|M)=\int\mathcal{L}({\rm
data}|\theta,M)p(\theta|M)\mathrm{d}\theta.$ (J3)

The evidence measures the marginal probability of the data, after averaging
over all possible parameter values of the model. The preferred model is then
the one that maximizes the evidence of the model. For snapshot modeling we
select the preferred model by computing the log-evidence in each snapshot, and
then sum the log-evidence across all snapshots. Note, that we are only able to
estimate the evidence for the Comrade pipeline. DPI/Variational inference
cannot estimate the evidence, but instead can compute an effective lower bound
to use as a proxy.

##### M-ring Order

To assess the impact of different model choices on the posterior samples, we
considered an $m$-ring model with 1–4 modes. The results for the $m$-ring
model from the Comrade pipeline are shown in Figure 45. We find that the trend
for the dipole moment phase is consistent across model specifications, albeit
the posteriors become more uncertain for the higher-order $m$-ring models. The
recovered total evidence for each model is shown in Figure 45. For April 6
$m=4$ is the preferred model with a log-evidence of $1499$. On April 7, the
$m=2$ $m$-ring is preferred with a log-evidence of $1667$. On both days the
overall trend of the position angle maintains stable for $m=1,2,3$, but the
distributions become noticeably wider.

For the DPI pipeline we find similar PA trends for April 6 and 7. Using the
evidence-lower-bound (ELBO) that is calculated as part of variational methods,
$m=3$ being preferred on April $6$ and $m=1$ on April 7. Furthermore, on April
7 the PA posterior for $m=3,4$ is very broad, becoming essentially
unconstrained.

Comparing the DPI results to Comrade we find that $m=1$ is preferred on April
7 and $m=3$ on April 6 according to the evidence lower bound. DPI fits closure
products which are equivalent to placing uniform priors on gains, meaning that
the data is less constraining. Therefore, it is not surprising that DPI
prefers simpler models compared to Comrade. To select a fiducial model across
days and bands we considered both the evidence, the relative impact on the
posterior, and impact on each pipeline. For these reasons we take $m=2$ for
DPI and $m=2$ and sometimes $m=3$ for Comrade. The $m=2$ model is preferred
for Comrade on April 7 and the $m=2$ distribution for DPI is similar but
broader than the $m=1$.

Figure 45: Comparing the $m$-ring results across orders 1, 2, 3, and 4 in both
Comrade (magenta) and DPI (blue). The impact of the different $m$-ring orders
on the position angle evolution for April 6 and 7 are shown in the middle and
bottom rows. The top row is the change in the log-evidence across $m$-ring
orders and days. The evidence lower bound produced by DPI is shown as a blue
upwards arrow.

## References

  * Akiyama et al. (2013) Akiyama, K., Takahashi, R., Honma, M., Oyama, T., & Kobayashi, H. 2013, PASJ, 65, 91
  * Akiyama et al. (2017a) Akiyama, K., Kuramochi, K., Ikeda, S., et al. 2017a, ApJ, 838, 1
  * Akiyama et al. (2017b) Akiyama, K., Ikeda, S., Pleau, M., et al. 2017b, AJ, 153, 159
  * Alberdi et al. (1993) Alberdi, A., Lara, L., Marcaide, J. M., et al. 1993, A&A, 277, L1
  * Arras et al. (2019) Arras, P., Frank, P., Leike, R., Westermann, R., & Enßlin, T. A. 2019, A&A, 627, A134
  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33
  * Astropy Collaboration et al. (2018) Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., et al. 2018, AJ, 156, 123
  * Baganoff et al. (2003) Baganoff, F. K., Maeda, Y., Morris, M., et al. 2003, ApJ, 591, 891
  * Balick & Brown (1974) Balick, B., & Brown, R. L. 1974, ApJ, 194, 265
  * Bardeen (1973) Bardeen, J. M. 1973, Les Astres Occlus, 215
  * Blackburn et al. (2020) Blackburn, L., Pesce, D. W., Johnson, M. D., et al. 2020, ApJ, 894, 31
  * Blackburn et al. (2019) Blackburn, L., Chan, C.-k., Crew, G. B., et al. 2019, ApJ, 882, 23
  * Bouman (2017) Bouman, K. L. 2017, PhD thesis, Massachusetts Institute of Technology
  * Bouman et al. (2018) Bouman, K. L., Johnson, M. D., Dalca, A. V., et al. 2018, IEEE Transactions on Computational Imaging, 4, 512
  * Bouman et al. (2016) Bouman, K. L., Johnson, M. D., Zoran, D., et al. 2016, in 2016 IEEE Conference on Computer Vision and Pattern Recognition (CVPR), 913
  * Bower & Backer (1998) Bower, G. C., & Backer, D. C. 1998, ApJ, 507, L117
  * Bower et al. (2004) Bower, G. C., Falcke, H., Herrnstein, R. M., et al. 2004, Science, 304, 704
  * Bower et al. (2014) Bower, G. C., Markoff, S., Brunthaler, A., et al. 2014, ApJ, 790, 1
  * Bower et al. (2015) Bower, G. C., Markoff, S., Dexter, J., et al. 2015, ApJ, 802, 69
  * Bower et al. (2018) Bower, G. C., Broderick, A., Dexter, J., et al. 2018, ApJ, 868, 101
  * Bower et al. (2019) Bower, G. C., Dexter, J., Asada, K., et al. 2019, ApJ, 881, L2
  * Boyce et al. (2019) Boyce, H., Haggard, D., Witzel, G., et al. 2019, ApJ, 871, 161
  * Brinkerink et al. (2016) Brinkerink, C. D., Müller, C., Falcke, H., et al. 2016, MNRAS, 462, 1382, (B16)
  * Broderick et al. (2022) Broderick, A., et al. 2022, in prep
  * Broderick et al. (2020a) Broderick, A. E., Gold, R., Karami, M., et al. 2020a, ApJ, 897, 139
  * Broderick et al. (2020b) Broderick, A. E., Pesce, D. W., Tiede, P., Pu, H.-Y., & Gold, R. 2020b, ApJ, 898, 9
  * Broderick et al. (2016) Broderick, A. E., Fish, V. L., Johnson, M. D., et al. 2016, ApJ, 820, 137
  * Carpenter et al. (2017) Carpenter, B., Gelman, A., Hoffman, M. D., et al. 2017, Journal of Statistical Software, 76, 1
  * Chael (2019) Chael, A. A. 2019, PhD thesis, Harvard University, Graduate School of Arts & Sciences
  * Chael et al. (2018) Chael, A. A., Johnson, M. D., Bouman, K. L., et al. 2018, ApJ, 857, 23
  * Chael et al. (2016) Chael, A. A., Johnson, M. D., Narayan, R., et al. 2016, ApJ, 829, 11
  * Cho et al. (2022) Cho, I., Zhao, G.-Y., Kawashima, T., et al. 2022, ApJ, 926, 108
  * Clark (1980) Clark, B. G. 1980, A&A, 89, 377
  * Cornwell et al. (1999) Cornwell, T., Braun, R., & Briggs, D. S. 1999, in Astronomical Society of the Pacific Conference Series, Vol. 180, Synthesis Imaging in Radio Astronomy II, ed. G. B. Taylor, C. L. Carilli, & R. A. Perley, 151
  * Cornwell (2008) Cornwell, T. J. 2008, IEEE Journal of Selected Topics in Signal Processing, 2, 793
  * Davies et al. (1976) Davies, R. D., Walsh, D., & Booth, R. S. 1976, MNRAS, 177, 319
  * Dexter et al. (2014) Dexter, J., Kelly, B., Bower, G. C., et al. 2014, MNRAS, 442, 2797
  * Do et al. (2019) Do, T., Hees, A., Ghez, A., et al. 2019, Science, 365, 664
  * Doeleman et al. (2009) Doeleman, S., Agol, E., Backer, D., et al. 2009, in Astronomy, Vol. 2010, astro2010: The Astronomy and Astrophysics Decadal Survey
  * Doeleman et al. (2008) Doeleman, S. S., Weintroub, J., Rogers, A. E. E., et al. 2008, Nature, 455, 78
  * Eckart & Genzel (1997) Eckart, A., & Genzel, R. 1997, MNRAS, 284, 576
  * Event Horizon Telescope Collaboration et al. (2019a) Event Horizon Telescope Collaboration, Akiyama, K., Alberdi, A., et al. 2019a, ApJ, 875, L1, (M87∗ Paper I)
  * Event Horizon Telescope Collaboration et al. (2019b) —. 2019b, ApJ, 875, L2, (M87∗ Paper II)
  * Event Horizon Telescope Collaboration et al. (2019c) —. 2019c, ApJ, 875, L3, (M87∗ Paper III)
  * Event Horizon Telescope Collaboration et al. (2019d) —. 2019d, ApJ, 875, L4, (M87∗ Paper IV)
  * Event Horizon Telescope Collaboration et al. (2019e) —. 2019e, ApJ, 875, L5, (M87∗ Paper V)
  * Event Horizon Telescope Collaboration et al. (2019f) —. 2019f, ApJ, 875, L6, (M87∗ Paper VI)
  * Event Horizon Telescope Collaboration et al. (2021a) —. 2021a, ApJ, 910, L12, (M87∗ Paper VII)
  * Event Horizon Telescope Collaboration et al. (2021b) —. 2021b, ApJ, 910, L13, (M87∗ Paper VIII)
  * Event Horizon Telescope Collaboration et al. (2022a) —. 2022a, submitted to ApJ, (Paper I)
  * Event Horizon Telescope Collaboration et al. (2022b) —. 2022b, submitted to ApJ, (Paper II)
  * Event Horizon Telescope Collaboration et al. (2022c) —. 2022c, submitted to ApJ, (Paper III)
  * Event Horizon Telescope Collaboration et al. (2022d) —. 2022d, submitted to ApJ, (Paper IV)
  * Event Horizon Telescope Collaboration et al. (2022e) —. 2022e, submitted to ApJ, (Paper V)
  * Event Horizon Telescope Collaboration et al. (2022f) —. 2022f, submitted to ApJ, (Paper VI)
  * Falcke et al. (1998) Falcke, H., Goss, W. M., Matsuo, H., et al. 1998, ApJ, 499, 731
  * Falcke et al. (2009) Falcke, H., Markoff, S., & Bower, G. C. 2009, A&A, 496, 77
  * Falcke et al. (2000) Falcke, H., Melia, F., & Agol, E. 2000, ApJ, 528, L13
  * Farah et al. (2021) Farah, J. R., Galison, P. L., & EHT Collaboration. 2021, ApJ, 62
  * Fish et al. (2011) Fish, V. L., Doeleman, S. S., Beaudoin, C., et al. 2011, ApJ, 727, L36
  * Fish et al. (2014) Fish, V. L., Johnson, M. D., Lu, R.-S., et al. 2014, ApJ, 795, 134
  * Fish et al. (2016) Fish, V. L., Johnson, M. D., Doeleman, S. S., et al. 2016, ApJ, 820, 90
  * Frail et al. (1994) Frail, D. A., Diamond, P. J., Cordes, J. M., & van Langevelde, H. J. 1994, ApJ, 427, L43
  * Gebhardt et al. (2011) Gebhardt, K., Adams, J., Richstone, D., et al. 2011, ApJ, 729, 119
  * Genzel et al. (2010) Genzel, R., Eisenhauer, F., & Gillessen, S. 2010, Reviews of Modern Physics, 82, 3121
  * Genzel et al. (2003) Genzel, R., Schödel, R., Ott, T., et al. 2003, Nature, 425, 934
  * Georgiev et al. (2021) Georgiev, B., A, A., B, B., & C, C. 2021, in prep., 00, 0
  * Ghez et al. (1998) Ghez, A. M., Klein, B. L., Morris, M., & Becklin, E. E. 1998, ApJ, 509, 678
  * Ghez et al. (2004) Ghez, A. M., Wright, S. A., Matthews, K., et al. 2004, ApJ, 601, L159
  * Ghez et al. (2008) Ghez, A. M., Salim, S., Weinberg, N. N., et al. 2008, ApJ, 689, 1044
  * Gillessen et al. (2009) Gillessen, S., Eisenhauer, F., Fritz, T. K., et al. 2009, ApJ, 707, L114
  * Gillessen et al. (2017) Gillessen, S., Plewa, P. M., Eisenhauer, F., et al. 2017, ApJ, 837, 30
  * Goddi et al. (2019) Goddi, C., Martí-Vidal, I., Messias, H., et al. 2019, PASP, 131, 075003
  * Goddi et al. (2021) —. 2021, ApJ, 910, L14
  * Goldreich & Sridhar (1995) Goldreich, P., & Sridhar, S. 1995, ApJ, 438, 763
  * Gravity Collaboration et al. (2018a) Gravity Collaboration, Abuter, R., Amorim, A., et al. 2018a, A&A, 618, L10
  * Gravity Collaboration et al. (2018b) —. 2018b, A&A, 615, L15
  * Gravity Collaboration et al. (2019) —. 2019, A&A, 625, L10
  * Gravity Collaboration et al. (2020) Gravity Collaboration, Bauböck, M., Dexter, J., et al. 2020, A&A, 635, A143
  * Gwinn et al. (2014) Gwinn, C. R., Kovalev, Y. Y., Johnson, M. D., & Soglasnov, V. A. 2014, ApJ, 794, L14
  * Harris et al. (2020) Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357
  * Hilbert (1917) Hilbert, D. 1917, Nachrichten von der Königlichen Gesellschaft der Wissenschaften zu Göttingen - Mathematisch-physikalische Klasse (Berlin: Weidmannsche Buchhandlung), 53
  * Högbom (1974) Högbom, J. A. 1974, A&AS, 15, 417
  * Honma et al. (2014) Honma, M., Akiyama, K., Uemura, M., & Ikeda, S. 2014, PASJ, 66, 95
  * Hora et al. (2014) Hora, J. L., Witzel, G., Ashby, M. L. N., et al. 2014, ApJ, 793, 120
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science & Engineering, 9, 90
  * Issaoun et al. (2019) Issaoun, S., Johnson, M. D., Blackburn, L., et al. 2019, ApJ, 871, 30
  * Issaoun et al. (2021) —. 2021, ApJ, 915, 99
  * Janssen et al. (2019) Janssen, M., Goddi, C., van Bemmel, I. M., et al. 2019, A&A, 626, A75
  * Johnson (2016) Johnson, M. D. 2016, ApJ, 833, 74
  * Johnson & Gwinn (2015) Johnson, M. D., & Gwinn, C. R. 2015, ApJ, 805, 180
  * Johnson & Narayan (2016) Johnson, M. D., & Narayan, R. 2016, ApJ, 826, 170
  * Johnson et al. (2015) Johnson, M. D., Fish, V. L., Doeleman, S. S., et al. 2015, Science, 350, 1242
  * Johnson et al. (2017) Johnson, M. D., Bouman, K. L., Blackburn, L., et al. 2017, The Astrophysical Journal, 850, 172
  * Johnson et al. (2018) Johnson, M. D., Narayan, R., Psaltis, D., et al. 2018, ApJ, 865, 104, (J18)
  * Johnson et al. (2020) Johnson, M. D., Lupsasca, A., Strominger, A., et al. 2020, Science Advances, 6, eaaz1310
  * Kim et al. (2020) Kim, J.-Y., Krichbaum, T. P., Broderick, A. E., et al. 2020, A&A, 640, A69
  * Krichbaum et al. (1993) Krichbaum, T. P., Zensus, J. A., Witzel, A., et al. 1993, A&A, 274, L37
  * Krichbaum et al. (1998) Krichbaum, T. P., Graham, D. A., Witzel, A., et al. 1998, A&A, 335, L106
  * Lartillot & Philippe (2006) Lartillot, N., & Philippe, H. 2006, Systematic Biology, 55, 195
  * Lee & Gammie (2021) Lee, D., & Gammie, C. F. 2021, ApJ, 906, 39
  * Levis et al. (2021) Levis, A., Lee, D., Tropp, J. A., Gammie, C. F., & Bouman, K. L. 2021, in Proceedings of the IEEE/CVF International Conference on Computer Vision, 2340
  * Lo et al. (1985) Lo, K. Y., Backer, D. C., Ekers, R. D., et al. 1985, Nature, 315, 124
  * Lo et al. (1998) Lo, K. Y., Shen, Z.-Q., Zhao, J.-H., & Ho, P. T. P. 1998, ApJ, 508, L61
  * Lu et al. (2011a) Lu, R.-S., Krichbaum, T. P., Eckart, A., et al. 2011a, A&A, 525, A76
  * Lu et al. (2011b) —. 2011b, A&A, 525, A76
  * Lu et al. (2018) Lu, R.-S., Krichbaum, T. P., Roy, A. L., et al. 2018, ApJ, 859, 60
  * Luminet (1979) Luminet, J.-P. 1979, A&A, 75, 228
  * Macquart & Bower (2006) Macquart, J.-P., & Bower, G. C. 2006, ApJ, 641, 302
  * Marcaide et al. (1999) Marcaide, J. M., Alberdi, A., Lara, L., Pérez-Torres, M. A., & Diamond, P. J. 1999, A&A, 343, 801
  * Marrone et al. (2006) Marrone, D. P., Moran, J. M., Zhao, J.-H., & Rao, R. 2006, ApJ, 640, 308
  * Martí-Vidal et al. (2016) Martí-Vidal, I., Vlemmings, W. H. T., & Muller, S. 2016, A&A, 593, A61
  * Matthews et al. (2018) Matthews, L. D., Crew, G. B., Doeleman, S. S., et al. 2018, PASP, 130, 015002
  * McMullin et al. (2007) McMullin, J. P., Waters, B., Schiebel, D., Young, W., & Golap, K. 2007, in Astronomical Society of the Pacific Conference Series, Vol. 376, Astronomical Data Analysis Software and Systems XVI, ed. R. A. Shaw, F. Hill, & D. J. Bell, 127
  * Menten et al. (1997) Menten, K. M., Reid, M. J., Eckart, A., & Genzel, R. 1997, ApJ, 475, L111
  * Meyer et al. (2008) Meyer, L., Do, T., Ghez, A., et al. 2008, ApJ, 688, L17
  * Miller-Jones et al. (2019) Miller-Jones, J. C., Tetarenko, A. J., Sivakoff, G. R., et al. 2019, Nature, 569, 374
  * Muñoz et al. (2012) Muñoz, D. J., Marrone, D. P., Moran, J. M., & Rao, R. 2012, ApJ, 745, 115
  * Narayan (1992) Narayan, R. 1992, Royal Society of London Philosophical Transactions Series A, 341, 151
  * Narayan & Nityananda (1986a) Narayan, R., & Nityananda, R. 1986a, ARA&A, 24, 127
  * Narayan & Nityananda (1986b) —. 1986b, ARA&A, 24, 127
  * Neilsen et al. (2013) Neilsen, J., Nowak, M. A., Gammie, C., et al. 2013, ApJ, 774, 42
  * Neilsen et al. (2015) Neilsen, J., Markoff, S., Nowak, M. A., et al. 2015, ApJ, 799, 199
  * Ortiz-León et al. (2016) Ortiz-León, G. N., Johnson, M. D., Doeleman, S. S., et al. 2016, ApJ, 824, 40, (O16)
  * Palumbo et al. (2019) Palumbo, D. C., Doeleman, S. S., Johnson, M. D., Bouman, K. L., & Chael, A. A. 2019, ApJ, 881, 62
  * pandas development team (2020) pandas development team, T. 2020, pandas-dev/pandas: Pandas
  * Pesce (2021) Pesce, D. W. 2021, AJ, 161, 178
  * Psaltis et al. (2018) Psaltis, D., Johnson, M., Narayan, R., et al. 2018, ArXiv e-prints, arXiv:1805.01242 [astro-ph.HE]
  * Rau (2012) Rau, U. 2012, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 8500, Image Reconstruction from Incomplete Data VII, ed. P. J. Bones, M. A. Fiddy, & R. P. Millane, 85000N
  * Rauch et al. (2016) Rauch, C., Ros, E., Krichbaum, T. P., et al. 2016, A&A, 587, A37
  * Reid (2009) Reid, M. J. 2009, International Journal of Modern Physics D, 18, 889
  * Reid & Brunthaler (2004) Reid, M. J., & Brunthaler, A. 2004, ApJ, 616, 872
  * Rickett (1990) Rickett, B. J. 1990, ARA&A, 28, 561
  * Roelofs et al. (2017) Roelofs, F., Johnson, M. D., Shiokawa, H., Doeleman, S. S., & Falcke, H. 2017, ApJ, 847, 55
  * Scherrer & Zhao (2020) Scherrer, C., & Zhao, T. 2020, in International Conference on Probabilistic Programming (Zenodo)
  * Schwab (1984) Schwab, F. R. 1984, AJ, 89, 1076
  * Schwarz (1978) Schwarz, U. J. 1978, A&A, 65, 345
  * Shen et al. (2005) Shen, Z.-Q., Lo, K. Y., Liang, M.-C., Ho, P. T. P., & Zhao, J.-H. 2005, Nature, 438, 62
  * Shepherd (2011) Shepherd, M. 2011, Difmap: Synthesis Imaging of Visibility Data, Astrophysics Source Code Library, ascl:1103.001
  * Shepherd (1997) Shepherd, M. C. 1997, in Astronomical Society of the Pacific Conference Series, Vol. 125, Astronomical Data Analysis Software and Systems VI, ed. G. Hunt & H. Payne, 77
  * Shepherd et al. (1995) Shepherd, M. C., Pearson, T. J., & Taylor, G. B. 1995, Bull. Am. Meteorol. Soc., 27, 903
  * Speagle (2020) Speagle, J. S. 2020, MNRAS, 493, 3132
  * Stewart et al. (2011) Stewart, I., Fenech, D., & Muxlow, T. 2011, Astronomy & Astrophysics, 535, A81
  * Sun & Bouman (2021) Sun, H., & Bouman, K. L. 2021, in AAAI (AAAI Press), 2628
  * Sun et al. (2022) Sun, H., Bouman, K. L., Tiede, P., et al. 2022, arXiv e-prints, arXiv:2201.08506
  * Syed et al. (2019) Syed, S., Bouchard-Côté, A., Deligiannidis, G., & Doucet, A. 2019, arXiv e-prints, arXiv:1905.02939
  * Thompson et al. (2017) Thompson, A. R., Moran, J. M., & Swenson, Jr., G. W. 2017, Interferometry and Synthesis in Radio Astronomy, 3rd Edition (Springer International Publishing)
  * Tiede (2022) Tiede, P. 2022, in prep.
  * Tiede et al. (2022) Tiede, P., Broderick, A. E., & Palumbo, D. C. M. 2022, ApJ, 925, 122
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261
  * von Laue (1921) von Laue, M. 1921, Die Relativitätstheorie, Vol. Vol. 1. Das Relativitätsprinzip der Lorentztransformation - Vol. 2. Die allgemeine Relativitätstheorie und Einsteins Lehre von der Schwerkraft (Braunschweig: Friedrich Vieweg & Sohn)
  * Wes McKinney (2010) Wes McKinney. 2010, in Proceedings of the 9th Python in Science Conference, ed. Stéfan van der Walt & Jarrod Millman, 56
  * Whitney et al. (2004) Whitney, A. R., Cappallo, R., Aldrich, W., et al. 2004, Radio Science, 39
  * Wiaux et al. (2009a) Wiaux, Y., Jacques, L., Puy, G., Scaife, A. M. M., & Vandergheynst, P. 2009a, MNRAS, 395, 1733
  * Wiaux et al. (2009b) Wiaux, Y., Puy, G., Boursier, Y., & Vandergheynst, P. 2009b, MNRAS, 400, 1029
  * Wielgus et al. (2022) Wielgus, M., Marchili, N., Martí-Vidal, I., et al. 2022, in prep., 00, 0
  * Wielgus et al. (2020) Wielgus, M., Palumbo, D., & Hamilton, L. 2020, LCG metric, https://github.com/wielgusm/mwtools/blob/master/LCG_metric.py
  * Witzel et al. (2018) Witzel, G., Martinez, G., Hora, J., et al. 2018, ApJ, 863, 15