InseeFr / disaggR

Two-Steps Benchmarks for Time Series Disaggregation (French Quarterly National Accounts methodology)
https://inseefr.github.io/disaggR/
Other
12 stars 6 forks source link

reuseBenchmark : problème si on réutilise le benchmark sur un indicateur dégfini sur une autre plage #64

Closed ThomasLaurentInsee closed 2 years ago

ThomasLaurentInsee commented 2 years ago

On observe un décalage annuel d'une constante lorsqu'on applique un benchmark estimé sur une série brute à une série CVS, lorsque la série CVS est définie sur une plage plus petite que la série brute (mais la plage définie par start et end.domain est bien comprise dans la plage de définition des deux séries).

exemple ci-dessous : il faut imaginer que indic02 est la série cvs, définie sur 2002-2019 alors que la série brute est définie sur 2000-2019, mais on calcule l'étalonnage sur 2002-2019.

manifestement, ça vient de la façon dont est géré l'indétermination entre la constante et le résidu.

compte <- construction
indic <- turnover

compte02 <- window(compte, 2002)
indic02  <- window(indic, 2002)

tsb_1 <- twoStepsBenchmark(indic02, compte02,
                           include.differenciation = TRUE,
                           include.rho = FALSE,
                           start.coeff.calc = 2002,
                           end.coeff.calc = 2019,
                           start.benchmark = 2002,
                           end.benchmark = 2019,
                           start.domain = 2002,
                           end.domain = 2019 )

tsb_1

tsb_2 <- twoStepsBenchmark(indic, compte02,
                           include.differenciation = TRUE,
                           include.rho = FALSE,
                           start.coeff.calc = 2002,
                           end.coeff.calc = 2019,
                           start.benchmark = 2002,
                           end.benchmark = 2019,
                           start.domain = 2002,
                           end.domain = 2019 )

tsb_2

smoothed.part(tsb_1)
smoothed.part(tsb_2)

# Les résidus ne sont pas égaux
smoothed.part(tsb_1)-smoothed.part(tsb_2)

# Si on applique le modèle 2 à l'indic02, on crée un écart en niveau égal à la différence entre les résidus.

tsb_3 <- reUseBenchmark(indic02, benchmark = tsb_2)

tsb_3 - tsb_2
ThomasLaurentInsee commented 2 years ago

Je crois avoir compris le problème :

Dans un modèle en différenciation, on calcule la tendance (constant) en référence à la plage de définition de l'indic, avant application de la plage domain (https://github.com/InseeFr/disaggR/blob/master/R/twoStepsBenchmark.R#L433). Donc si on change la plage de définition de l'indic sans changer l'indic, cela change la valeur de constante.

Il me semble qu'on ferait mieux de définir cette constante en référence à la plage de domaine.

ça fait un calcul un peu long, mais c'est pas non plus trop énorme :

constant <- (1L:NROW(hfserie)+(start(hfserie)[1]+(start(hfserie)[2]-1)/frequency(hfserie)-start.domain)*(frequency(hfserie)/(frequency(lfserie))))*(frequency(lfserie)/frequency(hfserie))^2

constant02 <- (1L:NROW(hfserie02)+(start(hfserie02)[1]+(start(hfserie02)[2]-1)/frequency(hfserie02)-start.domain)*(frequency(hfserie02)/(frequency(lfserie))))*(frequency(lfserie)/frequency(hfserie02))^2

# on fait un cbind de constant et hfserie L 448, ce qui revient à définir constant comme une ts commençant à start(hfserie)
ts(constant02, start=start(hfserie02), frequency = frequency(hfserie02)) - 
  ts(constant, start=start(hfserie), frequency = frequency(hfserie))
arnaud-feldmann commented 2 years ago

Hello @ThomasLaurentInsee Ce problème me semble assez mineur, et me semble être insuffisant pour motiver une complexification du calcul de la trend qui pourrait rendre le code moins lisible. (un lecteur non averti se demanderait pourquoi on détermine le niveau de manière si tarabiscotée et se dirait qu'il y a une importance au niveau de cette série dans le calcul. J'aime bien le 1L:n qui annonce directement la couleur -> le niveau n'a aucune importance)

Plus simplement, j'ai une suggestion assez proche de ce que tu m'as dit tout à l'heure au téléphone : modifier reUseBenchmark, et, de manière très faible, twoStepsBenchmark

Je ne veux pas faire appel à la routine interne _impl dans reUseBenchmark, je pense que c’est assez moche car ça brise l’encapsulation fonctionnelle du code, un peu comme si on modifiait un objet sans ses méthodes.

La meilleur solution me semble être de permettre une matrice qui contient le colname « constant » d’être donnée en arg à twoStepsBenchmark, qui le cas échéant bloque la génération de l'indicateur de constante. Dans la mesure où cet usage est pour l’instant bloqué par :

 if (anyDuplicated(colnames(hfserie))) stop("Invalid colnames for hfserie",
                                             call. = FALSE)

On ne cause aucune rupture de code (avoir "constant" dans les noms n'est pour l'instant pas possible. De plus, tester constant %in% colnames(hfserie) avant de rajouter l'indicateur de constante me semble assez économe.

ThomasLaurentInsee commented 2 years ago

Oui cette solution marche. Une autre solution est de faire de constant un argument de twostepbenchmark, et de ne générer la constante que quand cet argument est null. Avec une fonction constant(benchmark), et l'ajout de constant=constant(benchmark) dans l'appel à TSB de reUsebenchmark, on traiterait de manière proche outliers et constante.

arnaud-feldmann commented 2 years ago

https://github.com/InseeFr/disaggR/tree/proposition_reUseBenchmark tente de répondre à ce bug. (au passage j'ai vu un pb dans les tests de reUseBenchmark qui appliquaient une cvs multi en additive). Ce n'était pas très grave dès lors qu'on ne faisait pas attention à ton point

arnaud-feldmann commented 2 years ago

Les tests externes à reUseBenchmark ne sont pas changés à l'exception du colnames() == "constant" qui ne provoque plus d'erreur. Il faut quand même que je rajoute un test pour vérifier que si c'est la seule colonne, là erreur

arnaud-feldmann commented 2 years ago

si on replace ça par un argument masqué comme je le fais pour set.benchmark, il faut bien que ce soit masqué car ce genre de bidouillage n'a rien à faire dans l'UI principale accessible aux utilisateurs. Et pour masquer, il faut un appel à une fonction implicite à l'intérieur de la fonction

arnaud-feldmann commented 2 years ago

Je fais une proposition alternative toute à l'heure avec argument masqué si tu veux. À choisir la moins moche

arnaud-feldmann commented 2 years ago

Attends la proposition que j'ai faite ne marche pas si on enregistre un modèle, puis un an après on le réutilise (par exemple pour une campagne annuelle) puisqu'il n'étend pas la contante. Tu as peut-être raison, peut-être que modifier la constante est la solution la moins moche finalement. @ThomasLaurentInsee sinon il faut étendre la constante dans reUseBenchmark

ThomasLaurentInsee commented 2 years ago

Attends la proposition que j'ai faite ne marche pas si on enregistre un modèle, puis un an après on le réutilise (par exemple pour une campagne annuelle) puisqu'il n'étend pas la contante. Tu as peut-être raison, peut-être que modifier la constante est la solution la moins moche finalement. @ThomasLaurentInsee sinon il faut étendre la constante dans reUseBenchmark

Oui tu as raison, mais de toute façon reubenchamrk ne peut pas ête utilisé dans ce cas, car le domain doit être le même.

On peut définir la même constante qu'actuellement, avec en référence la date de benchmark (mieux que domain car on peut réutiliser un modèle sur un autre modèle, mais sur un autre benchmark ça n'a pas de sens), avec ce calcul :

constant <- (time(hfserie)-start.benchmark[1])*(frequency(lfserie)/frequency(hfserie))

La seule différence est qu'on commence par zéro et non 1/144. Je mets le [1] pour traiter le cas ou la date de benchmark est de longueur 2 (étalonnage Trim vers mensuel par exemple)

arnaud-feldmann commented 2 years ago

Je suis en train de faire un truc sur ta proposition première proposition, ça brise moins de tests que je pensais. Il faut fuir comme la peste time() pour 2 raisons : c'est une fonction lente pour pas grand chose, et aussi et surtout elle engendre des corruptions de double à cause de seq.int(xtsp[1L], xtsp[2L], length.out = n) qui répartit n'importe comment entre les 2 bornes en cas d'arrondi non-juste, et est très peu soluble avec un floor (par exemple) qui viendrait après. Mieux vaut un seq avec départ exact, période exacte, et ajustement de l'arrivée (c'est un comportement moins surprenant)

arnaud-feldmann commented 2 years ago

de plus il ne faut pas oublier que start.benchmark peuvent être des couples ou des singletons, il faut mieux ne pas les utiliser je pense. Je pars plutôt sur l'idée d'une trend générale qui prend t= 0 comme origine pour tous les étalonnages

ThomasLaurentInsee commented 2 years ago

Oui c'est peut-être le plus simple. Mais je pense que t=2000 ou 2020 est mieux pour que les smoothed_part n'est pas une trop grande constante dans le cas usuel.

Ou alors il faut dire que smoothed_part c'est le résidu plus la constante, ce qui le rend smoothed_part indépendant du choix de normalisation

arnaud-feldmann commented 2 years ago

C'était déjà le cas dans notre procédure sous sas, le niveau de la cale n'a pas d'importance et j'ai déjà prévenu dans la vignette. le choix arbitraire de niveaux n'a pas beaucoup d'importance car aucun des graphiques n'utilise la cale en niveaux. Si j'avais uni la constante avec la cale, ce qui est un peu trop tard dans la mesure où tout le package est conçu sur cette séparation stricte, ça aurait causé d'autres problèmes, bien plus graves qu'un choix de niveaux : impossibilité de séparer la contribution de la cale et celle de la trend. Il ne faut pas oublier que si la cale n'est pas séparée de la constante en niveaux, la séparation est importante en différences

arnaud-feldmann commented 2 years ago

on va sur l'idée d'une base à 2000 du coup ? 2000 c'est bien, ça évite d'avoir des cales à 2000, ça communique bien un côté arbitraire et ça laisse la cale suffisamment éloignée de zéro pour éviter d'inciter trop à donner un sens à un objet qui n'est pas trop fait pour en avoir.

ThomasLaurentInsee commented 2 years ago

Oui 2000 c'est bien.

Oui tu as raison qu'il faut garder les deux séparés pour les contributions, et le fonctionnement du package tourne autour de ça.

Après une petite fonction qui fait la somme des deux aurait son intérêt. Ca a du sens pour regarder la part du compte expliqué par l'indicateur par exemple (entre un pur lissage à 0 et un indicateur parfait à 100). Mais c'est secondaire et il faudrait lui trouver un nom.

arnaud-feldmann commented 2 years ago

Je fais ça calmement et sérieusement mercredi. A plus et merci pour tes remarques J'avais appelé ça "smoothed.part" parce que ya pas de mot pour dire exactement "cale" en anglais mais c'est vrai que c'est difficile de trouver des noms lol 100% explicites. Bonne soirée

arnaud-feldmann commented 2 years ago

@ThomasLaurentInsee C'est bon tu peux checker. Mais finalement j'ai repensé à ce que tu dis et, quitte à changer, peut-être que comme tu le dis, utiliser une base à 0 sur start.benchmark serait mieux ? Je ne sais pas. Dans tous les cas, si on le fait sur start.domain, en plus d'un switch length, il faudrait faire un floor low frequency (car start domain peut commencer par exemple en février). start.benchmark serait plus logique car c'est là que l'agrégation commence... Mais ça obscurcit peut-être le code pour pas grand chose, alors même que ces histoires de fenêtres (éventuellement un peu corrompues d'ailleurs) sont, dans la structure du programme, pensées pour être gérées plus tard. Qu'en penses-tu ? Est-ce satisfaisant avec le simple axiome "2000"?

arnaud-feldmann commented 2 years ago

https://github.com/InseeFr/disaggR/tree/proposition_reUseBenchmark_3 cette branche là est avec start.benchmark

ThomasLaurentInsee commented 2 years ago

Merci Arnaud, on a commencé à regarder la version '2000' : cela semble fonctionner et réduit les plus gros écarts que l'on observait sur la consommation (quelques 100 de M). Mais dans certains cas (un peu baroque) il semble y avoir encore des écarts (2-3 M) mais cela vient peut être de la façon dont on teste. On continue à regarder ça demain.

Je trouve que la version 2000 a le mérite de la simplicité.

arnaud-feldmann commented 2 years ago

@ThomasLaurentInsee et, pour faciliter la lecture, j'ai bien nommé la constante arbitrary_constant lol