ON THE RELIABILITY OF POLARIZATION ESTIMATION USING ROTATION MEASURE SYNTHESIS

It is often necessary to recover the Faraday dispersion function, F, by deconvolving $\tilde{F}$ from the RMTF, particularly if the source contains emission at multiple Faraday depths and if the S/N of the observations is high.⁴ To this end, we implemented a variation of the CLEAN algorithm (Högbom 1974; Heald et al. 2009) to handle the fact that both $\tilde{F}$ and the RMTF are complex valued. It searches for the highest peak in $|\tilde{F}|$ over a specified range of ϕ and removes a suitably scaled copy of the RMTF centered on its corresponding Faraday depth. The algorithm iterated over successive peaks until either a maximum number of allowed iterations was reached or the rms of the fluctuations in $\tilde{F}$ in the search region matched the rms of $\tilde{F}$ exterior to the search region. We used a loop gain of γ = 1.0 in the CLEAN algorithm. Lower loop gain values, as normally used in two-dimensional image deconvolution, were trialed, but it was found that since these stop recovering the polarization signal once it reaches the noise floor, solutions in which γ < 1 are used systematically underestimate the polarization amplitude by at least 1σ.

The RM synthesis code was supplemented by several additional features. Because the functions $\tilde{F}$ and R were computed on a grid of Faraday depths much finer than the width of R, it is sometimes desirable to group together CLEAN components that are very tightly clustered together and identify them as members of a single component. By grouping components in this way, we could obtain the χ² estimate for each independent component. We introduced a Faraday depth bunching tolerance, Δϕ, such that any two components with depths ϕ_j and ϕ_k were merged if |ϕ_j − ϕ_k| ⩽ Δϕ. A tolerance Δϕ = 5 rad m⁻² was found to be suitable for the data set examined in the present paper.

The code outputs a list of clean components ranked in order of amplitude, together with their corresponding Faraday depth and polarization vector ${\cal P}_i$ evaluated at the mean-square wavelength λ₀. The code also computed the reduced χ² associated with the polarization solution as a function of the number of clean components used. To be specific, for measurements with standard deviations σ_Q and σ_U in Q and U, respectively, and clean components with Faraday depths ϕ_c1, ϕ_c2, ... and associated polarization vectors ${{\cal P}_c}_1,{{\cal P}_c}_2, \ldots$, the reduced χ² associated with a fit to the data by including the first n clean components is

$$\begin{eqnarray} \chi ^2 &=& \frac{1}{2 N - 3 n} \sum _i^N \frac{[Q_i - Q({\phi _c}_1,\ldots, {\phi _c}_n;{{\cal P}_c}_1,\ldots,{{\cal P}_c}_n;\lambda _i) ]^2}{\sigma _Q^2}\nonumber\\ &&+ \frac{[U_i - U({\phi _c}_1,\ldots, {\phi _c}_n;{{\cal P}_c}_1,\ldots,{{\cal P}_c}_n;\lambda _i) ]^2}{\sigma _{U}^2}, \end{eqnarray} \tag{ 7 }$$

where

$$\begin{eqnarray} &&Q({\phi _c}_1,\ldots, {\phi _c}_n;{{\cal P}_c}_1,\ldots,{{\cal P}_c}_n;\lambda _k)\nonumber\\ && + i U({\phi _c}_1,\ldots, {\phi _c}_n;{{\cal P}_c}_1,\ldots,{{\cal P}_c}_n;\lambda _k)\nonumber\\ &&\quad= \sum _{j=1}^n {{\cal P}_c}_j e^{2 i \lambda _k^2 {\phi _c}_j}. \end{eqnarray} \tag{ 8 }$$

This χ² statistic was found to be an excellent estimator of the number of significant clean components required to model the polarization behavior of each source. It was also useful to identify cases in which RM synthesis solution differed significantly from the polarization behavior of the source; this highlighted cases in which the polarization is not noise-like but where its behavior is chiefly governed by a mechanism other than Faraday rotation (Section 5 discusses specific examples of this).

We performed RM synthesis on the 1005 sources in the Cen A field to investigate the reliability of the derived Faraday depths as a function of S/N. Several recent papers adopt an S/N of 7 as the threshold below which the results of RM synthesis are regarded as unreliable and are not reported (Feain et al. 2009; Law et al. 2011). It is therefore propitious to investigate whether this threshold represents a sharp barrier to the believability of RM synthesis results or whether results below this level only decline gradually in accuracy. The large size of our sample enables us to empirically derive the level of confidence associated with RM synthesis solutions over a large range of S/N.

Figures 1 and 2 show the distribution of derived RMs for varying S/N thresholds. For the purposes of this analysis, we have used only the single most significant Faraday depth derived by RM synthesis; the small number of sources in which there are multiple significant RMs has a negligible influence on the statistics (see Section 5). The signal was derived from the value of |F| at the most significant Faraday depth. Some RFI had to be flagged in the time-based data, resulting in a small variation in S/N between channels. This also varied from source to source. The noise, σ, is taken to be the average rms noise per channel for each source divided by the square root of the weighted number of channels used. To simplify the theoretical analysis in the following sections, we have assumed equal noise per channel, but it was then essential to correct for the effect of this non-uniform noise in the data by including the non-uniform weighting in the S/N estimates and in the estimates of the number of independent channels. We estimate this average rms noise per channel from the observations using the fact that the rms errors add in quadrature. Feain et al. (2009) provide specifics on the estimation of errors for the 8 MHz Q and U channels, which were interpolated from the correlator output. The effect of RFI flagging, which resulted in non-uniform noise varying from source to source, had to be included in the analysis. Even though these effects are small and subtle, it is essential to account for them when analyzing the polarization statistics.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

It is evident that all the observed RMs for high-S/N detections are in the range −150 to +50 rad m⁻², and we assume that the lower S/N detections are drawn from the same distribution. Polarization detections with RMs outside this range are most likely to be due to noise. At the lower S/N regime of interest here, the expected distribution is broadened by noise, so we increase the false RM range by 2σ and use the range −210 to +110 rad m⁻², which is appropriate down to signal strengths of 5σ. If Faraday depths associated with false detections are evenly distributed over the depth search range—as appears to be the case based on Figures 1 and 2—there is only a 4% chance that any given Faraday depth found by RM synthesis will fall in the range of expected real detections by accident. Table 1 gives the distribution of "true" and "false" polarization detections as a function of S/N.

Some RMs outside this range could be real (see Section 4.5), and while this does not significantly affect the statistical analysis, individual sources with high intrinsic RM are of interest (see Section 5). This population is rare in the flux density range of our sample since only one of the 309 sources detected above 6σ possesses an intrinsic RM magnitude larger than 200 rad m⁻².

The one high RM source in the 6  <  S/N  <  7 bin, 132446−460218, merits comment. This source has a clear double peak in the RM synthesis spectrum. The two peaks have nearly equal polarization: one is at 297 rad m⁻² but the second at −59 rad m⁻² does fall within the expected RM range. The S/N difference in amplitude is <10%. As discussed at the end of Section 5, the statistical treatment of a second RM component in the presence of a significant polarized signal is well beyond the analysis in this paper, and this may be a correct but rare detection.

The traditional method to determine the RM has been to determine the position angle of polarization in each frequency channel and then to make a scalar least-squares fit to the position angles. Figure 3 (top) shows the result of this analysis when applied to sources in the Cen A polarization catalog. When compared to the RM synthesis solution for the same sources (bottom panel in Figure 3), it can be seen that a small number of spurious high RM values are obtained. Only those sources with integrated polarization signals S/N > 7, as identified by RM synthesis, were used in the analysis. It is evident that there is a tail in the distribution of the Faraday depths derived by fitting polarization position angles. Of the 251 sources with S/N > 7 used in the analysis, 16 of these fall outside the range [ − 210, 110] rad m⁻². This can be explained by noting that with the lower S/N in each channel the least-squares fit for weak sources can find solutions that are randomly distributed over the RM space included in the least-squares solution. This high RM tail is reminiscent of the tail of high RMs sometimes seen in existing RM catalogs, e.g., see the comparison of the NVSS catalogs and the Kronberg & Newton–McGee compilation in Pshirkov et al. (2011). The small number of genuine high RMs seen in our data set (<1%) indicates that this class of radio source is rare in surveys at 1.4 GHz, presumably because sources with high intrinsic RM are more likely to be strongly depolarized.

Zoom In Zoom Out Reset image size

Taylor et al. (2009) report the RMs of 37,543 sources from the NVSS based on Very Large Array (VLA) snapshot images in two 42 MHz wide bands centered on 1364.9 MHz and 1435.1 MHz. It is opportune to investigate the accuracy of this catalog given that it currently represents one of the largest resources of RMs and is rapidly being adopted as a de facto standard in analyses of the Galactic and extragalactic magnetic fields (Harvey-Smith et al. 2011; Law et al. 2011; Oppermann et al. 2011; Govoni et al. 2010; McClure-Griffiths et al. 2010; Schnitzeler 2010; Stasyszyn et al. 2010).

We used the Cen A polarization catalog to investigate the accuracy of this technique used to derive NVSS RMs. The data were obtained at a comparable frequency, and we have collapsed them into two spectral channels with similar center frequencies of 1364 and 1436 MHz. As for the NVSS sample, the range of Faraday depths in our sample is sufficiently small that phase wrap ambiguities will not occur.⁵ The RM results for the two-band measurements are in full agreement with the RM synthesis results. The agreement with the two-band RM synthesis results is unsurprising given that the bands are contiguous.

Although RM synthesis outperforms traditional frequency-channel-based position-angle fitting algorithms such as the ones discussed above, it remains to be seen whether it makes the best use of the polarization information available. We address this by considering the output of a brute-force minimum least-squares search for detected polarization at any RM on all the sources, subject to the assumption that the source contains only a single polarized component that is assumed to have a flat spectrum. To be specific, we implemented a code that seeks the parameters P, χ₀, and ϕ (polarization amplitude, position angle, and Faraday depth) that minimize the penalty function,

$$\begin{eqnarray} \chi ^2 &=& \frac{1}{2 N - 3} \sum _j^N \frac{[Q_j - {\rm Re}[p_j(P,\chi _0,\phi) ]^2}{\sigma _Q^2}\nonumber\\ &&+ \frac{[U_j - {\rm Im}[p_j(P,\chi _0,\phi) ]^2}{\sigma _U^2}, \end{eqnarray} \tag{ 9 }$$

where the model for the polarization of the j th channel is

$$\begin{eqnarray} &&p_j(P,\chi _0,\phi) = P \exp \left[ 2 i \left(\chi _0 + \lambda _j^2 \phi\right)\right]. \end{eqnarray} \tag{ 10 }$$

Figure 4 compares the solutions from the RM synthesis with the brute-force least-squares fit approach. The two methods agree very well in both RM and polarization amplitude. At high S/N the solutions are identical. The increased scatter at lower S/N is not surprising as the two methods find different solutions, but the gap (along the Faraday depth axis) surrounding the identical solution indicates that both methods converge to the identical solution if the solutions are close. There is no systematic change in the ratio of signal amplitude between the two methods with S/N, and the average ratio for the two methods is 1.0.

Zoom In Zoom Out Reset image size

The behavior of the solutions for the sources with a "false" detection based on the RM has a number of undesirable properties. The brute-force solution has a much narrower distribution of RMs at low S/N (± 500 rad m⁻²). The amplitudes of the bad solutions are lower and the reduced χ² worse, indicating that it does not find the global minimum for noise-dominated signals.

In Section 4.2, it is demonstrated that the noise does fall as N^−1/2 for RM synthesis. The naive expectation would be that it falls similarly for traditional methods, since the estimated standard deviation in the slope of a regression line scales as N^−1/2. However, in practice this dependence is more complicated for traditional methods in a way that is hard to analyze; nπ phase ambiguities in the position angle introduce additional degrees of freedom into the fitting problem, and the fit error depends in detail on how the channels are spaced.

Here, we derive the distribution of polarization amplitudes recovered by RM synthesis from pure noise. The input measurements Q_j and U_j are therefore assumed to possess the properties of thermal noise and are independent and normally distributed with zero mean and noise per spectral channel σ.

Recall that RM synthesis de-rotates the channelized measurements of Q and U for a range of different trial RMs, ϕ. The resultant summed (complex) number is

$$\begin{eqnarray} {\cal P} &=& \frac{1}{N} \sum _j^N P_j e^{i \psi _j} e^{2 i \lambda _j^2 \phi }. \end{eqnarray} \tag{ 11 }$$

The polarization vector derived by RM synthesis is the value of ${\cal P}$ for which the choice of ϕ gives a maximum in the value of P. The task here is to determine the distribution of the amplitude of ${\cal P}$. If ϕ is given, it is straightforward to recover the polarization distribution. However, an additional complication arises because the RM synthesis algorithm is allowed the extra degree of freedom to chose any value of ϕ that maximizes P. The freedom to rotate the measured polarization vectors over a range of trial values of ϕ ensures that the polarization vector amplitude returned by RM synthesis is at least as large as that found in the ϕ = 0 case. This is because the algorithm may de-rotate the measured polarization vectors with any value of ϕ that improves the coherence of the sum in Equation (11).

There is a direct analogy between the noise characteristics of RM synthesis and those encountered in interferometric radio imaging, which involves deriving the distribution of noise amplitudes from (complex) visibilities (see Section 9.3 of Thompson et al. 2001). In both cases one is attempting to "wind up" complex valued quantities whose real and imaginary components are both normally distributed. In the case of RM synthesis, they are wound up for various values of Faraday depth. In the interferometric imaging case the complex-valued interferometric visibilities are searched across a range of fringe rotation values, corresponding to the possible position of the source. Thus, the two treatments are mathematically identical. We elaborate on this connection in the analysis below.

4.1.1. Constant ϕ

First consider the simple, well-known case in which we explicitly hold the RM at a fixed value. Since the input data are purely noise-like, the actual value of ϕ is unimportant: any noise-like polarization vector that is transformed by de-rotation to any nonzero value of ϕ still possesses the same statistical noise-like properties. Here, for the purpose of simplicity, we set ϕ = 0. The polarization found by RM synthesis is

$$\begin{eqnarray} &&{\cal P}_{\phi =0} = \frac{1}{N} \left[ \sum _j^N Q_j + i \sum _j^N U_j \right]. \end{eqnarray} \tag{ 12 }$$

The statistics of the sum of N independent normally distributed random variates is also normally distributed, with standard deviation $\sigma _{Q,U} \sqrt{N}$. Thus, S_Q = N⁻¹∑^N_jQ_j also follows a normal distribution with zero mean and standard deviation $\sigma _{Q,U}/ \sqrt{N}$. An identical relation holds for S_U = N⁻¹∑^N_jU_j. It is then straightforward to show that the probability distribution for the polarization amplitude when there is no signal present and ϕ = 0 is just a Rayleigh distribution,

$$\begin{eqnarray} &&\vspace*{-2pt} p_P(P)dP = \frac{P}{\Sigma ^2} \exp \left(- \frac{P^2}{2 \Sigma ^2} \right) d P, \end{eqnarray} \tag{ 13 }$$

where we define Σ² = σ²/N. This result is well known in the analogous interferometric imaging case pertaining to noise in very long baseline interferometry (VLBI) observations, where it corresponds to the statistics of the fringe amplitude for a specific fringe rotation when no signal is present (Thompson et al. 2001, p. 318).

4.1.2. Unconstrained ϕ

Now consider the general case in which RM synthesis performs a fit to the polarization data that is unconstrained in Faraday depth, ϕ. This consists of two separate steps. The first is to find the specific value of ϕ that maximizes the value of P for a given data set. The second step is to find the distribution of P corresponding to these maxima.

The square of the polarization amplitude is

$$\begin{eqnarray} \vspace*{-2pt} {\cal PP^*} &=& \sum _{k,m}^N P_k P_m e^{i (\psi _k - \psi _m)} e^{2 i \phi (\lambda _k^2 - \lambda _m^2)} \nonumber \\ &=& \sum _{k,m}^N P_k P_m \cos \left[ \psi _k - \psi _m + 2 \phi \left(\lambda _k^2 - \lambda _m^2\right)\right]. \end{eqnarray} \tag{ 14 }$$

The second line follows because the sine function has odd symmetry, so the imaginary part of every pair k, m in the sum is canceled by every pair m, k.

The RM ϕ should obey the following condition at the maximal value of P:

$$\begin{eqnarray} &&\vspace*{-2pt} \frac{\partial ({\cal PP^*})}{\partial \phi } = 0 = \sum _{k,m}^N 2 i \big(\lambda _k^2 - \lambda _m^2\big) P_k P_m e^{i (\psi _k - \psi _m)} e^{2 i \phi (\lambda _k^2 - \lambda _m^2)}.\nonumber\\ \end{eqnarray} \tag{ 15 }$$

Since the imaginary terms are antisymmetric with respect to exchange of the indices k and m, we see that the imaginary part of this sum is zero for any value of ϕ. So only the real part of the sum yields a nontrivial constraint on ϕ:

$$\begin{eqnarray} &&\vspace*{-2pt} 0 = \sum _{k,m}^N \big(\lambda _k^2 -\lambda _m^2\big) P_k P_m \sin \big[\psi _k - \psi _m + 2 \phi \big(\lambda _k^2 - \lambda _m^2\big) \big].\nonumber\\ \end{eqnarray} \tag{ 16 }$$

The foregoing relations show that the distribution of ${\cal P}$ exhibits a detailed dependence on the spectral coverage (i.e., the values of λ²_k). This indicates that the exact location and amplitude of the noise peak depend on the λ²_k, and that the noise distribution for any given spectral setup is best found empirically.

It is possible to construct a simple approximate analytic formulation independent of the detailed spectral coverage. This result has already been derived in the analogous case of the noise properties in VLBI images, where one is concerned with the statistics of the fringe amplitude across (Thompson et al. 2001, p. 211–319). We briefly summarize the argument here and place it in the context of RM synthesis.

The statistics of the polarization amplitude at any one value of ϕ obey a Rayleigh distribution. At any given value of ϕ the probability that we detect a polarization amplitude lower than some prescribed threshold amplitude, say, Z, is the cumulative probability distribution corresponding to Equation (13), namely, p(Z ⩽ P) = 1 − exp (− Z²/2Σ²). If the RM synthesis search is conducted over a large range of possible Faraday depths, −ϕ₀ < ϕ < ϕ₀, and the width of the amplitude of the RMTF is approximately Δϕ, there are N ≈ 2ϕ₀/Δϕ independent Faraday depth "pixels," or trials, over which the RM search is being performed. The probability that we do not detect a polarization amplitude above some amplitude Z after these N trials is p^N. The probability of finding $Z \le |{\cal P|}$ is just the cumulative probability distribution function of the amplitude of the polarization vector,

$$\begin{eqnarray} &&{\rm CDF}_{n} (Z) = \left[ 1 - \exp \left(- \frac{Z^2}{2 \Sigma ^2} \right) \right]^N, \end{eqnarray} \tag{ 17 }$$

where the subscript n emphasizes that the distribution is associated with a purely noise-like signal. The corresponding probability density function is

$$\begin{eqnarray} &&p_n(Z) = \frac{N Z}{\Sigma ^2} \left[ 1 - \exp \left(- \frac{Z^2}{2 \Sigma ^2} \right) \right]^{N-1} \exp \left(- \frac{Z^2}{2 \Sigma ^2} \right). \nonumber \\ \end{eqnarray} \tag{ 18 }$$

This expression describes the expected distribution of the polarization amplitude when the input data are purely noise-like and in which the RM is searched over approximately N-independent Faraday depth pixels. It is identical to the characteristics of the noise in a search for fringes from N-independent samples (Thompson et al. 2001, Equation (9.57)). The exact noise distribution pertaining to the Cen A data set is discussed at the conclusion of this section. We pre-empt this by noting that we find empirically that the noise distribution for the Cen A observations (shown in Figure 6) differs slightly from the analytic formula derived here.

The behavior of the noise distribution has a number of important implications for the statistics of polarization amplitudes determined by RM synthesis. Equation (18) shows that the distribution of polarization amplitudes derived from noise narrows as N increases (for constant Σ), and the peak of the distribution shifts to increasing values of P/Σ. The location of the peak is determined by the solution to the transcendental equation,

$$\begin{eqnarray} &&P^2 [ N -\exp (P^2/2 \Sigma ^2)] = \Sigma ^2 [ 1 - \exp (P^2/2 \Sigma ^2)]. \end{eqnarray} \tag{ 19 }$$

This has solutions P/Σ = 1.93, 2.24, 2.84, 3.07, and 3.54, respectively, for N = 5, 10, 50, 100, and 500.

Thus, the threshold S/N, P/Σ, for the believability of a polarization signal extracted using RM synthesis exhibits an important, albeit weak, dependence on the number of independent trials, N. The distribution peaks sharply on the low P/Σ side of the distribution, implying that one expects very few polarization amplitudes with S/N values less than 2.0 once N exceeds 10. For the analysis in this paper N = 30.

There is an important distinction between the polarization distribution obtained from a purely noise-like signal from RM synthesis and the polarization distribution derived when an intrinsic (nonzero) polarization signal is present. The distinction arises because, in the absence of an intrinsic polarization signal, the RM synthesis algorithm picks out the polarization vector associated with the Faraday depth that maximizes the length of the "Faraday-de-rotated" noise vector. However, this freedom to optimize the length of the polarization noise vector is removed provided that the intrinsic signal is sufficiently large that the RM synthesis code finds the highest amplitude at the correct Faraday depth of the intrinsic signal (i.e., there is no other Faraday depth at which a purely noise-like signal yields a higher polarization amplitude). In this instance the RM synthesis code is constrained to de-rotate the signal at the Faraday depth of the intrinsic polarization signal.

We compute the probability distribution of the polarization signal subject to the assumption that the intrinsic polarization is nonzero and that the RM synthesis code de-rotates the signal at the Faraday depth associated with the intrinsic signal.⁶ Suppose we have a nonzero intrinsic signal whose Faraday de-rotated channel measurements of the linear polarization are ${\cal P}_0=(Q_0, U_0)$. Associated with each spectral channel there is also a noise contribution. The noise follows a Gaussian distribution. It is also uncorrelated with the intrinsic signal and thus its statistical properties remain identical after de-rotation at the Faraday depth of the intrinsic signal. We write the de-rotated noise polarization contribution for the ith spectral channel as ΔQ_i and ΔU_i. Then the polarization vector derived by RM synthesis is

$$\begin{eqnarray} &&{\cal P} = (Q_0 + S_Q, U_0 + S_U ), \end{eqnarray} \tag{ 20 }$$

where S_Q = N⁻¹∑^N_iΔQ_i and S_U = N⁻¹∑^N_iΔU_i (as defined previously). The polarization vector, ${\cal P}$, therefore obeys the statistics of a random walk with a constant offset ${\cal P}_0=(Q_0,U_0)$. Thus, the polarization amplitude in the presence of an intrinsic signal follows the Rician distribution:

$$\begin{eqnarray} &&p_{s}(Z=P \vert P_0) = \frac{Z}{\Sigma ^2} \exp \left[ - \frac{Z^2 + P_0^2}{2 \Sigma ^2}\right] I_0 \left(\frac{Z | P_0|}{\Sigma ^2} \right), \end{eqnarray} \tag{ 21 }$$

where I₀(x) is the modified Bessel function of the first kind of order zero and the subscript s denotes that the distribution is valid only when a nonzero polarization signal is present so that the signal is de-rotated at the Faraday depth corresponding to the intrinsic signal, ${\cal P}_0$. The associated cumulative probability distribution function (i.e., the probability of Z ⩽ P) is

$$\begin{eqnarray} &&{\rm CDF}_s (Z \vert P_0) = 1 - Q_1 \left(\frac{ P_0}{\Sigma }, \frac{Z}{\Sigma } \right), \end{eqnarray} \tag{ 22 }$$

where Q_m(a, b) is Marcum's Q function. This result is analogous to the probability distribution of the interferometric fringe amplitude in which one measures a nonzero signal in the presence of noise (Thompson et al. 2001, Equation (9.37)).

Some remarks about the treatment of polarization noise bias are in order. The polarization amplitude distribution in Equation (21) demonstrates that the measured amplitude is, on average, biased by noise. In particular, one has

$$\begin{eqnarray} &&\langle P^2 \rangle = P_0^2 + 2 \Sigma ^2. \end{eqnarray} \tag{ 23 }$$

Thus, the appropriate noise bias to subtract from a single measurement of P² derived on the basis of RM synthesis is twice the channel-integrated noise variance, Σ². While correct on average, bias subtraction is clearly inexact for any given source; in as much as the noise can assume a range of values, so is there a range of values of intrinsic P₀ that are consistent with a given measured P. Our analysis shows that it is entirely inappropriate to simply remove a noise bias from the data and treat the resulting polarization amplitudes as if they were distributed with a Gaussian error distribution.

In the two foregoing subsections, we have derived the probability distributions of the polarization amplitude both when a signal is absent and when one is present and its amplitude is sufficiently large that it dominates the noise (i.e., provided that P₀ exceeds some threshold).

Here, we deduce the probability distribution of the polarization amplitude that is valid for all P₀. The key to this derivation is that there is a nonzero probability that the value of P drawn from the signal distribution p_s is lower than the value that would be drawn from the noise distribution p_n. Under this circumstance, the value of P picked by the RM synthesis algorithm is the greater of the two values, namely, that drawn from the noise distribution.

The probability of obtaining a polarization amplitude of Z = P is the sum of (1) the probability of drawing Z from the noise distribution and requiring that its value exceeds the amplitude drawn from the signal distribution and (2) the probability that Z was instead drawn from the signal distribution and requiring its value exceeds the amplitude drawn from the noise distribution. Thus, one obtains the combined probability distribution,

$$\begin{eqnarray} &&p_{\rm c}(Z = P \vert \, P_0) = p_n (Z) {\rm CDF}_s(Z | P_0) + p_s (Z | P_0) {\rm CDF}_n(Z). \nonumber\\ \end{eqnarray} \tag{ 24 }$$

This distribution is valid for any value of P₀. It obeys the expected property that at low Z and P₀ the distribution is dominated by the noise distribution, p_n(Z), because CDF_n(Z) is almost zero in this region, while it is dominated by p_s(Z|P₀) at high Z because p_n(Z) is almost zero in this domain.

By way of illustration, Figure 5 compares the combined probability distribution for some specific choices of P₀ against p_s and p_n.

Zoom In Zoom Out Reset image size

One can extend the foregoing logic to determine the likelihood that a signal found by RM synthesis at a given amplitude, Z, is a genuine signal. The probability that the signal has been de-rotated to the correct value of ϕ is

$$\begin{eqnarray} &&p_{\phi \hbox{ correct}}(Z | P_0) = \frac{p_s (Z | P_0) {\rm CDF}_n(Z)}{p_c (Z | P_0)}. \end{eqnarray} \tag{ 25 }$$

It is apparent that, even for small values of Z and P₀, there is a nonzero probability that the RM synthesis code deduces the correct Faraday depth of the signal.

This is borne out by the Cen A data: close inspection of the lower panel in Figure 1 indicates that even at S/N ∼3.5 the RM synthesis algorithm preferentially finds signals within the range of physically acceptable values of ϕ.

The foregoing formalism provides an obvious prescription to determine the distribution of intrinsic polarization amplitudes, N(P₀), based on the measured distribution of polarization amplitudes as deduced by RM synthesis, namely, N(P):

$$\begin{eqnarray} &&N(P) = \int d P_0 \,\, p(P \vert P_0) N (P_0). \end{eqnarray} \tag{ 26 }$$

The solution of Equation (26) to find N(P₀) is formally an ill-posed problem, but it is of a form that is often encountered in numerical optimization theory.

It is instructive to demonstrate how the solution of Equation (26) is readily amenable to numerical solution. The distribution of intrinsic polarization amplitudes, N(P₀), must be solved by sampling onto a grid with resolution, say, Δ|P|. We represent this distribution by the vector y_j. The corresponding distribution of measured polarization amplitudes is discretized onto a grid in P with identical resolution, and this is represented by the vector x_i. If one similarly discretizes p_c(Z | P₀) into a matrix whose i, jth element corresponds to the parameters Z_i and P_0j, with resolution ΔP in both dimensions, the integral in Equation (26) can be cast in the form

$$\begin{eqnarray} &&x_i = \sum _j A_{ij} y_j \Delta P. \end{eqnarray} \tag{ 27 }$$

This class of problems is amenable to solution by iterative deconvolution and expectation-maximization algorithms.

A simpler but more robust approach is to forward model the expected distribution $N(|{\cal P}_0|)$ to determine whether it is consistent with the observed amplitude distribution.

As an illustration of the foregoing formalism, we apply it to determine the distribution of intrinsic polarization amplitudes in the Cen A data set.

In order to model the effect of noise in this data set precisely, we determined the noise distribution p_n empirically based on the specific spectral coverage of the ATCA observations. The same RM synthesis code used to process the Cen A data was employed. We synthesized 50, 000 purely noise-like sources where, for each source, synthetic Q and U data for each spectral channel were drawn from a normal distribution of zero mean and unit standard deviation and assigned to spectral channels with frequencies identical to those used in the Cen A observations.

For each source the "signal" was extracted at the value of ϕ at which the amplitude in the Faraday dispersion spectrum was a maximum. The statistics of the polarization amplitude at a fixed Faraday depth, ϕ = 0, were also recorded. The resulting probability distributions are shown in Figure 6. There is excellent agreement between the distribution for the ϕ = 0 case and the Rayleigh distribution discussed in Section 4.1.1. The lower panel of Figure 6 shows the amplitude distribution determined from the synthetic noise data set when the RM synthesis algorithm is allowed to find a source anywhere over the Faraday depth range, [ − 4000, +4000] rad m⁻². One expects naively, on the basis of the width of the RMTF and the range of the search in Faraday depth, that there are N ≈ 29 independent Faraday depth trials over the search space. There is good agreement in the shape of the theoretically and empirically derived distributions for this value of N. However, it is evident that the theoretical distribution peaks at a lower S/N than the empirically derived distribution. This cannot be accounted for in terms of a change in N since, although this would shift the peak of the distribution to higher S/N, it would also narrow the distribution to a value lower than is observed. We note that, from a purely empirical viewpoint, the two distributions agree closely if the S/N axis is scaled such that the recovered S/N values from the simulated sources are systematically higher than those expected in the theoretical distribution by 8.5%. The distribution corresponding to this is also shown in Figure 6. However, the reason for such an S/N discrepancy is not obvious; the lower right panel of Figure 4 shows that there is no systematic bias in the S/N of sources recovered by the RM synthesis technique and a brute-force fit to the polarization data.

Zoom In Zoom Out Reset image size

One can directly compare the polarization amplitude distribution for the sources in the Cen A field with "false" RMs (see Section 3.1) against the distribution, p_n as determined by the synthetic pure-noise sources. This comparison is shown in Figure 7. In this plot we see excellent agreement between the polarization amplitude distribution of sources with false RMs and the expected noise distribution. However, there is an excess of sources in the region 4 < S/N < 5 that deserves further investigation. A likely explanation for this excess is that a number of these sources are mistakenly identified as false detections. It is possible that several contain real signal, but their S/N is sufficiently low that their Faraday depth is poorly localized by RM synthesis (i.e., they are only localized to about the ∼280 rad m⁻² width of the RM transfer function), thus placing the detection RM outside the nominal range [ − 210, 110] rad m⁻² considered to constitute a "good" detection. The excess of sources can be seen in Figure 2 as the concentration of sources in the 4 < S/N < 5 bin with RMs in the range [0, 500] rad m⁻². The results of Section 4.2 (and Figure 5) indicate that when some small signal is present the polarization amplitude distribution shifts to higher S/N. This suggests that some sources with supposedly "false" RMs do indeed contain intrinsic signal.

Zoom In Zoom Out Reset image size

A related statistic is p_ϕ correct, the expected fraction of sources whose Faraday depths are identified correctly as a function of P for a given P₀. We compute this quantity using the numerically derived distribution, p_c. Figure 8 displays this likelihood. This gives a direct estimate that the probability is correct at any given RM. At S/N ≈3 the probability of retrieving the correct Faraday depth is ≈50%, and it approaches 100% for S/N values between 5 and 6, depending on the value of the intrinsic polarization value, P₀.

Zoom In Zoom Out Reset image size

It is interesting to compare the differential source counts from the Cen A data in Table  1 against this curve, despite uncertainty in the intrinsic polarization amplitudes that contribute to the Cen A data. At S/N values of 3, 4, and 5, the measured fractions of good RMs are 26%, 57%, and 92%, respectively, while the predicted fractions at these S/N values are in the ranges 25%–57%, 61%–93%, and 94%–100%, respectively. In estimating the range of fractions, we have assumed that the S/N of the intrinsic polarization, P₀, is between 2.5 and 5.

4.5.1. Forward Modeling of Polarization Source Counts

It is beyond the scope of this paper to implement a complete iterative deconvolution scheme to recover the intrinsic polarization source counts, N(P₀), using Equation (26) in conjunction with p_c(Z | P₀).

We have instead developed a simple iterative forward-modeling scheme in which the intrinsic polarized source count distribution is taken to be of the form

$$\begin{eqnarray} N(P_0) = A \left\lbrace \begin{array}{@{}ll}(C/\Sigma)^{\alpha _1}(P_0 / C)^{\alpha _2}, & P_0 / \Sigma < C \\ (P_0 / \Sigma)^{\alpha _1}, & P_0 / \Sigma < C, \\ \end{array} \right. \end{eqnarray} \tag{ 28 }$$

which is a broken power-law distribution with index α₁ < 0 that turns over at an S/N of C and scales with an index −4 < α₂ < 4 below this S/N break, and where A is an overall normalization.⁷ This distribution was gridded onto the vector y with a resolution in ΔP/Σ = 0.1 and substituted into Equation (27) to find the predicted distribution of measured P values, which we label x_pred. This is compared against the actual distribution vector x_meas by computing the reduced χ² of the fit,

$$\begin{eqnarray} &&\chi ^2=\frac{1}{N_{\rm bins}-4} \sum _i^{N_{\rm bins}} \frac{|{\bf x}_{{\rm pred}_i} - {\bf x}_{{\rm meas}_i} |^2}{\sigma _i^2}, \end{eqnarray} \tag{ 29 }$$

where we take the errors to be σ²_i = x_measi (and set σ_i = 1 when $x_{{\rm meas}_i}=0$). The quantity χ² was searched for a minimum in terms of the parameters A, C, α₁, and α₂. The minimum occurs at χ² = 0.85 for (C, α₁, α₂, A) = (3.72, −2.17, −0.30, 2.60 × 10⁴), and the best-fitting solution is plotted over the measured polarization distribution in Figure 9. However, the minimum formally supports a range of indices in the range −2.9 ≲ α₁ ≲ −1.0, and we note that the values of α₁ and α₂ and C are all highly correlated. To illustrate this, Table 2 below lists the derived fit parameters for a range in the break S/N, C.

Zoom In Zoom Out Reset image size

Table 2. The Fit Quality and Associated Fit Parameters as a Function of Location in the Spectral Break, C

C	χ2	α1	α2	A
0.0	1.12	−0.96	–	119
1.0	0.93	−1.53	4	1.56 × 103
2.0	0.88	−1.90	1.01	8.09 × 103
3.0	0.86	−2.04	−0.05	1.48 × 104
3.72	0.85	−2.17	−0.30	2.60 × 104
4.0	0.85	−2.22	−0.37	3.27 × 104
5.0	0.87	−2.40	−0.54	7.39 × 104
6.0	0.91	−2.54	−0.66	1.45 × 105
7.0	0.95	−2.64	−0.74	2.50 × 105

Notes. This table illustrates the covariance between the break S/N and the spectral indices of the polarization amplitude count distribution.

Download table as:  ASCIITypeset image

The best-fitting index of the slope of the differential counts in P at high S/N, α₁ = −2.17, is very close to the value of −2.2 obtained by fitting just the stronger well-detected sources in the sample with P in the range 0.6–6 mJy. This indicates that the source counts in P continue with the same slope down to at least 0.35 mJy, and possibly lower. Note that this forward-modeling method does not require any assumptions about polarization amplitude bias at low S/N.

Following the technique used by Tucci et al. (2004) we can now estimate the change in fractional polarization as a function of flux density by comparing the source counts in I and P at the same source density. For both FIRST (White et al. 1997) and these Centaurus A field sources the differential counts in I have a slope of −1.8 in this flux range. The fractional polarization is 3.0% for the strongest sources (>200 mJy), and this increases to 3.5% for sources below 20 mJy. This can be compared with the mean fractional polarization of ATLSB sources (S_ν > 0.4 mJy) of 4.3%. Taylor et al. (2007) make a detailed analysis of polarization in the Elais N1 field. They have 786 sources with 83 polarization detections. We have 1005 sources and at least 346 polarization detections (5Σ) in the same flux density range. We agree well in the source density and polarized fraction at the higher flux density, but our fractional polarization for sources below 20 mJy (3.7%) is less than the increase to 4.8% seen by Taylor et al. (2007).