Project AMIGA: The Circumgalactic Medium of Andromeda^*

The science goals of our HST large program require estimating the spatial distributions of the kinematics and metal column densities of the M31 CGM gas within about 1.1 ${R}_{\mathrm{vir}}$ as a function of azimuthal angle and impact parameter. The search radius was selected based on our pilot study, where we detected M31 CGM gas up to ∼1.1 ${R}_{\mathrm{vir}}$, but essentially not beyond (LHW15) (a finding that we revisit in this paper with a larger archival sample; see below). With our HST program, we observed 18 QSOs at R ≲ 1.1 ${R}_{\mathrm{vir}}$ that were selected to span the M31 projected major axis, minor axis, and intermediate orientations. The sightlines do not sample the impact parameter space or azimuthal distribution at random. Instead, the sightlines were selected to probe the azimuthal variations systematically. The sample was also limited by a general lack of identified UV-bright active galactic nuclei (AGNs) behind the northern half of M31's CGM owing to higher foreground MW dust extinction near the plane of the MW disk. Combined with seven archival QSOs, these sightlines probe the CGM of M31 in azimuthal sectors spanning the major and minor axes with a radial sample of 7–10 QSOs in each ∼100 kpc bin in R.

In addition to target locations, the 18 QSOs were optimized to be the brightest available QSOs (to minimize exposure time) and to have the lowest available redshifts (in order to minimize the contamination from unrelated absorption from high-redshift absorbers). For targets with no existing UV spectra prior to our observations, we also required that the Galaxy Evolution Explorer near-UV (NUV) and far-UV (FUV) flux magnitudes are about the same to minimize the likelihood of an intervening Lyman limit system (LLS) with optical depth at the Lyman limit τ_LL > 2. An intervening LLS could absorb more or all of the QSO flux we would need to measure foreground absorption in M31. This strategy successfully kept QSOs with intervening LLS out of the sample.

As we discuss below and as detailed by LHW15, the MS crosses through the M31 region of the sky at radial velocities that can overlap with those of M31 (see also Nidever et al. 2008; Fox et al. 2014). To understand the extent of MS contamination and the extended gas around M31 beyond the virial radius, we also searched for targets beyond 1.1 ${R}_{\mathrm{vir}}$ with COS G130M and/or G160M data. This search identified another 18 QSOs at 1.1 ≲ R/${R}_{\mathrm{vir}}$ < 1.9 that met the data quality criteria for inclusion in the sample.²² Our final sample consists of 43 sightlines probing the CGM of M31 from 25 to 569 kpc; 25 of these probe the M31 CGM from 25 to 342 kpc, corresponding to 0.08–1.1 ${R}_{\mathrm{vir}}$. Figure 1 shows the locations of each QSO in the M31–M33 system (the filled circles being targets obtained as part of our HST program PID: 14268 and the open circles being QSOs with archival HST COS G130M/G160M data), and Table 1 lists the properties of our sample QSOs ordered by increasing projected distances from M31. In this table, we list the redshift of the QSOs (${z}_{\mathrm{em}}$), the J2000 right ascension (R.A.) and declination (decl.), the MS coordinate (${l}_{\mathrm{MS}}$, ${b}_{\mathrm{MS}}$; see Nidever et al. 2008 for the definition of this coordinate system), the radially (R) and Cartesian (X, Y) projected distances, the program identification of the HST program (PID), the COS grating used for the observations of the targets, and the signal-to-noise ratio (S/N) per COS resolution element of the COS spectra near the Si iii transition.

Zoom In Zoom Out Reset image size

To search for M31 CGM absorption and to determine the properties of the CGM gas, we use ions and atoms that have their wavelengths in the UV (see Section 2.4). Any transitions with λ > 1144 Å are in the HST COS bandpass. All the targets in our sample were observed with HST using the COS G130M grating (R_λ ≈ 17,000). All the targets observed as part of our new HST program were also observed with COS G160M, and all the targets but one within R < 1.1 ${R}_{\mathrm{vir}}$ have both G130M and G160M wavelength coverage.

We also searched for additional archival UV spectra in MAST, including the FUSE (R_λ ≈ 15,000) archive to complement the gas-phase diagnostics from the COS spectra with information from the O vi absorption. We use the FUSE observations for 11 targets with adequate S/N near O vi (i.e., ≳5): RX J0048.3+3941, IRAS F00040+4325, Mrk 352, PG 0052+251, Mrk 335, UGC 12163, PG 0026+129, Mrk 1502, NGC 7469, Mrk 304, and PG 2349−014 (only the first six targets in this list are at R ≲ 1.1 ${R}_{\mathrm{vir}}$). We did not consider FUSE data for quasars without COS observations because the available UV transitions in the FUV spectrum (O vi, C ii, C iii, Si ii, Fe ii) are either too weak or too contaminated to allow for a reliable identification of the individual velocity components in their absorption profiles.

There are also three targets (Mrk 335, UGC 12163, and NGC 7469) with HST STIS E140M (R_λ ≃ 46,500) observations that provide higher-resolution information.²³

Information on the design and performance of COS, STIS, and FUSE can be found in Green et al. (2012), Woodgate et al. (1998), and Moos et al. (2000), respectively. For the HST data, we use the pipeline-calibrated final data products available in MAST. The HST STIS E140M data have an accurate wavelength calibration, and the various exposure and echelle orders are combined into a single spectrum by interpolating the photon counts and errors onto a common grid, adding the photon counts and converting back to a flux.

The processing of the FUSE data is described in detail by Wakker et al. (2003) and Wakker (2006). In short, the spectra are calibrated using version 2.1 or version 2.4 of the FUSE calibration pipeline. The wavelength calibration of FUSE can suffer from stretches and misalignments. To correct for residual wavelength shifts, the central velocities of the MW interstellar lines are determined for each detector segment of each individual observation. The FUSE segments are then aligned with the interstellar velocities implied by the STIS E140M spectra or with the velocity of the strongest component seen in the 21 cm H i spectrum. Since the O vi absorption can be contaminated by H₂ absorption, we remove this contamination following the method described in Wakker (2006). This contamination can be removed fairly accurately with an uncertainty of about ±0.1 dex on the O vi column density (Wakker et al. 2003).

For the COS G130M and G160M spectra, the spectral lines in separate observations of the same target are not always aligned, with misalignments of up to ±40 km s⁻¹ that vary as a function of wavelength. This is a known issue that has been reported previously (e.g., Savage et al. 2014; Wakker et al. 2015). While the COS team has improved the wavelength solution, we find that this problem can still be present sometimes. Since accurate alignment is critical for studying multiple gas phases probed by different ions, and since there is no way to determine a priori which targets are affected, we uniformly apply the Wakker et al. (2015) methodology to co-add the different exposures of the COS data to ensure proper alignment of the absorption lines. In short, we identify the various strong ISM and intergalactic medium (IGM) weak lines and record the component structures and identify possible contamination of the ISM lines by IGM lines. We cross-correlate each line in each exposure, using a ∼3 Å wide region, and apply a shift as a function of wavelength to each spectrum. To determine the absolute wavelength calibration, we compare the velocity centroids of the Gaussian fits to the interstellar UV absorption lines (higher velocity absorption features being Gaussian fitted separately) and the H i emission observed from our 9' GBT H i survey (Paper I) or otherwise from 21 cm data from the Leiden/Argentine/Bonn (LAB) survey (Kalberla et al. 2005) or the Parkes Galactic All Sky Survey (GASS; Kalberla et al. 2010). The alignment is coupled with the line identification into an iterative process to simultaneously determine the most accurate alignment and line identification (see Section 2.3). To combine the aligned spectra, we add the total counts in each pixel and then convert back to flux, using the average flux/count ratio at each wavelength (see also Tripp et al. 2011; Tumlinson et al. 2011); the flux error is estimated from the Poisson noise implied by the total count rate.

We are interested in the velocity range −700 ≤ ${v}_{\mathrm{LSR}}$ ≤ −150 km s⁻¹, where absorption from the M31 CGM may occur (see Section 2.4 for the motivation of this velocity range). It is straightforward to identify M31 absorption or its absence in this predefined velocity range, but we must ensure either that there is no contamination from higher-redshift absorbers or, if there is, that we can correct for it.

For ions with multiple transitions, it is relatively simple to determine whether contamination is at play by comparing the column densities and the shapes of the velocity profiles of the available transitions. The profiles of atoms or ions with a single transition can be compared to other detected ions to check whether there is some obvious contamination in the single transition absorption. However, some contamination may still remain undetected if it directly coincides with the absorption under consideration. Furthermore, when only a single ion with a single transition is detected (Si iii λ1206 being the prime example), the only method that determines whether it is contaminated or not is to undertake a complete line identification of all absorption features in each QSO spectrum.

For the 18 targets in our large HST program, our instrument setup ensures that we have the complete wavelength coverage with no gap between 1140 and 1800 Å. As part of our target selection, we also favor QSOs at low redshift (44% are at ${z}_{\mathrm{em}}$ ≤ 0.1, 89% at ${z}_{\mathrm{em}}$ ≤ 0.3). This assures that Lyα remains in the observed wavelength range out to the redshift of the QSO (Lyα redshifts out the long end of the COS band at z = 0.48) and greatly reduces the contamination from EUV transitions in the COS bandpass. The combination of wavelength coverage and low QSO redshift ensures the most accurate line identification. At R < 351 kpc (i.e., ≲1.2 ${R}_{\mathrm{vir}}$), 93% have Lyα coverage down to z = ${z}_{\mathrm{em}}$ that remains in the observed wavelength range (one target has only observation of G130M, and another QSO is at z = 0.5; see Table 1). On the other hand, for the targets at R > 351 kpc, the wavelength coverage is not as complete over 1140–1800 Å (55% of the QSOs have only one COS grating—all but one have G130M, and four QSOs have ${z}_{\mathrm{em}}$ ≳ 0.48). We note that the QSOs of 6/10 G130M observations have ${z}_{\mathrm{em}}$ < 0.17, setting all the Lyα transitions within the COS G130M bandpass.

The overall line identification process is as follows. First, we mark all the ISM absorption features (i.e., any absorption that could arise from the MW or M31) and the velocity components (which is done as part of the overall alignment of the spectra; see Section 2.2). Local (approximate) continua are fitted near the absorption lines to estimate the equivalent widths (W_λ), and their ratios for ions with several transitions are checked to determine whether any are potentially contaminated. Next, we search for any absorption features at z = ${z}_{\mathrm{em}}$, again identifying any velocity component structures in the absorption. We then identify possible Lyα absorption and any other associated lines (other H i transitions and metal transitions) from the redshift of QSOs down to z = 0. In each case, if there are simultaneous detections of Lyα, Lyβ, and/or Lyγ (and weaker transitions), we check that the equivalent width ratios are consistent. If there are any transitions left unidentified, we check whether it could be O vi λλ1031, 1037, as this doublet can sometimes be detected without any accompanying H i (Tripp et al. 2008). Finally, we check that the alignment in each absorber with multiple detected absorption lines is correct, or whether it needs some additional adjustment.

In the region R ≲ 1.1 ${R}_{\mathrm{vir}}$ and for 84% of the sample at any R, we believe that the line identifications are reliable and accurate at the 98% confidence level. In the appendices, we provide some additional information regarding the line identification, in particular for the troublesome cases. We also make available in a machine-readable format the full line identification for all the targets listed in Table 1 (see Appendix A).

Our systematic search window for absorption that may be associated with the CGM of M31 is −700 ≤ ${v}_{\mathrm{LSR}}$ ≤ −150 km s⁻¹ (LHW15). The −700 km s⁻¹ cutoff corresponds to about −100 km s⁻¹ less than the most negative velocities from the rotation curve of M31 (∼−600 km s⁻¹; see Chemin et al. 2009). The −150 km s⁻¹ cutoff is set by the MW lines that dominate the absorption in the velocity range −150 ≲ ${v}_{\mathrm{LSR}}$ ≲ +50 km s⁻¹. At −100 ≲ ${v}_{\mathrm{LSR}}$ ≲ −50 km s⁻¹, the absorption is dominated by low- and intermediate-velocity clouds that are observed in and near the MW disk. Absorption from Galactic high-velocity clouds (HVCs) is seen to velocities ${v}_{\mathrm{LSR}}$ ∼ −150 km s⁻¹ toward distant Galactic halo stars in the general direction of M31 (Lehner & Howk 2011; Lehner et al. 2012, 2015). Since the M31 disk rotation velocities extend to about −150 km s⁻¹ in the northern tip of M31, there is a small window that is inaccessible for studying the CGM of M31 (see also Lehner et al. 2015 and Section 3.2).²⁴

To search for M31 CGM gas and determine its properties, we use the following atomic and ionic transitions: O i λ1302, C ii λλ1036, 1334, C iv λλ1548, 1550, Si ii λλ1190, 1193, 1260, 1304, 1526, Si iii λ1206, Si iv λλ1393, 1402, O vi λ1031, Fe ii λλ1144, 1608, and Al ii λ1670. We also report results (mostly upper limits on column densities) for N v λλ1238, 1242, N i λ1199 (N i λλ1200, 1201 being typically blended in the velocity range of interest, −700 ≤ ${v}_{\mathrm{LSR}}$ ≤ −150 km s⁻¹), P ii λ1301, S iii λ1190, and S ii λλ1250, 1253, 1259.

To determine the column densities and velocities of the absorption, we use the apparent optical depth (AOD) method (see Section 2.4.2), but in Appendix D we confront the AOD results with column densities estimated from Voigt profile fitting (PF; see also Section 2.4.3). As much as possible at COS resolution, we derive the properties of the absorption in individual components. Especially toward M31, this is important since along the same line of sight in the velocity window −700 ≤ ${v}_{\mathrm{LSR}}$ ≤ −150 km s⁻¹, there can be multiple origins of the gas (including the CGM of M31 or MS; see Figure 1 and LHW15) as we detail in Section 2.5. However, the first step to any analysis of the absorption imprinted on the QSO spectra is to model the QSO's continuum.

2.4.1. Continuum Placement

2.4.2. Velocity Components and AOD Analysis

The next step of the analysis is to determine the velocity components and integrate them to determine the average central velocities and column densities for each absorption feature. In Figure 2, we show an example of the normalized velocity profiles. In the online figure set, we provide a similar figure for each QSO in our sample. Although we systematically search for absorption in the full velocity range −700 ≤ ${v}_{\mathrm{LSR}}$ ≤ −150 km s⁻¹, the most negative velocity of detected absorption in our sample is ${v}_{\mathrm{LSR}}$ = −508 km s⁻¹; that is, we do not detect any M31 absorption in the range −700  ≲ ${v}_{\mathrm{LSR}}$ ≲ −510 km s⁻¹. In Figure 2, MW absorption at −100  ≲ ${v}_{\mathrm{LSR}}$ ≲ 100 km s⁻¹ is clearly seen in all species but N v. Absorption observed in the −510 ≤ ${v}_{\mathrm{LSR}}$ ≤ −150 km s⁻¹ that is not color-coded is produced by higher-redshift absorbers or other MW lines.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

View figure set (43 panels)

Tarball of figure set images

To estimate the column density in each observed component, we use the AOD method (Savage & Sembach 1991). In this method, the absorption profiles are converted into apparent optical depth per unit velocity, τ_a(v) = ln[F_c(v)/F_obs(v)], where F_c(v) and F_obs(v) are the modeled continuum and observed fluxes as a function of velocity. The AOD, τ_a(v), is related to the apparent column density per unit velocity, N_a(v), through the relation N_a(v) = 3.768 × 10¹⁴τ_a(v)/(fλ (Å)) cm⁻² (km s⁻¹)⁻¹, where f is the oscillator strength of the transition and λ is the wavelength in Å. The total column density is obtained by integrating the profile over the predefined velocity interval, $N\,={\int }_{{v}_{1}}^{{v}_{2}}{N}_{a}(v){dv}$, where [v₁, v₂] are the boundaries of the absorption. We estimate the line centroids with the first moment of the AOD ${v}_{a}\,=\int v{\tau }_{a}(v){dv}/\int {\tau }_{a}(v){dv}$ km s⁻¹. As part of this process, we also estimate the equivalent widths, which we use mainly to determine whether the absorption is detected at the ≥2σ level. In cases where the line is not detected at ≥2σ significance, we quote a 2σ upper limit on the column density, which is defined as twice the 1σ error derived for the column density assuming that the absorption line lies on the linear part of the curve of growth.

For features that are detected above the 2σ level, the estimated column densities are stored for further analysis. Since we have undertaken a full identification of the absorption features in each spectrum (see Section 2.3, Appendix A), we can reliably assess whether a given transition is contaminated using in particular the conflict plots described in the appendix (see Appendix B). If there is evidence of some line contamination and several transitions are available for this ion (e.g., Si ii, Si iv, C iv), we exclude it from our list.

We find that contamination affects the Si iii and C ii in the velocity range −700 ≤ ${v}_{\mathrm{LSR}}$ ≤ −150 km s⁻¹ in a few rare cases (six components of Si iii and three components of C ii λ1334).²⁵ For all but one of these contaminated Si iii components, we can correct the contamination because the interfering line is a Lyman series line from a higher redshift and the other H i transitions constrain the equivalent width of the contamination. The one case we cannot correct this way is the −340 km s⁻¹ component toward PHL 1226 (see also Appendix A), which is associated with the MS. In the footnote of Table 2, we list the ions that are found to be contaminated at some level. For any column density that is corrected for contamination, the typical correction error is about 0.05–0.10 dex depending on the level of contamination, as well as the S/Ns of the spectrum in that region.

Table 2.  Summary of the Results

Target	Ion	v1	v2	v	σv	$\mathrm{log}N$	${\sigma }_{\mathrm{log}N}^{1}$	${\sigma }_{\mathrm{log}N}^{2}$	fN	fMS
		(km s−1)	(km s−1)	(km s−1)	(km s−1)	[cm−2]
RX J0048.3+3941	Al ii	−480.0	−320.0	−397.2	19.6	11.92	0.15	0.23	0	−1
RX J0048.3+3941	C ii	−480.0	−320.0	−381.4	1.2	14.15	0.01	0.01	−2	−1
RX J0048.3+3941	C iv	−480.0	−320.0	−390.8	5.3	13.25	0.05	0.05	0	−1
RX J0048.3+3941	Fe ii	−480.0	−320.0	−432.9	25.2	13.40	0.15	0.23	0	−1
RX J0048.3+3941	N i	−480.0	−320.0	⋯	⋯	12.88	0.18	0.30	−1	−1
RX J0048.3+3941	N v	−480.0	−320.0	⋯	⋯	12.77	0.18	0.30	−1	−1
RX J0048.3+3941	O i	−480.0	−320.0	⋯	⋯	13.47	0.18	0.30	−1	−1
RX J0048.3+3941	O vi	−480.0	−320.0	−375.5	6.2	13.85	0.05	0.06	0	−1
RX J0048.3+3941	S ii	−480.0	−320.0	⋯	⋯	13.74	0.18	0.30	−1	−1
RX J0048.3+3941	Si ii	−480.0	−320.0	−374.0	2.4	12.93	0.02	0.02	0	−1
RX J0048.3+3941	Si iii	−480.0	−320.0	−376.0	1.6	13.06	0.02	0.02	−2	−1
RX J0048.3+3941	Si iv	−480.0	−320.0	⋯	⋯	12.39	0.18	0.30	−1	−1

Note. The velocities v₁ and v₂ correspond to the integration of the absorption component. The 1σ errors ${\sigma }_{\mathrm{log}N}^{1}$ and ${\sigma }_{\mathrm{log}N}^{2}$ are the positive and negative errors, respectively, on the logarithm of the column density. The flag f_N has the following definition: 0 = detection (not saturated or contaminated); −1 = upper limit; −2 = lower limit (due to saturation of the line). The flag f_MS has the following definition: 0 = not contaminated by the MS; −1 = contaminated by the MS (see Section 2.5).

Only a portion of this table is shown here to demonstrate its form and content. A machine-readable version of the full table is available.

Download table as:  DataTypeset image

The last step is to check for any unresolved saturation. When the absorption is clearly saturated (i.e., the flux level reaches zero flux in the core of the absorption), the line is automatically marked as saturated and a lower limit is assigned to the column density. In Section 2.5, we will show how we separate the MS from the M31 CGM absorption, but we note that only the Si iii components associated with the MS and the MW have their absorption reaching zero-flux level, not the components associated with the CGM of M31.

When the flux does not reach a zero-flux level, the procedure for checking saturation depends on the number of transitions for a given ion or atom. We first consider ions with several transitions (Si ii, C iv, Si iv, sometimes C ii) since they can provide information about the level of saturation for a given peak optical depth. For ions with several transitions, we compare the column densities with different fλ-values to determine whether there is a systematic decrease in the column density as fλ increases. If there is not, we estimate the average column density using all the available measurements and propagate the errors using a weighted mean. For the Si ii transitions, Si ii λ1526 shows no evidence for saturation when detected based on the comparison with stronger transitions, while Si ii λ1260 or λ1193 can be saturated if the peak optical τ_a ≳ 0.9. For doublets (e.g., C iv, Si iv), we systematically check whether the column densities of each transition agree within 1σ error; if they do not and the weak transition gives a higher value (and there is no contamination in the weaker transition), we correct for saturation following the procedure discussed in Lehner et al. (2018) (and see also Savage & Sembach 1991). For C iv and Si iv, there is rarely any evidence for saturation (we only correct once for saturation of C iv in the third component observed in the Mrk 352 spectrum; in that component the peak optical τ_a ∼ 0.9). For single strong transitions (in particular Si iii and often C ii), if the peak optical depth is τ_a ≥ 0.9, we conservatively flag the component as saturated and adopt a lower limit for that component. We adopt τ_a ≥ 0.9 as the threshold for saturation based on other ions with multiple transitions (in particular Si ii) where the absorption starts to show some saturation at this peak optical depth.

To estimate how the column density of silicon varies with R (which has a direct consequence for the CGM mass estimates derived from silicon in Sections 4.5 and 4.8), it is useful to assess the level of saturation of Si iii, which is the only silicon ion that cannot be directly corrected for saturation.²⁶ The lower limits of the Si iii components associated with the CGM of M31 are mostly observed at R ≲ 140 kpc (only two are observed at R > 140 kpc), but they do not reach zero-flux level; these components are conservatively marked as saturated because their peak apparent optical depth is τ_a > 0.9 (not because τ_a ≫ 2) and because the comparison between the different Si ii transitions shows in some cases evidence for saturation (see above). Hence, the true values of the column densities of these saturated components are most likely higher than the adopted lower-limit values but are very unlikely to be overestimated by a factor ≫3–4. We can estimate how large the saturation correction for Si iii might be using the strong Si ii lines (e.g., Si ii λ1193 or Si ii λ1260) compared to the weaker ones (e.g., Si ii λ1526). Going through the eight sightlines showing some saturation in the components of Si iii associated with the CGM of M31 (see Table 2), for all the targets beyond 50 kpc, the saturation correction is likely to be small, <0.10–0.15 dex, based on the fact that many show no evidence of saturation in Si ii λ1260 (when there is no contamination for this transition) or Si ii λ1193. On the other hand, for the two innermost targets, the saturation correction is at least 0.3 dex and possibly as large as 0.6 dex based on the column density comparison between saturated Si ii and weaker, unsaturated transitions. The latter would put N_Si ≃ 14.5 close to the maximum values derived with photoionization modeling in the COS-Halos sample (see Section 5.2). Therefore, for the components associated with the CGM of M31 at R > 50 kpc when we estimate the functional form of N_Si with R, we adopt an increase of 0.1 dex of the lower limits. For the two inner targets at R < 50 kpc, we explore how an increase of 0.3 and 0.6 dex affects the estimation of N_Si(R).

2.4.3. High-resolution Spectra and Profile Fitting Analysis

Prior to determining the properties of the gas associated with the CGM of M31, we need to identify that gas and distinguish it from the MW and the MS. We have already removed from our analysis any contamination from higher-redshift intervening absorbers and from the MW (defined as −150 ≲ ${v}_{\mathrm{LSR}}$ ≲ 100 km s⁻¹). However, as shown in Figure 1 and discussed in LHW15, the MS is another potentially large source of contamination: in the direction of M31, the velocities of the MS can overlap with those expected from the CGM of M31. The targets in our sample have MS longitudes and latitudes in the range −132° ≤ ${l}_{\mathrm{MS}}$ ≤ −86° and $-14^\circ \,\leqslant {b}_{\mathrm{MS}}\ \leqslant +41^\circ$. The H i 21 cm emission GBT survey by Nidever et al. (2010) finds that the MS extends to about ${l}_{\mathrm{MS}}$ ≃ −140°. Based on this and previous H i emission surveys, Nidever et al. (2008, 2010) found a relation between the observed LSR velocities of the MS and ${l}_{\mathrm{MS}}$ that can be used to assess contamination in our targeted sightlines based on their MS coordinates. Using Figure 7 of Nidever et al. (2010), we estimate the upper and lower boundaries of the H i velocity range as a function of ${l}_{\mathrm{MS}}$, which we show in Figure 3 by the curved colored area. The MS velocity decreases with decreasing ${l}_{\mathrm{MS}}$ up to ${l}_{\mathrm{MS}}$ ≃ −120°, where there is an inflection point where the MS LSR velocity increases. We note that the region beyond ${l}_{\mathrm{MS}}$ ≲ −135° is uncertain but cannot be larger than shown in Figure 3 (see also Nidever et al. 2010); however, this does not affect our survey since all our data are at ${l}_{\mathrm{MS}}$ ≳ −132°.

Zoom In Zoom Out Reset image size

We take a systematic approach to removing the MS contamination that does not reject entire sightlines based on their MS coordinates since not all velocity components may be contaminated even on sightlines close to the MS. In Figure 3, we show LSR velocity of the Si iii components as a function of the MS longitude. We choose Si iii, as this ion is the most sensitive to detect both weak and strong absorption and is readily observed in the physical conditions of the MS and M31 CGM (Fox et al. 2014; Lehner et al. 2015). We consider the individual components, as for a given sightline, several components can be observed falling in or outside the boundary region associated with the MS as illustrated in Figure 3. We find that 28/74 ≃ 38% of the detected Si iii components are within the MS boundary region shown in Figure 3. We note that changing the upper boundary by ±5 km s⁻¹ would change this number by about ±3%.

To our own sample, we also add data from two different surveys: the HST/COS MS survey by Fox et al. (2014) and the M31 dwarfs (McConnachie 2012; see also Section 3). For the MS survey, we restrict the sample to −150° ≤ ${l}_{\mathrm{MS}}$ ≤ −20°, i.e., overlapping with our sample but also including higher ${l}_{\mathrm{MS}}$ values while still avoiding the Magellanic Clouds region, where conditions may be different. The origin of the sample for the M31 dwarf galaxies is fully discussed in Section 3. The larger galaxy M33 is excluded here from that sample, as its large mass is not characteristic. The LSR velocities of the M31 dwarfs as a function of ${l}_{\mathrm{MS}}$ are plotted with a star symbol in Figure 3. For the MS survey, we select the LSR velocities of Si iii for the MS survey (note that these are average velocities that can include multiple components), which are shown with squares in Figure 3. Most (∼90%) of the squares fall between the two curves in Figure 3, confirming the likelihood that these sightlines probe the MS (although we emphasize that this test was not initially used by Fox et al. 2014 to determine the association with the MS).

The M31 dwarf galaxies are of course not contaminated by the MS as the gas may be. However, if we assume that the dwarfs' companions have a similar kinematic distribution, the frequency with which the dwarfs fall within the velocity range where MS contamination is likely gives us guidance as to how frequently CGM absorption might be flagged as contaminated. For ${l}_{\mathrm{MS}}$ ≳ −132° (the range probed by the background QSOs), only 9% (2/22) of the dwarfs are within the velocity region where MS contamination occurs. If the velocity distributions of the M31 dwarfs and M31 CGM gas are similar, this suggests that gas components flagged as MS material are highly likely to be MS gas. We note, however, that two additional dwarfs just miss being included in the MS velocity region; a small increase to our MS region would change the frequency of the dwarfs consistent with MS velocities to 18%.

Observations of H i 21 cm emission toward the QSOs observed with COS in the MS survey (Fox et al. 2014) and Project AMIGA (Howk et al. 2017) show only H i detections within $| {b}_{\mathrm{MS}}| \lesssim 11^\circ$. In the region defined by −150° ≤ ${l}_{\mathrm{MS}}$ ≤ −20°, the bulk of the H i 21 cm emission is observed within $| {b}_{\mathrm{MS}}| \lesssim 5^\circ$ (Nidever et al. 2010). We therefore expect the metal ionic column densities to have a strong absorption when $| {b}_{\mathrm{MS}}| \lesssim 10^\circ$ and a weaker absorption as $| {b}_{\mathrm{MS}}|$ increases. In Figure 4, we show the total column densities of Si iii for the velocity components from the Project AMIGA sample found within the MS boundary region shown in Figure 3, i.e., we added the column densities of the components that are likely associated with the MS. We also show in the same figure the results from the Fox et al. (2014) survey. Both data sets show the same behavior of the total Si iii column densities with $| {b}_{\mathrm{MS}}|$, an overall decrease in N_{Si iii} as $| {b}_{\mathrm{MS}}|$ increases. Treating the limits as values, combining the two samples, and using the Spearman rank order, the test confirms the visual impression that there is a strong monotonic anticorrelation between N_{Si iii} and $| {b}_{\mathrm{MS}}|$ with a correlation coefficient r_S = −0.72 and a p-value ≪ 0.1%.²⁷ There is a large scatter (about ±0.4 dex around the dotted line) at any ${b}_{\mathrm{MS}}$, making it difficult to determine whether any data points may not be associated with the MS (e.g., the three very low N_{Si iii} at $12^\circ \lt | {b}_{\mathrm{MS}}| \lt 18^\circ$ from our sample or the very high value at $| {b}_{\mathrm{MS}}| \sim 27^\circ$ from the Fox et al. 2014 sample).

Zoom In Zoom Out Reset image size

In Figure 5, we show the individual column densities of Si iii as a function of the impact parameter from M31 for the Project AMIGA sightlines, where we separate components associated with the MS from those that are not. Looking at Figures 1 and 4, we expect the strongest column densities associated with the MS to be at $| {b}_{\mathrm{MS}}| \lesssim 10^\circ$ and R ≳ 300 kpc, which is where they are located on Figure 5. We also expect a positive correlation between N_{Si iii} and R for the MS contaminated components, while for uncontaminated components we expect the opposite (see LHW15). Treating again limits as values, the Spearman rank order test demonstrates a strong monotonic correlation between N_{Si iii} and R (r_S = 0.68 with p ≪ 0.1%), while for uncontaminated components there is a strong monotonic anticorrelation (r_S = −0.57 with p ≪ 0.1%), in agreement with the expectations. Based on these results, it is therefore reasonable to consider any absorption components observed in the COS spectra within the MS boundary region defined in Figure 3 as most likely associated with the MS. We therefore flag any of these components (28 out of 74 components for Si iii) as contaminated by the MS, and they are not included further in our sample.

Zoom In Zoom Out Reset image size

Finally, we noted above that only a small fraction of the dwarfs are found in the MS contaminated region. While that fraction is small (9%), this could still suggest that in the MS contaminated region some of the absorption could be a blend between both MS and M31 CGM components. However, considering the uncontaminated velocities along sightlines in (29 components) and outside (17 components) the contaminated regions, with a p-value of 0.74 the Kolmogorov–Smirnov (K-S) comparison of the two samples cannot reject the null hypothesis that the distributions are the same. This strongly suggests that the correction from the MS contamination does not bias much the velocity distribution associated with the CGM of M31 (assuming that there is no strong change of the velocity with the azimuth Φ; as we explore this further in Sections 3.2 and 4.10, there is, however, no strong evidence of a velocity dependence with Φ).

So far we have used LSR velocity to characterize MW and MS contamination of gas in the M31 halo. However, as we now consider relative motions over 30° on the sky, we cannot simply subtract M31's systemic radius velocity to place these relative motions in the correct reference frame. Over such large sky areas, tangential motion must be accounted for because the "systemic" projected radial velocity of the M31 system changes with sightline. To eliminate the effects of "perspective motion," we follow Gilbert et al. (2018) (and see also Veljanoski et al. 2014) by first transforming the heliocentric velocity (v_⊙) into the Galactocentric frame, v_Gal, which removes any effects the solar motion could have on the kinematic analysis. We converted our measured radial velocities from the heliocentric to the Galactocentric frame using the relation from Courteau & van den Bergh (1999) with updated solar motions from McMillan (2011) and Schönrich et al. (2010):

$$\begin{eqnarray}\begin{array}{rcl}{v}_{\mathrm{Gal}} & = & {v}_{\odot }+251.24\,\sin (l)\cos (b)\\ & & +\ 11.1\,\cos (l)\cos (b)+7.25\,\sin (b),\end{array}\end{eqnarray} \tag{ 1 }$$

where (l, b) are the Galactic longitude and latitude of the object. To remove the bulk motion of M31 along the sightline to each object, we use the heliocentric systemic radial velocity for M31 of −301 km s⁻¹ (van der Marel & Guhathakurta 2008; Chemin et al. 2009), which is v_M31,r = −109 km s⁻¹ in the Galactocentric velocity frame. The systemic transverse velocity of M31 is v_M31,t = −17 km s⁻¹ in the direction on the sky given by the position angle θ_t = 287° (van der Marel et al. 2012). The removal of M31's motion from the sightline velocities resulting in peculiar line-of-sight velocities for each absorber or dwarf, v_M31, is then given by (van der Marel & Guhathakurta 2008)

$$\begin{eqnarray}\begin{array}{rcl}{v}_{{\rm{M}}31} & = & {v}_{\mathrm{Gal}}-{v}_{{\rm{M}}31,{\rm{r}}}\,\cos (\rho )\\ & & +\ {v}_{{\rm{M}}31,{\rm{t}}}\,\sin (\rho )\cos (\phi -{\theta }_{t}),\end{array}\end{eqnarray} \tag{ 2 }$$

where ρ is the angular separation between the center of M31 and the QSO or dwarf position and ϕ is the position angle of the QSO or dwarf with respect to M31's center. We note that the transverse term in Equation (2) is more uncertain (van der Marel & Guhathakurta 2008; Veljanoski et al. 2014), but its effect is also much smaller, and indeed including it or not would not quantitatively change the results; we opted to include that term in the velocity transformation. We apply these transformations to change the LSR velocities to heliocentric velocities to Galactocentric velocities to peculiar velocities for each component observed in absorption toward the QSOs and for each dwarf. With this transformation, an absorber or dwarf with no peculiar velocity relative to M31's bulk motion has v_M31 = 0 km s⁻¹, regardless of its position on the sky (Gilbert et al. 2018).

In Figure 6, we compare the M31 peculiar velocities of the absorbers using Si iii and dwarfs against the projected distance (see Section 3.1). In Figure 6, we also show the expected escape velocity, v_esc, as a function of R for a 1.3 × 10¹² M_⊙ point mass. We conservatively divide v_esc by $\sqrt{3}$ in that figure to account for remaining unconstrained projection effects. Nearly all the CGM gas traced by Si iii within ${R}_{\mathrm{vir}}$ is found at velocities consistent with being gravitationally bound, and this is true even at larger R for most of the absorbers. This finding also holds for most of the dwarf galaxies, and, as demonstrated by McConnachie (2012), it holds when the galaxies' 3D distances are used (i.e., using the actual distance of the dwarf galaxies, instead of the projected distances used in this work). Therefore, both the CGM gas and galaxies probed in our sample at both small and large R are consistent with being gravitationally bound to M31.

Zoom In Zoom Out Reset image size

Figure 6 also informs us that the dwarf satellite and CGM gas velocities overlap to a high degree but do not follow identical distributions. The mean and standard deviation of the M31 velocities for the dwarfs are +27.1 ± 109.5 and +25.3 ± 67.3 km s⁻¹ for the CGM (Si iii) gas. There is therefore a slight asymmetry favoring more positive peculiar motions. A simple two-sided K-S test of the two samples rejects the null hypothesis that the distributions are the same at 95% level confidence (p = 0.04). And indeed, while the two distributions overlap and the means are similar, the velocity dispersion of the dwarfs is larger than that of the QSO absorbers.²⁹ For the QSO absorbers, all the components but one have their M31 velocities in the interval −80  ≤v_M31 ≤ +160 km s⁻¹, but 9/32 (28%) of the dwarfs are outside that range. Four of the dwarfs are in the range +160 < v_M31 ≤ +210 km s⁻¹, a velocity interval that cannot be probed in absorption owing to foreground MW contamination. The other five dwarfs have v_M31 < −80 km s⁻¹, while only one out of 46 Si iii components (2%) have v_M31 < −80 km s⁻¹. Both the small fraction of dwarfs at v_M31 > +160 km s⁻¹ and v_M31 < −80 km s⁻¹ and the even smaller fraction of absorbers at v_M31 < −80 km s⁻¹ suggest that there is no important population of absorbers at the inaccessible velocities v_M31 > +160 km s⁻¹ (see also Section 2.5).

Using the information from Table 3, we cross-match the sample of dwarf galaxies and QSOs to determine the QSO sightlines that are passing within a dwarf's R₂₀₀ radius. There are 11 QSOs (with 58 Si iii components) within R₂₀₀ of 16 dwarfs. In Table 4, we summarize the results of this cross-match. Figure 7, we show the map of the QSOs and dwarf locations in our survey, where the M31 velocities of the Si iii components and dwarfs are color-coded on the same scale and the circles around each dwarf represent their ${R}_{200}^{{\rm{d}}{\rm{w}}{\rm{a}}{\rm{r}}{\rm{f}}}$ R₂₀₀ radius.

Zoom In Zoom Out Reset image size

Table 4 and Figure 7 show that several absorbers can be found within R₂₀₀ of several dwarfs when Si iii is used as the gas tracer. For example, the two components observed in Si iii toward Zw 535.012 are found within R₂₀₀ of six dwarf galaxies. In Table 4, we also list the escape velocity (v_esc) at the observed projected distance of the QSO relative to the dwarf, as well as the velocity separation between the QSO absorber and the dwarf ($\delta v\equiv | {v}_{{\rm{M}}31,\mathrm{Si}{\rm\small{III}}}-{v}_{{\rm{M}}31,\mathrm{dwarf}}|$). So far we have not considered the velocity separation δv between the dwarf and the absorber, but it is likely that if δv ≫ v_esc then the observed gas traced by the absorber is unlikely to be bound to the dwarf galaxy even if Δ_sep = R/R₂₀₀ < 1.

If we set $\delta v\lt {v}_{\mathrm{esc}}/\sqrt{3}$, then the sample of components would be reduced to 31 instead of 58. The sample is reduced still further down to 12 if the two most massive dwarfs (M32 and NGC 205) are removed from the sample, and down to 6 if the most massive dwarfs with M_h > 3.9 × 10¹⁰ M_⊙ are removed from the sample. Applying a cross-match where δv < v_esc and Δ_sep < 1 can reduce the degeneracy between different galaxies, especially if one excludes the four most massive galaxies. For example, RXS J0118.8+3836 is located at 0.40R₂₀₀ and 0.72R₂₀₀ from Andromeda XV and Andromeda XXIII, but only in the latter case is δv ≪ v_esc (and in the former case δv > v_esc), making the two components observed toward RXS J0118.8+3836 more likely associated with Andromeda XXIII.

Several sightlines therefore pass within Δ_sep < 1 of a dwarf galaxy and show a velocity absorption within the escape velocity. This gas could be gravitationally bound to the dwarf. However, there are also five absorbers where $\delta v\lt {v}_{\mathrm{esc}}/\sqrt{3}$, but the QSO is at 1 < Δ_sep ≤ 2 from the dwarf, i.e., the velocity separation is small, but the spatial projected separation makes it unlikely to be bound to the dwarf. Here the velocity match may be a coincidence or a result of the relative proximity of the dwarfs and QSOs in the CGM of M31 assuming that the gas and dwarfs both follow the same global velocity motion of the M31 CGM. As illustrated in Figure 7, there are, however, some dwarfs with Δ_sep ≲ 2 with a radial velocity very different from that observed in absorption toward the QSO or vice versa, implying that not all the dwarfs and gas velocities are tightly connected.

In summary, it is plausible that absorbers with Δ_sep < 1 and δv ≪ v_esc trace gas associated with the CGM of a dwarf, but we cannot confirm unambiguously this association. We inspected a number of gas properties (e.g., column densities, ionization levels, kinematics) but did not find any that can differentiate clearly a dwarf CGM origin from an M31 CGM origin. Nothing in the properties of the components found within Δ_sep < 1 of a dwarf and having δv ≪ v_esc makes them outliers. This is certainly not surprising since any association assumes that the dwarf galaxies have a rich gas CGM. Yet all the satellites listed in the cross-matched Table 4 are dSph galaxies, which are known to be neutral gas poor (Grebel et al. 2003). The dSph galaxies are also likely ionized gas deficient since the favored mechanism to strip their gas is ram pressure, a stripping mechanism efficient on both the neutral and ionized gas (Grebel et al. 2003; Mayer et al. 2006). Therefore, these galaxies are unlikely to have gas-rich CGM, and based on our observations, we do not find any persuasive evidence that gas associated with M31 satellites causes the absorption we see in the M31 CGM.

Radio observations have not detected any H i 21 cm emission toward any of the QSO targets in Project AMIGA down to a 5σ level of $\mathrm{log}{N}_{{\rm{H}}{\rm\small{I}}}$ ≳ 17.6 (Paper I); many sightlines could therefore have $\mathrm{log}{N}_{{\rm{H}}{\rm\small{I}}}$ ≪ 17.6. As a consequence of this, we cannot directly estimate the metallicities of the CGM in our sample. However, we have some weak detections of O i in four components at better than the 3σ level. Since O i and H i have nearly identical ionization potentials and are strongly coupled through charge exchange reactions (Field 1971), O i is an excellent proxy for H i, requiring no or very small ionization correction as long as the photoionization spectrum is not too hard (e.g., Lehner et al. 2003). Therefore, O i can be compared to the limit of H i to put a lower limit on the metallicity. The O i logarithmic column densities are in the range of 13.3−13.7 dex (see Table 2), with a mean of 13.5 dex. This implies [O i/H i] = $\mathrm{log}{({N}_{{\rm{O}}{\rm\small{I}}}/{N}_{{\rm{H}}{\rm\small{I}}})-\mathrm{log}({\rm{O}}/{\rm{H}})}_{\odot }\gtrsim -0.7$ or a metallicity Z ≳ 0.2Z_⊙. This lower limit, however, assumes that there is no beam dilution effect, i.e., we assume that the limit on the H i column density in the 2 kpc beam (at the distance of M31) would be the same as in a pencil beam observed in absorption. Any beam dilution would increase the limit on H i, and therefore the metallicity limit could be less stringent. We therefore caution the reader not to take this limit as a hard lower limit.

While the metallicity remains quite uncertain, from the relative abundances of detected ions we can assess the level of ionization, dust depletion, and nucleosynthetic history. For assessing depletions, we can compare refractory elements like Fe to less refractory elements like Si (e.g., Savage & Sembach 1996; Jenkins 2009). Fe ii and Si ii have similar ionization energies (8–16 eV), and their observed ratio should be minimally affected by differential ionization. Hence, the ratio [Fe ii/Si ii] = $\mathrm{log}{({N}_{\mathrm{Fe}{\rm\small{II}}}/{N}_{\mathrm{Si}{\rm\small{II}}})-\mathrm{log}(\mathrm{Fe}/\mathrm{Si})}_{\odot }$ traces dust depletion levels. Unfortunately (but perhaps not surprisingly), Fe ii is only detected in the sightline closest to M31, and in that sightline we derive [Fe ii/Si ii] = −0.13 ± 0.16. In 10 other sightlines, we place upper limits on that ratio where the two smallest upper limits imply [Fe ii/Si ii] ≲ 0, while all the others are above 0 dex. While the information is minimal, this still demonstrates that there is no evidence for significant dust depletion in the CGM of M31. As depletions get stronger in denser gas, it is perhaps not surprising that we find little evidence for it in a sample where the sightlines all have $\mathrm{log}{N}_{{\rm{H}}{\rm\small{I}}}$ ≲ 17.6 and low ions are not commonly detected. While we assume that dust would be the major factor to deplete Fe relative to Si, the lack of evidence for depletion of Fe also points to a negligible nucleosynthesis effect on that ratio that would produce a nonsolar α-particle (e.g., Si) enhancement relative to Fe (e.g., Welty et al. 1997).

Using ratios of elements with different nucleosynthetic origins, we can assess the chemical enrichment history of the M31 halo gas by measuring departures from a solar relative abundance ratio in elements of different nucleosynthetic origin. For instance, the $[{\rm{C}}/\alpha ]$ ratio should be sensitive to nucleosynthesis effects since there is a time lag between the production of α-elements and carbon (see, e.g., Cescutti et al. 2009; Mattsson 2010). This analysis would be complicated by large depletions, but as we have shown above, the Fe/Si ratios show little if any evidence of large depletions. As Fe is typically the most depleted element in these conditions (Savage & Sembach 1996; Welty et al. 1999; Jenkins 2009), we can reliably assume that $[{\rm{C}}/\alpha ]$ does not suffer large depletions and can therefore be used as a nucleosynthetic indicator. We use C ii/Si ii as a proxy for $[{\rm{C}}/\alpha ]$. However, we must also consider ionization effects since differential ionization can affect the C ii/Si ii ratio (C ii has a higher ionization energy range [12–25 eV] than Si ii [8–16 eV]). To assess this, we use the nine absorbers with detections of both Si ii and C ii to estimate [C ii/Si ii] = $\mathrm{log}{({N}_{{\rm{C}}{\rm\small{II}}}/{N}_{\mathrm{Si}{\rm\small{II}}})-\mathrm{log}({\rm{C}}/\mathrm{Si})}_{\odot }$. Since this subsample includes both detections and lower limits owing to saturation of C ii, we use a survival analysis where the four censored lower limits are included (Feigelson & Nelson 1985; Isobe et al. 1986). We find that the mean [C ii/Si ii] = 0.07 ± 0.09 (where the error is the error on the mean from the Kaplan–Meier estimator) and the 1σ dispersion is 0.19 dex. This ratio is consistent with a solar value, i.e., nonsolar nucleosynthesis effects are negligible. If nondetections of Si ii are included, the mean rises to [C ii/Si ii] = 0.52 ± 0.11, strongly indicating that ionization affects this ratio owing to photons ionizing Si ii into Si iii.

Therefore, based on the relative abundances of Fe and C to Si, there is no evidence for strong dust depletion or nonsolar nucleosynthesis effects in the CGM of M31. We emphasize that this does not mean that there is no dust in the CGM of M31, and indeed several studies have shown that the CGM of galaxies can have a substantial mass of dust (e.g., Ménard et al. 2010; Peek et al. 2015). However, its effect on elemental abundances must be smaller than in the dense regions of galaxies. The lack of nucleosynthesis effects on the abundance of Fe or C relative to Si strongly suggests that the overall metallicity of the gas is not extremely low, as enhancements of α-elements are seen in low-metallicity MW halo stars and in low-metallicity gas in CGM absorbers over a range of redshift. For a sample of H i-selected absorbers with 15 ≲ $\mathrm{log}{N}_{{\rm{H}}{\rm\small{I}}}$ ≲ 18 at 0.2 ≲ z ≲ 1, Lehner et al. (2019) found little correlation between $[{\rm{C}}/\alpha ]$ and the metallicity. However, in stars and H ii regions in the local universe, there is evidence of a trend between $[{\rm{C}}/\alpha ]$ and the metallicity where $[{\rm{C}}/\alpha ]$ ≃ −0.6 at $-2\lesssim \mathrm{log}Z/{Z}_{\odot }\lesssim -0.5$ and $[{\rm{C}}/\alpha ]$ ≃ 0 near solar metallicities (e.g., Akerman et al. 2004; Fabbian et al. 2010). Therefore, the metallicity of the M31 CGM could still be subsolar, but it is unlikely to be much below 1/3 Z_⊙. This is consistent with the rough metallicity estimate set in Section 4.1.

The ionization fraction of the CGM gas can be estimated directly by comparing the column densities of O i to those of Si ii, Si iii, and Si iv (e.g., Lehner et al. 2001; Zech et al. 2008). O i is an excellent proxy for neutral gas (see Section 4.1). Si ii is found in both neutral and ionized gas, and Si iii and Si iv arise only in ionized gas. O and Si are both α-elements with similar nucleosynthetic origins and have similar levels of dust depletion in the diffuse gas (Savage & Sembach 1996; Jenkins 2009). Therefore, if the ratio [O i/Si] = $\mathrm{log}{({N}_{{\rm{O}}{\rm\small{I}}}/{N}_{\mathrm{Si}})-\mathrm{log}({\rm{O}}/\mathrm{Si})}_{\odot }$ is subsolar, ionization is important in the M31 CGM.

To obtain the total Si column density, we use the individual ion columns listed in Table 2. In the case of nondetections of the Si ions or O i, we conservatively add the upper limits to the column densities. When there are lower limits present, we add the column densities using the lower limit values. When both detections and nondetections are present, we consider the two extreme possibilities, where we either set the column density of the nondetection to the upper limit value (i.e., the absorption is nearly detected—case 1) or neglect the upper limit (i.e., it is a true nondetection—case 2). For 28 targets, we can estimate the [O i/Si] ratio. Considering case 1, we find that the mean and dispersion are [O i/Si] < −0.95 ± 0.38 and the full range is [<−1.78, <−0.34], i.e., on average the gas is ionized at the >89% level. In case 2, the mean and dispersion are [O i/Si] < −0.74 ± 0.51, so that the ionization fraction is still >81% on average. These are upper limits because typically O i is not detected. However, even in the five cases where O i is detected, four of five are upper limits too because Si iii is saturated and hence only a lower limit on the column density of Si can be derived. In that case, [O i/Si] ranges from <−1.78 to −0.43 (or to <−0.85 if we remove the absorber where the O i absorption is just detected at the 2σ level), i.e., even when O i is detected to more than 3σ, the gas is still ionized at levels >86%–98%.

The combination of Si ii, Si iii, and Si iv allows us to probe gas within the ionization energies 8–45 eV, i.e., the bulk of the photoionized CGM of M31. The high ions, C iv and O vi, have ionization energies 48–85 eV and 114–138 eV, respectively, and are not included in the above calculation. The column density of H can be directly estimated in the ionization energy 8–45 eV range from the observations via $\mathrm{log}{N}_{{\rm{H}}}$ = $\mathrm{log}{N}_{\mathrm{Si}}-\mathrm{log}Z/{Z}_{\odot }$. As we show below, Si varies strongly with R with values $\mathrm{log}{N}_{\mathrm{Si}}\gtrsim 13.7$ at R ≲ 100 kpc and $\mathrm{log}{N}_{\mathrm{Si}}\lesssim 13.3$ at R ≳ 100 kpc, which implies N_H ≳ 1.5 × 10¹⁸ (Z/Z_⊙)⁻¹ cm⁻² and ≲0.6 × 10¹⁸ (Z/Z_⊙)⁻¹ cm⁻², respectively. For the high ions, a ionization correction needs to be added, and, e.g., for O vi, $\mathrm{log}{N}_{{\rm{H}}}$ = $\mathrm{log}{N}_{{\rm{O}}{\rm\small{VI}}}-\mathrm{log}Z/{Z}_{\odot }-\mathrm{log}{f}_{{\rm{O}}\,{\rm\small{VI}}}^{i}$, where ${f}_{{\rm{O}}\,{\rm\small{VI}}}^{i}\lesssim 0.2$ is the ionization fraction of O vi that peaks around 20% for any ionizing models (e.g., Gnat & Sternberg 2007; Oppenheimer & Schaye 2013a; Lehner et al. 2014). As discussed below, there is little variation of N_{O vi} with R, and it is always such that $\mathrm{log}{N}_{{\rm{O}}{\rm\small{VI}}}$ ≳ 14.4–14.9 within 300 kpc from M31, which implies N_H ≳ (2.5–8.1) × ${10}^{18}{(Z/{Z}_{\odot })}^{-1}$ cm⁻². Therefore, the CGM of M31 not only is mostly ionized (often at levels close to 100%) but also contains a substantial fraction of highly ionized gas, with even higher column densities than the weakly photoionized gas at R ≳ 100 kpc.

In Figure 8, we show the logarithmic (left) and linear (right) values of the total column densities of the components associated with M31 for C ii, Si ii, Si iii, Si iv, C iv, and O vi as a function of the projected distances from M31. Gray data points are upper limits, while blue data with upward-pointing arrows are lower limits owing to saturated absorption. Overall, the column densities decrease at higher impact parameter. As the ionization potentials of the ions increase, the decrease in the column densities becomes shallower; O vi is almost flat. These conclusions were already noted in LHW15, but now that the region from 50 to 350 kpc is filled with data, these trends are even more striking. However, our new sample shows also an additional feature: there is a remarkable change in the quantity and quality of absorption around R₂₀₀ ≃ 230 kpc. This change is especially visible for the low and intermediate ions, whereby high column densities of C ii, Si ii, Si iii, and Si iv are observed solely at R ≲ R₂₀₀. Low column densities C ii, Si ii, Si iii, and Si iv are observed at all R, but strong absorption is observed only at R ≲ R₂₀₀. The frequency of strong absorption is also larger at R ≲ 0.6R₂₀₀ than at larger R for all ions. While the higher ionization gas traced by C iv and O vi is also weaker beyond R₂₀₀, the changes are less extreme. For C ii, Si ii, Si iii, and Si iv, the difference between low and high column densities is a factor ≳5–10, while it drops to a factor of 2–4 for C iv, and possibly even less for O vi.

Zoom In Zoom Out Reset image size

In Figure 9, we show the logarithmic values of the column densities derived from the individual components for C ii, Si ii, Si iii, Si iv, C iv, and O vi as a function of the projected distances from M31. Similar trends are observed in Figure 8, but Figure 9 additionally shows that (1) more complex velocity structures (i.e., multiple velocity components) are predominantly observed at R ≲ R₂₀₀ and (2) factor ≳2–10 changes in the column densities are observed across multiple velocity components along a given sightline.

Zoom In Zoom Out Reset image size

With Si ii, Si iii, and Si iv, we can estimate the total column density of Si within the ionization energy range 8–45 eV without any ionization modeling. Gas in this range should constitute the bulk of the cool photoionized CGM of M31 (see Section 4.3). In Figure 10, we show the total column density of Si (estimated following Section 4.3) against the projected distance R from M31. The vertical ticked bars in Figure 10 indicate data with some upper limits, and the length of the vertical bar represents the range of N_Si values allowed between cases 1 and 2 (see Section 4.3).

Zoom In Zoom Out Reset image size

Figure 10 reinforces the conclusions observed from the individual low ions in Figures 8 and 9. Overall there is a decrease of the column density of Si at larger R. This decrease has a much stronger gradient in the inner region of the M31 CGM between R ≲ 25 kpc and about R ∼ 100–150 kpc than at R ≳ 150 kpc. N_Si changes by a factor >5–10 between about 25 and 150 kpc, while it changes by a factor ≲2 between 150 and 300 kpc. The scatter in N_Si is also larger in the inner regions of the CGM than beyond ≳120–150 kpc.

To model this overall trend (which is also useful to determine the baryon and metal content of the CGM; see Section 4.8), we consider three models, a hyperbolic (H) model, a single power-law (SPL) model, and a Gaussian process (GP) model. We refer the reader to Appendix F, where we fully explain the modeling process and how lower and upper limits are accounted for in the modeling. Figure 10 shows that these three models greatly overlap. The nonparametric GP model overlaps more with the SPL model than with the H model in the range 250 kpc ≲ R ≲ 400 kpc and at R < 90 kpc (especially for the high H model; see Figure 10). While there are some differences between these models (and we will explore in Section 4.8 how these affect the mass estimates of the CGM), they all further confirm the strong evolution of the column density of Si with R between ≲25 and 90–150 kpc and a much shallower evolution with R beyond 200 kpc. In Section 4.8, we use these models to constrain the metal and baryon masses of the cool CGM gas probed by Si ii, Si iii, and Si iv.

As noted in Section 4.4, the diagnostic ions behave differently with R in a way that reflects the underlying physical conditions. For example, Si ii has a high detection rate within R < 100 kpc, a sharp drop beyond R > 100 kpc, and a total absence at R ≳ 240 kpc (see Figures 8 and 9). On the other hand, Si iii and O vi are detected at all R traced by our survey. The column densities of Si iii, however, fall significantly with R, while O vi remains relatively flat. In this section, we further quantity the detection rates, or the covering factors, for each ion.

To calculate the covering factors of the low and high ions, we follow the methodology described in Paper I for H i by assuming a binomial distribution. We assess the likelihood function for values of the covering factor given the number of detections against the total sample, i.e., the number of targets within a given impact parameter range (see Cameron 2011). As demonstrated by Cameron (2011), the normalized likelihood function for calculating the Bayesian confidence intervals on a binomial distribution with a noninformative (uniform) prior follows a β-distribution.

In Figure 11, we show the cumulative covering factors (${f}_{{\rm{c}}}$) for the various ions, where each point represents the covering factor for all impact parameters less than the given value of R. The vertical bars are 68% confidence intervals. As discussed in Sections 4.4 and 4.5, the highest column densities are found exclusively at R ≤ R₂₀₀ (and for some ions at R ≲ 100–150 kpc) for all ions except O vi. For the covering factors, we therefore consider (1) the entire sample (most of the upper limits—nondetections—are at the level of lowest column densities of a detected absorption, so it is adequate to do that) and (2) the sample where we set a threshold column density (N_th) to be included in the sample. In the left panel of Figure 11 we show the first case, while in the right panel we focus on the strong absorbers only. For the Si ions, we use $\mathrm{log}{N}_{\mathrm{th}}=13;$ for the C ions, $\mathrm{log}{N}_{\mathrm{th}}=13.8;$ for O vi, $\mathrm{log}{N}_{\mathrm{th}}=14.6$. These threshold column densities are chosen to separate strong absorbers that are mostly observed within R ≲ R₂₀₀ from the weaker absorbers that are observed at any R (see Figure 8). We also show in the right panel of Figure 11 the results for the H i emission from Paper I.

Zoom In Zoom Out Reset image size

These results must be interpreted in light of the fact that the intrinsic strength of the diagnostic lines varies by ion. The oscillator strengths, fλ, of these transitions are listed in Table 5, along with the solar abundances of these elements. The optical depth scales as τ ∝ fλN (see Section 2.4.2), and fλ is a good representation of the strength of a given transition and N of the abundance. For the Si ions, Si iii has the strongest transition, a factor of 2.7–5.5 stronger than Si iv and a factor of 1.3–5.7 stronger than Si ii (the weaker Si ii λ1526 is sometimes used to better constrain the column density of Si ii if the absorption is strong). Si ii and Si iv have more comparable strength, which is also the case between C ii and C iv. Comparing between different species, (fλ)_{Si iii} ≃ 14.4(fλ)_{O vi}, but this is counterbalanced by oxygen being 15 times more abundant than silicon (and a similar conclusion applies comparing Si iii with C ii or C iv).

With that in mind, we first consider the left panel of Figure 11, showing the cumulative covering factors for detections at any column density. We fitted four low-degree polynomials to the data: Si iii, O vi, and treating in pairs C ii+C iv and Si ii+Si iv, which follow each other well. For C ii+C iv and Si ii+Si iv, we fit the mean covering factors of each ionic pair. For O vi, we only fitted data beyond 200 kpc owing to the smaller size sample (there are only three data points within 200 kpc and 11 in total; see Figure 8). It is striking how the cumulative covering factors of Si iii and O vi vary with R quite differently from each other and from the other ions. The cumulative covering factor of Si iii appears to increase with R, but overall it is essentially consistent with about a unity covering out to 360 kpc, where it starts to decreases. At any R, the covering factor of Si iii remains much higher than ${f}_{{\rm{c}}}$ of C ii+C iv or Si ii+Si iv. The cumulative covering factor of O vi monotonically increases with R up to R ∼ 569 kpc. In contrast, while the cumulative covering factors of C ii+C iv or Si ii+Si iv are offset from each other, they both monotonically decrease with R.

Turning to the right panel of Figure 11, where we show ${f}_{{\rm{c}}}$ with column density thresholds that change with species (see above), the relation between ${f}_{{\rm{c}}}$ and R is quite different. For all the ions, the cumulative covering factors monotonically decrease with increasing R. For C ii, C iv, Si ii, Si iii, and Si iv, the covering factors are essentially the same within 1σ, and the orange line in Figure 11 shows a second-degree polynomial fit to the mean values of ${f}_{{\rm{c}}}$ between these different ions. Ignoring data at R < 200 kpc owing to the small sample size, O vi has a similar evolution of ${f}_{{\rm{c}}}$ with R, but overall ${f}_{{\rm{c}}}$ is tentatively ∼1.5 times larger than for the other ions at any R.

The contrast between the two panels of Figure 11 and insight from Figures 8 and 9 strongly suggest that the CGM of M31 has three main populations of absorbers: (1) the strong absorbers that are found mostly at R ≲ 100–150 kpc (0.3–0.5 ${R}_{\mathrm{vir}}$) probing the denser regions and multiple gas phase (singly to highly ionized gas) of the CGM, (2) weak absorbers probing the diffuse CGM traced principally by Si iii (but also observed in higher ions and more rarely in C ii) that are found at any surveyed R but more frequent at R ≲ R₂₀₀, and (3) hotter, more diffuse CGM probed by O vi, O vi having the unique property compared to the ions that its column density remains largely invariant with the radius of the M31 CGM.

In Section 4.3, we show that the ratio of O i to Si ions provides a direct estimate of the ionization fraction of the CGM gas of M31. Using ratios of the main ions studied here (C ii, C iv, Si ii, Si iii, Si iv, O vi), we can further constrain the ionization and physical conditions in the CGM of M31 and how they may change with R. To estimate the ionic ratios, we consider the component analysis of the absorption profiles, i.e., we compare the column densities estimated over the same velocity range, noting, however, that coincident velocities do not necessarily mean that they probe the same gas, especially if their ionization potentials are quite different (such as for C ii and C iv). In Figure 12, we show the results for several ion ratios as a function of R.

Zoom In Zoom Out Reset image size

4.7.1. The Si ii/Si iii and Si iv/Si iii Ratios

4.7.2. The C ii/C iv Ratio

4.7.3. The C iv/Si iv Ratio

4.7.4. The C iv/O vi Ratio

With a better understanding of the column density variation with R, we can estimate with more confidence the metal and baryon masses of the M31 CGM than in LHW15, which had very little information between 50 and 300 kpc. The metal mass can be directly estimated from the column densities of the metal ions. With the silicon ions, we have information on its three dominant ionization stages in the T < 7 × 10⁴ K ionized gas (ionizing energies in the range of 8–45 eV; see Section 4.5), so we can obtain a direct measured metal mass without any major ionization corrections. Following LHW15 (see also Peeples et al. 2014), the metal mass of the cool photoionized CGM is

$$\begin{eqnarray*}&&{M}_{{\rm{Z}}}^{\mathrm{cool}}=2\pi \,{\mu }_{\mathrm{Si}}^{-1}\,{m}_{\mathrm{Si}}\,\int R\,{N}_{\mathrm{Si}}(R)\,{dR},\end{eqnarray*}$$

where μ_Si = 0.064 is the solar mass fraction of metals in silicon (i.e., $12+\mathrm{log}{(\mathrm{Si}/{\rm{H}})}_{\odot }=7.51$ and Z_⊙ = 0.0142 from Asplund et al. 2009), m_Si = 28m_p, and for N_Si(R) we use the H (Equation (F1)), SPL (Equation (F2)), and GP models that we determine in Section 4.5 and Appendix F (see Figure 10).

A direct method to estimate the total mass is to convert the total observed column density of Si to total hydrogen column density via N_H = N_{H i} + N_{H ii} = N_Si (Si/H)_⊙⁻¹ (Z/Z_⊙)⁻¹. The baryonic mass of the CGM of M31 is then

$$\begin{eqnarray*}&&{M}_{{\rm{g}}}^{\mathrm{cool}}=2\pi \,{m}_{{\rm{H}}}\,\mu \,{f}_{{\rm{c}}}\ {\left(\displaystyle \frac{\mathrm{Si}}{{\rm{H}}}\right)}_{\odot }^{-1}{\left(\displaystyle \frac{Z}{{Z}_{\odot }}\right)}^{-1}\,\int R\,{N}_{\mathrm{Si}}(R){dR},\end{eqnarray*}$$

where μ ≃ 1.4 (to correct for the presence of He), m_H = 1.67 × 10⁻²⁴ g is the hydrogen mass, ${f}_{{\rm{c}}}$ is the covering fraction (which is 1 over the considered radii), and $\mathrm{log}{(\mathrm{Si}/{\rm{H}})}_{\odot }\,=-4.49$ is the solar abundance of Si. Inserting the values for each parameter, ${M}_{{\rm{g}}}^{\mathrm{cool}}$ can be simply written in terms of ${M}_{{\rm{Z}}}^{\mathrm{cool}}$: ${M}_{{\rm{g}}}^{\mathrm{cool}}\simeq {10}^{2}{(Z/{Z}_{\odot })}^{-1}{M}_{{\rm{Z}}}^{\mathrm{cool}}$.

In Table 6, we summarize the estimated metal mass over different regions of the CGM for the three models of N_Si(R), within R₂₀₀ (first entry), within ${R}_{\mathrm{vir}}$ (second entry), within 1/2 ${R}_{\mathrm{vir}}$ (third entry), between 1/2 ${R}_{\mathrm{vir}}$ and ${R}_{\mathrm{vir}}$ (fourth entry), and within 360 kpc (fifth entry), which corresponds to the radius where at least one of the Si ions is always detected (beyond that, the number of detections drastically plummets). A key difference between the H/SPL models and the GP model is that the range of values for the H/SPL models is derived using the low (dotted) and high (dashed) curves in Figure 10, while for the GP models we actually use the confidence intervals from the low and high models (i.e., the top and bottom of the shaded blue curve in Figure 10). Hence, it is not surprising that the mass ranges for the GP model are larger. Nevertheless, there is a large overlap between the three models. As the GP results overlap with the other models and provide empirical confidence intervals, we adopt them for the remaining of the paper. Within ${R}_{\mathrm{vir}}$, the metal and cool gas masses are therefore (2.0 ± 0.5) × 10⁷ and 2 × 10⁹ (Z/Z_⊙)⁻¹ M_⊙, respectively. Owing to the new functional form of N_Si(R) and how the lower limits are treated, this explains the factor of 1.4 increase in the metal mass compared to that derived in LHW15.

These masses do not include the more highly ionized gas traced by O vi, C iv, or even higher (unobserved) ionization states. Even though the sample with O vi is smaller than C iv, we use O vi to probe the higher-ionization gas phase because, as we show above, the properties of O vi (column density and covering fraction as a function of R) are quite different from those of all the lower ions, including C iv, which behaves more like the other, lower ions. Furthermore, Lehner et al. (2011), using 1.5–3 km s⁻¹ resolution UV spectra, show that C iv can probe cool and hotter gas, while the profiles of N v and O vi are typically broad and more consistent with hotter gas. Since O vi is always detected and there is little evidence for variation with R (see Figure 8), we can simply use the mean column density $\mathrm{log}{N}_{{\rm{O}}{\rm\small{VI}}}$ = 14.46 ± 0.10 (error on the mean using the survival method for censoring) to estimate the baryon mass assuming a spherical distribution:

$$\begin{eqnarray*}&&{M}_{{\rm{g}}}^{\mathrm{warm}}=\pi {r}^{2}\,{m}_{{\rm{H}}}\,\mu \,{f}_{{\rm{c}}}\,\displaystyle \frac{{N}_{{\rm{O}}{\rm\small{VI}}}}{{f}_{{\rm{O}}\,{\rm\small{VI}}}^{i}}{\left(\displaystyle \frac{{\rm{O}}}{{\rm{H}}}\right)}_{\odot }^{-1}{\left(\displaystyle \frac{Z}{{Z}_{\odot }}\right)}^{-1},\end{eqnarray*}$$

where the O vi ionization fraction is ${f}_{{\rm{O}}\,{\rm\small{VI}}}^{i}\lesssim 0.2$ (see Section 4.3) and ${f}_{{\rm{c}}}$ = 1 for O vi at any R (see Figure 8). Within ${R}_{\mathrm{vir}}$, we find M^warm_g ≳ 9.3 × 10⁹(Z/Z_⊙)⁻¹ M_⊙ or ${M}_{{\rm{g}}}^{\mathrm{warm}}\gtrsim 4.4\,{M}_{{\rm{g}}}^{\mathrm{cool}}$ (assuming that the metallicity is about similar in the cool and warm gas phases). At R₂₀₀, we find ${M}_{{\rm{g}}}^{\mathrm{warm}}\gtrsim 5.5\times {10}^{9}{(Z/{Z}_{\odot })}^{-1}$ M_⊙. These are lower limits because the fraction of O vi could be much smaller than 20% and the metallicity of the cool or warm ionized gas is also likely to be less than solar (see below). In terms of metal mass in the highly ionized gas phase, this gives ${M}_{{\rm{g}}}^{\mathrm{warm}}\,\simeq {10}^{2}{(Z/{Z}_{\odot })}^{-1}{M}_{{\rm{Z}}}^{\mathrm{warm}}$ and ${M}_{{\rm{Z}}}^{\mathrm{warm}}\gtrsim 4.4\,{M}_{{\rm{Z}}}^{\mathrm{cool}}$. Since O vi is detected out to the maximum surveyed radius of 569 kpc, within that radius (i.e., 1.9 ${R}_{\mathrm{vir}}$) we find ${M}_{{\rm{g}}}^{\mathrm{warm}}\gtrsim 34\,\times {10}^{9}{(Z/{Z}_{\odot })}^{-1}$ M_⊙.

By combining both the cool and hot gas-phase masses, we can find the baryon mass for gas in the temperature range of ∼10⁴–10^5.5 K within ${R}_{\mathrm{vir}}$:

$$\begin{eqnarray*}\begin{array}{rcl}{M}_{{\rm{g}}} & = & {M}_{{\rm{g}}}^{\mathrm{cool}}+{M}_{{\rm{g}}}^{\mathrm{warm}}\\ & \gtrsim & \ 1.1\times {10}^{10}{\left(\displaystyle \frac{Z}{{Z}_{\odot }}\right)}^{-1}\,{M}_{\odot }.\end{array}\end{eqnarray*}$$

Within R₂₀₀, the total mass M_g ≳ 7.2 × 10⁹ M_⊙. As the stellar mass of M31 is about 10¹¹ M_⊙ (e.g., Geehan et al. 2006; Tamm et al. 2012), the mass of the diffuse weakly and highly ionized CGM of M31 within 1 ${R}_{\mathrm{vir}}$ is therefore at least 10% of the stellar mass of M31 and could be significantly larger than 10%.

This estimate does not take into account the hot (T ≳ 10⁶ K) coronal gas. The diffuse X-ray emission is observed to extend to about 30–70 kpc around a handful of massive, nonstarbursting galaxies (Anderson & Bregman 2011; Bregman et al. 2018) or in stacked images of galaxies (Anderson et al. 2013; Bregman et al. 2018), but beyond 50 kpc, the CGM is too diffuse to be traced with X-ray imaging, even though a large mass could be present. Using the results summarized recently by Bregman et al. (2018), the hot gas mass of spiral galaxy halos is in the range ${M}_{{\rm{g}}}^{\mathrm{hot}}\simeq (1$–10) × 10⁹ M_⊙ within 50 kpc (assuming or estimating metallicities in the range 0.13–0.5Z_⊙). For M31, ${M}_{{\rm{g}}}={M}_{{\rm{g}}}^{\mathrm{cool}}+{M}_{{\rm{g}}}^{\mathrm{warm}}\,\gtrsim 0.4\times {10}^{9}$ M_⊙ within 50 kpc. Extrapolating the X-ray results to ${R}_{\mathrm{vir}}$, Bregman et al. (2018) find masses of the hot X-ray gas similar to the stellar masses of these galaxies in the range ${M}_{{\rm{g}}}^{\mathrm{hot}}\simeq (1$–10) × 10¹¹ M_⊙. For the MW hot halo within 1 ${R}_{\mathrm{vir}}$, Gupta et al. (2017; but see also Gupta et al. 2012, 2014; Wang & Yao 2012; Henley et al. 2014) derive (3–10) × 10¹⁰ M_⊙, i.e., on the low side of the mass range listed in Bregman et al. (2018). The hot gas could therefore dominate the mass of the CGM of M31. There are, however, two caveats to that latter conclusion. First, if f_{O vi} ≪ 0.2, then M^warm_g could become much larger than the lower limits we give. Second, the metallicity of the hot X-ray gas ranges from 0.1 to 0.5Z_⊙ with a mean metallicity of 0.3Z_⊙ (Bregman et al. 2018; Gupta et al. 2017), while for the cooler gas we have conservatively adopted a solar abundance. If instead we adopt a 0.3Z_⊙ metallicity (consistent with the rough limits set in Sections 4.1, 4.2), then M_g ≃ 3.7 × 10¹⁰ M_⊙ within ${R}_{\mathrm{vir}}$, which is now comparable to the hot halo mass of the MW. If we adopt the average metallicity derived for the X-ray gas, then ${M}_{{\rm{g}}}^{\mathrm{cool}}+{M}_{{\rm{g}}}^{\mathrm{warm}}$ would be comparable to the hot gas mass if ${M}_{{\rm{g}}}^{\mathrm{hot}}\sim 5\times {10}^{10}$ M_⊙ within ${R}_{\mathrm{vir}}$ for M31. Depending on the true metallicities and the actual state of ionization, the cool and warm gas in the M31 halo could therefore contribute to a substantial enhancement of the total baryonic mass compared to our conservative assumptions.

Thus far, we have ignored the distribution of the targets in azimuthal angle (Φ) relative to the projected minor and major axes of M31, where different physical processes may occur. In Figure 13, we show the distribution of the column densities of each ion in the X-Y plane near M31, where the circles represent detections and downward-pointing triangles are nondetections. Multiple colors in a given circle indicate multiple components along that sightline for that ion. We also show the projected minor and major axes of M31 (dashed lines). The overall trends that are readily apparent from Figure 13 are the ones already described in the previous sections: (1) overall the column density decreases with increasing R, (2) the decrease in N is much stronger for low ions than high ions, and (3) Si iii and O vi are observed at all R probed, while singly ionized species tend to be more frequently observed at small impact parameters. This figure (and Figure 9) also reveals that absorption with two or more components is observed more frequently at R < 200 kpc: using Si iii, 64%–86% of the sightlines have at least two velocity components at R < 200 kpc, while this drops to 14%–31% at R > 200 kpc (68% confidence intervals using the Wilson score interval); similar results are found using the other ions. However, the complexity of the velocity profiles does not change with Φ.

Zoom In Zoom Out Reset image size

Considering various radius ranges (25–50 kpc, 50–100 kpc, etc.) up to 1 ${R}_{\mathrm{vir}}$, there is no indication that the column densities depend strongly on Φ. Considering Si iii first, it is equally detected along the projected major and minor axes and in between (wherever there is a sightline), and overall the strength of the absorption mostly depends on R, not Φ. Considering the other ions, they all show a mixture of detections and nondetections, and the nondetections (which are mostly beyond 50 kpc) are not preferentially observed along a certain axis or one of the regions shown in Figure 13. We therefore find no strong evidence of an azimuthal dependence in the column densities.

Beyond ≳1.1 ${R}_{\mathrm{vir}}$, the situation is different, with all but one detection (in C iv and O vi only) being near the southern projected major axis and about 52° west off near the X = 0 kpc axis. There is detection in this region of Si iii, C iv, Si iv, O vi, and also C ii. That is the main region where C ii is detected beyond 200 kpc. In contrast, between the X = 0 kpc axis and southern projected minor axis, the only region where there are several QSOs beyond ${R}_{\mathrm{vir}}$, there is no detection in any of the ions (excluding O vi because there are no FUSE observations in these directions). Although that direction is suspiciously in the direction of the MS, it is very unlikely to be additional contamination from the MS because (1) the velocities are inconsistent with the MS expectations in these directions (see Figure 3 and also Section 4.10), and (2) we see a decrease in the column densities as $| {b}_{\mathrm{MS}}|$ decreases (closer to the core of the MS), the opposite of expectations if the material is associated with the MS. Therefore, while at R < ${R}_{\mathrm{vir}}$ there is no apparent trend between N and Φ for any ions (although we keep in mind that the azimuthal information for O vi is minimal), most of the detections at R > ${R}_{\mathrm{vir}}$ are near the southern projected major axis and 52° west off of that axis.

The fact that the gas is observed mainly in a specific region of the CGM beyond ${R}_{\mathrm{vir}}$ suggests an IGM filament feeding the CGM of M31, as is observed in some cosmological simulations. In particular, Nuza et al. (2014) study the gas distribution in simulated recreations of MW and M31 using a constrained cosmological simulation of the Local Group from the Constrained Local UniversE Simulations (CLUES) project. In their Figures 3 and 6, they show different velocity and density projection maps where the central galaxy (M31 or MW) is edge-on. They find that some of the gas in the CGM can flow in a filament-like structure, coming from outside the virial radius all the way down to the galactic disk.

How the velocity field of the gas is distributed in R and Φ beyond 25–50 kpc is a key diagnostic of accretion and feedback. However, a statistical survey using one sightline per galaxy (such as COS-Halos) cannot address this problem because it observes many galaxies in an essential random mix of orientations and inclinations, which necessarily washes out any coherent velocity structures. An experiment like Project AMIGA is needed to access information about large-scale flows in a sizable sample of lines of sight for a single galaxy. The velocity information remains limited because we have only the (projected) radial velocity along pencil beams piercing the CGM at various R and Φ. Nevertheless, as we show below, some trends are apparent thanks to the large size of the sample. We use here the v_M31 peculiar velocities as defined by Equation (2). By definition, in the M31 velocity frame, an absorber with no peculiar velocity relative to M31's bulk motion has v_M31 = 0 km s⁻¹. In Section 3.2, we show that the M31 peculiar velocities of the absorbers seen toward the QSOs and the velocities of the M31 dwarf satellites largely overlap. We now review how the velocities of the absorbers are distributed in the CGM of M31 over the entire surveyed range of R.

In Figures 14 and 15, we show the distribution of the M31 peculiar velocities of the individual components identified for each ion and the column-density-weighted average velocities of each ion, respectively. Circles with several colors again indicate that the observed absorption appears in more than one component. Both Figures 14 and 15 demonstrate that in many cases there is some overlap in the velocities between lower ions (Si ii, C ii, Si iii) and higher ions (Si iv, C iv, O vi). This strongly implies that the CGM of M31 has multiple gas phases with overlapping kinematics when they are observed in projection (a property also readily observed from the normalized profiles shown in Figure 2 and in the figure set in the online Journal). There are also some rarer cases where there is no velocity correspondence in the velocities between Si iii and higher ions (see, e.g., near X ≃ −335, Y ≃ −95 kpc), indicating that the observed absorption in each ion is dominated by a single phase–that is, the components are likely to be distinct single-phase objects.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The full range of velocities associated with the CGM of M31 is −249  ≤ v_M31 ≤ +175 km s⁻¹ for Si iii, but for all the other ions it is −53 ≲ v_M31 ≤ +175 km s⁻¹. Furthermore, there is only one absorber/component of Si iii that has v_M31 = −249 km s⁻¹. We emphasize that the rarity of velocity v_M31 < −249 km s⁻¹ (corresponding to ${v}_{\mathrm{LSR}}$ < −510 km s⁻¹ in the direction of this sightline) is not an artifact since velocities below these values are not contaminated by any foreground gaseous features.

We show in Section 3.2 that the M31 dwarf satellites have a velocity dispersion that is larger (110 km s⁻¹ for the dwarfs vs. 68 km s⁻¹ for the Si iii absorbers) and the M31 dwarfs have some velocities in the velocity range contaminated by the MW and MS. While the CGM gas velocity field distribution may not follow that of the dwarf satellites, it remains plausible that some of the absorption from the extended region of the M31 CGM could be lost owing to contamination from the MW or MS, i.e., we may not be fully probing the entire velocity distribution of the M31 CGM. However, as discussed in Section 2.5, there is no evidence that the velocity distributions of the Si iii component in and outside the MS contamination zone are different (see also Figures 1 and 14), and hence it is quite possible that at least the MS contamination does not affect much the velocity distribution of the M31 CGM. With these caveats, we now proceed to describe the apparent trends of the velocity distribution in the CGM of M31.

From Figure 14, the first apparent property was already noted in the previous section: the velocity complexity (and hence full width) of the absorption profiles increases with decreasing R (see Section 4.9). Within R ≲ 100 kpc or ≲200 kpc, about 75% of the Si iii absorbers have at least two components (at the COS G130M−G160M resolution). This drops to about 33% at 200  < R ≲ 569 kpc.

The second property evident from either Figure 14 or Figure 15 is that the M31 peculiar velocities are larger at R ≲ 100 kpc than at higher R. Table 7 lists the average M31 velocities, their standard deviations, and their interquartile ranges (IQRs) for the individual components and averaged components in three samples: the full AMIGA set, the subset with R ≤ 100 kpc, and the subset with R > 100 kpc. From this table and for all the ions besides O vi, $\langle {v}_{{\rm{M}}31}\rangle =+90$ km s⁻¹ at R ≤ 100 kpc, while at R > 100 kpc, $\langle {v}_{{\rm{M}}31}\rangle =+20$ km s⁻¹, a factor of 4.5 smaller. There are only two data points for O vi, at R ≤ 100 kpc, but the average at R > 100 kpc is also $\langle {v}_{{\rm{M}}31}\rangle =+22$ km s⁻¹, following a similar pattern to that observed for the other ions. For all the ions but C ii, the velocity dispersions or IQRs are smaller at R ≤ 100 kpc than at R > 100 kpc.

Table 7.  Summary of the Velocities of the CGM of M31

Ion	Region	$\langle {v}_{{\rm{M}}31}\rangle$	IQR
		(km s−1)	(km s−1)
Individual Components
C ii	All R	31.3 ± 52.3	[−6.5, 68.0]
C ii	R ≤ 100 kpc	80.6 ± 39.9	[37.9, 101.1]
C ii	R > 100 kpc	4.7 ± 36.6	[−21.2, 24.9]
Si ii	All R	65.7 ± 53.1	[53.8, 95.2]
Si ii	R ≤ 100 kpc	90.6 ± 32.3	[67.2, 117.8]
Si ii	R > 100 kpc	22.3 ± 54.5	[−27.0, 76.8]
Si iii	All R	26.8 ± 45.5	[−10.1, 54.7]
Si iii	R ≤ 100 kpc	77.9 ± 34.1	[59.9, 95.9]
Si iii	R > 100 kpc	7.4 ± 32.3	[−12.4, 36.1]
Si iv	All R	39.5 ± 43.1	[8.3, 71.0]
Si iv	R ≤ 100 kpc	82.2 ± 7.5	[76.3, 88.4]
Si iv	R > 100 kpc	12.8 ± 33.7	[−0.6, 26.9]
C iv	All R	42.6 ± 47.3	[10.8, 70.2]
C iv	R ≤ 100 kpc	88.3 ± 23.8	[71.1, 102.1]
C iv	R > 100 kpc	13.5 ± 33.4	[4.9, 40.4]
O vi	All R	28.9 ± 30.8	[13.1, 43.0]
Averaged Components
C ii	All R	34.2 ± 58.9	[−23.0, 68.1]
C ii	R ≤ 100 kpc	67.0 ± 55.3	[22.7, 92.4]
C ii	R > 100 kpc	10.5 ± 49.3	[−33.2, 37.2]
Si ii	All R	68.1 ± 50.6	[54.2, 98.7]
Si ii	R ≤ 100 kpc	86.4 ± 35.1	[61.5, 111.2]
Si ii	R > 100 kpc	22.3 ± 54.5	[−27.0, 76.8]
Si iii	All R	25.3 ± 67.3	[−18.9, 66.9]
Si iii	R ≤ 100 kpc	69.3 ± 48.1	[43.8, 103.3]
Si iii	R > 100 kpc	4.1 ± 64.9	[−24.0, 47.7]
Si iv	All R	39.6 ± 55.4	[−2.5, 78.8]
Si iv	R ≤ 100 kpc	69.9 ± 56.2	[52.3, 112.4]
Si iv	R > 100 kpc	14.9 ± 40.4	[−22.9, 50.3]
C iv	All R	43.3 ± 55.7	[8.4, 84.4]
C iv	R ≤ 100 kpc	91.3 ± 33.7	[69.9, 107.4]
C iv	R > 100 kpc	14.0 ± 45.2	[−8.0, 52.8]
O vi	All R	27.9 ± 44.0	[−9.8, 68.5]

Note. Average, standard deviation, and IQR velocity range are listed for the individual components and averaged M31 velocity components.

Download table as:  ASCIITypeset image

The third property observed in Figure 14 or Figure 15 is that at R ≤ 100 kpc there is no evidence for negative M31 velocities, while at R > 100 kpc about 40% of the Si iii sample has blueshifted v_M31 velocities. This partially explains the previous result, but even if we consider the absolute velocities, $\langle | {v}_{{\rm{M}}31}| \rangle =+40$ km s⁻¹ at R > 100 kpc, implying $\langle | {v}_{{\rm{M}}31}(R\,\gt 100)| \rangle =0.44 \langle | {v}_{{\rm{M}}31}| (R\leqslant 100) \rangle$, i.e., in absolute terms or not, v_M31 is smaller at R > 100 kpc than at R ≤ 100 kpc. Therefore, at R > 100 kpc, not only are the peculiar velocities of the CGM gas less extreme, but they are also more uniformly distributed around the bulk motion of M31. At R < 100 kpc, the peculiar velocities of the CGM gas are more extreme and systematically redshifted relative to the bulk motion of M31.

The fourth property appears in Figure 7, where we compare v_M31 velocities of the M31 dwarfs and Si iii absorbers, which shows that overall the velocities of the satellites and the CGM gas do not follow each other. As noted in Section 3.3 (see also Table 4), some velocity components seen in absorption toward the QSOs are found with Δ_sep < 1 and have δv < v_esc,dwarf. However, the last two trends described above for the CGM gas are not observed for the dwarfs. First, contrasting with the M31 CGM gas, Figure 7 shows that both blue- and redshifted v_M31 velocities for the dwarfs are observed at any R and azimuth. Second, at R > 100 kpc and R ≤ 100 kpc, the average v_M31 velocities for the dwarfs are also similar, while for the CGM gas they are quite different (see above and Table 7). These differences strongly suggest that the velocity fields of the dwarfs and CGM gas are decoupled. We infer from this decoupling that (1) gas bound to satellites does not make a significant contribution to CGM gas observed in this way, and (2) the velocities of gas removed from satellites via tidal or ram pressure interactions, if it is present, become decoupled from the dwarf that brought it in (as one might expect from its definition as unbound to the satellites).

The fifth property is more readily apparent considering the average velocities shown in Figure 15, where, considering the CGM gas in different annuli, there is an apparent change in the sign of the average v_M31 velocities with on average a positive velocity at R < 200 kpc ($\langle {v}_{{\rm{M}}31}\rangle =+58.5\,\pm 43.5$ km s⁻¹ for Si iii), negative velocity at 200 < R < 300 kpc ($\langle {v}_{{\rm{M}}31}\rangle \,=\,-11.3\,\pm 24.2$ km s⁻¹ for Si iii), and again positive velocity at 300 < R < 400 kpc ($\langle {v}_{{\rm{M}}31}\rangle =+12.3\pm 21.8$ km s⁻¹ for Si iii). This is more evident with Si iii, where the sample of absorbers is larger, but taking the average velocities in the different annuli, the same pattern is observed for C ii, Si iii, Si iv, and C iv.

We emphasize again that the MS contamination does not really alter these properties and neither is the source of these properties. As shown in Figure 1 (see also Section 2.5), the MS contamination predominantly occurs in the region X < 0. There is no evidence that these properties change with Φ and in particular between the halves of the map X < 0 and X > 0 (see Section 2.5). Absorption occurring in the velocity range −50 ≲ v_M31 ≲ +150 km s⁻¹ is also not contaminated by the MS.

A key finding of Project AMIGA is the ubiquitous presence of metals throughout the CGM of M31. While the search for H i with $\mathrm{log}{N}_{{\rm{H}}{\rm\small{I}}}$ ≳ 17.5 in the CGM of M31 toward pointed radio observations has been unsuccessful in the current sample (Paper I and see Figure 1), the covering factor of Si iii (29 sightlines) is essentially 100% out to 1.2 ${R}_{\mathrm{vir}}$, while O vi associated with M31 is detected toward all 11 sightlines with FUSE data, all the way out to 1.9 ${R}_{\mathrm{vir}}$, the maximum radius of our survey (see Sections 4.4, 4.6). From the ionization range probed by Project AMIGA, we further show that Si iii and O vi are key probes of the diffuse gas (see Section 4.7). With information from Si ii, Si iii, and Si iv, we demonstrate that Si iii is the dominant ion in the ionizing energy range of 8–45 eV (see Section 4.7.1). The high covering factors of Si iii and O vi imply that these ions could be pervasively distributed, but it is quite possible also that small structures could lead to these large covering factors in projection.

The finding of ubiquitous metals in the CGM of M31 is a strong indication of ongoing and past gas outflows that ejected metals well beyond their formation site. Based on a specific SFR of SFR/M_⋆ = (5 ± 1) × 10⁻¹² yr⁻¹ (using the stellar mass M_⋆ and SFR from Geehan et al. 2006; Kang et al. 2009), M31 is not currently in an active star-forming episode. In fact, Williams et al. (2017) show that the bulk of star formation occurred in the first ∼6 billion years and the last strong episode happened over ∼2 billion years ago (see also Figure 6 in Telford et al. 2019 for a metal production model of M31). For a typical outflow with gas moving at 100 km s⁻¹, the gas would have traveled about 200 kpc in 2 Gyr, i.e., this is the maximum distance at which the last burst of star formation could have impacted the CGM of M31. Hence, many of the metals seen in the CGM of M31 must have also been ejected by previous star-forming episodes and/or stripped from its dwarfs and more massive companions. However, the most distant metals, those beyond ${R}_{\mathrm{vir}}$ (and especially given that they are predominantly found along a preferred direction from M31; see Section 4.9), may be supplied from elsewhere. For example, they could have originated in a Local Group medium bearing metals lost from other galaxies (e.g., the MW or dwarfs) or in a pre-polluted IGM filament.

In Section 4.8, we estimate that the mass of metals M^cool_Z = (2.0 ± 0.5) × 10⁷ M_⊙ within ${R}_{\mathrm{vir}}$ for the predominantly photoionized gas probed by Si ii, Si iii, and Si iv. For the gas probed by O vi, we find that ${M}_{{\rm{Z}}}^{\mathrm{warm}}\gt 4.4{M}_{{\rm{Z}}}^{\mathrm{cool}}\gtrsim 9\,\times {10}^{7}$ M_⊙ within ${R}_{\mathrm{vir}}$ (this is a lower limit because the ionization fraction of O vi is an upper limit; see Section 4.8). The sum of these two phases yields a lower limit to the CGM metal mass because the hotter phase probed by X-rays and metals bound in dust is not included. If the hot baryon mass of M31 is not too different from that estimated for the MW (see Section 4.8), then we expect ${M}_{{\rm{Z}}}^{\mathrm{hot}}\approx {M}_{{\rm{Z}}}^{\mathrm{warm}}$. The CLUES simulation of the Local Group estimates that the mass of the hot (>10⁵ K) gas is a factor of 3 larger than the cooler (<10⁵ K) gas (Nuza et al. 2014). The CGM dust mass remains quite uncertain but could be at the level of 5 × 10⁷ M_⊙ according to estimates around 0.1–1L* galaxies (Ménard et al. 2010; Peeples et al. 2014; Peek et al. 2015). Hence, a plausible lower limit on the total metal mass of the CGM of M31 out to ${R}_{\mathrm{vir}}$ is ${M}_{{\rm{Z}}}^{{\rm{C}}{\rm{G}}{\rm{M}}}\,\gtrsim 2.5\times {10}^{8}$ M_⊙, a factor of ∼2 times larger than the cool+warm gas-phase metal mass.

The stellar mass of M31 is (1.5 ± 0.2) × 10¹¹ M_⊙ (e.g., Williams et al. 2017). Using this result, Telford et al. (2019) estimated that the current metal mass in stars is 3.9 × 10⁸ M_⊙, i.e., about the same amount that is found in the entire CGM of M31 up to ${R}_{\mathrm{vir}}$. Telford et al. (2019) also estimated the metal mass of the gas in the disk of M31 to be around (0.8–3.2) × 10⁷ M_⊙, while Draine et al. (2014) estimated the dust mass in the disk to be around 5.4 × 10⁷ M_⊙, yielding a total metal mass in the disk of M31 of about ${M}_{{\rm{Z}}}^{{\rm{d}}{\rm{i}}{\rm{s}}{\rm{k}}}\simeq 5\,\times {10}^{8}$ M_⊙. Therefore, M31 has in its CGM within ${R}_{\mathrm{vir}}$ at least 50% of the present-day metal mass in its disk. As we show in Sections 4.4 and 4.8 and discuss above, metals are also found beyond ${R}_{\mathrm{vir}}$, especially in the more highly ionized phase traced by O vi (and even higher ions). These metals may have come from M31 or from other galaxies in the Local Group such as the MW or dwarf galaxies.

The Project AMIGA experiment is quite different from most of the surveys of the CGM of galaxies done so far. Outside the local universe, surveys of the CGM of galaxies involve assembling samples of CGM gas in aggregate by using one sightline per galaxy (see Section 1), in some nearby cases up to three to four sightlines (e.g., Bowen et al. 2016; Keeney et al. 2017). By assembling a sizable sample of absorbers associated with galaxies in a particular subpopulation (e.g., L*, sub-L*, passive or star-forming galaxies), one can then assess how the column densities change with radii around that kind of galaxy, and from this estimate average surface densities, mass budgets, etc., can then be evaluated. By contrast, Project AMIGA has assembled almost as many sightlines surrounding M31 as COS-Halos had for its full sample of 44 galaxies. We can now make a direct comparison between these two types of experiments. For this comparison, we use the COS-Halos survey of 0.3 < L/L* < 2 galaxies at z ≃ 0.2, which selected galaxies within about 160 kpc from the sightline (Tumlinson et al. 2011, 2013; Werk et al. 2013, 2014). The full mass range of the COS-Halos galaxies is quite large, 11.5 ≲ log M₂₀₀ ≲ 13.7, but most of the star-forming galaxies are in the range $11.5\lesssim \mathrm{log}{M}_{200}\lesssim 12.5$, and most of the passive quiescent galaxies have $13.0\lesssim \mathrm{log}{M}_{200}\lesssim 13.7$. As a reminder, M31 has $\mathrm{log}{M}_{200}=12.1$ (see Section 1).

5.2.1. Column Densities of Si versus R

In Figure 16, we show the total column densities of Si as a function of R/R₂₀₀ for the COS-Halos galaxies and M31. For the COS-Halos survey, each data point corresponds to an absorber at some impact parameter from a galaxy, while for M31, each data point is an absorber probing the CGM at a different impact parameter from the same galaxy. For COS-Halos, we consider two cases: (1) the column densities of Si estimated in a similar fashion as those for M31, and (2) the column densities of Si estimated from photoionization modeling. For case 2 we use the results from Werk et al. (2014; see also Prochaska et al. 2017b), which were used to determine the metal mass of the CGM of the COS-Halos galaxies in Peeples et al. (2014). For case 1, we use the results from Werk et al. (2013) and follow the procedure in Section 4.5 to estimate N_Si from the column densities of Si ii, Si iii, and Si iv. We require that all three Si ions are available, except in the cases where there are only lower limits for two ions (typically Si ii and Si iii) since in that case the resulting column density is a lower limit that encompasses any missing column density from the remaining Si ion (typically Si iv). The sample size in case 1 is 35, while in case 2 it is 33, with some overlap between the two subsamples. In Figure 16, we also show the modeled column density of Si as a function of R, in gray for COS-Halos and blue for M31 (we show the adopted GP model; see Section 4.5).

Zoom In Zoom Out Reset image size

A striking difference between the COS-Halos and M31 data immediately apparent from Figure 16 is that a large number of very high Si column densities are observed in COS-Halos at R/R₂₀₀ ≲ 0.3 that are absent in Project AMIGA. The reason that the COS-Halos Si column densities have higher lower limits than those of M31 is because the weak transitions of Si ii are saturated in COS-Halos, a situation not observed in M31 toward any of the sightlines—the lower limits of Si arise only because Si iii is saturated. These high Si column densities also correspond to the very strong N_{H i} ($\mathrm{log}{N}_{{\rm{H}}{\rm\small{I}}}$ ≳ 18) absorbers observed in COS-Halos, but again not in M31 (see Figure 5 and Section 4.1 in Paper I). Beam dilution could have affected somewhat the comparison between H i column densities from COS-Halos (H i absorption) and M31 (H i emission); for the metal ions, this is not an issue. Therefore, the higher frequency of saturated Si ii transitions in COS-Halos compared to M31 is a real effect, not an artifact.

Besides this difference, the estimated Si column densities from the observations in the COS-Halos and Project AMIGA surveys are distributed with a similar scatter at larger impact parameters (R/R₂₀₀ ≳ 0.4), where the gas is more ionized (see top panel in Figure 16). The photoionization-modeled COS-Halos Si displayed in the bottom panel has some higher values than observed in the top panel, but in the impact parameter region 0.4 ≲ R/R₂₀₀ ≲ 0.8 where they are observed, there are also several lower limits. Beyond R > 0.9R₂₀₀, there are no COS-Halos observations (owing to the design of the survey). The extrapolated model to the COS-Halos observations shown in gray in Figure 16 is a factor of 2–4 higher than the models of the Project AMIGA data shown in blue depending on R/R₂₀₀.

A likely explanation for the higher column density absorbers is that some of these COS-Halos absorbers could be fully or partly associated with a galaxy near to the initially targeted COS-Halos galaxies. In this case the gas would be more likely to show high columns of neutral and weakly ionized gas. Indeed, while the COS-Halos galaxies were selected to have no bright companions, that selection did not preclude fainter nearby companions such as dwarf satellites (see Tumlinson et al. 2013). Galaxy observation follow-up by Werk et al. (2012) found several L > 0.1L* galaxies within 160 kpc of the targeted COS-Halos galaxy. Comparing the results from other surveys of galaxies/absorbers at low redshift (Bowen et al. 2002; Stocke et al. 2013), Bregman et al. (2018) also noted a higher preponderance of high H i column density absorbers in the COS-Halos survey. However, the higher COS-Halos column densities at large radii could also be an effect of evolution in the typical CGM, as COS-Halos probed a slightly higher cosmological redshift. It is also possible that the M31 CGM is less rich in neutral gas at these radii than the typical L* galaxy at z ∼ 0.2, because of its star formation history or environment.

5.2.2. CGM Mass Comparison

A key discovery from Project AMIGA is that the properties of the CGM of M31 change with R. This is reminiscent of our earlier survey (LHW15), but the increase in the size sample has transformed some of the tentative results of our earlier survey into robust findings. In particular, the radius around R ∼ 100–150 kpc appears critical in view of several properties changing near this threshold radius:

• 1.  
  For any ions, the frequency of strong absorption is larger at R ≲ 100–150 kpc than at larger R.
• 2.  
  The column densities of Si and C ions change by a factor >5–10 between about 25 and 150 kpc, while they change only by a factor ≲2 between 150 and 300 kpc.
• 3.  
  The detection rate of singly ionized species (C ii, Si ii) is close to 100% at R < 150 kpc but sharply decreases beyond (see Figure 9), and therefore the gas has a more complex gas-phase structure at R < 150 kpc.
• 4.  
  The peculiar velocities of the CGM gas are more extreme and systematically redshifted relative to the bulk motion of M31 at R ≲ 100 kpc, while at R ≳ 100 kpc the peculiar velocities of the CGM gas are less extreme and more uniformly distributed around the bulk motion of M31.

There are also two other significant regions: (1) beyond R₂₀₀ ≃ 230 kpc the gas is becoming more ionized and more highly ionized than at lower R (e.g., there is a near total absence of Si ii absorption beyond R₂₀₀ [see Figure 9], or a higher C ii/C iv ratio on average at R ≳ R₂₀₀ than at lower R [see Section 4.7.2]), and (2) beyond 1.1 ${R}_{\mathrm{vir}}$ the gas is not detected in all the directions away from M31, as it is at smaller radii, but only in a cone near the southern projected major axis and about 52° west of the X = 0 kpc axis (see Section 4.9).

The overall picture that can be drawn out from these properties is that the inner regions of the CGM of M31 are more dynamic and complex, while the more diffuse regions at R ≳ 0.5 ${R}_{\mathrm{vir}}$ are more static and "simpler." Zoom-in cosmological simulations capture in more detail and more accurately the structures of the CGM than large-scale cosmological simulations thanks to their higher mass and spatial resolution. Below we use several results from zoom simulations to gain some insights on these observed changes with R. However, the results laid out in Section 4 also now provide a new test bed for zoom simulations, so that not only qualitative but also quantitative comparison can be undertaken. We note that most of the zoom simulations discussed here have only a single massive halo. However, according to the ELVIS simulations of Local Group analogs (Garrison-Kimmel et al. 2014), there should be no major difference at least within about ${R}_{\mathrm{vir}}$ for the distribution of the gas between isolated and paired galaxies.

5.3.1. Visualization and Origins of the CGM Variation

5.3.2. Quantitative Comparison in the CGM Variation between Observations and Simulations

Two simulations of M31-like galaxies in different environments at widely separated epochs show some similarity with some of the observed trends in the CGM of M31. We now take one step further by quantitatively comparing the column density variation of the different ions as a function of R in three different zoom-in cosmological simulations, two being led by members of the Project AMIGA team (FIRE and FOGGIE collaborations), and a zoom-in simulation from the Evolution and Assembly of GaLaxies and their Environments (EAGLE) simulation project (Crain et al. 2015; Schaye et al. 2015; Oppenheimer et al. 2018a).

(1) Comparison with FIRE-2 Zoom Simulations

We first compare our observations with column densities modeled using cosmological zoom-in simulations from the FIRE project.³¹ Details of the simulation setup and CGM modeling methods are presented in Ji et al. (2020). Briefly, the outputs analyzed here are FIRE-2 simulations evolved with the GIZMO code using the meshless finite mass (MFM) solver (Hopkins 2015). The simulations include a detailed model for stellar feedback including core-collapse and Type Ia supernovae, stellar winds from OB and AGB stars, photoionization, and radiation pressure (for details, see Hopkins et al. 2018). We focus on the "m12i" FIRE halo, which has a mass M_vir ≈ 1.2M₂₀₀ ≈ 1.2 × 10¹² M_⊙ at z = 0, which is comparable to the halo mass of M31. However, neither the SFR history nor the present-day SFR is similar. The "m12i" FIRE halo has a factor of 10–12 higher SFR (see Figure 3 in Hopkins et al. 2020) than the present-day SFR of M31 of 0.5 M_⊙ yr⁻¹ (e.g., Kang et al. 2009). We compare Project AMIGA to FIRE-2 simulations with two different sets of physical ingredients. The "MHD" run includes magnetic fields, anisotropic thermal conduction, and viscosity, and the "CR" run includes all these processes plus the "full physics" treatment of stellar cosmic rays (CRs). The CR simulation assumes a diffusion coefficient ${\kappa }_{| | }=3\,\times \,{10}^{29}$ cm⁻² s⁻¹, which was calibrated to be consistent with observational constraints from γ−ray emission of the MW and some other nearby galaxies (Chan et al. 2019; Hopkins et al. 2020). Ji et al. (2020) showed that CRs can potentially provide a large or even dominant nonthermal fraction of the total pressure support in the CGM of low-redshift ∼L* galaxies. As a result, in the fiducial CR run analyzed here, the volume-filling CGM is much cooler (∼10⁴–10⁵ K) and is thus photoionized in regions where in the run without CRs it is filled with hot gas that is more collisionally ionized.

The column densities are generated as discussed in Ji et al. (2020) using the Trident code (Hummels et al. 2017) to model the ion populations and density projections. For the photoionization modeling, a hybrid treatment combining the FG09 (Faucher-Giguère et al. 2009) and HM12 (Haardt & Madau 2012) UV background models is used.³²

In Figure 17, we compare the ion column densities from FIRE-2 simulations with observationally derived total column densities around M31 as a function of R/R₂₀₀. The green and orange curves show the median simulated column densities for the MHD and CR runs, respectively, while the shaded regions show the full range of columns for all sightlines at a given impact parameter (the lowest values are truncated to match the scales that are adequate for the observations; see Ji et al. 2020 for the full range of values). The CR run produces higher column densities and better agreement with observations than the MHD run for all ions presented. The much higher column densities of low/intermediate ions (C ii, Si ii, Si iii, and Si iv) in the CR run are due to the more volume-filling and uniform cool phase, which produces higher median values of ion column densities and smaller variations across different sightlines. In contrast, in the MHD run the cool phase is pressure confined by the hot phase to compact and dense regions, leading to smaller median columns but larger scatter for the low and intermediate ions. We note, however, that even in the CR runs the predicted column densities are lower than observations at the larger impact parameters R ≳ 0.5R₂₀₀. This might be due to insufficient resolution to resolve fine-scale structure in outer halos, or it may indicate that feedback effects are more important at large radii than in the present simulations. This difference is quite notable owing to the fact that the star formation of the "m12i" galaxy has been continuous with an SFR in the range of 5–20 M_⊙ yr⁻¹ (Hopkins et al. 2020) over the past ∼8 billion years, while M31 had only a continuous SFR around 6–8 M_⊙ yr⁻¹ over its first 5 billion years, while over the past 8 billion years it had only two short bursts of star formation about 4 and 2 billion years ago (Williams et al. 2017). While there are some clear discrepancies, the simulations also follow some similar trends to those seen in the empirical results detailed in Section 4: (1) the simulated column densities of the low ions decrease more rapidly with R than the high ions; (2) O vi is observed beyond 1.7R₂₀₀, where there is no substantial amount of low/intermediate ions; and (3) a larger scatter is observed in the column densities of the low and intermediate ions than in O vi.

Zoom In Zoom Out Reset image size

In the FIRE-2 simulations, both collisional ionization and photoionization can contribute significantly to the simulated O vi columns, typically with an increasing contribution from photoionization with increasing impact parameter, driven by decreasing gas densities. In the MHD run, most of the O vi in the inner halo (R ≲ 0.5R₂₀₀) is produced by collisional ionization, but photoionization can dominate at larger impact parameters. In the CR run, collisional ionization and photoionization contribute comparably to the O vi mass at radii 50 < R < 200 kpc (Ji et al. 2020). The actual origins of the CGM in terms of gas flows in FIRE-2 simulations without magnetic fields or CRs were analyzed in Hafen et al. (2020), although the results are expected to be similar for simulations with MHD only. In these simulations, O vi exists as part of a well-mixed hot halo, with contributions from all the primary channels of CGM mass growth: IGM accretion, wind, and contributions from satellite halos (reminiscent of the Eris2 simulations; see above and Shen et al. 2013). The metals responsible for O vi absorption originate primarily in winds, but IGM accretion may contribute a large fraction of total gas mass traced by O vi since the halo is well mixed and IGM accretion contributes ≳60% of the total CGM mass. In the simulations, the hot halo gas persists in the CGM for billions of years, and the gas that leaves the CGM does so primarily by accreting onto the central galaxy (Hafen et al. 2019).

(2) Comparison with FOGGIE Simulations

We also compare the observed total column densities to the MW-mass "Tempest" (M₂₀₀ ≈ 4.2 × 10¹¹ M_⊙) halo from the FOGGIE simulations,³³ which has a halo mass of M₂₀₀ ≈ 4.2 × 10¹¹ M_⊙ (Peeples et al. 2019). Again we employ the software Trident (Hummels et al. 2017) to generate the ion populations and density projections, limiting potential sources of discrepancy. We use the z = 0 output (see Zheng et al. 2020 for simulation details), but because of the size difference between M31 and the Tempest galaxy, we again scale all distances by R₂₀₀ (R₂₀₀ = 159 kpc for the simulated halo compared to 230 kpc for M31). The only "feedback" included in this FOGGIE run is thermal explosion-driven supernova outflows. While this feedback is limited in scope compared to FIRE, FOGGIE achieves higher mass resolution than FIRE-2 by using a "forced refinement" scheme that applies a fixed computational cell size of ∼381 h⁻¹ pc within a moving cube centered on the galaxy that is ∼200 h⁻¹ ckpc on a side. This refinement scheme enforces constant spatial resolution on the CGM, resulting in a variable and very small mass resolution in the low-density gas, with typical cell masses of (≲1–100 M_⊙). The individual small-scale structures that contribute to the observed absorption profiles can therefore be resolved. These small-scale structures that become only apparent in high-resolution simulations are hosts to a significant amount of cool gas, enhancing the column densities in especially the low ionization state of the gas (Peeples et al. 2019; Corlies et al. 2020; see also Hummels et al. 2019; Rhodin et al. 2019; van de Voort et al. 2019).

As for the FIRE-2 simulations, we compare the Project AMIGA total column densities to FOGGIE because in the simulation we do not (yet) separate individual components, but look at the projected column densities through the halo. We note that the CGM is not necessarily self-similar, so some differences between the simulation predictions and M31 observations at rescaled impact parameter could be due to the halo mass difference. This is especially so since the halo mass range M_h ≈ 3 × 10¹¹–10¹² M_⊙ corresponds to the expected transition between cold and hot accretion (e.g., Birnboim & Dekel 2003; Keres et al. 2003; Faucher-Giguère et al. 2011; Stern et al. 2020).

In Figure 18, we compare the simulated and observed column densities for each ion probed by our survey. The pink and green shaded areas are the data points from the simulation (with and without satellite contribution, respectively) and show the total column density in projection through the halo. The scatter in the simulated data points comes from variation in the structures along the mock sightline, and most of the scatter is in fact below 10¹¹ cm⁻². The peaks in the column densities are due to small satellites in the halo, which enhance primarily the low-ion column densities. We show the green points to highlight the difference between the mock column densities with and without satellites. For the high-ionization lines the difference is negligible, while the difference in the low ions is significant.

Zoom In Zoom Out Reset image size

Overall, the metal-line column densities are systematically lower than in the observations at any R. Only at R ≲ 0.3R₂₀₀ is there some overlap for the singly ionized species between the FOGGIE simulation and observations. However, the discrepancy is particularly striking for Si iii and the high ions. This can be understood by the current feedback implementation in FOGGIE, which does not expel enough metals from the stellar disk into the CGM (Hamilton-Campos et al. 2020) to be consistent with known galactic metal budgets (Peeples et al. 2014). This effect is expected to be stronger for the high ions than the low ions, due to the additional heating and ionization of the CGM that would be expected from stronger feedback, and indeed the discrepancy between the simulation and observations is larger for the high ions (and Si iii) than for the singly ionized species. However, while the absolute scale of the column densities is off, there are also some similarities between the simulation and observations in the behavior of the relative scale of the column density profiles with R: (1) the column densities of the low ions drop more rapidly with R than the high ions; (2) despite the inadequate feedback in the current simulations, the O vi-bearing gas (and C iv to a lesser extent) is observed well beyond R₂₀₀; and (3) a larger scatter is observed in the column densities of the low and intermediate ions than in O vi. It is striking that the overall slope of the O vi profile resembles the observations but at significantly lower absolute column density. In the FOGGIE simulation, the low ions tracing mainly dense, cool gas are preferentially found in the disk or satellites, while the hotter gas traced by the higher ions is more homogeneously distributed in the halo.

(3) Comparison with EAGLE Simulations

Finally, we compare our results with the EAGLE zoom-in simulations (EAGLE Recal-L025N0752 high-resolution volume) discussed in length in Oppenheimer et al. (2018a). The EAGLE simulations have successfully reproduced a variety of galaxy observables (e.g., Crain et al. 2015; Schaye et al. 2015) and achieved "broad but imperfect" agreement with some of the extant CGM observations (e.g., Turner et al. 2016; Oppenheimer et al. 2018a; Rahmati & Oppenheimer 2018; Lehner et al. 2019; Wotta et al. 2019).

Oppenheimer et al. (2018a) aimed to directly study the multiphase CGM traced by low metal ions and to compare with the COS-Halos survey (see Section 5.2). As such, they explored the circumgalactic metal content traced by the same ions explored in Project AMIGA in the CGM galaxies with masses that comprise that of M31. Overall Oppenheimer et al. find agreement between the simulated and COS-Halos samples for Si ii, Si iii, Si iv, and C ii within a factor of two or so and larger disagreement with O vi, where the column density is systematically lower. With Project AMIGA, we can directly compare the results with one of the EAGLE galaxies that has a mass very close to M31 and also compare the column densities beyond 160 kpc, the maximum radius of the COS-Halos survey (Tumlinson et al. 2013; Werk et al. 2013). We refer the reader to Oppenheimer et al. (2016), Rahmati & Oppenheimer (2018), and Oppenheimer et al. (2018a) for more details on the EAGLE zoom-in simulations, which employ the SpecExBin code (Oppenheimer & Davé 2006) for ion modeling and column density projections. We also refer the reader to Figure 1 in Oppenheimer et al. (2018a), where in the middle column they show the column density map for galaxy halo mass of $\mathrm{log}{M}_{200}=12.2$ at z ≃ 0.2, which qualitatively shows similar trends described in Section 5.3.1.

In Figure 19, we compare the EAGLE and observed column densities as a function of the impact parameter out to ${R}_{\mathrm{vir}}$. As in the previous two figures, the blue and gray circles are detections and nondetections in the halo of M31, respectively. The green curve in each panel represents the mean column density for each ion as a function of the impact parameter for the EAGLE galaxy with $\mathrm{log}{M}_{200}=12.1$ at z = 0. In contrast to FIRE-2 or FOGGIE simulations, the EAGLE simulations appear to produce a better agreement between N and R for low and intermediate ions (Si ii, Si iii, Si iv), and C iv out to larger impact parameters. However, as already noted in Oppenheimer et al. (2018a), this agreement is offset by producing too much column density for the low and intermediate ions at small impact parameters (see, e.g., Si ii, which is not affected by lower limits and is clearly overproduced at R ≲ 80 kpc). The flat profile of O vi, with very little dependence on R, is similar to the observations and other models, but overall the EAGLE O vi column densities are a factor of 0.2–0.6 dex smaller than observed. Oppenheimer et al. (2018a) (and also Oppenheimer et al. 2016) already noted that issue from their comparison with the COS-Halos galaxies (see also Section 5.2), requiring additional source(s) of ionization for the O vi such as AGN flickering (Oppenheimer & Schaye 2013b; Oppenheimer et al. 2018b) or possibly CRs as shown for the FIRE-2 simulations (see Ji et al. 2020 and above). While the results are shown only to ${R}_{\mathrm{vir}}$, as in the other simulations and M31, O vi is also observed well beyond ${R}_{\mathrm{vir}}$ in the EAGLE simulations (see Figure 1 in Oppenheimer et al. 2018a).

Zoom In Zoom Out Reset image size

5.3.3. Insights from the Observation/Simulation Comparison

Based on our findings, it is likely that the MW has not only an extended hot CGM (Gupta et al. 2014, 2017; see also Henley & Shelton 2010; Miller & Bregman 2013) but also an extended CGM of cool (Si ii, Si iii, Si iv) and warm–hot (C iv, O vi) gas that extends all the way to about 300 kpc (${R}_{\mathrm{vir}}$), and even farther away for the O vi. In fact, the MW and M31 O vi CGMs most likely already overlap as can be seen, e.g., in the CLUES simulations of the Local Group (Nuza et al. 2014) since the distance between M31 and MW is only 752 kpc.

The multiple gas-phase MW halo has largely been studied using HVCs because the velocities of these clouds are high enough to separate them from the disk absorption (e.g., Wakker & van Woerden 1997; Putman et al. 2012; Richter et al. 2017). However, we emphasize that the large majority of HVCs, including the predominantly ionized HVCs, are not at hundreds of kiloparsecs from the MW, but most of them are within 15–20 kpc from the Sun (e.g., Wakker 2001; Thom et al. 2008; Wakker et al. 2008; Lehner & Howk 2011; Lehner et al. 2012), i.e., in a radius not even explored by Project AMIGA and many other surveys of the galaxy CGM at higher redshifts (e.g., Werk et al. 2013; Liang & Chen 2014; Borthakur et al. 2016; Burchett et al. 2016). Only the MS allows us to probe the interaction between the MW and the Magellanic Clouds in the CGM of the MW out to about 50–100 kpc (e.g., D'Onghia & Fox 2016). The results from Project AMIGA strongly suggest that the CGM of the MW is hidden in the low-velocity absorption arising from its disk (see also Zheng et al. 2015). To complicate the matter, the column densities of the low ions, intermediate ions, and C iv drop substantially beyond 100–150 kpc (see, e.g., Figures 8 and 10). Owing to its strength and little dependence on R, O vi is among the best ultraviolet diagnostics of the extended CGM (see also the recent FOGGIE simulation results by Zheng et al. 2020).