Reputation: 349
I have two plots displaying supply and demand, and one plot in which I have subtracted the demand from the supply to show the resulting asymmetry. I would like to shade the area between the x-axis and the negative part of the asymmetry, to show the extent of the deficit.
I currently use the following code:
plot.asymmetry <- ggplot(data=df.overview.month,
aes(x=Date.Time, y=Asymmetry)) +
geom_area(data=subset(df.overview.month, Asymmetry < 0),
aes(x=Date.Time, y=Asymmetry))
However - as could be expected - this does not shade the area between geom_line and the x-axis, but only between negative values of the asymmetry data, which is something else entirely, as shown in the resulting graph:
Is there any way to overcome this problem?
/Edit: some example data:
time.initial <- as.POSIXct("2010-12-31 23:00:00", tz="GMT")
Date.Time<-vector()
for(i in 1:24) {
Date.Time[i] <- time.initial + i*3600
}
Demand<-vector()
for(i in 0:23) {
Demand[i+1] <- 155 + 20*sin((pi/12)*i - (pi/2)) + 10*sin((pi/4380)*i + (pi/2))
}
Supply<-vector()
for(i in 0:23) {
Supply[i+1] <- 165 + 5*sin((pi/4380)*i - (pi/2)) + rnorm(1, mean=0, sd=0.20*165)
}
df.overview.month <- data.frame(Date.Time, Demand, Supply, Asymmetry=Supply-Demand)
Upvotes: 4
Views: 4841
Reputation: 77096
Below is some code ported from Matlab to calculate the intersection between segments. If you apply it between the x axis (fixed) and each pair of successive points, you get a list of new coordinates that indicate the crossing points between your geom_line
and the x axis. From this it's an easy step to shade the relevant polygons. Note that I haven't properly tested the ported Matlab code.
## Ported from Matlab to R
## Copyright (c) 2010, U. Murat Erdem
## All rights reserved.
## http://www.mathworks.com/matlabcentral/fileexchange/27205
lineSegmentIntersect <- function(XY1, XY2){
n_rows_1 <- nrow(XY1)
n_cols_1 <- ncol(XY1)
n_rows_2 <- nrow(XY2)
n_cols_2 <- ncol(XY2)
stopifnot(n_cols_1 == 4 && n_cols_2 == 4)
nc <- n_rows_1 * n_rows_2
X1 <- matrix(XY1[,1], nrow=nc, ncol=1)
X2 <- matrix(XY1[,3], nrow=nc, ncol=1)
Y1 <- matrix(XY1[,2], nrow=nc, ncol=1)
Y2 <- matrix(XY1[,4], nrow=nc, ncol=1)
XY2 <- t(XY2)
X3 <- matrix(XY2[1,], nrow=nc, ncol=1)
X4 <- matrix(XY2[3,], nrow=nc, ncol=1)
Y3 <- matrix(XY2[2,], nrow=nc, ncol=1)
Y4 <- matrix(XY2[4,], nrow=nc, ncol=1)
X4_X3 <- X4-X3
Y1_Y3 <- Y1-Y3
Y4_Y3 <- Y4-Y3
X1_X3 <- X1-X3
X2_X1 <- X2-X1
Y2_Y1 <- Y2-Y1
numerator_a <- X4_X3 * Y1_Y3 - Y4_Y3 * X1_X3
numerator_b <- X2_X1 * Y1_Y3 - Y2_Y1 * X1_X3
denominator <- Y4_Y3 * X2_X1 - X4_X3 * Y2_Y1
u_a <- numerator_a / denominator
u_b <- numerator_b / denominator
INT_X <- X1 + X2_X1 * u_a
INT_Y <- Y1 + Y2_Y1 * u_a
INT_B <- (u_a >= 0) & (u_a <= 1) & (u_b >= 0) & (u_b <= 1)
PAR_B <- denominator == 0
COINC_B <- (numerator_a == 0 & numerator_b == 0 & PAR_B)
data.frame(x=INT_X[INT_B], y=INT_Y[INT_B])
}
set.seed(123)
x <- sort(runif(50, -10, 10))
y <- jitter(sin(x), a=2)
n <- length(x)
xy1 <- matrix(c(-10, 0, 10, 0), ncol=4)
xy2 <- cbind(x[-n], y[-n], x[-1], y[-1])
test <- lineSegmentIntersect(xy1, xy2)
library(ggplot2)
d <- data.frame(x=x, y=y)
d2 <- rbind(d, test)
d2 <- subset(d2[order(d2$x), ], y <=0)
p <- qplot(x, y, data=d, geom="path")
p + geom_ribbon(data=d2, aes(ymin = 0, ymax = y), fill="red")
Upvotes: 5
Reputation: 60944
What about this as inspiration. Now you only need to add additional data points where the asymmetry is equal to zero (like @baptiste suggested). I create a new column which is NA
when the asymmetry is above zero, in this way no geom_ribbon will be drawn there. Just subsetting the data will not lead to the required plot.
df.overview.month$Assym_ribbon = ifelse(df.overview.month$Asymmetry > 0,
NA,
df.overview.month$Asymmetry)
ggplot(aes(x = Date.Time, y = Asymmetry),
data = df.overview.month) +
geom_line() +
geom_ribbon(aes(ymin = 0, ymax = Assym_ribbon),
data = , fill = "red")
Some additional notes about the way you constructed your example. The most important one is that R is vectorized. For example:
set.seed(1)
Supply<-vector()
for(i in 0:23) {
Supply[i+1] <- 165 +
5*sin((pi/4380)*i -
(pi/2)) +
rnorm(1, mean=0, sd=0.20*165)
}
is equivalent to:
set.seed(1)
i = 0:23
Supply_vec <- 165 + 5*sin((pi/4380)*i -
(pi/2)) +
rnorm(length(i), mean=0, sd=0.20*165)
> all.equal(Supply_vec, Supply)
[1] TRUE
In this case the reduction in code is modest, but in other (more realistic) settings using vectorization will save you dozens of lines of code.
Upvotes: 3