2 Functional data in a finite population

Hervé Cardot, Alain Dessertaine, Camelia Goga, Étienne Josserand and Pauline Lardin

Previous | Next

Consider a finite population $U=\left\{\mathrm{1,...,}N\right\}$ of size N and assume that for each unit $k$ in the population $U,$ we can observe the deterministic curve $Y_k={\left(Y_k\left(t\right)\right)}_{t\in \left[\mathrm{0,}T\right]}$. The objective is to estimate the mean curve of the population, which is defined for any instant $t\in \left[\mathrm{0,}T\right]\mathrm{,}$ by

$$\mu \left(t\right)=\frac{1}{N}{\displaystyle \sum _{k\in U}}{Y}_k\left(t\right).$$

Let $s$ be a sample of fixed size $n,$ selected randomly in $U,$ according to a sampling design $p\left(\mathrm{.}\right)\mathrm{.}$ Let ${\pi }_k=\mathrm{Pr}\left(k\in s\right)$ and ${\pi }_{kl}=\mathrm{Pr}\left(k\&l\in s\right)$ be the first- and second- order inclusion probabilities respectively. Assume that ${\pi }_k>0$ for any unit $k$ in population $U\mathrm{.}$

The mean curve $\mu$ is estimated using the Horvitz-Thompson estimator (Cardot et al. 2010) as follows:

$$\widehat{\mu }\left(t\right)=\frac{1}{N}{\displaystyle \sum _{k\in s}}\frac{Y_k\left(t\right)}{{\pi }_k}=\frac{1}{N}{\displaystyle \sum _{k\in U}}\frac{Y_k\left(t\right)}{{\pi }_k}{1}_{k\in s}\mathrm{,}t\in \left[\mathrm{0,}T\right]\mathrm{,}\left(2.1\right)$$

where $1_{k\in s}$ is the indicator that unit $k$ belongs to the sample $s.$ For each instant $t\in \left[\mathrm{0,}T\right]\mathrm{,}$ the estimator $\widehat{\mu }\left(t\right)$ is unbiased for $\mu \left(t\right)\mathrm{,}$ meaning that $E\left(\widehat{\mu }\left(t\right)\right)=\mu \left(t\right)$ where the expectation is considered in relation to the sampling design.

Generally, the trajectories $Y_k\left(t\right)$ are not observed continuously for $t\in \left[\mathrm{0,}T\right]$ but only for a set of $D$ measurement instants $0=t_1<t_2<\dots <t_D=T\mathrm{.}$  In functional data analysis, a classical strategy is to interpolate or smooth discretized trajectories to obtain objects that are truly functions (Ramsay and Silverman 2005). This also makes it possible to deal with curves whose measurement instants are not identical. In the context of surveys, Cardot and Josserand (2011) studied linear interpolation where there is no measurement error at the discretized points, while Cardot et al. (2013) examined smoothing procedures. If there are enough discretization points and the trajectories are fairly regular (but not necessarily derivable), the approximation error due to smoothing or interpolation is negligible in relation to the sampling error. We subsequently assume that the trajectories are observed at any point $t$ of the interval $\left[\mathrm{0,}T\right]\mathrm{.}$

The Horvitz-Thompson covariance function $\gamma \left(r,t\right)=\text{cov}\left(\widehat{\mu }\left(r\right)\mathrm{,}\widehat{\mu }\left(t\right)\right)$ is given by

$$\gamma \left(r\mathrm{,}t\right)=\frac{1}{N^2}{\displaystyle \sum _{k\in U}}{\displaystyle \sum _{l\in U}}{\Delta }_{kl}\frac{Y_k\left(r\right)}{{\pi }_k}\frac{Y_l\left(t\right)}{{\pi }_l}$$

for any $\left(r\mathrm{,}t\right)\in \left[\mathrm{0,}T\right]\times \left[\mathrm{0,}T\right]$ and ${\Delta }_{kl}={\pi }_{kl}-{\pi }_k{\pi }_l\mathrm{.}$ If we assume that the second-order probabilities of inclusion satisfy ${\pi }_{kl}>\mathrm{0,}$ an unbiased estimator of $\gamma \left(r,t\right)$ is given by the Horvitz-Thompson unbiased estimator of the variance,

$$\widehat{\gamma }\left(r\mathrm{,}t\right)=\frac{1}{N^2}{\displaystyle \sum _{k\in s}}{\displaystyle \sum _{l\in s}\frac{{\Delta }_{kl}}{{\pi }_{kl}}\frac{Y_k\left(r\right)}{{\pi }_k}\frac{Y_l\left(t\right)}{{\pi }_l}}\left(2.2\right)$$

for any $\left(r\mathrm{,}t\right)\in \left[\mathrm{0,}T\right]\times \left[\mathrm{0,}T\right]\mathrm{.}$

2.1 Using auxiliary information for estimating the mean trajectory

It is well known that using auxiliary information that effectively explains the variable of interest can greatly improve the precision of the Horvitz-Thompson estimator. In the case of the EDF data, the outside temperature or the type of contract could probably be useful auxiliary variables. A stratification based on geographic position would also yield estimates for different regions. In this study, we have as an auxiliary variable the total electricity consumption for the previous week. We assume that this variable (a real one) is observed for all units in the population.

In this section, we present the Horvitz-Thompson estimator for the mean curve as well as an estimate of the covariance function of this estimator, both for a stratified design using simple random sampling without replacement (SRSWOR) in each stratum, denoted hereafter as STRAT, and for PPS sampling without replacement, which will be denoted as $\pi ps$. We also consider an estimator of the mean curve, assisted by a functional linear model.

2.1.1 Stratified sampling with SRSWOR in each stratum (STRAT)

The population $U$ is assumed to be stratified into a fixed number $H$ of strata $U_1\mathrm{,}\dots \mathrm{,}{U}_H$ of sizes $N_1\mathrm{,}\dots \mathrm{,}{N}_H\mathrm{.}$ Within each stratum $U_h,$ a sample $s_h$ of size $n_h$ is drawn according to an SRSWOR design.

We denote ${\mu }_h\left(t\right)={\displaystyle {\sum }_{k\in U_h}{Y}_k\left(t\right)/N_h}$, for $t\in \left[\mathrm{0,}T\right]$, the mean curve in each stratum, and ${\widehat{\mu }}_h\left(t\right)={\displaystyle {\sum }_{k\in s_h}{Y}_k\left(t\right)/n_h}$, its estimate. The estimator of the mean curve $\mu$ is then defined by

$${\widehat{\mu }}_{\text{strat}}\left(t\right)=\frac{1}{N}{\displaystyle \sum _{h=1}^H{N}_h}{\widehat{\mu }}_h\left(t\right)={\displaystyle \sum _{h=1}^H\frac{N_h}{N}}\left(\frac{1}{n_h}{\displaystyle \sum _{k\in s_h}{Y}_k\left(t\right)}\right)\mathrm{,}t\in \left[\mathrm{0,}T\right]\mathrm{.}\left(2.3\right)$$

The Horvitz-Thompson estimator of the covariance function $\gamma$ is then

$${\widehat{\gamma }}_{strat}\left(r\mathrm{,}t\right)=\frac{1}{N^2}{\displaystyle \sum _{h=1}^H}{N}_h^2\left(\frac{1}{n_h}-\frac{1}{N_h}\right)S_{Y\left(r\right)Y\left(t\right)\mathrm{,}{s}_h}r\mathrm{,}t\in \left[\mathrm{0,}T\right]\mathrm{,}\left(2.4\right)$$

where

$$S_{Y\left(r\right)Y\left(t\right)\mathrm{,}{s}_h}=\frac{1}{n_h-1}{\displaystyle \sum _{k\in s_h}}\left(Y_k\left(r\right)-{\widehat{\mu }}_h\left(r\right)\right)\left(Y_k\left(t\right)-{\widehat{\mu }}_h\left(t\right)\right)$$

is the estimator of the covariance function $S_{Y\left(r\right)Y\left(t\right)\mathrm{,}{U}_h}$ in stratum $h.$ For $r=t\in \left[\mathrm{0,}T\right]\mathrm{,}$ we obtain the estimator of the variance function as follows:

$${\widehat{\gamma }}_{strat}\left(r\right)=\frac{1}{N^2}{\displaystyle \sum _{h=1}^H}{N}_h^2\left(\frac{1}{n_h}-\frac{1}{N_h}\right)S_{Y\left(r\right)\mathrm{,}{s}_h}^2\mathrm{,}$$

where

$$S_{Y\left(r\right)\mathrm{,}{s}_h}^2=\frac{1}{n_h-1}{\displaystyle \sum _{k\in s_h}}{\left(Y_k\left(r\right)-{\widehat{\mu }}_h\left(r\right)\right)}^2$$

is the estimator of the variance $S_{Y\left(r\right)\mathrm{,}{U}_h}^2$ in stratum $h\mathrm{.}$ Cardot and Josserand (2011) propose an extension, in the functional framework, of Neyman's optimal allocation. When the sizes $n_h$ of the samples $s_h$ verify

$$n_h=n\frac{N_h\sqrt{{\displaystyle {\int }_0^T}{S}_{Y\left(r\right)\mathrm{,}{U}_h}^{2}dr}}{{\displaystyle \sum _{h=1}^H{N}_h}\sqrt{{\displaystyle {\int }_0^T}{S}_{Y\left(r\right)\mathrm{,}{U}_h}^{2}dr}}\mathrm{,}h=\mathrm{1,...,}H\mathrm{,}\left(2.5\right)$$

the integrated variance, ${\displaystyle {\int }_0^T}{\widehat{\gamma }}_{\text{strat}}\left(t\right)dt,$ of the stratified estimator is minimized. This allocation is similar to the one obtained in a multivariate context by Cochran (1977). By replacing the variable $Y$ by another variable $X$ that is known for the entire population and is highly correlated with the variable of interest, we obtain an allocation that can be described as $x­$ optimal.

Note 2.1 For $H=\mathrm{1,}$ we obtain the simple random design without replacement (SRSWOR), and the mean curve $\mu \left(t\right)$ is estimated by

$${\widehat{\mu }}_{\text{srswor}}\left(t\right)=\frac{1}{n}{\displaystyle \sum _{k\in s}{Y}_k\left(t\right)\mathrm{,}}t\in \left[\mathrm{0,}T\right]\mathrm{.}\left(2.6\right)$$

The estimator of the covariance function defined in (2.2) is then

$${\widehat{\gamma }}_{\text{srswor}}\left(r,t\right)=\left(\frac{1}{n}-\frac{1}{N}\right)S_{Y\left(r\right)Y\left(t\right)\mathrm{,}s}\mathrm{.}\left(2.7\right)$$

2.1.2  PPS sampling without replacement ($\pi ps$)

PPS sampling designs with or without replacement are often used in practice because they are more effective than equal probability designs when the variable of interest is basically proportional to an auxiliary variable $X$ that has strictly positive values.

In the case of samples of fixed size $n$ drawn without replacement, it is possible to give the equivalent of the formula of Yates and Grundy (1953) and Sen (1953). The covariance function of $\widehat{\mu }$ verifies

$$\gamma \left(r,t\right)=-\frac{1}{2}\frac{1}{N^2}{\displaystyle \sum _{k\in U}}{\displaystyle \sum _{l\in U\mathrm{,}l\ne k}}\left({\pi }_{kl}-{\pi }_k{\pi }_l\right)\left(\frac{Y_k\left(r\right)}{{\pi }_k}-\frac{Y_l\left(r\right)}{{\pi }_l}\right)\left(\frac{Y_k\left(t\right)}{{\pi }_k}-\frac{Y_l\left(t\right)}{{\pi }_l}\right)\mathrm{,}r\mathrm{,}t\in \left[\mathrm{0,}T\right]\mathrm{.}\left(2.8\right)$$            Assume that the values $x_k$ of variable $X$ are known for all units $k$ in the population. It is then possible to define the inclusion probabilities as follows:

$${\pi }_k=n\frac{x_k}{{\displaystyle \sum _{k\in U}}{x}_k}\mathrm{.}$$

Methods have been proposed in the literature for the case ${\pi }_k>1$ (Särndal et al. 1992).

Second-order inclusion probabilities are generally very difficult to calculate for $\pi ps$ designs, and therefore Formula (2.2) cannot be used. However, there is a simple asymptotic approximation of the variance, which was proposed by Hájek (1964) and which entails only first-order inclusion probabilities. This approximation proves to be very effective when the sample is large and the entropy of the sampling design is close to maximum entropy. To select sample $s$ with inclusion probabilities ${\pi }_k,$ the cube algorithm (Deville and Tillé 2004) balanced on the variable $\pi ={\left({\pi }_k\right)}_{k\in U}$ can be used. Deville and Tillé (2005) show that for this particular sampling design, the Hàjek formula is highly effective for estimating the variance of a total or a mean. This formula for approximating the variance can also be used for the covariance, which is then estimated by

$${\widehat{\gamma }}_{\pi \text{ps}}\left(r,t\right)=\frac{1}{N^2}{\displaystyle \sum _{k\in s}}\left(1-{\pi }_k\right)\left(\frac{Y_k\left(r\right)}{{\pi }_k}-\widehat{R}\left(r\right)\right)\left(\frac{Y_k\left(t\right)}{{\pi }_k}-\widehat{R}\left(t\right)\right)\mathrm{,}r\mathrm{,}t\in \left[\mathrm{0,}T\right]\mathrm{,}\left(2.9\right)$$

where

$$\widehat{R}\left(t\right)=\frac{{\displaystyle \sum _{k\in s}}\frac{Y_k\left(t\right)}{{\pi }_k}\left(1-{\pi }_k\right)}{{\displaystyle \sum _{k\in s}}\left(1-{\pi }_k\right)}\text{.}$$

We also used the systematic sampling with unequal probabilities proposed by Madow (1949), since it is simple to use. Unfortunately, it is difficult to estimate the variance for this type of design, and we will therefore not use it to construct confidence bands.

2.2 The model-assisted estimator

Consider $p$ real auxiliary variables $X_1\mathrm{,}\dots \mathrm{,}{X}_p$ and let $x_{kj}$ be the value of the variable $X_j$ for the $k^{th}$ individual. Let $x_k=\left(x_{k1}\mathrm{,...,}{x}_{kp}\right)\prime$ denote the vector containing the values of $p$ auxiliary variables measured on the $k^{th}$ individual. We consider that the relationship between the variable of interest and the auxiliary variables is modeled by the following superpopulation model

$$\xi :Y_k\left(t\right)={x^{\prime }}_k\beta \left(t\right)+{\epsilon }_{kt}\mathrm{,}t\in \left[\mathrm{0,}T\right]\left(2.10\right)$$

with

$$E_{\xi }\left({\epsilon }_{kt}\right)=\mathrm{0,}{E}_{\xi }\left({\epsilon }_{kt}{\epsilon }_{l{t}^{\prime }}\right)=0\text{ for }k\ne l\text{ and }{E}_{\xi }\left({\epsilon }_{kt}{\epsilon }_{k{t}^{\prime }}\right)={\sigma }_{t{t}^{\prime }}^2\text{ for }k=l.$$

This model is an immediate generalization of the functional linear model proposed by Faraway (1997) to several auxiliary variables.

The estimate of $\beta$ based on the model $\xi$ and the sampling design $p\left(\mathrm{.}\right)$ is given by

$$\widehat{\beta }\left(t\right)={\left({\displaystyle \sum _{k\in s}}\frac{x_k{x^{\prime }}_k}{{\pi }_k}\right)}^{-1}{\displaystyle \sum _{k\in s}}\frac{x_k{Y}_k\left(t\right)}{{\pi }_k}\mathrm{,}t\in \left[\mathrm{0,}T\right]\mathrm{.}\left(2.11\right)$$

Note that the sampling weights do not depend on the time $t\in \left[\mathrm{0,}T\right]\mathrm{.}$ Let ${\widehat{Y}}_k\left(t\right)={x^{\prime }}_k\widehat{\beta }\left(t\right)$ be the estimator based on the sampling design for the prediction of $Y_k\left(t\right)$ under the model $\xi$. By direct analogy with the univariate case (Särndal et al. 1992), we finally obtain the following estimator for the mean, for $t\in \left[\mathrm{0,}T\right]\mathrm{,}$

$$\begin{array}{c}{\widehat{\mu }}_{MA}\left(t\right)=\frac{1}{N}{\displaystyle \sum _{k\in s}{\widehat{Y}}_k\left(t\right)}-\frac{1}{N}{\displaystyle \sum _{k\in s}}\frac{\left({\widehat{Y}}_k\left(t\right)-Y_k\left(t\right)\right)}{{\pi }_k}\left(2.12\right)\\ =\frac{1}{N}{\displaystyle \sum _{k\in U}}\frac{Y_k\left(t\right)-{x^{\prime }}_k\widehat{\beta }\left(t\right)}{{\pi }_k}+\frac{1}{N}\left({\displaystyle \sum _{k\in U}{x^{\prime }}_k}\right)\widehat{\beta }\left(t\right).\end{array}$$

If the $\xi$ contains the constant variable 1, then the estimator becomes

$${\widehat{\mu }}_{MA}\left(t\right)=\frac{1}{N}{\displaystyle \sum _{k\in U}{\widehat{Y}}_k\left(t\right)\mathrm{,}}t\in \left[\mathrm{0,}T\right]\mathrm{.}\left(2.13\right)$$

For fixed $r$ and $t$, the asymptotic covariance of ${\widehat{\mu }}_{MA}\left(r\right)$ and ${\widehat{\mu }}_{MA}\left(t\right)$ can be calculated according to the classical residual technique (Särndal et al. 1992),

$${\gamma }_{MA}\left(r,t\right)\simeq \frac{1}{N^2}{\displaystyle \sum _{k\in U}}{\displaystyle \sum _{l\in U}}\left({\pi }_{kl}-{\pi }_k{\pi }_l\right)\frac{\left(Y_k\left(r\right)-{\tilde{Y}}_k\left(r\right)\right)}{{\pi }_k}\frac{\left(Y_l\left(t\right)-{\tilde{Y}}_l\left(t\right)\right)}{{\pi }_l}\mathrm{,}\left(2.14\right)$$

where ${\tilde{Y}}_k\left(r\right)={x^{\prime }}_k\tilde{\beta }\left(t\right)$ is the prediction of $Y_k\left(t\right)$ under the superpopulation model and $\tilde{\beta }\left(t\right)={\left({\displaystyle {\sum }_U}{x}_k{x^{\prime }}_k\right)}^{-1}\left({\displaystyle {\sum }_U}{x}_k{Y}_k\left(t\right)\right)$ is the estimate of $\beta$ at the level of the population and $r\mathrm{,}t\in \left[\mathrm{0,}T\right].$  Cardot, Goga and Lardin (2013) show that this result remains valid uniformly in $r\mathrm{,}t\in \left[\mathrm{0,}T\right]\mathrm{.}$

As an estimator of the covariance function ${\gamma }_{MA}\left(r\mathrm{,}t\right)$, we propose the Horvitz-Thompson estimator of asymptotic covariance given by (2.14) where $\tilde{\beta }\left(t\right)$ is replaced by its estimator $\widehat{\beta }\left(t\right)$ based on the sampling design,

$${\widehat{\gamma }}_{MA}\left(r,t\right)=\frac{1}{N^2}{\displaystyle \sum _{k\mathrm{,}l\in s}}\frac{{\pi }_{kl}-{\pi }_k{\pi }_l}{{\pi }_{kl}}\frac{\left(Y_k\left(r\right)-{\widehat{Y}}_k\left(r\right)\right)}{{\pi }_k}\frac{\left(Y_l\left(t\right)-{\widehat{Y}}_l\left(t\right)\right)}{{\pi }_l}\mathrm{,}r\mathrm{,}t\in \left[\mathrm{0,}T\right]\mathrm{.}\left(2.15\right)$$

Note 2.2 It is entirely possible to consider a superpopulation model $\xi$ that is more general than the linear model proposed here. Estimation techniques based on smoothing by B-splines (Goga and Ruiz-Gazen 2012) can then also be considered. In our study, the relationship between consumption at instant $t$ and the mean consumption for the previous week is almost linear (cf. Figure 4.1), which justifies not using these non-parametric approaches.

Previous | Next