Repaso 1

library(fpp3)

── Attaching packages ────────────────────────────────────────────── fpp3 0.5 ──

✔ tibble      3.2.1     ✔ tsibble     1.1.3
✔ dplyr       1.1.4     ✔ tsibbledata 0.4.1
✔ tidyr       1.3.0     ✔ feasts      0.3.1
✔ lubridate   1.9.3     ✔ fable       0.3.3
✔ ggplot2     3.4.4     ✔ fabletools  0.3.4

── Conflicts ───────────────────────────────────────────────── fpp3_conflicts ──
✖ lubridate::date()    masks base::date()
✖ dplyr::filter()      masks stats::filter()
✖ tsibble::intersect() masks base::intersect()
✖ tsibble::interval()  masks lubridate::interval()
✖ dplyr::lag()         masks stats::lag()
✖ tsibble::setdiff()   masks base::setdiff()
✖ tsibble::union()     masks base::union()

library(tidyverse)

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ forcats 1.0.0     ✔ readr   2.1.5
✔ purrr   1.0.2     ✔ stringr 1.5.1

── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter()     masks stats::filter()
✖ tsibble::interval() masks lubridate::interval()
✖ dplyr::lag()        masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

aus_livestock |> 
  distinct(Animal)

# A tibble: 7 × 1
  Animal                    
  <fct>                     
1 Bulls, bullocks and steers
2 Calves                    
3 Cattle (excl. calves)     
4 Cows and heifers          
5 Lambs                     
6 Pigs                      
7 Sheep                     

ovejas <- aus_livestock |> 
  filter(Animal == "Sheep",
         State == "New South Wales")

ovejas

# A tsibble: 558 x 4 [1M]
# Key:       Animal, State [1]
      Month Animal State            Count
      <mth> <fct>  <fct>            <dbl>
 1 1972 Jul Sheep  New South Wales 669400
 2 1972 Aug Sheep  New South Wales 581100
 3 1972 Sep Sheep  New South Wales 468100
 4 1972 Oct Sheep  New South Wales 515300
 5 1972 Nov Sheep  New South Wales 564500
 6 1972 Dec Sheep  New South Wales 565400
 7 1973 Jan Sheep  New South Wales 630300
 8 1973 Feb Sheep  New South Wales 528900
 9 1973 Mar Sheep  New South Wales 594500
10 1973 Apr Sheep  New South Wales 443000
# ℹ 548 more rows

ovejas |> 
  autoplot(Count)

Train/Test splits

ovejas_train <- ovejas |> 
  filter_index(. ~ "Dec. 2012")

ovejas_train

# A tsibble: 486 x 4 [1M]
# Key:       Animal, State [1]
      Month Animal State            Count
      <mth> <fct>  <fct>            <dbl>
 1 1972 Jul Sheep  New South Wales 669400
 2 1972 Aug Sheep  New South Wales 581100
 3 1972 Sep Sheep  New South Wales 468100
 4 1972 Oct Sheep  New South Wales 515300
 5 1972 Nov Sheep  New South Wales 564500
 6 1972 Dec Sheep  New South Wales 565400
 7 1973 Jan Sheep  New South Wales 630300
 8 1973 Feb Sheep  New South Wales 528900
 9 1973 Mar Sheep  New South Wales 594500
10 1973 Apr Sheep  New South Wales 443000
# ℹ 476 more rows

Visualization

Time plot

ovejas_train |> 
  autoplot(Count)

lambda <- ovejas_train |> 
  features(Count, features = guerrero) |> 
  pull(lambda_guerrero)
lambda

[1] 0.5806898

ovejas_train |> 
  autoplot(box_cox(Count, lambda))

Season plot

ovejas_train |> 
  gg_season(Count)

ovejas_train |> 
  gg_season(box_cox(Count, lambda))

ovejas_train |> 
  features(Count, feat_stl)

# A tibble: 1 × 11
  Animal State          trend_strength seasonal_strength_year seasonal_peak_year
  <fct>  <fct>                   <dbl>                  <dbl>              <dbl>
1 Sheep  New South Wal…          0.914                  0.584                  9
# ℹ 6 more variables: seasonal_trough_year <dbl>, spikiness <dbl>,
#   linearity <dbl>, curvature <dbl>, stl_e_acf1 <dbl>, stl_e_acf10 <dbl>

STL Decomposition

ovejas_train |> 
  model(
    stl = STL(box_cox(Count, lambda), robust = TRUE)
  ) |> 
  components() |> 
  autoplot()

Model estimation

ovejas_fit <- ovejas_train |> 
  model(
    snaive_bc = SNAIVE(box_cox(Count, lambda)),
    snaive = SNAIVE(Count),
    combi_bc = decomposition_model(
      STL(box_cox(Count, lambda), robust = TRUE),
      SNAIVE(season_year),
      RW(season_adjust ~ drift())
    ),
    combi = decomposition_model(
      STL(Count, robust = TRUE),
      SNAIVE(season_year),
      RW(season_adjust ~ drift())
    )
  )
ovejas_fit

# A mable: 1 x 6
# Key:     Animal, State [1]
  Animal State           snaive_bc   snaive                  combi_bc
  <fct>  <fct>             <model>  <model>                   <model>
1 Sheep  New South Wales  <SNAIVE> <SNAIVE> <STL decomposition model>
# ℹ 1 more variable: combi <model>

Autocorrelation diagnostics

ovejas_fit |> 
  select(snaive_bc) |> 
  gg_tsresiduals()

Warning: Removed 12 rows containing missing values (`geom_line()`).

Warning: Removed 12 rows containing missing values (`geom_point()`).

Warning: Removed 12 rows containing non-finite values (`stat_bin()`).

ovejas_fit |> 
  select(snaive) |> 
  gg_tsresiduals()

Warning: Removed 12 rows containing missing values (`geom_line()`).

Warning: Removed 12 rows containing missing values (`geom_point()`).

Warning: Removed 12 rows containing non-finite values (`stat_bin()`).

ovejas_fit |> 
  select(combi) |> 
  gg_tsresiduals()

Warning: Removed 12 rows containing missing values (`geom_line()`).

Warning: Removed 12 rows containing missing values (`geom_point()`).

Warning: Removed 12 rows containing non-finite values (`stat_bin()`).

ovejas_fit |> 
  select(combi_bc) |> 
  gg_tsresiduals()

Warning: Removed 12 rows containing missing values (`geom_line()`).

Warning: Removed 12 rows containing missing values (`geom_point()`).

Warning: Removed 12 rows containing non-finite values (`stat_bin()`).

ovejas_fit |> 
  augment() |> 
  features(.innov, ljung_box, dof = 0, lag = 24)

# A tibble: 4 × 5
  Animal State           .model    lb_stat lb_pvalue
  <fct>  <fct>           <chr>       <dbl>     <dbl>
1 Sheep  New South Wales combi        145.         0
2 Sheep  New South Wales combi_bc     138.         0
3 Sheep  New South Wales snaive      1021.         0
4 Sheep  New South Wales snaive_bc   1098.         0

accuracy(ovejas_fit) |> 
  arrange(MAPE)

# A tibble: 4 × 12
  Animal State      .model .type      ME   RMSE    MAE     MPE  MAPE  MASE RMSSE
  <fct>  <fct>      <chr>  <chr>   <dbl>  <dbl>  <dbl>   <dbl> <dbl> <dbl> <dbl>
1 Sheep  New South… combi… Trai…    451. 4.70e4 34632.  -1.16   12.1 0.458 0.457
2 Sheep  New South… combi  Trai…    713. 4.73e4 34649.  -0.863  12.3 0.459 0.459
3 Sheep  New South… snaiv… Trai… -10469. 1.03e5 75551. -11.7    30.7 1     1    
4 Sheep  New South… snaive Trai… -10469. 1.03e5 75551. -11.7    30.7 1     1    
# ℹ 1 more variable: ACF1 <dbl>

Forecast

ovejas_fc <- ovejas_fit |> 
  forecast(h = "6 years")

ovejas_fc

# A fable: 288 x 6 [1M]
# Key:     Animal, State, .model [4]
   Animal State           .model       Month              Count   .mean
   <fct>  <fct>           <chr>        <mth>             <dist>   <dbl>
 1 Sheep  New South Wales snaive_bc 2013 Jan t(N(1435, 255334)) 115655.
 2 Sheep  New South Wales snaive_bc 2013 Feb t(N(1304, 255334))  99577.
 3 Sheep  New South Wales snaive_bc 2013 Mar t(N(1251, 255334))  93376.
 4 Sheep  New South Wales snaive_bc 2013 Apr t(N(1157, 255334))  82864.
 5 Sheep  New South Wales snaive_bc 2013 May t(N(1293, 255334))  98198.
 6 Sheep  New South Wales snaive_bc 2013 Jun t(N(1020, 255334))  68775.
 7 Sheep  New South Wales snaive_bc 2013 Jul t(N(1163, 255334))  83551.
 8 Sheep  New South Wales snaive_bc 2013 Aug t(N(1355, 255334)) 105589.
 9 Sheep  New South Wales snaive_bc 2013 Sep t(N(1708, 255334)) 152765.
10 Sheep  New South Wales snaive_bc 2013 Oct t(N(2119, 255334)) 217309.
# ℹ 278 more rows

ovejas_fc |> 
  autoplot(ovejas |> filter_index("Jan. 2003" ~ .), level = NULL)

ovejas_fc |> 
  accuracy(ovejas) |> 
  arrange(MAPE)

# A tibble: 4 × 12
  .model  Animal State .type      ME   RMSE    MAE   MPE  MAPE  MASE RMSSE  ACF1
  <chr>   <fct>  <fct> <chr>   <dbl>  <dbl>  <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 combi_… Sheep  New … Test  -10324. 47074. 40622. -11.3  25.8 0.538 0.457 0.759
2 snaive… Sheep  New … Test   36872. 65520. 54675.  16.9  29.5 0.724 0.637 0.744
3 snaive  Sheep  New … Test   65819. 85025. 73106.  34.0  39.1 0.968 0.826 0.731
4 combi   Sheep  New … Test   80445. 94525. 80445.  42.2  42.2 1.06  0.918 0.797

ovejas |> 
  model(
    combi_bc = decomposition_model(
      STL(box_cox(Count, lambda), robust = TRUE),
      SNAIVE(season_year),
      RW(season_adjust ~ drift())
    )
  ) |> 
  
  forecast(h = "6 years") |> 
  autoplot(ovejas)

Quarto

Quarto enables you to weave together content and executable code into a finished document. To learn more about Quarto see https://quarto.org.

Running Code

When you click the Render button a document will be generated that includes both content and the output of embedded code. You can embed code like this:

1 + 1

[1] 2

You can add options to executable code like this

[1] 4

The echo: false option disables the printing of code (only output is displayed).