The star that moves: radial velocity

In January 1995 Michel Mayor and Didier Queloz computed a first orbit for the star 51 Pegasi. Its velocity along our line of sight rose and fell by about 59 metres per second every 4.23 days. They had been watching it only since the previous September, as one of 142 Sun-like stars in a survey with the ELODIE spectrograph at Haute-Provence, and the period was so short that they did not trust it. They tested the orbit over eight consecutive nights in July and eight more in September, checked whether starspots or pulsations could fake the signal, and found that the shapes of the spectral lines did not change with the orbital period. The circular solution from 35 measurements gave $P=4.2293\pm0.0011$ days, $K=59\pm3$ m/s and a minimum mass of $0.47$ Jupiter masses. The paper reached Nature on 29 August and appeared on 23 November 1995, after two other teams had confirmed the signal in mid-October (Mayor & Queloz 1995).

That number, 59 m/s, frames this chapter. A Jupiter-mass planet in a four-day orbit pulls its star around at the speed of a fast car. The Earth pulls the Sun around at about 9 cm/s, roughly 650 times less: a slow crawl. Thirty years of work have taken the best instruments from about 15 m/s to well under 1 m/s. Along the way the problem changed character. A modern Doppler spectrograph can repeat itself within a night at the level of 10 cm/s, but the star’s surface boils, rotates and goes through magnetic cycles that move its spectral lines by metres per second. Finding an Earth twin by radial velocity (RV) is now less a measurement problem than a separation problem: which part of a tiny spectral change is the whole star moving, and which part is the star changing?

The computation also has an unusual shape. An exposure on a modern instrument records millions of detector pixels. Extraction compresses them into one velocity. A planet search then works on a few hundred such velocities. Most of the information that could distinguish a planet from a stellar surface effect is either kept or thrown away at that first compression, and most of the recent progress, and most of the remaining room, lies there.

From Struve to 51 Pegasi

A proposal ahead of its time

In October 1952 Otto Struve published a two-page note in The Observatory proposing “a project of high-precision stellar radial velocity work”. His argument was that nature need not copy the Solar System. A planet of Jupiter’s mass at 1/50 au would orbit in about a day and move its star by about ±0.2 km/s, which the coudé spectrographs of the time might detect. Such a planet would also eclipse its star by about 0.02 magnitudes if the orbit were edge-on (Struve 1952). This was a calculation, not a discovery. Its insight was to look where the signal would be largest rather than where the Solar System says planets ought to be. Forty-three years later 51 Pegasi b turned out to be almost exactly the kind of planet Struve had described.

Putting the reference into the starlight

The first sustained search at the metres-per-second level came from Bruce Campbell and Gordon Walker in Canada. Their trick was to pass the starlight through a cell of hydrogen fluoride gas before it entered the spectrograph. The gas imprinted its own absorption lines on the stellar spectrum, and because both sets of lines then travelled through the same optics, drifts of the spectrograph moved them together. The velocity of the star relative to the gas could be measured even though the spectrograph itself was unstable. Gerhard Herzberg had advised on the choice of gas. Working at the 3.6-m Canada–France–Hawaii Telescope from 1980, Campbell, Walker and Stephenson Yang reported in 1988 a mean external error of 13 m/s for sixteen stars over six years. Several stars showed slow trends, and γ Cephei showed a residual signal of about 25 m/s with a 2.7-year period on top of its binary motion. They did not have the baselines or the activity diagnostics to call these signals planets, and they did not (Campbell, Walker & Yang 1988).

A year later David Latham, Tsevi Mazeh, Robert Stefanik, Mayor and Gilbert Burki found an unseen companion to HD 114762 while monitoring stars intended as velocity standards. Their 280 observations from four instruments had a combined residual scatter of 0.42 km/s, comparable to the signal itself, yet a period search found a coherent 84-day modulation, and an orbital fit gave a minimum mass of about 11 Jupiter masses. The title called it “a probable brown dwarf” (Latham et al. 1989). The detection-theory lesson is the one developed in Chapter 1: a coherent signal smaller than the error on any single measurement becomes convincing when many measurements are combined along the right template.

From fitting the instability to removing it

Two schools then developed. Geoffrey Marcy, R. Paul Butler and colleagues replaced hydrogen fluoride with iodine and modelled, exposure by exposure, how the spectrograph blurred the iodine lines. Their 1996 paper reached about 3 m/s from minutes to a year with the Lick and Keck spectrographs, which were general-purpose instruments not built for stability (Butler et al. 1996). The Geneva school instead built instruments so stable that there was little left to model. HARPS, on ESO’s 3.6-m telescope at La Silla, saw first light on 11 February 2003 with resolving power $R\simeq115{,}000$, coverage from 380 to 690 nm and a simultaneous thorium–argon reference (Mayor et al. 2003). ESPRESSO carried the approach to the 8-m VLT, with first light on 27 November 2017 and better than 25 cm/s within a night during commissioning (Pepe et al. 2021). NEID on the WIYN 3.5-m telescope and KPF at Keck continue the line. KPF has been offline since a vacuum-pump failure announced in April 2026, a reminder that a published precision is not a promise of current operation (Keck Observatory 2026).

The figure gathers the numbers the papers quote. Each one is a different quantity, and the chapter will keep returning to that point: an external scatter over years, an instrument accuracy, a photon-limited error, a within-night repeatability and a formal single-exposure error are not interchangeable. The most recent points sit within a factor of two or three of the Earth’s signal. None of them is a demonstration that an Earth twin can be found.

The signal

Reflex motion and what the mass means

A star and a planet orbit their common centre of mass. For non-interacting planets on Keplerian orbits the star’s line-of-sight velocity at time $t_i$ is

where $\gamma_j$ is a zero point for instrument or era $j$, $\dot\gamma$ a linear trend, and each planet $p$ has semi-amplitude $K_p$, eccentricity $e_p$ and argument of periastron $\omega_p$. The true anomaly $\nu$ follows from the mean anomaly through Kepler’s equation,

with $T_p$ the time of periastron. The semi-amplitude is

which for a planet much lighter than its star becomes

(Cumming 2004, Eq. 30). Three numbers from this formula set the scale. Jupiter, with a period of 11.86 years, moves the Sun by $28.4\times11.86^{-1/3}\approx12.5$ m/s. The Earth, with 1/318 of Jupiter’s mass and a one-year period, gives $28.4/318\approx0.089$ m/s. A planet of $0.47M_J$ in a 4.23-day orbit around a solar-mass star gives $28.4\times0.47\times(4.23/365.25)^{-1/3}\approx59$ m/s, the 51 Pegasi signal. The $P^{-1/3}$ factor and the $M_\star^{-2/3}$ factor explain why RV surveys found close-in giants and planets around small stars first.

Only $m_p\sin i$ appears, where $i$ is the inclination of the orbit to the sky. The data measure the mass function $f(m)=PK^3(1-e^2)^{3/2}/(2\pi G)=m_p^3\sin^3 i/(M_\star+m_p)^2$, so an RV planet comes with a minimum mass unless something else fixes $i$: a transit (Chapter 30), astrometry (Chapter 34) or dynamical interactions (Chapter 31). Interacting planets in resonance can also break the independent sum of Keplerians, and then an N-body model is needed.

A shift of one part in three hundred million

A Doppler shift is a translation in the logarithm of wavelength. With $x=\ln\lambda$ and a small velocity $v$, $\delta x\simeq v/c$, so 1 m/s is a fractional shift of $3.3\times10^{-9}$. At $R=115{,}000$ one resolution element spans $c/R\approx2.6$ km/s, so a 1 m/s shift moves each line by about 1/2,600 of its own width (a derived number). No single line can be centred that well. The measurement works only because thousands of lines are combined and because the photons from each line are very many.

Precise work also has to be careful about what “velocity” means. The observatory moves around the Sun at up to 30 km/s and spins with the Earth, so the measured shift is dominated by our own motion. Jason Wright and Jason Eastman showed that removing it at the 1 cm/s level requires multiplying redshifts, $1+z_{\rm corr}=(1+z_{\rm meas})(1+z_B)$, rather than adding velocities, together with accurate times, observatory coordinates, stellar proper motion and photon-weighted exposure midpoints (Wright & Eastman 2014). The parallel with Solar-System searches is direct. Just as the apparent motion of a distant body is mostly the reflection of the Earth’s orbit, the raw velocity of a star is mostly the reflection of the Earth’s velocity. In both cases the observer’s motion must be removed exactly before the faint signal of interest can be searched for, and in both cases errors in that correction come back with a period of one year.

Where the Doppler information lives

A useful forward model for the counts $m_{iq}$ in detector pixel $q$ of exposure $i$ is

where $S$ is the stellar spectrum shifted by $\delta_i$, $T_i$ the transmission of the Earth’s atmosphere, $I$ an optional absorption cell, $C_i$ the continuum and throughput, $L_i$ the instrumental line-spread function and $b$ the background. The order matters: the star and the atmosphere multiply before the instrument blurs them, and deep telluric absorption removes photons that no later correction can restore (Artigau et al. 2022; Zechmeister et al. 2018; Bedell et al. 2019).

Let $g_q=\partial m_q/\partial v$ be the change in pixel $q$ per unit velocity and $\sigma_q$ its noise. Linearising around a template, the best velocity shift and its limiting error are

This is the Connes–Bouchy information argument as developed by Artigau et al. (2022). It says something simple and important: velocity information lives on the sloping wings of lines. Flat continuum carries none, and neither does the flat bottom of a saturated line. More lines help in proportion to their slopes and their photons, not their number. Rapid rotation and low resolving power flatten the slopes and remove information.

The same algebra handles nuisance parameters. Suppose continuum tilts, line-width changes or other deformations are allowed through a design matrix $B$, with pixel covariance $\Sigma$. Profiling them out leaves the Doppler information

This Schur complement, a standard result of linear estimation, is the most important equation in the chapter. If an allowed deformation of the spectrum is parallel to the Doppler derivative $g$, the information about velocity drops to zero and no optimiser can recover it. A stellar process that shifts lines without changing their shape is, spectroscopically, a planet. The only escapes are other observables, other wavelengths with a different response, or time.

From photons to velocities

Calibration

A production pipeline first does the work of any spectrograph reduction: bias and dark subtraction, scattered light, tracing and optimal extraction of echelle orders, flat-fielding, cosmic-ray rejection. It then fits the map from pixel to wavelength and tracks its drift. NEID’s public algorithm documentation is unusually complete and is a good guide to reproducing these stages (NEID DRP documentation).

Three kinds of wavelength reference are used. An absorption cell (hydrogen fluoride, then iodine) shares the optical path with the star and constrains the instrumental profile, but it absorbs photons and blends its lines with the star’s. A thorium–argon lamp in a second fibre supplies laboratory emission lines and monitors drift, but its lines are unevenly spaced and the lamp ages. A laser frequency comb supplies modes at exactly known frequencies $f_n=f_{\rm CEO}+nf_{\rm rep}$, a dense and regular ruler. Known frequencies do not guarantee perfect centroids: detector response, asymmetric line-spread functions and unstable comb intensity still matter (Fischer et al. 2016). NEID’s comb modes are 20 GHz apart. NEID’s documentation notes that unstable comb flux can degrade RV precision by 5 to 10 m/s if bad calibrations propagate, which is why the pipeline checks every frame and every mode.

Iodine forward modelling

Butler and colleagues model an observed spectrum as $I_{\rm obs}(\lambda)\simeq\{I_\star(\lambda/D_\star)\,I_{\rm I_2}(\lambda)\}\otimes\mathrm{PSF}$, a Doppler-shifted stellar template times the iodine transmission, convolved with the instrumental profile. The iodine atlas comes from a Fourier-transform spectrometer. The stellar template comes from deconvolving an iodine-free observation of the same star. Each exposure is then cut into chunks of about 2 Å, and in each chunk the wavelength solution, the velocity and ten profile parameters (thirteen in all) are fitted to data oversampled fourfold. Their τ Ceti example used 704 chunks over the 5,000–6,200 Å iodine region, and a signal-to-noise ratio of 200 gave about 3 m/s (Butler et al. 1996).

The method fits the instrument’s instability instead of removing it, and its cost is a large number of small nonlinear fits. For $N$ exposures, $C_h$ chunks of $p$ samples, $k$ parameters and $I$ iterations, Gauss–Newton costs roughly $O(NC_hI(pk^2+k^3))$ (an operation count, not a measured runtime). Chunks and exposures are independent and parallelise trivially. The accuracy is limited by the template deconvolution and by degeneracies between profile and velocity, and a GPU accelerates a wrong profile model just as well as a right one.

Cross-correlation with a mask

ELODIE’s approach, inherited by HARPS and most fibre-fed instruments, multiplies the spectrum by a binary or weighted mask of line positions and slides it over trial velocities $u$,

The result is a mean line profile whose minimum gives the velocity, and whose width, contrast and bisector (the locus of midpoints across the profile at each depth) give shape diagnostics. The cost is $O(UL)$ for $U$ trial velocities and $L$ mask lines, tiny compared with extraction. A CCF works from one exposure, is easy to inspect for a second stellar spectrum, and has been extraordinarily productive. It is not a sufficient statistic: a binary mask discards gradient information, and one average profile mixes lines with different sensitivity to activity and contamination.

Template matching and line-by-line velocities

Template matching uses the star itself as its own reference. HARPS-TERRA (Anglada-Escudé & Butler 2012) and SERVAL (Zechmeister et al. 2018) coadd all spectra of a star into a high signal-to-noise template, then fit each exposure for a shift and a low-order continuum,

and iterate between template and velocities. SERVAL evaluates the template at the observed pixel wavelengths rather than resampling the data, scans velocity on a coarse grid (100 m/s by default) and refines the minimum with a parabola. It also returns a chromatic index (how the velocity depends on wavelength) and a differential line width, two diagnostics that help later to recognise stellar activity.

Étienne Artigau and colleagues pushed the idea to many small domains. Their line-by-line (LBL) method applies the derivative estimator above within each spectral segment between neighbouring template maxima, so that each “line” yields a velocity $\widehat{\Delta v}_\ell$ and an error $s_\ell$, and then combines them with a mixture model that down-weights inconsistent lines rather than letting one damaged feature drag the mean. On SPIRou near-infrared spectra of Barnard’s star the RMS fell from 7.53 m/s with a CCF to 2.59 m/s with LBL over 2.7 years (Artigau et al. 2022). The cost is $O(NP)$ per iteration for $P$ pixels, and the line-by-epoch velocity matrix it leaves behind, $O(NL)$ numbers, is itself a valuable data product.

A rare same-data comparison now exists. The NEID Earth Twin Survey extracted the same spectra of three stars three ways, with the bad-column fix applied and instrumental eras separated:

SERVAL’s formal errors are about 40% smaller than the corrected CCF’s on all three stars. These are internal uncertainties, not the scatter from stellar variability, and they rank implementations on these stars only (Gupta et al. 2026).

Line selection is a tempting further step. Pandey and colleagues fitted 2,809 vetted line–order pairs in 383 daily NEID solar spectra and greedily removed lines that most reduced the scatter. Their best single run reached 1.122 m/s with 24 lines, against 2.012 m/s for the full list. But a control that selected neighbouring lines reached 1.510 m/s, and the authors themselves raise overfitting as a possible explanation for part of the gain (Pandey et al. 2026). Choosing lines on the same time series used for detection is not the same as reducing independent noise. The selection must be validated on separate time blocks and with injected planets.

Joint stellar and telluric models

Megan Bedell, David Hogg and colleagues’ wobble treats each log-flux spectrum as the sum of a stellar component, translated by the stellar velocity, and a telluric component, fixed in the observatory frame, each with optional low-rank variability (Bedell et al. 2019). Because the two components move differently through the year, they can be separated without a model of the atmosphere. This is the same principle that lets difference imaging separate a fixed sky from a moving asteroid (Chapter 10): two signals with different motions are identifiable if the observations sample enough relative motion, and degenerate if they do not. The published TensorFlow implementation processed 91 HARPS spectra of 51 Pegasi over 72 orders in about an hour on a desktop Mac. That is not a GPU benchmark.

Two recent results argue for keeping more than one number per exposure. Langford and colleagues fitted a shared orbit to order-by-order NEID velocities and reported minimum-mass uncertainties 1.5 to 6.8 times smaller than from wavelength-averaged velocities for three known systems (Langford et al. 2025). The generative model Æstra II, discussed below, makes the same point more strongly for activity. A scalar velocity is a lossy summary, and what it loses is exactly what distinguishes a planet from the star.

Searching the velocities

Likelihood and the price of a fitted noise level

After extraction the data are $N$ velocities $v_i$ with formal errors $\sigma_i$. Write the residuals $\boldsymbol r=\boldsymbol v-\boldsymbol\mu(\theta)$ and a covariance $C_{ij}=(\sigma_i^2+s_{j(i)}^2)\delta_{ij}+k(t_i,t_j)$, where $s_j$ is an extra white “jitter” per instrument and $k$ a correlated component. Then

When the noise parameters are fitted, the determinant term is essential. Without it the likelihood rewards inflating the noise until any residual looks acceptable (Fulton et al. 2018).

How small a signal can a survey expect to find? For well-sampled circular orbits much shorter than the baseline, Andrew Cumming derived an approximate 50% detection threshold

where $F$ is the false-alarm probability allowed for the whole search and $M$ the effective number of independent frequencies (Cumming 2004). With $\sigma=1$ m/s, $N=100$, $M=10^4$ and $F=0.01$ this gives $K_0\approx0.1\times\sqrt{4\times13.8}\approx0.74$ m/s (an illustrative calculation). Two lessons follow. The threshold falls only as $\sqrt N$. And the trial factor enters only logarithmically, so a search over a wider period range costs little sensitivity. Beyond the observing baseline the picture changes completely: Cumming found that the amplitude needed for 50% completeness grows roughly as $P$, and for 99% completeness as $P^2$, because an orbit sampled over only part of its period is partly absorbed by the fitted zero points and trends.

The generalised Lomb–Scargle periodogram

The basic search fits, at each trial angular frequency $\omega$,

with weights $1/\sigma_i^2$. Mathias Zechmeister and Martin Kürster’s generalised Lomb–Scargle (GLS) periodogram gives analytic weighted sums for this fit and shows why the offset $c$ must be fitted at every frequency: subtracting the sample mean beforehand is not equivalent when the phase coverage is uneven (Zechmeister & Kürster 2009). Instrument offsets, trends or activity indicators can be added as further columns of the baseline. With diagonal noise the scan costs $O(NM)$ for $M$ frequencies. With a fixed dense covariance there is an $O(N^3)$ factorisation once, then $O(MN^2)$ whitening, a cost that structured covariances (below) reduce enormously.

Two details trip up newcomers. The natural width of a periodogram peak is about $1/T$ for baseline $T$, but the grid must be finer than that to avoid missing the top of a peak, and the number of grid points is not the number of independent trials. And irregular sampling has no single Nyquist frequency. The observing window

has strong peaks at one day, one sidereal day and one year, and each true signal produces copies (aliases) at its frequency plus and minus those of the window. A star observed once per night has a daily alias for every planet. The highest peak is a candidate, not an orbit. Competing aliases must be fitted jointly or tested with new observations. Ambiguity from nightly sampling is familiar from the linking of sparse asteroid detections (Chapter 12), where the same diurnal cadence makes several orbits fit the same points.

How significant is the highest peak?

For a single fixed frequency with known Gaussian errors, $z=(\chi_0^2-\chi_f^2)/2$ has an exponential tail. The highest of many correlated peaks does not. Roman Baluev applied the Rice method for the maximum of a random process and obtained, in the high-threshold, weak-alias regime,

where $\Delta f$ is the searched frequency band (Baluev 2008). Jean-Baptiste Delisle, Nathan Hara and Damien Ségransan extended the calculation to a known correlated-noise covariance and showed how strongly the assumed noise changes the answer (Delisle et al. 2020). Shuffling residuals to build a null distribution destroys their temporal correlation and gives invalid false-alarm probabilities. The honest procedure is to simulate the fitted null model and rerun the whole selection, including every choice made along the way.

Keplerian periodograms: a small template bank

A sinusoid is the circular-orbit limit. For an eccentric orbit at fixed $(P,e,T_p)$ the model is linear in the coefficients of $\cos\nu$, $\sin\nu$ and the offset (the term $eK\cos\omega$ is absorbed into the offset), so these can be solved analytically while the nonlinear parameters are scanned (Zechmeister & Kürster 2009):

fit baseline: instrument offsets and noise terms
for each trial period P:
    for each eccentricity e and periastron time T_p:
        solve Kepler's equation at every observation time
        solve the linear amplitude/offset coefficients
        record the likelihood improvement
    refine the best seeds
refit surviving candidates jointly; compare aliases and the null

This is a matched-filter template bank in the sense of Chapter 38, and the structure is the same as in the Solar-System searches of Part V. A shift-and-stack search (Chapter 17) scans the nonlinear parameters of a trajectory and profiles the one linear parameter, the flux, analytically. A Keplerian periodogram scans the nonlinear orbital elements and profiles the linear amplitudes analytically. The difference is scale. Shift-and-stack must evaluate billions of trajectories against gigapixels; an RV star has a few hundred numbers, and a direct bank costs $O(ME_gH_gNJ)$ for $E_g$ eccentricity seeds, $H_g$ phase seeds and $J$ Kepler-solver iterations, which a laptop handles. The hard part is not the arithmetic but the bank’s edges. High eccentricity concentrates the signal into a short periastron passage that sparse sampling can miss, and Cumming found that even noiseless tests missed about 10–20% of $e=0.9$ orbits through optimisation failure.

Real searches add planets one at a time. After each detection the previously found signals and all shared nuisance parameters must be refitted when a new one is tested, because irreversible “prewhitening” can lock onto an alias and suppress the true companion. A current example is the RVSearch pipeline in the NEID survey: frequency spacing $1/(2\pi T)$, periods from two days to twice the baseline, and a threshold at an empirically calibrated 0.1% false-alarm probability on a $\Delta$BIC periodogram (Gupta et al. 2026).

Searching for several frequencies at once

Two methods avoid committing to one frequency before considering the others. Hara, Gwenaël Boué, Jacques Laskar and Alexandre Correia solve a sparse-recovery problem,

where the columns of $A$ are sinusoids on a fine frequency grid and $B$ holds offsets and trends. The solution is a sparse spectrum that often suppresses the forest of aliases a GLS periodogram shows (Hara et al. 2017). A dense dictionary costs $O(NM)$ in storage and per iteration. On the datasets of the 2018 evidence challenge an improved implementation took 1.5 seconds on average, or 4.6 seconds including significance assessment (Nelson et al. 2020).

Baluev’s FREDEC (“frequency decomposer”) takes the combinatorial route. It starts from a deliberately generous pool of candidate frequencies and tests which subsets survive together, using analytic multi-frequency significance estimates. A pool of $q$ frequencies has $2^q-1$ subsets, so cheap preliminary tests prune the tree before expensive fits. The CUDA implementation reached about 30 times overall acceleration in double precision on an NVIDIA Tesla C2075 against single-threaded CPU code, with hundreds-fold gains on the trigonometric sums. A 150-fold single-precision result came with explicit warnings about numerical instability. Practical pools were limited to about 25–27 candidates, and the method could return several viable combinations without choosing among them (Baluev 2013). This is an early demonstration that order-of-magnitude GPU gains in RV search are real, and an equally early warning that precision and combinatorics do not go away.

Bayesian evidence and counting planets

The cleanest question is: how many planets do the data support? For a model with $n$ planets the evidence is

which penalises parameter volume the data do not need. Priors on period, eccentricity, amplitude and noise are part of the answer, not a formality. Nested sampling turns the integral into a sequence of likelihood-constrained prior volumes. The kima package uses diffusive nested sampling and treats the number of planets as a parameter, returning orbits and model probabilities together (Faria et al. 2018). Its 2026 extension fits Gaia epoch astrometry jointly with velocities, which can break the $\sin i$ degeneracy (Baycroft et al. 2026).

Evidence integrals are hard to compute reliably. In a 2018 challenge Benjamin Nelson and colleagues gave teams six simulated datasets of 200 velocities each, with common likelihoods and priors. The spread of the teams’ evidence estimates grew from a factor of about three for zero planets to about ten for one, $10^2$–$10^3$ for two, and more than $10^4$ for three, and many internal error estimates were too optimistic (Nelson et al. 2020). The difficulty is exploration, not evaluation: missed modes and badly estimated mode volumes, not slow likelihoods.

The star fights back

Three impostors

Stellar activity has humbled RV astronomy repeatedly, and three cases show the range of the problem.

In 2001 Queloz and colleagues reported that HD 166435 showed an apparently planetary velocity period of 3.7987 days with ELODIE. Photometry and Ca II H&K emission showed magnetic activity, the line bisectors correlated with velocity, and the brightness was a quarter of a cycle out of phase with the velocity, as expected for a rotating spotted star. Although an orbit fitted the data over more than two years, the signal stayed coherent for only about 30 days at a time. The title of the paper was “No planet for HD 166435” (Queloz et al. 2001).

The transiting planet CoRoT-7b showed how much an inferred mass depends on the activity model. Different teams’ filtering, prewhitening and spot models gave different masses. Raphaëlle Haywood, Andrew Collier Cameron, Queloz and colleagues obtained simultaneous CoRoT photometry and HARPS velocities in January 2012 and fitted a model in which the photometry informed a Gaussian process for the velocities. They found $4.73\pm0.95\,M_\oplus$ for CoRoT-7b, did not confirm a proposed third planet, and placed the stellar rotation at $23.81\pm0.03$ days (Haywood et al. 2014).

The sharpest warning was α Centauri B b. The claimed planet had a semi-amplitude of $0.51\pm0.04$ m/s and a period of 3.24 days in 459 HARPS observations over four years. Vinesh Rajpaul, Suzanne Aigrain and Stephen Roberts showed that the observing window plus the published detrending procedure could manufacture an apparently significant signal at that period from planet-free synthetic data sampled at the same times (Rajpaul, Aigrain & Roberts 2015). The lesson is that a pipeline, not just a star, can make a planet. A large Bayes factor between two inadequate models did not help.

Where activity comes from

A planet translates the whole spectrum. The star changes the spectrum it emits. To first order in a log-wavelength coordinate $x$, a common velocity is $S(x-v/c)\simeq S(x)-(v/c)S'(x)$, so a planet acts only along the derivative direction $S'$. Stellar processes act along other directions, but nothing makes those directions perpendicular to $S'$, and the component along $S'$ is precisely what a velocity estimator reports as motion (Rajpaul et al. 2015; Collier Cameron et al. 2021). The main processes are:

• Spots and faculae. As the star rotates, dark spots and bright faculae change the brightness weighting of the approaching and receding hemispheres, distorting the disk-integrated line. Their growth, decay and differential rotation make the signal quasi-periodic.
• Suppressed convective blueshift. On a quiet surface, hot rising gas is brighter than the cool sinking gas, so the net line is blueshifted. Magnetic regions suppress convection and reduce the blueshift. For slowly rotating Sun-like stars this can dominate the activity signal even when spots are photometrically weak (Haywood et al. 2014; Rescigno & Al Moulla 2025).
• Pressure modes. Stochastically excited acoustic oscillations with periods of a few minutes produce short-timescale velocity noise.
• Granulation and supergranulation. Convection cells of different sizes produce correlated velocity noise on timescales of minutes to days, including substantial low-frequency power (Meunier & Lagrange 2019).
• Magnetic cycles. Over years the coverage of active regions changes the mean convective suppression, adding a slow drift.

The Sun is the test case for all of this, because its centre-of-mass motion is known exactly. HARPS-N has observed it through a small solar telescope since 2015. The 2025 release contains 109,466 curated spectra over a decade, with median photon uncertainty 0.28 m/s, median daily RMS 0.49 m/s and an RMS over the decade of 2.95 m/s before removing the magnetic-cycle signature (Dumusque et al. 2025). The dominant noise changes with timescale.

Averaging what can be averaged

A finite exposure is a filter. For a rectangular exposure of length $T$, a velocity oscillation at frequency $f$ is reduced in amplitude by $\operatorname{sinc}(\pi fT)$, so the residual variance is

William Chaplin and colleagues integrated realistic oscillation spectra and found the exposure times needed to bring residual p-mode noise below 0.1 m/s: 3.9 minutes for τ Ceti, 3.6 for α Cen B and 9.7 for the Sun, with hotter or evolved stars needing over 100 minutes (Chaplin et al. 2019). These figures concern oscillations alone.

Granulation and supergranulation cannot be averaged out within a night. Niamh O’Sullivan, Aigrain and colleagues fitted the HARPS-N solar data across a magnetic cycle and found that supergranulation contributes 0.75 to 0.83 m/s, more than twice the granulation, with a characteristic timescale of about one to one and a half days that changes with activity (O’Sullivan et al. 2025):

A noise process with a timescale of a day and an amplitude of 0.8 m/s is ten times larger than the Earth’s signal and cannot be beaten by longer exposures. It can only be beaten by many independent nights, which makes scheduling part of the algorithm. Under the assumptions of the same paper, cycling through ten stars with five-minute exposures and five-minute overheads measures supergranulation timescales to 30% in 23 nights, while staring at one star for an hour at a time takes 40 nights to do the same for six stars. Spreading observations in time beats concentrating them, despite the overheads.

The ideal arithmetic is sobering. For a coherent circular signal, independent errors $\sigma$ and uniform phase coverage, the amplitude error is $\sigma_K\simeq\sigma\sqrt{2/N}$. A $5\sigma$ detection of $K=0.09$ m/s then needs about 6,200 independent measurements at $\sigma=1$ m/s, or about 560 at $\sigma=0.3$ m/s (illustrative calculations from this formula). Real surveys have unknown periods, seasonal gaps and correlated noise, and need more.

Regression on indicators, and the $FF'$ model

The cheapest mitigation regresses velocities on activity indicators measured from the same spectra: Ca II H&K emission, Hα, the CCF’s width, contrast and bisector span. It works when the relation is linear and in phase, which it often is not: velocity and indicators can be shifted in phase or trace different mixtures of surface processes. Suzanne Aigrain, Frédéric Pont and Shay Zucker derived a physically motivated alternative from photometry. If $f(t)$ is the fractional flux deficit caused by spots and $f_{\max}$ its maximum, the rotational and convective velocity perturbations are approximately

(Aigrain, Pont & Zucker 2012). This "$FF'$" model predicts velocities from a light curve with two free scaling quantities, and needs a smooth estimate of the photometric derivative.

Gaussian processes

A Gaussian process (GP) models correlated variability as a random function with a covariance kernel. The usual quasi-periodic choice is

with $\ell$ an evolution timescale, $P_{\rm rot}$ the rotation period and $w$ a harmonic-complexity parameter. Rajpaul, Aigrain, Michael Osborne, Steven Reece and Roberts turned this into a joint model. A single latent process $g(t)$ and its derivative generate the velocity and the indicators,

with planets only in the velocity mean $m_{\rm RV}$ (Rajpaul et al. 2015). The covariance between channels follows from derivatives of the kernel, so a chromospheric index that tracks $g$ constrains the part of the velocity that tracks $g$. What remains after conditioning on the indicators is a Schur complement,

An indicator uncorrelated with the troublesome process cannot remove it, however precisely it is measured.

A dense GP costs $O(n^3)$ per likelihood evaluation for $n$ total observations, which once forced teams in the 2017 RV Fitting Challenge into approximations (Dumusque et al. 2017). Daniel Foreman-Mackey, Eric Agol and colleagues’ celerite restricts the kernel to sums of damped complex exponentials, interpretable as stochastically driven oscillators, and evaluates the likelihood in $O(nr^2)$ for structured rank $r$ (Foreman-Mackey et al. 2017). Delisle, Nicolas Unger, Hara and Ségransan’s s+leaf 2 extends this to several observables with derivative couplings and different observing calendars. On 246 HARPS epochs of HD 13808 with three channels each (738 measurements) it reproduced the two-planet solution with more than a hundredfold reduction in compute, using a structured approximation to the quasi-periodic kernel (Delisle et al. 2022). Any new accelerated method should be compared against this, not against a dense Cholesky factorisation.

A GP is flexible enough to absorb a planet near the rotation period, and inflexible enough, with the wrong kernel, to leave correlated residuals that look like planets. Marginalising over its parameters and testing with injections are needed in both directions.

Separating shift from shape in the spectrum

The most promising methods return to the spectrum. Collier Cameron, Eric Ford, Sahar Shahaf and colleagues’ SCALPELS takes the autocorrelation of each CCF, $A_i(\Delta)=\int C_i(v)C_i(v+\Delta)\,dv$, which on an infinite domain is unchanged by a pure translation but responds to changes of shape. A singular-value decomposition of the autocorrelations across epochs yields temporal shape vectors $U$, and the velocities split into a shape-driven part $UU^T\boldsymbol v$ and a shift part $(I-UU^T)\boldsymbol v$. On 853 daily solar observations over five years, injected $K=0.40$ m/s signals were recovered with amplitude precision about 0.066 m/s (Collier Cameron et al. 2021). Translation invariance of the features does not make the subsequent regression planet-proof: a planet’s time series can correlate with $U$ by chance, so the planet must be fitted simultaneously.

Michaël Cretignier and colleagues’ YARARA corrects tellurics, interference fringes, detector stitching and ghosts in spectral time series, then extracts line-by-line velocities. Version 2 applies principal component analysis to those velocities after removing the common mean, on the idea that a planet moves all lines together while activity and instruments move them differently (Cretignier et al. 2021; Cretignier et al. 2023). Across 20 HARPS targets it found common instrumental components and about 0.9 m/s long-baseline RMS for several stars. Its injection test is the clearest demonstration in the literature of how cleaning can eat planets. A 122-day, 3 m/s signal fell to 1.50 m/s when the nuisance basis was fitted without a planet term, and a spurious annual peak appeared. Fitting the candidate Keplerians jointly with the basis restored injected 2 and 3 m/s signals to 1.99 and 2.89 m/s.

Joseph Salzer, Jessi Cisewski-Kehe, Ford and colleagues fitted a sparse fixed-effects model in which every NEID solar line has its own offset and its own response to measured shape changes, around a common epoch velocity. The RMS fell from 1.919 m/s without shape terms to 0.575 m/s with them, and grouped cross-validation gave a more honest 0.580–0.612 m/s (Salzer et al. 2025). Federica Rescigno and Khaled Al Moulla grouped lines by formation temperature and found that segments formed at 4,000–4,750 K gave the lowest dispersion, while hotter segments traced suppressed convective blueshift strongly; residuals still did not go below about 0.5 m/s (Rescigno & Al Moulla 2025).

Learned models

Machine learning entered through the CCF. Zoe de Beurs, Andrew Vanderburg, Christopher Shallue and colleagues trained linear, dense and convolutional models to predict activity velocity from the shape of centred HARPS-N solar CCFs, using no timestamps. On held-out solar data the test-set standard deviation fell from 1.736 to 0.967 m/s with the CNN. A frequently quoted 1.753 to 1.039 m/s is a different statistic, half the 16th–84th percentile width over the full dataset (de Beurs et al. 2022). Training took 5–20 minutes per model on a single CPU. A 2026 extension adding chromatic CCFs, photometry and magnetic flux recovered an injected 0.40 m/s, 400-day signal confidently, detected a 0.30 m/s signal weakly, and missed 0.10 and 0.20 m/s (McWilliam et al. 2026).

The strongest recent result works on full spectra. Yan Liang, Joshua Winn, Peter Melchior and colleagues’ Æstra II learns a generative model of 29,352 NEID solar spectra, each with 170,923 wavelength bins in 42 orders, using 20 latent activity coordinates. Its spectral encoder is shown to cancel a common translation in its aligned features, so it is trained once and the injected planets enter only the downstream activity-and-planet decomposition (Liang et al. 2026). The authors injected 500 single planets with periods 2.5–400 days and amplitudes 0.1–0.7 m/s:

At thresholds tuned to give zero spurious detections on this suite, Æstra II recovered 238 planets and a CCF-plus-indicator baseline 9. The baseline needed a nominal false-alarm probability below $10^{-19}$ to reach zero spurious detections, while the conventional 0.1% threshold accepted 322 wrong top candidates. That last number says more about white-noise false-alarm probabilities in the presence of activity than about any neural network. The qualifications matter equally: all 500 injections share one realisation of solar noise, a “recovery” allows 20% error in period, and the threshold was chosen on the evaluation suite itself. For 500 independent trials with no failures, the “rule of three” gives an approximate 95% upper limit on the failure rate of $3/500=0.006$, but these trials are neither independent nor held out, so even that bound does not apply.

Other 2026 models point in the same direction. CANSTAR orthogonalises CCF shift and distortion directions and feeds their time sequences to an attention network trained on simulations, reducing residual RMS to 52.5% of the original for ε Eridani (Blanco-Pozo et al. 2026). doppleriann compresses temperature-grouped spectral shells into 9×9 maps and reports reliable recovery at 25 cm/s under five-fold cross-validation, but a more conservative 0.35 m/s under chronological hold-out (Gómez-Vargas et al. 2026). The gap between those two numbers is the gap between random and time-blocked validation, and it is the difference that matters for a real star. Learned methods are discussed more generally in Chapter 22.

Tellurics

The Earth’s atmosphere imprints absorption lines that stay roughly fixed in the observatory frame while the stellar lines move with the barycentric velocity. Their depth varies with humidity and airmass, so uncorrected tellurics produce velocity errors that correlate with the annual geometry. Because $(S\,T)*\mathcal L\ne(S*\mathcal L)(T*\mathcal L)$, dividing an observed spectrum by a separately convolved telluric model is not exact. StellarSpectraObservationFitting (SSOF) models variable stellar and telluric components in linear flux with a wavelength-dependent line-spread function, and reports photon-limited errors 2 to 6 times smaller than the NEID CCF pipeline. That is an extraction gain, not a demonstrated activity-corrected detection gain (Gilbertson et al. 2026). Controlled simulations by Sharon Wang and colleagues found optical telluric residuals that can approach 1–10 cm/s under favourable assumptions, but metre-per-second residuals in the red and near-infrared (Wang et al. 2022).

Calibrating the whole procedure

Every flexible step, from line selection to GP kernels to neural networks, can both remove planets and invent them. The only defence is an injection–recovery experiment that reruns the complete pipeline:

freeze data cuts, model choices, search range and acceptance rule
reserve independent calibration and evaluation sets (time-blocked)
for each stellar state and timestamp realisation:
    generate planet-free nulls with plausible correlated noise
    inject planets with random phases and declared P, K, e
    inject into spectra or CCFs, before any Doppler-sensitive step
    rerun preprocessing, nuisance fitting, search and selection
    match candidates with declared tolerances, including aliases
tabulate completeness, spurious detections and parameter bias

Injecting into scalar velocities tests only the steps after extraction. The 2017 RV Fitting Challenge defined the figure of merit $K/N\equiv(K/\mathrm{RV_{rms}})\sqrt{N_{\rm obs}}$ and found that many methods recovered 80–90% of planets above $K/N\approx7.5$, with performance falling below that (Dumusque et al. 2017). The 2022 EXPRES Stellar Signals Project compared 22 methods on common observations of four stars and found that none consistently reached sub-metre-per-second residuals, and that methods often disagreed about the corrected velocities (Zhao et al. 2022).

Landmark results

Sources: as cited in the text, plus Faria et al. (2022) and Figueira et al. (2025). The ESPRESSO study of Figueira and colleagues is instructive. It reached below 10 cm/s on short timescales and about 40 cm/s over 3.5 years, found no conclusive ESPRESSO-only planets around three of its targets, and for τ Ceti reported sensitivity to about 1.7 Earth masses below 100 days and 2–5 Earth masses in the habitable zone. For context, radial velocity is the second most productive discovery method in the NASA Exoplanet Archive, with 1,204 of 6,451 confirmed planets as of October 2026; Chapter 35 gives the full census.

Limits and open problems

Precision is not one number. A photon-information bound, a formal extraction error, calibration stability, within-night repeatability, long-term residual scatter and injection-recovery completeness are different quantities. A method can improve one and worsen another, and papers that quote one are often read as quoting another.

Identifiability is physical. A surface process that produces a pure translation with no shape signature cannot be distinguished from orbital motion by spectroscopy at one instant. Time coherence, other wavelengths and other observables can help. No architecture can recover a distinguishing signal that is absent from the data. Deep telluric absorption and finite photon counts are equally hard limits.

Most quoted “floors” are conditional. A residual of half a metre per second is not a theorem. It reflects the calibration, the tracers available, the compression, the kernel and the phase coverage of one study. The temperature-GP study, for example, identifies unmodelled variability at its floor.

Long periods are the hardest. Instrument offsets absorb information, and an orbit sampled over only part of its period is partly fitted away by trends. Supergranulation is neither white nor stationary across a cycle. The habitable zones of Sun-like stars sit at periods of a year, where barycentric systematics and telluric errors also live.

Solar success need not transfer. The Sun offers extraordinary cadence, signal-to-noise and known dynamics. Other stars differ in inclination, spectral type, rotation and instrument history. A model trained on the Sun can learn calibration epochs or sampling patterns. Random train/test splits can share an active region across the split, while time-blocked evaluation tests the realistic case.

Significance is conditional. A 0.1% false-alarm probability rejects a specified noise-only model, not every non-planetary explanation. No result reviewed for this chapter demonstrates routine discovery of a strict Earth–Sun twin around a broad population of Sun-like stars.

Where the room is

The arithmetic of an RV planet search is small. A few hundred velocities and a bank of Keplerian templates fit comfortably on a laptop, and published GPU gains such as FREDEC’s 30× already exist. Accelerating the periodogram alone will not find Earth twins. The room is elsewhere, in four places.

Better compression before the search. Near a fixed template, the Doppler numerator $g^T\Sigma^{-1}r$ and denominator $g^T\Sigma^{-1}g$ carry the local shift information. Adding projections onto continuum, width, asymmetry and telluric modes, with their cross-information matrix, gives a small set of statistics that is sufficient for that local linear model and interpretable. This is the same move as reducing images to $\psi$ and $\phi$ arrays before shift-and-stack (Chapter 17): reduce once, keep the information, search cheaply. These reductions are batched matrix operations, partitioned by order, segment and epoch, and map naturally onto GPUs (Chapter 40). The test is whether injected planets survive realistic stellar variability, not whether the kernel reaches a high fraction of peak throughput.

Translation-equivariant spectral operators. Æstra II’s main computational idea is amortisation: a spectral representation learned once and provably insensitive to a common translation can be reused for every injected or candidate planet. That proof must cover interpolation, masks, weights, continuum fits and template construction, or the reuse quietly changes the experiment. A general library of such operators, with validated covariance propagation, would let many groups do this safely. The data volume is real: Æstra II’s arrays alone hold about $5.0\times10^9$ flux samples, 20 GB in float32 and 40 GB with weights (derived from its published dimensions), which calls for streaming orders and accumulating sufficient statistics rather than holding everything on one device.

Structure before scale. Spectral derivatives, line-depth responses, formation temperatures and observer-frame telluric modes give compact coordinates. Randomised SVDs, sparse fixed-effects regressions and batched least squares on these coordinates may beat a larger network on the same information. On the time axis, structured GPs such as s+leaf 2 already give two orders of magnitude, and batching many posterior states, stars, null realisations and injections across a GPU is the natural next step. Each recurrence is sequential in time, but the independent problems are many.

Exploration, not just evaluation, for evidence. Evidence failures in the 2018 challenge came from missed modes and wrong mode volumes. Parallel independent nested-sampling runs, importance samplers seeded from several aliases and amortised proposals learned across simulated stars can improve both throughput and reliability, provided the final estimate is corrected to the true posterior. A neural score is not a marginal likelihood.

Co-design the schedule. The information about a planet’s amplitude, after marginalising nuisance parameters, is a Schur complement that depends on when the star is observed. Optimising it directly, and adaptively revisiting to separate oscillations, granulation, rotation and aliases, can be worth more than additional exposures.

Two cautions run through all of these. First, numerical precision: subtracting nearly equal spectra, accumulating weighted sums and solving poorly conditioned nuisance systems expose rounding, and FREDEC’s single-precision instability is direct evidence. Centre times and velocities, scale design matrices, and accept reduced precision only after comparing velocity bias and recovery against float64. Second, the benchmark: a claimed tenfold improvement should mean tenfold lower cost at the same recovery curve and the same false-positive rate, with the entire adaptive pipeline rerun under the null, on public spectra with frozen calibration eras and masks. The research for this chapter found no such end-to-end benchmark on modern accelerators.

The hard limits are finite photons, missing photons, incomplete phase coverage and nuisance signals that are indistinguishable within the measurements. The engineering limits are calibration fidelity, detector modelling, data movement, redundant trial work and inadequate posterior exploration. Algorithms can move the second set a long way. They cannot move the first. Chapter 42 places these in the context of the whole book.

References

1. Struve, O. (1952). Proposal for a project of high-precision stellar radial velocity work. The Observatory 72, 199. ADS 1952Obs....72..199S
2. Campbell, B., Walker, G. A. H. & Yang, S. (1988). A search for substellar companions to solar-type stars. ApJ 331, 902. doi:10.1086/166608
3. Latham, D. W., Mazeh, T., Stefanik, R. P., Mayor, M. & Burki, G. (1989). The unseen companion of HD 114762: a probable brown dwarf. Nature 339, 38. doi:10.1038/339038a0
4. Mayor, M. & Queloz, D. (1995). A Jupiter-mass companion to a solar-type star. Nature 378, 355. doi:10.1038/378355a0
5. Baranne, A., Queloz, D., Mayor, M. et al. (1996). ELODIE: a spectrograph for accurate radial velocity measurements. A&AS 119, 373. ADS 1996A&AS..119..373B
6. Butler, R. P., Marcy, G. W., Williams, E., McCarthy, C., Dosanjh, P. & Vogt, S. S. (1996). Attaining Doppler precision of 3 m/s. PASP 108, 500. ADS 1996PASP..108..500B
7. Queloz, D., Henry, G. W., Sivan, J. P. et al. (2001). No planet for HD 166435. A&A. arXiv:astro-ph/0109491
8. Mayor, M., Pepe, F., Queloz, D. et al. (2003). Setting new standards with HARPS. The Messenger 114, 20. ESO Messenger
9. Cumming, A. (2004). Detectability of extrasolar planets in radial velocity surveys. MNRAS. arXiv:astro-ph/0408470
10. Baluev, R. V. (2008). Assessing the statistical significance of periodogram peaks. MNRAS. arXiv:0711.0330
11. Zechmeister, M. & Kürster, M. (2009). The generalised Lomb–Scargle periodogram. A&A. arXiv:0901.2573
12. Aigrain, S., Pont, F. & Zucker, S. (2012). A simple method to estimate radial velocity variations due to stellar activity using photometry. MNRAS. arXiv:1110.1034
13. Anglada-Escudé, G. & Butler, R. P. (2012). The HARPS-TERRA project I. ApJS. arXiv:1202.2570
14. Baluev, R. V. (2013). Detecting multiple periodicities in observational data with the multifrequency periodogram. II. Frequency Decomposer. Astronomy and Computing. arXiv:1309.0100
15. Haywood, R. D., Collier Cameron, A., Queloz, D. et al. (2014). Planets and stellar activity: hide and seek in the CoRoT-7 system. MNRAS. arXiv:1407.1044
16. Wright, J. T. & Eastman, J. D. (2014). Barycentric corrections at 1 cm/s for precise Doppler velocities. PASP. arXiv:1409.4774
17. Rajpaul, V., Aigrain, S., Osborne, M. A., Reece, S. & Roberts, S. (2015). A Gaussian process framework for modelling stellar activity signals in radial velocity data. MNRAS. arXiv:1506.07304
18. Rajpaul, V., Aigrain, S. & Roberts, S. (2015). Ghost in the time series: no planet for Alpha Cen B. MNRAS Letters. arXiv:1510.05598
19. Fischer, D. A. et al. (2016). State of the field: extreme precision radial velocities. PASP. arXiv:1602.07939
20. Dumusque, X., Borsa, F., Damasso, M. et al. (2017). Radial-velocity fitting challenge II. A&A. arXiv:1609.03674
21. Hara, N. C., Boué, G., Laskar, J. & Correia, A. C. M. (2017). Radial velocity data analysis with compressed sensing techniques. MNRAS. arXiv:1609.01519
22. Foreman-Mackey, D., Agol, E., Ambikasaran, S. & Angus, R. (2017). Fast and scalable Gaussian process modeling with applications to astronomical time series. AJ. arXiv:1703.09710
23. Zechmeister, M. et al. (2018). Spectrum radial velocity analyser (SERVAL). A&A. arXiv:1710.10114
24. Fulton, B. J., Petigura, E. A., Blunt, S. & Sinukoff, E. (2018). RadVel: the radial velocity modeling toolkit. PASP. arXiv:1801.01947
25. Faria, J. P., Santos, N. C., Figueira, P. & Brewer, B. J. (2018). kima: exoplanet detection in radial velocities. JOSS 3, 487. arXiv:1806.08305
26. Bedell, M., Hogg, D. W., Foreman-Mackey, D., Montet, B. T. & Luger, R. (2019). wobble: a data-driven analysis technique for time-series stellar spectra. AJ. arXiv:1901.00503
27. Chaplin, W. J., Cegla, H. M., Watson, C. A. et al. (2019). Filtering solar-like oscillations for exoplanet detection in radial velocity observations. AJ. arXiv:1903.00657
28. Meunier, N. & Lagrange, A.-M. (2019). Unexpectedly strong effect of supergranulation on the detectability of Earth twins orbiting Sun-like stars with radial velocities. A&A. arXiv:1904.09089
29. Delisle, J.-B., Hara, N. & Ségransan, D. (2020). Efficient modeling of correlated noise I. Statistical significance of periodogram peaks. A&A. arXiv:2001.10319
30. Nelson, B. E. et al. (2020). Quantifying the Bayesian evidence for a planet in radial velocity data. AJ. arXiv:1806.04683
31. Pepe, F. et al. (2021). ESPRESSO@VLT: on-sky performance and first results. A&A. arXiv:2010.00316
32. Collier Cameron, A., Ford, E. B., Shahaf, S. et al. (2021). Separating planetary reflex Doppler shifts from stellar variability in the wavelength domain. MNRAS. arXiv:2011.00018
33. Cretignier, M., Dumusque, X., Hara, N. C. et al. (2021). YARARA: significant improvement in RV precision through post-processing of spectral time series. A&A. arXiv:2106.07301
34. Delisle, J.-B., Unger, N., Hara, N. C. & Ségransan, D. (2022). Efficient modeling of correlated noise III. Scalable methods for jointly modeling several observables' time series with Gaussian processes. A&A. arXiv:2201.02440
35. de Beurs, Z. L., Vanderburg, A., Shallue, C. J. et al. (2022). Identifying exoplanets with deep learning IV. Removing stellar activity signals from radial velocity measurements using neural networks. AJ. arXiv:2011.00003
36. Faria, J. P. et al. (2022). A candidate short-period sub-Earth orbiting Proxima Centauri. A&A. arXiv:2202.05188
37. Artigau, É. et al. (2022). Line-by-line velocity measurements: an outlier-resistant method for precision velocimetry. AJ 164, 84. arXiv:2207.13524
38. Wang, S. X., Latouf, N., Plavchan, P. et al. (2022). Characterizing and mitigating the impact of telluric absorption in precise radial velocities. AJ. arXiv:2206.07287
39. Zhao, L. L., Fischer, D. A., Ford, E. B. et al. (2022). The EXPRES Stellar Signals Project II. State of the field in disentangling photospheric velocities. AJ. arXiv:2201.10639
40. Cretignier, M., Dumusque, X., Aigrain, S. et al. (2023). YARARA V2: reaching sub-m/s precision over a decade using PCA on line-by-line radial velocities. A&A. arXiv:2308.11812
41. Rubenzahl, R. A. et al. (2024). KPF confirms a polar orbit for KELT-18 b. AJ. arXiv:2407.21196
42. Rescigno, F. & Al Moulla, K. (2025). Gaussian process regression of temperature-dependent radial velocities. MNRAS. arXiv:2501.02959
43. Salzer, J., Cisewski-Kehe, J., Ford, E. B. et al. (2025). Searching for low-mass exoplanets amid stellar variability with a fixed effects linear model of line-by-line shape changes. AJ. arXiv:2502.11930
44. O'Sullivan, N. K., Aigrain, S., Cretignier, M. et al. (2025). Measuring the Sun's radial velocity variability due to supergranulation over a magnetic cycle. MNRAS. arXiv:2506.23693
45. Figueira, P. et al. (2025). A comprehensive study on radial velocity signals using ESPRESSO: pushing precision to the 10 cm/s level. A&A. arXiv:2507.07514
46. Langford, Z. et al. (2025). Order-by-order modeling of exoplanet radial velocity data. arXiv:2510.16139
47. Dumusque, X., Al Moulla, K., Cretignier, M. et al. (2025). A decade of solar high-fidelity spectroscopy and precise radial velocities from HARPS-N. A&A. arXiv:2510.27635
48. Gilbertson, C., Ford, E. B., Halverson, S. et al. (2026). Data-driven modeling of telluric features and stellar variability with StellarSpectraObservationFitting.jl. AJ. arXiv:2408.17289
49. McWilliam, N., de Beurs, Z. L., Vanderburg, A. et al. (2026). Identifying exoplanets with deep learning VI. Enhancing neural network mitigation of stellar activity RV signals with additional metrics. AJ. arXiv:2602.17760
50. Blanco-Pozo, J., Perger, M., Anglada-Escudé, G. et al. (2026). Mitigating stellar radial velocity jitter using orthogonal activity indices and a time-aware neural network. A&A. arXiv:2605.04587
51. Liang, Y., Winn, J. N., Melchior, P., Lu, S. & Tran, Q. H. (2026). AESTRA II: generative spectral modeling of the Sun as a star for precise radial velocities. arXiv:2606.13574
52. Gómez-Vargas, I., Dumusque, X., Zhao, Y. et al. (2026). Modeling Doppler shifts in radial-velocity data with deep learning toward Earth-mass exoplanet detection. A&A. arXiv:2606.18464
53. Baycroft, T. A., Faria, J. P. & Delisle, J.-B. (2026). Bayesian analysis of Gaia epoch astrometry and radial velocities with kima. A&A. arXiv:2606.24132
54. Pandey, K. et al. (2026). Improving the precision of line-by-line radial velocities: a data-driven iterative algorithm for spectral line selection. arXiv:2606.29552
55. Gupta, A. F. et al. (2026). The NEID Earth Twin Survey VI. Discovery of a low-mass planet orbiting HD 126053 and RV signals near the rotation periods of HD 168009 and HD 10780. arXiv:2609.32097
56. NEID team (2026). NEID data reduction pipeline: algorithms. Living documentation, accessed 2026-10-16. neid.ipac.caltech.edu
57. Keck Observatory (2026). KPF instrument status. Accessed 2026-10-16. www2.keck.hawaii.edu