Wavelet Coefficient Estimation Under Random Design

Why Random Design Changes a Wavelet Coefficient

In random-design regression, observations arrive at covariates sampled from a distribution, rather than at a fixed grid. If some regions are sampled more often than others, averaging responses times a wavelet evaluates an integral weighted by the sampling density. But the usual wavelet coefficient of a regression function is an integral with respect to Lebesgue measure. These are different targets unless the design is uniform or a special cancellation occurs.

This distinction matters even when the wavelet basis is orthonormal and the regression function is smooth. The basis determines the coefficient to be estimated; the design determines how observations sample the function. We will compare the direct empirical average with an inverse-density-corrected average, and derive the variance of each.

Model and Coefficient Estimators

Let $(X_i,Y_i)$, $i=1,\ldots,n$, be independent and identically distributed observations on $[0,1]$. Assume $X_i$ has density $g$ with respect to Lebesgue measure, and write the regression model as $Y_i=f(X_i)+\varepsilon_i$, where $E[\varepsilon_i\mid X_i]=0$ and $\operatorname{Var}(\varepsilon_i\mid X_i=x)=\sigma^2(x)$. Thus $E[Y_i\mid X_i=x]=f(x)$.

Let $\psi_{jk}$ be a boundary-corrected wavelet basis function at resolution $j$ and location $k$, normalized so that $\int_0^1\psi_{jk}(x)^2\,dx=1$. Its Lebesgue wavelet coefficient is $\beta_{jk}=\int_0^1 f(x)\psi_{jk}(x)\,dx$.

Definition: The unweighted empirical coefficient and the density-corrected empirical coefficient are, respectively, $$ \widehat{\beta}^{\,U}_{jk} =\frac{1}{n}\sum_{i=1}^n Y_i\psi_{jk}(X_i), \qquad \widehat{\beta}^{\,C}_{jk} =\frac{1}{n}\sum_{i=1}^n \frac{Y_i\psi_{jk}(X_i)}{g(X_i)}. $$ The corrected estimator assumes that $g$ is known and positive wherever the coefficient is being estimated. The superscripts $U$ and $C$ distinguish the unweighted and corrected estimators.

The unweighted estimator averages under the distribution of $X$, so its population target is $\gamma_{jk}=\int_0^1 f(x)\psi_{jk}(x)g(x)\,dx$. The corrected estimator uses the factor $1/g(X_i)$ to recover the Lebesgue integral defining $\beta_{jk}$.

Worked Example: A Constant Signal Can Have a Nonzero Unweighted Coefficient

Use the Haar wavelet $\psi(x)=1$ on $[0,1/2)$ and $\psi(x)=-1$ on $[1/2,1]$. It has squared integral $1$. Let $f(x)=1$, and let the design density be $g(x)=3/2$ on the first half and $g(x)=1/2$ on the second half. This density integrates to $(3/2)(1/2)+(1/2)(1/2)=1$, and is bounded above and below by positive constants.

The Lebesgue coefficient is $$ \beta=\int_0^{1/2}1\,dx-\int_{1/2}^1 1\,dx =\frac12-\frac12=0. $$ In contrast, the unweighted estimator has expectation $$ E[Y\psi(X)] =\int_0^{1/2}\frac32\,dx-\int_{1/2}^1\frac12\,dx =\frac34-\frac14=\frac12. $$ The density-corrected expectation is $$ E\!\left[\frac{Y\psi(X)}{g(X)}\right] =\int_0^1 f(x)\psi(x)\,dx=0. $$ Thus the nonzero expectation of the unweighted average reflects how often each half is sampled, not a nonzero Lebesgue wavelet coefficient.

Expectations and Exact Variances

Write $q(x)=E[Y^2\mid X=x]=f(x)^2+\sigma^2(x)$. The following result gives the targets and variances directly. In addition to the model assumptions above, suppose the integrals displayed below are finite; for example, bounded $q$ and a density bounded away from zero are sufficient for the corrected variance.

Proposition: For a fixed coefficient $(j,k)$, $$ E[\widehat{\beta}^{\,U}_{jk}]=\gamma_{jk}, \qquad \operatorname{Var}(\widehat{\beta}^{\,U}_{jk}) =\frac{1}{n}\left\{ \int_0^1 q(x)\psi_{jk}(x)^2g(x)\,dx-\gamma_{jk}^2 \right\}, $$ and $$ E[\widehat{\beta}^{\,C}_{jk}]=\beta_{jk}, \qquad \operatorname{Var}(\widehat{\beta}^{\,C}_{jk}) =\frac{1}{n}\left\{ \int_0^1\frac{q(x)\psi_{jk}(x)^2}{g(x)}\,dx-\beta_{jk}^2 \right\}. $$

Proof. For one observation, conditional expectation gives $E[Y\mid X=x]=f(x)$. Therefore $$ E[Y\psi_{jk}(X)] =\int_0^1 f(x)\psi_{jk}(x)g(x)\,dx =\gamma_{jk}, $$ whereas $$ E\!\left[\frac{Y\psi_{jk}(X)}{g(X)}\right] =\int_0^1\frac{f(x)\psi_{jk}(x)}{g(x)}g(x)\,dx =\beta_{jk}. $$ For the unweighted summand, its second moment is $$ E[Y^2\psi_{jk}(X)^2] =\int_0^1q(x)\psi_{jk}(x)^2g(x)\,dx. $$ For the corrected summand, it is $$ E\!\left[\frac{Y^2\psi_{jk}(X)^2}{g(X)^2}\right] =\int_0^1\frac{q(x)\psi_{jk}(x)^2}{g(x)}\,dx. $$ In each case, variance equals second moment minus squared expectation. Since the $n$ summands are independent and identically distributed, the variance of their average is the variance of one summand divided by $n$. Substituting the expectations and second moments proves both formulas. $\square$

The proposition separates two issues. The unweighted estimator can have bias for $\beta_{jk}$, equal to $\gamma_{jk}-\beta_{jk}$, but it is unbiased for $\gamma_{jk}$. The corrected estimator is unbiased for the desired Lebesgue coefficient, but its second moment contains $1/g(x)$. Regions with low sampling density can therefore contribute more to its variance.

Worked Example: Exact Variances Under a Piecewise Design

Continue with the Haar wavelet and density from the first example, but now take $f(x)=x$ and $Y=X+\varepsilon$, where $E[\varepsilon\mid X]=0$ and $\operatorname{Var}(\varepsilon\mid X)=1$. Then $q(x)=x^2+1$.

First calculate the target: $$ \beta=\int_0^{1/2}x\,dx-\int_{1/2}^1x\,dx =\frac18-\frac38=-\frac14. $$ The unweighted target is $$ \gamma=\frac32\int_0^{1/2}x\,dx-\frac12\int_{1/2}^1x\,dx =\frac{3}{16}-\frac{3}{16}=0. $$ So the unweighted estimator is biased for $\beta$, while the corrected one has expectation $-1/4$.

For the unweighted estimator, use $\int_0^{1/2}x^2\,dx=1/24$ and $\int_{1/2}^1x^2\,dx=7/24$. Its one-observation second moment is $$ \int_0^1(x^2+1)g(x)\,dx =\frac32\left(\frac{1}{24}+\frac12\right) +\frac12\left(\frac{7}{24}+\frac12\right) =\frac{29}{24}. $$ Since $\gamma=0$, the variance of the average is $29/(24n)$.

For the corrected estimator, the one-observation second moment is $$ \int_0^1\frac{x^2+1}{g(x)}\,dx =\frac23\left(\frac{1}{24}+\frac12\right) +2\left(\frac{7}{24}+\frac12\right) =\frac{13}{36}+\frac{19}{12} =\frac{35}{18}. $$ Subtracting $\beta^2=1/16$ gives the one-observation variance $35/18-1/16=271/144$. Thus $\operatorname{Var}(\widehat{\beta}^{\,C})=271/(144n)$. The calculation illustrates that correction changes both the target and the variance; neither estimator is automatically preferable without specifying the target and considering the resulting uncertainty.

Variance Bounds and Wavelet Resolution

A useful general bound follows from the density assumptions and wavelet normalization. Assume that, for constants $0<g_{\min}\leq g_{\max}<\infty$ and $M<\infty$, $g_{\min}\leq g(x)\leq g_{\max}$ and $q(x)\leq M$ almost everywhere. These assumptions ensure finite variances for both estimators.

Result: If $\int_0^1\psi_{jk}(x)^2\,dx=1$, then $$ \operatorname{Var}(\widehat{\beta}^{\,U}_{jk}) \leq \frac{M g_{\max}}{n}, \qquad \operatorname{Var}(\widehat{\beta}^{\,C}_{jk}) \leq \frac{M}{n g_{\min}}. $$

Proof. In each exact variance formula, the subtracted squared expectation is nonnegative. For the unweighted estimator, $$ \operatorname{Var}(\widehat{\beta}^{\,U}_{jk}) \leq \frac{1}{n}\int_0^1q(x)\psi_{jk}(x)^2g(x)\,dx \leq \frac{M g_{\max}}{n}\int_0^1\psi_{jk}(x)^2\,dx =\frac{M g_{\max}}{n}. $$ For the corrected estimator, the lower density bound implies $1/g(x)\leq1/g_{\min}$, so $$ \operatorname{Var}(\widehat{\beta}^{\,C}_{jk}) \leq \frac{1}{n}\int_0^1\frac{q(x)\psi_{jk}(x)^2}{g(x)}\,dx \leq \frac{M}{n g_{\min}}\int_0^1\psi_{jk}(x)^2\,dx =\frac{M}{n g_{\min}}. $$ This proves both bounds. $\square$

Notably, these bounds do not increase with $j$ for an individual $L^2$-normalized coefficient. A wavelet at a fine resolution has a narrower support and a larger amplitude, but the integral of its square remains one; the variance calculation uses that integral. The number of coefficients at a resolution does grow with resolution, so the total variability across many coefficients can grow even though each coefficient separately has an $O(1/n)$ bound. That is a different calculation from the variance of one coefficient.

Worked Example: Heteroskedastic Noise and the Cost of Correction

Use the same Haar wavelet and piecewise design density, take $f(x)=0$, and let the conditional error variance be $1$ on the first half and $4$ on the second half. Both coefficient estimators have expectation zero because the regression function is zero.

For the unweighted estimator, the one-observation second moment is $$ \int_0^1\sigma^2(x)\psi(x)^2g(x)\,dx =\frac12(1)\frac32+\frac12(4)\frac12 =\frac34+1=\frac74. $$ Its variance is therefore $7/(4n)$.

For the corrected estimator, the one-observation second moment is $$ \int_0^1\frac{\sigma^2(x)\psi(x)^2}{g(x)}\,dx =\frac12\frac{1}{3/2}+\frac12\frac{4}{1/2} =\frac13+4=\frac{13}{3}. $$ Its variance is $13/(3n)$. Here inverse-density weighting is required to target Lebesgue coefficients in general, but it gives greater variance for this coefficient. The low-density half has both a larger noise variance and a larger inverse-density weight.

When Each Estimator Is Appropriate

If the target is the standard wavelet coefficient $\beta_{jk}$, density correction gives an unbiased estimator when $g$ is known and bounded away from zero. Without correction, the empirical average targets $\gamma_{jk}$, a different integral. If the design is uniform on $[0,1]$, then $g=1$, and the two estimators coincide. More generally, they can agree for a particular coefficient if $\int f(x)\psi_{jk}(x)(g(x)-1)\,dx=0$, even when the design is nonuniform.

The lower density bound is substantive. If $g$ is very small on a region where the wavelet has nonzero values, division by $g(X_i)$ can make the corrected estimator highly variable. If $g$ vanishes there, the displayed correction is not defined, and that region cannot be recovered by this inverse-density estimator from the observed design alone. The upper and lower bounds control this issue; they are not merely algebraic conveniences.

Common Mistake

Pitfall: “The wavelets are orthonormal, so $n^{-1}\sum_iY_i\psi_{jk}(X_i)$ must estimate $\int f\psi_{jk}$.” Orthonormality is with respect to Lebesgue measure; the expectation of the empirical average is $\int f\psi_{jk}g$, because $X_i$ is sampled with density $g$. To estimate the Lebesgue coefficient under nonuniform random design, use the inverse-density factor and check that $g$ is positive and that the resulting variance is finite.

Check Your Understanding

Work through these questions using the definitions and results above.

What is the expectation of the unweighted empirical coefficient, and under what design density does it equal the Lebesgue coefficient for every integrable regression function?
Why does the density-corrected variance involve $1/g(x)$, rather than $g(x)$?
For $q(x)\leq M$, $g(x)\geq g_{\min}>0$, and an $L^2$-normalized wavelet, state the variance bound for the corrected coefficient.
Explain why an individual wavelet coefficient can have variance $O(1/n)$ even at a fine resolution.
What can happen to the corrected estimator’s variance when the design density is very small on part of the wavelet’s support?