{'fraction_non_alphanumeric': 0.04426869350862778, 'fraction_numerical': 0.010271158586688579, 'mean_word_length': 4.135196905766526, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 10, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
An improved,"phase-relaxed"F-statistic for gravitational-wave data analysis

Rapidly rotating, slightly non-axisymmetric neutron stars emit nearly periodic gravitational waves (GWs), quite possibly at levels detectable by ground-based GW interferometers. We refer to these sources as"GW pulsars". For any given sky position and frequency evolution, the F-statistic is the optimal (frequentist) statistic for the detection of GW pulsars. However, in"all-sky"searches for previously unknown GW pulsars, it would be computationally intractable to calculate the (fully coherent) F-statistic at every point of a (suitably fine) grid covering the parameter space: the number of gridpoints is many orders of magnitude too large for that. Here we introduce a"phase-relaxed"F-statistic, which we denote F_pr, for incoherently combining the results of fully coherent searches over short time intervals. We estimate (very roughly) that for realistic searches, our F_pr is ~10-15% more sensitive than the"semi-coherent"F-statistic that is currently used. Moreover, as a byproduct of computing F_pr, one obtains a rough determination of the time-evolving phase offset between one's template and the true signal imbedded in the detector noise. Almost all the ingredients that go into calculating F_pr are already implemented in LAL, so we expect that relatively little additional effort would be required to develop a search code that uses F_pr.


I. INTRODUCTION
The F -statistic is the optimal frequentist statistic for the detection of nearly monochromatic gravitational waves (GWs) from a neutron star with known (or assumed) sky location and frequency evolution [1]. The basic idea behind the F -statistic is simply this: for any given sky location and frequency evolution, the set of possible GW signals forms a four-dimensional (real) vector space [1,2]. The four basis vectors are the two quadratures (sin and cos) of each of the two polarization bases, + and ×. Because the set is a vector space (not just a 4-d manifold), it is computationally trivial to maximize the likelihood function over this set; F is the maximized log-likelihood. For Gaussian noise, the probability distribution function (pdf) of F is also particularly simple: a (perhaps non-central) χ 2 distribution with 4 degrees of freedom. In their original paper by Jaranowski, Krolak & Schutz [1] (hereinafter referred to JKS), the F -statistic was derived only for the case of a single GW detector and a single GW pulsar. Cutler & Schutz [3] showed how the F -statistic can be generalized in a straightforward manner to the cases of 1) a network of detectors noise curves, and 2) an entire collection of known sources.
In practice, searches for nearly-monochromatic GWs are separated into a few different types, depending on how much is know about the source. The different types of searches can have vastly different computational requirements. While the search for a GW counterpart to a known radio pulsar is trivial in terms of computational burden, "all-sky" searches for GW pulsars with no known counterpart (and hence unknown frequency and frequency derivatives) are currently limited by the available computational power. That is, we could dig deeper into the existing data sets if we possessed either larger computational resources or more efficient algorithms. In this paper we demonstrate a way of significantly improving on the existing algorithms.
Currently, the most sensitive all-sky searches are based on the following idea [4]. For a GW pulsar with unknown frequency evolution (i.e., unknown f ,ḟ , etc.), computational power required for an optimal (i.e, fully coherent, F -based) search grows as a high-power of the total observation time, T . Therefore in practice one divides T into some number N (typically of order 10 2 ) short intervals of duration ∆T = T /N , performs an optimal, coherent search one each short interval, and then "adds up" the power from the all subintervals. More specifically, the current method is to calculate a "semi-coherent" detection statistic F sc , defined as the sum of the F values from each of the short intervals: Because calcuting each F i involves a maximization over 4 free parameters, the pdf for 2F sc is a χ 2 distribution with 4N d.o.f. But this is far more parameters than are actually needed to describe the physical system! Consider the very first interval. The imbedded GW signal is described by 4 parameters: (h 0 , ι, ψ, Φ). But the triplet (h 0 , ι, ψ) are the same for all N intervals. All that changes from interval to interval is the overall phase Φ. Assuming maximum ignorance of the signal's phase evolution, one therefore needs N − 1 additional phases to fully describe the signal. So the GW signal is fully described by only N + 3 parameters. We define F pr to be the maximized log-likelihood on this N +3-dimensional space. Our main aims in this paper are (i) to demonstrate an efficient al-gorithm for calculating F pr and (ii) to and to illustrate its superiority as a detection statistic (superior in the sense of improved ROC curves, where ROC stands for "Receiver Operating Characteristics").
We note that our work is rather similar in spirtit to a recent paper by Dergachev [6], though we believe that we have advanced the idea considerably further. We note, too, a recent paper by Pletsch [7], which addresses the followng, rather different weakness in the detection statistic Eq. (1.1). To understand the problem, let (t 0 , t 1 , · · · , t N ) be the boundary points of the N short intervals, and consider any one such boundary time, say t 50 . In Eq. (1.1), data sampled ony slightly earlier than t 50 gets combined incoherently with data sampled only slightly later than t 50 , which clearly exacts a price in sensitivity. Pletsch overcomes this problem (which stems from the rather arbitrary choice of boundary points) using his "sliding window" technique [7]. We suspect that relatively simple refinements of the detection statistic that we develop in this paper could also capitalize on Pletsch's basic insight, but we leave such refinements to future work.
The plan of this paper is as follows. In Sec. II we briefly review the fully coherent F -statistic, partly to establish notation. We generally try to align our notation with that of JKS, to ease comparison with their work. We also review the (currently used) semi-coherent F -statistic and its properties.
As further motivation for our work, in Sec. IV we consider a "warm-up" problem that we can easily treat analytically and that qualitatively has much in common with our actual problem. In Sec. V, in order to illustrate both the use of F pr , we present numerical results for one example search. For the same search, we also investigate the relative power/sensitivty of F pr versus F sc . Our conclusions are summarized in Sec. VI. This paper represents our "first-cut" analysis of the F pr statistic. There remains significant follow-up work to better elucidate the properties of F pr and to implement it in realistic, hierarchical searches. This future work is also summarized in Sec. VI.

II. REVIEW OF SIGNAL PROCESSING FOR GW PULSARS
In this section we review the rudiments of signal processing that we will require, partly to fix notation. We also review both the coherent and semi-coherent versions of the F -statistic. For simplicity of exposition, in this paper we will restrict the case of a single detector, and assume that the detector noise is stationary. The extension to the more realistic case of multiple detectors with slowly changing noise spectra is completely straightforward.

A. Mathematics of signal processing
We begin by reviewing the basic mathematics of signal processing. For more details, we refer the reader to Thorne (1987), Finn & Chernoff (1993), and/or Cutler & Flanagan (1994) [2,8,9].
Assuming that the noise is stationary and Gassuian, the noise spectral density S h (f ) determines a natural inner product (. . . | . . .) on the vector space of all detector outputs x(t): wherex(f ) andỹ(f ) denote the Fourier transforms of x(t) and y(t), and S h (f ) is the single-sided spectral density of the noise. In terms of this inner product, the probability distribution function (pdf) for the noise n(t) takes the form where N is a normalization constant. Using · · · to denote "expectation value" (over many realizations of the noise), it follows from Eq. (2.2) that In this paper, we will be concerned with waveforms h(t) that are nearly monochromatic (here meaning that their frequencies f (t) are slowly varying). In this case their inner product is equally simple in the time domain. Taking the measurement time interval to be 0 to T , we have

4)
B. The fully coherent F-statistic Next we briefly review the use of the coherent Fstatistic in GW pulsar searches. For more details we refer the reader to Cutler & Schutz (2005) [3]. Consider a nearly monochromatic GW signal from an individual source with known sky location and known frequency evolution f (t). The GW signal is then characterized by four remaining unknowns: an overall amplitude A (equivalent to the combination h 0 sin ζ sin 2 θ in the notation of JKS), two angles ι and ψ that characterize the waves' polarization (equivalent to determining the direction of the NS's spin axis), and an overall phase Φ.
The GW signal h(t) registered by the detectgor depends nonlinearly on ι, ψ, Φ, but, crucially, one can make a simple change of variables-to λ 1 , λ 2 , λ 3 , λ 4 -such that dependence of h(t) is linear in these new variables: where the four basis waveforms h α a (t) are defined by (2.6) Here Φ(t) is the waveform phase at the detector: where f (t ′ ) is the measured GW frequency at the detector at time t ′ . The measured frequency includes the Doppler effect from the detector's motion relative to the source, as well as Einstein and Shapiro delays associated with the Earth's orbit around the Sun. When the GW pulsar is in a binary, then f (t ′ ) also includes the Roemer, Einstein, and Shapiro delays associated with that binary orbit. ( We emphasize that the known-pulsar searches described here do not require that the GW pulsar be isolated, but just that there exists an accurate timing model for the emitted waves.) The F + (t) and F × (t) terms in Eq. (2.6) are the beam-pattern functions describe the detector's response to the + and × polarizations, respectively. We note that the exact form of F + (t) and F × (t) depends on one's convention for decomposing the waveform into "plus" and "cross" polarizations; a oneparameter family of choices is possible, corresponding to the freedom to rotate the axes around the line of sight. JKS follow the conventions of Bonazzola & Gourgoulhon [10]. Next we define the 4 × 4 matrix Γ ab by Because both the observation time and 1 day (the timescale on which the F α +,× (t) vary) are vastly larger than the period of the sought-for GWs (typically 10 −2 − 10 −3 s), we can replace cos 2 Φ(t), sin 2 Φ(t), and cosΦ(t)sinΦ(t) by their time-averages: Then 2F , which is defined to be twice the log of the maximized likelihood ratio, is just Using 2F as one's detection statistic satisfies the Neyman-Pearson criterion for an optimum test: it minimizes the false dismissal (FD) probability for any given false alarm (FA) probability. Writing x = n + h, and plugging into Eq. (2.12), we find where we have used Eq. (2.3) and the fact that (h | n) = 0. More generally, it is easy to show that y ≡ 2F follows a χ 2 distribution with 4 degrees of freedom (d.o.f) and non-centrality parameter ρ 2 ≡ (h|h): (2.14) As pointed out by JKS, if we use the following complexified variables, and D ≡ AB − C 2 . (Note that the A,B,C terms defined here are, in the single-detector case, larger than the A,B,C terms in JKS by a factor of the observation time T .) C. The "semi-coherent" F-statistic As mentioned above, the current method of incoherently combining the coherent results from successive intervals is just to sum of the F -statistics from all the intervals: It also easy to show that y ≡ 2F sc follows a χ 2 distribution with 4N degrees of freedom: where the non-centrality parameter ρ 2 = N i=1 ρ 2 i . In the cases of interest to us, 4N will generally be large, and then the χ 2 distribution with 4N d.o.f. can often be approximated as a Gaussian. Let y ≡ 2F , and let ρ 2 tot . For example, using this approximation (and the fact that for a Gaussian P (y), events 2.326σ above the mean occur 1% of the time), we see that the threshhold value y th that yields a 1% FA probability is We are now ready to define F pr . Basically, F pr coincides with the full matched-filtering SN R 2 , under the assumption that the manifold of waveforms is N + 3dimensional (i.e., 4 parameters for the first segment, and N − 1 for the relative phase offsets of the remaining segments). What makes F pr useful in practice is that we have also found a simple and efficient method for calculating it.

A. Motivation and definition
We begin by defining complex basis functions H + and H × by This complex representation is especially convenient for our purposes because H + and H × both transform very simply under an overall phase shift in Φ(t): under Φ(t) → Φ(t) + δ, H + and H × transform as H + (t) → e −iδ H + (t) and H × (t) → e −iδ H × (t). (Note that the minus sign in the exponent in the term e −iδ stems from the minus signs in the definitions of of H + and H × in Eq.( 3.1). For these complexified signals, our usual inner product becomes a Hermitian one; for nearly monochromatic signals near frequency f , this Hermitian inner product is given simply by Next we define Γ αβ by where α and β run over +, ×. It follows immediately that for GW data x(t), (twice the) F -statistic is given by

4)
Now imagine breaking up the full integration time T into N intervals of duration ∆T i , for i = 1, 2, · · · , N . (We expect that in practice the ∆T i will generally be of approximately the same length, but this is not required.) Next define x i (t) to be the restriction of x(t) to the i th interval; i.e, x i (t) = x(t) for t in the i th interval, and x i (t) = 0 for t outside the i th interval. Then clearly we have To motivate our definition of F pr , recall that if we had practically limitless computer power at our disposal, then the most sensitive search would be a coherent matchedfilter search over a fine grid covering the entire GWpulsar parameter space. However for "blind" GW pulsar searches (i.e., searches for GW pulsars whose sky location and/or time-changing frequency are unknown), maintaining phase coherence between the template signal and true imbedded signal, over timescales of months to years, would require an extremely fine grid on parameter space, and (one easily shows) many of orders of magnitude more computing power than is realistic [11] The basic idea behind "semi-coherent" searches is to employ a detection statistic that is less sensitive to phase decoherence across the whole observation time, which allows one to use a much coarser grid on parameter space. In effect, one sacrifices some sensitivity in the interest of computational practicality. For our phaserelaxed F -statistic, the idea is that the search-template signal should remain approximately in phase with the true, imbedded signal in each interval ∆T i -up to some constant phase "offset δ i -but that the δ i should be allowed to vary from interval to interval. That is, we replace Finally, we define (twice) F pr to be the rhs of (3.6), maximized over all phase-offsets δ i : The whole point of developing alternatives to the fully coherent F -statistic is to save on computational cost, so for our phase-relaxed F -statistic to be useful, we need a reasonably efficient way of maximizing over the δ i . In this section we demonstrate one efficient method. We demonstrate only the simplest version of this method, which we regard as basically an "existence proof" that efficient methods do exist. It should be clear by the end of this section that there are many variations on our basic method by which one might attempt to improve its efficiency, but we defer such improvements to later work.
Our method is as follows. We can simplify the appearance of the equations by defining and defining v to be the following N -dimensional vector formed out of the phase offsets: v ≡ (e iδ1 , e iδ2 , · · · e iδn ) . (3.9) Note that K ji is Hermitian (i.e., K ij = K ji * ), and that we can now re-write Eq. (3.6) as Of course, K ji is completely determined by the two complex templates H α and their inner products with the data, while our goal is to find the v i that maximize v * j K ji v i , subject to the N constraints that v i v i * = 1∀i. (We emphasize that i is not summed over in these constraints.) Put another way: v lies on the unit N-torus (i.e, the unit circle cross itself N times.) Naturally, we employ the method of Lagrange multipliers to maximize v * j K ji v i on this constraint surface. Since there are N constraints, we obtain N equations with N (real) Langrange multipliers λ j : (3.11) We emphasize that Eq. (3.11) is not an eigenvalue equation, since in general the N values λ j will all be different. Next we find it convenient to introduce a projection operator P operating on C N .
Let w = (c 1 e iδ1 , c 2 e iδ2 , · · · , c N e iδN ), where the c i are all real. Then P is defined by P w = (e iδ1 , e iδ2 , · · · , e iδN ) . (3.12) I.e., the operator P takes any vector in C N and projects it down onto the unit torus. (Note that P is not a linear operator, but it is true that P 2 = P .) Then Eq. (3.11) is clearly equivalent to the requirement that Hence the solution v is a fixed point of the operator P K. In fact, numerical experience shows that it is an attractive fixed point. That is, let v 0 be some initial guess, and then operate on it repeatedly with P K. Define (P K) 2 ≡ (P K)(P K), (P K) 3 ≡ (P K)(P K)(P K), etc. Then for v 0 sufficiently close to the true solution v 0 , we find that as m increases. In practice, we find that the convergence is quite rapid, and that the initial guess v 0 need not be particulary close to the solution v 0 . In numerical experiments (in many thousands of cases, and covering a large range of N ) we found that the following initial guess always led to converge of the iterated sequence. For each segment ∆T i , it trivial to calculate the fully coherent F i and the corresponding best-fit parameters for that segment alone: (A i , ι i , ψ i Φ i ). Then we take as our initial guess v 0 = (e iΦ1 , e iΦ2 , · · · e iΦn ) ; (3.15) i.e., the initial guess for the phase offset in each segment is the best-fit offset for that segment by itself.

IV. ANALYTIC RESULTS FOR A RELATED, WARM-UP PROBLEM
It is common sense that when one goes to solve some problem numerically, it is useful to have analytical results with which to compare it-ideally for a special case of the true problem, or, failing that, for some qualitatively similar problem. In this section we derive analytic results for the following case: Consider a vector space of waveforms that is completely described by 2 parameters per interval-so 2N parameters in all, where N is large-and consider two different searches: one search that maximizes the fit over those 2N parameters, and another, less efficient search, that begins with a 4N -dimensional vector space (in which the true, 2N -dim vector space lies), whose detection statistic is the maximized log-likelihood on the 4N -dimensional space. That is, our two detection statistics are the 2N -and 4N -dimensional F -statistics, which in this section we will denote F 2N and F 4N .
For each search, there is a threshold value ρ th such that the signal is detectable with F A = 0.01 and F D = 0.5. We can solve both problems at the same time, by considering the general M-dimensional search. Then the expectation value of F M is F M = M + ρ 2 and its standard deviation is σ M = (2M + 4ρ 2 ) 1/2 . For large M , the χ 2 function approaches a Gaussian, so we will approximate the pdf of F M as a Gaussian with this mean and standard deviation. Then the threshold for detection with F A = 0.01 is Therefore using the correct statistic allows one to see sources 19% farther away. For comparison with results in the next section, we also plot in Fig. 1 the FA vs. FD curves for the two statistics, for a range of ρ values.

V. NUMERICAL RESULTS FOR ONE EXAMPLE SEARCH
To illustrate the utility of our F pr statistic, in this section we present results for a simple, one-parameter family of examples (where the varied parameter is the strength of the embedded GW signal), and we compare the the effectiveness of F pr and F sc .
We fix the number of intervals at N = 100, and evaluate search effectiveness for signals with a range of total ρ 2 ≡ N i=1 ρ 2 i . We will eventually consider a range of ρ, but for now imagine ρ as fixed. For simplicity, in this example we will consider a case where the ρ i are the same for all i, so ρ 2 i = ρ 2 /N , and where the Γ i ab matrices are also the same for all i: Γ i ++ = 3, Γ i +× = 1 = Γ i ×+ and Γ i ×× = 1 for all i. We decompose the measured signal x i into waveform plus noise, For each i, filtering the data with H + and H x produces two complex numbers: c i + ≡ x i | H + and c i × ≡ x i | H × . Clearly, the measured signals c i α can be decomposed as We simulate the noise piece m i α by taking random draws of (pairs of) complex numbers from a Gaussian distribu-tion with covariance matrix Again for simplicity, we will consider a case where the g α i are the same for each i, modulo a random, complex phase factor. Our particular (and rather arbitrary) choice is where the ϕ i are random phases drawn uniformly from [0, 2π). One easily checks that The inclusion of the e iϕ i terms reflects our goal of modeling a case where frequency evolutions of the tempate and the true signal are so mismatched that their relative phases jump significantly and randomly from one interval to the next. Choosing the ϕ i randomly corresponds to the "worst-case scenario", where the trueversus-template phase offsets show no pattern. In practice, we expect that the situation will often be much more favorable for searches: i.e., the phase offsets might very often be well fit by some low-order polynomial in time. In a later paper we plan to investigate the extent to which the time-evolution of the offsets can be fit by a few parameters, and how that information can be exploited to speed up other parts of the search.
Given one simulated data set c i α (200 complex numbers), we compute 2F pr . We repeat for 10000 data sets to determine the distribution of 2F pr , and calculate its mean 2F pr (ρ) , standard deviation σ pr (ρ), skewness, and kurtosis. In practice, we find that the skewness and kurtosis are relatively small (as might be expected, since our N is large), so for the rest of this section we will approximate p(2F pr ; ρ 2 ) as simply a Gaussian with our measured mean and standard deviation. Given these distributions it is completely straightforward to determine the false alarm probability F A for any threshold value 2F th pr , and to calculate the false dismissal probability F D for any pair of 2F th pr and ρ 2 . We expect that, in practical searches, F th pr will find its main use in hierarchical search algorithms, in which a very coarse search at relatively low threshold identifies candidates for further examination, and these are winnowed down in successive stages [12,13]. In this context, one generally wants a fairly small FD rate (< 1%, say), so as not to lose any events, and strongly prefers a very low FA probability, to reduce the computational cost of follow-ups. With this application in mind, in Figs. 2 and 3. we plot FA as a function of FD, for several values of ρ 2 .
How much does employing F pr increase the sensitity of a blind GW pulsar search? Obtaining a useful and accurate answer to this question is much more complicated than it might initially seem, since the most sensitive known search algorithms for GW pulsars are hierarchical searches [12,13]. These searches involve involve several stages, with successive stages "ruling out" an ever-increasing fraction of the parameter space. We imagine that the most sensitive search-at some fixed, realistic computational cost-might use F pr for only some of its stages. And while calculating F pr clearly requires more floating point operations than calculating F sc , it is premature to compare these costs in a detailed way, since (i) there has as yet been no attempt to speed up our iterative relaxation scheme, and (ii) the extra information that comes with the phase offsets δ i can presumably be used to speed up other parts of the search. Despite these difficulties we can obtain a rough estimate of the sensitivity improvement from using ρ 2 pr , as follows. We have seen that given some detection statistic S, one naturally obtains a map from the total ρ 2 to curves in the F A−F D plane. Let C(S; ρ 2 ) denote that curve. Then we can look for pairs ρ 2 pr and ρ 2 sc such that C(F pr ; ρ 2 pr ) lies close to C(F sc ; ρ 2 sc ). Three such pairs are shown in Fig. 4. We see that, in the most relevant portion of F D − F A plane, C(F pr ; 140) lies close to C(F sc ; 170), C(F pr ; 160) lies close to C(F sc ; 200), and C(F pr ; 190) lies close to C(F sc ; 250). Thus, based on this example, we might estimate that using F pr rather than F sc affords an increase in sensitivity of ∼ 20 − 30% in ρ 2 , or ∼ 10 − 15% in ρ. Clearly, to obtain a more reliable estimate we should perform a Monte Carlo simulation (based on random locations of the source on the sky, and random orientations of the GW pulsar's spin axis). We plan to do this in follow-up work.
Comparing Fig. 1 to Fig. 2, we see that the sensitivity gain from replacing the detection statistic F sc with F pr is qualitatively similar to the gain from F 4N → F 2N , but that the latter gain is greater (at least for our one numerical example, and for total SNR ∼ 10). That may seem surprising, since in the former case we are eliminating 3N − 3 redundant parameters, while in the latter case we are eliminating only 2N redundant parameters. We conjecture that the main reason that the replacement F sc → F pr "buys us less" in sensitivity is the following. In calculating F 2N , the noise contributions from different intervals i still get combined incoherently. However in F pr , the maximization over the phase offsets δ i allows the noise contributions to combine coherently. Indeed for ρ i ∼ < 1, the offsets δ i are determined more by the noise than by the imbedded waveform. However we leave a thorough investigation of this effect to future work.

VI. SUMMARY AND FUTURE WORK
In this paper we defined a "phase-relaxed" F -statistic, denoted F pr , to be used in cases where the total observation time is sufficiently long that straightforward calculation of the fully coherent F -statistic over the relevant parameter space is computationally intractable. The calculation of F pr takes as input the results from coherent searches over N shorter time intervals. Our F pr coincides with the fully coherent F -statistic under the approximation that the phase offsets between template and imbedded signal are treated as an additional N −1 independent parameters. We also demonstrated one efficient, iterative method for calculating F pr . We regard our iterative method as an "existence proof" for efficient algorithms. In future work we intend to explore variations on our basic method that we suspect would lead to substantial improvements in computational cost. We illustrated the use of F pr in one simple family of examples, in which the sensitivity improvement (compared to F sc ) was shown to be ∼ 10 − 15%.
Our example was based on the "worst-case" assumption that where the phase offsets δ i are completely random. In follow-on work we intend to examine the more realistic case where the δ i are a (reasonably low-order) polynomial in time, and we plan to calculate the increased sensitivity based on a very large number of cases, in Monte Carlo fashion.
Other follow-on projects that we intend to work on include i) the development of new vetoes [14] for instrumental artifacts, since, e.g., for true GW pulsars the parameters (h 0 , ι, ψ) calculated from the first half of the observation should be consistent with those calculated from the second half; ii) the use of the δ i to quickly converge on improved estimates for the GW pulsar's frequency and spindown parameters, and iii) the optimal use of F pr in multi-stage, hierarchical searches.
Finally, while our primary interest in this paper has been the application of the F pr -statistic to GW pulsars, it has not escaped our notice that the same idea and formalism can be applied, with only trivial modifications, to searches for (quasi-circular, non-precessing) inspiraling binaries. With binaries, we expect F pr to be most useful in those cases where the observed GW signal has a very large number of cycles; e.g., searches for neutron-star binaries by proposed GW detectors that have reasonable sensitivity at ∼ 1 Hz, such as the Einstein Telescope, Decigo, or the Big Bang Observer.