Author Identifiers in Scholarly Repositories
Simeon Warner

Cornell Information Science and

Cornell University Library

Ithaca, NY 14850, USA

<EMAIL_ADDRESS>
Submitted: 2009-10-09

## Abstract

Bibliometric and usage-based analyses and tools highlight the value of
information about scholarship contained within the network of authors,
articles and usage data. Less progress has been made on populating and using
the author side of this network than the article side, in part because of the
difficulty of unambiguously identifying authors. I briefly review a sample of
author identifier schemes, and consider use in scholarly repositories. I then
describe preliminary work at arXiv to implement public author identifiers,
services based on them, and plans to make this information useful beyond the
boundaries of arXiv.

##  1 Context

In an ideal scholarly communication system there would be tools to browse,
navigate, make recommendations and assess influence based on the complete
graph of all actors (people, collaborations, institutions) and all
communication artifacts (articles, comments, blog posts, usage
data111Logically usage data would be links between actors and artifacts.
However, for historical, cultural and practical reasons most usage data is
treated as anonymous even though co-usage information may be extracted.). As a
shorthand I will call this complete graph the publication network. Contained
within it are the familiar citation, usage, co-authorship, and co-citation
graphs. In recent bibliometric and usage-based work, significant progress has
been made with the artifact part of this graph (see, for example the work of
the MESUR project [3]). Much less progress has been made with the actor part
of the graph, in part because it is much harder to unambiguously identify
authors than articles.

Lastname, Initial | Count
---|---
Zhang, Y | 100
Lee, J | 97
Wang, Y | 89
Wang, J | 84
Chen, Y | 77
Kim, J | 77
Wang, X | 76
Lee, S | 74
Kim, S | 69
Liu, Y | 69
Table 1:  Most frequently occurring ${lastname,initial}$ pairs in arXiv user
accounts. There may be a few duplicate accounts but this indicates that nearly
100 different people named “Zhang, Y” have created user accounts at arXiv (as
of May 2009).

Consider table 1 which shows the most frequently occurring
${lastname,initial}$ pairs in arXiv user accounts. This illustrates one facet
of the name disambiguation problem, namely that there are many authors with
the same name. This is compounded by inconsistent spellings, use of initials
or full first names, and even name changes. Within a single repository such as
arXiv it is not usually possible to accurately answer the question “show me
all the articles by this Zhang, Y”. In recent years there has been
considerable work on unsupervised and supervised author name disambiguation
using many different heuristic, machine learning and clustering techniques,
and many different properties including co-authorship, citations and
subjects/topics. While much better than naive approaches, these techniques are
still far from perfect.

In a recent Nature Correspondence, Raf Aerts asked “If it is possible to have
DOIs for objects (or, so they say, enough IPv6 addresses for every molecule on
Earth), why is it so difficult to implement DAIs [Digital Author Identifiers]
for authors?” [1]. Raf had earlier hinted at part of the answer by pointing
out that he has more than one identifier in Scopus [6]. As we have already
discussed, it is difficult to mine existing data to disambiguate references to
authors. The more fundamental part of the answer is that it is much easier to
create DOIs for articles when the one owner for an article creates the one DOI
for it and presents it with the article (ignoring the issue of multiple
versions of articles). As authors, we are not owned by a single authority and
even if an identifier were created for us at birth by the appropriate
government, there would be significant privacy concerns about using it for
everything. Consider, for example, concerns over the uses and misuses of
social security numbers in the USA. While we want to link a single author’s
works together, do we want that identity to immediately link us to all other
digital information about the private life of the individual?

##  2 Author Identifiers

To illustrate the diversity of currently used author identifiers, table 2
shows several example schemes used in the scholarly domain. A more detailed
inventory is provided on the repinf wiki [5]. The OpenID and ISNI schemes are
not limited to the scholarly domain. OpenID is aimed primarily at
authentication, however, if it continues to see growing acceptance it may well
be a useful open system that repositories could use. It is not clear whether
ISNI will develop into a widely used system. The largest efforts to create
author identifiers specifically for the scholarly domain, Scopus Author
Identifiers and ResearcherID, come from commercial entities and are clearly
motivated by the desire to provide improved services based upon then. It is
not clear how open the interfaces based on these identifiers will be, or what
data about them will be openly available.

Scheme |  Example |  Scope |  Authority
---|---|---|---
OpenID |  http://samruby.myopenid.com/ |  People |  Distributed, anyone supporting the protocol, relies upon DNS
ISNI (International Standard Name Identifier) |  ISNI 1422 4586 3573 0476 |  People |  Draft ISO standard requiring central DB operated by proposed International Agency
Scopus Author Id |  7103063073 |  Academia |  Elsevier
ResearcherID |  A-1637-2009 |  Academia |  Thomson Reuters
Digital Author Id |  info:eu-repo/dai/nl/304825271 |  Dutch Researchers |  Dutch Universities and Research Institutes
RePEc Author Service |  pzi1 |  Economics |  RePEc
arXiv Author Id |  http://arxiv.org/a/warner_s_1 |  arXiv.org |  arXiv.org
Table 2: A sample of identifier schemes used for scholarly author identifiers

The three other examples in table 2 illustrate decreasing scopes. In the
Netherlands the Digital Author Id (DAI) “is a unique national number assigned
to every author who has been appointed to a position at a Dutch university or
research institute or has some other relevant connection with one of these
organizations” [7]. The DAI provides a join point for data in different
repositories and enables services based on this combined data (e.g. NARCIS).
The RePEc scheme identifies authors and is used to link their publications
together within the RePEc system serving the economic community. The
AuthorClaim project aims to extend the RePEc model, using the same software
infrastructure, to the entire academic domain. Finally, the arXiv author
identifier, described below, is local to a single repository.

The arguments above and the understanding that there are many different
interests in, and uses for, author/person identities, suggest that there will
be many different systems and multiple identities for each author. In the
scholarly communication domain there will be a patchwork of overlapping
publication networks. Unless one system grows to dominate, the different
patches of publication network will identify authors using different
identifiers. However, it will be vastly easier to match multiple identifiers
for each author than to disambiguate multiple authors with the same name. A
significant aid will be the addition of assertions that link identities in
different networks (e.g. Author A2 in Repository 1 is the same person as
Author A4 in Repository 2) as illustrated in figure 1. This linking
information might be expressed either via the Semantic Web or in repository
metadata. The ability to record foreign identifiers in the author record
within a repository, and the ability to match articles between repositories,
will allow the joining of data across repository boundaries.

Figure 1: If authors A2 and A3 in repository 1 are identified with authors A4
and A5 in repository 2 respectively, then the components of the publication
graph from the two repositories can be joined. The P# nodes indicate papers.

##  3 Author Identifiers at arXiv

There are a significant number of physicists for whom all articles, or at
least all recent articles, are available on arXiv. It is not uncommon to find
web homepages with a link to arXiv author search in place of a bibliography —
why maintain the information in a second place when arXiv will do it
automatically? Fielded author search has been used in this way for many years
and has exactly the same problems of author disambiguation as text-based
efforts to build the publication network.

With the introduction of user accounts, arXiv, like many other repositories,
started to collect data on which user made each submission and whether he or
she claimed to be an author. This start to building authority records was
augmented by attempts to retrospectively associate older papers with users
based on email address matching, and the introduction of facilities by which
users could “claim ownership” of existing submissions. Use of the claim
ownership facility was motivated through the introduction of an endorsement
system222http://arxiv.org/help/endorsement where users must be known as
authors of a certain number of papers in order to endorse new users. Various
heuristics are used to limit what papers can be claimed automatically and so
far these have proved adequate to avoid incorrect claims being automatically
accepted.

Demanding identification of all authors at submission time was considered
impractical. For articles with one or two authors identification would not be
too burdensome, but for papers with 10 or even 2500 authors333Articles from
high-energy physics collaborations often have many authors. See, for example,
the recent ATLAS collaboration paper http://arxiv.org/abs/0901.0512 with
$>2500$ authors it is clearly impractical. A solution that uses arXiv
administrator effort to deal with each article is also impractical because
just two administrators handle all user queries relating to arXiv’s 58,000
submissions/year — most submissions must be entirely automated. We thus
decided on an approach that will create useful services based on a public
author identifier which we internally link to our user records. We hope that
by providing useful services our users will be motivated to further improve
the authority records on which these services depend.

###  3.1 Author URI and Services

We have opted for a web-centric approach using Linked Data [2] style HTTP
access. Each arXiv author identifier is a unique URI (e.g.
http://arxiv.org/a/warner_s_1) which supports HTTP content-negotiation. These
URIs are designed to be human copyable and are based on an ASCII dumb-down of
the author name. By default, or if selected via content-negotiation headers,
the arXiv author URI redirects to an HTML page listing all arXiv publications
authored by the given individual based on our user records. An example HTML
page is shown in figure 2. In cases where the author-article associations are
complete this facility already solves the problem of name collision in arXiv
author search and so provides a more reliable link than our text-based author
search. To allow the data to be used by other applications or to allow display
or monitoring with a feed reader, the list of articles associated with an
author id is also available as an Atom feed. Figure 3 shows that same data as
figure 2 but rendered from the Atom feed. As of September 2009 we see about
300 different author id URIs being accessed per week to return HTML pages or
Atom feeds.

Figure 2: HTML screen returned when an arXiv author id is accessed and HTTP
content negotiation results in HTML, or HTML is explicitly requested by
appending .html to the author id (http://arxiv.org/a/kurtz_m_1.html) Figure 3:
Web browser (Firefox) rendering of the Atom feed returned when an arXiv author
id is accessed and HTTP content negotiation results in Atom, or Atom is
explicitly requested by appending .atom to the author id
(http://arxiv.org/a/kurtz_m_1.atom)

A list of articles on the arXiv site is still one click away from the user’s
homepage. We thus provide JavaScript code, which we call the myarticles
widget, that a user may include in their personal homepage to dynamically
include an up-to-date publication list from arXiv. Various formatting options
are provided and the content may be styled using CSS. Figure 4 shows two
screen shots from early adopters of the myarticles widget. This facility is
based upon a content-negotiated request for an Atom representation of the
arXiv author id resource which results in a machine readable Atom feed of
paper information (in the same format as the arXiv
API444http://arxiv.org/api).

Figure 4: Two early adopter examples of use of the myarticles widget. The
example on the left shows the “Google ads” formatting option used to produce a
compact display in the lower left corner of the browser window. The example on
the right shows the “arXiv list” formatting option which picks up local
stylesheet information. In both cases the data from arXiv is seamlessly
embedded in the user’s homepage.

arXiv’s second use of arXiv author ids is to leverage this automatically
generated and updated list of publications to lower the effort required to
integrate arXiv papers into social networking sites. Facebook was chosen as
the first site to work with but the OpenSocial API is also being investigated.
Once the arXiv Facebook application has been told the association between a
user’s Facebook account and their arXiv author identifier, a list of
publications is immediately available as either a panel or a tab on their
Facebook profile as shown in figure 5. All title, author list, abstract and
linking information is automatically imported from arXiv. New or old
publications may be reported in the user’s feed, with optional comments, and
thus show up in friends’ news feeds. This application was released in March
2009 and sees approximately 600 users per week as of September 2009. Use is
steady and increasing, but there is not rapid adoption. We continue to
experiment with new facilities and modes of interaction.

Figure 5: Example screenshots from the arXiv Facebook application. The window
on the right shows a list of article titles and authors that is automatically
generated because Joshua Erlich enabled the arXiv application from with his
Facebook account and told it the association with his arXiv author id. Any new
arXiv article owned by him will automatically be added to this list.

###  3.2 Helping to Build the Publication Network

arXiv is making the multiple-identifier problem one identifier worse by
creating arXiv specific identifiers. Deduping articles is a key problem in
bibliometrics, and we don’t want to create a similar deduping problem with
author identifiers. The OpenID scheme explicitly caters for multiple
identifiers for a single person, and even for multiple identifiers for each
persona a single person might use. Facilities to express and leverage multiple
identities a described in the Yadis/XRDS document [8].

At arXiv we can go some way to addressing this problem by augmenting the Atom
format machine readable authorship information arXiv exposes with
correspondences between author identifiers in different schemes. The issue
then is how to encourage authors to supply and update alternative identifiers
associated with their account. Again we believe that the solution will be to
build useful services, such a links to other systems, that depend upon this
data.

Another option to encourage use of arXiv data on the relationships between
authors, papers authored and identifiers in other schemes is to expose it in
RDF. This is a good application for OAI-ORE Resource Maps [4] and we intend to
provide OAI-ORE resource maps as another representation available from the
arXiv author id. An example showing alternate name, article information and
how identifiers in other schemes can be exposed is illustrated in figure 6.

Figure 6: Possible OAI-ORE resource map representing the aggregation of papers
by an author Ang Lee. The node ReM is the Resource Map which describes the
aggregation Agg that has the author id as its URI
(http://arxiv.org/a/lee_a_1). The aggregation includes the three papers (P1,
P2, P3) authored by this Ang Lee. Through the ore:similarTo relation we also
indicate two other related resources: the identities in other author id
schemes. The DAI is a URI and so can be related directly as Id2. The
ResearcherID is a string and must therefore be related via an additional node
Id1.

##  4 Conclusions

There is growing interest in accurate author identification based on explicit
author identifiers. The many different commercial and non-commercial parties
have varying motives and goals, and so are adopting different solutions. It
seems likely that there will continue to be many different systems and
multiple identities for each author. While not perfect, this situation will
greatly assist the assembly of publication network data linking authors and
articles which will facilitate bibliometric analyses and support new discovery
tools that span multiple repositories. The implementation of author
identifiers at arXiv, and of services to promote their use has been described
to illustrate one approach at the repository level. Early use is encouraging
but it remains to be seen how quickly the use of author identification is
adopted and accepted by the scholarly community.

##  5 Acknowledgements

I am pleased to acknowledge contributions from Nathan Woody (Facebook and
JavaScript interface for arXiv), Thorsten Schwander and Paul Ginsparg. This
work is supported by Microsoft through a Technical Computing Initiative (TCI)
Grant. This paper is based on a presentation given at Open Repositories 2009
on 18 May 2009 but with updated usage data through September 2009.

##  6 Note added in proof

There has been considerable recent activity around author identifiers, the
most notable effort being ORCID (Open Researcher Contributor Identification
Initiative, http://orcid.securesites.net/) which has significant commercial
and community participation. Both Elsevier and Thomson Reuters are
participants in ORCID and the resulting system may replace both Scopus Author
Id and ResearcherID with a more open and more broadly adopted scheme.

## References

  * [1] Raf Aerts.  Digital identifiers work for articles, so why not for authors?  Nature, 453, 2008.  http://dx.doi.org/10.1038/453979b.
  * [2] Chris Bizer, Richard Cyganiak, and Tom Heath.  How to Publish Linked Data on the Web, 2007.  http://www4.wiwiss.fu-berlin.de/bizer/pub/LinkedDataTutorial/20070727/.
  * [3] Johan Bollen, Herbert Van de Sompel, and Marko A. Rodriguez.  Towards Usage-based Impact Metrics: First Results from the MESUR Project.  Proceedings of the Joint Conference on Digital Libraries 2008 (JCDL08), June 16, 2008, Pittsburgh, Pennsylvania, USA, 2008.  http://arxiv.org/abs/0804.3791.
  * [4] Open Archives Initiative – Object Reuse and Exchange Specifications, 2007.  http://www.openarchives.org/ore/toc.
  * [5] repinf — Author identification wiki.  http://repinf.pbworks.com/Author-identification.  International Repositories Infrastructure wiki initially created to support the March 2009 International Repositories Infrastructure Workshop in Amsterdam, NL. Last accessed 2009-10-09.
  * [6] Scopus Author Identifier.  http://www.info.scopus.com/.  The main Scopus service http://www.scopus.com/scopus/home.url and the details of author identifier are accessible from http://tinyurl.com/yta9m7.
  * [7] SURFfoundation.  Digital Author Identifier (DAI).  http://www.surffoundation.nl/smartsite.dws?ch=eng&id=13480.
  * [8] Yadis 1.0: The Identity and Accountability Foundation for Web 2.0, March 2006\.  http://yadis.org/papers/yadis-v1.0.pdf.  Page listing all versions is http://yadis.org/wiki/Yadis_Documents.


# A gap in the densities of small planets orbiting M dwarfs: rigorous
statistical confirmation using the open-source code RhoPop

J.G. Schulze School of Earth Sciences, The Ohio State University, 125 South Oval Mall, Columbus OH, 43210, USA  Ji Wang (王吉) Department of Astronomy, The Ohio State University, 140 West 18th Avenue, Columbus Ohio, 43210, USA  J. A. Johnson Department of Astronomy, The Ohio State University, 140 West 18th Avenue, Columbus Ohio, 43210, USA  B. S. Gaudi Department of Astronomy, The Ohio State University, 140 West 18th Avenue, Columbus Ohio, 43210, USA  R. Rodriguez Martinez Center for Astrophysics | Harvard & Smithsonian, 60 Garden St, Cambridge, MA 02138, USA  C.T. Unterborn Southwest Research Institute, 6220 Culebra Rd, San Antonio, TX 78238, USA  W.R. Panero School of Earth Sciences, The Ohio State University, 125 South Oval Mall, Columbus OH, 43210, USA  Division of Earth Sciences, U.S. National Science Foundation, Alexandria, VA, USA

(Received 2023 October 16; Revised 2024 January 16; Accepted 2024 January 30)

###### Abstract

Using mass-radius-composition models, small planets ($\mathrm{R}\lesssim
2\mathrm{R_{\oplus}}$) are typically classified into three types: iron-rich,
nominally Earth-like, and those with solid/liquid water and/or atmosphere.
These classes are generally expected to be variations within a compositional
continuum. Recently, however, Luque & Pallé (2022) observed that potentially
Earth-like planets around M dwarfs are separated from a lower-density
population by a density gap. Meanwhile, the results of Adibekyan et al. (2021)
hint that iron-rich planets around FGK stars are also a distinct population.
It therefore remains unclear whether small planets represent a continuum or
multiple distinct populations. Differentiating the nature of these populations
will help constrain potential formation mechanisms. We present the RhoPop
software for identifying small-planet populations. RhoPop employs mixture
models in a hierarchical framework and a nested sampler for parameter and
evidence estimates. Using RhoPop, we confirm the two populations of Luque &
Pallé (2022) with $>4\sigma$ significance. The intrinsic scatter in the Earth-
like subpopulation is roughly half that expected based on stellar abundance
variations in local FGK stars, perhaps implying M dwarfs have a smaller spread
in the major rock-building elements (Fe, Mg, Si) than FGK stars. We apply
RhoPop to the Adibekyan et al. (2021) sample and find no evidence of more than
one population. We estimate the sample size required to resolve a population
of planets with Mercury-like compositions from those with Earth-like
compositions for various mass-radius precisions. Only 16 planets are needed
when $\sigma_{M_{p}}=5\%$ and $\sigma_{R_{p}}=1\%$. At $\sigma_{M_{p}}=10\%$
and $\sigma_{R_{p}}=2.5\%$, however, over 154 planets are needed, an order of
magnitude increase.

Extrasolar rocky planets, Ocean planets, Open source software

††journal: PSJ

##  1 Introduction

Small planets ($\mathrm{R}\lesssim 2\mathrm{R_{\oplus}}$) are among the most
common in the Galaxy. Standard formation models predict that the relative
amounts of the major rock-building elements (Fe/Mg and Si/Mg) in small planets
should reflect that of their host stars (e.g., Bond et al., 2010a, b; Elser et
al., 2012; Johnson et al., 2012; Thiabaud et al., 2015). The inferred
compositions of observed small exoplanets (Schulze et al., 2021; Adibekyan et
al., 2021; Unterborn et al., 2023) and rocky exoplanetary material (Bonsor et
al., 2021) are consistent with this host-planet compositional connection at
current observational precisions except in extreme cases, e.g., Kepler-107 c
(Bonomo et al., 2019) and HD 137496 b (Azevedo Silva et al., 2022). In the
Solar System, Earth and Mars have compositions that are consistent with the
solar abundances of the major rock-building elements (e.g., Wanke & Dreibus,
1994; Lodders, 2003; Wang et al., 2019; Yoshizaki & McDonough, 2020). The
composition of Venus is unknown but expected to be consistent with Earth
(e.g., Zharkov, 1983; Aitta, 2012; Dumoulin et al., 2017). Mercury, in
contrast, has a 200-400% overabundance of Fe relative to the Sun, indicative
of a fundamentally distinct formation regime and/or subsequent evolution from
the other small planets in the Solar System (e.g., Cameron, 1985; Benz et al.,
1988; Wurm et al., 2013; Johansen & Dorn, 2022). If we, therefore, assume that
small planets are barren and rocky, then inconsistencies between inferred
planet composition and host composition indicate deviations from standard
planet formation models/theories.

Recently, Luque & Pallé (2022) observed a bimodal density distribution for 27
small planets in orbit around M dwarf hosts. The authors interpret this
bimodality as two compositionally distinct subpopulations, one with Earth-like
compositions and one with 50% water and 50% rock by mass, but do not provide
robust statistical evidence to support this claim or a comparison with host
abundances of the major rock-building elements. There is a paucity of known M
dwarf abundances of the major rock-building elements as measuring M dwarf
abundances is a non-trivial problem and an active area of research (e.g.,
Allard et al., 2000; Souto et al., 2022). As such, it is unclear whether there
are indeed two subpopulations of planets around M dwarfs and to what degree
these planets deviate from the host-planet compositional connection. Adibekyan
et al. (2021) show that 22 likely rocky small planets in orbit around FGK
stars have compositions correlated with their hosts’ but not with a 1:1
relationship. The authors identify a class of 5 planets with a Mercury-like
overabundance of iron separated from the general star-planet compositional
link by a visually apparent, but not statistically conclusive, gap in density.
Barros et al. (2022) identify two additional planets around HD 23472 that
appear to belong to the same iron-enriched population identified in Adibekyan
et al. (2021). Together these studies suggest that the small-planet density
distribution may not be a continuum. Instead, planet formation mechanisms
might create up to three distinct classes of small planets: (1) those that
follow the star-planet compositional link, (2) those enriched in
water/atmosphere and/or depleted in iron, and (3) those enriched in iron
and/or depleted in silicates.

Statistical host-planet composition comparisons like that of Adibekyan et al.
(2021) provide a powerful tool for confirming iron-enriched/silicate-depleted
or water-enriched/iron-depleted small planets (e.g., Santerne et al., 2018;
Schulze et al., 2021; Rodríguez Martínez et al., 2022; Azevedo Silva et al.,
2022; Unterborn et al., 2023). There is, however, a paucity of planets with
known densities and host abundances of Fe, Mg, and Si, and no single catalog
lists all of these quantities simultaneously. Combining the NASA Exoplanet
Archive111https://exoplanetarchive.ipac.caltech.edu/index.html with the
Hypatia Catalog stellar abundance
database222https://www.hypatiacatalog.com/hypatia (Hinkel et al., 2014), we
estimate there are approximately 43 planets with measured densities and host
Fe, Mg, and Si abundances. Of these, only 12 are small planets. Direct host-
planet comparisons are therefore limited to single planets or small samples
(e.g., 11, 22, and 7 in Schulze et al., 2021; Adibekyan et al., 2021;
Unterborn et al., 2023, respectively).

To overcome this lack of host-star abundances, some studies instead compare
inferred small-planet compositions with a large sample of stellar abundances
of the major rock-building elements (Scora et al., 2020; Plotnykov & Valencia,
2020; Unterborn et al., 2023). Unterborn et al. (2023) used the Hypatia
Catalog stellar abundance database (Hinkel et al., 2014) to quantify a
nominally rocky planet zone (NRPZ): the expected 99.7% range of rocky planet
densities owing to the intrinsic spread in stellar compositions. Where the
density and corresponding uncertainty of a small planet exclude it from the
NRPZ, its bulk composition likely deviates from that predicted through
standard formation models. Excluded planets whose densities exceed that of the
NRPZ are then interpreted as having an iron enrichment and/or silicate
depletion. Excluded planets with densities less than the NRPZ are interpreted
as being iron-depleted, water-rich, or volatile-rich. Applying this technique
to a large sample of small planets, Unterborn et al. (2023) find 21 planets
that are excluded from the NRPZ with $\geq 90\%$ confidence: 8 are likely
iron-enriched and 13 are iron-depleted, water-rich, or volatile-rich.

It remains an open question whether the densities of small planets reflect a
continuous variation in composition or multiple, distinct subpopulations.
Robust identification of these subpopulations, or lack thereof, is an
imperative first step to better constraining their formation and subsequent
evolution. Direct host-planet compositional comparisons are a powerful method
for identifying planets that deviate from standard formation models but are
sample-size-limited, making it difficult to resolve distinct subpopulations of
small planets. This limitation is often circumvented by instead comparing
planet compositions with a large sample of stellar abundances. While the
latter has been used to identify a number of individual planets that are
inconsistent with the NRPZ, no study has investigated whether these planets
form a continuum with the NRPZ or are members of fundamentally different
subpopulations.

In this work, we present the open-source RhoPop
333https://zenodo.org/records/10368148 code for robust identification of
compositionally distinct populations of small planets. RhoPop models the
density distribution of an input planet sample as a Gaussian mixture of up to
three compositionally-distinct subpopulations in a hierarchical Bayesian
framework. Similar models have been used for exoplanet population analyses,
but small planets are treated as a single population (Chen & Kipping, 2017) or
assume a fixed composition for rocky materials (Neil & Rogers, 2020; Neil et
al., 2022). RhoPop uses compositional grids that span from pure-iron to pure-
water planets that are normalized relative to the average composition of the
NRPZ, so that model results are easily comparable with the expected
compositional range of small water-/volatile-free planets. We describe the
inner workings of our software in Sect. 2. We use RhoPop to validate the
results of Luque & Pallé (2022) and Adibekyan et al. (2021) in Sect. 3.1 and
3.2, respectively. In Sect. 4, we use RhoPop to estimate the minimum number of
planets needed to identify a population of Mercury-like planets from a
population of Earth-like planets for various mass and radius uncertainties.

##  2 Methods

###  2.1 Density-mass-composition grids

We use the open-source ExoPlex444https://github.com/CaymanUnterborn/ExoPlex
mass-radius-composition calculator (Unterborn et al., 2023) to build density-
mass grids from pure water worlds to pure-iron planets. ExoPlex finds the
radius for a user-input planet mass and bulk composition by self-consistently
solving the five coupled differential equations: the mass within a sphere,
hydrostatic equilibrium, adiabatic temperature profile, Gauss’s law of gravity
in one dimension, and the thermally dependent equation of state. We assume a
pure-Fe core and adopt ExoPlex’s default equation of state (EOS) for liquid
iron (Anderson & Ahrens, 1994). ExoPlex couples the thermodynamic database of
Stixrude & Lithgow-Bertelloni (2011) with the Perple_X thermodynamic
equilibrium software (Connolly, 2009) to find the equilibrium mineralogy and
corresponding material density throughout the mantle. For outer water layers,
ExoPlex adopts the Seafreeze (Journaux et al., 2020) software for liquid
water, Ice Ih, II, III, V, or VI. For Ice VI, ExoPlex couples the isothermal
EOS of Journaux et al. (2020) with the empirical equations of Asahara et al.
(2010) and Fei et al. (1993).

We build isocomposition rock and liquid/solid water grids over a mass range of
0.05 to 10 $M_{\oplus}$ in increments of
$\Delta\log_{10}(\mathrm{M/M_{\oplus}})\simeq 0.023$. For all planets, we
assume the silicate mantle is Fe-free and fix the mantle molar ratios to Si/Mg
= 0.79, Al/Mg = 0.07, and Ca/Mg = 0.09 per Unterborn et al. (2023) as the
Fe/Mg ratio is the primary control on density for water-free rocky planets
(Dorn et al., 2015; Unterborn et al., 2016). For pure-rock compositions, we
vary $\log_{10}(\mathrm{Fe/Mg})$ from -0.2 to 1.5 (Fe/Mg = 0.01-30) in
increments of $\Delta\log_{10}(\mathrm{Fe/Mg})\simeq 0.035$. Our range of
Fe/Mg corresponds to a core mass fraction (CMF) range of 0.006 to 0.95. For
water-rich compositions, we build a grid of WMF = 0.001, 0.005, 0.01, 0.025,
and then in increments of $\Delta\mathrm{WMF}=0.025$ from 0.025 to 1.0. We fix
the interior rocky portion of these planets to have Fe/Mg = 0.71 corresponding
to the average Hypatia values per Unterborn et al. (2023). For compositions
between our computed grid lines, we linearly interpolate between the two
nearest compositional neighbors to approximate the intermediate isocomposition
line.

We note that there is an overlap between our rock grid and water grid between
WMF = 0 and 0.1-0.15, depending on mass, and Fe/Mg = 0 - 0.71 (CMF = 0-0.29),
meaning there are two grid solutions for a given mass and radius in this
regime: one with water and one without. This is a result of a well-known
degeneracy when inferring planet composition from density alone (e.g., Rogers,
2015; Dorn et al., 2015; Unterborn et al., 2016) which persists despite our
simplifications of a pure-Fe core and Fe-free silicate mantle. We circumvent
this degeneracy by using a reduced rock grid from Fe/Mg = 0.71-30 (CMF =
0.29-0.95). For planet densities lower than a water-free planet with Fe/Mg =
0.71, RhoPop switches to the water grid. Where our models return a WMF =
0.0-0.15, however, we also report the corresponding best fit Fe/Mg and CMF for
a water-free planet.

Figure 1: Water grid (cyan) and reduced rock grid (black). The $\rho_{scaled}$
parameter corresponds to the density at mass of a water- and volatile-free
rocky planet with a liquid pure-iron core and Fe-free silicate mantle with
Fe/Mg = 0.71, Si/Mg = 0.79, Ca/Mg = 0.09, and Al/Mg = 0.07.

Following the work of Luque & Pallé (2022), we normalize all planet densities
to a scaled density parameter denoted as $\rho_{scaled}$. The scaled density
parameter is a function of planet mass and a reference composition. It gives
the density a planet of a given mass would have with the reference
composition. Both our work and Luque & Pallé (2022) assume volatile-/water-
free compositions for reference. Where Luque & Pallé (2022) uses a nominally
‘Earth-like’ CMF = 0.325, however, we assume the average CMF = 0.29 of the
Hypatia Catalog corresponding to Fe/Mg = 0.71, Si/Mg = 0.79, Ca/Mg = 0.09, and
Al/Mg = 0.07. That is, $\rho_{scaled}(M)=\rho(M,\mathrm{CMF=0.29})$. We choose
this value so any planet populations inferred with RhoPop may be directly
compared with the nominally rocky planet zone. We can then define a general
density ratio parameter $\rho_{ratio}\equiv\rho/\rho_{scaled}$. For the
ExoPlex-derived density-mass-grids, this is a function of mass and composition

$\rho_{ratio}=\rho(M,\mu)/\rho_{scaled}(M)$ (1)

where $\mu$ is the composition in terms of a CMF from 0.29 to 1.0 or a WMF
from 0 to 1.0. Our water grid and reduced rock grid are shown relative to
$\rho_{scaled}$ in Fig. 1.

###  2.2 Hierarchical Bayesian Model Framework

Inspired by the population insights of Chen & Kipping (2017) and Neil et al.
(2022), RhoPop employs Gaussian mixture models in a hierarchical Bayesian
framework. Hierarchical Bayesian models (HBM) employ two sets of parameters,
local and hyper. We treat each observed data point as a Gaussian random
variable centered on its true value with a standard deviation equal to the
observational uncertainty. Since the true values can not be known a priori,
these are treated as free model parameters called ‘local parameters.’ The
parameters that govern the model that predicts the local parameters are the
‘hyperparameters’. Here the hyperparameters are those that describe the
underlying populations of planets.

The local parameters of our HBM are the true planet masses, $M_{t}$, and
density ratios, $\rho_{ratio,t}$. Observed planet mass, $M_{\mathrm{ob}}$, is
modeled as a random Gaussian variable centered on the true mass, $M_{t}$, with
a standard deviation equal to the observed mass uncertainty,
$\sigma_{M_{\mathrm{ob}}}$. That is,

$M_{\mathrm{ob}}^{i}\sim\mathcal{N}\left(M_{t}^{i},\sigma_{M_{\mathrm{ob}}^{i}}\right)$
(2)

where the superscript $i$ denotes the $i$-th planet in the sample. For planets
with asymmetric observational mass uncertainties, we take the average of the
upper and lower uncertainties to calculate $\sigma_{M_{\mathrm{ob}}^{i}}$.
Evaluating the right-hand side of Eq. 2 at $M_{\mathrm{ob}}^{i}$ gives the
likelihood function for $M_{t}^{i}$, which we denote as $\mathcal{L}_{m}^{i}$,
i.e.,

$\mathcal{L}_{m}^{i}=\mathcal{N}\left(M_{\mathrm{ob}}^{i}\bigg{|}M_{t}^{i},\sigma_{M_{\mathrm{ob}}^{i}}\right).$
(3)

We link the local parameters to the hyperparameters that describe the planet
populations using a Gaussian mixture model (GMM). Gaussian mixture models
represent some global population as a combination of $N_{c}$ normally
distributed subpopulations. Each subpopulation has a mixing weight, $w$, which
is its fractional size relative to the global population. The sum of the
mixing weights over all $N_{c}$ components must sum to unity, i.e.,

$\sum_{j=1}^{N_{c}}w^{j}=1$ (4)

where we use the superscript $j$ to denote the $j$-th subpopulation. We assume
that each $j$ component of our Gaussian mixture models is centered on
composition $\mu_{comp}^{j}$ with an intrinsic scatter of density ratios
$\sigma_{\rho_{ratio,pop}^{j}}$. Within a given population,
$\sigma_{\rho_{ratio,pop}^{j}}$ represents the variation that cannot be
accounted for by the observational uncertainties of the planets that belong to
it. We then model the local parameter $\rho_{ratio,t}^{i}$ as a Gaussian
mixture of all $N_{c}$ components

$\rho_{ratio,t}^{i}\sim\sum_{j=1}^{N_{c}}w^{j}\times\mathcal{N}\left(f(M_{t}^{i},\mu_{comp}^{j}),\sigma_{\rho_{ratio,pop}^{j}}\right).$
(5)

The argument of the Gaussian in Eq. 5, $f(M_{t}^{i},\mu_{comp}^{j})$, is
calculated directly from the compositional grids (Sect. 2.1) and represents
the predicted density ratio at $M_{t}^{i}$ given a central composition of
population $j$, $\mu_{comp}^{j}$. We then model the observed density ratio of
the $i$-th planet, $\rho_{ratio,\mathrm{ob}}^{i}$, as being centered on
$\rho_{ratio,t}^{i}$ with an error term corresponding to its observed
uncertainty,

$\begin{split}\rho_{ratio,\mathrm{ob}}^{i}\sim&\sum_{j=1}^{N_{c}}\left[w^{j}\times\mathcal{N}\left(f\left(M_{t}^{i},\mu_{comp}^{j}\right),\sigma_{\rho_{ratio,pop}^{j}}\right)+\mathcal{N}\left(0,\sigma_{\rho_{ratio,\mathrm{ob}}^{i}}\right)\right]\\\
&=\sum_{j=1}^{N_{c}}w^{j}\times\mathcal{N}\left(f\left(M_{t}^{i},\mu_{comp}^{j}\right),\sqrt{\left(\sigma_{\rho_{ratio,pop}^{j}}\right)^{2}+\left(\sigma_{\rho_{ratio,\mathrm{ob}}^{i}}\right)^{2}}\right).\end{split}$
(6)

As with Eq. 2, evaluating the right-hand side of Eq. 6 at
$\rho_{ratio,\mathrm{ob}}^{i}$ gives the likelihood function for
$\rho_{ratio,t}^{i}$,

$\mathcal{L}_{\rho}^{i}=\sum_{j=1}^{N_{c}}w^{j}\times\mathcal{N}\left(\rho_{ratio,\mathrm{ob}}^{i}\bigg{|}f\left(M_{t}^{i},\mu_{comp}^{j}\right),\sqrt{\left(\sigma_{\rho_{ratio,pop}^{j}}\right)^{2}+\left(\sigma_{\rho_{ratio,\mathrm{ob}}^{i}}\right)^{2}}\right).$
(7)

Finally, combing Equations 3 and 7, taking the natural logarithm, and summing
over all $N_{p}$ planets gives the log-likelihood function for the entire
model as a function of the hyperparameters,

$\begin{split}\ln\mathcal{L}&=\sum_{i=1}^{N_{p}}\left[\ln\mathcal{L}_{\rho}^{i}+\ln\mathcal{L}_{m}^{i}\right]\\\
&=\sum_{i=1}^{N_{p}}\bigg{[}\ln\left(\sum_{j=1}^{N_{c}}w^{j}\times\mathcal{N}\left(\rho_{ratio,\mathrm{ob}}^{i}\bigg{|}f\left(M_{t}^{i},\mu_{comp}^{j}\right),\sqrt{\left(\sigma_{\rho_{ratio,pop}^{j}}\right)^{2}+\left(\sigma_{\rho_{ratio,\mathrm{ob}}^{i}}\right)^{2}}\right)\right)\\\
&+\ln\left(\mathcal{N}\left(M_{\mathrm{ob}}^{i}\bigg{|}M_{t}^{i},\sigma_{M_{\mathrm{ob}}^{i}}\right)\right)\bigg{]}.\end{split}$
(8)

Our HBM framework is summarized schematically in Figure 2 and all parameters
are summarized in Table 1. In this work, we compare models through Bayesian
model evidence $\mathcal{Z}$: the integral of the likelihood (Eq. 8) over the
entire parameter prior space. Mathematically,

$P(\mathbf{D}|M)\equiv\mathcal{Z}=\int_{\Omega_{\mathbf{\Theta}}}\mathcal{L}(\mathbf{D}|\mathbf{\Theta},\mathcal{M})Pr(\mathbf{\Theta}|\mathcal{M})d\mathbf{\Theta},$
(9)

where $\mathbf{\Theta}$ represents a generic set of local parameters and
hyperparameters, $Pr$ are the corresponding priors, and $\mathcal{M}$ is the
model. We choose wide uninformative uniform priors on all hyperparameters to
remain agnostic when validating previous results in Sections 3.1 and 3.2. We
summarize all hyperparameter priors in Table 2. RhoPop handles up to three
populations, i.e., $N_{c}=1$, 2, or 3. RhoPop makes no a priori assumptions
about the central compositions of the populations $\mu_{comp}^{j}$ other than
$\mu^{1}_{comp}\geq\mu^{2}_{comp}$ for $N_{c}=2$ and
$\mu^{1}_{comp}\geq\mu^{2}_{comp}\geq\mu^{3}_{comp}$ for $N_{c}=3$. In the
3-population scenario, we allow $w^{1}$ and $w^{2}$ to vary from 0 to 1.
However, if the sum of these 2 is greater than one, the likelihood function
returns a $-\inf$, and that set of parameters is rejected since Eq. 4 must be
satisfied. We assume priors of
$M_{t}^{i}\sim\mathcal{N}(M_{\mathrm{ob}}^{i},\sigma_{M_{\mathrm{ob}}^{i}})$
for all $i$ in $N_{p}$.

Figure 2: HBM framework. Hyperparameters are shown in blue, local parameters in white, and observed parameters in gray. The observed and local parameters always have $N_{p}$ members corresponding to the number of planets in the sample. $N_{c}=1,2,$ or 3 when the 1-, 2-, and 3-population model is considered, respectively. Observed masses and density ratios are modeled as Gaussian random variables centered on their unknown true values (local parameters) with a standard deviation equal to their observational uncertainties (Eq. 2 and 6, respectively). The local and hyperparameters are linked through the mixture model defined in Eq. 5. Parameter | Parameter Type | Description
---|---|---
$\mu_{comp}^{j}$ | hyper |  Central composition of population $j$ in WMF or CMF
$\sigma_{\rho_{ratio,pop}^{j}}$ | hyper |  intrinsic 1$\sigma$ density ratio scatter of population $j$
$w^{j}$ | hyper |  Mixing weight of population $j$
$\rho_{ratio,t}^{i}$ | local |  True density ratio of planet $i$
$M_{t}^{i}$ | local |  True mass of planet $i$
$M_{\mathrm{ob}}^{i}$ | data |  Central observed mass of planet $i$
$\sigma_{M_{\mathrm{ob}}^{i}}$ | data |  $1\sigma$ uncertainty in $M_{\mathrm{ob}}^{i}$
$\rho_{ratio,\mathrm{ob}}^{i}$ | data |  Central observed density ratio of planet $i$
$\sigma\rho_{ratio,\mathrm{ob}}^{i}$ | data |  $1\sigma$ uncertainty in $\rho_{ratio,\mathrm{ob}}^{i}$
Table 1: Summary of parameter definitions.

In practice, calculating $\mathcal{Z}$ analytically is not tractable and it
must be computed numerically. As such, RhoPop uses the
dynesty555https://github.com/joshspeagle/dynesty dynamic nested sampling
software package with the default weighting scheme which is optimized to
simultaneously estimate the parameter posteriors and model evidence. We point
the reader to Speagle (2020) for more details.

####  2.2.1 Model selection

We adopt the methods outlined in Benneke & Seager (2013) which use Bayes
factors to compare models of increasing complexity. A general Bayes factor
$B_{ij}$ is defined as the ratio of the Bayesian evidence of model $i$ to the
Bayesian evidence of model $j$, i.e.,
$B_{ij}\equiv\mathcal{Z}_{i}/\mathcal{Z}_{j}$. In terms of hypothesis testing,
the denominator represents the evidence for the null hypothesis and the
numerator is the alternative hypothesis. A value of $B_{ij}<1$ suggests the
data favor model $j$ over model $i$ i.e., the data are consistent with the
null hypothesis. A $B_{ij}$ from 1-3 corresponds to inconclusive support for
the more complex model $i$, and the null hypothesis cannot be rejected. A
$B_{ij}$ from 3-12, 12-150, and $>150$ correspond to weak, moderate, and
strong support, respectively, for the alternative hypothesis. For a given
sample in this work, we first run the 1- and 2-population models and compare
via $B_{21}=\mathcal{Z}_{2}/\mathcal{Z}_{1}$ where $\mathcal{Z}_{1}$ and
$\mathcal{Z}_{2}$ correspond to the evidence for the 1- and 2-population
models, respectively. Where there is significant support for the 2-population
model ($B_{21}>150$), we then run the 3-population model and compare it with
the 2-population model via $B_{32}$. Benneke & Seager (2013), and references
therein, provide formulae to convert Bayes factors to the frequentist
$p$-value (see their Eq. 10) and corresponding sigma significance $n_{\sigma}$
(see their Eq. 11). For all Bayes factors in this work, we use these formulae
to calculate and report the equivalent $p$-values and $n_{\sigma}$
significance levels for the frequentist reader.

Param | $N_{c}=1$ | $N_{c}=2$ | $N_{c}=3$
---|---|---|---
$\mu_{comp}^{1}$ |  $\mathcal{U}$[1.0 WMF, 0.95 CMF] |  $\mathcal{U}$[1.0 WMF, 0.95 CMF] |  $\mathcal{U}$[1.0 WMF, 0.95 CMF]
$\sigma_{\rho_{ratio,pop}^{1}}$ |  $\mathcal{U}$[0,2] |  $\mathcal{U}$[0,2] |  $\mathcal{U}$[0,2]
$w^{1}$ | 1 |  $\mathcal{U}$[0,1] |  $\mathcal{U}$[0,1]
$\mu_{comp}^{2}$ |  | $\mathcal{U}[0,1]\times\mu_{comp}^{1}$ | $\mathcal{U}[0,1]\times\mu_{comp}^{1}$
$\sigma_{\rho_{ratio,pop}^{2}}$ |  |  $\mathcal{U}$[0,2] |  $\mathcal{U}$[0,2]
$w^{2}$ |  |  1-$w^{1}$ |  $\mathcal{U}$[0,1]
$\mu_{comp}^{3}$ |  |  | $\mathcal{U}[0,1]\times\mu_{comp}^{2}$
$\sigma_{\rho_{ratio,pop}^{3}}$ |  |  |  $\mathcal{U}$[0,2]
$w^{3}$ |  |  |  1$-w^{1}-w^{2}$
Table 2: Hyperparameter priors. The default prior on $\mu_{comp}^{1}$ spans
the entirety of our water and reduced rock grids.

##  3 Validating Previous Results

###  3.1 Luque & Pallé (2022) sample

The full sample of Luque & Pallé (2022) contains 34 planets with mass and
radius precisions of better than 25% and 8%, respectively, and radii less than
4 $R_{\oplus}$. The authors identify 7 planets as mini-Neptunes. We do not
consider these planets further as they require an extended H/He envelope to
explain and do not meet our definition of a small planet. We instead focus on
the remaining 27 planets from which Luque & Pallé (2022) observe a bimodal
density distribution: a higher density population with 21 members and a lower
density population with 6 members. The authors interpret the higher-density
population as being purely rocky in composition and the lower-density
population as having a rocky interior surrounded by a 50% H2O layer by mass.

We run this 27-planet sample through RhoPop with $N_{c}=1$ and 2, i.e., 1- and
2-population models. We visualize our results in Figure 3 and summarize the
hyperparameter posterior statistics in Table 3. We find the 2-population model
is strongly preferred with a Bayes factor $B_{21}=1042$. In frequentist
parlance, this corresponds to a $p$-value of $3.4\times 10^{-5}$ or a
4.1$\sigma$ detection of 2-populations. For good measure, we run this sample
using the 3-population model. We find $B_{32}=0.26$ meaning the 2-population
model is indeed favored over the 3-population model.

As further validation, our 2-population model identifies the same 21 planets
as Luque & Pallé (2022) belonging to a higher density population. We find this
population is best fit by a central WMF of 0.02. This WMF falls within the
degenerate zone described in Sect. 1. As such, we apply a $\chi^{2}$
minimization with these 21 planets over our full rock grid to find the
equivalent CMF assuming the planets are water-free. We find a CMF = 0.28
(Fe/Mg = 0.67), meaning this population is either slightly wet or has a $\sim
6\%$ Fe-depletion relative to the average volatile-/water-free NRPZ
composition of CMF = 0.29 (Fe/Mg = 0.71).

We find the same 6 planets that Luque & Pallé (2022) identify as being 50%
water worlds as indeed belonging to a separate lower-density population. We,
however, find this population is best fit by a central WMF of 0.78. While this
WMF value is substantially more than that predicted by Luque & Pallé (2022),
this is primarily due to the differences in the underlying material equations
of state used in this work versus Luque & Pallé (2022). Both populations are
found to have very little intrinsic scatter. Even at 10$M_{\oplus}$, where the
difference in density ratios is minimized over the mass range considered, the
two populations are separated by more than 12x their mutual intrinsic scatter
parameters, meaning there is a robust density gap between these two small-
planet populations.

| $N_{c}=1$ | $N_{c}=2$
---|---|---
Param | Med | 68% CI | Med | 68% CI
$\mu_{comp}^{1}$ | 0.09 WMF | [0.14 WMF, 0.06 WMF] | 0.02 WMF | [0.03 WMF, 0.017 WMF]
$\sigma_{\rho_{ratio,pop}^{1}}$ | 0.21 | [0.17, 0.25] | 0.02 | [0.005, 0.04]
$w^{1}$ | 1 |  | 0.76 | [0.67, 0.83]
$\mu_{comp}^{2}$ |  |  | 0.78 WMF | [0.87 WMF, 0.69 WMF]
$\sigma_{\rho_{ratio,pop}^{2}}$ |  |  | 0.03 | [0.008, 0.07]
$w^{2}$ |  |  |  1-$w^{1}$ |
|  | $\ln\mathcal{Z}_{1}=7.20\pm 0.24$ |  | $\ln\mathcal{Z}_{2}=14.15\pm 0.32$

$B_{21}$ |  $p$-value | Sigma | Interpretation
---|---|---|---
1042 |  3.4$\times 10^{-5}$ |  4.1$\sigma$ |  Strong evidence for two populations

Table 3: Summary of parameter posteriors, model evidence, and model selection
parameters for our 1- and 2-population model runs with the Luque & Pallé
(2022) sample. CI = credible interval

.

Figure 3: Density ratio as a function of mass for our $N_{c}=1$ (left) and
$N_{c}=2$ (right) models applied to the Luque & Pallé (2022) sample of small
transiting planets in orbit around M-dwarfs. We color code the highest-density
population as green and the lowest-density population as purple. The planets
are assigned the same color as the population to which they most likely
belong. The solid lines are the isocomposition lines for each
$\mu_{comp}^{j}$. The bands are the corresponding intrinsic density ratio
scatter parameter $\sigma_{\rho_{ratio,pop}}^{j}$. For all $\mu_{comp}^{j}$
and $\sigma_{\rho_{ratio,pop}}^{j}$ we use the median values given in Table 3.

###  3.2 Adibekyan et al. (2021) sample

The Adibekyan et al. (2021) sample includes 22 likely volatile-free rocky
planets with $M_{p}<10M_{\oplus}$ and mass and radius uncertainties $<30\%$
around FGK stars. The authors identify a visually apparent, but not
statistically significant, subpopulation of 5 iron-rich planets. We run our 1-
and 2-population models with this sample and find $B_{21}=0.39$, meaning there
is no support for the more complex 2-population over the 1-population model.
Thus, the single population model shown in Fig. 4 is our preferred model for
this data set. Given the lack of support for the 2-population model for this
data set, we do not consider the 3-population scenario. The results for the
2-population and preferred 1-population model are summarized in Table 4.

| $N_{c}=1$ | $N_{c}=2$
---|---|---
Param | Med | 68% CI | Med | 68% CI
$\mu_{comp}^{1}$ | 0.31 CMF | [0.02 WMF, 0.36 CMF] | 0.36 CMF | [0.30 CMF, 0.53 CMF]
$\sigma_{\rho_{ratio,pop}^{1}}$ | 0.16 | [0.11, 0.22] | 0.17 | [0.06, 0.69]
$w^{1}$ | 1 |  | 0.64 | [0.07, 0.94]
$\mu_{comp}^{2}$ |  |  | 0.05 WMF | [0.53 WMF, 0.33 CMF]
$\sigma_{\rho_{ratio,pop}^{2}}$ |  |  | 0.19 | [0.09, 0.77]
$w^{2}$ |  |  |  1-$w^{1}$ |
|  | $\ln\mathcal{Z}_{1}=-25.991\pm 0.223$ |  | $\ln\mathcal{Z}_{2}=-26.93\pm 0.21$

$B_{21}$ |  $p$-value | Sigma | Interpretation
---|---|---|---
0.39 | 1.0 | 0 | No evidence for more than one population

Table 4: Summary of parameter posteriors, model evidence, and model selection
parameters for our 1- and 2-population model runs with the Adibekyan et al.
(2021) sample. CI = credible interval

.

Figure 4: Density ratio as a function of mass for the preferred $N_{c}=1$
model (left) and the $N_{c}=2$ model (right)

applied to the Adibekyan et al. (2021) sample. We color code the highest-
density population as green and the lowest-density population as purple. The
planets are assigned the same color as the population to which they most
likely belong. The solid lines are the isocomposition lines for each
$\mu_{comp}^{j}$. The bands are the corresponding intrinsic density ratio
scatter parameter $\sigma_{\rho_{ratio,pop}}^{j}$. For all $\mu_{comp}^{j}$
and $\sigma_{\rho_{ratio,pop}}^{j}$ we use the median values given in Table 4.

##  4 Minimum number of planets required to resolve a population of iron-rich
planets

Here we assume the 5 potentially iron-rich planets identified by Adibekyan et
al. (2021) do indeed represent a distinct higher-density population from a
lower-density population comprised of the remaining 17 planets and that our
null finding in Sect. 3.2 is a result of observational uncertainty alone. We
then estimate the minimum number of planets needed to resolve a population of
iron-enriched planets from a lower-density population of volatile-/water-free
rocky planets. Evidencing the Solar System scaled-density dichotomy between
Mercury and the other rocky planets, we assume the lower-density population
has an Earth-like composition, i.e., $\mu_{comp}^{\mathrm{el}}=0.325$ CMF, and
the higher density population has a Mercury-like core mass fraction of
$\mu_{comp}^{\mathrm{ml}}=0.7$ CMF. The latter choice is further supported by
the fact that $\chi^{2}$ is minimized between the 5 Fe-rich planets of
Adibekyan et al. (2021) and our full rock grid when CMF$\sim$ 0.73. For a
sense of scale, the Hypatia catalog has Fe/Mg = $0.79\pm 0.18$ corresponding
to CMF = $0.29\pm 0.05$, meaning this population is $\sim 8\sigma$ above the
central CMF of Hypatia. Per the Adibekyan et al. (2021) sample, Fe-rich
planets represent 0.23 of the total population. For simplicity, we round this
up to a quarter and set the mixing weight of our Mercury-like population to
$w^{\mathrm{ml}}=0.25$. The corresponding mixing weight for the Earth-like
population is then $w^{\mathrm{el}}=0.75$.

For each sample, we randomly draw $N_{p}$ planet masses from
$M_{p}\sim\mathcal{U}[1,10]M_{\oplus}$. We then assume $N_{ml}=w_{ml}\times
N_{p}$ of these are Mercury-like and $N_{el}=w_{el}\times N_{p}$ are Earth-
like in composition. We then calculate the density corresponding to
$\mu_{comp}^{\mathrm{ml}}=0.7$ CMF and $\mu_{comp}^{\mathrm{el}}=0.325$ CMF at
mass for each Mercury-like and Earth-like planet, respectively, and then
calculate the corresponding planet radius, $R_{p}$. Finally, we resample each
mass and radius from the Gaussian distributions
$\sim\mathcal{N}(M_{p},\sigma_{M_{p}})$ and
$\sim\mathcal{N}(R_{p},\sigma_{R_{p}})$, respectively, to simulate
observational uncertainties. We consider three sets of mass-radius
uncertainties: (1) $\sigma_{M_{p}}=10\%$ and $\sigma_{R_{p}}=2.5\%$, (2)
$\sigma_{M_{p}}=6\%$ and $\sigma_{R_{p}}=2\%$, and (3) $\sigma_{M_{p}}=5\%$
and $\sigma_{R_{p}}=1\%$. For each set of mass-radius uncertainties, we
calculate $B_{21}$ and $n_{\sigma}$ as a function of $N_{p}$.

We find that the number of planets required to achieve a strong detection of
both populations is highly dependent on the average mass and radius
uncertainties. For the near best-case scenario of $\sigma_{M_{p}}=5\%$ and
$\sigma_{R_{p}}=1\%$, we estimate only $\sim 12$ planets are needed to detect
both populations with moderate support ($B_{21}>12$ and $n_{\sigma}>2.7$) and
$\sim 16$ for detection with strong support ($B_{21}>150$ and
$n_{\sigma}>3.6$). For the most conservative, but still precise, case of
$\sigma_{M_{p}}=10\%$ and $\sigma_{R_{p}}=2.5\%$, a moderate detection
requires approximately 117-118 planets, and a strong detection requires at
least $\sim$ 154 planets, i.e., an almost factor of 10 increase in $N_{p}$
from the near best-case scenario. For comparison, the average mass and radius
uncertainties for the Adibekyan et al. (2021) sample are $14\%$ and $5\%$,
respectively. As such, it is not surprising that they do not find a conclusive
and distinct population of Fe-rich planets as it would require significantly
more than 154 planets to do so at their average precisions. We also show our
model results for the synthetic 20-planet sample with $\sigma_{M_{p}}=5\%$ and
$\sigma_{R_{p}}=1\%$ in detail in Appendix 8 to demonstrate RhoPop’s ability
to recover known hyperparameters.

Figure 5: $B_{21}$ (left vertical axis) and $n_{\sigma}$ (right vertical axis)
as a function of $N_{p}$ for average observational uncertainties of
$\sigma_{M_{p}}=5\%$ and $\sigma_{R_{p}}=1\%$ (cyan), $\sigma_{M_{p}}=6\%$ and
$\sigma_{R_{p}}=2\%$ (purple), and $\sigma_{M_{p}}=10\%$ and
$\sigma_{R_{p}}=2.5\%$ (red).

##  5 Discussion

###  5.1 Comparisons with the nominally rocky planet zone

We compare our results for the Luque & Pallé (2022) (Sect. 3.1) and Adibekyan
et al. (2021) (Sect. 3.2) samples with the nominally rocky planet zone
(Unterborn et al., 2023). From Eq. 24 and 27 of Unterborn et al. (2023), we
estimate that the NRPZ covers a density ratio range of 0.85-1.1. We reiterate
that the NRPZ represents the expected 3$\sigma$ compositional range of
volatile-free and water-free rocky planets. The preferred single population
model for the Adibekyan et al. (2021) sample yields a 1$\sigma$ density ratio
range that is similar in width to the NRPZ, meaning the 3$\sigma$ range for
this sample is $\sim 3\times$ larger than the NRPZ. This result is not
surprising given the number of exceptions to this simplistic picture of planet
formation we outline in Sect. 1.

Our 1-population model for the Luque & Pallé (2022) sample finds a central
composition of 0.09 WMF with an intrinsic 1$\sigma$ scatter of
$\sigma_{\rho_{pop}}=0.21$. This leads to a 1$\sigma$ density range of $\sim
0.6-1.0$ which is inconsistent with the NRPZ. This result is not surprising
since we expected to find 2 populations a priori given the results of Luque &
Pallé (2022), although we did not inject this expectation into our models.
This tension with the NRPZ is reduced for the strongly preferred 2-population
model. We find a higher density population with a central composition of WMF =
0.02 and $\sigma_{\rho_{pop}}=0.02$. The 2% water mass fraction isocomposition
line corresponds to a density ratio of $\sim 0.95$ and the 3$\sigma$ range in
density ratios is $\sim 0.89-1.01$. This density ratio range falls entirely
within the NRPZ meaning that this population of planets can be explained via
volatile-free rocky compositions, i.e., we don’t need to invoke water to
explain the densities of these planets. More notable, however, is the fact
that the 3$\sigma$ range of this population is approximately half the width of
the NRPZ suggesting the compositions of these planets may be less diverse than
stellar abundances predict. We note the important caveat that the NRPZ is
estimated using Fe, Mg, and Si abundances from the Hypatia Catalog database
(Hinkel et al., 2014) and is thus strongly biased towards FGK stars as there
is a lack of M-type stars with abundance measurements, especially those with
Fe, Mg, and Si measurements. This suggests that the Luque & Pallé (2022)
sample has a smaller spread in the major rock-building elements than FGK
stars. In fact, M dwarfs that host planets in general may have a smaller
spread. Testing this possibility, however, requires more M-dwarf abundances of
Fe, Mg, and Si to build an M-dwarf NRPZ, i.e., an M-NRPZ. While TOI-561 b in
the Adibekyan et al. (2021) sample has a $\rho_{ratio}\sim 0.5$, consistent
with the lower-density population of Luque & Pallé (2022), the lack of
evidence for a water-rich population around FGK stars further suggests that
small-planet compositions and formation processes may be dependent on host-
star type.

###  5.2 The TRAPPIST-1 planets

The primary goal of Sect. 3.1 is to validate the results of Luque & Pallé
(2022) using their sample as published. We note, however, that the TRAPPIST-1
system accounts for roughly a quarter of the LP22 sample. Further, the seven
TRAPPIST-1 planets have some of the smallest error bars with an average
density uncertainty of 5.6% compared to the remainder of the sample which has
an average density uncertainty of 19%. Together this suggests our detection of
2-populations of M dwarf planets may be driven in large part by the TRAPPIST-1
system. To test this, we removed the TRAPPIST-1 planets from the LP22 sample
and reran the 1- and 2-population models. Without these seven planets, the
detection of 2-populations of planets is reduced from 4.1$\sigma$ to
1.8$\sigma$, meaning the TRAPPIST-1 planets do indeed play a large role in the
detection of 2-populations of small planets around M dwarfs. Further, for the
$N_{c}=2$ scenario, we find a median value of $\mu^{1}_{comp}=0.327$ CMF,
nearly identical to the CMF of Earth, in contrast to $\mu^{1}_{comp}=0.02$
WMF, when the TRAPPIST-1 planets are included. This begs the question as to
whether volatile-/water-free rocky planets around M dwarfs are generally
Earth-like in composition and the TRAPPIST-1 planets are anomalously iron-
depleted or wet. A larger sample of small M dwarf planets with density
precisions comparable to the TRAPPIST-1 planets and/or an M-NRPZ are needed to
fully investigate this question.

###  5.3 Alternative explanations for water worlds

We note that our water grid assumes condensed water/ice phases. The
equilibrium temperatures of the five LP22 planets in the lower density
population range from $\sim 1.5-4.2\times$ Earth’s. As such, it is likely
substantial fractions of their water exist in the vapor and super-critical
phases. Our choice of condensed water/ice then overestimates the density of
these water layers, and, therefore, our WMF estimates serve as an upper limit.
We chose to work in terms of condensed water/ice for the most direct
comparison with the interpretations made in Luque & Pallé (2022), but a future
iteration of RhoPop will incorporate the super-critical and vapor phases of
water and planet equilibrium temperature.

We also note that the population of planets interpreted as water-worlds in
Luque & Pallé (2022) and this work is not uniquely explained as such. Rogers
et al. (2023) show that these planets are equally well explained as sub-
Neptunes that have retained a small fraction of their primordial H/He
envelopes and provide some observational evidence to support this
interpretation. It is also plausible such planets have outer gaseous envelopes
comprised of species heavier than H/He (e.g., Piaulet et al., 2023). In short,
the source M dwarf planet density gap remains to be fully explored. In a
future update, we will expand RhoPop to account for these various scenarios
and investigate whether the data favors one scenario over the others.

###  5.4 Caveats to estimates on the minimum number of planets needed to
resolve a Fe-rich population

In Sect. 4, we estimated that approximately 154 planets with average
$\sigma_{M_{p}}=10\%$ and $\sigma_{R_{p}}=2.5\%$ are needed to resolve a
population of planets with Mercury-like compositions from a population with
nominally Earth-like compositions. Interestingly, there are 154 small planets
on the NASA Exoplanet Archive666as of 07 AUG 2023 with density uncertainties
less than 100%. The average and median M-R uncertainties of these planets,
however, are $\sigma_{M_{p}}=27\%$ and $\sigma_{R_{p}}=6\%$ and
$\sigma_{M_{p}}=19\%$ and $\sigma_{R_{p}}=5\%$, respectively. These
uncertainties are $\gtrsim 2\times$ the precisions required to resolve a Fe-
rich population with high confidence for $N_{p}=154$. As such, even if such a
population exists, we do not expect it to be resolvable with the current
sample of small planets. We note, however, that our analysis assumes there is
1 Fe-rich planet for every 3 nominally Earth-like planets ($w^{ml}=0.25$)
based on the work of Adibekyan et al. (2021). If Fe-rich planets are more
abundant relative to Earth-like planets than 1:3, then the required sample
size to resolve this population will be reduced and vice versa. Further, we
assume a constant bulk composition for the iron-rich population of
$\mu_{comp}^{ml}=0.7$. While we use Mercury and the 5 potentially Fe-rich
planets from Adibekyan et al. (2021) to justify this choice, the average CMF
can be higher or lower which will act to make this population more or less
easily resolvable at a given mass and radius precision, respectively. We will
consider a wider range of mixing weights and compositions for the Fe-rich
population in future work.

##  6 Conclusions

We present the open-source software RhoPop for identifying distinct
populations of small planets and use it to validate the recent results of
Luque & Pallé (2022) and Adibekyan et al. (2021). While we show that a
distinct population of Fe-rich planets is likely unresolvable with the current
sample of small planets, RhoPop is designed to provide the community with a
ready-to-go tool for identifying such populations as the number and precisions
of small planets continue to increase. The next release of RhoPop will include
pre-loaded compositional grids with atmosphere layers comprised of various
molecular species. We will use these additional grids to investigate whether
the population of planets interpreted by Luque & Pallé (2022) as 50% water
worlds are better explained as having an atmosphere. In future work, we will
apply RhoPop to the entire available sample of small planets and expand our
analysis on the minimum number of planets required to resolve a population of
iron-rich planets.

##  7 Acknowledgements

W.R.P’s effort is based upon work while serving at the National Science
Foundation and was initially funded by NSF grant EAR-1724693. Work by B.S.G.
was supported by the Thomas Jefferson Chair for Space Exploration endowment
from The Ohio State University. C.T.U acknowledges support under grants
NNX15AD53G and 80NSSC23K0267 from the National Aeronautics and Space
Administration (NASA). The results reported herein benefited from
collaborations and/or information exchange within NASA’s Nexus for Exoplanet
System Science (NExSS) research coordination network sponsored by NASA’s
Science Mission Directorate. This work is supported by the National Science
Foundation under Grant No. 2143400.

## References

  * Adibekyan et al. (2021) Adibekyan, V., Dorn, C., Sousa, S. G., et al. 2021, Science, 374, 330, doi: 10.1126/science.abg8794
  * Aitta (2012) Aitta, A. 2012, Icarus, 218, 967, doi: 10.1016/j.icarus.2012.01.007
  * Allard et al. (2000) Allard, F., Hauschildt, P. H., & Schwenke, D. 2000, ApJ, 540, 1005, doi: 10.1086/309366
  * Anderson & Ahrens (1994) Anderson, W. W., & Ahrens, T. J. 1994, J. Geophys. Res., 99, 4273, doi: 10.1029/93JB03158
  * Asahara et al. (2010) Asahara, Y., Hirose, K., Ohishi, Y., Hirao, N., & Murakami, M. 2010, Earth and Planetary Science Letters, 299, 474, doi: https://doi.org/10.1016/j.epsl.2010.09.037
  * Azevedo Silva et al. (2022) Azevedo Silva, T., Demangeon, O. D. S., Barros, S. C. C., et al. 2022, A&A, 657, A68, doi: 10.1051/0004-6361/202141520
  * Barros et al. (2022) Barros, S. C. C., Demangeon, O. D. S., Alibert, Y., et al. 2022, A&A, 665, A154, doi: 10.1051/0004-6361/202244293
  * Benneke & Seager (2013) Benneke, B., & Seager, S. 2013, ApJ, 778, 153, doi: 10.1088/0004-637X/778/2/153
  * Benz et al. (1988) Benz, W., Slattery, W. L., & Cameron, A. 1988, Icarus, 74, 516, doi: https://doi.org/10.1016/0019-1035(88)90118-2
  * Bond et al. (2010a) Bond, J. C., Lauretta, D. S., & O’Brien, D. P. 2010a, Icarus, 205, 321, doi: 10.1016/j.icarus.2009.07.037
  * Bond et al. (2010b) Bond, J. C., O’Brien, D. P., & Lauretta, D. S. 2010b, ApJ, 715, 1050, doi: 10.1088/0004-637X/715/2/1050
  * Bonomo et al. (2019) Bonomo, A. S., Zeng, L., Damasso, M., et al. 2019, Nature Astronomy, 3, 416, doi: 10.1038/s41550-018-0684-9
  * Bonsor et al. (2021) Bonsor, A., Jofré, P., Shorttle, O., et al. 2021, MNRAS, 503, 1877, doi: 10.1093/mnras/stab370
  * Cameron (1985) Cameron, A. 1985, Icarus, 64, 285, doi: https://doi.org/10.1016/0019-1035(85)90091-0
  * Chen & Kipping (2017) Chen, J., & Kipping, D. 2017, ApJ, 834, 17, doi: 10.3847/1538-4357/834/1/17
  * Connolly (2009) Connolly, J. A. D. 2009, Geochemistry, Geophysics, Geosystems, 10, doi: https://doi.org/10.1029/2009GC002540
  * Dorn et al. (2015) Dorn, C., Khan, A., Heng, K., et al. 2015, A&A, 577, A83, doi: 10.1051/0004-6361/201424915
  * Dumoulin et al. (2017) Dumoulin, C., Tobie, G., Verhoeven, O., Rosenblatt, P., & Rambaux, N. 2017, Journal of Geophysical Research (Planets), 122, 1338, doi: 10.1002/2016JE005249
  * Elser et al. (2012) Elser, S., Meyer, M. R., & Moore, B. 2012, Icarus, 221, 859, doi: 10.1016/j.icarus.2012.09.016
  * Fei et al. (1993) Fei, Y., Mao, H., & Hemley, R. J. 1993, The Journal of Chemical Physics, 99, 5369, doi: 10.1063/1.465980
  * Hinkel et al. (2014) Hinkel, N. R., Timmes, F. X., Young, P. A., Pagano, M. D., & Turnbull, M. C. 2014, AJ, 148, 54, doi: 10.1088/0004-6256/148/3/54
  * Johansen & Dorn (2022) Johansen, A., & Dorn, C. 2022, A&A, 662, A19, doi: 10.1051/0004-6361/202243480
  * Johnson et al. (2012) Johnson, T. V., Mousis, O., Lunine, J. I., & Madhusudhan, N. 2012, ApJ, 757, 192, doi: 10.1088/0004-637X/757/2/192
  * Journaux et al. (2020) Journaux, B., Brown, J. M., Pakhomova, A., et al. 2020, Journal of Geophysical Research: Planets, 125, e2019JE006176, doi: https://doi.org/10.1029/2019JE006176
  * Lodders (2003) Lodders, K. 2003, ApJ, 591, 1220, doi: 10.1086/375492
  * Luque & Pallé (2022) Luque, R., & Pallé, E. 2022, Science, 377, 1211, doi: 10.1126/science.abl7164
  * Neil et al. (2022) Neil, A. R., Liston, J., & Rogers, L. A. 2022, ApJ, 933, 63, doi: 10.3847/1538-4357/ac609b
  * Neil & Rogers (2020) Neil, A. R., & Rogers, L. A. 2020, The Astrophysical Journal, 891, 12, doi: 10.3847/1538-4357/ab6a92
  * Piaulet et al. (2023) Piaulet, C., Benneke, B., Almenara, J. M., et al. 2023, Nature Astronomy, 7, 206, doi: 10.1038/s41550-022-01835-4
  * Plotnykov & Valencia (2020) Plotnykov, M., & Valencia, D. 2020, Monthly Notices of the Royal Astronomical Society, 499, 932, doi: 10.1093/mnras/staa2615
  * Rodríguez Martínez et al. (2022) Rodríguez Martínez, R., Gaudi, B. S., Schulze, J. G., et al. 2022, arXiv e-prints, arXiv:2208.07883, doi: 10.48550/arXiv.2208.07883
  * Rogers et al. (2023) Rogers, J. G., Schlichting, H. E., & Owen, J. E. 2023, ApJ, 947, L19, doi: 10.3847/2041-8213/acc86f
  * Rogers (2015) Rogers, L. A. 2015, ApJ, 801, 41, doi: 10.1088/0004-637X/801/1/41
  * Santerne et al. (2018) Santerne, A., Brugger, B., Armstrong, D. J., et al. 2018, Nature Astronomy, 2, 393, doi: 10.1038/s41550-018-0420-5
  * Schulze et al. (2021) Schulze, J. G., Wang, J., Johnson, J. A., et al. 2021, The Planetary Science Journal, 2, 113, doi: 10.3847/psj/abcaa8
  * Scora et al. (2020) Scora, J., Valencia, D., Morbidelli, A., & Jacobson, S. 2020, MNRAS, 493, 4910, doi: 10.1093/mnras/staa568
  * Souto et al. (2022) Souto, D., Cunha, K., Smith, V. V., et al. 2022, ApJ, 927, 123, doi: 10.3847/1538-4357/ac4891
  * Speagle (2020) Speagle, J. S. 2020, MNRAS, 493, 3132, doi: 10.1093/mnras/staa278
  * Stixrude & Lithgow-Bertelloni (2011) Stixrude, L., & Lithgow-Bertelloni, C. 2011, Geophysical Journal International, 184, 1180, doi: 10.1111/j.1365-246X.2010.04890.x
  * Thiabaud et al. (2015) Thiabaud, A., Marboeuf, U., Alibert, Y., Leya, I., & Mezger, K. 2015, A&A, 580, A30, doi: 10.1051/0004-6361/201525963
  * Unterborn et al. (2023) Unterborn, C. T., Desch, S. J., Haldemann, J., et al. 2023, ApJ, 944, 42, doi: 10.3847/1538-4357/acaa3b
  * Unterborn et al. (2016) Unterborn, C. T., Dismukes, E. E., & Panero, W. R. 2016, ApJ, 819, 32, doi: 10.3847/0004-637X/819/1/32
  * Wang et al. (2019) Wang, H. S., Lineweaver, C. H., & Ireland, T. R. 2019, Icarus, 328, 287, doi: 10.1016/j.icarus.2019.03.018
  * Wanke & Dreibus (1994) Wanke, H., & Dreibus, G. 1994, Philosophical Transactions of the Royal Society of London Series A, 349, 285, doi: 10.1098/rsta.1994.0132
  * Wurm et al. (2013) Wurm, G., Trieloff, M., & Rauer, H. 2013, ApJ, 769, 78, doi: 10.1088/0004-637X/769/1/78
  * Yoshizaki & McDonough (2020) Yoshizaki, T., & McDonough, W. F. 2020, Geochim. Cosmochim. Acta, 273, 137, doi: 10.1016/j.gca.2020.01.011
  * Zharkov (1983) Zharkov, V. N. 1983, Moon and Planets, 29, 139, doi: 10.1007/BF00928322

## Appendix

##  8 Synthetic planet sample example

Here we show the detailed results of RhoPop applied to a synthetic 20-planet
sample where the hyperparameters are known a priori to demonstrate their
successful recovery. To reiterate, we sample 5 planets from the Mercury-like
population with $\mu_{comp}^{ml}=0.7$ CMF and intrinsic scatter
$\sigma_{\rho_{ratio,pop}^{1}}$ = 0.0. Similarly, we sample 15 planets from an
Earth-like population with $\mu_{comp}^{el}=0.325$ CMF and intrinsic scatter
$\sigma_{\rho_{ratio,pop}^{2}}$ = 0.0. This gives $w^{el}=0.75$ and
$w^{ml}=1-w^{el}=0.25$. We use the same hyperparameter priors given in Table 2
and assign 5% uncertainties in $M_{\mathrm{ob}}$ and 1% in $R_{\mathrm{ob}}$
to each planet.

We apply the 1- and 2-population RhoPop models to this sample and find,
unsurprisingly, that the 2-population model is favored over the 1-population
model at the $>4\sigma$ level. The results are summarized in Table 5 and
visualized in Fig. 6. The 1-population model is naturally a poor fit for this
sample and recovers a mean composition of CMF = 0.42 which is $\simeq
w^{el}\times\mu_{comp}^{el}+w^{ml}\times\mu_{comp}^{ml}$. None of the true
hyperparameter values are contained within the 68% CIs of the estimated
$\mu_{comp}$ for the 1-population model. The 2-population model, however,
recovers accurate values for both the Earth-like and Mercury-like populations.
RhoPop finds $\mu_{comp}^{1}=0.71^{+0.03}_{-0.04}$ CMF corresponding to the
Mercury-like population we fixed to have $\mu_{comp}^{ml}=0.7$ CMF. Similarly,
RhoPop finds a lower density population with
$\mu_{comp}^{2}=0.32^{+0.02}_{-0.02}$ CMF corresponding to the Earth-like
population we fixed to have $\mu_{comp}^{el}=0.325$ CMF. Further, this model
finds $w^{1}=0.27^{-0.09}_{+0.1}$ which is consistent with the injected value
of $w^{ml}=0.25$. We show the hyperparameter posteriors for the 1- and
2-population models in Fig. 7 and Fig. 8, respectively.

| $N_{c}=1$ | $N_{c}=2$
---|---|---
Param | Med | 68% CI | Med | 68% CI
$\mu_{comp}^{1}$ | 0.44 CMF | [0.39 CMF, 0.48 CMF] | 0.71 CMF | [0.67 CMF, 0.74 CMF]
$\sigma_{\rho_{ratio,pop}^{1}}$ | 0.20 | [0.16, 0.24] | 0.06 | [0.02, 0.17]
$w^{1}$ | 1 |  | 0.27 | [0.18, 0.37]
$\mu_{comp}^{2}$ |  |  | 0.32 CMF | [0.30 CMF, 0.34 CMF]
$\sigma_{\rho_{ratio,pop}^{2}}$ |  |  | 0.02 | [0.005, 0.04]
$w^{2}$ |  |  |  1-$w^{1}$ |
|  | $\ln\mathcal{Z}_{1}=-3.34\pm 0.20$ |  | $\ln\mathcal{Z}_{2}=3.79\pm 0.27$

$B_{21}$ |  $p$-value | Sigma | Interpretation
---|---|---|---
1250 |  2.8$\times 10^{-5}$ |  4.2$\sigma$ |  Strong evidence for two populations

Table 5: Summary of parameter posteriors, model evidence, and model selection
parameters for our 1- and 2-population model for a synthetic sample of 20
planets with Mercury-like and Earth-like subpopulations. CI = credible
interval

.

Figure 6: Summary of results for a synthetic 20-planet sample of planets with
two-populations: a Mercury-like one with $\mu_{comp}^{ml}=0.7$ CMF and
$w^{ml}=0.25$ and an Earth-like population with $\mu^{el}_{comp}=0.325$ CMF
and $w^{el}=0.75$.

Figure 7: hyperparameter posterior corner plot for the 1-population model
applied to our synthetic sample. The +/- values correspond to the 95% credible

intervals. This figure was generated with the dynesty dynamic nested sampling
software package (Speagle, 2020).

Figure 8: hyperparameter posterior corner plot for the 2-population model
applied to our synthetic sample. The +/- values correspond to the 95% credible

intervals. This figure was generated with the dynesty dynamic nested sampling
software package (Speagle, 2020).