[R] Multiple plots on one device using stl

[Ben Madin]

G'day Bill,

On 18/05/2011, at 10:36 AM, <Bill.Venables at csiro.au> wrote:

> If you 
> 
> ?plot.stl
> 
> you will see that that the second argument, set.pars, is a list of argument settings for par(), including a (variable) default setting for mfrow.  I.e. plot.stl overrides your external setting (which will also override any layout() setting).  
> 
> It looks like to override it back, you may need to take charge yourself.  Here is the gist of a partial way round it, perhaps, sort of. ...
> 
> ################
> par(mar = c(0, 6, 0, 6), oma = c(6, 0, 4, 0),
>   tck = -0.01, mfcol = c(4, 2))
> 
> plot(stl.1, set.pars = list())
> plot(stl.2, set.pars = list())
> ################

Perfect - thank you very much - I had played with the set.pars, but it never occurred to me to send it a null list (both times)

> Things might be a bit more flexible if you were to use xyplot.stl in the latticeExtra package.  With lattice the operations of making the plot and displaying it are more cleanly separated.

I'll look into it - thanks for the tip.


cheers

Ben




> I am looking at monthly reports, and have three series of monthly data from 2007 to 2009. I would like to show the season decomposition of these two series side by side on the one device, however using plot doesn't seem to respect any use of layout(matrix(1:3, ncol=3)) or par(mfcol=c(1,3)).
> 
> I'm guessing that this means that the plot(stl) perhaps uses them, but I can't find anywhere the / a plot.stl() - ie, I can't work out where the plot() call is going? 
> 
> I've attached a small example of data and some probably overly verbose code. There are only two sets of example data, not three as mentioned above.
> 
> # load the data
> vals.1 <-
> structure(list(year = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 
> 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
> 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), .Label = c("2007", 
> "2008", "2009"), class = c("ordered", "factor")), month = structure(c(1L, 
> 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 1L, 2L, 3L, 4L, 
> 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 
> 8L, 9L, 10L, 11L, 12L), .Label = c("Jan", "Feb", "Mar", "Apr", 
> "May", "Jun", "Jul", "Aug", "Sep", "Oct", "Nov", "Dec"), class = c("ordered", 
> "factor")), count = c(105, 100, 64, 44, 49, 65, 88, 90, 99, 92, 
> 93, 88, 212, 146, 96, 131, 220, 143, 137, 138, 395, 362, 349, 
> 222, 294, 268, 298, 310, 426, 348, 287, 101, 66, 112, 105, 4)), .Names = c("year", 
> "month", "count"), row.names = c("1", "2", "3", "4", "5", "6", 
> "7", "8", "9", "10", "11", "12", "13", "14", "15", "16", "17", 
> "18", "19", "20", "21", "22", "23", "24", "25", "26", "27", "28", 
> "29", "30", "31", "32", "33", "34", "35", "36"), class = "data.frame")
> vals.2 <-
> structure(list(year = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 
> 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
> 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), .Label = c("2007", 
> "2008", "2009"), class = c("ordered", "factor")), month = structure(c(1L, 
> 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 1L, 2L, 3L, 4L, 
> 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 1L, 2L, 3L, 4L, 5L, 6L, 7L, 
> 8L, 9L, 10L, 11L), .Label = c("Jan", "Feb", "Mar", "Apr", "May", 
> "Jun", "Jul", "Aug", "Sep", "Oct", "Nov", "Dec"), class = c("ordered", 
> "factor")), count = c(45, 34, 17, 6, 7, 16, 12, 11, 14, 17, 11, 
> 20, 27, 12, 10, 14, 22, 23, 92, 144, 385, 274, 320, 252, 240, 
> 146, 222, 142, 122, 117, 163, 89, 51, 89, 108)), .Names = c("year", 
> "month", "count"), row.names = c("1", "2", "3", "4", "5", "6", 
> "7", "8", "9", "10", "11", "12", "13", "14", "15", "16", "17", 
> "18", "19", "20", "21", "22", "23", "24", "25", "26", "27", "28", 
> "29", "30", "31", "32", "33", "34", "35"), class = "data.frame")
> 
> # create the time series
> series.1 <- ts(vals.1$count, start=c(2007,1), freq=12)
> series.2 <- ts(vals.2$count, start=c(2007,1), freq=12)
> 
> # apply the seasonal decomposition
> stl.1 <- stl(series.1, "per", robust=TRUE)
> stl.2 <- stl(series.2, "per", robust=TRUE)
> 
> # set up the device for side by side display
> layout(matrix(1:2, ncol=2))
> 
> # little check
> layout.show(2)
> 
> # plot the first series and second series
> plot(stl.1, labels = c('Count by month','Seasonal Component','Trend','Remainder'), main='Series.1 Decomposition')
> plot(stl.2, labels = c('Count by month','Seasonal Component','Trend','Remainder'), main='Series.2 Decomposition')
> 
> # hmm, used whole device twice, try again
> 
> par(mfcol=c(1,2))
> 
> # and now the second
> 
> plot(stl.1, labels = c('Count by month','Seasonal Component','Trend','Remainder'), main='Series.1 Decomposition')
> plot(stl.2, labels = c('Count by month','Seasonal Component','Trend','Remainder'), main='Series.2 Decomposition')
> 
> # oh what about split.screen()
> 
> split.screen(c(1,2))
> screen(1)
> # now plot
> plot(stl.1, labels = c('Count by month','Seasonal Component','Trend','Remainder'), main='Series.1 Decomposition') 
> 
> # something wrong with the plot, not seeing original series at the top.
> screen(2)
> 
> # Error in par(split.screens[[n]]) : parameter "j" in "mfg" is out of range




-- 

Ben Madin

AusVet Animal Health Services
P.O. Box 5467
Broome   WA   6725
Australia

AusVet's website:  http://www.ausvet.com.au