In spacetime: Classes and Methods for Spatio-Temporal Data

Spatial and spatio-temporal objects in google charts

This vignette shows how Google charts can be used to display some spatio-temporal data sets used in the JSS paper on spacetime.

Google charts are interactive graphs in web pages that use Google proprietary code. R package googleVis converts R data.frame objects into Google charts. This vignette uses knitr and markdown to create a web page from the R-markdown source file. It was inspired by, and copies small sections of, the corresponding googleVis vignette with knitr examples.

Set the googleVis options first to change the behaviour of plot.gvis, so that only the chart component of the HTML file is written into the output file.

library(googleVis)
op <- options(gvis.plot.tag='chart')

Following plot statements for gvis objects will automatically return the HTML required for the 'knitted' output.

Geo Charts

Geo charts work with country or region (state, administrative regions, e.g. NUTS1, formally ISO 3166-2) data. We will try to make the DE_NUTS1 data in spacetime ready for this.

ISOCodes data set from CRAN

ISOcodes are available from CRAN, including a mapping to German states:

library(spacetime)
data(air) # loads rural and DE_NUTS1
rural = STFDF(stations, dates, data.frame(PM10 = as.vector(air)))
library(ISOcodes)
data("ISO_3166_2")
## State names are already in German
ISO_3166_2_DE <- ISO_3166_2[grep("DE-", ISO_3166_2$Code), ]
plot(
  gvisTable(ISO_3166_2_DE)
  )

We load the two data sets rural and DE_NUTS1 from the spacetime package and add the German state names.

ISO_3166_2_DE <- ISO_3166_2_DE[order(ISO_3166_2_DE$Name),]
DE_NUTS1$name <- ISO_3166_2_DE$Name

Plotting Shape_Area, a variable present for all regions in DE_NUTS1:

library(googleVis)
## Create list with options for Geo Chart to be used
geoChartDE <- list(region="DE", 
                   resolution="provinces",
                   legend="{numberFormat:'#,###.00'}") 
plot(
  gvisGeoChart(DE_NUTS1@data, locationvar = "name", 
               colorvar = "Shape_Area",
               options=geoChartDE)
  )

reveals that Brandenburg is not recognized! We will now try with the ISO 3166-2 codes:

DE_NUTS1$code = ISO_3166_2_DE$Code
plot(
  gvisGeoChart(DE_NUTS1@data, locationvar = "code", 
               colorvar = "Shape_Area",
               options=geoChartDE)
  )

which reveals that we now have all regions displayed.

An air quality example Geo chart

We can compute yearly average PM10 concentration over each of the states and present a map with a table next to it:

DE_NUTS1.years = STF(DE_NUTS1, as.Date(c("2008-01-01", "2009-01-01")))
agg = aggregate(rural[,"2008::2009"], DE_NUTS1.years, mean, na.rm=TRUE)
d = agg[,1]@data # select time step one, take attr table of resulting SpatialPolygonsDataFrame object
d$code = ISO_3166_2_DE$Code # add region codes
d$PM10 = round(d$PM10, 2) # avoid 12-digit numbers printed when hovering over

M <- gvisGeoChart(na.omit(d), locationvar = "code", 
                 colorvar = "PM10",
                 options=c(geoChartDE, height=400)) # drop NA values for Geo chart

## Add German state names for table
d$code <- ISO_3166_2_DE$Name
Tbl <- gvisTable(d, options=list(height=400), 
                 formats=list(PM10="#,###.0"))
plot(gvisMerge(M, Tbl, horizontal=TRUE))

The white states received no value and correspond to NA values in the R object d; they were omitted by the na.omit call.

Irish wind station mean and standard deviations

The Irish wind data has 12 stations, shown here on a map, colour denoting mean wind speed, symbols size its standard deviation:

library(sp)
data("wind", package = "gstat")
wind.loc$y = as.numeric(char2dms(as.character(wind.loc[["Latitude"]])))
wind.loc$x = as.numeric(char2dms(as.character(wind.loc[["Longitude"]])))
wind.loc$yx = paste(wind.loc$y, wind.loc$x, sep=":")
m = round(apply(wind,2,mean), 3)
sd = round(apply(wind,2,sd), 3)
wind.loc$mean = m[match(wind.loc$Code, names(m))]
wind.loc$sd = sd[match(wind.loc$Code, names(sd))]
plot( gvisGeoChart(wind.loc, "yx", "mean", "sd", "Station",
                   options = list(region = "IE", 
                                  legend="{numberFormat:'#.00'}"))
      )

Time lines

Two time lines are shown that have interesting interaction possibilities. For each of them, a subset of the full data sets (Irish wind and German air quality) are given, as the full data sets take rather long to process (Rmd to html) and makes the browser slow when loading the graph and interacting with it.

Iris Wind in Time Annotation

We select the last three years of the wind data set, to not overload the browser:

# select:
wind = wind[wind$year > 76,] # select three years, 1976-1978
time = ISOdate(wind$year+1900, wind$month, wind$day)
wind.stack = stack(wind[,4:15])
names(wind.stack) = c("wind_speed", "station")
wind.stack$time = rep(time, 12)
plot(
  gvisAnnotationChart(wind.stack, "time", "wind_speed", "station",
                        options = list(width = "1000px", 
                                       height = "500px"))
)

Irish wind in a Line Chart

The Line Chart also supports time, and allows for identifying separate lines, when clicking the line or the station. It does not allow zooming, and does not work well for many observations or many lines.

wind$time = time
plot(
  gvisLineChart(wind[1:200,], "time", names(wind)[4:9])
  )

Line charts deal well with gaps in time series:

wind$time = time
wind[4:5, 7:9] = NA
plot(
  gvisLineChart(wind[1:10,], "time", names(wind)[4:9])
  )

Rural PM10 data from spacetime: Annotation Chart

The PM10 data from the ``spacetime vignette on spatio-temporal overlay and aggregation'' contain missing values. These are connected by straight lines in the Annotated Time Line, rendering this visualisation useless.

In contrast, Annotation Charts deal with gaps (missing values) in time series, and show them as discontinuities in the lines:

library(spacetime)
data(air)
d = as(rural[1:10,"2001"], "data.frame") # daily, 2001
plot(
  gvisAnnotationChart(d, "time", "PM10", "sp.ID",
                      options = list(width = "1000px", 
                                     height = "500px"))
  )

The result is still a bit flaky; just try it for the full time series.

Motion Chart

The following example shows a motion chart of the Produc data used in the spacetime manual, but does not convert the data. Time is simply year (integer).

data("Produc", package = "plm")
plot(
  gvisMotionChart(Produc, idvar="state", timevar="year")
)

(Please note that the Motion Chart is only displayed when hosted on a web server, or if placed in a directory which has been added to the trusted sources in the [security settings of Macromedia] (http://www.macromedia.com/support/documentation/en/flashplayer/help/settings_manager04.html). See the googleVis package vignette for more details. )

## Set options back to original options
options(op)

spacetime documentation built on April 6, 2023, 1:09 a.m.