A fast test for the identification and confirmation of massive black hole binaries

We present a new observational test to identify massive black hole binaries in large multi-epoch spectroscopical catalogues and to probe the real nature of already proposed binary candidates. The test is tailored for binaries with separations large enough to allow each black hole to retain its own broad line region. In this limit the fast AGN variability typically observed over months cannot be associated to the much longer binary period, and it is assumed (as for the case of single black holes) to be the consequence of the evolution of the innermost regions of the two accretion discs. A simple analysis of the cross-correlation between different parts of individual broad emission lines can therefore identify the presence of two massive black holes whose continua vary independently of each other. Our analysis indicates that to be less affected by the noise in the spectra the broad lines should be divided in two close-to-equal-flux parts. This ensures that in the single massive black hole scenario the cross-correlation will always be high. With monitoring campaigns similar to those performed for reverberation mapping studies, on the other way, a binary can show any value of the cross-correlation and can therefore be distinguished from a standard AGN. The new test can be performed over timescales orders of magnitude shorter than the alternative tests already discussed in literature, and can be a powerful complement to the massive black hole binary search strategies already in place.


Introduction
Gravitationally bound massive black hole binaries (MBHBs) are the outcome of galaxy mergers, and were originally predicted by Begelman et al. (1980).They are among the loudest sources of gravitational waves (GWs) detectable with current pulsar timing array (PTA) campaigns (Verbiest et al. 2016) and by future space-based gravitational wave interferometers (e.g.LISA Amaro-Seoane et al. 2017, 2023).The GW background generated by MBHBs, and the rate of detectable binary coalescence, are still quite uncertain.This is partly a consequence of the uncertainties on the efficiency of MBH pairing during the early stages of galaxy mergers (Fiacconi et al. 2013;del Valle et al. 2015;Tamburello et al. 2017;Souza Lima et al. 2017;Bortolas et al. 2020Bortolas et al. , 2022)).The electromagnetic (EM) identification of a sample of already bound MBH binaries would therefore greatly reduce the uncertainties on the signals observable with both PTA and LISA.
To date, unequivocal observational confirmation of an MBHB has not yet been found.The only spatially resolved MBHB candidate is 0402+379 (Rodriguez et al. 2009;Burke-Spolaor 2011), which has two flat-spectrum radio cores with a projected separation of ≈ 7 pc.
⋆ massimo.dotti@unimib.itAt smaller, spatially unresolved scales, MBHBs have been consistently searched for either through peculiar spectral features (Tsalmantza et al. 2011;Eracleous et al. 2012;Ju et al. 2013;Shen et al. 2013;Wang et al. 2017) or through photometric variability (Valtonen et al. 2008;Ackermann et al. 2015;Graham et al. 2015;Li et al. 2016;Charisi et al. 2016;Sandrinelli et al. 2016Sandrinelli et al. , 2018;;Severgnini et al. 2018;Li et al. 2019;Liu et al. 2019;Chen et al. 2020).These two search strategies target MBHBs at different separations.The spectroscopic approach assumes that the broad line region (BLR) of at least one of the binary components does not extend beyond its MBH Roche radius.Therefore, by sharing the same kinematics as the MBH, the broad emission lines (BELs) are shifted in frequency with respect to the host galaxy rest frame and evolve in time over a binary orbital period τ orb .The constraint on the size of the BLR can be translated into a constraint on the minimum binary separation, which is of the order of ∼ 0.1 pc for a 10 8 M ⊙ MBH accreting at one-tenth of its Eddington ratio (see section 2 for more details).
At smaller separations, the BLR is either truncated by the time-dependent binary potential (Montuori et al. 2011) or shared by both MBHs and comoving with the MBHB centre of mass.Theoretical studies predict periodic variability on timescales of ∼ τ orb , either associated with the periodic fueling from the circumbinary material (e.g.Hayasaki et al. 2008) or with the Doppler boosting of the emission for very close binaries (D'Orazio et al. 2015).At close separations, such timescales can be as short as ∼ < 1 yr, allowing observational searches in current and future multi-epoch observations (see Graham et al. 2015;Charisi et al. 2016;Liu et al. 2019;Chen et al. 2020, and references therein).
Nevertheless, a major problem is that all of the abovementioned features used to identify MBHB candidates can have alternative explanations (see Dotti et al. 2023, for an overview of the possible alternative interpretations).Theoretical prediction of clear and unique observational features associated with MBHBs is therefore needed in order to test the actual binary nature of both small-(selected through their variability) and largescale (spectroscopically identified) MBHB candidates.Some tests have recently been proposed for the small-scale MBHB candidates, including the possibility of periodic gravitational lensing from the companion of the active component of the MBHB (for binaries observed nearly edge-on, see e.g.D 'Orazio & Di Stefano 2018;Davelaar & Haiman 2022b,a) and periodic evolution of the polarisation fraction and angle (Dotti et al. 2022).
For MBHBs at larger separations (e.g.Gaskell 1996;Eracleous et al. 1997), a straightforward test consists in observing the expected Doppler drift in the BEL profiles according to the motion of the active component (or components, if both MBHs are active) of the binary with a τ orb period (e.g.Gaskell 1996;Eracleous et al. 1997).A conclusive test would require following the system spectroscopically for ∼ τ orb .However, at these scales, τ orb can be as high as ∼ > 10 − 100 yr (see section 2), making such a test (dubbed slow periodicity test (SPT) hereafter) challenging and unpractical for some of the candidates (see e.g.Eracleous et al. 2012;Decarli et al. 2013;Runnoe et al. 2017).
In this paper, we propose a new test for the same largeseparation binary candidates based on the assumption that, at such large separations, the short-timescale ( ∼ < 1 day) intrinsic accretion variability is unrelated to the binary period, and that, when both MBHs are active, their variability patterns are uncorrelated.Each BLR reverberates to the varying continuum of its MBH (e.g.Blandford & McKee 1982).As we argue in the following, dividing the BELs into two components (one 'red' and one 'blue', which together account for the total flux in the BEL) and cross-correlating the time-evolution of their fluxes can provide a means to identify MBHBs.Indeed, in the binary scenario, each observed BEL would be composed into two different broad lines, each with a different velocity offset with respect to the galaxy rest frame, and with the two components varying independently.In this scenario, the cross-correlation between the blue and red parts of the BELs can be significantly smaller with respect to the case of an AGN powered by a single MBH, as long as the velocity offset is not negligible (i.e. each component associated to one of the MBHs contributes significantly to one of the blue or red parts of the global BEL).This implies that our test becomes decreasingly effective for larger and larger binary separations, with a cross-correlation between the two different parts of the BEL getting closer and closer to the high values expected for single MBHSs, and failing to detect any true binary for close-to-zero relative line-of-sight velocity between the two MBHs (see section 4).However, while other proposed spectroscopical tests are expected to fail when the velocity shift between the two components is equal to their width (e.g.Shen & Loeb 2010), in section 4 we demonstrate that our test can still identify binaries at larger separations and smaller velocity shifts, including binaries with periods of up to ∼ 1000 yr.
Our proposed method can be thought of as a simplified version of the tests proposed by Wang et al. (2018); Songsheng et al. (2020), in which a catalogue of two-dimensional transfer functions (TFs) for different binaries is constructed and compared to the full TF of a candidate binary.While our method is clearly less sensitive to the details of the BEL profile variability, it provides fast and quantitative verification of the binary nature of the candidate, even when the quality of the data is suboptimal for constraining the details of the TF.
The clearest advantage of our new test (dubbed fast uncorrelated variability test; FUVT) is that it can be performed over timescales comparable to the typical duration of the reverberation mapping (RM) campaigns (e.g.Bentz et al. 2009b), which is usually smaller than 1 yr (down to weeks, and orders of magnitude smaller than τ orb ).A second advantage is that, while FUVT can be applied to already identified MBHB candidates with shifted asymmetric or double-peaked BELs, it can also be used to identify new binaries with apparently 'standard' BELs in large multi-epoch samples of AGNs.
This paper is organised as follows: in § 2 we describe the timescales required by FUVT and its range of applicability in terms of intrinsic MBHB properties; in § 3 we describe how the FUVT is structured, and we present the results of its application to a real sample of type I AGN (sample I hereafter) studied using RM; in § 4 we explain how we constructed a mock catalogue of MBHBs (sample II) starting from sample I, and present the results of the application of FUVT to the new binary sample; we then conclude with § 5, summarising the main results of our study, describing the advantages and disadvantages of FUVT with respect to SPT, and discussing how future observations can further strengthen the test results.

Analytical estimates of FUVT range of applicability
As mentioned in Sect. 1, both SPT and FUVT assume that both the BLRs are bound to and comoving with the individual components of the binary.The BLR radius depends on the BEL we are focusing on and on the luminosities of the accreting MBHs.
For the broad Hβ line used in the following, the BLR radius is (Bentz et al. 2009a) where λL λ,5100 is the monocromatic luminosity of the AGN continuum at 5100 Å.Assuming that the bolometric luminosity is L bol ≈ 9 λL λ,5100 (Kaspi et al. 2000), we can express eq. 1 as a function of the individual MBH mass M and its Eddington ratio It is required that this radius be smaller than the Roche lobe radius of each individual MBH (Montuori et al. 2011).For the test to work, the fluxes from the two accretion discs (and therefore the BLR radii) have to be comparable.For this reason, the minimum separation between the two MBHs is set by the Roche lobe of the secondary MBH, which for circular binaries is (Eggleton 1983) where a is the separation between the two MBHs and q = M 2 /M 1 is the ratio between the secondary and the primary masses.
The minimum separation of the binary and, assuming circular Keplerian orbits, the minimum value of τ orb for which both SPT and FUVT are applicable is obtained by equating eq. 2 and eq.3: where f Edd,2 is the Eddington ratio of the secondary.It is interesting to note that τ orb,min has a stronger dependence on f Edd,2 than on M 2 .Therefore, at fixed luminosity, shorter minimum periods are expected for higher MBH masses.Similarly to the spectroscopic MBHB candidates proposed to date (e.g.Tsalmantza et al. 2011;Eracleous et al. 2012), secondary MBHs of 10 8 M ⊙ could be in binaries with minimum periods of as small as ∼ 27 yr (see also Decarli et al. 2013;Runnoe et al. 2017) if accreting at f Edd,2 = 0.01, that is, with the same bolometric luminosity as an MBH of 10 6 M ⊙ accreting at its Eddington limit.However, we note that for disc-like BLRs, the region within which circular orbits around a single MBH are stable is sizeably smaller (by a factor of ≈ 4 − 5) than the Roche lobe (e.g.Eggleton 1983;Runnoe et al. 2015).On top of this consideration, gas at distances of > R B−Hβ contribute to the central part of the broad Hβ line, and, if the binary separation is too small, such a contribution would be lost and the line would acquire a 'boxy' profile1 .In the following, we consider binaries with periods corresponding to BLR sizes similar to the secondary Roche lobe, as well as binaries with significantly larger periods, for which a stable disc-like BLR can survive around each of the two MBHs (see section 4).
Finally, the typical timescale for the Hβ reverberation to the change of the continuum can be estimated as τ B−Hβ = R B−Hβ /c (where c is the speed of light).This timescale was evaluated using RM and ranges from days to months depending on the luminosity of the accreting MBH.However, RM studies typically require longer durations -even for ∼ 10 6 M ⊙ MBHs-due to the need for sufficiently long baselines and sufficiently varying continua to allow for the test.As a consequence, the typical timescales needed to probe the variability both in the continuum and in the BELs are of the order of weeks to months (Bentz et al. 2009b).

FUVT description and results on a control sample of single MBHs
We apply FUVT to nine AGN for which τ B−Hβ (and therefore masses) were obtained by Bentz et al. (2009b): Mrk 142, SBS 1116+583A, Arp 151, Mrk 1310, Mrk202, NGC 4253, NGC 4748, NGC 5548, and NGC 6814.We assume, as a working hypothesis, that the whole sample is composed of only bona fide single MBHs, and use it as a comparison for the mock binary sample described in the following section.Nevertheless, we stress that at least one of the AGN in sample I (NGC 5548) has been proposed as a binary candidate (Li et al. 2016) based on the variability of the periodic modulation of its broad Hβ line.For all AGN, Bentz et al. (2009b) estimated delays between the continuum and the broad Hβ of between 2 and 7 days and masses in the range of 1 − 7 × 10 6 M ⊙ , with the exception of NGC 6814, which has a mass ≈ 1.85 × 10 7 M ⊙ , and NGC 5548, which has a mass of 8.2 × 10 7 M ⊙ (see Bentz et al. 2009b, for additional details on the properties of the single AGNs).
For our analysis, we used scaled spectra released by these latter authors after the application of a renormalisation used to set the flux of all spectra to a consistent scale.For every AGN, we removed the mean spectrum (averaged over all the observations) and worked on the variable part of the spectrum only, which is dominated by the AGN continua and the BELs.We then removed the AGN continuum by fitting a straight line to the spectrum near the Hβ BEL using the same wavelength intervals as those used in Bentz et al. (2009b). 2or each AGN, we then computed the root mean square (RMS) spectrum.We used it to identify seven λ dividing the RMS broad Hβ into eight parts of equal flux3 .For each λ in each spectrum, we divided the Hβ BEL into two components, a red one and a blue one (r-Hβ and b-Hβ, for wavelengths longer and shorter than λ, respectively).We then computed the light curves of the r-Hβ and b-Hβ components, and performed the crosscorrelation between the two, allowing for a small time-shift τ shift in the [−10 day, 10 day] interval to take into account different reverberation times for different parts of the broad lines (due to e.g. the inflowing or outflowing BLR dynamics; e.g.Bentz et al. 2009b).More specifically, we follow the standard RM practice of first using one light curve and interpolating the other, before swapping the two and repeating the exercise.The average crosscorrelation is then used, and the uncertainties on the the crosscorrelation and τ shift are estimated using the public Python version (Sun et al. 2018) 4 of the Monte Carlo code discussed in Peterson et al. (1998Peterson et al. ( , 2004)).
For each AGN and each λ, we measure max − CCF( λ), that is, the maximum value of the r-Hβ-b-Hβ cross-correlation, and its uncertainties.The max − CCF( λ) obtained following this procedure are shown in figure 1 for NGC 4748.The peak of max − CCF( λ) (CCF hereafter)5 corresponds to the λ threshold dividing the broad Hβ line into two parts of equal flux, while the cross-correlation decreases significantly (and its uncertainties increase) when moving λ towards the line wings, which are more affected by the noise in the spectra.More generally, in the sample, the values of CCF are comparable to or slightly higher than the peak of the cross-correlations between the whole broad Hβ and the photometric B and V light curves shown in Bentz et al. (2009b).The trends discussed for NGC 4748 are common for all nine of the AGN we examined: CCF is always found in the bulk of the line (hereby defined as the region in between the third and the fifth values of λ, i.e. when the least luminous part of the BEL has at least three-eighths of the total BEL flux), and it is always ∼ > 0.75. 6.

Application of FUVT to mock MBHBs
We apply FUVT to mock equal-mass MBHBs generated by duplicating the temporal series of spectra of each of the nine single MBHs discussed above.One of the copies is shifted in time by half of the observational period (assuming periodic boundary conditions, as done in standard reverberation mapping studies).As the magnitude of the time-shift τ shift in the test is constrained to be smaller than 10 days, such a shift makes the evolution of the two series of spectra independent.The original and shifted light curves (blue and red lines) of the whole Hβ for NGC 4748, together with the sum of the two (green line), mimicking the light curve of a MBHB, are shown in figure 2. 7The two series of spectra are then shifted in frequency by the same amount towards higher and lower wavelengths, respectively.The wavelength shifts are chosen in order to mimic circular MBHBs with periods of τ orb = 50, 100, 300, and 1000 yr, observed close to edge-on at the orbital phase that would maximise the Doppler effect.We chose such a configuration to check whether or not the test succeeds in identifying binaries in the most favourable scenario.Different configurations would result in smaller shifts between the two independent components of the BEL, and could in principle result in a missed detection.We note that, while all four periods fulfill the criterion set in equation 4, τ orb = 50, 100 would correspond to separations at which a sizable part of a disc-like BLR would be unstable.Furthermore, for those periods, the outer parts of the BLR would be removed, resulting in a more 'boxy' profile of each BEL component (see e.g.Nguyen et al. 2019), while only a few of the AGN discussed in section 3 have a boxy broad Hβ.For these reasons, hereafter we refer to binaries with τ orb = 300 and 1000 yr as 'solid' binaries, while keeping the τ orb = 50 and 100 yr binaries to highlight how the performance of our test scales with the orbital period and to account for possible alternative geometries and dynamics of the BLR.
We then add the two series up (adding the original errors in quadrature for each wavelength bin), obtaining four new mock MBHBs for each original MBH, for a grand total of 36 mock MBHBs.The upper panel of figure 3 shows an example of a mock MBHB spectrum based on the observations of NGC 4748 -without the application of any time shift, for clarity purposes-and for a binary with τ orb = 50 yr.In the lower panel the resulting line profile starting from the same spectrum are shown for four different mock binary periods.
We then apply FUVT to these new sets of mock spectra.The resulting profile of max − CCF( λ) as a function of λ for the four mock MBHBs constructed starting from NGC 4748 are shown in figure 4.
The four mock binaries show clearly different max−CCF( λ) profiles with respect to the single MBH case; they tend to have significantly smaller values of max − CCF( λ) than the single MBH cases close to the bulk of the broad line: the shorter the period, the larger the frequency separation between the binary component contributions, and the lower the max − CCF( λ).The highest values of the mock MBHB cross-correlation are instead found close to the line wings, where a correlation might be found due to the low signal-to-noise ratio of one of the two sides of the line.The anti-correlation observable for central values of λ is due to the shape of the broad Hβ light curve, which features a single peak appearing after a dimmer state (see figure 2 and the following discussion in this section).Such an anti-correlation should not be considered a solid feature associated with binarity; indeed, it is not observed in some of the other systems.We note that even the mock binary with a 1000 yr period shows a clear minimum in the max − CCF( λ) profile, with the two parts of the line being almost completely uncorrelated, regardless of the relatively small velocity shift between the two components as- sociated with the two MBHs.In this case, the relative velocity between the two MBHs is ≈ 500 km/s (equivalently, the shift of one of the peaks with respect to the centroid is ≈250 km/s), and the FWHM of the Hβ line is ≈2000 km/s (σ ≈ 1000 km/s) in the mean spectrum for NGC4748.If we consider the varying component only, looking at the RMS spectrum for the same system, we get FWHM≈1200 km/s (σ ≈650 km/s)8 .While the velocity shift is smaller than the FWHM, our test still manages to identify the mock as a binary.
More generally, four other sets of mock binaries based on other observed AGN show trends that are qualitatively similar to the one shown in figure 4; three sets have values of CCF that can be higher or lower than that of the single MBH case, depending on the choice of orbital period; and one set has lower CCF compared to the corresponding single MBH, which nevertheless remains relatively high, up to ≈ 0.75.An example mock If instead the Hβ light curve, for example, shows two peaks and two troughs, both equally spaced and of similar duration, the time shift will not significantly reduce the maximum crosscorrelation.We stress that such a variety of behaviours is due to the limited timescale of the spectroscopic monitoring campaigns.In principle, in the binary scenario, a campaign of arbitrary timescale would have a collection of BEL profiles (one per pointing) that are contributed from the two completely uncorrelated components (1 and 2) associated to the two MBHs.In this case, when cross-correlating the blue and red parts of the total BEL, the cross terms in the cross-correlation would be equal to zero 10 .The cross correlation between the two sides of the BEL would then read: where σ(x) is the standard deviation of the x quantity.In the simplifying case, in which the time delay maximising the crosscorrelation is the same for the two components, as is the case for our mock binaries, we obtain From eqs. 5 and 6, we can derive some interesting trends and limits: when the shifts of the two components (due to the velocity of the two MBHs along the line of side) tend to zero, and the line is split approximately in half, the denominator of the σ ratios tends to 2 σ(r 1 ) σ(b 1 ) ≈ 2 σ(r 2 ) σ(b 2 ), and the crosscorrelation in eqs. 5 and 6 tends to the arithmetic average of the cross-correlation between the blue and red parts of the two single components.If these cross-correlations are strong (as expected for a single MBH-BLR system; see section 3), the total crosscorrelation will be equally strong, as expected when there is zero velocity shift between the two components.On the other hand, the maximum cross-correlation in eq.6 will become lower and lower (while remaining positive) as the two components 1 and 2 become more shifted from the reference frame of the galaxy and contribute more asymmetrically to the red and blue parts of the whole line.Unfortunately, in many cases (such as those considered in this paper), the length of the RM campaigns will not be sufficient to ensure that the two components are completely uncorrelated, leading to the variety of outcomes mentioned above.
The dependence of the profile of max − CCF( λ) for binaries on the specific realisation of the variability observed during the RM campaign has important implications for the identification of the studied systems as MBHBs.In following section, we further discuss this point.
In principle, with a sufficiently large population of single MBHs studied with RM in order to constrain the distribution of CCF in the line bulk (i.e.measured for one of the three central λ), we could translate the values of CCF of our mock binaries into a probability of not belonging to the 'standard' population.To date, this distribution is largely unconstrained.As an example, in figure 6 we show the distributions of CCF for the observed single MBHs discussed in section 3 (blue histogram) and the mock binaries (red histograms).The two distributions are clearly different.Starting from the observed CCF for the nine AGNs, we construct a model of the CCF distribution for 'typical' single MBH, which we then use to infer the probability that mocks belong to the single MBH population.Our model is based on the assumption that (i) the population of single MBHs has a well-defined value of CCF, and that (ii) the red and blue parts of the Hβ flux data follow a bivariate Gaussian distribution, with the aforementioned, unknown correlation coefficient CCF.The latter allows us to exploit the analytical form of the sampling distribution of the CCF derived by Fisher (1928) given a set of observed CCF i .We infer a posterior distribution for CCF using a nested sampling algorithm, which we use as the model for the expected distribution of measured CCF for single MBHs.The cumulative probability of observing a given value of CCF (again restricting our search to the bulk of the BELs) is shown as a black curve in figure 6.The vertical blue and red ticks highlight the values of CCF measured for all the observed single and mock binaries in the bulk of the BEL.Mocks with lower values of CCF have lower probabilities, with 24 (out of 36) mocks having a probability of < 10 −3 .The values of CCF and of the associated probabilities are reported in Table 1. 11  11 We stress that the somewhat low probability values for two of the single MBHs (Mrk142 and NGC5548) should be considered with caution,

Discussion
We present a novel method (FUVT) to search for MBHBs and to test the binary hypothesis in already identified MBHB candidates.The method is tailored to seek large separation binaries, when the two MBHs can retain their own BLR, with orbital periods that can be ∼ > 50 yr.FUVT assumes that the two MBHs are at sufficiently large distances to ensure that the short-term variability of the two BLRs is uncorrelated.A simple correlation test between the long-and short-wavelength parts of BELs can therefore identify binaries if the maximum cross-correlation (once the BELs have been split into two comparable flux components) is sufficiently small compared to the high level of correlation expected and observed for single MBHs.
As opposed to the already proposed SPT, which requires follow-up campaigns with a minimum duration of ∼ τ orb in order to track the expected evolution of the orbital velocity of one active component of the binary, FUVT works on significantly shorter observational campaigns.The expected timescales of the test are from weeks to a few years in the worst-case scenario, provided that a sufficiently frequent time coverage (typical of RM studies) of the candidate spectra is achieved.
In addition to the obvious advantage of its short duration, FUVT has a number of specific advantages and disadvantages with respect to SPT, making the two procedures complementary: -FUVT performs best when the two MBHs in the binary have similar BEL luminosities, on average.According to theoretical studies (e.g.Roedig et al. 2012;Duffell et al. 2020, and references therein), secondary MBHs are expected to experience higher accretion rates, making FUVT particularly relevant either for systems of close to equal mass (and equal Eddington ratio) or for systems in which the secondary luminosity is close to its Eddington limit.If the flux coming from one of the two BLRs is too small compared to the other, FUVT could fail in distinguishing the subdominant contribution from the observational noise.
-If however the two MBHs contribute similarly to the BELs, the previously proposed tests (and SPT in particular) could miss real binaries, because at large separations, doublepeaked profiles are not always detectable (see e.g.fig 3) or could be misinterpreted as the indication of a disc-like BLR, as in the standard interpretation for the broad double-peaked emitters (e.g.Eracleous & Halpern 1994).Even if a binary is selected as a candidate at a time when one MBH was significantly brighter than the other, SPT could fail in observing a long-term frequency shift consistent with the orbital evolution of a binary if the other component undergoes a sizable rebrightening, changing the shape (and centroid) of the line.Such an evolution could result in the erroneous dismissal of real binaries if only SPT is performed, while the binary nature of the system would be easily identifiable by FUVT.-As discussed in section 4, the stochastic nature of the shortterm variability of each MBH accretion disc can occasionally result in strong cross-correlations even if the observed AGN hosts a real MBHB.If FUVT is applied to large spectroscopic catalogues to identify new MBHB candidates, a fraction of the MBHBs in the data are expected not to be missed.It is difficult to estimate a solid missed detection fraction because of the small number of single MBHs that we use to characterize the 'control sample'.However, taking the probabilities we estimated (with a priori assumptions on the because the functional form of the probability distribution was decided a priori for this exercise.The case on the left represents a class of objects in which the max cross-correlation in the bulk of the broad line can be higher or lower than that of the single MBH depending on the assumed period.The right panel shows the only case in which the single MBH has a value of CCF of higher than any of its sibling mock binaries, but the cross-correlation of the mocks can be quite high, up to ∼ 0.75 for the shortest period.Please note that the y axis changes from panel to panel in order to emphasise the differences between the different cases. 1.0 0.5 0.0 0.5 1.0 shape of p(CC), which will be verified when enough data become available), and assuming a 'missed' alarm threshold of p ≥ 0.01 (corresponding to CC ≥ 0.70), we get 10 systems out of 36 (≈ 28% of missed binaries), including binaries with a period shorter than 300 yr, or 9 out of 18 when considering only larger 'solid' binaries (see discussion in section 2).Such shortcomings can be circumvented when FUVT is applied to data taken from a significantly longer campaign with equally frequent observations, motivated for example by the objective to independently identify a source as a promising MBHB candidate.In this case, if the first observational campaign finds a high correlation in the red and blue part of the BELs, the test should be considered inconclusive, and a new campaign should be performed.Only after a few campaigns (for a total of a few years, in the worst-case scenario) can the binary hypothesis be ruled out (unless a low correlation period is found, in which case the binary nature is indeed confirmed).However, we note that the missed detection fraction estimated here cannot be used to infer the statistics of the global MBHB population, because in this first method paper we are considering only the idealised edge-on and Doppler-maximising configuration for equal-mass binaries.
A broader analysis of the binary parameter space, including a statistical estimate of the minimum number of observational campaigns needed to disprove the MBHB scenario with at a given confidence level, is postponed to a future study.-Finally, we stress that all the currently spectroscopically selected MBHB candidates have masses of ∼ > 10 8 M ⊙ (e.g.Tsalmantza et al. 2011;Eracleous et al. 2012), and no candidates in the 10 5 −10 7 M ⊙ range of interest for the future gravitational wave interferometer LISA have yet been identified.This might be due to a selection effect specific of the traditional spectroscopic search.Indeed, in order to be selected as MBHB candidates, the BELs are typically required to be shifted by ∼ > 1000 km s −1 (Tsalmantza et al. 2011;Eracleous et al. 2012) with respect to the host rest frame (traced by the narrow emission lines) in order to prevent the inclusion of single MBHs with slightly asymmetric BLRs.This threshold can be higher than the maximum velocity for which the secondary of a MBHB can retain its own BLR12 : v 2 ≈ 480 km s −1 × M 2 10 6 M ⊙ 0.24 f −0.26 Edd f (q) −0.5 , where f (q) = q 1/3 (1 + q) 0.6q 2/3 + ln (1 + q 1/3 ) , for small secondary masses and binaries of close to equal mass, unless the Eddington ratio is small.For example, M 2 = 10 6 M ⊙ and q = 1 would require f Edd ∼ < 0.01 to have a maximum secondary velocity of v 2 ≈ 1000 km s −1 .For such small masses and low f Edd , collecting a spectrum with a sufficiently high signal-to-noise ratio to perform the test might be prohibitively challenging.Nevertheless, FUVT could still succeed in identifying a MBHB (as it does for three mock MBHBs with periods of 1000 yr), if the random fluctuations of the two MBHs do not appear correlated by chance (see the discussion above).
As a final note, we stress that the procedure we discuss could, in principle, identify false-positive MBHB if the structure of the BLR around a single MBH is sufficiently complex, with strong asymmetries, including for example long-lived spiral waves (e.g.Storchi-Bergmann et al. 2003) or hot spots (e.g.Fries et al. 2023).If a large and representative distribution of CCF were available for bona fide symmetric single MBHs, FUVT could identify subpopulations of outliers.Their deviation from a reference CCF distribution would provide an estimate of the probability that such systems do not belong to the 'normal' single MBH distribution (as exemplified in section 4), allowing the selection of candidates to be further analysed to discriminate between the two (asymmetric single BLR and MBHB) scenarios.We speculate that the inclusion of the information on the crosscorrelation of the two (blue and red) parts of the BELs with the observed continuum could inform us about the most plausible scenario.A study quantifying false positives, that is, AGN with single asymmetric BELs, and the development of tests to identify them are ongoing.Nevertheless, we believe that currently the most limiting factor is the paucity of systems for which RM studies have been performed.Indeed, in order to estimate the probability that our mock MBHBs are outliers of the single MBH distribution, in section 4 we have to assume an a priori CCF distribution for single MBHs, as nine measurements are barely enough to characterise the distribution.However, the situation will improve significantly thanks to future large RM campaigns (e.g. the black hole mapper program of the SDSS-V, Almeida et al. 2023), adding approximately 1000 systems to the current sample, allowing far better determination of the distribution of CCF for single MBHs.With a large enough sample, one could even begin to identify subpopulations within the single MBH class, ultimately leading to a more accurate identification of binary MBHs candidates.An alternative improvement could be achieved with a dedicated RM campaign on double-peaked emitters 13 .Such a study could provide an additional test as to the nature of individual double-peaked emitters (see also Eracleous et al. 1997;Liu et al. 2016, for different tests) and, even if all the systems prove to be single MBHs with disk-like BELs (e.g.Eracleous & Halpern 1994), could serve as a comparison for future candidates.

Fig. 1 .
Fig. 1.NGC 4748: Maximum value of the cross correlation between the blue and red parts of the broad Hβ line as a function of the dividing wavelength λ.The peak of the cross-correlation is obtained when the line is divided into two parts of equal flux, which is when the signal-tonoise ratio is the highest for both the red and blue sides.

Fig. 2 .
Fig. 2. NGC 4748: Original light curve of the entire Hβ line (blue line), the same light curve shifted by half of the observational period (with periodic boundary conditions, red line), and the sum of the two (rigidly shifted towards lower fluxes for visualisation purposes, green line), mimicking a MBHB whose components evolve independently on short timescales.

Fig. 3 .
Fig. 3. NGC 4748.Upper panel: Broad Hβ component of a mock spectrum shifted redwards (red line) and bluewards (blue line) to mimic the orbital velocity of an equal-mass MBHB with orbital period τ orb = 50 yr.The green line shows the composite line associated to the mock MBHB.Lower panel: Composite broad Hβ lines for mock MBHBs with τ orb = 50, 100, 300, and 1000 yr (green, red, blue, and magenta lines, respectively).The Doppler shifts are evaluated assuming an edgeon circular binary with the MBH velocities perfectly aligned with the line of sight.

Fig. 4 .
Fig. 4. Same as figure 1, but for the four mock MBHBs obtained from NGC 4748.The green, red, cyan, and magenta points with uncertainties are for binaries with τ orb = 50, 100, 300, and 1000 yr, respectively.The original values and uncertainties of NGC 4748 are shown in blue for comparison.'Solid' binaries are shown with dashed lines, while dotted lines show shorter-period binaries.

Fig. 5 .
Fig.5.Same as figure4but for NGC 5548 and the associated mock binaries (left panel) and for Mrk 202 (right panel).The case on the left represents a class of objects in which the max cross-correlation in the bulk of the broad line can be higher or lower than that of the single MBH depending on the assumed period.The right panel shows the only case in which the single MBH has a value of CCF of higher than any of its sibling mock binaries, but the cross-correlation of the mocks can be quite high, up to ∼ 0.75 for the shortest period.Please note that the y axis changes from panel to panel in order to emphasise the differences between the different cases.
Fig.6.Normalised cumulative distribution of the observed single MBHs (blue histogram) and of our mock MBHBs (red histrogram).The solid black line shows the median cumulative probability that a system will have a CCF of lower than a given value (obtained from the single MBH data under the assumptions specified in the main text).The teal band indicates the 90% credible region for the cumulative distribution of CCF for single MBHs.The vertical red (blue) ticks highlight the CCF values for all of the 36 mock binaries (9 single MBHs).The numerical values of all of the CCF and cumulative probabilities are listed in table 1.