Impact of the PSR J0740+6620 radius constraint on the properties of high-density matter

X-ray pulse profile modeling of PSR J0740+6620, the most massive known pulsar, with data from the NICER and XMM-Newton observatories recently led to a measurement of its radius. We investigate this measurement's implications for the neutron star equation of state (EoS), employing a nonparametric EoS model based on Gaussian processes and combining information from other x-ray, radio and gravitational-wave observations of neutron stars. Our analysis mildly disfavors EoSs that support a disconnected hybrid star branch in the mass-radius relation, a proxy for strong phase transitions, with a Bayes factor of $6.9$. For such EoSs, the transition mass from the hadronic to the hybrid branch is constrained to lie outside ($1,2$) $M_{\odot}$. We also find that the conformal sound-speed bound is violated inside neutron star cores, which implies that the core matter is strongly interacting. The squared sound speed reaches a maximum of $0.75^{+0.25}_{-0.24}\, c^2$ at $3.60^{+2.25}_{-1.89}$ times nuclear saturation density at 90% credibility. Since all but the gravitational-wave observations prefer a relatively stiff EoS, PSR J0740+6620's central density is only $3.57^{+1.3}_{-1.3}$ times nuclear saturation, limiting the density range probed by observations of cold, nonrotating neutron stars in $\beta$-equilibrium.


I. INTRODUCTION
The properties and composition of matter at the highest densities achieved in neutron star (NS) cores remain uncertain [1][2][3][4]. The main observational constraints on the equation of state (EoS) of NS matter at densities 3 ρ nuc , where ρ nuc = 2.8 × 10 14 g/cm 3 is the nuclear saturation density, come from radio measurements of the masses of the heaviest known pulsars [5][6][7][8][9]. These observations place the maximum nonspinning NS mass above 2 M , which limits the softness of the high-density EoS and tends to decrease the probability of exotic degrees of freedom that reduce the pressure within NS matter.
Recently, a measurement of the radius of the 2.08 M pulsar PSR J0740+6620 [8,9] using x-ray data from NICER and XMM-Newton was reported by two independent analyses [55,56]. This radius constraint presents a rare glimpse of the properties of the most massive NSs, and a golden opportunity to obtain observational information about the maximum NS mass, M max , as well as potential phase transitions in NS cores. In the context of the preferred hotspot model in each analysis, [55] finds 13.7 +2.6 −1.5 km and [56] obtains 12.4 +1.3 −1.0 km for J0740+6620's radius (medians and symmetric 68% credible intervals). For context, the inference reported in [50] predicts the radius of 2.08 M NSs to be 12.08 +0. 79 −0.98 km at the 68% confidence level.
Observations of the most massive NSs, such as J0740+6620, have important implications beyond the EoS. They inform the NS mass distribution [57][58][59][60][61][62][63][64], the classification of the heaviest NS candidates observed with GWs [30], our understanding of the proposed mass gap between NSs and black holes [65,66], and the characteristics of NS merger remnants [67] that influence electromagnetic counterpart emission [68,69]. The properties of the high-density EoS are also connected to the properties at other density scales through correlations shaped by causality considerations [70,71].
To determine the implications of J0740+6620's radius measurement for NS matter, we employ a nonparametric model for the NS EoS based on Gaussian processes (GPs), which offers us the flexibility of an analysis that (i) is not tightly linked to specific nuclear models, (ii) can account for phase transitions, including strong firstorder phase transitions that result in disconnected stable branches in the mass-radius relation, and (iii) is not subject to the systematic errors that arise with parametrized EoS families described by a finite set of parameters. Additionally, the nonparametric EoS model allows us to probe a wider range of intra-density correlations in the EoS than parametric models, something especially relevant for the current data set, which targets a wide range of NS densities [72].
We find that the new J0740+6620 observation pushes the inferred radii and maximum mass for NSs to larger values: we obtain R 1.4 = 12.56 +1.00 −1.07 km for the radius of a 1.4 M NS and M max = 2.21 +0. 31 −0.21 M for the maximum nonrotating NS mass (we quote medians and 90% highest-probability-density credible regions unless otherwise noted). Despite significant statistical uncertainties, the inferred NS radii are consistent with being equal over a broad mass range, with a radius difference of ∆R ≡ R 2.0 − R 1.4 = −0.12 +0. 83 −0.85 km between 2.0 M and 1.4 M NSs. This conclusion rules out a large reduction in the radius for massive NSs, a feature that is sometimes characteristic of strong phase transitions in the mass regime of typical NSs. Further, we find that EoSs with at least one disconnected hybrid star branch in their mass-radius relation are disfavored compared to those with a single stable branch by a factor of approximately 6.9. This supports the current consensus that all dense-matter observations can be accommodated by a standard hadronic EoS, although the possibility of a phase transition remains viable; only the strongest firstorder phase transitions produce more than one stable sequence of compact stars. If, on the other hand, the massradius relation has multiple stable branches, the heaviest star on the first stable branch is either 1 M or ∼ 2 M . Our results disfavor a transition mass in the intermedi-ate mass regime, suggesting that either all NSs observed to date contain exotic cores, or virtually all are purely hadronic.
We also find support for a violation of the conjectured conformal bound on the sound speed c s in NS matter, c 2 s ≤ c 2 /3 [73][74][75], where c is the speed of light. Such a violation indicates that the sound speed does not rise monotonically to the perturbative QCD limit (c 2 s → c 2 /3) at asymptotically high densities [11] and signals the presence of strongly interacting matter in NS cores [50]. The stiff high-density EoS required by the massive pulsar observations already put the conformal bound in jeopardy [11,59,[75][76][77], but the softer low-density behavior favored by GWs and the NICER radius measurements help reach a Bayes factor of 1000 ± 340 (mean and standard deviation from Monte Carlo uncertainty), securely in favor of a violation. We infer that c 2 s reaches a maximum of 0.75 +0 Our results are comparable to other analyses of the new J0740+6620 radius measurement. Reference [55] examined the pressure-density relation, the NS radius, and M max using the same set of observational data as we do but did not comment on the possibility of phase transitions in the EoS. They adopted three different models for the EoS(including a simple, more restricted implementation of a GP) which each yielded different but overlapping constraints on the EoS. EoS models informed by chiral effective field theory (χEFT) at low densities and GW170817's electromagnetic counterpart were considered in Refs. [78,79]; the latter analysis also disfavors EoSs with strong first-order phase transitions, while the former compared two parametric EoS models, finding some model-dependence in their results. A hybrid nuclear parameterization and piecewise polytrope EoS model was employed in [19], which also accounted for the recent PREX-II measurement of the neutron skin of 208 Pb [17]. Compared to these studies, our less restrictive treatment of the EoS model broadly results in both qualitative and quantitative agreement. Nonetheless, it allows us more freedom to investigate the consequences of the J0740+6620 radius measurement for NS matter microphysics, including phase transitions, the conformal sound-speed bound, and the inferred stiffness of the EoS.
The remainder of the paper describes the details and results of our analysis. In Sec. II we briefly describe the methodology we employ as well as the relevant data sets. In Sec. III we present the results of our inference for macroscopic NS properties. In Sec. IV we discuss the constraints that can be placed on microscopic EoS properties in terms of the sound speed in NS matter and phase transitions. We conclude and discuss other studies of J0740+6620's EoS implications in the literature in more detail in Sec. V.

II. EQUATION OF STATE INFERENCE
Our analysis methodology closely follows that of [50]; here we briefly summarize the main features and discuss the updated treatment of J0740+6620.

A. Hierarchical inference
In order to combine information from multiple data sets that include statistical uncertainties, we use hierarchical inference [80]. The relevant formalism and equations are described in detail in Sec. III B of [50]. The marginal likelihood of each observation (for example, a GW tidal measurement) for a given EoS model is obtained by marginalizing over the relevant parameters for individual events (in the GW case, the binary masses and tidal parameters) assuming some prior distribution (i.e., population model) for the nuisance parameters (in the GW case, the binary masses). Similar to [50], we assume a fixed population for all observations given the relatively low number of observations to date. This simplification also makes the EoS likelihood independent of selection effects [81]. However, as the size of each data set increases (for example through the observation of additional GW signals), we will need to simultaneously fit the population in order to avoid biases in the EoS inference [62,82].
In the absence of knowledge of the true compact object mass distribution, we choose a uniform population model that extends beyond the maximum mass of all EoSs we consider. For a given EoS model, we further assume that all objects with m ≤ M max are NSs. In other words, we assume that it is the EoS, and not the astrophysical formation mechanism, that limits the maximum NS mass. Then, for observations of objects known a priori to be NSs (such as J0740+6620, but unlike the components of the GW events), the normalization of the mass prior mildly penalizes EoSs that predict a maximum mass larger than all observed NS masses. This is an Occam penalty that favors EoSs that occupy a smaller prior volume and do not predict unobserved data in the form of very massive NSs, all else being equal. If, instead, we truncated the NS mass distribution below M maxe.g., because we had knowledge of an astrophysical process that limits the maximum NS mass-all EoSs with M max greater than the largest population mass would be assigned equal marginal likelihood. However, such a choice would have to be accompanied by an arbitrary choice of the truncation mass, given our lack of prior knowledge about the upper limit of the astrophysical NS mass distribution.
The distinction between these scenarios is important for any analysis of J0740+6620, given its high mass. We employ a uniform mass distribution with a lower limit of 0.5 M ; hence, an EoS with M max = 3 M is disfavored in our inference compared to an EoS with M max = 2.5 M by a factor of (3 − 0.5)/(2.5 − 0.5) = 1.25. In the results presented in later sections, for example Fig. 2, this contributes to the fact that the tail of our M max posterior is slightly tighter than the prior. More details and a quantitative assessment of the effect of the mass prior are given in the Appendix.

B. Nonparametric EoS model
The procedure outlined above requires a model that describes the NS EoS and can be used to compute all relevant macroscopic NS properties, such as masses, radii, and tidal deformabilities [83][84][85]. Following [50], we use a nonparametric EoS model constructed through GPs conditioned on existing dense-matter EoS models available in the literature; see [50,86,87] for more details. While the GP never assumes a specific functional form for the EoS, unlike parametric analyses, it does assume probabilistic knowledge about correlations within the EoS. Each GP is constructed with different hyperparameters that specify a covariance kernel, which in turn controls the scale and strength of correlations between the sound speed at different pressures. The specific model employed here is described in detail in [87]. It is constructed as a mixture of ∼ 150 individual GPs with a broad set of hyperparameters, allowing us to probe a wide range of EoS models with different intra-density correlations.
The nuclear models on which the process is conditioned contain EoSs with different degrees of freedom, including purely hadronic, hyperonic, and quark models. We intentionally condition only loosely on these models, resulting in the process termed model-agnostic in [87]. 1 As a result, our EoS prior contains a large variety of EoS behavior, including phase transitions at different density scales and of different strengths; see for example Fig. 1 of [87]. This nonparametric approach offers two further advantages over more traditional parametric models [88][89][90]: it avoids (i) systematic errors and (ii) strong (and perhaps opaque) intra-density correlations [72] that arise from restricting the EoS to a specific functional form with finite parameters, which will inevitably not match the correct EoS.
To that end, Ref. [55] also employed a GP EoS model, citing the same benefits we point out here. However, [55] used a single GP with a single set of hyperparameters (compared to ∼ 150 GPs we consider) and chose those hyperparameters to approximate the variability observed within tabulated EoSs from the CompOSE database [91]. Therefore, the GP prior explored in [55] is more reminiscent of the model-informed prior considered in [86,87] than the model-agnostic prior considered here and in [14,18,50,86,87]. In fact, the hyperparameters used in [55] assume less variance (smaller σ) and stronger correlations between pressures (larger l) than any of the allowed hyperparameters within our hyperprior (see [87] for more details). Our results, therefore, intentionally explore broader ranges of possible EoS behavior and intradensity correlations than [55], particularly at high densities where the GP model in [55] forces the sound speed to approach the speed of light a priori. Their more closely tailored GP design may explain why Fig. 10 of [55] shows that their GP analysis leads to more stringent EoS constraints than parametric EoS inferences.
One notable difference compared to [50] is that J0740+6620 now appears in the list of both radio and x-ray observations. As described in [55,56], the measured mass of J0740+6620 is still dominated by the radio observations [9]. The most recent mass estimate of 2.08 +0.07 −0.07 M is slightly lower than the originally reported value of 2.14 +0. 10 −0.09 M [8] (68% confidence level), making it more consistent with other Galactic NS mass measurements [61]. To avoid double-counting, we include J0740+6620 either through its updated mass estimate in the radio list or through its mass and radius estimate in the x-ray list. The difference gives an estimate of the impact of the radius constraint alone on the NS EoS.
When treating J0740+6620 as either a radio or an xray observation, we explicitly account for the normalization of the mass prior in Eqs. (9) and (11) of [50], in accordance with our choice of fixed population model. 2 For the J0740+6620 x-ray data, we use either the NICER+XMM samples from [55,97] or the ST-U samples from [56,98]. Both sets of samples already incorporate the updated mass estimate from [9], though [55] inflates the uncertainty in this measurement by ±0.02 M out of concern for systematic uncertainties. For our analysis, we choose to revert back to the published result from [9] and remove 2 Unlike in [50], where we assumed the population of NICER targets ended at masses well below Mmax, we assume the population of NICER targets extends well beyond Mmax and include the proper normalization for the mass prior for both J0030+0451 and J0740+6620.
the additional uncertainty of 0.04 M . In practice, we use the posterior samples from [55,97] but weight each sample by N (2.08 M , 0.07 M )/N (2.08 M , 0.09 M ), the ratio of the inferred mass estimate from [9] to the inflated mass estimate used in [55]. This allows for a more direct comparison between the results of [55] and [56]. We find a negligible effect on our results when we repeat our analysis with the increased mass uncertainty. Table I summarizes the mass and radius data we use for J0740+6620.  Unlike J0030+0451, the two independent analyses of J0740+6620 arrive at slightly different values for its radius, even if one accounts for their different priors (flat in mass-radius [56] vs. flat in mass-compactness [55]) and their different treatments of the uncertainty in the mass estimate from [9]. Accounting for the prior differences increases the discrepancy between the two results, as the flat-in-compactness prior disfavors large radii. Miller et al. [55] use the nominal XMM-Newton calibration uncertainty, while Riley et al. [56] use a larger uncertainty. The main effect of the XMM-Newton data is to provide an estimate of the pulsar count rate, which aids in the determination of the relative modulation depth of the x-ray pulse profile, which is essential for placing an upper limit on J0740+6620's compactness. Consequently, the larger calibration uncertainty of [56] results in a less stringent lower bound on the radius. We focus on results based on the analysis in [55], since it uses the nominal calibration uncertainty, although we provide select comparisons to the results of [56].
Nonetheless, we stress that hierarchical EoS constraints are unaffected by the choice of prior for the J0740+6620 radius measurement; any discrepancies are solely due to systematic differences between the two analyses, such as the choice of XMM-Newton calibration uncertainty or issues of convergence within sampling algorithms (see the discussion in Sec. 4.6 of [55]).

III. NEUTRON STAR MASS AND RADIUS
We apply our analysis to the combined radio, GW, and x-ray data and present the resulting constraints for macroscopic NS properties, notably masses, radii, and tidal deformabilities. In what follows, whenever we refer Constraints on the NS mass-radius relation. Shaded regions enclose the 90% symmetric credible intervals for the radius for each value of the mass. The left panel shows the effect of the J0740+6620 radius constraint by comparing the prior (black), and results with (without) the J0740+6620 radius in blue (turquoise). The right panel presents cumulative constraints on the mass-radius relation as each type of data set is analyzed. In black we again show the prior. The turquoise region shows the posterior after including the mass measurement of the two heavy pulsars (including the updated J0740+6620 mass estimate). The green region correspond to constraints obtained after adding the GW data. Finally, the blue region correspond to constraints after further adding the J0030+0451 and J0740+6620 mass and radius constraints from NICER. In the last case we remove the J0740+6620 mass constraint from the list of radio constraints so as to avoid double-counting.
to results without the J0740+6620 radius measurement, we still use its updated mass estimate from [9] within the inference. Unless otherwise stated, all results make use of the Miller et al. [55] mass and radius constraints without the inflated mass uncertainty.
We infer the NS mass-radius relation shown in Fig. 1, which plots the 90% symmetric credible region for R as a function of m. 3 The left panel focuses on the effect of the new J0740+6620 radius measurement: it tightens the 90% credible constraint on the radius from the low side by 0.57 km at 1.4 M and 0.71 km at 2.0 M . The right panel shows cumulative constraints on the massradius relation as the different data sets (radio, GW, x-ray) are added one at a time. As discussed in [50], the radio and x-ray observations tend to drive the lower bound on the NS radius, while the GW data and causality set the upper bound. This is because the GW measurements mainly constrain Λ ∼ (R/m) 6 from above, while the x-ray measurements primarily set an upper bound on the compactness m/R and therefore a lower bound on R. The J0740+6620 radius measurement reinforces this picture of complementary constraints.
One-and two-dimensional marginalized priors and posteriors for various macroscopic and microscopic parameters are given in Fig. 2, while Table II presents medians and 90% highest-probability-density credible regions for these and other quantities of interest. Like the left panel of Fig. 1, we compare the prior and posterior with and without the J0740+6620 radius constraint. However, we no longer restrict the prior to EoSs that support stable NSs at a given mass scale: it includes EoSs with M max significantly smaller than 1 M , such that, for example, the prior on Λ 1.4 peaks at the black hole value of zero. This distinction is less relevant for the posterior, as the data significantly disfavor EoSs that do not support M max 2 M . Nonetheless, it explains the shape of some priors in Fig. 2.
On the whole, we find that the J0740+6620 radius constraint increases support for stiffer EoSs with larger radii and tidal deformabilities. Our inferred M max is also slightly increased. This is because J0740+6620's radius is no smaller than that of a lower-mass NS, indicating that the turning point in the mass-radius relation occurs above the pulsar's mass. As discussed above, the tail of the M max posterior is slightly lower than its prior. This is driven by two factors: first the bound on R 1.4 provided by GW170817 that limits M max via causality considerations, and second, our assumption that the maximum NS possible is determined by the EoS and not NS formation mechanisms, resulting in EoSs that predict very heavy (and unobserved) NSs being disfavored. An upper limit on M max 2.2 − 2.6 M has been proposed by assuming that the electromagnetic counterpart to GW170817 suggests that the merger remnant collapsed to a BH shortly after merger [68,[99][100][101][102][103]. We do not employ this upper limit here (nor any other information from the GW170817 counterpart), and thus our inferred M max extends to higher values. Indeed the data sets we use can only stringently constrain M max from below. The effect of folding in such an upper limit is demonstrated in [79].
Based on Fig. 2, we also see that R 2.0 is more strongly correlated with the pressure at 2 ρ nuc than at 6 ρ nuc [49,104]. Additionally, J0740+6620's radius measurement from [55] eliminates the bimodality in the posterior on Λ 1.4 [47], now favoring the (initially subdom- The general trend in favor of stiffer EoSs also increases the lower bound of the 90% highest-probability-density credible region for Λ 1.4 (respectively, Λ 2.0 ) from 168 (7) to 265 (14). Setting the tidal deformability equal to this lower limit, we can obtain a conservative estimate of the signal-to-noise ratio (SNR) required for a GW observation to confidently detect tidal effects, i.e., bound Λ away from zero. The measurement uncertainty in Λ was ∼ 700 at an SNR of ∼ 33 for GW170817 [47]. Assuming that this measurement is typical and that uncertainties scale inversely with the SNR [105], a back-of-the-envelope estimate suggests that tidal effects can be measured to within 265 (14) for a binary with masses of 1.4 M (2.0 M ) with SNR of 44 (770). The threshold SNR for Λ 1.4 is within reach of current advanced detectors [106], although the SNR for Λ 2.0 will require next-generation detectors, consistent with the findings of [107].
The full EoS inference also allows us to obtain an updated radius estimate for J0740+6620 informed by all the data, as plotted in Fig. 3. We find 12.41 +0.93 −1.16 km at the 90% level, compared to 13.24 +2.25 −1.93 km when using only the J0740+6620 x-ray data conditioned on our nonparametric EoS model. For reference, the J0740+6620 measurement from [55] is 14.30 +4.33 −2.97 km at 90% credibility when adjusted to remove the 0.04 M systematic error estimate and intrinsic flat-in-compactness prior. The radius uncertainty for J0740+6620 at the 90% level is reduced by 3.12 km by conditioning on our EoS prior and further by 2.09 km when additionally including all our astrophysical data. Most of this improvement comes from the exclusion of large radii due to two reasons: (i) the EoS prior model favors realistic EoSs and a radius below ∼ 17 km, see prior in Fig. 1, and (ii) the GW data are inconsistent with large radii above ∼ 13 − 14 km. The updated radius estimate is consistent with the constraint of 12.28 +0. 60 −0.68 km from [55] (68% level) after conditioning on other data and their EoS prior.
We also investigate how our results change if we use the J0740+6620 data from [56] in place of the data from [55]. The two sets of inferred NS properties are compared in Table II. Because of their more conservative treatment of calibration error, the Riley et al. [56] data place a less constraining lower bound on J0740+6620's radius and therefore result in a more modest shift towards stiff EoSs. Out of the ∼ 0.8 km difference between the lower bounds of the 68% credible intervals on the pulsar's radius obtained by the two analyses, [55] attributes 0.55 km to the calibration difference and choices of prior boundaries. Our hierarchical analysis is immune to the prior difference, and after conditioning on all the observational data we find an overall difference of 0.4 km (respectively, 0.29 km) in the lower bound of the 90% credible interval on R 1.4 (R 2.0 ) due to other systematic differences between [55] and [56]

IV. PROPERTIES OF DENSE MATTER
We now turn our attention to the properties of dense matter and examine the implications of the J0740+6620 radius constraint. We begin in Fig. 4 with the inferred pressure-density relation. In the left panel, we show the effect of the new J0740+6620 radius constraint: it restricts the low-pressure side of the EoS at densities of 2-3 ρ nuc . This is comparable, but a bit lower, than the central density of J0740+6620, denoted by the magenta contours. In the right panel, we show the cumulative constraints that result from adding the different data sets sequentially. Red contours here denote the central pressure-density posterior for the maximum-mass NS.
The central density of J0740+6620 is 10.0 +3.5 −3.6 × 10 14 g/cm 3 ∼ 3.57 +1.3 −1.3 ρ nuc , as inferred from all available data under our EoS model. The relatively low inferred central density for a ∼ 2 M NS is indicative of a relatively stiff EoS at densities ∼ 1-2 ρ nuc ; see, e.g., Table  III of [108] for a comparison between two representative hadronic models. However, our analysis intentionally does not closely follow specific nuclear theoretic predictions. At low densities (up to ∼ 2 ρ nuc ), theoretical predictions from χEFT may place an upper limit on the pressure, which would tend to increase the central density of J0740+6620, although the most recent measurement of the neutron skin thickness of 208 Pb [17]   ρ ≤ ρc(Mmax)), and therefore the exact density range over which we maximize depends on the EoS. Columns correspond to the prior, the posterior with only the two heavy pulsars, and the posterior with and without the radius constraint from J0740+6620. Results with only the two heavy pulsars and without the radius constraint include the updated mass measurement of J0740+6620.
The column with the pulsar-only posterior is similar to the second column of Table IV in [50]. We include it as it roughly corresponds to the assumption that all objects up to ∼ 2 M are NSs (as opposed to our prior in some cases). We also present results based on both the Miller et al. [55] and the Riley et al. [56] analyses.  Fig. 1 but for the pressure-density relation. In the left panel, magenta contours give the 50% and 90% level of the central pressure-density posterior for J0740+6620 inferred from all available data. In the right panel, red contours give the 50% and 90% level of the central pressure-density posterior for the maximum-mass NS. a relatively stiff EoS below and around ρ nuc ; see [18] for more discussion. 5 We further investigate the NS central densities in Fig. 5, which shows the mass-central density posterior inferred using all the data. The central density of the maximum-mass NS is 1.5 +0. 3 −0.4 ×10 15 g/cm 3 ∼ 5.4 +1.1 −1.4 ρ nuc , corresponding to the maximum matter density that can be probed with observations of cold, nonspinning NSs. Table II also gives the central densities for NSs of 1.4 M and 2.0 M . In general, we can understand the trends in the central densities within the same context as Figs. 1  and 4. Typically, the central density remains low (stiff EoS) until masses are 2 M . Beyond this limit, set primarily by J0740+6620, the EoS can soften appreciably and the central density can increase considerably. Indeed, the density range explored by NSs above 2 M could be a factor of two times larger than what is explored by canonical 1.4 M stars. High-mass NSs may yet have surprises in store for future measurements.

A. Speed of sound
We examine the speed of sound inside NSs in Figs. 6 and 7. Figure 6 shows the speed of sound squared (c 2 s ) as a function of density with and without the J0740+6620 radius measurement. Already in [50], we concluded that the conformal limit of c 2 s /c 2 = 1/3 is likely violated inside NSs, primarily due to the combination of a soft lowdensity and a stiff high-density EoS and in agreement with [11,59,[75][76][77]. We here find that the lower limit on the J0740+6620 radius agrees with this picture and pushes the maximum of the marginal 90% lower limit for c 2 s to lower densities. In other words, the pressure needs to increase more rapidly at even lower densities in order to accommodate the relatively large radius of J0740+6620. The red contours corresponds to the central speed of sound and central density of the maximum-mass NS, again bounding the densities that can be probed observationally. The central speed of sound is essentially unconstrained, which means that, for some EoSs, the speed of sound sharply decreases after it reaches its maximum value. Figure 7 shows the maximum c 2 s inside NSs and the density at which it is reached. For each EoS, we maximize c 2 s over all densities smaller than the central density of the maximum-mass stellar configuration (i.e., a different range for each EoS). Comparing the posterior to the conformal limit, we again find that the latter is violated inside NSs with a maximum c 2 s /c 2 of 0.75 +0.25 −0.24 achieved 5 Fig. 2 of [109] depicts the central densities obtained by extrapolating the realistic two-and three-nucleon interactions predicted by microscopic theory to higher densities. The central values of pressure at around 2 ρnuc inferred from our analysis (see Table II) point to the stiffest EoS compatible with low-density chiral effective-field-theory (χEFT) [12,13,78,79,104].
at a density of 1.01 +6.3 −5.3 × 10 14 g/cm 3 (3.60 +2.25 −1.89 ρ nuc ). Compared to results without J0740+6620, the maximum speed of sound is slightly lower and occurs at slightly lower densities, as also seen in Fig. 6. This behavior was also observed for J0030+0451 [50]. Since both J0030+0451 and J0740+6620 data place a lower limit on the NS radius we interpret the reduced value for the maximum speed of sound as follows: the preference for a stiffer EoS at ∼ 2ρ nuc means that the stiff EoS at ∼ 5ρ nuc can be achieved with a milder pressure-density slope and thus a smaller speed of sound. The strongest support for a large speed of sound comes from the combination of GW and heavy pulsar data that point to a soft lowdensity and stiff high-density EoS respectively, thus necessitating a steep slope in between. Figure 7 also shows our prior on the maximum c 2 s and, even though it is consistent with the conformal limit, it certainly disfavors it.
To further assess the impact of data on the conformal limit in relation to the prior, Table III compares the evidence for EoSs that violate the conformal limit (max c 2 s > c 2 /3) with those that obey the conformal limit within nonrotating NSs through the corresponding Bayes factor: We find strong support that the conformal limit is violated: B c 2 s >c 2 /3 c 2 s ≤c 2 /3 10 3 . Although our prior is consistent with EoSs that obey the conformal limit, it includes relatively few realizations that do so. As such, our Bayes factors are subject to sizeable sampling uncertainty from the finite number of Monte Carlo samples we employ, making it hard to conclude whether support for the violation of the conformal limit increases or decreases due to J0740+6620. Nonetheless, we recover large Bayes factors, even considering this sampling uncertainty, and typically find that B This measures how well each type of EoS is able to fit the observed data, and Table III shows that EoSs that violate the conformal limit are typically favored over those that obey it by between a factor of 40-110.

B. Strong first-order phase transitions
We now turn our attention to the implications of J0740+6620 for strong phase transitions. Figure 8 compares the pressure-density posterior inferred with EoSs that support different numbers of stable branches in the mass-radius relation, used here as a proxy for strong   (3) and (4)]. We report point estimates and standard deviation from Monte Carlo sampling uncertainty. Data sets are labeled in the same way as in Table II. phase transitions. While strong first-order phase transitions can lead to EoSs with multiple stable branches and possibly even "twin stars", i.e., stars with roughly the same mass but very different radii [110], the converse is not necessarily true. Only the strongest phase transitions lead to disconnected branches, and so what follows concerns only the most extreme phase transitions. Typically, strong phase transitions and multiple branches result in a large decrease in the radius between subsequent branches [111][112][113][114]. The lower limit on the radius of J0740+6620 constrains such a sudden decrease.
Indeed, we find that the full posterior is similar (though not identical) to the one obtained from restricting to EoSs with only a single stable branch. This suggests that the full posterior marginalized over the number of branches is dominated by EoSs with a single stable branch, though this is also true of the prior. Table III reports the evidence ratios for EoSs with different numbers of stable branches: We find a Bayes factor of 6.9 in favor of a single stable branch, compared to < 5 without the J0740+6620 radius measurement. Astrophysical data generally disfavor the existence of multiple stable branches, driven primarily by the requirement that the EoS supports ∼ 2 M stars. As expected, the lower limit on the J0740+6620 further reduces the evidence for multiple stable branches. However, even the most extreme preference for a single stable branch only suggest a Bayes factor of 8.
Just as with the EoSs that obey vs. violate the conformal limit, we also report ratios in the maximum likelihood observed with EoSs that have a single stable branch vs. those with multiple stable branches Similar to B n>1 n=1 , we find a preference for EoSs with a single stable branch, but it is small (at most a factor of 2). Previous work reported B n>1 n=1 additionally conditioned on the existence of massive pulsars a priori [50], equivalent to dividing any B n>1 n=1 by the result using only the massive pulsar observations. This amounts to examining whether the GW and x-ray data are consistent with multiple stable branches, after we have already assumed the existence of ∼ 2 M stars. If we follow suit, One-and two-dimensional marginalized prior and posterior of the maximum c 2 s /c 2 encountered inside the NS and the density at which this happens. We show the prior in black and the posterior with (without) the J0740+6620 radius measurement in blue (turquoise). The vertical line denotes the conformal limit of c 2 s /c 2 = 1/3.
we obtain B n>1 n=1 ∼ 1.2 conditioning on the existence of massive PSRs a priori and including the x-ray observations of J0740+6620, compared to ∼ 1.8 reported in [50]. As such, we again find that x-ray observations of J0740+6620 lower the evidence in favor of multiple stable branches. Our conclusions are generally consistent with those reported in [79] (B n>1 n=1 ∼ 0.2), although our results disfavor multiple stable branches slightly more strongly.
A direct comparison is difficult as [79] do not quote uncertainties in their estimates. However, the observed differences could easily be due to priors (e.g., our priors allow for more model freedom and therefore contain more EoSs with multiple stable branches that are not forced a priori to support massive stars) or by how exactly phase transitions are defined (here we define them in terms of stable branches).
Several caveats should be kept in mind when interpreting our Bayes factors. Most importantly, it is well documented that Bayes factors are affected by the prior coverage of each model under consideration, particularly if they span regions of parameter space without any support a posteriori. That is to say, the marginal likelihoods that appear in, e.g., Eqs. (1) and (3) are averages of the likelihood over each prior; if priors span large regions of parameter space with small likelihood values, their marginal evidence will be smaller even if they achieve the same maximum likelihood (match the data just as well) as other, more compact priors. 6 Indeed, this is why we additionally report max L c 2 s >c 2 /3 c 2 s ≤c 2 /3 and max L n>1 n=1 . In particular, differences in the prior support are thought to be a driving factor behind the Bayes factors' apparent preference for EoSs with a single stable branch (multiplebranch EoSs span a broader range of behavior, comparing max L n>1 n=1 and B n>1 n=1 in Table III) as well as the preference for χEFT models over more agnostic EoS priors [14,115]. While this type of Occam factor is desirable in many cases (see, e.g., discussion in Sec. II A), one needs to take care when drawing conclusions based on such effects. Although not guaranteed, we generally find that prior choices of this kind shift our Bayes factor by only a factor of a few, typically much smaller than the variability due to different realizations of experimental noise [116], which can be as large as factors of O(10). It is therefore always prudent to check both the priors and posteriors for the behavior in question, for example checking both Fig. 7 and Table III when considering the conformal limit.
In the case of an EoS with multiple stable branches, we find a pressure-density envelope that is morphologically similar to the one in Fig. 4 of [50] (obtained without the J0740+6620 radius data). These plots show that, if the EoS has multiple stable branches, then the pressure is higher below nuclear saturation and lower at 2-3 ρ nuc , hinting towards a phase transition in this density regime and suggesting that all observed NSs may already contain an exotic core. Besides such a low-density phase transition, another possibility is a phase transition at higher densities. Such an effect is expected to lead to a reduction in the radius of no more than ∼ 3 km [104,108,113,114] for the most massive NSs compared to R 1.4 , which would have been undetectable before the J0740+6620 radius lower limit.
To further explore the implications of a sudden decrease in the radius, in Fig. 9 we plot the posterior for the radius difference ∆R ≡ R 2.0 −R 1.4 and the maximum c 2 s , broken down again by the number of stable branches.
A large negative value for ∆R suggests a strong phase transition at low densities, a scenario tightly constrained by the lower limit on the J0740+6620 radius. We find ∆R = −0.12 +0.83 −0.85 km, consistent with zero although with large uncertainties, as also demonstrated in [55,56]. This effectively rules only out the most extreme case of phase transitions that lead to a 2 km decrease in radii [114] but still remains consistent with milder or smooth phase transitions [108,117]. In the case of multiple stable branches, we find that the maximum c 2 s is higher than the single-branch case, though this does not seem to affect ∆R.
The larger speed of sound is consistent with previous work that suggests that, in the case of sharp phase transitions, the post-transition speed of sound in general needs to be larger in order to compensate for the intrinsic softening induced by the phase transition [113]. Indeed, as previously studied in [104], the absolute bounds on NS radii assuming an EoS with a constant sound speed at high densities is very sensitive to the assumed value of max c 2 s /c 2 . The lower (upper) radius bound decreases (increases) as max c 2 s /c 2 increases. This is in agreement with Fig. 9. The absolute lower bound on NS radii derived in [104] corresponds to the most negative value of ∆R induced by the strongest possible phase transition (limited by c 2 s ≤ c 2 at high densities) compatible with M max . On the other hand, a positive ∆R suggests weaker phase transitions, progressing towards the absolute upper bound on NS radii. We also note that for various physical models of hadronic matter (with or without a smooth crossover to exotic matter), ∆R −1.5 km is typical, although a few exhibit an increase from R 1.4 to R 2.0 [108,118,119].
To further explore this, in Fig. 10 we plot the posterior for the transition mass M t , defined as the largest mass of the first stable EoS branch, and the transition density ρ t , defined as the central density at the transition mass, and select macroscopic quantities. Fig. 10 also considers only EoSs with multiple stable branches. We plot R 1.4 , R 2.0 , and M max , which roughly represent the main observables from GWs, the two NICER pulsars, and the radio mass observations. We find that, if the EoS has multiple stable branches, the transition from the first branch probably happens for masses 1 M or ∼ 2 M . The corresponding transition density is 2.2 ρ nuc or ∼ 4.5 ρ nuc . Highdensity phase transitions would be the most challenging to detect, as they could result in small changes in the radius and thus be indistinguishable from EoSs without a phase transition [111].
The posteriors also indicate that transition mass M t and the radius difference ∆R are anticorrelated. If the transition mass is very low, then the entire star is mostly composed of exotic matter. As expected for quark stars, we find that ∆R is closer to zero and can even be positive, i.e., the most massive star is bigger (as expected for self-bound configurations). This is similar to the behavior of the two brown curves in Fig. 1 of [114]. As the transition mass increases, ∆R becomes more nega-tive. This is similar to the purple and red curves from Fig. 1 of [114] that result in stars that are hadronic in the outer layers but possess a large quark core. If future GW detections place further upper limits on R 1.4 , then large negative values for ∆R will be further constrained, thus pushing M t even lower.
The current data disfavor phase transitions that lead to multiple stable branches occurring in the mass range ∼ (1-2) M , suggesting that the majority of NSs we observe belong in a single branch: if the true EoS has multiple stable branches, either all sub-2 M contain exotic material or none do. The two-dimensional M t -M max plot reveals that this is due to the requirement that M max 2 M , which disfavors M t ∼ 1.5 M a priori and "splits" the M t posterior into two modes [112]. This behavior is expected, for example, from Fig. 3 of [113] which shows that an M max measurement constrains the intermediate values of the transition pressure. We leave extraction of further characteristics of the phase transition (such as the transition strength) and EoSs with phase transitions that do not lead to multiple stable branches to future work [120].
Finally, Fig. 11 shows distributions of the same variates as Fig. 2, but separates the EoSs with one and multiple stable branches. We find that all posteriors are consistent with each other, though EoSs with multiple stable branches are on average consistent with a softer EoS at low densities around 1-2 ρ nuc and a stiffer high-density EoS than single-branch EoSs. This is expressed through a slightly higher maximum mass and pressure at 6 ρ nuc , but slightly lower radii, tidal parameters, and pressure at 2 ρ nuc . The trend towards a stiffer high-density EoS agrees with the maximum speed of sound of Fig. 9. Similarly, the softer low-density EoS agrees with the pressure-density curves of Fig. 8.
Overall, we find that the data mildly disfavor multiple stable branches, though they do not rule out their presence. However, if the true EoS indeed has multiple branches, then this would suggest that some extra softening in the EoS, ostensibly due to a phase transition, has already taken place at densities below ∼ 2.2 ρ nuc . Currently, of all the available astronomical data sets, the GW data dominate the upper limits on the stiffness or pressure of EoS around 2 ρ nuc . Should further GW observations continue to push in the same direction, then the evidence for the presence of a strong phase transition in the relevant density region could be strengthened. To this end, astrophysical observations have limited constraining power at very low densities [14,18], and improved theoretical calculations or terrestrial experiments will likely determine whether the pressure is small or large near ρ nuc (see, e.g., [17,20]).

V. DISCUSSION
In summary, the new radius measurement for J0740+6620 refines our inference of the EoS by tightening the constraint on the pressure at densities ∼ 2-3 ρ nuc . Like NICER's previous observation of J0030+0451, this constraint comes mainly from the soft side, as the xray pulse-profiling available to date primarily bounds NS radii, tides, and pressures from below. We infer that all observed NSs have the same radius to within ∼ 2 km. This picture is consistent with other recent studies of J0740+6620 [55,78,79]. Our analysis draws three further principal conclusions: (i) the sound speed in NS cores very likely exceeds the conformal bound; (ii) the lack of a large radius difference between high-and low-mass NSs renders the existence of a separate stable branch in the mass-radius relation less likely; and (iii) the stiff EoS around 2 ρ nuc implied by the ensemble of observations results in a relatively low central density of 3.57 +1.3 −1.3 ρ nuc for J0740+6620, capping the density range that astronomical observations of nonrotating NSs can probe to date. However, the fact that the radius of J0740+6620 is comparable to R 1.4 suggests that J0740+6620 at ∼ 2.08 M might not be at the turning-point of the mass-radius curve and more massive NSs are possible; the inferred central density of the maximum-mass NS is 5.4 +1.1 −1.4 ρ nuc . Our main results are based on the J0740+6620 massradius constraint from [55], mainly due to the fact that this analysis uses the nominal relative NICER/XMM-Newton calibration uncertainty. Nonetheless, we find broadly consistent results when using the data from [56] instead. The larger calibration uncertainty assumed by [56] results in a weaker lower bound on the radius of J0740+6620, and after conditioning on all the observational data this translates to a 0.4 km difference in the lower limit of the 90% credible interval we extract for R 1.4 . Our conclusions about strong phase transitions and the violation of the conformal sound-speed bound are unaltered when the data from [56] is used. Prior differences in the two analyses (flat-in-compactness [55] vs. flat-in-radius [56]) do not affect results within the hierarchical inference formalism.
A direct numerical comparison between our results and [19,55,78,79] must be done with care due to the different data sets used and other assumptions. For example, [78,79] include GW170817 counterpart models, which we omit here due to concerns about systematic errors, and they assume a priori that GW190425 was a binary neutron star merger, which informs the M max inference because of its large primary mass. Nonetheless, with those caveats in mind, we can compare posterior constraints on the radius of a 1.4 M NS and the maximum NS mass. Ref. [55]  and more tightly constrained than our corresponding estimates, using a χEFT-informed parametric EoS model. Besides the aforementioned caveats, this difference can be partly attributed to the fact that Ref. [79] reports the radius with respect to a flat prior, whereas we report it, like all our constraints, with respect to the prior informed by our nonparametric EoS model. Both of these results refer to the J0740+6620 data from [55] and can therefore be compared to the second-last column in Table II ing a nuclear parameterization for the EoS with a piecewise polytrope extension. These numbers can be compared to the last column in Table II as they are based on the J0740+6620 data from [56]. The results from [78] in particular match our inferred values very closely, though our uncertainties are broader, which we attribute to the larger model freedom inherent in our GP prior. All these results are further broadly consistent with radius estimates from x-ray observations of NSs in low-mass xray binaries during quiescent or bursts phases [121][122][123][124][125], though these are subject to considerable modeling uncertainties.
This comparison of our results with the existing literature [55,56,78,79] brings forward the issue of model dependence in EoS constraints obtained from observations, experiments, and calculations that span many orders of magnitude in density. By design, the GP EoS prior used in [55] does not allow for as much model freedom as our model-agnostic process due to the strong intra-density correlations it assumes a priori. This is especially true at high densities. Another approach to nonparametric inference is to use neural networks, as in [126], though the model constructed in that study deliberately seeks to closely reproduce the behavior of a handful of tabulated EoSs from the literature. In this sense, the nonparametric models used in [126] and [55] are more analogous to the model-informed GP prior from [87] that makes relatively strong prior assumptions about correlations within the EoS. Parametric EoS models, such as piecewise polytropes [88], the spectral decomposition [89], and the speed-of-sound parameterization [11,90], impose even more restrictive assumptions on EoS morphology by virtue of specifying the functional form of the EoS with a finite number of parameters to describe an infinitedimensional function space. Examples of such model dependence are given in Fig. 10 of [55] and the variation between the two models presented in [78].
These considerations pose the problem of the degree to which EoS constraints are driven by the data, rather than by correlations between different densities imposed by the EoS model. Under that light, it is interesting to consider the effect of folding nuclear theoretic calculations into the inference of the EoS. Figure 4 shows that the J0740+6620 radius measurement does not inform the EoS below ρ nuc , something also confirmed by [55]. References [14,18] further show that our GP EoS prior is designed with no strong correlations between low-density information and high-density physics by explicitly showing that the same results are obtained at high densities regardless of whether the EoS is conditioned on χEFT at low densities or not. As such, we do not expect the J0740+6620 radius data to offer new insights about χEFT predictions within its regime of applicability, i.e.
1 − 2ρ nuc , nor do we expect χEFT predictions to influence our conclusions about NS matter at high densities. In contrast, the parametric EoS inference in [78] is sensitive to the χEFT calculations they condition on up to 1.1 ρ nuc even at the highest densities probed. Figure 7 of [78] shows that the NS radii and pressures they infer with both of their parametric EoS models have some dependence on which χEFT calculation is assumed. This suggests that statements both about the validity of nuclear calculations based on astrophysical data and the inference of NS properties after assuming a specific low-density calculation must take care to avoid introducing unwanted systematic modeling assumptions through the choice of high-density EoS representation.
In addition to χEFT and other theoretical models, several terrestrial experiments probe the EoS at densities up to ρ nuc . In particular, the PREX collaboration recently measured the neutron skin thickness of lead R 208 Pb skin [17], which is tightly correlated with the density dependence of the nuclear symmetry energy (the difference in the energy per particle for matter that contains only neutrons and matter that contains an equal number of neutrons and protons) and therefore the pressure at ρ nuc [16,127,128]. Using a model-agnostic nonparametric analysis similar to ours, Ref. [18] found no strong correlation between the results of several low-density experiments and high-density NS observables. Claims to the contrary [16,19], therefore, are driven by specific modeling assumptions, which may not be justified. Nevertheless, the large R 208 Pb skin reported in [17] suggests a relatively stiff EoS at low densities, although there are other low-energy experiments (e.g., [20,129,130]) and alternative interpretations of the data [131] that favor softer EoSs. 7 A stiff EoS at low densities may increase the evidence in favor of multiple stable branches, but we expect the effect to be small with current experimental uncertainties. Given the slight tension between nuclear experiments and the fact that additional constraints at low densities will not strongly influence our conclusions from J0740+6620's radius measurement, we omit nuclear experimental data from our current analysis and leave such investigations to future work.
Nevertheless, the growing number of constraints on the NS EoS is progressively sharpening our picture of dense matter. The radius measurement for J0740+6620 is a reminder of how different observations, experiments, and theoretical calculations complement each other by targeting different density scales inside NSs. Joint analyses of this ensemble of data require models for the EoS that span many orders of magnitude in pressure and density. As a result, it is important to understand how different EoS models, both parametric and nonparametric, correlate different densities to distinguish data-driven features from those driven by the prior. The nonparametric model we use is deliberately constructed to emphasize flexibility in EoS morphology and impose few correlations between high and low densities besides those dictated by the physical requirements of causality and thermodynamic stability. The intra-density correlations introduced by different parametric and nonparametric EoS models will be investigated in quantitative detail in upcoming work [72].
has made use of data, software and/or web tools obtained from the Gravitational Wave Open Science Center (https://www.gw-openscience.org), a service of LIGO Laboratory, the LIGO Scientific Collaboration and the Virgo Collaboration. Virgo is funded by the French Centre National de Recherche Scientifique (CNRS), the Italian Istituto Nazionale della Fisica Nucleare (INFN) and the Dutch Nikhef, with contributions by Polish and Hungarian institutes. This material is based upon work supported by NSF's LIGO Laboratory which is a major facility fully funded by the National Science Foundation. The authors are grateful for computational resources provided by the LIGO Laboratory and supported by National Science Foundation Grants PHY-0757058 and PHY-0823459.

APPENDIX: CHOOSING A MASS PRIOR
Substantial uncertainty persists in the distribution of compact-object masses, including the question of whether the NS and BH mass distributions overlap. However, any hierarchical analysis of the EoS needs to assume a compact-object mass distribution to specify the mass prior for a given EoS. 8 The full framework was laid out in Sec. III B of [50] where the likelihood for an EoS model ε was written in terms of the compact-object mass distribution P (m|ε). In the current study and [50], we assume a uniform mass distribution, such that the mass prior takes the general form However, a choice still needs to be made for the lower and upper limits, M lower and M upper . We set M lower = 0.5 M , but M upper is subject to three possible assumptions.
Assumption 1: The compact object might not be a NS The first option accounts for the possibility that the compact object in question is a BH, rather than a NS. If we do not have definite prior knowledge that it is a NS, our mass prior does not require its mass to be below the maximum NS mass. In this case, M upper is the maximum formation mass for the relevant type of compact object. If m > M max (ε), with M max (ε) the TOV mass of EoS ε, the radius and tidal deformability of the compact object are set to their Schwarzschild BH values, 2Gm/c 2 and 0, respectively. If m ≤ M max (ε), the compact object is a NS, and its properties are set by the EoS ε. The mass prior itself is independent of the EoS ε in this scenario. This is the assumption we employ for both GW170817 and GW190425. While it is the most agnostic assumption possible, it is clearly erroneous for compact objects detected as pulsars, such as J0740+6620.
Assumption 2: The astrophysical formation mechanism limits the maximum possible mass In the second scenario, the compact object under consideration is assumed to be a NS, but astrophysical formation mechanisms (for example, supernovae) are known a priori not to produce NSs above a certain mass, M pop . This upper limit might be comparable to M max (ε) for some EoSs. In this case, M upper = min (M pop , M max (ε)) .
While it is plausible that there may be an astrophysical limit to the mass of NSs, in practice this assumption comes at the expense of a completely arbitrary choice for M pop , given the current state of compact-object population knowledge. As such, we do not employ it for any of the compact objects analyzed in the main body of the paper. Nonetheless, Fig. 12 compares the impact of one possible choice of M pop against the other two assumptions.

Assumption 3: The EoS limits the maximum possible mass
Under the third assumption, the compact object under consideration is known to be a NS, and astrophysical formation mechanisms can produce NSs as heavy as the EoS can support. Thus, M upper = M max (ε). In this case, the prior depends on the EoS both through the upper limit (which rejects any masses above M max ) and the normalization (which constitutes an Occam penalty against EoSs that predict masses larger than have been observed). Stated differently, an EoS with M max (ε) slightly above the most massive known pulsar will be favored compared to an EoS that predicts the existence of much more massive NSs that have not been observed. We employ this assumption for all pulsars in our main study, including both radio and x-ray observations. To the best of our understanding, the same assumption is employed in [19,55,78].
The advantage of this assumption is that it does not rely on an explicit choice of M pop . The disadvantage is that the lack of observations of more massive NSs is attributed to (and therefore informs) the EoS, while other factors (astrophysical conditions and selection effects) are ignored, even as potential higher-mass NS candidates have been identified [132]. A simultaneous inference of the compact-object population and the EoS would obviate the need for choosing between Assumptions 2 and 3; instead, it would select the appropriate case as a function of the population model realization within the inference. This is possible because Assumption 3 is really just a special case of Assumption 2 in which M pop ≥ M max (ε) for all viable EoSs.
In order to quantitatively assess the impact of assumptions about the mass prior on J0740+6620, we repeat the main analysis with the alternative assumptions. In Assumption 1, we effectively assume that J0740+6620 could be a BH. In Assumption 2, we arbitrarily select M pop = 2.3 M , motivated by the approximate upper limit of the inferred J0740+6620 mass posterior [9]. In Assumption 3 (same as the main body of the text), we assume M pop = 3.0 M , which is larger than M max for the vast majority of EoSs in our prior. In Fig. 12, we plot the 2-dimensional and 1-dimensional M max − R 1.4 marginalized prior and posterior under these three assumptions. We find that the different mass prior choices only affect the inferred value of the maximum mass. Even then, the effect is small compared to current statistical uncertainties. Quantities determined at lower density scales, such as R 1.4 , are essentially unaffected.