Starlight

A prior the survey measures through its own noise

Whether a weak periodic signal is a planet depends on how large planets are expected to be, and the best source for that expectation is the population the surveys themselves have measured. But the population's small end is built from the marginal detections whose reality depends on the prior. Breaking the circle means measuring the population and the planets together, with the noise and the threshold written into the model — and doing it that way changes the answers by factors of two where they matter most.

Assumes Periodograms, Reflex velocity and Occurrence rates.

A radial-velocity signal near the noise can be judged in two ways. The first asks how often noise alone would produce a peak as tall, which gives a false-alarm probability that depends only on the noise. The second asks how much more probable the data are with a planet than without one, which gives a Bayes factor — and the Bayes factor depends on what amplitude a planet would be expected to have. A prior that allows planets of any size up to large amplitudes spreads its belief thinly and penalises a small signal; a prior concentrated at small amplitudes favours it. The same data can give odds of several to one either way.

The obvious remedy is to make the prior a measurement. The distribution of planet masses and periods has been measured by the radial-velocity surveys themselves, and converted to semi-amplitudes for a given star it is exactly the prior the Bayes factor needs: how large are planets, around stars like this, actually found to be? The difficulty is that the population’s small end — the part that matters for a marginal signal — was measured from marginal signals. The prior for judging a weak detection comes from weak detections, whose reality was judged with some earlier prior. This essay is about what that circle does to the numbers and how it is broken.

An amplitude is biased before it is selected

The first effect is present even before a threshold is applied, and it comes from what an amplitude is. At a known period, a radial-velocity signal is fitted as a sine and a cosine, each with its own noise; the amplitude is the length of the vector they make. A length cannot be negative, and noise added to both components makes the vector longer on average, not shorter.

What an amplitude looks like when it is measured. The distribution of the amplitude fitted at a known period, in units of its noise, for true amplitudes of 0, 1, 2, 4 (a Rice distribution). An amplitude is the length of a vector with two noisy components, so it cannot be negative and its noise does not cancel: a true amplitude of zero is measured, on average, as 1.25, and one of 1 as 1.55. Only well above the noise, at 4, is the mean close to the truth (4.13). Every amplitude measured near the noise is therefore biased upward before any threshold is applied — and a threshold, which keeps the high draws and discards the low, biases the kept ones further.
Fig. 1 The distribution of the amplitude fitted at a known period, in units of its noise, for true amplitudes of 0, 1, 2 and 4 (Rice distributions). A true amplitude of zero is measured on average as 1.25; one of 1 as 1.55; one of 4 as 4.13. Amplitudes near the noise are biased upward before any threshold is applied.

A planet with no signal at all — a star with nothing orbiting at that period — is measured, on average, with an amplitude of 1.25 times the noise. A planet with an amplitude equal to the noise is measured as 1.55. The distributions are Rice distributions, the same that describe the amplitude of any sinusoid measured through Gaussian noise, and they approach the truth only when the signal is several times the noise. Every catalogue of amplitudes contains this bias for its weakest entries, and it is quite separate from the question of whether those entries are planets.

The same bias in every length

The upward bias of a noisy amplitude is not peculiar to radial velocities. It afflicts every quantity that is the length of a vector estimated from noisy components. An orbit’s eccentricity is one: it is the length of a vector whose two components are fitted from the data, and an eccentricity fitted to noisy data never comes back zero, even for a perfectly circular orbit — the classic form of the problem, noticed in spectroscopic binaries half a century ago. The degree of polarisation of faint light is another, and so is the amplitude of any periodic signal found in a time series, from a pulsating star’s to the reflex motion of a planet’s host.

The remedies are the same in every case. A single measurement can be corrected only with a prior on the true values, and the prior is only as good as the population it comes from; a population can be estimated without bias only by modelling the measurement process, not by averaging the measured values. The eccentricity distribution of exoplanets has been revised for exactly this reason: the fraction of apparently eccentric orbits among planets with small amplitudes was inflated by the same arithmetic as their amplitudes, and hierarchical fits that model the bias find the small planets’ orbits closer to circular than the raw catalogues suggest.

The threshold keeps the lucky ones

A detection threshold selects the upper part of each distribution. A planet whose measured amplitude clears the threshold is, near the threshold, a planet whose noise happened to lengthen its amplitude; the ones whose noise shortened it fell below and were not counted. Since the population contains more small planets than large ones, more are promoted across the threshold than demoted, and the detections just above it are systematically too large — the same arithmetic that inflates marginal transit depths, and every flux-limited count near its limit.

Amplitudes measured near the noise, against the amplitudes that were there. A simulated survey of 6000 stars, each with one planet whose true radial-velocity semi-amplitude K is drawn from a population rising towards small amplitudes (dN/d ln K ∝ K^−0.6, from 0.2 to 20 m/s), each measured with an amplitude uncertainty of 0.6 m/s; 1500 are drawn. Across: the true amplitude; up: the amplitude fitted at the planet's period; both logarithmic. The dashed diagonal is a perfect measurement; the horizontal line is the detection threshold, 4 times the uncertainty. Points above it are the detections (filled). Well above the threshold they sit on the diagonal; near it they sit above it — the detections within one unit of the threshold have measured amplitudes 1.33 times their true ones on average — and below it the measured amplitudes stop falling at all, because an amplitude is a length and noise alone gives it about 0.75 m/s.
Fig. 2 A simulated survey of 6,000 stars, each with a planet whose true semi-amplitude is drawn from a population rising towards small amplitudes (dN/d ln K ∝ K^−0.6), measured with an uncertainty of 0.6 m/s. Across: true amplitude; up: measured; both logarithmic. Detections (filled) lie above the 4σ threshold (dashed). Near the threshold they sit above the diagonal — on average 1.33 times their true amplitude — and below it the measured amplitudes level off at the noise.

In the simulated survey, the planets detected within one unit of the threshold have measured amplitudes a third larger than their true ones, on average. Since the minimum mass of a planet is proportional to its semi-amplitude, a minimum mass that is itself a lower bound is, for a marginal detection, inflated by the same third before the unknown inclination is even considered. Well above the threshold the detections sit on the diagonal and the bias disappears. Below it the measured amplitudes level off near the noise’s own contribution, whatever the true amplitude, and the undetected planets are indistinguishable from stars without planets.

What a population prior does to one planet

If the population were known, each marginal planet’s amplitude could be corrected by asking what true amplitudes, weighted by how common they are, would produce the measured one. That is the posterior mean under the population prior.

The amplitude a population prior returns for each measured one. The posterior mean of a planet's true amplitude given its measured amplitude, both in units of the noise, when the prior is the population's own amplitude distribution, dN/d ln K ∝ K^−β, for β = 0, 0.6, 1.2. The dashed diagonal is the measured value taken at face value. At the detection threshold of 4, the posterior means are 3.52, 3.25, 2.88: the steeper the population — the more small planets there are to have been scattered up — the more a marginal measurement is pulled down. Well above the noise the curves meet the diagonal and the prior stops mattering. A population prior does to a marginal detection what the survey's own statistics say should be done to it; the difficulty is that the population is what the survey is trying to measure.
Fig. 3 The posterior mean of a planet’s true amplitude against its measured amplitude, in units of the noise, under population priors dN/d ln K ∝ K^−β for β = 0, 0.6 and 1.2. At the 4σ threshold the posterior means are 3.52, 3.25 and 2.88: the steeper the population, the more a marginal measurement is pulled down. Well above the noise the curves meet the diagonal.

The correction is a shrinkage: the measured amplitude is pulled towards the part of the population where true amplitudes are common, which, for a population rising towards small planets, is downward. At the threshold, a measured amplitude of four times the noise corresponds to a true amplitude of about 3.3 for the population drawn — a fifth smaller — and to 2.9 for a steeper population. Well above the noise the prior no longer matters and the measurement stands on its own. This is the individual-planet side of the population prior: it does to each marginal detection what the survey’s statistics say should be done to it. It also shows why the circle matters. The correction depends on the steepness of the population, and the steepness is measured from the detections being corrected.

One marginal planet, corrected

The size of the correction for a single planet is easiest to see with numbers attached. Suppose a star’s radial velocities, a hundred measurements with a metre per second of scatter each, show a signal at some period with a fitted semi-amplitude of 0.6 metres per second, where the amplitude’s uncertainty is about 0.14 — a signal at a little over four times its noise, just past a typical threshold. Taken at face value, around a Sun-like star at a period of a few weeks, that is a minimum mass of a couple of Earth masses. Under a population prior of the steepness drawn, the posterior mean amplitude is about a fifth lower, and so is the minimum mass. The planet is still probably there — its evidence is enough to clear the threshold — but its most probable mass is smaller than the number quoted, and its uncertainty, once the prior’s own uncertainty is included, is larger.

For a single planet that is a modest adjustment. For a catalogue it is not, because the adjustment is always in the same direction. A list of marginal detections quoted at their fitted amplitudes is a list of planets all somewhat too massive, and a mass–radius relation, an occurrence rate or a comparison with formation models built on it inherits the shift. The correction is also what changes how many planets a star appears to have: a system whose third planet sits at the threshold looks like a two-planet system in some analyses and a three-planet one in others, and the population of multiplicities measured from such systems depends on how marginal members were treated.

The slope of a population read from its detections

The natural way to measure the population is to fit a power law to the detected planets’ amplitudes, above the threshold where the survey is complete. It is also wrong, in a predictable direction.

The population slope read from detections, naively and by forward modelling. The power-law index β of the amplitude distribution recovered from the detected planets of simulated surveys (threshold 4σ, σ = 0.6 m/s), against the index that was put in. Open circles: a power law fitted directly to the measured amplitudes of the detections, as if they were exact and complete above the threshold — 0.00 for 0, 0.40 for 0.3, 0.72 for 0.6, 1.04 for 0.9, 1.49 for 1.2. Filled: the index for which a simulated survey, with the same noise and the same threshold, produces detections whose amplitudes are distributed like the observed ones — 0.00, 0.34, 0.64, 0.86, 1.22. The naive fit comes out steeper than the truth, and more so the steeper the population: planets scattered up from below the threshold crowd the amplitudes just above it, so the detections look richer in small planets than the population is. Forward modelling recovers the input to within its statistical error, because it puts the noise and the selection into the model instead of pretending they are absent.
Fig. 4 The power-law index recovered from the detections of simulated surveys (threshold 4σ, σ = 0.6 m/s) against the index put in. Open circles: a law fitted directly to the measured amplitudes above the threshold — 0.40 for 0.3, 0.72 for 0.6, 1.04 for 0.9, 1.49 for 1.2. Filled: the index for which a simulated survey with the same noise and threshold produces detections distributed like the observed ones — 0.34, 0.64, 0.86, 1.22.

The direct fit comes out steeper than the truth, and the error grows with the steepness. Planets scattered up from below the threshold crowd the amplitudes just above it, so the detected sample looks richer in small planets than the population really is, and the fitted law rises too fast towards the threshold. For a flat population the effect vanishes, because there are as many planets above the threshold to be scattered down as below to be scattered up; for the steep populations actually observed, it is a third of the index or more.

The repair is to fit not the detections but the survey. For each trial population, simulate the survey — the same stars, the same noise, the same threshold, the same fitting — and compare the simulated detections with the real ones; the population whose simulated detections look like the observed ones is the estimate. The forward model puts the noise and the selection where they belong, inside the model, instead of treating the measured amplitudes as exact and the sample as complete. In the simulated surveys it recovers the input index to within its statistical error at every steepness drawn.

The count below the threshold

The quantity surveys most want is not the slope but what the slope implies for planets too small to detect — Earth-mass planets, whose semi-amplitudes around Sun-like stars are a few tenths of a metre per second, below the reach of all but the most precise instruments.

How far off the count below the threshold comes out. The fraction of stars with a planet of semi-amplitude between 0.3 and 1 m/s — all of them below the detection threshold of 2.4 m/s — as predicted by extrapolating a power law from the planets detected at more than twice the threshold, divided by the true fraction, on a logarithmic scale, for simulated surveys with population indices β = 0.3, 0.6, 0.9, 1.2. With the law fitted naively to the detected amplitudes, the predictions are 1.32, 1.11, 2.31, 2.37 times the truth; with the law found by forward modelling, 1.13, 0.89, 1.25, 1.15. A small error in the slope, measured over the factor of two or three in amplitude where planets are detected, is multiplied by the lever arm of an extrapolation over a further factor of three to eight, and the planets a survey most wants to count — the ones smaller than it can see — are where its method matters most.
Fig. 5 The fraction of stars with a planet of semi-amplitude 0.3–1 m/s, all below the 2.4 m/s threshold, predicted by extrapolating a power law from planets detected at more than twice the threshold, divided by the true fraction, for populations with β = 0.3, 0.6, 0.9 and 1.2. Naive fits give 1.32, 1.11, 2.31 and 2.37 times the truth; forward-modelled ones 1.13, 0.89, 1.25 and 1.15.

Extrapolating the naive law over the factor of three to eight in amplitude between the detected planets and the band of interest multiplies its error in the slope. For the steeper populations the naive prediction of how many small planets there are is more than twice the truth; the forward-modelled one is within a quarter. That is the planets that were not seen counted by extrapolation, and the method of the extrapolation matters more than the number of detections it rests on — on top of the part of any such rate that is a definition, the choice of which box of masses and periods is being counted. It is also exactly the prior that a Bayes factor for a new marginal signal would use. A survey that measured its population naively would assign every new small signal a prior that expects too many small planets, and would be too ready to believe in them.

Why the naive estimate lasted

The direct fit to detected amplitudes was the standard method for years, and not through carelessness. Its errors are invisible in the data it uses. A power law fitted to detections describes the detections well; residuals from the fit look like noise; the fitted slope has a small statistical error when there are many detections. Nothing inside the analysis signals that the slope is biased, because the bias lives in the relation between the detections and the planets that were not detected, and the naive method never models the latter. Only a simulation of the survey, with known inputs, shows the gap — the same lesson the luminosity functions of galaxies taught, where the scatter of measured brightnesses about true ones steepened the bright end of the count until it was modelled.

The forward model’s cost is that it must be told everything about the survey: which stars were observed, how often, with what noise, with what search and what threshold, and how candidates were vetted. A published catalogue of detections does not contain that information, which is why population analyses by groups outside a survey have tended to use the naive method, and why the most reliable population estimates have come from the surveys themselves, which could simulate their own pipelines. That dependence is itself a reason for surveys to publish their non-detections and their sensitivity alongside their planets.

Measuring the population and the planets together

The circle is broken by not separating the two steps. A hierarchical model treats the population’s parameters — its slope, its normalisation, its dependence on period and on the star — as unknowns alongside every planet’s amplitude, period and phase, and fits them all to all the data at once: every star’s radial velocities, including the stars with no detection, whose velocities say that any planet they have is small. The population informs each planet’s amplitude through the prior; each planet’s data inform the population through its likelihood; and nothing is thresholded, because a marginal signal enters with exactly the weight its data give it, neither promoted to a detection nor discarded.

The result is a population that is consistent with the individual planets and a set of individual planets consistent with the population, and the marginal signals are handled by their probabilities rather than by a cut. It is computationally heavy — tens of thousands of stars, each with its own period search — and it depends on modelling the noise correctly for every star, which is its real weak point. But it answers the question the Bayes factor left open: the prior for a marginal planet is the population, and the population is measured, not guessed, with the marginal planets contributing only what they are worth.

What the simulations leave out

The simulated surveys make the problem cleaner than real ones in three ways. Each star has exactly one planet at a known period, so the amplitude is fitted at the right period; a real signal has to be found first, and a search over many periods raises the threshold and adds false detections whose amplitudes are pure noise, which bias the population further towards small amplitudes. The noise is white and identical for every star; real radial-velocity noise is dominated at the metre-per-second level by the star’s own surface — granulation, oscillations, spots rotating across the disc — which is correlated in time, different for every star, and can imitate a planet at the star’s rotation period and its harmonics. And the population is a single power law in amplitude, when the real one depends on period and on the star’s mass and composition, so that a survey’s population prior for one star is a slice through a many-dimensional distribution.

None of these changes the direction of the effects drawn — amplitudes near the noise are inflated, detections near a threshold more so, and a law fitted to detections is too steep — but each changes their size, and the correlated stellar noise changes it most, because a signal whose noise is not the assumed noise has a detection significance that is not the computed one.

Why the number at the bottom matters

The occurrence of Earth-like planets around Sun-like stars is one of the handful of numbers that determine how to design the next generation of instruments and missions. Direct imaging of such planets requires knowing how many there are within reach, and radial-velocity campaigns designed to find them need to know how many stars to observe and for how long. The number currently comes partly from transit surveys, with their own extrapolations, and partly from radial-velocity surveys whose smallest detections sit at their thresholds. The figures say that the method of extrapolation can change that number by a factor of two, independently of the data, and that forward-modelled estimates are the ones to use — which is the practice the recent analyses have adopted, and the reason their answers differ from the earliest ones.

Still open: a prior for planets nobody has measured

The instruments now reaching ten centimetres per second in principle are aimed at planets whose amplitudes are below any population anyone has measured by radial velocity: Earth analogues around Sun-like stars, at nine centimetres per second. For those, there is no population to use as a prior except one extrapolated from larger planets — or one borrowed from transit surveys, which measure radii rather than masses and relate to amplitudes only through a mass–radius relation with its own scatter. The first detections at that level will be judged with a prior that is itself an extrapolation, and the population measured from them will inherit whatever that prior assumed. Whether the hierarchical approach can be carried into a regime where every signal is marginal — where the population and the planets must be measured from data none of which, alone, would count as a detection — is the question the next decade of radial-velocity surveys will have to answer.

About the same objects

Not linked from either essay — found by the objects both name.

What links here

Essays that link to this one from their own argument.

The objects this essay names

Each one links to every other essay that touches it.

Bayes factorEddington biasHierarchical modelOccurrence ratePeriodogramPriorRadial-velocitySelection effect