--- title: Essais et développement de la méthode de calage sur marges keywords: fastai sidebar: home_sidebar nb_path: "notebooks/calmar/calmar.ipynb" ---
{% raw %}
{% endraw %}

Imports

{% raw %}
import unittest

tc = unittest.TestCase()
{% endraw %} {% raw %}
import pandas as pd
import seaborn as sns
from leximpact_socio_fisca_simu_etat.config import Configuration
from matplotlib import pyplot as plt

config = Configuration(project_folder="leximpact-prepare-data")
{% endraw %}

Code issu de survey_manager

{% raw %}
{% endraw %} {% raw %}

linear[source]

linear(u)

Args: u:

Returns:

{% endraw %} {% raw %}
{% endraw %} {% raw %}

linear_prime[source]

linear_prime(u)

Args: u:

Returns:

{% endraw %} {% raw %}
{% endraw %} {% raw %}

raking_ratio[source]

raking_ratio(u)

Args: u:

Returns:

{% endraw %} {% raw %}
{% endraw %} {% raw %}

raking_ratio_prime[source]

raking_ratio_prime(u)

Args: u:

Returns:

{% endraw %} {% raw %}
{% endraw %} {% raw %}

logit[source]

logit(u, low, up)

Args: u: low: up:

Returns:

{% endraw %} {% raw %}
{% endraw %} {% raw %}

logit_prime[source]

logit_prime(u, low, up)

Args: u: low: up:

Returns:

{% endraw %} {% raw %}
{% endraw %} {% raw %}

build_dummies_dict[source]

build_dummies_dict(data)

Args: data:

Returns:

{% endraw %} {% raw %}
{% endraw %} {% raw %}

calmar[source]

calmar(data_in, margins, initial_weight='wprm_init', method='linear', lo=None, up=None, use_proportions=False, xtol=1.49012e-08, maxfev=256)

Calibrates weights to satisfy margins constraints

Args: data_in (pd.DataFrame): The observations data margins (dict): Margins is a dictionnary containing for each variable as key the following values

  - a scalar for numeric variables
  - a dictionnary with categories as key and populations as values
  - eventually a key named `total_population` with value the total population. If absent it is initialized to the actual total population
initial_weight (str, optional): Initial weight variable. Defaults to 'wprm_init'.
method (str, optional): Calibration method. Should be 'linear', 'raking ratio' or 'logit'. Defaults to 'linear'.
lo (float, optional): Lower bound on weights ratio. Mandatory when using logit method. Should be < 1. Defaults to None.
up (float, optional): Upper bound on weights ratio. Mandatory when using logit method. Should be > 1. Defaults to None.
use_proportions (bool, optional): When True use proportions if total population from margins doesn't match total population. Defaults to False.
xtol (float, optional): Relative precision on lagrangian multipliers.  Defaults to 1.49012e-08 (fsolve xtol).
maxfev (int, optional): Maximum number of function evaluation. Defaults to 256.

Raises: Exception: [description] Exception: [description] Exception: [description]

Returns: np.array: Margins adjusting weights float: Lagrangian parameter dict: Updated margins

{% endraw %} {% raw %}
{% endraw %} {% raw %}

check_calmar[source]

check_calmar(data_in, margins, initial_weight='wprm_init', pondfin_out=None, lambdasol=None, margins_new_dict=None)

Args: data_in: margins: initial_weight: (Default value = 'wprm_init') pondfin_out: (Default value = None) lambdasol: (Default value = None) margins_new_dict: (Default value = None)

Returns:

{% endraw %} {% raw %}
{% endraw %}

Essais

Génération d'un échantillon

{% raw %}
erfs_03 = pd.read_hdf(
    config.get("DATA_OUT") + "03_erfs_rfr_cal_ind" + config.get("YEAR_ERFS") + ".h5"
)
erfs_03.columns
erfs_03.tail()
{% endraw %} {% raw %}
dataset = erfs_03[
    [
        "idfoy",
        "quifoy",
        "wprm",
        "salaire_de_base",
        "chomage_brut",
        "retraite_brute",
        "f4ba",
        "rfr",
    ]
]

# Pour tester sur une partie de la base
# dataset = dataset.sample(n=1000, random_state=53)
print("On a ", len(dataset), "individus dans cet échantillon")

# Ici, wprm sont les poids initiaux dits "de sondage"
dataset.rename(columns={"wprm": "wprm_init"}, inplace=True)
wprm_init = dataset["wprm_init"]

dataset.head()
{% endraw %} {% raw %}
margins_avant = {}

margins_avant.update(
    {
        "salaire_de_base": (dataset["wprm_init"] * dataset["salaire_de_base"]).sum(),
        "chomage_brut": (dataset["wprm_init"] * dataset["chomage_brut"]).sum(),
        "retraite_brute": (dataset["wprm_init"] * dataset["retraite_brute"]).sum(),
        "f4ba": (dataset["wprm_init"] * dataset["f4ba"]).sum(),
    }
)

margins_avant
{% endraw %} {% raw %}
tc.assertEqual(dataset.isna().sum().sum(), 0)

# Si c'est le cas, on pourra les enlever:
# dataset = dataset.fillna(0) ou dataset = dataset.dropna(0)
{% endraw %} {% raw %}
sns.pairplot(
    dataset[["salaire_de_base", "chomage_brut", "retraite_brute", "rfr"]],
    diag_kind="kde",
)
{% endraw %}

TODO: voir si on fait une normalization

Calage sur marges

Marges

{% raw %}
margins = {}
margins.update(
    {
        "salaire_de_base": 650855163531.0,  # Salaire imposable POTE 2019
        "chomage_brut": 0.97
        * (
            dataset["wprm_init"] * dataset["chomage_brut"]
        ).sum(),  # On ne connait pas le chiffre, donc on met un ecart pour calibrer
        "retraite_brute": 307254581479.0,  # Retraite POTE 2019
        "f4ba": 1.34 * (dataset["wprm_init"] * dataset["f4ba"]).sum(),
    }
)  # On ne connait pas le chiffre, donc on met un ecart pour calibrer

# On ajoute une limite de population
margins.update({"total_population": 38549926})
margins
{% endraw %}

Calage

Calibrates weights to satisfy margins constraints

Args:

  • data_in (pd.DataFrame): The observations data
  • margins (dict): Margins is a dictionnary containing for each variable as key the following values
    • a scalar for numeric variables
    • a dictionnary with categories as key and populations as values
    • eventually a key named total_population with value the total population. If absent it is initialized to the actual total population
  • initial_weight (str, optional): Initial weight variable. Defaults to 'wprm_init'.
  • method (str, optional): Calibration method. Should be 'linear', 'raking ratio' or 'logit'. Defaults to 'linear'.
  • lo (float, optional): Lower bound on weights ratio. Mandatory when using logit method. Should be < 1. Defaults to None.
  • up (float, optional): Upper bound on weights ratio. Mandatory when using logit method. Should be > 1. Defaults to None.
  • use_proportions (bool, optional): When True use proportions if total population from margins doesn't match total population. Defaults to False.
  • xtol (float, optional): Relative precision on lagrangian multipliers. Defaults to 1.49012e-08 (fsolve xtol).
  • maxfev (int, optional): Maximum number of function evaluation. Defaults to 256.
{% raw %}
pondfin_out, lambdasol, margins_new_dict = calmar(
    dataset,
    margins,
    initial_weight="wprm_init",
    method="raking ratio",
    lo=0.33,
    up=3,
    use_proportions=False,
    xtol=1.49012e-08,
    maxfev=256,
)

# pondfin_out
# lambdasol
# margins_new_dict
{% endraw %}

Analyse

On va maintenant regarder l'évolution des distributions avant et après la pondération

{% raw %}
variables = ["salaire_de_base", "chomage_brut", "retraite_brute"]
{% endraw %} {% raw %}
def compute_dataset_pondere(dataset, wprm, type_cal):

    data_pond = dataset.copy()
    for variable in dataset.columns:
        data_pond[variable] = data_pond[variable] * wprm

    data_pond["wprm"] = wprm
    data_pond["type"] = type_cal
    return data_pond
{% endraw %} {% raw %}
def analyse_calage(dataset, pondfin_out):

    # Resultat sur les marges
    check_calmar(
        dataset,
        margins,
        initial_weight="wprm_init",
        pondfin_out=pondfin_out,
        lambdasol=lambdasol,
        margins_new_dict=margins_new_dict,
    )

    # Données pondérées après calage
    data_cal_pond = compute_dataset_pondere(dataset, pondfin_out, "cal")

    # Analyse
    sns.pairplot(
        data_cal_pond,
        diag_kind="kde",
    )
    plt.show()
{% endraw %} {% raw %}
analyse_calage(dataset[variables], pondfin_out)
{% endraw %}

Comparer les méthodes de calibration

Differents calages

Calibrates weights to satisfy margins constraints

Args:

  • data_in (pd.DataFrame): The observations data
  • margins (dict): Margins is a dictionnary containing for each variable as key the following values
    • a scalar for numeric variables
    • a dictionnary with categories as key and populations as values
    • eventually a key named total_population with value the total population. If absent it is initialized to the actual total population
  • initial_weight (str, optional): Initial weight variable. Defaults to 'wprm_init'.
  • method (str, optional): Calibration method. Should be 'linear', 'raking ratio' or 'logit'. Defaults to 'linear'.
  • lo (float, optional): Lower bound on weights ratio. Mandatory when using logit method. Should be < 1. Defaults to None.
  • up (float, optional): Upper bound on weights ratio. Mandatory when using logit method. Should be > 1. Defaults to None.
  • use_proportions (bool, optional): When True use proportions if total population from margins doesn't match total population. Defaults to False.
  • xtol (float, optional): Relative precision on lagrangian multipliers. Defaults to 1.49012e-08 (fsolve xtol).
  • maxfev (int, optional): Maximum number of function evaluation. Defaults to 256.
{% raw %}
pondfin_lin, lambdasol, margins_new_dict = calmar(
    dataset,
    margins,
    initial_weight="wprm_init",
    method="linear",
    lo=0,
    up=3,
    use_proportions=False,
    xtol=1.49012e-08,
    maxfev=256,
)

analyse_calage(dataset[variables], pondfin_lin)
{% endraw %} {% raw %}
pondfin_rr, lambdasol, margins_new_dict = calmar(
    dataset,
    margins,
    initial_weight="wprm_init",
    method="raking ratio",
    lo=0,
    up=3,
    use_proportions=False,
    xtol=1.49012e-08,
    maxfev=256,
)

analyse_calage(dataset[variables], pondfin_rr)
{% endraw %} {% raw %}
pondfin_log, lambdasol, margins_new_dict = calmar(
    dataset,
    margins,
    initial_weight="wprm_init",
    method="logit",
    lo=0,
    up=3,
    use_proportions=False,
    xtol=1.49012e-08,
    maxfev=256,
)

analyse_calage(dataset[variables], pondfin_log)
{% endraw %}

Check des distributions selon la calibration

{% raw %}
def matrice_violin(dataset, results):

    # Initial
    data_init_pond = compute_dataset_pondere(
        dataset[variables], dataset["wprm_init"], "init"
    )
    matrix = data_init_pond.copy()
    # Pour chaque calibration on concatene dans une matrice
    for cal in results:
        data_cal_pond = compute_dataset_pondere(dataset[variables], results[cal], cal)
        matrix = pd.concat([matrix, data_cal_pond])  # .reset_index()

    return matrix
{% endraw %} {% raw %}
results = {}
results.update(
    {"linear": pondfin_lin, "raking_ratio": pondfin_rr, "logit": pondfin_log}
)

len(dataset)
matrix = matrice_violin(dataset, results)
len(matrix)
matrix.head()
matrix.tail()
{% endraw %} {% raw %}
for var in variables:
    sns.violinplot(data=dataset, y=var)
    plt.show()
{% endraw %} {% raw %}
for var in variables:
    ax = sns.violinplot(data=matrix, y=var, x="type", split=True)
    plt.show()
{% endraw %}

Linear Calibration

On voit bien qu'elle détonne. Je suis étonnée des résultats négatifs...

{% raw %}
min(results["linear"])
{% endraw %} {% raw %}
results2 = results.copy()
results2.pop("linear")
results2
{% endraw %} {% raw %}
matrix2 = matrice_violin(dataset, results2)

# On trace les distributions APRES
for var in variables:
    ax = sns.violinplot(data=matrix2, y=var, x="type", split=True)
    plt.show()
{% endraw %}

Analyse de la prédiction

Marges

{% raw %}
margins
{% endraw %} {% raw %}
for var in variables:
    print(var)
    tot_rr = (dataset[var] * results2["raking_ratio"]).sum()
    tot_log = (dataset[var] * results2["logit"]).sum()

    error_rr = abs(tot_rr - margins[var]) / margins[var]
    error_log = abs(tot_log - margins[var]) / margins[var]

    print("Total initial", (dataset[var] * dataset["wprm_init"]).sum())
    print("Total Raking Ratio", tot_rr, ", erreur : ", error_rr)
    print("Total Logit", tot_log, ", erreur : ", error_log, "\n")
{% endraw %}

Prédictions

{% raw %}
tot_rr = (dataset["rfr"] * results2["raking_ratio"]).sum()
tot_log = (dataset["rfr"] * results2["logit"]).sum()
print("Total initial", (dataset["rfr"] * dataset["wprm_init"]).sum())
print("Total Raking Ratio", tot_rr)
print("Total Logit", tot_log)
{% endraw %}

Observation de distributions pondérées

Ici on met en place pour tracer les distributions pondérées car les précédentes ne représentent rien:

[x, wprm] = [5, 10] -> aujourd'hui on plot une personne de salaire 50, alors qu'on devrait plotter 10 personnes de salaire 5 !

Pour représenter ça on a 2 options:

1 - On génère une matrice avec 'wprm_k' personnes pour chaque ligne k (précis)

2 - On fait des histogrames en adaptant la population totale (rapide)

{% raw %}
dataset.head()
{% endraw %} {% raw %}
# results.pop("linear")
{% endraw %} {% raw %}
def calmar_hist(dataset, variables, wprm_init, results, min_=None, max_=None):

    plot_dict = {}
    for var in variables:
        if min_ is None:
            min_ = 0
        if max_ is None:
            max_ = max(dataset[var])

        fig, ax = plt.subplots(figsize=(15, 5))
        plt.hist(
            dataset[var],
            bins=100,
            weights=wprm_init,
            range=(min_, max_),
            histtype="step",
            alpha=0.5,
            linewidth=2,
        )
        legend = ["Initial"]

        for method in results.keys():
            plt.hist(
                dataset[var],
                bins=100,
                weights=results.get(method),
                range=(min_, max_),
                histtype="step",
                alpha=0.5,
            )
            legend.append(method)

        # Mise en forme
        plt.legend(legend)
        plt.title("Comparaison des méthodes de calage sur marges pour le " + var)
        plt.xlabel(var)
        plt.ylabel("Total de population")

        plot_dict.update({f"Calage_{var}": fig})
        plt.savefig(f"Calage_{var}.png")

    return plot_dict
{% endraw %}

Histogrammes comparatifs

Distribution globale

{% raw %}
plot_dict = calmar_hist(dataset, variables, wprm_init, results)
{% endraw %}

Distribution centrale

{% raw %}
plot_dict_cent = calmar_hist(dataset, variables, wprm_init, results, min_=500_000)
{% endraw %}

Distribution hauts revenus

{% raw %}
plot_dict_hr = calmar_hist(dataset, variables, wprm_init, results, min_=2_000_000)
{% endraw %}

Histogrammes comparatifs normés

On a du mal à visualiser les évolutions donc je normalise les paramètres

{% raw %}
def normalizer(dataset, variables):

    norm = {}
    for var in variables:
        mean_ = dataset[var].mean()
        std_ = dataset[var].std()

        dataset[var] = (dataset[var] - mean_) / std_

        norm.update({f"norm_{var}": [mean_, std_]})
        print(var, dataset[var].min(), dataset[var].max())

    return norm, dataset
{% endraw %} {% raw %}
norm, datanorm = normalizer(dataset, variables)
{% endraw %} {% raw %}
plot_dict = calmar_hist(datanorm, variables, wprm_init, results, min_=0)
{% endraw %}

Distribution hauts revenus

{% raw %}
plot_dict = calmar_hist(datanorm, variables, wprm_init, results, min_=90)
{% endraw %}

Conclusion

Il semble que la méthode de calmar nous fasse disparaître les hauts revenus ? On regarde les poids

{% raw %}
var_to_check = "salaire_de_base"
analyse = dataset[[var_to_check, "wprm_init"]]
analyse["logit"] = pondfin_log
analyse["raking_ration"] = pondfin_rr
analyse["linear"] = pondfin_lin

analyse.sort_values(by=var_to_check, inplace=True)
analyse
{% endraw %} {% raw %}
dataset
{% endraw %}

Next steps:

  • Ajouter des hauts revenus AVANT la calibration
  • Voir si un neural network ne permettrait pas une meilleure calibration
  • Est-ce qu'un calage par buckets permettrait d'affiner la chose ? Mais alors il faut des buckets qui regroupent les gens pour toutes les variables confondues...

Reunion Mahdi - 2022-01-26

Checks calmar:

  • ecart poids de calage vs poids de sondage
  • distribution des vars

Calage ERFS 2018 - POTE 2018 -> attention, ici on tire forcement la base dans une direction Vieillissement 2019 -> est-ce que c'est pertinent ? Insertion donnees Rk 2019 Calage avec POTE 2019 -> ca permet de tester la pertinence de nos procedures de vieillissement Vieillissement -> 2021

Vieillissement: il faut adapter la distribution: parfois on a une nouvelle distribution (on adapte les var de l'inflation, le total de pop, mais aussi le nb de gens qui ont/n'ont pas le chomage/la retraite) -> on peut ajouter une cible (=indicatrice) "au chomage" en binaire pour caler sur cette nouvelle marge

Separer l'exercice:

  • il faut avoir la vision economique de ce qu'on fait: quelles sont les indicatrices à utiliser (pour caler et/ou pour vieillir) ? On peut introduire des variables auxilliaires pour representer une distribution (i.e. plus d'infos que juste la somme totale)
  • vieillissement != amelioration de la base
  • il faut definir AVANT chaque etape les mesures de succes de cette operation

Etat de l'art (IPP)

  • un travail sur la base initiale (appariemment et agregats)
  • vieillissement : essentiellement (?) de l'inflation (les != types)
  • voir pour recaler la pyramide des ages

Pour faire du calage en utilisant les buckets de POTE:

  • on introduit des variables auxilliaires qui representent un bucket et on cale sur ces marges là

TO DO

  • Faire tout en raking ratio et voir comment on parameterise la fonction (ecart de poids min et max)
  • Lister nos variables: salaire brut, impot, prestations familiales
  • Haut revenus : voir ce qu'il faut ajouter