Healpy fit_dipole: Cosmic Dipole Estimator

• healpy.fit_dipole is a least-squares estimator that fits a monopole+dipole model to HEALPix pixelized number-count maps.
• It accurately recovers the cosmic dipole vector in its native Cartesian form but exhibits increased variance with reduced sky coverage.
• Transforming the unbiased Cartesian estimates to spherical coordinates can introduce apparent biases in amplitude and direction.

Searching arXiv for the specified paper and closely related context. arxiv_search.query({"search_query":"ti:\"Cosmic Multipoles in Galaxy Surveys II: Comparing Different Methods in Assessing the Cosmic Dipole\" OR id:(Mittal et al., 25 Sep 2025)","start":0,"max_results":5,"sort_by":"submittedDate","sort_order":"descending"}) healpy.fit_dipole denotes, in the comparative framework of "Cosmic Multipoles in Galaxy Surveys II: Comparing Different Methods in Assessing the Cosmic Dipole" (Mittal et al., 25 Sep 2025), a monopole+dipole least-squares estimator applied to a pixelized sky map of source number densities. In that treatment, the algorithm is used to determine the number count dipole from cosmological surveys and is analyzed as a representative estimator-based alternative to a Bayesian dipole analysis. The central result is that the estimator is unbiased in its native Cartesian parameterization of the dipole vector, but that its variance becomes large for low sky coverage, so practical reliability degrades in that regime. The same study further shows that apparent biases in dipole amplitude and direction can arise when outputs are transformed to spherical parameters and averaged there, even when the underlying Cartesian estimator remains unbiased (Mittal et al., 25 Sep 2025).

1. Problem setting and fitted sky model

The algorithm is applied to a HEALPix map of pixelized number-count density,

$$\mathcal N_i \equiv \text{number density / counts in pixel } i,$$

with pixel $i$ centered on the unit-vector direction

$\hat{\mathbf p}_i = (x_i,y_i,z_i).$

In the simulations used for comparison, the sky is pixelized with nside = 64, giving

$n_{\rm pix} = 49152$

equal-area pixels. The quantity of interest is the dipole modulation of the number-count map, represented internally by four coefficients,

$(A_0, A_1, A_2, A_3),$

where $A_0$ is the monopole or mean number density and $(A_1,A_2,A_3)$ are the Cartesian dipole components of the modulated count field (Mittal et al., 25 Sep 2025).

The assumed sky model is a linear monopole+dipole field on the sphere,

$$\mathcal N_i \approx A_0 + A_1 x_i + A_2 y_i + A_3 z_i,$$

or equivalently

$$\mathcal N(\hat{\mathbf p}) = A_0 \left[1 + \mathbf{\mathcal D}\cdot \hat{\mathbf p}\right].$$

The inferred dimensionless dipole vector is

$$\mathbf{\mathcal D} = \left(\frac{A_1}{A_0}, \frac{A_2}{A_0}, \frac{A_3}{A_0}\right) \equiv (\mathcal D_x,\mathcal D_y,\mathcal D_z).$$

This formulation places the estimator squarely within the standard first-order dipolar modulation model for number-count dipoles.

The paper also situates the fitted field in the broader cosmological dipole problem. The expected kinematic number-count modulation is quoted as

$i$0

where $i$1 is the line-of-sight unit vector, $i$2 is the observer velocity, $i$3 is the speed of light, $i$4 is the slope of cumulative counts $i$5, and $i$6 is the spectral index defined by $i$7. In this setting, the increase in discordance between the number count dipole and the CMB's kinematic dipole motivates a careful comparison of methods for dipole determination.

2. Least-squares formulation and regression structure

In the study, fit_dipole is explicitly interpreted as a least-squares estimator and, under Gaussian residual assumptions with uninformative priors, as the maximum of a Gaussian likelihood and the posterior mode. Its objective function is

$i$8

with the sum taken over unmasked pixels only. The stationarity equations are written as

$i$9

where $\hat{\mathbf p}_i = (x_i,y_i,z_i).$0. This is the normal-equations form of an ordinary linear least-squares problem (Mittal et al., 25 Sep 2025).

Appendix A rewrites the fit in matrix form as

$\hat{\mathbf p}_i = (x_i,y_i,z_i).$1

with

$\hat{\mathbf p}_i = (x_i,y_i,z_i).$2

$\hat{\mathbf p}_i = (x_i,y_i,z_i).$3

and

$\hat{\mathbf p}_i = (x_i,y_i,z_i).$4

This makes the connection to linear regression explicit: the fitted basis is $\hat{\mathbf p}_i = (x_i,y_i,z_i).$5, with design-matrix rows $\hat{\mathbf p}_i = (x_i,y_i,z_i).$6.

For a full, sufficiently symmetric sky, the paper notes the approximations

$\hat{\mathbf p}_i = (x_i,y_i,z_i).$7

so that $\hat{\mathbf p}_i = (x_i,y_i,z_i).$8 becomes approximately diagonal. In that case,

$\hat{\mathbf p}_i = (x_i,y_i,z_i).$9

Substituting the exact dipole template recovers the input dipole vector in the noiseless full-sky case, which clarifies why the method is unbiased in that limit.

3. Parameterization, coordinate transformation, and simulation setup

The template map used in the simulations is

$n_{\rm pix} = 49152$0

where $n_{\rm pix} = 49152$1 is the mean counts per pixel, $n_{\rm pix} = 49152$2 is the dipole amplitude, and $n_{\rm pix} = 49152$3 is the angular separation between pixel $n_{\rm pix} = 49152$4 and the dipole direction. Writing the dipole direction as

$n_{\rm pix} = 49152$5

and the pixel direction as

$n_{\rm pix} = 49152$6

one has

$n_{\rm pix} = 49152$7

so

$n_{\rm pix} = 49152$8

This identifies

$n_{\rm pix} = 49152$9

The natural estimator outputs are therefore the Cartesian dipole components $(A_0, A_1, A_2, A_3),$0 rather than spherical amplitude-direction coordinates (Mittal et al., 25 Sep 2025).

The corresponding spherical parameters are

$(A_0, A_1, A_2, A_3),$1

$(A_0, A_1, A_2, A_3),$2

with

$(A_0, A_1, A_2, A_3),$3

The paper stresses that Cartesian space is the correct parameterization in which to study bias and covariance, because the mapping to $(A_0, A_1, A_2, A_3),$4 is nonlinear.

The main benchmark simulations adopt the following settings: nside = 64, $(A_0, A_1, A_2, A_3),$5 pixels, mean count per pixel

$(A_0, A_1, A_2, A_3),$6

fiducial dipole amplitude

$(A_0, A_1, A_2, A_3),$7

and fiducial direction

$(A_0, A_1, A_2, A_3),$8

corresponding to the CMB dipole direction. Noisy realizations are generated by Poisson sampling,

$(A_0, A_1, A_2, A_3),$9

The fitted maps are therefore pixelized number-count maps contaminated primarily by Poisson shot noise, even though the least-squares objective corresponds to Gaussian residual assumptions.

4. Masking, coverage dependence, and frequentist behavior

The estimator is run on masked HEALPix maps by setting masked pixels to UNSEEN. The study tests Galactic plane masks with

$A_0$0

equatorial polar-cap masks with radii

$A_0$1

and a discontinuous survey mask comprising many small patches below declination $A_0$2. For the discontinuous case, the mean count per pixel is increased to

$A_0$3

specifically to study sparse or discontinuous sky coverage at high number density. The estimator analysis uses

$A_0$4

mock catalogs for each mask scenario (Mittal et al., 25 Sep 2025).

The central frequentist result is that healpy.fit_dipole is unbiased in its native Cartesian dipole-vector parameterization. The paper states that

$A_0$5

is consistent with the true dipole vector across different noise levels, different sky coverages, and different mask shapes tested in the analysis. In the noiseless case, the appendix argues that the estimator recovers the true dipole exactly, with negligible variance, even with masks.

That conclusion is sharply qualified by the estimator variance. For sufficient sky coverage, the variance is small and the recovered amplitudes and directions are tightly centered on the truth. As sky coverage decreases, the variance grows strongly. The abstract emphasizes that, although the estimator gives unbiased results regardless of noise levels and sky coverage, low sky coverage leads to large variance and therefore practical unreliability. For Galactic masks, the discussion highlights a critical regime around

$A_0$6

where the estimator variance is large. For equatorial caps, large $A_0$7 values produce amplitude distributions shifted high and a bimodal direction distribution clustering near the two equatorial poles.

The discontinuous high-density benchmark is an important counterexample to a simplistic geometry-based interpretation. In that case, estimator outputs remain low-variance, both averaging methods agree with the truth, and the estimator is effectively unbiased and reliable. This suggests that sky geometry alone is not fatal; the decisive quantity is the effective information content, especially sky coverage together with count density.

5. Apparent bias in spherical coordinates

A recurrent misconception addressed in the paper is the claim that fit_dipole is intrinsically biased because average recovered amplitudes and directions can drift away from the injected truth. The study argues that this effect is largely a parameterization artifact. If one converts each realization from Cartesian dipole components to spherical parameters $A_0$8 and then averages amplitude and direction directly, apparent biases emerge: amplitudes shift high, and directions can shift toward structures associated with the mask boundaries, such as the Galactic equator or the equatorial poles (Mittal et al., 25 Sep 2025).

The explanation is the nonlinearity of the transformation

$A_0$9

A roughly Gaussian cloud in Cartesian space becomes distorted in spherical coordinates, so in general

$(A_1,A_2,A_3)$0

for nonlinear $(A_1,A_2,A_3)$1. The paper explicitly uses

$(A_1,A_2,A_3)$2

to explain why high-amplitude solutions at low latitudes can still average to the correct $(A_1,A_2,A_3)$3-component. The underlying estimator can therefore remain unbiased in Cartesian space even when averages in $(A_1,A_2,A_3)$4 appear displaced.

This distinction governs uncertainty propagation as well. The stated procedure is to estimate covariance in Cartesian components first and then transform to spherical coordinates via a Jacobian if needed. Directly reading uncertainty from transformed histograms in $(A_1,A_2,A_3)$5 does not correctly represent estimator uncertainty. In that sense, the paper treats Cartesian dipole-vector space as the native inferential space and spherical coordinates as a derived, potentially misleading summary representation.

6. Relation to Bayesian dipole analysis and methodological limits

The Bayesian comparison is based on the dipole model

$(A_1,A_2,A_3)$6

normalized over the observed sky as

$(A_1,A_2,A_3)$7

The likelihood is

$(A_1,A_2,A_3)$8

with $(A_1,A_2,A_3)$9 denoting a monopole-only model $$\mathcal N_i \approx A_0 + A_1 x_i + A_2 y_i + A_3 z_i,$$0 and $$\mathcal N_i \approx A_0 + A_1 x_i + A_2 y_i + A_3 z_i,$$1 a monopole+dipole model $$\mathcal N_i \approx A_0 + A_1 x_i + A_2 y_i + A_3 z_i,$$2. The Bayesian approach returns a posterior distribution over $$\mathcal N_i \approx A_0 + A_1 x_i + A_2 y_i + A_3 z_i,$$3 and permits model comparison through

$$\mathcal N_i \approx A_0 + A_1 x_i + A_2 y_i + A_3 z_i,$$4

with

$$\mathcal N_i \approx A_0 + A_1 x_i + A_2 y_i + A_3 z_i,$$5

quoted as “strong support” for a dipole (Mittal et al., 25 Sep 2025).

The comparison clarifies the role of fit_dipole. The estimator solves for a single best-fit vector, is naturally Cartesian, and does not provide built-in model comparison; uncertainty must be assessed by simulations or covariance estimation. By contrast, the Bayesian method uses an explicit likelihood for the count distribution on pixels, returns full posterior distributions, and can become inconclusive when the data are weak. This difference is most consequential at low sky coverage. There, fit_dipole still has unbiased mean in Cartesian space but can yield many high-amplitude realizations because of large variance, whereas Bayesian posteriors broaden and Bayes factors can fall to near zero or negative values. The study emphasizes that such inconclusive Bayesian outcomes safeguard against incorrect conclusions.

The paper also identifies several limitations and failure modes of the estimator. Low sky coverage is the dominant one, because it broadens amplitude distributions, distorts direction distributions, and undermines the reliability of single-realization estimates. Incomplete sky and higher multipoles can also induce power leakage from higher multipoles into the dipole amplitude, with worsening contamination as the masked fraction increases. Finally, fit_dipole is built for a single dipole fit to one map. A joint analysis across multiple catalogs would require a generalized likelihood or different estimator, written abstractly as

$$\mathcal N_i \approx A_0 + A_1 x_i + A_2 y_i + A_3 z_i,$$6

Accordingly, the paper’s practical recommendation is conditional rather than absolute: healpy.fit_dipole is acceptable and useful when sky coverage is sufficient and the analysis is conducted in Cartesian space, but Bayesian analysis is preferable when sky coverage is poor, when model comparison is required, or when a true multi-sample inference problem is being posed.