Add standard deviation as grey smoothed bands to lattice xyplot

Question

It's quite a time I am trying to resolve this question and I can't figure it out... I am sure it's quite simple. I would like to add the +/- MM standard deviation (column SD.MM) to my MM line as a line above/below it and the same for the KK one. Here is my data

data <- read.table(textConnection(
    "Day         Period  MM      KK     SD.MM   SD.KK
1   Tuesday     1   20.0    131.0   1.74    11.40
2   Tuesday     2   76.8    203.0   10.60   28.10
3   Tuesday     3   26.7    118.0   13.00   57.50
4   Wednesday   1   33.8    143.0   2.64    11.20
5   Wednesday   2   74.1    232.0   10.30   32.60
6   Wednesday   3   34.0    130.0   17.60   67.00
7   Thursday    1   46.4    203.0   3.95    16.80
8   Thursday    2   59.5    165.0   8.26    23.10
9   Thursday    3   32.1    120.0   16.50   61.90
10  Friday      1   24.6    98.8    1.96    7.95
11  Friday      2   139.0   367.0   18.20   47.90
12  Friday      3   91.4    216.0   20.90   49.40
13  Saturday    1   30.6    158.0   2.42    11.60
14  Saturday    2   78.2    295.0   8.02    33.40
15  Saturday    3   90.7    310.0   51.60   176.00
16  Sunday      1   28.0    115.0   2.47    9.91
17  Sunday      2   80.4    232.0   10.90   30.40
18  Sunday      3   42.5    211.0   19.70   97.80"),
header=TRUE,as.is=TRUE)

Here the code for producing the image below it.

library(lattice)
library(gridExtra)
a <- xyplot(MM~Period|Day,data=data,layout=c(6,1), type="o", 
            xlab="Period", 
            ylab="Loads",
            ylab.right = "Loads",
            main ="Daily trend",
            index.cond=list(c(5,6,4, 1:3)),
            ylim = c(0,150),
            scales=list(x=list(at=seq(1,3,1),labels=c("B","F","A")),
                        y=(list(alternating=3))),
            key = list(text = list("MM"), points = list(pch=1, col="steelblue2"),
                       text = list("KK"), points = list(pch=1, col="violet"),
                       text = list("d-KK = 306"), lines = list(lty=4, lwd=2, col="green4"),
                       text = list("u-KK = 312"), lines = list(lty=4, lwd=2, col="red"),
                       text = list("s-KK = 41"), lines = list(lty=4, lwd=2, col="gold"),
                       border = F
            )
)

b <- xyplot(KK~Period|Day,data=data,layout=c(6,1), type="o", col="violet",
            xlab="Period", 
            ylab="Loads",
            ylab.right = "Loads",
            index.cond=list(c(5,6,4, 1:3)), 
            ylim = c(0,400),
            scales=list(x=list(at=seq(1,3,1),labels=c("B","F","A")),
                        y=(list(alternating=3))),
            panel = function(x,y,...) {
              panel.abline(h=306,lty = 4, lwd=2, col="green4")
              panel.abline(h = 312, lty = 4, lwd=2, col="red")
              panel.abline(h = 41, lty = 4, lwd=2, col="gold")
              panel.xyplot(x,y,...)
            }
)
grid.arrange(a,b, nrow=2)

I would like something resulting like the grey band on the xyplot here R - Lattice xyplot - How do you add error bars to groups and summary lines?.

Any help ? Thanks in advance !

Answer 1

Its a tricky one with creating a gray smooth lines - your best bid is to play with panel.smoother from latticeExtra package, but it has a lot of weird parameters so its not easy to hack it out of box.

However, doing typical standard deviation bars from your data is easy with xYplot from Hmisc package.

First, you have to define lower and upper limits for your standard deviation bars based on your corresponding data. I just use your standard deviation symmetrically on both sides of the values, but you can divide it by two depending on how you calculated it:

data$MM_sd_hi <- data$MM+data$SD.MM  
data$MM_sd_low <- data$MM-data$SD.MM  
data$KK_sd_hi <- data$KK+data$SD.KK  
data$KK_sd_low <- data$KK-data$SD.KK  

Than you add those limits within Cbind function inside xYplot. You need to increase your ylim a bit more to accommodate this, because xYplot does something weird with scales. Here is the code:

library(Hmisc)
library(gridExtra)
a <- xYplot(Cbind(MM, MM_sd_low, MM_sd_hi)~Period | Day, data=data, layout=c(6,1), 
            type="o", xlab="Period",ylab="Loads", ylab.right = "Loads", 
            main ="Daily trend", index.cond=list(c(5,6,4, 1:3)), ylim = c(0,170),
            scales=list(x=list(at=seq(1,3,1),labels=c("B","F","A")),
                        y=(list(alternating=3))),
            key = list(text = list("MM"), points = list(pch=1, col="steelblue2"),
                       text = list("KK"), points = list(pch=1, col="violet"),
                       text = list("d-KK = 306"), lines = list(lty=4, lwd=2, col="green4"),
                       text = list("u-KK = 312"), lines = list(lty=4, lwd=2, col="red"),
                       text = list("s-KK = 41"), lines = list(lty=4, lwd=2, col="gold"),
                       border = F)
           )   

b <- xYplot(Cbind(KK, KK_sd_low, KK_sd_hi)~Period|Day,data=data,layout=c(6,1), 
            type="o", col="violet", xlab="Period", ylab="Loads", ylab.right = "Loads",
            index.cond=list(c(5,6,4, 1:3)), ylim = c(0,500),
            scales=list(x=list(at=seq(1,3,1),labels=c("B","F","A")),
                        y=(list(alternating=3))),
            panel = function(x,y,...) {
              panel.abline(h=306,lty = 4, lwd=2, col="green4")
              panel.abline(h = 312, lty = 4, lwd=2, col="red")
              panel.abline(h = 41, lty = 4, lwd=2, col="gold")
              panel.xYplot(x,y,...)
            }
)
grid.arrange(a,b, nrow=2)