EURASIP Journal on Applied Signal Processing 2005:15, 2413–2425 c ○ 2005 Hindawi Publishing Corporation Adapted Method for Separating Kinetic SZ Signal from Primary CMB Fluctuations

In this first attempt to extract a map of the kinetic Sunyaev-Zel'dovich (KSZ) temperature fluctuations from the cosmic microwave background (CMB) anisotropies, we use a method which is based on simple and minimal assumptions. We first focus on the intrinsic limitations of the method due to the cosmological signal itself. We demonstrate using simulated maps that the KSZ reconstructed maps are in quite good agreement with the original input signal with a correlation coefficient between original and reconstructed maps of on average, and an error on the standard deviation of the reconstructed KSZ map of only % on average. To achieve these results, our method is based on the fact that some first-step component separation provides us with (i) a map of Compton parameters for the thermal Sunyaev-Zel'dovich (TSZ) effect of galaxy clusters, and (ii) a map of temperature fluctuations which is the sum of primary CMB and KSZ signals. Our method takes benefit from the spatial correlation between KSZ and TSZ effects which are both due to the same galaxy clusters. This correlation allows us to use the TSZ map as a spatial template in order to mask, in the map, the pixels where the clusters must have imprinted an SZ fluctuation. In practice, a series of TSZ thresholds is defined and for each threshold, we estimate the corresponding KSZ signal by interpolating the CMB fluctuations on the masked pixels. The series of estimated KSZ maps is finally used to reconstruct the KSZ map through the minimisation of a criterion taking into account two statistical properties of the KSZ signal (KSZ dominates over primary anisotropies at small scales, KSZ fluctuations are non-Gaussian distributed). We show that the results are quite sensitive to the effect of beam convolution, especially for large beams, and to the corruption by instrumental noise.


INTRODUCTION
The cosmic microwave background (CMB) temperature anisotropies field encloses so-called primary anisotropies, directly related to the initial density fluctuations at early stages of the universe, and so-called secondary anisotropies generated after matter and radiation decoupled. The secondary anisotropies arise from the interaction of the CMB photons with gravitational potential wells or with ionised gas along their way towards us. More "local" contributions to the CMB signal are due to foreground emissions from our galaxy and from distant galaxies. One of the major goals of observational cosmology is to use the CMB anisotropies to probe the cosmological model mainly through cosmological parameter estimation. This is already performed by a number of groups using ground-based and balloon-borne experiments such as TOCO [1], BOOMERanG [2], MAXIMA [3], DASI [4], and Archeops [5], which achieved a firm detection of the so-called "first peak" in the CMB anisotropy angular power spectrum at the degree scale. This detection was recently confirmed by the WMAP satellite [6]. A series of small scale CMB experiments (e.g., VSA [7], CBI [8], ACBAR [9]) showed rather conclusive evidence for a second and a third peak. The positions, heights, and widths of these features in the angular power spectrum already give us a good idea of the cosmological model.
It is clear however that such constraints necessitate the "cleanest" possible cosmological signal, or in other words they need the best possible monitoring of contaminating signals such as secondary anisotropies or foreground emissions. This is the objective of the component separation for CMB observations. In most cases, the different contributing signals to CMB anisotropies exhibit both different frequency (ν) and spectral (in Fourier or spherical harmonic domains) dependences. As a consequence CMB experiments often observe at several frequencies to be able to separate the different astrophysical signals. Numerous methods were adapted and developed for the CMB problem like the Wiener filtering, maximum entropy, independent component analysis, and so forth [10,11,12,13,14]. All gave very satisfactory results and showed clearly that one can extract the CMB primary signal from the observed mixture. Obviously the success of the separation methods greatly depend on how different from each other the signals, in the observed mixture, are. Signals sharing with the primary anisotropies the same frequency dependence and/or the same power spectrum will be badly detected or even undetected by the separation methods mentioned above.
We present here a new method optimised to extract, from the primary signal, temperature fluctuations which have the same frequency dependence and are associated with a major secondary effect, namely, the kinetic Sunyaev-Zel'dovich (KSZ) effect (for details see [15]). Our method is based on a two-step strategy in which we derive the best estimate of the CMB signal on masked pixels by interpolation, and then we deduce the best estimate of the KSZ map by minimisation. We show that we retrieve the KSZ signal, in the best possible way, in terms of its amplitude, its distribution and its power spectrum provided we use (i) a well-chosen spatial template for the masked pixels, and (ii) adapted signal processing techniques for both interpolating and minimising. For the first point, we use the TSZ maps as a template since both TSZ and KSZ are associated with the same objects (clusters of galaxies). For the second point, namely, the signal processing techniques, there are several issues to address in order to optimise the extraction of the KSZ signal from the primary CMB. We have thus to make sure at each step that we use adapted techniques. At the interpolation on the masked pixels defined by the template, we can obviously use several methods like for example constrained 2D realisations of the underlying CMB (sensitive to our knowledge of the CMB through the confidence intervals on the cosmological parameters), textures [16] (sensitive to the morphological description of the processes), or simply cubic B-spline methods. The latter, which we use in the present study, is a classical and robust method giving reliable results especially when we set proper continuity and boundary conditions. In order to separate the primary CMB signal from KSZ, several methods can be used such as principal component analysis or least square minimisation. Here, we choose to minimise over statistical criteria. At the minimisation step, it is thus important to use the tools which emphasize the different statistical characteristics of the signal (power, non-Gaussian signatures). Among the numerous possible tools used to exhibit non-Gaussian signatures (higher-order moments in real space, biand tri-spectrum (e.g., [17]), Minkowski functionals (e.g., [18]), higher criticism (e.g., [19]), the multiscale transforms seem to be the most satisfactory and will thus be used in the following. In the same spirit, we use the most sensitive non-Gaussian estimator among the coefficients in a biorthogonal wavelet transform, namely, the diagonal coefficients. We test and develop our method on a set of numerical simulations that will be described in Section 2. In Section 3, we detail the methodology adopted to separate KSZ from primary CMB signal. We focus on the way our method is intrinsically limited by pure cosmological signals primary CMB and SZ effect. We perform some sensitivity tests in Section 4 and explore the effects of beam convolution and instrumental white noise on our results. Finally, we discuss our results in Section 5.

THE ASTROPHYSICAL SIGNALS
Among all secondary anisotropies, the dominant contribution to CMB signal comes from the Sunyaev-Zel'dovich (SZ) effect [20,21] which represents the inverse Compton scattering of CMB photons by free electrons in ionised and hot intracluster gas. This so-called thermal SZ (TSZ) effect, whose amplitude is given by the Compton parameter y, is the integral of the pressure along the line of sight. The inverse Compton effect moves the CMB photons from the lower to the higher frequencies of the spectrum. This results in a peculiar spectral signature with a brightness decrement at long wavelengths and an increment at short wavelengths ( Figure 1, solid line). When the galaxy cluster moves with respect to the CMB rest frame, with a peculiar radial velocity v r , a Doppler shift called the kinetic SZ (KSZ) effect generates a temperature anisotropies with the same spectral signature as the primary CMB fluctuations, at first order ( Figure 1, dashed line).
The importance of the SZ effect for cosmology has been recognized very early (see reviews by [22,23]). It is a powerful tool to detect high-redshift galaxy clusters since it is redshift independent. In combination with X-ray observations, it can be used to determine the Hubble constant and probe the intracluster gas distribution. Moreover, the KSZ effect may be one of the best ways of measuring the cluster peculiar velocities by combining thermal and kinetic effects [21]. Formally, the KSZ can be distinguished from the TSZ effect due to their different frequency dependences. The KSZ intensity reaches its maximum at ∼218 GHz, just where the TSZ intensity is zero (Figure 1). In practice, very few measurements of the peculiar velocities were attempted [24,25].
With usual component separation techniques it has been shown that the TSZ signal can be extracted from the CMB rather easily (its frequency dependence is quite different from a black-body spectrum) while the KSZ component remains indistinguishable from the primary CMB anisotropies due to same frequency distribution. In early works, [26] used an optimal filtering (Wiener), with a spatial filter derived from X-ray observations of galaxy clusters, that minimises confusion with CMB, and [27] used a matched filter optimised on simulated data and independent of the underlying CMB model.
We simulate (512 × 512 pixels) primary CMB (∆T/ T) RMS CMB = σ CMB = 1.9 ×10 −5 ), TSZ (mean σ y = 1.17 ×10 −5 , i.e., mean (∆T/T) RMS TSZ = −2.34 ×10 −5 at 2 mm) and KSZ (mean (∆T/T) RMS KSZ = σ KSZ = 1.85 ×10 −6 ) maps with a pixelsize of 1.5 arcmin. A precise description of the SZ simulations is given in [28]. The CMB signal is a Gaussian field whose power spectrum is computed from an inflationary flat, low matter density, model. The KSZ effect induces temperature fluctuations given by δ KSZ T = (∆T/T) KSZ = −(v r /c)τ, with c and τ the velocity of light and the cluster Thomson optical depth. The primary CMB and the KSZ anisotropies having the same spectral shape (at first order), we construct maps of thermodynamic temperature fluctuations, δ T , by adding the two signals δ T = (∆T/T) KSZ + (∆T/T) CMB . We are thus left with two simulated datasets of pure cosmological signals, one consisting of the temperature fluctuation maps (CMB+KSZ) and the other consisting of the Compton parameter maps, y, for the TSZ effect. For our study, we restrict the analysis to 15 simulated maps which span a representative range of amplitudes for all signals.
Note that in "real life," it is the multifrequency CMB experiments that can provide us, after classical component separation, with two sets of maps. One contains both CMB and KSZ temperature fluctuations, as they are indistinguishable, and the second consists of Compton parameter maps associated with the TSZ effect.

METHOD FOR SEPARATING KSZ FROM CMB SIGNAL
From the two available types of maps (y maps for TSZ and δ T maps for CMB + KSZ). Our goal is to obtain the best possible estimate of the KSZ buried in the dominant CMB signal. We benefit for this from the fact that TSZ and KSZ features are spatially correlated (e.g., [29,30]). The spatial correlation simply means that both effects are due to galaxy clusters. Therefore, where TSZ signal is present, so are KSZ fluctuations regardless of their signs or amplitudes. Conversely, where the TSZ fluctuations are absent, so are the KSZ fluctuations and the signal at that position in the δ T map is associated with the CMB primary anisotropies only. Note however that clusters at rest with respect to CMB (v r = 0) will have no KSZ signal.
From this simple statement we build a two-step component separation strategy in which (i) we first derive the best estimate of the CMB map by interpolation, and (ii) consequently deduce the best estimate of the KSZ map by minimisation.

Estimating primary CMB anisotropies
The basic idea in order to estimate the primary CMB anisotropies is to use the TSZ map as a mask to select in the δ T map the pixels where only primary anisotropies are present. These pixels contain no TSZ fluctuations in the y map. The rest of the pixels in the δ T map are masked pixels. We then interpolate the δ T signal on these masked pixels with the constraint that pixels where the signal is associated with only primary CMB, keep their values after the interpolation. We therefore obtain an estimated primary CMB map. It is worth noting that y map is an observable quantity that is rather easy to obtain from multifrequency observations due to its spectral signature. This is what makes it useful for the mask definition.
Formally, the KSZ map can then be estimated simply by computing the difference between the original unmasked δ T map and the primary CMB map obtained from the interpolation.

Interpolation of the masked pixels
The reconstruction of the KSZ map depends on the performances of the interpolation. We use the method described in [31] and consider the problem of the minimisation of a general criterion written as where f is an input image, u is the desired solution, w ≥ 0 is a map of space-varying weights, and d x and d y are the horizontal and vertical gradient operators, respectively. The second space-invariant term in (1) is a membrane spline regulariser; the amount of smoothness is controlled by the parameter λ.
Taking the partial derivative of (1) with respect to u, we find that u is the solution of the differential equation where W is the diagonal weight matrix, f w = W f the weighted data vector, L is the discrete Laplacian operator, and A = W + λL a symmetric definite matrix. The inversion of (2) is achieved using a multigrid technique [32]. Typically, we need two V-cycles with two iterations in the smoothing Gauss-Seidel part of the algorithm to reach a residual of the order of 10 −6 . The interpolation of the primary CMB map can be achieved by setting the weights to zero where the data are missing, that is in the masked pixels, and to one elsewhere and by resolving (2). The value of λ then determines the tightness of the fit at the known data points (unmasked pixels), while the surface u is interpolated such that the value of the Laplacian of u is zero elsewhere. In the present work, we impose a low value for λ so that the recovered values at the known data points are equal to the original values. This criterion can be relaxed to take into account corruption of the data by additive white noise [31]. In this case, the optimum regularisation parameter λ can be defined as where σ 2 is the variance of the noise and E( f · L f ) denotes an estimate of the correlation between the noisy image f and its Laplacian L f . In the case of nonwhite noise, the optimal regularisation parameter λ may be determined from the data using cross-validation methods [33], or from a given measurement model of the signal + noise [34]. The performances of the interpolation are improved if the values of the Laplacian of u at the missing data points are nonzero. Moreover, the values are set such that the first and second derivatives of the interpolated signal are continuous throughout the interval. These continuity conditions characterise the cubic B-spline functions which are known for their simplicity and their performances in terms of signal reconstruction [35,36]. In practice, these conditions imply that the source term f w in (2) is modified to impose nonzero values at the points where the weights are set to zero (i.e., masked pixels). An equivalent way to solve (2) with the above-mentioned conditions is to replace the Laplacian operator L by the quadratic operator L 2 .
Obviously other interpolation methods can be proposed and used to estimate the CMB data in the masked pixels. We could for example improve the interpolation by using textures [16]. The latter account for the morphological properties of the signal. Such method is limited by our knowledge of these characteristics. We could also think of using constrained 2D realisations of the CMB to obtain the values in the mask. This method is simple; however, it suffers from the precision to which the CMB power spectrum is estimated, or in other words the precision on the cosmological parameters used for the realisation.

Defining the masked pixels
We now define the mask, that is, how we select the missing data points. Besides the pixels that actually contain no galaxy clusters, that is, no SZ contributions, we fix a threshold value for the TSZ amplitude below which TSZ signal is considered too small to be detected. The corresponding pixels in the δ T maps are then associated only with primary CMB signal. On the contrary, above this threshold, pixels in the δ T map are considered to be the missing data points (masked pixels) that we want to interpolate. The number and location of the missing data depend on the threshold. The choice of this threshold has thus important consequences on the quality of the interpolation. When the threshold is high, the number of missing data is small and the interpolated surface is good but the selection retains only the clusters with the highest TSZ and misses the majority of clusters. In this case, we expect a low correlation coefficient between the retrieved KSZ map (obtained by subtracting the interpolated CMB map from the δ T map) and the original KSZ map. When the threshold is low, we take into account a large fraction of clusters, but the interpolated surfaces are large and the quality of the interpolation suffers from that. Moreover, the characteristic scale of the interpolated surfaces becomes of the order of the CMB fluctuations leading to "confusion effects." Since it is difficult to choose one single optimal TSZ threshold, we retrieve a set of interpolated CMB maps corresponding to a set of TSZ threshold values. The later are defined as follows: we compute the cumulative distribution function of the TSZ values in the y map and we search for the values corresponding to 15%-95% of the total number of pixels (with a step of 5%). This gives us a set of 17 threshold values. All pixels in the TSZ map that have y parameters above the thresholds are identified as missing data points in the simulated δ T map, that is, the mask.

Results
For each of the 15 simulated maps of our datasets, we obtain 17 TSZ thresholds, and thus 17 masked δ T maps. We interpolate the missing data points to recover the primary CMB signal in the masked regions. We evaluate 17 associated KSZ maps by subtracting the interpolated primary CMB maps from the total δ T map.
We compute for each of the 17 KSZ estimated maps the correlation coefficient between the original input KSZ map and the 17 estimated KSZ maps. The correlation coefficients are plotted as a function of the standard deviation of the estimated KSZ map for each of the 17 threshold values ( Figure 2). The diamonds and the dashed line represent the case where the interpolation is such that the Laplacian values are set to zero, and the triangles and the solid line are for the case in which the Laplacian values are nonzero. Figure 2a shows our best recovery case in terms of correlation coefficient. Figure 2b is for our worst case.
From Figure 2, we see that the correlation coefficient between the original and the estimated KSZ maps is higher when the Laplacian values are nonzero than when they are set to zero. This is especially true for the maps with low standard deviations. The improvement due to the biharmonic operator is of the order of 20% in our worst case (Figure 2b). We will therefore use, in the following, the L 2 operator as it gives better interpolation. In addition, we see that the correlation coefficient increases when the standard deviation of the estimated KSZ map increases (i.e., TSZ threshold decreases). The correlation coefficient reaches a maximum value and then decreases for the highest KSZ standard deviations (i.e., the lowest TSZ thresholds).

Reconstructing the KSZ map
From the previous step, we obtained a set of 17 estimated KSZ maps. Now, we search for a method that gives us either the reconstructed KSZ map which is the closest to the original KSZ signal or the combination of 17 KSZ maps giving the best estimate of the original KSZ map. We compare the reconstructed maps with the original KSZ maps. This allows us to calibrate our method and thus provides us with the intrinsic limitations of the reconstructing methods.

Method
We test the decorrelation by principal component analysis (PCA). The PCA gives us a reconstructed KSZ signal which is rather close to the original. The correlation coefficient, averaged over the 15 maps, between reconstructed and original KSZ reaches 0.73. However, the standard deviation is on average smaller by almost 50% than the original. This is not satisfactory. We can also search for a linear combination of the 17 estimated KSZ maps that is the closest to the original KSZ in the sense of least squares. This minimisation is done using a standard singular value decomposition (SVD). The average correlation coefficient (over the 15 simulated input maps) between the original KSZ map and the reconstructed map is 0.8, slightly higher than the PCA result.
However, in this case, the standard deviations of the reconstructed maps are lower than the original KSZ signal by almost 25% on average. Furthermore, the results of the SVD least square minimisation depend on the set of estimated maps that are used which is clearly undesirable. The two previous attempts being not quite satisfactory, we thus need as much map-independent results as possible. We must identify a criterion, to minimise on, which should ideally give the largest possible correlation coefficient and the reconstructed KSZ maps with the closest possible standard deviations to those of the original KSZ signal. Moreover, a good minimisation criterion would characterise the KSZ signal only, excluding the primary CMB signatures. We have identified two properties of the KSZ fluctuations that fulfill this definition.
(i) The KSZ signal dominates primary CMB at high wave numbers (small angular scales). (ii) The KSZ effect is a highly non-Gaussian process contrary to primary CMB which is a Gaussian process.
The KSZ effect is due to galaxy clusters whose typical sizes are a few to a few tens of arcmin. As a result, SZ anisotropies produced either by KSZ or TSZ intervene at small angular scales where they show a maximum amplitude (Figure 3, dashed line). At those scales, primary CMB anisotropies are severely damped and the angular power spectrum decreases sharply (Figure 3, solid line). Therefore, at small scales, both the power and the statistical properties of the total δ T signal should be those of the dominant signal, that is, KSZ effect. In order to focus on the KSZ signal and also to enhance the signal-to-noise ratio, we perform multiscale wavelet decomposition of the δ T map. The above-mentioned properties remain true in the wavelet domain as it was first recognised by [37] and applied by [38]. Thus, the statistical properties of the wavelet coefficients at the lowest decomposition scale (3 arcmin) reflect the properties of the SZ effect only.
We use the decimated biorthogonal wavelet transform which decomposes a signal s as follows: where φ and ψ are, respectively, the scaling and the wavelet functions. J is the number of resolutions used in the decomposition, w j the wavelet coefficients (or details) at scale j, and c J a smooth version of s ( j = 1 corresponds to the finest scale, highest frequencies). The two-dimensional algorithm gives three wavelet subimages at each decomposition scale. Within this choice, the wavelet analysis provides us with the wavelet coefficients associated with diagonal, vertical, and horizontal details of the analysed map.
Using this tool, we have demonstrated [39,40] that the excess kurtosis of the wavelet coefficients in a biorthogonal decomposition allows us to discriminate between a Gaussian primary CMB signal and a non-Gaussian process like SZ effect better than with an orthogonal wavelet decomposition. Moreover, we have shown that diagonal details are  the most sensitive to non-Gaussian signatures (recently confirmed and explained in [41]). We therefore choose to use the diagonal details in a biorthogonal wavelet decomposition at the smallest decomposition scale to obtain the best results.
In Table 1, we compare, using the 9/7 biorthogonal filter bank [42] for the worst and best cases, the statistical properties of the diagonal details of KSZ maps and CMB + KSZ maps at the first decomposition scale (3 arcmin). We also give the values for the primary CMB maps. As expected, we note that the wavelet coefficients for KSZ and KSZ + CMB share the same statistical properties and are quite different from those of the primary CMB alone. This confirms that KSZ signal dominates over primary CMB in wavelet domain (same standard deviation means same power, c.f. Figure 3 in real space), and that non-Gaussian signatures in the KSZ + CMB maps are associated with the KSZ effect (same skewness and excess kurtosis) alone at the smallest decomposition scale (3 arcmin). Consequently, we can confidently minimise on the statistical properties of power and non-Gaussianity at the smallest decomposition scale. In practice, we choose the following criterion minimised over the 17 estimated KSZ maps (for each of the 15 maps of our dataset): where w 0 is the distribution of diagonal wavelet coefficients for the known δ T map (KSZ + CMB) and w is the distribution of diagonal wavelet coefficients for the desired solution map (the reconstructed map). M 2 and M 4 are respectively the second and the fourth moments of the wavelet coefficients. This criterion takes into account both the energy content of the coefficients, through second moment, and the non-Gaussian character, through fourth moment. We have chosen the fourth moment because it is the one for which the KSZ signal is the most sensitive to non-Gaussianity. Clearly, we might also include the third moment of the wavelet coefficients to the criterion. This would be needed in particular if we were dealing with a "skewed" signal (e.g., weak lensing signal). Taking the fourth moment in the minimisation criterion allows us in turn to focus on the reconstruction of KSZ maps excluding any skewed signal that might contribute at small scales. In addition to the conditions of power and non-Gaussian character, we make use in the minimisation process, of a nice property of the wavelet transform, which is that it preserves the spatial information. Thus instead of minimising over all wavelet coefficients of the data map (w 0 in (5)), we can minimise only over those corresponding to clusters. This enhances the non-Gaussian character and reduces the influence of other possible non-Gaussian processes that could affect the anisotropy map δ T . Note the correlation coefficient between original and reconstructed KSZ maps of ∼0.9 and the agreement between total power P real and P est .

Results
In Figure 4, we present the standard deviations of the 15 original simulated KSZ maps (triangles) and of the 15 reconstructed KSZ maps (squares) obtained by the abovementioned minimisation technique. The agreement is good even for the lowest standard deviations with an error only of the order of ∼5%. This is much smaller than what was obtained from the PCA method (∼50%) or from the least square minimisation (∼25%). Furthermore, the mean value (over the 15 original maps) of the correlation coefficient between the original and the reconstructed KSZ maps is 0.78. Multipole l (c) Figure 6: The same as Figure 5. This is our worst reconstruction case which corresponds to the original KSZ map with the lowest standard deviation. Note the low correlation coefficient 0.62.
The quality of the KSZ map reconstruction can be observed in Figures 5 and 6 which display, for our best and worst cases respectively, the histograms of the temperature fluctuations and the power spectra of both original (solid line) and reconstructed (dashed line) KSZ maps as well as the ratio of these two power spectra. Note that the ratio is close to one over a large range of multipoles (wave number in the spherical harmonic decomposition) even in the domain where the primary CMB dominates the KSZ signal (see Figure 3). We also notice the correlation coefficient between original and reconstructed KSZ maps which reaches ∼0.9 in our best case and 0.62 in our worst case. The comparison between standard deviations of original and reconstructed maps σ real and σ est also gives a global indication on how well the method works. The method allows us to obtain such results because we are able to estimate correctly the amplitude of the KSZ signal for most clusters together with their angular separation, as well as the amplitude of the background (primary CMB). This is nicely exhibited by the superposition of the cuts across reconstructed (dashed line) and original (solid line) KSZ maps, for the best and worst cases (Figures 7 and 8, resp.). The method partially fails to find broad KSZ features due to their important level of confusion with primary CMB fluctuations. Moreover, since the minimisation process is an overall procedure, relatively large features (i.e., of the order of 10 −5 in absolute ∆T/T) are occasionally poorly recovered.

SENSITIVITY TESTS
We have shown in the previous section that statistical minimisation with a well-chosen criterion gives very good recon- structions of the KSZ original maps, both in terms of correlation coefficient, power spectrum and pixel distribution. We now investigate some of the effects that can affect our results.

Amplitude of the input KSZ signal
The previous results were obtained in a specific model which predicts the amplitude of the KSZ signal and thus its ratio to primary CMB anisotropies. Obviously the KSZ amplitude can vary for many physical reasons (number of clusters, distribution of velocities, etc.). It is thus important to test what is the performances of our separation method are in response to different mixing ratios. For illustration, we take one KSZ map and add it to the same primary CMB map. The standard deviation of the KSZ signal is reduced while the CMB standard deviation is kept the same (i.e., we reduce the KSZ contribution to the 5δ T map). We reduce the standard deviation following a geometrical progression σ i = σ 0 i = 0, 6 and σ 0 = 2.5 ×10 −7 . The highest standard deviation is then σ max = 2.0 ×10 −6 which is a typical value for our dataset. At the same time, we test the sensitivity of our method to the wavelet base by comparing results obtained using two different biorthogonal wavelet bases, the 9/7 tap filter and the 6/10 tap filter [43].
The results for this test are displayed in Table 2. We first notice that results do not depend much on the wavelet basis. As expected, the quality of the reconstruction (in terms of correlation coefficient) increases with the standard deviation of the original KSZ map from 0.5 to ∼0.8. The smallest coefficients are obtained for very-low-standard deviations (< 10 −6 ). It is worth noting that decreasing the KSZ amplitude by a factor 2 (i.e., a factor 4 in power) still gives a reasonably good correlation coefficient. At the same time, the standard deviation of the reconstructed KSZ map is very close to the original even when the input KSZ signal is decreased by one order of magnitude in terms of standard deviation. This is illustrated in Figure 9 by the reconstructed power spectrum of the KSZ map.

Beam convolution
Our separation method is based on two steps; the first is the interpolation and the second is the minimisation. Obviously, when the sky is observed by an instrument, the δ T map suffers from beam dilution, which means that the signal is damped at the typical scale of the beam size. The same is true for the y map for which the damping can be even more severe since the signal is mainly at small scales. As a consequence, the definition of the mask based on the TSZ template and used in the interpolation is also affected by beam dilution. We expect that this reduces the quality of the interpolation and in turn that of the 17 estimated KSZ maps. Moreover, the minimisation criterion is based on two properties of the KSZ signal (non-Gaussian character and excess of power) as compared to the primary CMB, which are mostly true at small scales. When the δ T map is convolved by the beam instrument, the contribution from KSZ signal is reduced affecting also the statistical minimisation criterion.
All these effects depend on the size of the beam. The smallest the beam is, the less affected the recovered signal is. For a beam-size of 1.5 arc-minute (like that of some planned SZ experiments), there should be no effect on our results since our minimum resolution is 1.5 arc-minute.
To illustrate the effect of a larger beam, we have convolved our observed maps (y and δ T maps) by a Gaussian-shaped beam (for simplicity) with a size of 3 arcmin. We find that the reconstructed KSZ map is not satisfying neither in terms of the correlation coefficients (mean coefficient of 0.59), nor in terms of the average amplitudes (standard deviations of the 15 reconstructed maps are typically 40% smaller than the original), nor in terms of the power spectrum. We show in Figure 10 a cut across a reconstructed KSZ map (which is not our best case) and its original counterpart. We note that only the largest amplitude features are reconstructed but with amplitudes which are lower than the original. As expected, we find that the results get worse for larger beam sizes.
One way to improve our results in the case of large beam experiments might be to use a minimisation criterion in (5) based on other wavelet decomposition scales which should be less affected by the beam dilution. For example, the second smallest decomposition scale could be used in the case of a 3 arc-minute beam. At that scale, the non-Gaussian character of the KSZ signal is indeed still preserved (see [39]), however the power is no more dominated by KSZ but rather by the CMB primary signal. More adapted criteria should then be investigated, but they will likely require more "a priori" knowledge of both KSZ and primary CMB signals.

Noise
We illustrate possible effects of noise on our separation method by adding to the observed δ T and y maps a white noise at the pixel size whose RMS amplitude in terms of temperature fluctuation is 2 × 10 −6 . This corresponds to a noise level of about 6 µK which is the typical noise of most future SZ experiments. We note that the RMS noise level is of the order of the mean standard deviation of the original input KSZ signal. It is twice as large as the standard deviation of some KSZ maps. This not only modifies the amplitude of the fluctuations in the δ T at a given position in the map, but also significantly modifies the position of the maxima and shape of the fluctuations associated with the KSZ signal. As a consequence, the spatial correlation between TSZ and KSZ is decreased and the reconstructed KSZ signal is different from the input map (see Figure 11). The correlation coefficients between the original and reconstructed KSZ maps are obviously very low in this case with values ranging between 0.24 and 0.54. . Note the excess of near-zero values in the histogram of the estimated map (logarithmic scale). Note also the very low correlation coefficient 0.48. This is for the worst case.
As noted in Section 3.1.1, the white noise can in principle be accounted for at the interpolation stage in the regularisation parameter. Such possibilities should have to be addressed.

DISCUSSION
We presented a method for separating the KSZ signal from primary CMB anisotropies based on two steps: (1) interpolation and (2) reconstruction. In our case this corresponds to the interpolation of a correlated noise (the CMB). The KSZ reconstruction is based on a set of KSZ estimated maps obtained with a choice of TSZ thresholds (from the cumulative distribution of the pixels in the TSZ template map), more sophisticated methods optimising the series of TSZ thresholds can be proposed. Using the set of KSZ estimated maps, we can investigate several methods to reconstruct the final KSZ maps. We tested a decorrelation-based approach using the PCA. The decorrelation is a blind method whose advantage is that no a priori criteria are needed to obtain the KSZ map. However, the resulting maps are of low quality in terms of standard deviation. More sophisticated methods such as the independent component analysis [12,44,45] can be used but the results obtained need to be rescaled using external constraints.
In our study, we choose to use a reconstruction method based on a minimisation technique. We propose a minimisation criterion taking into account statistical properties of the KSZ signal: (i) KSZ dominates over primary anisotropies at small angular scales, and (ii) KSZ fluctuations follow a non-Gaussian distribution. We use the excess kurtosis of the diagonal wavelet coefficients to characterise the non-Gaussian signatures of the KSZ effect. The minimisation method gives reconstructed KSZ maps that are in quite good agreement with the original signal with an average correlation coefficient between original and reconstructed KSZ maps of 0.78, and an error of 5% on the standard deviation of reconstructed KSZ maps. The KSZ reconstruction through minimisation depends on the minimisation criteria and therefore on our knowledge of the signals. The available CMB data seem to agree on the fact that primary CMB anisotropies are Gaussian distributed at least at small scales [46,47,48]; see [49] for large scales. The KSZ effect is dominant at small scales since it is associated with galaxy clusters. We have tested our results against the relative amplitude of KSZ to primary signal. We find satisfactory results even when KSZ is twice as small (in RMS) as predicted.
The results above are for the case where only the two signals CMB and SZ are taken into account, which allows us to investigate the intrinsic limitations of the method. Additional astrophysical contributions should be partly treated in a first-step component separation (from which we obtain the observed y and δ T maps). For example if some contribution from the TSZ signal remains in the δ T map, it will act as a correlated noise. We can account for it at the interpolation stage with the additional constraint that the skewness should be zero (which is the case for the primary CMB anisotropies), or in the minimisation procedure using a generalised criterion including the skewness as well as the excess kurtosis. In the present study, we have tested for the presence of an instrumental white noise at the pixel scale with 6 µK RMS amplitude. We find that such noise level affects the KSZ map reconstruction making it difficult to recover the KSZ signal buried in CMB. In theory, instrumental noise can be taken into account in the interpolation step by relaxing the parameter λ. Another way to deal with noise is to minimise not on the non-Gaussian character of the KSZ, but rather on the statistical properties of the remainder (i.e., CMB+noise+other components) at scales where CMB dominates. We should then obtain an estimate of all the components except KSZ that can then be subtracted from the total signal. These methods will need to be investigated in the future.
Another key element of our separation method is the use of a spatial template. The choice of a spatial template is an important issue since it is used to define the mask and hence the interpolated regions. The template should then be the closest possible to the signal. In our case, the optimal choice is the TSZ signal itself as it allows us to evaluate the temperature fluctuations associated with KSZ in the map without resorting to the knowledge or the measurement of cluster parameters. However, the beam dilution caused by observation suppresses the signal at small scales and can significantly affect the results especially for large beam sizes. (Note that our previous results are equivalent to a 1.5 arc-minute beamsize.) One way around the problem is to resort to multiscale minimisation criteria at the reconstruction step; we will investigate this question in the future. Olivier Forni is a planetologist at the Institut d'Astrophysique Spatiale in Orsay (France). His main research activities deal with the evolution of the planets and satellites of the solar system. Recently he has been working on the statistical properties of the cosmic microwave background (CMB) and of the secondary anisotropies by means of multiscale transform analysis. He also got involved in component separation techniques in order to improve the detection of low-power signatures and to analyse hyperspectral infrared data on Mars. Nabila Aghanim is a cosmologist at the Institut d'Astrophysique Spatiale in Orsay (France). Her main research interests are cosmic microwave background (CMB) and large-scale structure. She has been working during the last ten years on the statistical characterisation of CMB temperature anisotropies through power spectrum analyses and higher-order moments of wavelet coefficients. She naturally got invloved and interested in signal processing techniques in order to improve the detection of low signal-to-noise ratios such as those associated with secondary anisotropies and separate them from the primary signal.