7. Conclusions

Piero Demetrio Falorsi and Paolo Righi

Previous

The paper proposes a new approach for defining the optimal inclusion probabilities in various survey contexts, which are characterized by the need to disseminate survey estimates of prefixed accuracy, for a multiplicity of both variables and domains of interest.

This paper’s main contribution is the practical computation of these probabilities by means of a new algorithm, which is suitable for a general multi-way sampling design in which the standard stratified sampling represents a special case. The proposed approach, the algorithm and the final computation are domain- and variable-driven.

In our framework, the domain membership indicator variables are assumed to be known, while the variables of interest are not known. The procedure is, then, applied on the predicted values of the characteristics of interest via a superpopulation model, and the algorithm enables taking into account model uncertainty; this reflects the non-knowledge of the values of variables of interest. Using the Anticipated Variance as the measure of the estimators’ precision, this approach overcomes the limits of the standard algorithms for the sample allocation, in which the variables of interest driving the solution are assumed to be known.

The proposed algorithm exploits standard procedure, but does present some computational innovations which may be useful for dealing with the complexity deriving from the fact that the Anticipated Variances are implicit functions of the inclusion probabilities. The algorithm was tested on simulated and real survey data, to evaluate its performance and properties. The results of a small set of experiments are presented here. They confirm an improvement, in terms of efficiency, of the sampling strategy. A natural generalization of the case examined here may be developed by considering, as known during the design planning stage, the indicators of the domains and other quantitative independent variables. We note that the Anticipated Variance considering only the domain indicators is larger than the Anticipated Variance of this more general case. Thus, our solution represents an upper (and somehow robust) boundary solution in the design phase. Furthermore, the algorithmic solution can be easily adapted to this more general situation.

Acknowledgements

This research was funded by the partnership of the Global Strategy to improve Agricultural and Rural Statistics: http://www.fao.org/economic/ess/ess-capacity/ess-strategy/en/.

Appendix

Appendix A1

AV of the HT estimator

Let us consider the residual ${\eta }_{\left(dr\right)k}$ as expressed by equation (3.5), and replace the term $y_{rk}$ with ${\tilde{y}}_{rk}+u_{rk},$ thus obtaining

$${\eta }_{\left(dr\right)k}=\left({\tilde{y}}_{rk}+u_{rk}\right){\gamma }_{dk}-{\pi }_k{{\delta }^{\prime }}_k{\left[A\left(\pi \right)\right]}^{-1}{\displaystyle {\sum }_{j\in U}{\pi }_j{\delta }_j\left({\tilde{y}}_{rj}+u_{rj}\right){\gamma }_{dj}\left(1/{\pi }_j-1\right)}.(\text{A}1.1)$$

The weighted least predictions of ${\tilde{y}}_{rk}{\gamma }_{dk}$ and $u_{rk}{\gamma }_{dk},$ with predictors ${\pi }_k{\delta }_k$ and weights $1/{\pi }_k-1,$ are

$${\widehat{\tilde{y}}}_{\left(dr\right)k}={\pi }_k{a}_{\left(dr\right)k}(\text{A}1.2)$$

and

$${\widehat{u}}_{\left(dr\right)k}={\pi }_k{{\delta }^{\prime }}_k{\left[A\left(\pi \right)\right]}^{-1}{\displaystyle {\sum }_{j\in U}{\pi }_j{\delta }_j{u}_{rj}{\gamma }_{dj}}\left(1/{\pi }_j-1\right),(\text{A}1.3)$$

with

$$a_{\left(dr\right)k}\left(\pi \right)={{\delta }^{\prime }}_k{\left[A\left(\pi \right)\right]}^{-1}{\displaystyle {\sum }_{j\in U}{\pi }_j{\delta }_j{\tilde{y}}_{rj}{\gamma }_{dk}}\left(1/{\pi }_j-1\right).(\text{A}1.4)$$

Using the formulae (A1.2) and (A1.3), the expression (A1.1) may be reformulated as ${\eta }_{\left(dr\right)k}=\left({\tilde{y}}_{rk}+u_{rk}\right){\gamma }_{dk}-\left[{\widehat{\tilde{y}}}_{\left(dr\right)k}+{\widehat{u}}_{\left(dr\right)k}\right].$ Therefore, the model expectation of ${\eta }_{\left(dr\right)k}^2$ is

$$E_M\left({\eta }_{\left(dr\right)k}^2\right)={\left({\tilde{y}}_{rk}{\gamma }_{dk}-{\widehat{\tilde{y}}}_{\left(dr\right)k}\right)}^2+E_M\left[{\left(u_{rk}{\gamma }_{dk}-{\widehat{u}}_{\left(dr\right)k}\right)}^2\right]+\text{Mean zero terms,}(\text{A}1.5)$$

because $E_M\left(u_{rk}\right)=0.$ Furthermore,

$$E_M\left[{\left(u_{rk}{\gamma }_{dk}-{\widehat{u}}_{\left(dr\right)k}\right)}^2\right]={\sigma }_{rk}^2{\gamma }_{dk}+E_M{\left({\widehat{u}}_{\left(dr\right)k}\right)}^2-2E_M\left(u_{rk}{\gamma }_{dk},{\widehat{u}}_{\left(dr\right)k}\right),(\text{A}1.6)$$

where $E_M\left(u_{rk}{\gamma }_{dk}{\widehat{u}}_{\left(dr\right)k}\right)={\pi }_k{b}_{\left(dr\right)k}\left(\pi \right)$ and $E_M{\left({\widehat{u}}_{\left(dr\right)k}\right)}^2={\pi }_k^2{c}_{\left(dr\right)k}\left(\pi \right),$ with

$$b_{\left(dr\right)k}\left(\pi \right)={{\delta }^{\prime }}_k{\left[A\left(\pi \right)\right]}^{-1}{\delta }_k{\sigma }_{rk}^2{\gamma }_{dk}\left(1-{\pi }_k\right)(\text{A}1.7)$$

and

$$c_{\left(dr\right)k}\left(\pi \right)={{\delta }^{\prime }}_k{\left[A\left(\pi \right)\right]}^{-1}\left[{\displaystyle {\sum }_{j\in U}{\delta }_j{{\delta }^{\prime }}_j{\sigma }_{rj}^2{\gamma }_{dj}{\left(1-{\pi }_j\right)}^2}\right]{\left[A\left(\pi \right)\right]}^{-1}{\delta }_k.(\text{A}1.8)$$

Expression (4.5) is easily derived by plugging expressions from (A1.2) to (A1.8) into equation (4.3).

Appendix A2

Convergence of the algorithm

The optimization problem (5.1) is solved by two nested fixed point iterations. Given an unknown vector $x$ of dimension $q,$ the fixed point iteration chooses an initial guess ${}^{\left(0\right)}x.$ Then, it computes subsequent iterates by ${}^{\left(\tau +1\right)}x=g\left({}^{\left(\tau \right)}x\right),$ with $\tau =1,2,\dots ,$ with $g(\cdot )$ being a system of $q$ updating equations. The multivariate function $g$ has a fixed point in a domain $Q\subseteq {\Re }^q$ if $g$ maps $Q$ in $Q.$ Let $J_g\left(x\right)$ be the Jacobian matrix of first partial derivate of $g$ evaluated at $x,$ if there exists a constant $\rho <1$ such that, in some natural matrix norm, $\Vert J_g\left(x\right)\Vert \le \rho ,x\in Q,$ $g$ has a unique fixed point $x^{\ast }\in Q,$  and the fixed point iteration is guaranteed to converge to $x^{\ast }$  for any initial guess chosen in $Q.$ As regards the proposed algorithm, the convergence of the IL and OL is obtained when the terms ${}^{\left(\alpha \tau \right)}\text{A}{\text{AV}}_{3\left(dr\right)}$ converge to the fixed point. This means that the vectors ${}^{\left(\alpha \right)}\pi$ and ${}^{\left(\alpha \tau \right)}\pi$  do not change in the OL and IL iterations. The demonstration below considers the method proposed by Chromy (1987) to solve the LCSP of system (5.7), and makes use of some reasonable assumptions: (1)  ${\widehat{u}}_{\left(dr\right)k}\cong 0;$ (2) $\left[N/\left(N-H\right)\right]\cong 1;$ (3) ${\widehat{\tilde{y}}}_{rk}\cong {\tilde{y}}_{rk};$ (4) ${}^{\left(\alpha \right)}\pi {}_k\cong {}^{\left(\alpha \tau \right)}\Delta {}^{\left(\alpha \tau \right)}\pi {}_k$ with $0<{}^{\left(\alpha \tau \right)}\Delta \le 1;$ (5) $c_k\cong \overline{c}.$ Assumption (1) corresponds to the upward approximation of the Anticipated Variance, given in Remark 4.1, and implies that $b_{\left(dr\right)k}\left({}^{\left(\alpha \right)}\pi \right)=c_{\left(dr\right)k}\left({}^{\left(\alpha \right)}\pi \right)=0.$ Assumption (3) implies that $a_{\left(dr\right)k}\left({}^{\left(\alpha \right)}\pi \right){\tilde{y}}_{rk}{\gamma }_{dk}\cong {\tilde{y}}_{rk}^2{\gamma }_{dk}/{}^{\left(\alpha \right)}\pi {}_k.$ Assumption (4) states that the structure of the inclusion probabilities remains roughly constant in the different IL iterations. The assumption becomes reasonable considering that the updating equation A2.2 below (of a given inclusion probability) is essentially determined by the variance threshold that requires the largest sample size. It is plausible to hypothesize that this threshold remains more or less the same in the subsequent IL iterations of a given OL.

Proof of convergence of the Inner Loop. By reformulating expression (4.6) in accordance with the assumptions from (1) to (4),

$${}^{\left(\alpha \tau +1\right)}\text{A}{\text{AV}}_{3\left(dr\right)}={\displaystyle {\sum }_{k\in U}\left[\left(\frac{1}{{}^{\left(\alpha \tau +1\right)}\pi {}_k}-1\right)\left(2\frac{{\tilde{y}}_{rk}^2{\gamma }_{dk}}{{}^{\left(\alpha \tau +1\right)}\Delta }-\frac{{\tilde{y}}_{rk}^2{\gamma }_{dk}}{{}^{\left(\alpha \tau +1\right)}\Delta {}^2}\right)\right]}.(\text{A}2.1)$$

Considering in problem (5.7) that the ${}^{\left(a\tau \right)}\text{A}{\text{AV}}_{3\left(dr\right)}$ values are fixed, each value of the vector ${}^{\left(\alpha \tau +1\right)}\pi$ is obtained as a solution of the LCSP with the Chromy algorithm. Denote with $\alpha \tau v*$ the iteration of the Chromy algorithm into which it converges, where ${}^{\left(\alpha \tau v*+1\right)}\pi \cong {}^{\left(\alpha \tau v*\right)}\pi .$ Then, the IL updates the generic probability in accordance with the expression

$${}^{\left(\alpha \tau +1\right)}\pi {}_k={\left[{\displaystyle {\sum }_{\left(dr\right)}{}^{\left(\alpha \tau v*+1\right)}\varphi {}_{\left(dr\right)}}\frac{\left({\tilde{y}}_{rk}^2+{\sigma }_{rk}^2\right){\gamma }_{dk}}{\overline{c}}\right]}^{1/2},(\text{A}2.2)$$

where the right-hand term represents the updating formula of the Chromy algorithm, and ${\displaystyle \sum {}_{\left(dr\right)}}$ stands for ${\displaystyle {\sum }_{d=1}^D{\displaystyle \sum {}_{r=1}^R}},$ and ${}^{\left(\alpha \tau v*+1\right)}\varphi {}_{\left(dr\right)}$ is the generalized Lagrange multiplier, where

$$\begin{array}{lll}{}^{\left(\alpha \tau v*+1\right)}\varphi {}_{\left(dr\right)}\hfill & =\hfill & {}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(dr\right)}{\left[\frac{{}^{\left(\alpha \tau v*\right)}V{}_{\left(dr\right)}}{{\stackrel{⃛}{V}}_{\left(dr\right)}+{}^{\left(\alpha \tau \right)}\text{A}{\text{AV}}_{3\left(dr\right)}}\right]}^2,\hfill \\ {}^{\left(\alpha \tau v*\right)}V{}_{\left(dr\right)}\hfill & =\hfill & {\displaystyle {\sum }_{k\in U}\frac{\left({\tilde{y}}_{rk}^2+{\sigma }_{rk}^2\right){\gamma }_{dk}}{{}^{\left(\alpha \tau v*\right)}\pi {}_k}}\hfill \end{array}(\text{A}2.3)$$

and

$${\stackrel{⃛}{V}}_{\left(dr\right)}={\overline{V}}_{\left(dr\right)}+{\displaystyle {\sum }_{k\in U}\left({\tilde{y}}_{rk}^2+{\sigma }_{rk}^2\right){\gamma }_{dk}}.$$

The Kuhn-Tucker theory states that ${}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(dr\right)}\left[{}^{\left(\alpha \tau v*\right)}V{}_{\left(dr\right)}-\left({\stackrel{⃛}{V}}_{\left(dr\right)}+{}^{\left(\alpha \tau \right)}\text{A}{\text{V}}_{3\left(dr\right)}\right)\right]=0;$ therefore, ${}^{\left(\alpha \tau v*+1\right)}\varphi {}_{\left(dr\right)}={}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(dr\right)}$ and ${}^{\left(\alpha \tau v*+1\right)}\varphi {}_{\left(dr\right)}>0$  iff ${}^{\left(\alpha \tau v*\right)}V{}_{\left(dr\right)}/\left({\stackrel{⃛}{V}}_{\left(dr\right)}+{}^{\left(\alpha \tau \right)}\text{A}{\text{V}}_{3\left(dr\right)}\right)=1.$ Chromy asserts that few ${}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(dr\right)}$ $\left(\text{for }r=1,\dots ,R;d=1,\dots ,D\right)$ are larger than zero, and that in most cases, only one value is strictly positive. Denoting with ${}^{\left(\alpha \tau \right)}AA{V}_3={\left({}^{\left(\alpha \tau \right)}\text{A}{\text{AV}}_{3\left(11\right)},\dots ,{}^{\left(\alpha \tau \right)}\text{A}{\text{AV}}_{3\left(1R\right)},\dots ,{}^{\left(\alpha \tau \right)}\text{A}{\text{AV}}_{3\left(DR\right)}\right)}^{\prime },$ we define ${}^{\left(\alpha \tau +1\right)}AA{V}_3=g\left({}^{\left(\alpha \tau \right)}AA{V}_{\text{3}}\right)$ as the system of $D\times R$ updating equations where the generic $\left(\overline{dr}\right)$ equation of the system

$$\begin{array}{lll}{g}_{\left(\overline{dr}\right)}\left({}^{\left(\alpha \tau \right)}AA{V}_3\right)\hfill & \cong \hfill & {\displaystyle {\sum }_{k\in U}\left(2\frac{{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}{{}^{\left(\alpha \tau +1\right)}\Delta }-\frac{{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}{{}^{\left(\alpha \tau +1\right)}\Delta {}^2}\right)}\hfill \\ \hfill & \times \hfill & \left\{{\left[{\displaystyle {\sum }_{\left(dr\right)}{}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(dr\right)}}{\left[\frac{{}^{\left(\alpha \tau v*\right)}V{}_{\left(dr\right)}}{{\stackrel{⃛}{V}}_{\left(dr\right)}+{}^{\left(\alpha \tau \right)}\text{A}{\text{AV}}_{3\left(dr\right)}}\right]}^2\frac{\left({\tilde{y}}_{rk}^2+{\sigma }_{rk}^2\right){\gamma }_{dk}}{\overline{c}}\right]}^{-1/2}-1\right\},(\text{A}2.4)\hfill \end{array}$$

is obtained by plugging expression (A2.2) into (A2.1). If the convergence is obtained, then in the last iteration, ${}^{\left(\alpha \tau +1\right)}AA{V}_3\cong {}^{\left(\alpha \tau \right)}AA{V}_{\text{3}}.$  The function of equation (A2.4) is continuous and differentiable. Moreover, it maps onto the interval of the possible values of ${\text{AAV}}_{3\left(dr\right)}.$ Then, the IL converges if the following condition is fulfilled:

$$\Vert J_g\left(AA{V}_3\right)\Vert \le 1.(\text{A}2.5)$$

The Jacobian matrix is positive semi-defined, and a well-known result states that $\text{trace }\left(J_g{J^{\prime }}_g\right)\le \text{trace }{\left(J_g\right)}^2.$ By considering the Frobenius norm ${\Vert J_g\Vert }_F=\sqrt{\text{trace }\left(J_g{J^{\prime }}_g\right)},$ it is ${\Vert J_g\Vert }_F\le \text{trace }\left(J_g\right).$ Thus we can take into account the trace of the Jacobian matrix to verify condition (A2.5). Let ${g^{\prime }}_{\left(\overline{dr}\right)}=\partial g_{\left(\overline{dr}\right)}\left({}^{\left(\alpha \tau -1\right)}AA{V}_{3\left(dr\right)}/\partial {}^{\left(\alpha \tau -1\right)}\text{A}{\text{AV}}_{3\left(\overline{dr}\right)}\right)$ be the $\left(\overline{dr}\right)$ element of the diagonal of $J_g\left(AA{V}_3\right).$ Using the Kuhn-Tucker condition ${}^{\left(\alpha \tau v*\right)}V{}_{\left(dr\right)}/\left({\stackrel{⃛}{V}}_{\left(dr\right)}+{}^{\left(\alpha \tau \right)}\text{A}{\text{V}}_{3\left(dr\right)}\right)=1,$

$$\begin{array}{lll}{g^{\prime }}_{\left(\overline{dr}\right)}\hfill & =\hfill & {\displaystyle {\sum }_{k\in U}\left(2\frac{{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}{{}^{\left(\alpha \tau +1\right)}\Delta }-\frac{{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}{{}^{\left(\alpha \tau +1\right)}\Delta {}^2}\right){\left[{\displaystyle {\sum }_{\left(dr\right)}{}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(dr\right)}\frac{\left({\tilde{y}}_{rk}^2+{\sigma }_{rk}^2\right){\gamma }_{dk}}{\overline{c}}}\right]}^{-3/2}}\hfill \\ \hfill & \times \hfill & {}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}\frac{1}{{}^{\left(\alpha \tau v*\right)}V{}_{\left(\overline{dr}\right)}}\frac{\left({\tilde{y}}_{\overline{r}k}^2+{\sigma }_{\overline{r}k}^2\right){\gamma }_{\overline{d}k}}{\overline{c}}.\hfill \end{array}$$

Since many ${}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}=0$ (Chromy 1987), the respective ${g^{\prime }}_{\left(\overline{dr}\right)}$ is null. When ${}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}>0,$ then

$$\begin{array}{lll}{g^{\prime }}_{\left(\overline{dr}\right)}\hfill & \le \hfill & {\displaystyle {\sum }_{k\in U}\left(2\frac{{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}{{}^{\left(\alpha \tau +1\right)}\Delta }-\frac{{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}{{}^{\left(\alpha \tau +1\right)}\Delta {}^2}\right){\left[{}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}\frac{\left({\tilde{y}}_{\overline{r}k}^2+{\sigma }_{\overline{r}k}^2\right){\gamma }_{\overline{d}k}}{\overline{c}}\right]}^{-3/2}\times {}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}\frac{1}{{}^{\left(\alpha \tau v*\right)}V{}_{\left(\overline{dr}\right)}}\frac{\left({\tilde{y}}_{\overline{r}k}^2+{\sigma }_{\overline{r}k}^2\right){\gamma }_{\overline{d}k}}{\overline{c}}}\hfill \\ \hfill & =\hfill & {\displaystyle {\sum }_{k\in U}\left(2\frac{{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}{{}^{\left(\alpha \tau +1\right)}\Delta }-\frac{{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}{{}^{\left(\alpha \tau +1\right)}\Delta {}^2}\right)\frac{1}{\sqrt{{}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}\frac{\left({\tilde{y}}_{\overline{r}k}^2+{\sigma }_{\overline{r}k}^2\right){\gamma }_{\overline{d}k}}{\overline{c}}}{}^{\left(\alpha \tau v*\right)}V{}_{\left(\overline{dr}\right)}}}\hfill \\ \hfill & \le \hfill & {\displaystyle {\sum }_{k\in U}\frac{\frac{{\tilde{y}}_{\overline{r}k}{\gamma }_{\overline{d}k}}{{}^{\left(\alpha \tau +1\right)}\Delta }\left(2-\frac{1}{{}^{\left(\alpha \tau +1\right)}\Delta }\right)}{\sqrt{\overline{c}{}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}{\gamma }_{\overline{d}k}}{}^{\left(\alpha \tau v*\right)}V{}_{\left(\overline{dr}\right)}}}<<1.\hfill \end{array}$$

Therefore, the $\text{trace}\left(J_g\right)$ should be less than 1.

Proof of convergence of the Outer Loop. Let ${}^{\left(\alpha \tau +1\right)}\pi$ be the fixed point solution of the IL; then, the OL updates the vector ${}^{\left(\alpha \right)}\pi$ with ${}^{\left(\alpha +1\right)}\pi ={}^{\left(\alpha \tau +1\right)}\pi .$ Under conditions (1), (2) and (3),

$${}^{\left(\alpha +1\right)}\text{A}{\text{AV}}_{3\left(dr\right)}={\displaystyle {\sum }_{k\in U}\left(\frac{1}{{}^{\left(\alpha \tau +1\right)}\pi {}_k}-1\right){\tilde{y}}_{rk}^2{\gamma }_{dk}.}(\text{A}2.6)$$

Plugging expression (A2.2) into formula (A2.6) when the IL converges, the system of $D\times R$ updating equations of ${}^{\left(\alpha +1\right)}AA{V}_3$ is given by ${}^{\left(\alpha +1\right)}AA{V}_3=j\left({}^{\left(\alpha \tau \right)}AA{V}_{\text{3}}\right),$ where the generic equation of $j$ is

$$\begin{array}{lll}{}^{\left(\alpha +1\right)}\text{A}{\text{AV}}_{3\left(dr\right)}\hfill & =\hfill & j_{\left(\overline{dr}\right)}\left({}^{\left(\alpha \tau \right)}AA{V}_{\text{3}}\right)\hfill \\ \hfill & =\hfill & {\displaystyle {\sum }_{k\in U}{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}\left({\left[{\displaystyle {\sum }_{\left(dr\right)}{}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(dr\right)}}{\left[\frac{{}^{\left(\alpha \tau v*\right)}V{}_{\left(dr\right)}}{{\stackrel{⃛}{V}}_{\left(\overline{d}r\right)}+{}^{\left(\alpha \tau \right)}\text{A}{\text{AV}}_{3\left(\overline{d}r\right)}}\right]}^2\frac{\left({\tilde{y}}_{rk}^2+{\sigma }_{rk}^2\right){\gamma }_{dk}}{\overline{c}}\right]}^{-1/2}-1\right).(\text{A}2.7)\hfill \end{array}$$

Denoting with ${}^{\left(\alpha \right)}AA{V}_{\text{3}}={}^{\left(\alpha \tau =0\right)}AA{V}_{\text{3}},$ the system j may be expressed in a recursive form

$${}^{\left(\alpha +1\right)}AA{V}_3\cong j\left(g\left({}^{\left(\alpha \tau -1\right)}AA{V}_{\text{3}}\right)\right)=j\left(g\left(g\left(\mathrm{.....}g\left({}^{\left(\alpha \tau =0\right)}AA{V}_{\text{3}}\right)\right)\right)\right)=f\left({}^{\left(\alpha \right)}AA{V}_{\text{3}}\right),$$

with $f(\cdot )=j\left(g\left(g\left(\mathrm{.....}g(\cdot )\right)\right)\right)$ as the system of $D\times R$ updating equations of ${}^{\left(\alpha +1\right)}AA{V}_3,$ with respect to the previous values of the OL, ${}^{\left(\alpha \right)}AA{V}_3.$ To demonstrate the convergence of OL, it is necessary to demonstrate that the Jacobian norm $\Vert J_f\left(AA{V}_3\right)\Vert$ is lower than 1. Using standard results of matrix algebra,

$$\Vert J_f\left(AA{V}_3\right)\Vert \le \Vert J_j\left({}^{\left(\alpha \tau \right)}AA{V}_3\right)\Vert \times \Vert J_g\left({}^{\left(\alpha \tau -1\right)}AA{V}_3\right)\Vert \times \dots \times \Vert J_g\left({}^{\left(\alpha \tau =0\right)}AA{V}_3\right)\Vert ,$$

in which the generic norm $\Vert J_g(\cdot )\Vert$ is lesser than 1 (see the IL proof of convergence). Let ${j^{\prime }}_{\left(\overline{dr}\right)}$ be the $\left(\overline{dr}\right)$ element on the diagonal of $J_j\left({}^{\left(\alpha \tau \right)}AA{V}_3\right).$ It is

$$\begin{array}{lll}{j^{\prime }}_{\left(\overline{dr}\right)}\hfill & =\hfill & {\displaystyle {\sum }_{k\in U}{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}{\left[{\displaystyle {\sum }_{\left(dr\right)}{}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(dr\right)}}\frac{\left({\tilde{y}}_{rk}^2+{\sigma }_{rk}^2\right){\gamma }_{dk}}{\overline{c}}\right]}^{-3/2}\hfill \\ \hfill & \times \hfill & {}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}\frac{1}{{}^{\left(\alpha \tau v*\right)}V{}_{\left(\overline{dr}\right)}}\frac{\left({\tilde{y}}_{\overline{r}k}^2+{\sigma }_{\overline{r}k}^2\right){\gamma }_{\overline{d}k}}{\overline{c}}.(\text{A}2.8)\hfill \end{array}$$

Therefore, we have

$$\begin{array}{lll}{j^{\prime }}_{\left(\overline{dr}\right)}\hfill & \le \hfill & {\displaystyle {\sum }_{k\in U}{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}{\left[{}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}\frac{\left({\tilde{y}}_{\overline{r}k}^2+{\sigma }_{\overline{r}k}^2\right){\gamma }_{\overline{d}k}}{\overline{c}}\right]}^{-3/2}{}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}\frac{1}{{}^{\left(\alpha \tau v*\right)}V{}_{\left(\overline{dr}\right)}}\frac{\left({\tilde{y}}_{\overline{r}k}^2+{\sigma }_{\overline{r}k}^2\right){\gamma }_{\overline{d}k}}{\overline{c}}\hfill \\ \hfill & =\hfill & \frac{1}{{}^{\left(\alpha \tau v*\right)}V{}_{\left(\overline{dr}\right)}}{\displaystyle {\sum }_{k\in U}{\tilde{y}}_{\overline{r}k}^2{\gamma }_{\overline{d}k}}{\left[{}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}\frac{\left({\tilde{y}}_{\overline{r}k}^2+{\sigma }_{\overline{r}k}^2\right){\gamma }_{\overline{d}k}}{\overline{c}}\right]}^{-1/2}.\hfill \end{array}$$

The following inequality holds

$${j^{\prime }}_{\left(\overline{dr}\right)}<\frac{{\displaystyle {\sum }_{k\in U}{\tilde{y}}_{\overline{r}k}{\gamma }_{\overline{d}k}}}{\sqrt{\overline{c}{}^{\left(\alpha \tau v*\right)}\varphi {}_{\left(\overline{dr}\right)}}{}^{\left(\alpha \tau v*\right)}V{}_{\left(\overline{dr}\right)}}<<1.$$

Consequently, the norm $\Vert J_j\left({}^{\left(\alpha \tau \right)}AA{V}_3\right)\Vert <1,$  and therefore the OL converges.

Annexe A3

Proof that the approximation of Remark 4.1 is upward

Since ${\widehat{u}}_{\left(dr\right)k}$ is the weighted least square prediction of $u_{rk}{\gamma }_{dk},$ by using a different value of the ${\widehat{u}}_{\left(dr\right)k},$ such as ${\widehat{u}}_{\left(dr\right)k}=0,$ we obtain

$${\displaystyle {\sum }_{k\in U}\left(1/{\pi }_k-1\right)E_M\left[{\left(u_{rk}{\gamma }_{dk}-{\widehat{u}}_{\left(dr\right)k}\right)}^2\right]}\le {\displaystyle {\sum }_{k\in U}\left(1/{\pi }_k-1\right)E_M\left[{\left(u_{rk}{\gamma }_{dk}-0\right)}^2\right]},$$

where $E_M\left[{\left(u_{rk}{\gamma }_{dk}-0\right)}^2\right]={\sigma }_{rk}^2{\gamma }_{dk}.$ Replacing the terms $E_M\left[{\left(u_{rk}{\gamma }_{dk}-{\widehat{u}}_{\left(dr\right)k}\right)}^2\right]$ with ${\sigma }_{rk}^2{\gamma }_{dk}$ in expression (A1.5), the AAV (4.3) is inflated. The approximation ${\widehat{u}}_{\left(dr\right)k}=0$ implies that $b_{\left(dr\right)k}\left(\pi \right)=c_{\left(dr\right)k}\left(\pi \right)=0.$ Finally, we emphasize that in most cases, the upward is slight, since the ${\widehat{u}}_{\left(dr\right)k}$ are obtained by the $z_k$ variables that generally have a very low predictive power for the $u_{rk}{\gamma }_{dk}$  values (see Section 4). In these situations ${\widehat{u}}_{\left(dr\right)k}\cong \left(1/N\right){\displaystyle {\sum }_{k\in U}{u}_{rk}{\gamma }_{dk}\cong 0}.$ So $E_M\left(u_{rk}{\gamma }_{dk}{\widehat{u}}_{\left(dr\right)k}\right)\cong 0$ and $E_M{\left({\widehat{u}}_{\left(dr\right)k}\right)}^2\cong 0.$

Annexe A4

Proof of expression (4.7)

In this case, each ${\delta }_k$ vector has $H-1$ zero elements and 1 element equal to 1 (corresponding to the planned population to which the unit $k$ belongs). Given the input values, the optimization procedure ${\pi }_k={\pi }_h$ for $k\in U_h.$ Under the above assumption, ${\left[A\left(\pi \right)\right]}^{-1}$ is a diagonal matrix with the $h{h}^{\text{th}}$ element given by ${\left[A_{hh}\left(\pi \right)\right]}^{-1}={\left[N_h{\pi }_h^2\left(1/{\pi }_h-1\right)\right]}^{-1}.$ Considering ${\tilde{y}}_{rk}={\overline{Y}}_{rh},$  expressions (A1.2) and (A1.3) can be reformulated as, respectively,

$${\widehat{\tilde{y}}}_{\left(dr\right)k}={\pi }_h{{\delta }^{\prime }}_k{\left[A\left(\pi \right)\right]}^{-1}{N}_h{\pi }_h\left(1/{\pi }_h-1\right){\overline{Y}}_{rh}={\overline{Y}}_{rh}.(\text{A}4.1)$$

$${\widehat{u}}_{\left(dr\right)k}={\pi }_h{{\delta }^{\prime }}_k{\left[A\left(\pi \right)\right]}^{-1}{\pi }_h\left(1/{\pi }_h-1\right){\displaystyle {\sum }_{j\in U}{u}_{rj}}={\left({\pi }_h{N}_h\right)}^{-1}{\displaystyle {\sum }_{j\in U_h}{u}_{rj}},(\text{A}4.2)$$

but ${\displaystyle {\sum }_{j\in U_h}{u}_{rj}}=0$ as the sum of the residual of a regression model.

Using the formulae (A4.1) and (A4.2), expression (4.5) is given by

$$\begin{array}{lll}\text{AAV}\left({\widehat{t}}_{\left(dr\right)}\right)\hfill & =\hfill & \left[N/\left(N-H\right)\right]{\displaystyle {\sum }_h\left(\frac{1}{{\pi }_h}-1\right){\displaystyle {\sum }_{k\in U_h}{E}_M{\left(u_{rk}{\gamma }_{dk}\right)}^2}}\hfill \\ \hfill & =\hfill & \left[N/\left(N-H\right)\right]{\displaystyle {\sum }_{d=1}^D{\displaystyle {\sum }_{h\in H_d}{\sigma }_{rh}^2}{N}_h\left(N_h/n_h-1\right),}\hfill \end{array}$$

since ${\pi }_h=n_h/N_h$  and expression (4.7) may be obtained.

References

Bethel, J. (1989). Sample allocation in multivariate surveys. Survey Methodology, 15, 1, 47-57.

Boyd, S., and Vanderberg, L. (2004). Convex Optimization. Cambridge University Press.

Breidt, F.J., and Chauvet, G. (2011). Improved variance estimation for balanced samples drawn via the cube method. Journal of Statistical Planning and Inference, 141, 479-487.

Chauvet, G., Bonnéry, D. and Deville, J.-C. (2011). Optimal inclusion probabilities for balanced sampling. Journal of Statistical Planning and Inference, 141, 984-994.

Choudhry, G.H., Rao, J.N.K. and Hidiroglou, M.A. (2012). On sample allocation for efficient domain estimation. Survey Methodology, 18, 1, 23-29.

Chromy, J. (1987). Design optimization with multiple objectives. Proceedings of the Survey Research Methods Section, American Statistical Association, 194-199.

Cochran, W.G. (1977). Sampling Techniques. New York: John Wiley & Sons, Inc.

Deville, J.-C., and Tillé, Y. (2004). Efficient balanced sampling: The cube method. Biometrika, 91, 893-912.

Deville, J.-C., and Tillé, Y. (2005). Variance approximation under balanced sampling, Journal of Statistical Planning and Inference, 128, 569-591.

Dykstra R. and Wollan P. (1987). Finding I-projections subject to a finite set of linear inequality constraints, Applied Statistics, 36, 377-383.

Ernst, L.R. (1989). Further applications of linear programming to sampling problems. Proceedings of the Survey Research Methods Section, American Statistical Association, 625-631.

Falorsi, P.D., and Righi, P. (2008). A balanced sampling approach for multi-way stratification designs for small area estimation. Survey Methodology, 34, 2, 223-234.

Falorsi, P.D., Orsini, D. and Righi, P. (2006). Balanced and coordinated sampling designs for small domain estimation. Statistics in Transition, 7, 1173-1198.

Gonzalez, J.M., and Eltinge, J.L. (2010). Optimal survey design: A review. Section on Survey Research Methods – JSM 2010, October.

Isaki, C.T., and Fuller, W.A. (1982). Survey design under a regression superpopulation model. Journal of the American Statistical Association, 77, 89-96.

Khan, M.G.M., Mati, T. and Ahsan, M.J. (2010). An optimal multivariate stratified sampling design using auxiliary information: An integer solution using goal programming approach. Journal of Official Statistics, 26, 695-708.

Kokan, A., and Khan, S. (1967). Optimum allocation in multivariate surveys: An analytical solution. Journal of the Royal Statistical Society, Series B, 29, 115-125.

Lu, W., and Sitter, R.R. (2002). Multi-way stratification by linear programming made practical. Survey Methodology, 28, 2, 199-207.

Nedyalkova, D., and Tillé, Y. (2008). Optimal sampling and estimation strategies under the linear model. Biometrika, 95, 521-537.

Tillé, Y. (2006). Sampling Algorithms. Springer-Verlag, New York.

Tillé, Y., and Favre, A.-C. (2005). Optimal allocation in balanced sampling. Statistics and Probability Letters, 74, 31-37.

Winkler, W.E. (2001). Multi-way survey stratification and sampling. Research Report Series, Statistics #2001-01. Statistical Research Division U.S. Bureau of the Census Washington D.C. 20233.

Previous