DATA 624 - Predictive Analystics - (Assignment 1)

EXercise 2.10.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.
• What is the time interval of each series?
• Use autoplot() to produce a time plot of each series.
• For the last plot, modify the axis labels and title.

data("aus_production")
data("pelt")
data("gafa_stock")
data("vic_elec")

help("aus_production")

help("pelt")

help("gafa_stock")

help("vic_elec")

# Bricks from aus_production
aus_production %>% filter(!is.na(Bricks)) -> bricks_data
head(aus_production, 10)%>% kable() %>% 
  kable_styling(bootstrap_options = "striped", font_size = 12,) %>% 
  scroll_box(height = "100%", width = "100%", fixed_thead = T)

# Lynx from pelt
head(pelt, 10)%>% kable() %>% 
  kable_styling(bootstrap_options = "striped", font_size = 12,) %>% 
  scroll_box(height = "100%", width = "100%", fixed_thead = T)

# Close from gafa_stock
head(gafa_stock, 10)%>% kable() %>% 
  kable_styling(bootstrap_options = "striped", font_size = 12,) %>% 
  scroll_box(height = "100%", width = "100%", fixed_thead = T)

# Demand from vic_elec
head(vic_elec, 10) %>% kable() %>% 
  kable_styling(bootstrap_options = "striped", font_size = 12,) %>% 
  scroll_box(height = "100%", width = "100%", fixed_thead = T)

# Bricks data
autoplot(bricks_data, Bricks) +
  ggtitle("Bricks Production Over Time") +
  xlab("Year") +
  ylab("Bricks (Millions)")

# Lynx data
autoplot(pelt, Lynx) +
  ggtitle("Lynx Trapping Over Time") +
  xlab("Year") +
  ylab("Lynx pelts")

# Close data
autoplot(gafa_stock, Close) +
  ggtitle("GAFA Stock Close Price Over Time") +
  xlab("Date") +
  ylab("Close Price")

# Demand data
autoplot(vic_elec, Demand) +
  ggtitle("Electricity Demand") +
  xlab("Date") +
  ylab("Electricity Demand (MW)") +
  labs(title = "Electricity Demand Over Time",
       x = "Time",
       y = "Demand (Megawatts)")

EXercise 2.10.2

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

# Inspect the gafa_stock dataset
head(gafa_stock, 10) %>% kable() %>% 
  kable_styling(bootstrap_options = "striped", font_size = 12,) %>% 
  scroll_box(height = "100%", width = "100%", fixed_thead = T)

gafa_stock %>%
  group_by(Symbol) %>% 
  filter(Close == max(Close)) %>% # Use filter() to find the days that correspond to the peak Close price
  arrange(Symbol)

## # A tsibble: 4 x 8 [!]
## # Key:       Symbol [4]
## # Groups:    Symbol [4]
##   Symbol Date        Open  High   Low Close Adj_Close   Volume
##   <chr>  <date>     <dbl> <dbl> <dbl> <dbl>     <dbl>    <dbl>
## 1 AAPL   2018-10-03  230.  233.  230.  232.      230. 28654800
## 2 AMZN   2018-09-04 2026. 2050. 2013  2040.     2040.  5721100
## 3 FB     2018-07-25  216.  219.  214.  218.      218. 58954200
## 4 GOOG   2018-07-26 1251  1270. 1249. 1268.     1268.  2405600

EXercise 2.10.3

Download the file tute1.csv from the book website, open it in Excel (or some other spreadsheet application), and review its contents. You should find four columns of information. Columns B through D each contain a quarterly series, labelled Sales, AdBudget and GDP. Sales contains the quarterly sales for a small company over the period 1981-2005. AdBudget is the advertising budget and GDP is the gross domestic product. All series have been adjusted for inflation.

1. You can read the data into R with the following script:

#Load the Data 
tute1 <- readr::read_csv("tute1.csv")
View(tute1)

2. Convert the data to time series

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

3. 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 you don’t include facet_grid().

When Facet_grid() is removed, the three time series (Sales, AdBudget, and GDP) are plotted on the same graph. Depending on the scale of each series, this makes the plot crowded and harder to interpret. Facet_grid() helps to separate the plots, making it easier to analyze each series independently. In the previous code all three time series graphs are are plotted on individual graphs.

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

EXercise 2.10.4

The USgas package contains data on the demand for natural gas in the US.

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

# Create a tsibble from us_total with year as the index and state as the key
usgas_tsibble <- us_total |>
  as_tsibble(index = year, key = state)

head(usgas_tsibble, 10) %>% kable() %>% 
  kable_styling(bootstrap_options = "striped", font_size = 12,) %>% 
  scroll_box(height = "100%", width = "100%", fixed_thead = T)

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

# Filter New England states
new_england_states <- c("Maine", "Vermont", "New Hampshire", "Massachusetts", "Connecticut", "Rhode Island")

new_england_data <- usgas_tsibble |>
  filter(state %in% new_england_states)

# Plot the annual natural gas consumption by state
new_england_data |>
  ggplot(aes(x = year, y = y, colour = state)) +
  geom_line() +
  ggtitle("Annual Natural Gas Consumption by State (New England)") +
  xlab("Year") +
  ylab("Natural Gas Consumption (Million Cubic Feet)") +
  theme_minimal()

EXercise 2.10.5

Download tourism.xlsx from the book website and read it into R using readxl::read_excel().

tourism_data <- readxl::read_excel("tourism.xlsx")

head(tourism, 10) %>% kable() %>% 
  kable_styling(bootstrap_options = "striped", font_size = 12,) %>% 
  scroll_box(height = "100%", width = "100%", fixed_thead = T)

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

tourism_tsibble <- tourism_data |>
  mutate(Quarter = yearquarter(Quarter)) |>
  as_tsibble(index = Quarter, key = c(Region, State, Purpose))

# View the tsibble
head(tourism_tsibble, 10) %>% kable() %>% 
  kable_styling(bootstrap_options = "striped", font_size = 12,) %>% 
  scroll_box(height = "100%", width = "100%", fixed_thead = T)

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

max_avg_trips <- tourism_tsibble |>
  group_by(Region, Purpose) |>
  summarize(avg_trips = mean(Trips, na.rm = TRUE)) |>
  ungroup() |>
  slice_max(avg_trips)

head(max_avg_trips, 10)%>% kable() %>% 
  kable_styling(bootstrap_options = "striped", font_size = 12,) %>% 
  scroll_box(height = "100%", width = "100%", fixed_thead = T)

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

state_total_trips <- tourism_tsibble |>
  index_by(Quarter) |>
  group_by(State) |>
  summarize(total_trips = sum(Trips, na.rm = TRUE)) |>
  ungroup() |>
  as_tsibble(index = Quarter, key = State)

# View the new tsibble
head(state_total_trips, 10)%>% kable() %>% 
  kable_styling(bootstrap_options = "striped", font_size = 12,) %>% 
  scroll_box(height = "100%", width = "100%", fixed_thead = T)

EXercise 2.10.8

Use the following graphics functions: autoplot(), gg_season(), gg_subseries(), gg_lag(), ACF() and 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.

us_employment |>
  filter(Title == 'Total Private') |>
  autoplot(Employed) +
  ggtitle("autoplot")

us_employment |>
  filter(Title == 'Total Private') |>
  gg_season(Employed) +
  ggtitle("gg_season Plot")

 us_employment |>
  filter(Title == 'Total Private') |>
  gg_subseries(Employed) +
  ggtitle("gg_subseries Plot") 

 us_employment |>
  filter(Title == 'Total Private') |>
  gg_lag(Employed) +
  ggtitle("gg_lag Plot") 

 us_employment |>
  filter(Title == 'Total Private') |>
  ACF(Employed) |>
  autoplot() +
  ggtitle("ACF Plot")

 # Bricks from aus_production
aus_production |>
  autoplot(Bricks) +
  ggtitle("Bricks Production in Australia")

## Warning: Removed 20 rows containing missing values or values outside the scale range
## (`geom_line()`).

aus_production |>
  gg_season(Bricks) +
  ggtitle("gg_season - Bricks Production Plot")

## Warning: Removed 20 rows containing missing values or values outside the scale range
## (`geom_line()`).

aus_production |>
  gg_subseries(Bricks) +
  ggtitle("gg_subseries - Bricks Production Plot")

## Warning: Removed 5 rows containing missing values or values outside the scale range
## (`geom_line()`).

aus_production |>
  gg_lag(Bricks, geom='point') +
  ggtitle("gg_lag - Bricks Production Plot")

## Warning: Removed 20 rows containing missing values (gg_lag).

aus_production |>
  ACF(Bricks) |>
  autoplot() +
  ggtitle("ACF - Bricks Production Plot")

# Hare from pelt
autoplot(pelt, Hare) +
  ggtitle("Hare Pelts")

gg_subseries(pelt, Hare) +
  ggtitle("gg_subseries - Hare Pelts Plot")

gg_lag(pelt, Hare, geom='point') +
  ggtitle("gg_lag - Hare Pelts Plot")

pelt %>% ACF(Hare) 

## # A tsibble: 19 x 2 [1Y]
##         lag     acf
##    <cf_lag>   <dbl>
##  1       1Y  0.658 
##  2       2Y  0.214 
##  3       3Y -0.155 
##  4       4Y -0.401 
##  5       5Y -0.493 
##  6       6Y -0.401 
##  7       7Y -0.168 
##  8       8Y  0.113 
##  9       9Y  0.307 
## 10      10Y  0.340 
## 11      11Y  0.296 
## 12      12Y  0.206 
## 13      13Y  0.0372
## 14      14Y -0.153 
## 15      15Y -0.285 
## 16      16Y -0.295 
## 17      17Y -0.202 
## 18      18Y -0.0676
## 19      19Y  0.0956

PBS |>
  filter(ATC2 == 'H02') |>
  autoplot(Cost) +
  ggtitle("H02 Drug Cost from PBS")

PBS |>
  filter(ATC2 == 'H02') |>
  gg_season(Cost) +
  ggtitle("gg_season - H02 Drug Cost from PBS")

PBS |>
  filter(ATC2 == 'H02') |>
  gg_subseries(Cost) +
  ggtitle("gg_subseries - H02 Drug Cost from PBS")

PBS |>
  ACF(Cost) |>
  ggtitle("ACF - H02 Drug Cost from PBS")

## $title
## # A tsibble: 7,690 x 6 [1M]
## # Key:       Concession, Type, ATC1, ATC2 [336]
##    Concession   Type        ATC1  ATC2       lag     acf
##    <chr>        <chr>       <chr> <chr> <cf_lag>   <dbl>
##  1 Concessional Co-payments A     A01         1M  0.737 
##  2 Concessional Co-payments A     A01         2M  0.466 
##  3 Concessional Co-payments A     A01         3M  0.239 
##  4 Concessional Co-payments A     A01         4M  0.0575
##  5 Concessional Co-payments A     A01         5M  0.0120
##  6 Concessional Co-payments A     A01         6M -0.0325
##  7 Concessional Co-payments A     A01         7M -0.0219
##  8 Concessional Co-payments A     A01         8M  0.0389
##  9 Concessional Co-payments A     A01         9M  0.180 
## 10 Concessional Co-payments A     A01        10M  0.373 
## # ℹ 7,680 more rows
## 
## $subtitle
## [1] "ACF - H02 Drug Cost from PBS"
## 
## attr(,"class")
## [1] "labels"

us_gasoline |>
  autoplot(Barrels) +
  ggtitle("US Gasoline Barrels Production")

us_gasoline |>
  gg_season(Barrels) +
  ggtitle("gg_season - US Gasoline Barrels Production")

us_gasoline |>
  gg_subseries(Barrels) +
  ggtitle("gg_subseries - US Gasoline Barrels Production")

us_gasoline |>
  gg_lag(Barrels, geom='point') +
  ggtitle("gg_lag - US Gasoline Barrels Production")

us_gasoline |>
  ACF(Barrels) |>
  ggtitle("ACF - US Gasoline Barrels Production")

## $title
## # A tsibble: 31 x 2 [1W]
##         lag   acf
##    <cf_lag> <dbl>
##  1       1W 0.893
##  2       2W 0.882
##  3       3W 0.873
##  4       4W 0.866
##  5       5W 0.847
##  6       6W 0.844
##  7       7W 0.832
##  8       8W 0.831
##  9       9W 0.822
## 10      10W 0.808
## # ℹ 21 more rows
## 
## $subtitle
## [1] "ACF - US Gasoline Barrels Production"
## 
## attr(,"class")
## [1] "labels"