title: TDI-infinity: time-delay interferometry without delays
authors: Michele Vallisneri Jean-Baptiste Bayle Stanislav Babak Antoine Petiteau
text: 
  
  
  
  
  
  Jet Propulsion Laboratory, California Institute of Technology, Pasadena CA 91109 USA
  
  
  AstroParticule et Cosmologie (APC), Université de Paris/CNRS, 75013 Paris, France
  
  
  
  
  The space-based gravitational-wave observatory LISA relies on a form of synthetic interferometry (time-delay interferometry, or TDI) where the otherwise overwhelming laser phase noise is canceled by linear combinations of appropriately delayed phase measurements. These observables grow in length and complexity as the realistic features of the LISA orbits are taken into account.
  In this paper we outline an implicit formulation of TDI where we write the LISA likelihood directly in terms of the basic phase measurements, and we marginalize over the laser phase noises in the limit of infinite laser-noise variance. Equivalently, we rely on TDI observables that are defined numerically (rather than algebraically) from a discrete-filter representation of the laser propagation delays.
  Our method generalizes to any time dependence of the armlengths; it simplifies the modeling of gravitational-wave signals; and it allows a straightforward treatment of data gaps and missing measurements.  
  
  
  
  
  
  TDI-∞: time-delay interferometry without delays
      Antoine Petiteau
      December 30, 2023
  ===============================================
  
  
  
  
  
  § INTRODUCTION
  
  
  Interferometry is not indispensable to the experiments that seek to detect gravitational waves (GWs) by monitoring the displacement of freely falling test masses.
  Sensitivity is set by disturbances to free fall (acceleration noise) and by the precision of the distance measurement (position noise).
  Interferometry becomes crucial when the ruler by which distance is measured (typically, the wavelength of an infrared laser) is not sufficiently stable at the GW frequencies of interest, so its fluctuations must be canceled out interferometrically.
  Such is the case of the space GW observatory LISA <cit.>, in which laser frequency noise is several orders of magnitude larger than acceleration and position noise. LISA is a strange kind of interferometer, where the laser-noise-canceling interferometric observables are not realized physically, but reconstituted in post-processing from the set of one- or two-way phase measurements between the pairs of spacecraft in the constellation.
  
  This reconstitution is known as time-delay interferometry (TDI, <cit.>) because the phase measurements are delayed by multiples of the LISA armlengths before they are combined.
  While this design is ingenious, and indeed seminal to the LISA concept, it is inconvenient for data analysis. GWs are completely buried in the laser-noise-dominated phase measurements, so both the phase data and the theoretical GW templates must undergo a time-domain transformation, which is computationally costly and time dependent (because the LISA armlengths are changing continuously).
  TDI compounds the difficulties of data reduction: for instance, gaps in the phase measurements are replicated multiple times across the TDI time series <cit.>; clock noise requires a complicated subtraction procedure <cit.>;
  
  stretching LISA armlengths couple noisily to the interpolation of the delays <cit.>; and more.
  
  In this paper we propose that the LISA phase measurements can be analyzed directly for the purpose of GW detection and parameter estimation, without transforming them explicitly into TDI observables with analytical forms derived a priori.
  Equivalently, the TDI observables can be computed numerically from the LISA armlengths and plugged directly into the calculation of the likelihood, the essential ingredient of GW data analysis.
  The mathematical counterparts of these qualitative statements are the formulation of a joint probability density for the phase measurements and laser noises, which is marginalized with respect to the latter to yield the likelihood used in data analysis; and the definition of TDI observables as the null-space basis vectors of the design matrix that models the delayed appearance of the laser noises in the phase measurements.
  We refer to these vectors as “TDI-∞” observables, since they cancel laser noise for any time dependence of the LISA armlengths, whereas “first-generation” TDI is limited to constant armlengths, “second-generation” TDI to linearly evolving armlengths with sufficiently small rates, and so on.
  
  
  
  § TOY PROBLEM AND DISCRETIZED REPRESENTATION
  
  
  
  
  
  In this paper we describe our proposed scheme in the context of a representative toy model of the LISA measurements (see Fig. <ref>). We consider a single laser c(t), propagated along arms 1 and 2 (with lengths that may be evolving with time), and reflected back by perfect mirrors; phase measurements y_1(t) and y_2(t) are performed at the origin, separately for each arm. Thus the measurement can be written as 
  
  
      y_1(t)    = c(t - ℓ_1(t)) - c(t) + n_1(t) 
      
              y_2(t)    = c(t - ℓ_2(t)) - c(t) + n_2(t)
  
  
  where ℓ_1,2(t) are the roundtrip flight times along the two arms for light pulses that are received at time t, and n_1,2(t) represent measurement noises. In terms of one-way armlengths we have
  ℓ_1(t) = L_1(t) + L_1(t - L_1(t)), with L_1 the incoming flight time along arm 1, and L_1 the outgoing flight time; the two will be different if the mirror is moving with respect to the origin. Crucially, we assume c(t) ≫ n_1,2(t). We do not model gravitational waves, but they would appear in both y_1 and y_2 with appropriate delays and geometric projections.
  
  In practical measurements, all continuous time series will be sampled discretely with sufficiently high cadence, so in what follows we adopt the language and notation of linear algebra. Doing so is appropriate also for laser noise, under the assumption that interferometric signals are filtered so that the Nyquist criterion is satisfied by the sampling.
  
  It is convenient to combine the two measurements y_1,2 and their noises n_1,2 into vectors 𝐲 and 𝐧, so we write
  
  
      𝐲 = 𝖬 𝐜 + 𝐧;
  
  
  here 𝐜 is the (discretized) laser noise time series, and 𝖬 is a design matrix that models the delayed finite differences of Eq. (<ref>) by way of fractional-delay finite-impulse-response filters.
  These very filters will be used in the post-processing of the LISA data to delay the interferometric measurements as required in TDI (see below).
  Therefore the approximation that we make by writing Eq. (<ref>) as a discrete equation is already implicitly accepted in standard usage.
  
  If we assume (without loss of generality) that the laser noises are switched on instantaneously at time t=0, and that the delays ℓ_1 and ℓ_2 are constant multiples 2 Δ t and 3 Δ t of the basic sample cadence, the application of the design matrix would look like
  
  
      (
          [ y_1(t_0); y_2(t_0); y_1(t_1); y_2(t_1); y_1(t_2); y_2(t_2); y_1(t_3); y_2(t_3); y_1(t_4); y_2(t_4);        ⋮ ])
      =
      (
          [ -1  0  0  0  0  ⋯; -1  0  0  0  0  ⋯;  0 -1  0  0  0  ⋯;  0 -1  0  0  0  ⋯;  1  0 -1  0  0  ⋯;  0  0 -1  0  0  ⋯;  0  1  0 -1  0  ⋯;  1  0  0 -1  0  ⋯;  0  0  1  0 -1  ⋯;  0  1  0  0 -1  ⋯;  ⋮  ⋮  ⋮  ⋮  ⋮  ⋱ ])
      ·(
          [ c(t_0); c(t_1); c(t_2); c(t_3); c(t_4);      ⋮ ]),
  
  
  where t_k = k Δ t, and where we have interleaved y_1 and y_2 measurements.
  If we obtain measurements at n epochs, then 𝐜 is an n-vector, 𝐲 a 2n-vector, and 𝖫 a 2n × n matrix.
  
  Fractional delays would spread out the leftmost 1s into the appropriate filter masks. In this paper we shall use delay filters based on Lagrange interpolation: that is, the m-point filter mask follows from approximating f(δ t), with 0 < δ t < 1, by evaluating the (m - 1)-order interpolating polynomial with nodes at f(-m/2+1), f(-m/2+2), …, f(-1), f(0), f(1), …, f(m/2). (In this illustration, for simplicity we have set the cadence Δ t equal to 1.) Filters with ⌊δ t ⌋≠ 0 are obtained by first shifting the nodes by that integer part, then evaluating the interpolating polynomial at δ t - ⌊δ t ⌋. These filters have the property of maximal flatness in the frequency domain at f = 0; we always use them with even m, so that fractionally delayed quantities are continuous across δ t = 1/2.
  
  
  
  Figure <ref> shows a graphical representation of the design matrix, this time for linearly evolving, non-integer delays implemented with fractional-delay masks of length m = 6.
  Odd and even rows, corresponding to y_1 and y_2, are shown as thicker and lighter lines respectively.
  The diagonal pattern of dips, common to both y_1(t) and y_2(t), corresponds to the “direct” -c(t) terms.
  The patterns below the diagonal correspond to delayed c(t - ℓ_1,2(t)) terms, as realized by way of fractional-delay filter masks, and they are seen to shift with changing ℓ_1,2(t).
  
  
  
  § CLASSICAL TDI
  
  
  In the classical TDI approach <cit.>, one derives laser-noise-free TDI observables written as linear combinations of delayed measurements y_1,2. In our toy model there is one single such observable, which we may identify with the standard Michelson combination M if ℓ_1,2 are equal and constant:
  
  
      M(t) = y_1(t) - y_2(t);
  
  
  with “first-generation TDI” X if ℓ_1,2 are unequal and constant:
  
  
      X(t) =    (y_1(t) + y_2(t - ℓ_1(t))) - 
         (y_2(t) + y_1(t - ℓ_2(t)));
  
  
  and with “second-generation TDI” X_1 if ℓ_1,2 are unequal and mildly evolving:
  
  
      X_1(t) =
      
      
      
      
      
      
      
      
                      (y_1(t) + y_2(t_,1) + y_2(t_,12) + y_1(t_,122)) - 
         (y_2(t) + y_1(t_,2) + y_1(t_,21) + y_2(t_,211)),
  
  
  where t_,1 = t - ℓ_1(t), t_,12 = t - ℓ_1(t) - ℓ_2(t - ℓ_1(t)), and so on. By inserting Eq. (<ref>) into Eqs. (<ref>)–(<ref>) one can verify that the laser noises cancel in pairs. This is trivial for M. In each line of the equations for X and X_1, the delayed-laser term of each measurement cancels the direct term on the next one; also, the direct terms of y_1(t) and y_2(t) cancel out, as do the delayed terms of the last measurements. This last cancellation is only approximate for X_1 when the delays are unequal and evolving: Taylor-expanding all laser noises, one sees that TDI cancels noise terms that are linear in ℓ̇_1,2, but not those of higher orders (terms that are O(ℓ̇^2_1,2), O(ℓ̈_1,2), and so on).
  Nevertheless, second-generation TDI is sufficient to reduce laser noise to levels compatible with the LISA requirements, as shown by experiments <cit.> and analytic and numerical studies <cit.>.
  
  This procedure has a beautiful geometric formulation in terms of synthesized interferometric paths <cit.>. It can also be formalized algebraically in terms of polynomial syzygies <cit.>.
  
  Last, it can be recast as an application of principal component analysis <cit.>—a formalism closely related to ours, and discussed further below.
  
  Going back to our linear-algebraic notation, we represent a TDI observable evaluated at times t_0, t_1, …as the vector
  
  
      𝐨 = 𝖳𝐲,
  
  
  with 𝖳 an n × 2n matrix that encodes the delays of Eqs. (<ref>)–(<ref>) by way of fractional-delay filters. Laser-noise cancellation then corresponds to
  
  
      𝖳𝖬≃ 0,
  
  
  where the cancellation is exact for M (or X) with constant and equal (unequal) ℓ_1,2, and approximate but very accurate for X_1 with evolving and unequal ℓ_1,2.
  
  GW searches and source parameter estimation proceed from the evaluation of the likelihood of the data as a function of GW parameters θ. In terms of the TDI vector 𝐨, we obtain the likelihood by postulating that the measurement noises n_1,2(t) are independent Gaussian processes with zero mean and covariance  N(t',t”), and by equating the probability of observing the measurement residual Δ𝐨 = 𝖳 (𝐲 - 𝐲_GW(θ)) to the sampling probability of the noise, appropriately mapped from 𝐲 to 𝐨:
  
  
      log p(Δ𝐨 = 𝖳Δ𝐲|θ) = 
      
      
      = -1/2Δ𝐨^† (𝖳𝖭𝖳^†)^-1Δ𝐨
      -1/2log |2 π𝖳𝖭𝖳^†|,
  
  
  where
  
  
      𝖭_(ai)(bj) = ⟨ y_a(t_i) y_b(t_j) ⟩ = δ_ab N(t_i, t_j).
  
  
  In the classical treatments of TDI, one would usually compute the spectral density S_X(f) of the TDI observable as a function of acceleration and position noise in each element of LISA, and then write the log likelihood as 
  
  
      -2  Re∫ΔX̃^*(f)  ΔX̃(f)/S_X(f) d f;
  
  
  this equation is exactly equivalent to Eq. (<ref>), where 𝖳𝖭𝖳^† plays the role of S_X(f) in the discretized time domain.
  
  
  
  § INTRODUCING TDI INFINITY
  
  
  Instead of formulating the TDI observables algebraically or geometrically by matching direct and delayed noise terms, resulting in equations similar to (<ref>)–(<ref>), we take the approach of defining the set of discretized TDI vectors by solving the matrix equation
  
  
      𝖳𝖬 = 0
  
  
  for 𝖳 given the design matrix 𝖬, which is determined by the LISA orbits and by the accurate times at which the y_1,2 are sampled. (See also Ref. <cit.>, which derives the TDI observables by solving Eq. (<ref>) in the frequency domain, under the assumption that the armlengths are constant.)
  
  The solution 𝖳̂ to Eq. (<ref>), which is unique up to affine transformations, provides a basis for the null space of 𝖬^†: any vector 𝐨̂ within the null space solves the equation 𝖬^†𝐨̂ = 0. Correspondingly, each row 𝐨̂_k of 𝖳̂ can be dotted into an observed vector 𝐲 to generate a laser-noise–free observation o_k. We refer to the rows 𝐨̂_k as TDI-∞ observables; by construction, they cancel laser noise for any time dependence of the light-propagation delays. Given that 𝖬 is a 2 n × n matrix of rank n, we obtain n linearly independent TDI-∞ observables.
  
  It should be clear from our theoretical development so far that 𝖳̂ can be used with Eq. (<ref>) to evaluate the TDI likelihood directly from the interferometric measurements 𝐲, without the additional step of computing time-delayed combinations of measurements and GW templates.
  Furthermore, the solution of Eq. (<ref>) and the computation of the inverse covariance matrix 𝖪^-1≡ (𝖳̂𝖭𝖳̂^†)^-1 can be performed offline, before the repeated evaluation of the likelihood in a search or parameter-estimation scheme. The online steps are the TDI-∞ projection Δ𝐨̂ = 𝖳̂ (𝐲 - 𝐲_GW(θ)) and the kernel product -1/2Δ𝐨̂^†𝖪^-1Δ𝐨̂.
  
  We further motivate our proposal by demonstrating that, in the LISA-appropriate limit of large laser noise, the TDI-∞ likelihood is equivalent to the likelihood written from first principles for the y measurements. That is, we can derive TDI-∞ from a complete generative model of the LISA measurements, without need to model laser-noise subtraction explicitly.
  
  Representing c(t) as a Gaussian process with mean zero and covariance function C(t',t”) (with  𝖢_ij≡ C(t_i, t_j)), we write the likelihood of 𝐜 and of the observed residuals Δ𝐲 = 𝐲 - 𝐲_GW(θ) as
  
  
      p(Δ𝐲,𝐜|θ) = |2 π𝖭|^-1/2 e^-1/2 (Δ𝐲(θ) - 𝖬𝐜)^†𝖭^-1 (Δ𝐲(θ) - 𝖬𝐜);
  
  
  integrating this likelihood with respect to 𝐜, after multiplying by their prior p(𝐜) = |2 π𝖢|^-1/2 e^-1/2𝐜^†𝖢^-1𝐜, yields the marginalized log likelihood <cit.>
  
  
      log p(Δ𝐲|θ) =    -1/2Δ𝐲^†(θ)
      (𝖭 + 𝖬𝖢𝖬^†)^-1Δ𝐲(θ) 
         - 1/2log |2 π (𝖭 + 𝖬𝖢𝖬^†)|.
  
  
  The marginalization can be seen as a probabilistic version of solving for the lasers, and then propagating the uncertainty of the solution to the remaining degrees of freedom. In Eq. (<ref>), the augmentation of the covariance matrix 𝖭 by 𝖬𝖢𝖬^†  has the effect of downweighting (or, in the limit c(t) ≫ n_1,2(t), completely projecting out) the linear combinations of the 𝐲 in which the laser noises are dominant.
  
  While Eq. (<ref>) could be used directly for GW applications, doing so carries the risk of losing numerical precision, possibly catastrophically.
  The reason is that for LISA the 𝐲 will always be strongly dominated by the laser noise 𝐜; while the specific form of the covariance matrix will (in effect) select the 𝐜-orthogonal components of the 𝐲, that projection will involve the dangerous cancellation of very large numbers.
  
  We can instead rely on Eqs. (<ref>) and (<ref>), which we show to be equivalent to Eq. (<ref>) in the limit of overwhelming laser noise.
  To realize that limit, we take 𝖢 = σ1 with σ→∞, and write the inverse Gaussian-process kernel of Eq. (<ref>) using the singular value decomposition (SVD) 𝖬 = 𝖴𝖲𝖵^†:
  
  
  
      (𝖭 + 𝖬𝖢𝖬^†)^-1    = 
      (𝖭 + σ𝖴𝖲𝖲^†𝖴^†)^-1
         =
      𝖴 (𝖴^†𝖭𝖴 + σ𝖲𝖲^†)^-1𝖴^†,
  
  
  where the second equality follows by inserting factors 𝖴𝖴^† = 𝖨 and shifting the 2n × 2n orthogonal matrix 𝖴 outside the inverse. We then refactor the second line of Eq. (<ref>) as a block-matrix product, subdividing the columns of 𝖴 as (𝖤, 𝖥), where 𝖤 spans the range of 𝖬 and 𝖥^† the null space of 𝖬^†:
  
  
      ( 𝖤 𝖥)
      ( [ 𝖤^†𝖭𝖤 + σ𝖲𝖲^†         𝖤^†𝖭𝖥;         𝖥^†𝖭𝖤         𝖥^†𝖭𝖥 ])^-1([ 𝖤^†; 𝖥^† ]).
  
  
  Using the block inverse formula, we find that all blocks are O(σ^-1) and disappear in the limit σ→∞, except for the bottom right block 𝖥^†𝖭𝖥 (the Schur complement). Thus
  
  
      lim_σ→∞log p(Δ𝐲|θ) = 
       -1/2Δ𝐲^†(θ) 𝖥 (𝖥^†𝖭𝖥)^-1𝖥^†Δ𝐲(θ) -1/2log |2 π 𝖥^†𝖭𝖥|.
  
  
  This limiting procedure is similar in spirit and mathematical detail to the marginalization over timing-model corrections in the time-domain analysis of pulsar-timing-array data <cit.>; in that case, as here, the degrees of freedom with very large variance are effectively projected out of the data vector.
  
  Now, 𝖥 is a 2n × n orthogonal matrix such that 𝖬^†𝖥 = 𝖥^†𝖬 = 0; given that the TDI-∞ matrix 𝖳̂ is full rank and that 𝖳̂𝖬 = 0, there must exist an invertible but not necessarily orthogonal matrix 𝖠 such that 𝖥^† = 𝖠𝖳. Inserting this representation in Eq. (<ref>) reproduces Eq. (<ref>), modulo an additive factor that does not depend on 𝖭.
  
  In addition to demonstrating the large-𝐜 equivalence of Eqs. (<ref>) and (<ref>), this derivation suggests that the numerical instability of Eq. (<ref>) is resolved in Eq. (<ref>), since the large components proportional to 𝖢 drop out of Eq. (<ref>), while the 𝖳̂ projection cancels out the large laser-noise contributions to the 𝐲 measurements. The projection does require sufficient measurement precision and linearity, but no more so than the computation of the delayed combinations of classical TDI.
  
  In Ref. <cit.>, Romano and Woan identify the TDI observables with the small-eigenvalue eigenvectors of the 𝐲 covariance matrix (𝖭 + 𝖬𝖢𝖬^† in our notation), and emphasize that its singular value decomposition factorizes the 𝐲 likelihood into a TDI term (a sufficient statistic for astrophysical inference), and a laser-dominated term (useful for laser-noise monitoring but not GW detection). They also recover the classical TDI expressions by analyzing the covariance matrix for integer-Δ t laser delays.
  In the limit of large laser noise, Romano and Woan's approach is equivalent to the null-space formulation discussed here: indeed, Eqs. (<ref>) and (<ref>) describe how the SVD of 𝖬 induces the factorization of the marginalized likelihood.
  Baghi and colleagues <cit.> perform the Romano–Woan eigenvector decomposition in the frequency domain, and work with the resulting 𝐲 likelihood to simultaneously fit the GW source parameters, the LISA armlengths, and the components of the covariance matrix.
  
  
  
  § THE OBSERVABLES OF TDI INFINITY
  
  
  
  
  
  The standard linear-algebra approach to computing a basis for the null space of a matrix consists of factorizing it by SVD and then selecting the rows of the right factor that correspond to the null singular values.
  These rows are orthogonal by construction, and in general they are dense across the matrix.
  This means that TDI-∞ vectors obtained from the SVD are nonlocal: they span the length of the data, instead of being restricted to a few multiples of ℓ_1,2 as the classical TDI observables.
  The left panel of Fig. <ref> shows such SVD vectors for the evolving-ℓ_1,2 design matrix of Fig. <ref>.
  The particular SVD implementation used here ( from LAPACK <cit.>) results in some diagonal structure in the second half of the plotted timespan, but little regularity overall.
  
  Dense TDI-∞ vectors do not necessarily lead to loss of precision, but they certainly have other disadvantages. They obfuscate the time dependence of GW signals and instrument noise; they make it hard to analyze data in chunks; and they guarantee that the offline and online phases of likelihood evaluation have maximum computational complexities O(n^3) and O(n^2) (where 2 n is the length of the vector 𝐲).
  
  Fortunately, applied mathematicians have developed algorithms that generate banded basis matrices for the null space of sparse banded matrices such as 𝖬^†. 
  One such algorithm is the turnback method, originally suggested by Topçu <cit.> in the context of the matrix force method for linear elastic analysis, and further refined in Refs. <cit.>. The turnback method begins with the standard LU decomposition <cit.> of the sparse matrix, followed by a number of triangular factorizations of its submatrices. The resulting basis vectors are not orthogonal, but they are concentrated along the diagonal.
  
  We experimented with the method using an implementation kindly provided by Thuy Van Dang and Keck Voon Ling <cit.>. For integer delays, turnback basis vectors reproduce exactly the observable X of first-generation TDI (Eq. (<ref>)); for constant fractional delays, they have the same bandwidth as X (the larger of the ℓ_1 and ℓ_2 plus half of the filter delay width), but with smoother structure; for evolving delays, as shown in the center panel of Fig. <ref>, they have bandwidth comparable to X rather than to X_1 (which has extent ≳ 3 ×ℓ_1,2), again with smoother structure.
  These statements discount the first few vectors at the top of the center panel, which have trivial structure because c(t) = 0 for t < t_0, so delayed laser-noise terms do not affect the y_1,2 until t_7.
  
  We note that in a realistic data-reduction scenario the bandwidth is likely to be dominated by the length m of the fractional-delay filters, which will span several LISA armlengths to achieve the required interpolation accuracy.
  While the offline phase of likelihood evaluation has again complexity O(n^3), as required by the turnback algorithm, the banded structure of 𝖳 and therefore 𝖳^†𝖭𝖳 may allow the optimization of the online phase to complexities lower than O(n^2).
  
  In the right panel of Fig. <ref> we demonstrate that the TDI-∞ approach automatically takes data gaps into account. Here we have modeled missing phase measurements at four epochs by removing eight rows of the 2 n × n design matrix. The solution of Eq. (<ref>) with the turnback method yields (n - 8) TDI-∞ vectors that combine the measurements around the gap, shown as the blank vertical band in the plot. Notably, the ℓ_1,2(t) used here allow for two vectors that bridge the disruption. Near the gap, the bandwidth of the observables increases by ∼ 50%. By contrast, in this example the second-generation TDI observable X_1 would be unavailable for more than 20 epochs, since it requires phase observations spanning ≳ 3 ×ℓ_1,2. This advantage washes out for longer gaps.
  
  
  
  § DISCUSSION
  
  
  By way of a toy model of LISA interferometry, we have offered a proof of principle that the LISA GW data analysis can be formulated and performed directly in terms of the phase measurements, without recourse to the analytical observables of classical TDI.
  This approach leads to the numerically defined observables of TDI-∞, which cancel laser noise for any time dependence of the armlengths, and which can be conveniently time-localized to bandwidths comparable to or smaller than those of second-generation TDI.
  The scheme has several additional advantages:
  
  
  
    * There is no need to select a set of analytical TDI observables, model their power-spectral densities, and track their data quality;
  
    * GW theoretical templates can be computed directly for the simpler phase measurements rather than the more complicated TDI observables, or even for the basic GW strain polarizations, and then projected to the phase measurements;
  
    * Measurement gaps are handled automatically and gracefully, including the shift between one, two, and three independent combinations when four, five, and six LISA laser links are available;
  
    * The link to the computation of matrix null-space bases, a linear-algebra problem with many practical applications, raises the possibility of adopting new sophisticated algorithms <cit.>, including parallelized or streaming variants suited to GPUs.
  
  
  While our toy model is extremely idealized and therefore limited,  we believe these advantages warrant a detailed investigation of the numerical implementation of TDI-∞ and of its implications for the LISA system, which we leave for future work.
  
  MV is thankful to John Baker for posing the challenge of distributing the LISA data as the set of the “𝐲”, and to Alvin Chua for helpful interactions and conversations.
  MV is especially grateful to CEA/IPhT and APC for hospitality during a 2019 sabbatical that led to the inception of this work.
  MV was supported by the JPL RTD program, and JBB by a NASA postdoctoral fellowship administered by USRA.
  Part of this research was carried out at the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration (80NM0018D0004).
  Copyright 2020. All rights reserved.
  
  apsrev4-1
  
  
  