Mini-TP : Calage sur marges

Auteur·rice

Paul Géhin

Date de publication

24 avr. 2026

🍰 Exercice 1 : Estimation du nombre de boulangeries dans le Nord. 🍞

Total

Un statisticien souhaite déterminer le nombre de boulangeries dans le département du Nord. Pour faire cela, il souhaite effectuer un sondage dans lequel il tirera des communes parmi les 638 appartenant au département et demandera le nombre total de boulangeries dans la commune.

Dans un premier temps, le statisticien décide d’utiliser un plan aléatoire simple sans remise de taille 60.

Les résultats sont disponibles dans le fichier suivant : ech_srs.csv.

Nous disposons cette fois d’informations auxiliaires disponibles pour l’ensemble de communes du Nord dans la table var_aux.csv.

  1. Décrivez la population, la variable d’intérêt et la fonction d’intérêt. Dans la suite, nous noterons \(\mathcal{U}\), la population et \(y_k\) la valeur de la variable d’intérêt pour l’individu \(k\).

La population est l’ensemble des communes du Nord, la variable d’intérêt (pour une commune) est le nombre de boulangeries au sein de la commune et la fonction d’intérêt est le nombre total de boulangeries dans la commune.

  1. On rappelle les résultats obtenus avec l’estimateur d’Horvitz-Thompson.
library("sampling")

ech <- read.csv("data/ech_srs.csv", sep = ",", colClasses = c("COG" = "character"))

#Vérification de la table
str(ech)

n <- 60
N <- 638

#Calcul des probabilités d'inclusion d'ordre 1
pik_srs <- rep(n/N, n)
tot_boul_srs <- HTestimator(y = ech$BOULANGERIES, pik = pik_srs)

cat("L'estimation du total basée sur l'estimateur d'Horvitz-Thompson avec un échantillon tiré selon un SRS donne : ", tot_boul_srs, " boulangeries.")

#Calcul des probabilités d'inclusion d'ordre 2
pikl <- matrix(n*(n-1)/(N*(N-1)), ncol = n, nrow = n)
# Par contre il faut remplacer les termes diagonaux car pi_kk = pi_k
diag(pikl) <- n/N

#Methode 1 : estimateur de Horvitz-Thompson de la variance de l'estimateur du total de HT
estim_var_tot_boul_srs <- varHT(y = ech$BOULANGERIES, pikl, method = "1")
#Methode 2 : estimateur de Sen-Yates-Grundy de la variance de l'estimateur du total de HT
varHT(y = ech$BOULANGERIES, pikl, method = "2")

#Même résultat car SRS.
niveau <- 0.90
alpha <- 1 - niveau

quantile_norm <- qnorm(p = 1- (alpha/2))

#Borne inférieure :
tot_boul_srs - quantile_norm * sqrt(estim_var_tot_boul_srs)

#Borne supérieure :
tot_boul_srs + quantile_norm * sqrt(estim_var_tot_boul_srs)
  1. Calculez les totaux des variables auxiliaires.
var_cal <- read.csv("data/var_aux.csv", sep = ",", colClasses = c("COG" = "character"))

tot <- colSums(var_cal[, c("boucherie", "csp_plus","pop")])

print(tot)
  1. Calculez l’estimateur calé en utilisant la méthode linéaire en utilisant respectivement comme variable de calage : boucherie, csp_plus, pop puis les trois variables en même temps. Calculez les variances associées et comparez au cas où l’estimateur est non calé.
ech <- read.csv("data/ech_srs.csv", sep = ",", colClasses = c("COG" = "character"))
ech <- merge(ech, var_cal, by = "COG", sort = FALSE)


n <- 60
N <- 638
ech$prob_inclu <- n/N

g_bouch <- calib(ech$boucherie, d=1/ech$prob_inclu, tot["boucherie"] , method = "linear")

#Calcul de la matrice des pikl du SRS(n,N)
pikl <- matrix((n*(n-1))/(N*(N-1)), n, n)
diag(pikl) <- n/N
 crossprod(g_bouch/ech$prob_inclu, ech$boucherie)

tot["boucherie"]
#Calcul de la variance

#Pour la variable boucherie
g_bouch <- calib(ech$boucherie, d=1/ech$prob_inclu, tot["boucherie"] , method = "linear")
tot_cal_csp <- crossprod(g_bouch/ech$prob_inclu, ech$BOULANGERIES)
calibev(ech$BOULANGERIES, ech$boucherie, tot["boucherie"], pikl, 1/ech$prob_inclu, g_bouch, with = TRUE)

#Pour la var csp
g_csp <- calib(ech$csp_plus, d=1/ech$prob_inclu, tot["csp_plus"] , method = "linear")
tot_cal_csp <- crossprod(g_csp/ech$prob_inclu, ech$BOULANGERIES)
calibev(ech$BOULANGERIES, ech$csp_plus, tot["csp_plus"], pikl, 1/ech$prob_inclu, g_csp, with = TRUE)

## Vérification que tous les poids sont plus grands que 0 (sachant que les poids initiaux plus grand que 0)
all(g_csp > 0)
#Pour la variable pop

g_pop <- calib(ech$pop, d=1/ech$prob_inclu, tot["pop"] , method = "linear")
tot_cal_pop <- crossprod(g_pop/ech$prob_inclu, ech$BOULANGERIES)
calibev(ech$BOULANGERIES, ech$pop, tot["pop"], pikl, 1/ech$prob_inclu, g_pop, with = TRUE)
## Vérification que tous les poids sont plus grands que 0 (sachant que les poids initiaux plus grand que 0)
all(g_pop > 0)

#Utilisation des trois variables initiales.

g_trois <- calib(ech[, c("pop","csp_plus", "boucherie")], d=1/ech$prob_inclu, tot[c("pop","csp_plus", "boucherie")] ,
 method = "linear")
tot_cal_trois <- crossprod(g_trois/ech$prob_inclu, ech$BOULANGERIES)
calibev(ech$BOULANGERIES, ech[, c("pop","csp_plus", "boucherie")],  tot[c("pop","csp_plus", "boucherie")] ,
 pikl, 1/ech$prob_inclu, g_trois, with = TRUE)