'update'

[talweg.git] / reports / report.gj

-----
<h2>Introduction</h2>

J'ai fait quelques essais dans deux configurations pour la méthode "Neighbors"
(la seule dont on a parlé, incorporant désormais la "variante Bruno/Michel").

 * avec simtype="mix" et raccordement simple ("Zero") dans le cas "non local", i.e. on va
   chercher des voisins n'importe où du moment qu'ils correspondent à deux jours consécutifs sans
   valeurs manquantes.
 * avec simtype="endo" et raccordement "Neighbor" dans le cas "local" : voisins de même niveau de
   pollution et même saison.

J'ai systématiquement comparé à une approche naïve : la moyenne des lendemains des jours
"similaires" dans tout le passé, ainsi qu'à la persistence -- reproduisant le jour courant ou
allant chercher le futur similaire une semaine avant.

Ensuite j'affiche les erreurs, quelques courbes prévues/mesurées, quelques filaments puis les
histogrammes de quelques poids. Concernant les graphes de filaments, la moitié gauche du graphe
correspond aux jours similaires au jour courant, tandis que la moitié droite affiche les
lendemains : ce sont donc les voisinages tels qu'utilisés dans l'algorithme.

<%
list_titles = ['Pollution par chauffage', 'Pollution par épandage', 'Semaine non polluée']
list_indices = ['indices_ch', 'indices_ep', 'indices_np']
%>
-----r
library(talweg)

P = ${P} #instant de prévision
H = ${H} #horizon (en heures)

ts_data = read.csv(system.file("extdata","pm10_mesures_H_loc_report.csv",package="talweg"))
exo_data = read.csv(system.file("extdata","meteo_extra_noNAs.csv",package="talweg"))
# NOTE: 'GMT' because DST gaps are filled and multiple values merged in above dataset.
# Prediction from P+1 to P+H included.
data = getData(ts_data, exo_data, input_tz = "GMT", working_tz="GMT", predict_at=P)

indices_ch = seq(as.Date("2015-01-18"),as.Date("2015-01-24"),"days")
indices_ep = seq(as.Date("2015-03-15"),as.Date("2015-03-21"),"days")
indices_np = seq(as.Date("2015-04-26"),as.Date("2015-05-02"),"days")

% for i in range(3):
-----
<h2 style="color:blue;font-size:2em">${list_titles[i]}</h2>
-----r
p_n = computeForecast(data, ${list_indices[i]}, "Neighbors", "Zero", horizon=H,
	simtype="mix", local=FALSE)
p_l = computeForecast(data, ${list_indices[i]}, "Neighbors", "Neighbors", horizon=H,
	simtype="endo", local=TRUE)
p_a = computeForecast(data, ${list_indices[i]}, "Average", "Zero", horizon=H)
p_p = computeForecast(data, ${list_indices[i]}, "Persistence", "Zero", horizon=H,
	same_day=${'TRUE' if loop.index < 2 else 'FALSE'})
-----r
e_n = computeError(data, p_n, H)
e_l = computeError(data, p_nl, H)
e_a = computeError(data, p_a, H)
e_p = computeError(data, p_p, H)
options(repr.plot.width=9, repr.plot.height=7)
plotError(list(e_n, e_p, e_a, e_l), cols=c(1,2,colors()[258], 4))

# Noir: Neighbors non-local, bleu: Neighbors local, vert: moyenne, rouge: persistence

i_np = which.min(e_n$abs$indices)
i_p = which.max(e_n$abs$indices)
-----r
options(repr.plot.width=9, repr.plot.height=4)
par(mfrow=c(1,2))

plotPredReal(data, p_n, i_np); title(paste("PredReal non-loc day",i_np))
plotPredReal(data, p_n, i_p); title(paste("PredReal non-loc day",i_p))

plotPredReal(data, p_l, i_np); title(paste("PredReal loc day",i_np))
plotPredReal(data, p_l, i_p); title(paste("PredReal loc day",i_p))

plotPredReal(data, p_a, i_np); title(paste("PredReal avg day",i_np))
plotPredReal(data, p_a, i_p); title(paste("PredReal avg day",i_p))

# Bleu: prévue, noir: réalisée
-----r
par(mfrow=c(1,2))
f_np_n = computeFilaments(data, p_n, i_np, plot=TRUE); title(paste("Filaments non-loc day",i_np))
f_p_n = computeFilaments(data, p_n, i_p, plot=TRUE); title(paste("Filaments non-loc day",i_p))

f_np_l = computeFilaments(data, p_l, i_np, plot=TRUE); title(paste("Filaments loc day",i_np))
f_p_l = computeFilaments(data, p_l, i_p, plot=TRUE); title(paste("Filaments loc day",i_p))
-----r
par(mfrow=c(1,2))
plotFilamentsBox(data, f_np_n); title(paste("FilBox non-loc day",i_np))
plotFilamentsBox(data, f_p_n); title(paste("FilBox non-loc day",i_p))

# Generally too few neighbors:
#plotFilamentsBox(data, f_np_l); title(paste("FilBox loc day",i_np))
#plotFilamentsBox(data, f_p_l); title(paste("FilBox loc day",i_p))
-----r
par(mfrow=c(1,2))
plotRelVar(data, f_np_n); title(paste("StdDev non-loc day",i_np))
plotRelVar(data, f_p_n); title(paste("StdDev non-loc day",i_p))

plotRelVar(data, f_np_l); title(paste("StdDev loc day",i_np))
plotRelVar(data, f_p_l); title(paste("StdDev loc day",i_p))

# Variabilité globale en rouge ; sur les 60 voisins (+ lendemains) en noir
-----r
par(mfrow=c(1,2))
plotSimils(p_n, i_np); title(paste("Weights non-loc day",i_np))
plotSimils(p_n, i_p); title(paste("Weights non-loc day",i_p))

plotSimils(p_l, i_np); title(paste("Weights loc day",i_np))
plotSimils(p_l, i_p); title(paste("Weights loc day",i_p))

# - pollué à gauche, + pollué à droite
-----r
# Fenêtres sélectionnées dans ]0,7] / non-loc à gauche, loc à droite
p_n$getParams(i_np)$window
p_n$getParams(i_p)$window

p_l$getParams(i_np)$window
p_l$getParams(i_p)$window
% endfor