Exercise 2.1

Explore the following four time series: Bricks from aus_production, Lynx from pelt, Close from gafa_stock, Demand from vic_elec.

Use ? (or help()) to find out about the data in each series.

?aus_production
?pelt
?gafa_stock
?vic_elec

From the help pages: aus_production is quarterly manufacturing output (Beer, Tobacco, Bricks, Cement, Electricity, Gas) from the Australian Bureau of Statistics. pelt is Hudson Bay Company fur trading records for Snowshoe Hare and Canadian Lynx, 1845-1935. gafa_stock is 2014-2018 daily stock prices for Google, Amazon, Facebook and Apple from Yahoo Finance. vic_elec is half-hourly electricity demand for Victoria, Australia, along with Melbourne temperature and a public holiday indicator.

What is the time interval of each series?

Printing each tsibble shows the index and its interval directly.

aus_production
## # A tsibble: 218 x 7 [1Q]
##    Quarter  Beer Tobacco Bricks Cement Electricity   Gas
##      <qtr> <dbl>   <dbl>  <dbl>  <dbl>       <dbl> <dbl>
##  1 1956 Q1   284    5225    189    465        3923     5
##  2 1956 Q2   213    5178    204    532        4436     6
##  3 1956 Q3   227    5297    208    561        4806     7
##  4 1956 Q4   308    5681    197    570        4418     6
##  5 1957 Q1   262    5577    187    529        4339     5
##  6 1957 Q2   228    5651    214    604        4811     7
##  7 1957 Q3   236    5317    227    603        5259     7
##  8 1957 Q4   320    6152    222    582        4735     6
##  9 1958 Q1   272    5758    199    554        4608     5
## 10 1958 Q2   233    5641    229    620        5196     7
## # ℹ 208 more rows
pelt
## # A tsibble: 91 x 3 [1Y]
##     Year  Hare  Lynx
##    <dbl> <dbl> <dbl>
##  1  1845 19580 30090
##  2  1846 19600 45150
##  3  1847 19610 49150
##  4  1848 11990 39520
##  5  1849 28040 21230
##  6  1850 58000  8420
##  7  1851 74600  5560
##  8  1852 75090  5080
##  9  1853 88480 10170
## 10  1854 61280 19600
## # ℹ 81 more rows
gafa_stock
## # A tsibble: 5,032 x 8 [!]
## # Key:       Symbol [4]
##    Symbol Date        Open  High   Low Close Adj_Close    Volume
##    <chr>  <date>     <dbl> <dbl> <dbl> <dbl>     <dbl>     <dbl>
##  1 AAPL   2014-01-02  79.4  79.6  78.9  79.0      67.0  58671200
##  2 AAPL   2014-01-03  79.0  79.1  77.2  77.3      65.5  98116900
##  3 AAPL   2014-01-06  76.8  78.1  76.2  77.7      65.9 103152700
##  4 AAPL   2014-01-07  77.8  78.0  76.8  77.1      65.4  79302300
##  5 AAPL   2014-01-08  77.0  77.9  77.0  77.6      65.8  64632400
##  6 AAPL   2014-01-09  78.1  78.1  76.5  76.6      65.0  69787200
##  7 AAPL   2014-01-10  77.1  77.3  75.9  76.1      64.5  76244000
##  8 AAPL   2014-01-13  75.7  77.5  75.7  76.5      64.9  94623200
##  9 AAPL   2014-01-14  76.9  78.1  76.8  78.1      66.1  83140400
## 10 AAPL   2014-01-15  79.1  80.0  78.8  79.6      67.5  97909700
## # ℹ 5,022 more rows
vic_elec
## # A tsibble: 52,608 x 5 [30m] <Australia/Melbourne>
##    Time                Demand Temperature Date       Holiday
##    <dttm>               <dbl>       <dbl> <date>     <lgl>  
##  1 2012-01-01 00:00:00  4383.        21.4 2012-01-01 TRUE   
##  2 2012-01-01 00:30:00  4263.        21.0 2012-01-01 TRUE   
##  3 2012-01-01 01:00:00  4049.        20.7 2012-01-01 TRUE   
##  4 2012-01-01 01:30:00  3878.        20.6 2012-01-01 TRUE   
##  5 2012-01-01 02:00:00  4036.        20.4 2012-01-01 TRUE   
##  6 2012-01-01 02:30:00  3866.        20.2 2012-01-01 TRUE   
##  7 2012-01-01 03:00:00  3694.        20.1 2012-01-01 TRUE   
##  8 2012-01-01 03:30:00  3562.        19.6 2012-01-01 TRUE   
##  9 2012-01-01 04:00:00  3433.        19.1 2012-01-01 TRUE   
## 10 2012-01-01 04:30:00  3359.        19.0 2012-01-01 TRUE   
## # ℹ 52,598 more rows

aus_production (Bricks) is quarterly [1Q], pelt (Lynx) is annual [1Y], gafa_stock (Close) shows [!] - irregular, since only trading days are included (weekends and holidays are skipped) - and vic_elec (Demand) is half-hourly [30m].

Use autoplot() to produce a time plot of each series.

aus_production |> autoplot(Bricks)

pelt |> autoplot(Lynx)

gafa_stock |> autoplot(Close)

vic_elec |> autoplot(Demand)

For the last plot, modify the axis labels and title.

vic_elec |>
  autoplot(Demand) +
  labs(title = "Half-Hourly Electricity Demand, Victoria",
       x = "Time", y = "Demand (MW)")

Exercise 2.2

Use filter() to find what days corresponded to the peak closing price for each of the four stocks in gafa_stock.

# group by Symbol first so max(Close) is computed per company, not across all four
gafa_stock |>
  group_by(Symbol) |>
  filter(Close == max(Close)) |>
  select(Symbol, Date, Close)
## # A tsibble: 4 x 3 [!]
## # Key:       Symbol [4]
## # Groups:    Symbol [4]
##   Symbol Date       Close
##   <chr>  <date>     <dbl>
## 1 AAPL   2018-10-03  232.
## 2 AMZN   2018-09-04 2040.
## 3 FB     2018-07-25  218.
## 4 GOOG   2018-07-26 1268.

AAPL peaked on 2018-10-03, AMZN on 2018-09-04, FB on 2018-07-25, and GOOG on 2018-07-26.

Exercise 2.3

Download tute1.csv from the book website, open it, and read into R with read_csv().

tute1 <- readr::read_csv("tute1.csv")
View(tute1)

View() opens an interactive data viewer tab in RStudio, which doesn’t work when knitting to a document, so for this PDF I read the data the same way and just print it instead.

tute1 <- readr::read_csv("tute1.csv")
head(tute1, n=10)
## # A tibble: 10 × 4
##    Quarter    Sales AdBudget   GDP
##    <date>     <dbl>    <dbl> <dbl>
##  1 1981-03-01 1020.     659.  252.
##  2 1981-06-01  889.     589   291.
##  3 1981-09-01  795      512.  291.
##  4 1981-12-01 1004.     614.  292.
##  5 1982-03-01 1058.     647.  279.
##  6 1982-06-01  944.     602   254 
##  7 1982-09-01  778.     531.  296.
##  8 1982-12-01  932.     608.  272.
##  9 1983-03-01  996.     638.  260.
## 10 1983-06-01  908.     582.  280.

Convert the data to time series.

mytimeseries <- tute1 |>
  mutate(Quarter = yearquarter(Quarter)) |>
  as_tsibble(index = Quarter)

Construct time series plots of each of the three series.

mytimeseries |>
  pivot_longer(-Quarter) |>
  ggplot(aes(x = Quarter, y = value, colour = name)) +
  geom_line() +
  facet_grid(name ~ ., scales = "free_y")

Check what happens when the facet_grid line is omitted.

mytimeseries |>
  pivot_longer(-Quarter) |>
  ggplot(aes(x = Quarter, y = value, colour = name)) +
  geom_line()

Without facet_grid, all three series get squeezed onto one shared y-axis. GDP is on a much bigger scale than Sales and AdBudget, so those two series are compressed near the bottom, making their patterns difficult to distinguish. facet_grid(scales = "free_y") gives each series its own axis so all three are actually readable.

Exercise 2.4

Install the USgas package.

install.packages("USgas")

Create a tsibble from us_total with year as the index and state as the key.

us_total_ts <- us_total |>
  as_tsibble(index = year, key = state)
us_total_ts
## # A tsibble: 1,266 x 3 [1Y]
## # Key:       state [53]
##     year state        y
##    <int> <chr>    <int>
##  1  1997 Alabama 324158
##  2  1998 Alabama 329134
##  3  1999 Alabama 337270
##  4  2000 Alabama 353614
##  5  2001 Alabama 332693
##  6  2002 Alabama 379343
##  7  2003 Alabama 350345
##  8  2004 Alabama 382367
##  9  2005 Alabama 353156
## 10  2006 Alabama 391093
## # ℹ 1,256 more rows

Plot the annual natural gas consumption by state for the New England area (Maine, Vermont, New Hampshire, Massachusetts, Connecticut, Rhode Island).

new_england <- c("Maine", "Vermont", "New Hampshire",
                  "Massachusetts", "Connecticut", "Rhode Island")

us_total_ts |>
  filter(state %in% new_england) |>
  autoplot(y) +
  labs(title = "Annual Natural Gas Consumption, New England",
       x = "Year", y = "Consumption (million cubic feet)")

Each line above is one state’s annual gas consumption, so this allows comparison of growth and scale across the six New England states.

Exercise 2.5

Download tourism.xlsx from the book website and read it into R.

tourism_raw <- readxl::read_excel("tourism.xlsx")
tourism_raw
## # A tibble: 24,320 × 5
##    Quarter    Region   State           Purpose  Trips
##    <chr>      <chr>    <chr>           <chr>    <dbl>
##  1 1998-01-01 Adelaide South Australia Business  135.
##  2 1998-04-01 Adelaide South Australia Business  110.
##  3 1998-07-01 Adelaide South Australia Business  166.
##  4 1998-10-01 Adelaide South Australia Business  127.
##  5 1999-01-01 Adelaide South Australia Business  137.
##  6 1999-04-01 Adelaide South Australia Business  200.
##  7 1999-07-01 Adelaide South Australia Business  169.
##  8 1999-10-01 Adelaide South Australia Business  134.
##  9 2000-01-01 Adelaide South Australia Business  154.
## 10 2000-04-01 Adelaide South Australia Business  169.
## # ℹ 24,310 more rows

Create a tsibble which is identical to the tourism tsibble from the tsibble package.

tourism_mine <- tourism_raw |>
  mutate(Quarter = yearquarter(Quarter)) |>
  as_tsibble(index = Quarter, key = c(Region, State, Purpose))
tourism_mine
## # A tsibble: 24,320 x 5 [1Q]
## # Key:       Region, State, Purpose [304]
##    Quarter Region   State           Purpose  Trips
##      <qtr> <chr>    <chr>           <chr>    <dbl>
##  1 1998 Q1 Adelaide South Australia Business  135.
##  2 1998 Q2 Adelaide South Australia Business  110.
##  3 1998 Q3 Adelaide South Australia Business  166.
##  4 1998 Q4 Adelaide South Australia Business  127.
##  5 1999 Q1 Adelaide South Australia Business  137.
##  6 1999 Q2 Adelaide South Australia Business  200.
##  7 1999 Q3 Adelaide South Australia Business  169.
##  8 1999 Q4 Adelaide South Australia Business  134.
##  9 2000 Q1 Adelaide South Australia Business  154.
## 10 2000 Q2 Adelaide South Australia Business  169.
## # ℹ 24,310 more rows

Find what combination of Region and Purpose had the maximum number of overnight trips on average.

# Convert to a tibble so Quarter is no longer an implicit grouping variable.
# This produces one average for each Region-Purpose combination.
tourism_mine |>
  as_tibble() |>
  group_by(Region, Purpose) |>
  summarise(avg_trips = mean(Trips)) |>
  arrange(desc(avg_trips))
## # A tibble: 304 × 3
## # Groups:   Region [76]
##    Region          Purpose  avg_trips
##    <chr>           <chr>        <dbl>
##  1 Sydney          Visiting      747.
##  2 Melbourne       Visiting      619.
##  3 Sydney          Business      602.
##  4 North Coast NSW Holiday       588.
##  5 Sydney          Holiday       550.
##  6 Gold Coast      Holiday       528.
##  7 Melbourne       Holiday       507.
##  8 South Coast     Holiday       495.
##  9 Brisbane        Visiting      493.
## 10 Melbourne       Business      478.
## # ℹ 294 more rows

The Sydney-Visiting combination had the highest average number of overnight trips.

Create a new tsibble which combines the Purposes and Regions, and just has total trips by State.

# Summarise in a tibble so Region and Purpose are collapsed into state totals,
# then recreate the tsibble with State as the key.
tourism_states <- tourism_mine |>
  as_tibble() |>
  group_by(State, Quarter) |>
  summarise(Trips = sum(Trips), .groups = "drop") |>
  as_tsibble(index = Quarter, key = State)
tourism_states
## # A tsibble: 640 x 3 [1Q]
## # Key:       State [8]
##    State Quarter Trips
##    <chr>   <qtr> <dbl>
##  1 ACT   1998 Q1  551.
##  2 ACT   1998 Q2  416.
##  3 ACT   1998 Q3  436.
##  4 ACT   1998 Q4  450.
##  5 ACT   1999 Q1  379.
##  6 ACT   1999 Q2  558.
##  7 ACT   1999 Q3  449.
##  8 ACT   1999 Q4  595.
##  9 ACT   2000 Q1  600.
## 10 ACT   2000 Q2  557.
## # ℹ 630 more rows

Exercise 2.8

Use the graphics functions autoplot(), gg_season(), gg_subseries(), gg_lag() and ACF() to explore features from the following time series: “Total Private” Employed from us_employment, Bricks from aus_production, Hare from pelt, “H02” Cost from PBS, and Barrels from us_gasoline. Can I spot any seasonality, cyclicity and trend? What do I learn about the series? What can I say about the seasonal patterns? Can I identify any unusual years?

Total private employment (us_employment)

us_total_private <- us_employment |>
  filter(Title == "Total Private")

us_total_private |> autoplot(Employed)

us_total_private |> gg_season(Employed)

us_total_private |> ACF(Employed) |> autoplot()

us_total_private |>
  gg_subseries(Employed) +
  facet_wrap(vars(id), ncol = 4)

us_total_private |>
  gg_lag(Employed, geom = "point") +
  facet_wrap(vars(.lag), ncol = 3) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5))

Total private employment has a strong long-term upward trend, interrupted by sharp declines in the early 1980s, early 1990s, early 2000s, and especially 2008-09. These declines do not repeat on a fixed annual schedule, so they are more consistent with cyclical movement than seasonality. Monthly seasonality is mild relative to the trend. The near-diagonal lag plots and slowly declining ACF show that employment values remain highly correlated across successive months, as expected for a strongly trended series.

Bricks (aus_production)

aus_production |> autoplot(Bricks)

aus_production |> gg_season(Bricks)

aus_production |> gg_subseries(Bricks)

aus_production |>
  gg_lag(Bricks, geom = "point") +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5))

aus_production |> ACF(Bricks) |> autoplot()

Brick production increased through the 1960s and 1970s, peaked around 1980-81, and then declined overall with substantial irregular variation. Production recovered to a smaller secondary high around 1989 before continuing its longer-run decline. Bricks data end after 2005 Q2, even though the overall aus_production tsibble continues through 2010 Q2 for other variables.

Quarterly seasonality is clear: Q3 is generally highest and Q1 lowest. The lag-4 plot is the tightest relationship because it compares the same quarter one year apart, and the ACF has small peaks at multiples of four. Early 1983 is the most unusual period, with the largest quarter-to-quarter fall followed immediately by the largest increase.

Hare (pelt)

pelt |> autoplot(Hare)

pelt |> ACF(Hare) |> autoplot()

pelt |>
  gg_lag(Hare, geom = "point") +
  facet_wrap(vars(.lag), ncol = 3) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5))

Because pelt is annual data, gg_season() and gg_subseries() do not provide a useful subannual seasonal view and are omitted. The series shows recurring boom-bust cycles roughly 8-10 years apart rather than a long-term trend. The oscillating ACF, with alternating positive and negative values, supports this cyclic pattern. The 1863-64 peak of about 150,000 pelts is unusually high relative to most other peaks.

Cost (PBS, ATC2 == “H02”)

h02 <- PBS |>
  filter(ATC2 == "H02") |>
  select(Month, Concession, Type, Cost) |>
  summarise(TotalC = sum(Cost)) |>
  mutate(Cost = TotalC / 1e6)

h02 |> autoplot(Cost)

h02 |> gg_season(Cost)

h02 |> gg_subseries(Cost) +
  facet_wrap(vars(id), ncol = 4)

h02 |> gg_lag(Cost, geom = "point")

h02 |> ACF(Cost, lag_max = 48) |> autoplot()

H02 (corticosteroid) drug cost has a steady upward trend from 1991 through the mid-2000s, and a very regular yearly pattern: cost jumps in January, drops to its lowest point in February, then climbs gradually back up all year until the next January jump. This is the same Pharmaceutical Benefits Scheme safety-net effect that resets each calendar year, which the book already mentions for the A10 series. The ACF out to 48 lags shows clear peaks every 12 lags, confirming the yearly seasonality on top of the overall upward trend. 2008 looks unusually low, but the data for that year only runs through June, so it’s a partial year rather than a real drop.

US gasoline (us_gasoline)

us_gasoline |> autoplot(Barrels)

us_gasoline |> gg_season(Barrels)

us_gasoline |> ACF(Barrels) |> autoplot()

us_gasoline |>
  gg_lag(Barrels, geom = "point") +
  theme(legend.position = "none")

us_gasoline |>
  gg_subseries(Barrels) +
  facet_wrap(vars(id))

us_gasoline is weekly data, so its natural seasonal period is 52 weeks - the number of panels a subseries plot draws always matches the data’s own time interval, and weekly data always gets a 52-week season. I arranged the 52 panels in a grid for readability. Comparing panels that are far apart - the first few weeks of January against the panels around midyear - the blue mean line sits visibly higher midyear than in January, but that difference is only noticeable comparing panels directly; it does not stand out scanning the grid casually, since neighboring weeks barely differ from each other.

US gasoline supplied increased gradually from the early 1990s through the mid-2000s, then flattened and declined somewhat through the late 2000s and 2010s. Seasonality here is mild - not obvious scanning the season plot, since the strong upward trend dominates it, and only visible in the subseries grid by comparing panels far apart rather than at a glance. Late August and September 2005 stand out as an unusual dip on the time plot, with a sharp drop over a few weeks. The timing is consistent with disruption from Hurricane Katrina, but the plot alone cannot establish the cause. The apparent fall in 2017 reflects data available only through week 3, not a real decline.