A generalization of inverse probability weighting
Section 2. The generalized inverse probability estimator

• Previous
• Next

Estimators in this paper utilize a positive definite matrix $\Sigma \in$ ${ℝ}^{N\times N}$. A matrix formulation of the estimators will therefore be useful. For a vector of interest $y\prime =\left(y_1,y_2,\dots ,y_N\right)$  and $1_{N\times 1}$  a vector of ones, the inverse probability estimator of the total $\theta ={\displaystyle {\sum }_{i=1}^N{y}_i}=y\prime 1_{N\times 1}$  can be written

$$\begin{array}{ll}{\widehat{\theta }}_{\text{IP}}\hfill & ={\displaystyle \sum _{i=1}^N\frac{{\delta }_i{y}_i}{{\pi }_i}}\hfill \\ \hfill & =y\prime {\Delta }_s{\left(E\left({\Delta }_s\right)\right)}^{-1}{1}_{N\times 1},(2.1)\hfill \end{array}$$

where ${\pi }_i=E\left({\delta }_i\right)$ is assumed greater than 0 for $i=1,2,\dots ,N,$ and ${\Delta }_s$ is the $N\times N$ diagonal matrix of the ${\delta }_i.$

The generalization of the inverse probability estimator relies on the Moore-Penrose inverse of a matrix $M,$ denoted $M^{†}.$ The Moore-Penrose inverse is unique and always exists; it is equal to the ordinary inverse if the latter exists. A precise definition and properties of the Moore-Penrose inverse can be found in Ben-Israel and Greville (2002). In particular, it can be verified that ${\Delta }_s^{†}={\Delta }_s.$ Since it is also true that ${\Delta }_s^2={\Delta }_s,$ if $I\in$ ${ℝ}^{N\times N}$ is the identity matrix, the inverse probability estimator can be written

$$\begin{array}{ll}{\widehat{\theta }}_{\text{IP}}\hfill & =y\prime {\Delta }_s{\left(E\left({\Delta }_s\right)\right)}^{-1}{1}_{N\times 1}\hfill \\ \hfill & =y\prime {\left({\Delta }_{s}I{\Delta }_s\right)}^{†}{\left(E{\left({\Delta }_{s}I{\Delta }_s\right)}^{†}\right)}^{-1}{1}_{N\times 1}.(2.2)\hfill \end{array}$$

If in (2.2), the identity matrix is replaced by any $N\times N$ positive definite matrix $\Sigma ,$ one obtains the generalized inverse probability estimator or the generalized Horvitz-Thompson estimator,

$${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right)=y\prime {\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}{\left(E{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}\right)}^{-1}{1}_{N\times 1}.(2.3)$$

In the phrase “inverse probability”, the matrix $E{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}$ is now the “probability” and ${\left(E{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}\right)}^{-1}$ is the new “inverse probability”. The ordinary inverse probability estimator is simply a special case of ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right),$ which can be obtained by choosing $\Sigma =I.$ As will be seen in c) below, one now has a family of unbiased estimators, ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right),$ parameterized by $\Sigma .$

2.1   Notes on the generalized inverse probability estimator

1. Although the vector $y$ appears in the estimator, only the sampled units affect the estimator’s value. This is because ${\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}={\Delta }_s{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†},$ thus (2.3) could have been written

$${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right)={\left({\Delta }_{s}y\right)}^{\prime }{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}{\left(E{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}\right)}^{-1}{1}_{N\times 1}.(2.4)$$

The proof of this and of many other results stated here may be found in Théberge (2017). The $N\times 1$ vector $w_{s\text{GIP}}\left(\Sigma \right)={\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}{\left(E{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}\right)}^{-1}{1}_{N\times 1}$ gives the weights of ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right),$ and all the units not in sample have a weight of zero.

2. The matrix $E{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}$ is invertible under the assumptions that ${\pi }_i=E\left({\delta }_i\right)$ is greater than zero for $i=1,2,\dots ,N$ and that $\Sigma$ is positive definite. Thus, (2.3) is well defined.
3. By taking the expectation of (2.3), one immediately sees that ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right)$ is unbiased for estimating $\theta =y\prime 1_{N\times 1}.$ This is true for any positive definite $\Sigma .$ A poor choice of $\Sigma$ may mean an estimator with a high variance, but it does not cause a bias.
4. Often, there is no closed-form formula for $E{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†},$ but for single stage sampling plans at least, it can be easily approximated. One simply takes the average of a large number of values of ${\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†},$ each computed for a different sample obtained with the sampling plan. The computation does not require the knowledge of any of the variables of interest. It is a “desk exercise” in the sense that it does not require contacting the units. It can even be carried out before the actual sample is selected.
5. It is well known that for a total estimator utilizing a regression vector $\beta ,$ $\widehat{T}\left(\beta \right),$ is asymptotically equivalent in terms of bias and variance to the estimator $\widehat{T}\left({\widehat{\beta }}_s\right)$ where ${\widehat{\beta }}_s$ is an estimator that converges in probability to $\beta .$ Similarly, ${\widehat{\theta }}_{\text{GIP}}\left({\widehat{\Sigma }}_s\right)$ has the same asymptotic bias and variance as ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right)$ if the positive definite matrix ${\widehat{\Sigma }}_s$ converges in probability to the positive definite matrix $\Sigma .$ In essence, if the sample size is sufficiently large, the error introduced by estimating $\Sigma$ by ${\widehat{\Sigma }}_s$ is negligible compared to the error in ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right)$ due to the sampling of units. All asymptotic results in this paper assume that the sampling plan is non informative (see, for example Cassel, Särndal and Wretman, 1977).
6. When $\Sigma =I,$ then ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right)$ reduces to the ordinary inverse probability estimator, ${\widehat{\theta }}_{\text{IP}},$ as given in (2.1). This is the justification for referring to ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right)$ as the generalized inverse probability estimator or the generalized Horvitz-Thompson estimator. It will be seen later, why this particular unbiased extension of the ordinary Horvitz-Thompson estimator is of interest.
7. An arbitrary symmetric positive definite matrix $\Sigma$ may contain up to $N\left(N+1\right)/2$ distinct parameters. It is not feasible to specify so many values. If the sample $s$ is utilized to estimate those parameters, the task of estimating $N\left(N+1\right)/2$ parameters from $n<N$ observations is clearly impossible. A simpler choice must be used. The simplest choice utilizes $\Sigma =I,$ as seen in f). There are other choices that have a reasonable number of parameters. One example is given in Section 6.
8. For estimating a domain total $y\prime c$ where $c={\left(c_1,\dots ,c_i,\dots ,c_N\right)}^{\prime }$ is a vector of known constants with $c_i=1$ or 0 depending on whether unit $i$ is in the domain or not, it suffices to replace (2.3), which is for estimating the population total, with $y\prime {\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}{\left(E{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}\right)}^{-1}c.$ The weight vector ${\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}{\left(E{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}\right)}^{-1}c$ varies with each domain described by $c;$ however the weight matrix, ${\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}{\left(E{\left({\Delta }_s\Sigma {\Delta }_s\right)}^{†}\right)}^{-1},$ does not depend on the domain. There are $N-n$ rows of this matrix that are nil. Even though there are potentially $nN$ elements of the weight matrix that are non zero, post-multiplication by $c$ will give the weight vector for any domain described by $c.$
9. One possibility for the matrix $\Sigma$ is one where all the diagonal elements are the same, and all the off-diagonal elements are the same. In this way, all the units are the same with respect to $\Sigma .$ However, if all units are the same with respect to the sampling plan, for example simple random sampling or Bernoulli sampling, and if all units are the same with respect to the parameter estimated, for example a total or an average for all units, then by symmetry, every sampled unit will have the same weight. Since both ${\widehat{\theta }}_{\text{IP}}$ and ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right)$ are unbiased, both estimators will have the same weights. Nonetheless, for domain parameters, because some units are in the domain and some not, the symmetry argument no longer holds and the value of the off-diagonal elements of $\Sigma$ may make a difference in ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right).$
10. By setting $y=1_{N\times 1}$ in ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right),$ the estimator simply becomes the sum of all the weights of the sampled units and the parameter estimated becomes $1_{1\times N}{1}_{N\times 1}=N,$ the known total number of units. However, the sum of the weights does not necessarily equal $N.$ This does not bode well for the variance of ${\widehat{\theta }}_{\text{GIP}}\left(\Sigma \right).$ To fix this, calibration can be used. Calibration was introduced by Deville and Särndal (1992). At its simplest, it would consist of scaling the inverse probability weights, generalized or not, by a common factor so that the resulting final weights do add up to $N.$ Even for the ordinary inverse probability estimator, for some sampling plans, the sum of the design weights does not necessarily equal $N,$ and here too, the solution lies in calibration. The subject of calibration is examined in the next section.

• Previous
• Article main page
• Next