X-Git-Url: https://git.auder.net/images/pieces/%22%20%20V.getPpath%28board%5Bi%5D%5Bj%5D%29%20%20%20%22.svg?a=blobdiff_plain;f=pkg%2FR%2Fz_plot.R;fp=pkg%2FR%2Fz_plot.R;h=cf36ed69e08a191c7a011adbb47f4187e92e8cf6;hb=357b8be00388d07c04c10d4bf7f503fd947185d2;hp=0000000000000000000000000000000000000000;hpb=a961f8a15492bf4b18d24bc117358d1f412dd078;p=aggexp.git

diff --git a/pkg/R/z_plot.R b/pkg/R/z_plot.R
new file mode 100644
index 0000000..cf36ed6
--- /dev/null
+++ b/pkg/R/z_plot.R
@@ -0,0 +1,149 @@
+#' @include z_plotHelper.R
+
+#' @title Plot forecasts/observations
+#'
+#' @description Plot the measures at one station versus all experts forecasts.
+#'
+#' @param r Output of \code{\link{runAlgorithm}}.
+#' @param station Name or index of the station to consider. Default: the first one
+#' @param interval Time interval for the plot. Default: all time range.
+#' @param experts Subset of experts for the plot. Default: all experts.
+#' @param ... Additional arguments to be passed to graphics::plot method.
+#'
+#' @export
+plotCurves = function(r, station=1, interval=1:(nrow(r$data)/length(r$stations)), experts=r$experts, cols=rainbow(length(experts)), ...)
+{
+	if (is.character(station))
+		station = match(station, r$stations)
+	if (is.numeric(experts))
+		experts = r$experts[experts]
+
+	XY = subset(r$data[interval,], subset = (Station == station), select = c(experts,"Measure"))
+	indices = getNoNAindices(XY)
+	XY = XY[indices,]
+	X = as.matrix(XY[,names(XY) %in% experts])
+	Y = XY[,"Measure"]
+
+	yRange = range(XY)
+	par(mar=c(5,4.5,1,1), cex=1.5)
+	for (i in 1:length(experts))
+	{
+		plot(X[,i],ylim=yRange,type="l",lty="dotted",col=cols[i],xlab="",ylab="",xaxt="n",yaxt="n", lwd=2, ...)
+		par(new=TRUE)
+	}
+	plot(Y, type="l", ylim=yRange, xlab="", ylab="", lwd=2, cex.axis=1.5, ...)
+	title(xlab="Time",ylab="Forecasts / Measures", cex.lab=1.6)
+	legend("topright", title="Historical PM10",lwd=c(2,1),lty=c("solid","dotted"),horiz=TRUE,legend=c("Measures","Forecasts"))
+}
+
+#' @title Plot error
+#'
+#' @description Plot the absolute error over time at one station.
+#'
+#' @param r Output of \code{\link{runAlgorithm}}.
+#' @param station Name or index of the station to consider. Default: the first one
+#' @param start First index to consider (too much variability in early errors)
+#' @param noNA TRUE to show only errors associated with full lines (old behavior)
+#' @param ... Additional arguments to be passed to graphics::plot method.
+#'
+#' @export
+plotError = function(r, station=1, start=1, noNA=TRUE, ...)
+{
+	if (is.character(station))
+		station = match(station, r$stations)
+
+	XY = subset(r$data, subset = (Station == station), select = c(r$experts,"Measure","Prediction"))
+	Y = XY[,"Measure"]
+	hatY = XY[,"Prediction"]
+	indices = !is.na(Y) & !is.na(hatY)
+	if (noNA)
+	{
+		X = XY[,names(XY) %in% r$experts]
+		indices = indices & getNoNAindices(X)
+	}
+	Y = Y[indices]
+	hatY = hatY[indices]
+
+	error = abs(Y - hatY)
+	par(mar=c(5,4.5,1,1), cex=1.5)
+	plot(error, type="l", xaxt="n", xlab="Time",ylab="L1 error", cex.lab=1.6, cex.axis=1.5, ...)
+	axis(side=1, at=(seq(from=start,to=length(Y),by=30) - start), labels=seq(from=start,to=length(Y),by=30), cex.axis=1.5)
+}
+
+#' @title Plot regret
+#'
+#' @description Plot the regret over time at one station.
+#'
+#' @param r Output of \code{\link{runAlgorithm}}.
+#' @param vs Linear weights to compare with. Can be obtained by the \code{getBestXXX} methods, or by any other mean.
+#' @param station Name or index of the station to consider. Default: the first one
+#' @param start First index to consider (too much variability in early errors)
+#' @param ... Additional arguments to be passed to graphics::plot method.
+#'
+#' @export
+plotRegret = function(r, vs, station=1, start=1, ...)
+{
+	if (is.character(station))
+		station = match(station, r$stations)
+
+	XY = subset(r$data, subset = (Station == station), select = c(r$experts,"Measure","Prediction"))
+	X = XY[,names(XY) %in% r$experts]
+	Y = XY[,"Measure"]
+	hatY = XY[,"Prediction"]
+
+	indices = !is.na(Y) & !is.na(hatY) & getNoNAindices(X)
+	X = as.matrix(X[indices,])
+	Y = Y[indices]
+	hatY = hatY[indices]
+
+	error2 = abs(Y - hatY)^2
+	vsError2 = abs(Y - X %*% vs)^2
+	cumErr2 = cumsum(error2) / seq_along(error2)
+	cumVsErr2 = cumsum(vsError2) / seq_along(vsError2)
+	regret = cumErr2 - cumVsErr2
+
+	par(mar=c(5,4.5,1,1), cex=1.5)
+	plot(regret, type="l", xaxt="n", xlab="Time", ylab="Regret", cex.lab=1.6, cex.axis=1.5, ...)
+	abline(a=0., b=0., col=2)
+	axis(side=1, at=(seq(from=start,to=length(Y),by=30) - start), labels=seq(from=start,to=length(Y),by=30), cex.axis=1.5)
+}
+
+#' @title Plot predicted/expected cloud
+#'
+#' @description Plot the cloud of forecasts/observations + statistical indicators.
+#'
+#' @param r Output of \code{\link{runAlgorithm}}.
+#' @param thresh Threshold to consider for alerts (usually 30 or 50)
+#' @param hintThresh thresholds to draw on the plot to help visualization. Often \code{c(30,50,80)}
+#' @param station Name or index of the station to consider. Default: the first one
+#' @param noNA TRUE to show only errors associated with full lines (old behavior)
+#' @param ... Additional arguments to be passed to graphics::plot method.
+#'
+#' @export
+plotCloud = function(r, thresh=30, hintThresh=c(30,50,80), station=1, noNA=TRUE, ...)
+{
+	if (is.character(station))
+		station = match(station, r$stations)
+
+	XY = subset(r$data, subset = (Station == station), select = c(r$experts,"Measure","Prediction"))
+	Y = XY[,"Measure"]
+	hatY = XY[,"Prediction"]
+	indices = !is.na(Y) & !is.na(hatY)
+	if (noNA)
+	{
+		X = XY[,names(XY) %in% r$experts]
+		indices = indices & getNoNAindices(X)
+	}
+	Y = Y[indices]
+	hatY = hatY[indices]
+
+	indics = getIndicators(r, thresh, station, noNA)
+
+	par(mar=c(5,5,3,2), cex=1.5)
+	plot(Y, hatY, xlab="Measured PM10", ylab="Predicted PM10",
+		cex.lab=1.6, cex.axis=1.5, xlim=c(0,120), ylim=c(0,120), ...)
+	abline(0,1,h=hintThresh,v=hintThresh,col=2,lwd=2)
+#	legend("topleft",legend=c(paste("EV   ",indics$EV),paste("RMSE ",indics$RMSE)),cex=1.2)
+	legend("topleft",legend=paste("RMSE ",indics$RMSE))
+	legend("bottomright",legend=c(paste("TS ",indics$TS)))
+}