fix window bounds
[talweg.git] / reports / report.gj
CommitLineData
63ff1ecb 1-----
ff5df8e3 2<h2>Introduction</h2>
63ff1ecb
BA
3
4J'ai fait quelques essais dans différentes configurations pour la méthode "Neighbors"
ee8b1b4e
BA
5(la seule dont on a parlé) et sa variante récente appelée pour l'instant "Neighbors2",
6avec simtype="mix" : deux types de similarités prises en compte, puis multiplication des poids.
7Pour Neighbors on prédit le saut (par la moyenne pondérée des sauts passés), et pour Neighbors2
8on n'effectue aucun raccordement (prévision directe).
63ff1ecb 9
ee8b1b4e
BA
10J'ai systématiquement comparé à une approche naïve : la moyenne des lendemains des jours
11"similaires" dans tout le passé, ainsi qu'à la persistence -- reproduisant le jour courant ou
12allant chercher le futur similaire une semaine avant.
63ff1ecb
BA
13
14Ensuite j'affiche les erreurs, quelques courbes prévues/mesurées, quelques filaments puis les
15histogrammes de quelques poids. Concernant les graphes de filaments, la moitié gauche du graphe
16correspond aux jours similaires au jour courant, tandis que la moitié droite affiche les
17lendemains : ce sont donc les voisinages tels qu'utilisés dans l'algorithme.
18
19<%
20list_titles = ['Pollution par chauffage', 'Pollution par épandage', 'Semaine non polluée']
21list_indices = ['indices_ch', 'indices_ep', 'indices_np']
22%>
63ff1ecb 23-----r
63ff1ecb
BA
24library(talweg)
25
d09b09b0
BA
26P = ${P} #instant de prévision
27H = ${H} #horizon (en heures)
28
63ff1ecb
BA
29ts_data = read.csv(system.file("extdata","pm10_mesures_H_loc_report.csv",package="talweg"))
30exo_data = read.csv(system.file("extdata","meteo_extra_noNAs.csv",package="talweg"))
ee8b1b4e
BA
31# NOTE: 'GMT' because DST gaps are filled and multiple values merged in above dataset.
32# Prediction from P+1 to P+H included.
33data = getData(ts_data, exo_data, input_tz = "GMT", working_tz="GMT", predict_at=P)
63ff1ecb
BA
34
35indices_ch = seq(as.Date("2015-01-18"),as.Date("2015-01-24"),"days")
36indices_ep = seq(as.Date("2015-03-15"),as.Date("2015-03-21"),"days")
37indices_np = seq(as.Date("2015-04-26"),as.Date("2015-05-02"),"days")
d4841a3f 38
ff5df8e3 39% for i in range(3):
63ff1ecb 40-----
ff5df8e3 41<h2 style="color:blue;font-size:2em">${list_titles[i]}</h2>
63ff1ecb 42-----r
ee8b1b4e
BA
43p_nn = computeForecast(data, ${list_indices[i]}, "Neighbors", "Neighbors", horizon=H)
44p_nn2 = computeForecast(data, ${list_indices[i]}, "Neighbors2", "Zero", horizon=H)
45p_az = computeForecast(data, ${list_indices[i]}, "Average", "Zero", horizon=H)
46p_pz = computeForecast(data, ${list_indices[i]}, "Persistence", "Zero", horizon=H, same_day=${'TRUE' if loop.index < 2 else 'FALSE'})
63ff1ecb 47-----r
ee8b1b4e
BA
48e_nn = computeError(data, p_nn, H)
49e_nn2 = computeError(data, p_nn2, H)
d09b09b0
BA
50e_az = computeError(data, p_az, H)
51e_pz = computeError(data, p_pz, H)
63ff1ecb 52options(repr.plot.width=9, repr.plot.height=7)
ee8b1b4e 53plotError(list(e_nn, e_pz, e_az, e_nn2), cols=c(1,2,colors()[258], 4))
63ff1ecb 54
ee8b1b4e 55# Noir: Neighbors, bleu: Neighbors2, vert: moyenne, rouge: persistence
63ff1ecb 56
ee8b1b4e
BA
57i_np = which.min(e_nn$abs$indices)
58i_p = which.max(e_nn$abs$indices)
63ff1ecb
BA
59-----r
60options(repr.plot.width=9, repr.plot.height=4)
61par(mfrow=c(1,2))
62
ee8b1b4e
BA
63plotPredReal(data, p_nn, i_np); title(paste("PredReal nn day",i_np))
64plotPredReal(data, p_nn2, i_p); title(paste("PredReal nn day",i_p))
63ff1ecb 65
ee8b1b4e
BA
66plotPredReal(data, p_nn2, i_np); title(paste("PredReal nn2 day",i_np))
67plotPredReal(data, p_nn2, i_p); title(paste("PredReal nn2 day",i_p))
63ff1ecb
BA
68
69plotPredReal(data, p_az, i_np); title(paste("PredReal az day",i_np))
70plotPredReal(data, p_az, i_p); title(paste("PredReal az day",i_p))
71
ff5df8e3 72# Bleu: prévue, noir: réalisée
63ff1ecb
BA
73-----r
74par(mfrow=c(1,2))
ee8b1b4e
BA
75f_np = computeFilaments(data, p_nn, i_np, plot=TRUE); title(paste("Filaments nn day",i_np))
76f_p = computeFilaments(data, p_nn, i_p, plot=TRUE); title(paste("Filaments nn day",i_p))
63ff1ecb 77
ee8b1b4e
BA
78f_np2 = computeFilaments(data, p_nn2, i_np, plot=TRUE); title(paste("Filaments nn2 day",i_np))
79f_p2 = computeFilaments(data, p_nn2, i_p, plot=TRUE); title(paste("Filaments nn2 day",i_p))
63ff1ecb
BA
80-----r
81par(mfrow=c(1,2))
ee8b1b4e
BA
82plotFilamentsBox(data, f_np); title(paste("FilBox nn day",i_np))
83plotFilamentsBox(data, f_p); title(paste("FilBox nn day",i_p))
63ff1ecb 84
ea5c7e56
BA
85# Generally too few neighbors:
86#plotFilamentsBox(data, f_np2); title(paste("FilBox nn2 day",i_np))
87#plotFilamentsBox(data, f_p2); title(paste("FilBox nn2 day",i_p))
63ff1ecb
BA
88-----r
89par(mfrow=c(1,2))
ee8b1b4e
BA
90plotRelVar(data, f_np); title(paste("StdDev nn day",i_np))
91plotRelVar(data, f_p); title(paste("StdDev nn day",i_p))
63ff1ecb 92
ee8b1b4e
BA
93plotRelVar(data, f_np2); title(paste("StdDev nn2 day",i_np))
94plotRelVar(data, f_p2); title(paste("StdDev nn2 day",i_p))
63ff1ecb 95
ff5df8e3 96# Variabilité globale en rouge ; sur les 60 voisins (+ lendemains) en noir
63ff1ecb
BA
97-----r
98par(mfrow=c(1,2))
ee8b1b4e
BA
99plotSimils(p_nn, i_np); title(paste("Weights nn day",i_np))
100plotSimils(p_nn, i_p); title(paste("Weights nn day",i_p))
63ff1ecb 101
ee8b1b4e
BA
102plotSimils(p_nn2, i_np); title(paste("Weights nn2 day",i_np))
103plotSimils(p_nn2, i_p); title(paste("Weights nn2 day",i_p))
63ff1ecb 104
ff5df8e3 105# - pollué à gauche, + pollué à droite
63ff1ecb 106-----r
ee8b1b4e
BA
107# Fenêtres sélectionnées dans ]0,7] / nn à gauche, nn2 à droite
108p_nn$getParams(i_np)$window
109p_nn$getParams(i_p)$window
63ff1ecb 110
ee8b1b4e
BA
111p_nn2$getParams(i_np)$window
112p_nn2$getParams(i_p)$window
63ff1ecb 113% endfor