Skip to content

Commit 6d09160

Browse files
ericward-noaaericward-noaa
committed
update
1 parent cd6497a commit 6d09160

File tree

5 files changed

+83
-63
lines changed

5 files changed

+83
-63
lines changed

‎Lectures/Week 6/lec_12_multivariate_GAMs.Rmd‎

Lines changed: 21 additions & 17 deletions
Original file line numberDiff line numberDiff line change
@@ -173,7 +173,7 @@ fit_S <- gam(value ~
173173

174174
## Model S fits
175175

176-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
176+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
177177
long_seal$pred_S <- predict(fit_S, long_seal)
178178
179179
ggplot(long_seal, aes(Year, pred_S)) +
@@ -232,7 +232,7 @@ MARSS model where each time series is a different trend
232232
* `y ~ s(t, k = 12)`
233233
* low uncertainty, poor fit
234234

235-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
235+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
236236
library(mgcv)
237237
set.seed(42)
238238
@@ -273,7 +273,7 @@ ggplot(df_pred, aes(x, fit)) +
273273
* `y ~ s(t, k = 50)`
274274
* good fit, high uncertainty
275275

276-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
276+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
277277
set.seed(42)
278278
279279
# Simulate a seasonal time series
@@ -359,7 +359,7 @@ $$ E[Y_{j,t}] = g^{-1} \left( \beta_0 + \sum_{i=1}^{I} s(x_{i,j,t}) + \sum_{m=1}
359359
* `fit <- mvgam(y ~ s(time, k = 50),
360360
data = df, family = "gaussian")`
361361

362-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
362+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
363363
set.seed(42)
364364
365365
# Simulate a seasonal time series
@@ -433,7 +433,7 @@ plot_mvgam_smooth(fit_GS)
433433

434434
## Including dynamic components: coastal upwelling
435435

436-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
436+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
437437
cuti <- read.csv("CUTI_monthly.csv")
438438
439439
library(dplyr)
@@ -456,14 +456,18 @@ cuti_long <- cuti_subset %>%
456456
cuti_long$latitude <- as.factor(cuti_long$latitude)
457457
cuti_long$time <- sort(rep(1:447, 9))
458458
cuti_long$series <- cuti_long$latitude
459+
```
460+
461+
```{r echo=TRUE, message = TRUE, results = 'asis'}
459462
head(cuti_long)
460463
```
461464

465+
462466
## Including dynamic components: coastal upwelling
463467

464468
* for speed, we'll filter the data to be 2023--2024 (train) and use 2025 for testing
465469

466-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
470+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
467471
cuti2324 <- dplyr::filter(cuti_long, year %in% 2023:2024)
468472
ggplot(cuti2324, aes(year+month/12, cuti)) +
469473
geom_line()+
@@ -505,7 +509,7 @@ plot_mvgam_smooth(fit_trend, smooth = "time")
505509

506510
* Predictions (last 3 points are forecasts)
507511

508-
```{r echo=FALSE, fig.height=4, fig.width=7}
512+
```{r echo=FALSE, fig.height=3, fig.width=7}
509513
newdata <- dplyr::filter(cuti_long, year > 2022)
510514
newdata$pred <- predict(fit_trend, newdata = newdata)[,"Estimate"]
511515
@@ -518,7 +522,7 @@ ggplot(newdata, aes(time, pred)) +
518522
## Including dynamic components: coastal upwelling
519523

520524
* Possibly some problematic residuals (autocorrelated)
521-
```{r echo=TRUE, fig.height=4, fig.width=7}
525+
```{r echo=TRUE, fig.height=3, fig.width=7}
522526
plot(fit_trend)
523527
```
524528

@@ -544,29 +548,29 @@ fit_trend_rw <- mvgam(cuti ~ s(time, bs="cr") +
544548
## Including dynamic components: coastal upwelling
545549

546550
* ACF now looks a little better (and QQ plot too)
547-
```{r echo=TRUE, fig.height=4, fig.width=7}
551+
```{r echo=TRUE, fig.height=3, fig.width=7}
548552
plot(fit_trend_rw)
549553
```
550554

551555
## Including dynamic components: coastal upwelling
552556

553-
```{r echo=TRUE, fig.height=4, fig.width=7}
557+
```{r echo=TRUE, fig.height=3, fig.width=7}
554558
plot_mvgam_smooth(fit_trend_rw)
555559
```
556560

557561
## Including dynamic components: coastal upwelling
558562

559563
* Each of the trends (here for Latitude 31) is ~ flat
560564

561-
```{r echo=TRUE, fig.height=4, fig.width=7}
565+
```{r echo=TRUE, fig.height=3, fig.width=7}
562566
plot_mvgam_trend(fit_trend_rw)
563567
```
564568

565569
## Including dynamic components: coastal upwelling
566570

567571
* Some fits to forecasts are better -- but training data actually seems worse
568572

569-
```{r echo=FALSE, fig.height=4, fig.width=7}
573+
```{r echo=FALSE, fig.height=3, fig.width=7}
570574
newdata <- dplyr::filter(cuti_long, year > 2022)
571575
newdata$pred <- predict(fit_trend_rw, newdata = newdata)[,"Estimate"]
572576
@@ -582,7 +586,7 @@ ggplot(newdata, aes(time, pred)) +
582586

583587
* Let's drop out the main smooth now and let that be explained by the trend
584588

585-
```{r echo=TRUE, fig.height=4, fig.width=7}
589+
```{r echo=TRUE, fig.height=3, fig.width=7}
586590
fit_trend_ar1 <- mvgam(cuti ~ s(time, bs="cr") +
587591
s(latitude, k=9, bs="re") +
588592
s(month, k = 12, bs="cc"),
@@ -595,7 +599,7 @@ fit_trend_ar1 <- mvgam(cuti ~ s(time, bs="cr") +
595599
## Including dynamic components: coastal upwelling
596600

597601
* Trends more flexible than the random walk appear to fit data a little better
598-
```{r echo=FALSE, fig.height=4, fig.width=7}
602+
```{r echo=FALSE, fig.height=3, fig.width=7}
599603
newdata <- dplyr::filter(cuti_long, year > 2022)
600604
newdata$pred <- predict(fit_trend_ar1, newdata = newdata)[,"Estimate"]
601605
@@ -611,7 +615,7 @@ ggplot(newdata, aes(time, pred)) +
611615

612616
* Instead, model trends as dynamic factors (n = 2)
613617

614-
```{r echo=TRUE, fig.height=4, fig.width=7}
618+
```{r echo=TRUE, fig.height=3, fig.width=6.5}
615619
fit_trend_fac <- mvgam(cuti ~ s(time, bs="cr") +
616620
s(latitude, k=9, bs="re") +
617621
s(month, k = 12, bs="cc"),
@@ -631,7 +635,7 @@ fit_trend_fac <- mvgam(cuti ~ s(time, bs="cr") +
631635
* By default mvdlm models these as AR1
632636

633637
* DFA model with AR1 trend + seasonal smooth
634-
```{r echo=TRUE, fig.height=4, fig.width=7}
638+
```{r echo=TRUE, fig.height=3, fig.width=6.5}
635639
plot(fit_trend_fac, type = 'factors')
636640
```
637641

@@ -641,7 +645,7 @@ plot(fit_trend_fac, type = 'factors')
641645

642646
* Factors for this application probably not flexible enough
643647

644-
```{r echo=FALSE, fig.height=4, fig.width=7}
648+
```{r echo=FALSE, fig.height=3, fig.width=6.5}
645649
newdata <- dplyr::filter(cuti_long, year > 2022)
646650
newdata$pred <- predict(fit_trend_fac, newdata = newdata)[,"Estimate"]
647651

‎Lectures/Week 6/lec_12_multivariate_GAMs.html‎

Lines changed: 20 additions & 14 deletions
Large diffs are not rendered by default.

‎docs/Lectures/Week 6/lec_12_multivariate_GAMs.Rmd‎

Lines changed: 21 additions & 17 deletions
Original file line numberDiff line numberDiff line change
@@ -173,7 +173,7 @@ fit_S <- gam(value ~
173173

174174
## Model S fits
175175

176-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
176+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
177177
long_seal$pred_S <- predict(fit_S, long_seal)
178178
179179
ggplot(long_seal, aes(Year, pred_S)) +
@@ -232,7 +232,7 @@ MARSS model where each time series is a different trend
232232
* `y ~ s(t, k = 12)`
233233
* low uncertainty, poor fit
234234

235-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
235+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
236236
library(mgcv)
237237
set.seed(42)
238238
@@ -273,7 +273,7 @@ ggplot(df_pred, aes(x, fit)) +
273273
* `y ~ s(t, k = 50)`
274274
* good fit, high uncertainty
275275

276-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
276+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
277277
set.seed(42)
278278
279279
# Simulate a seasonal time series
@@ -359,7 +359,7 @@ $$ E[Y_{j,t}] = g^{-1} \left( \beta_0 + \sum_{i=1}^{I} s(x_{i,j,t}) + \sum_{m=1}
359359
* `fit <- mvgam(y ~ s(time, k = 50),
360360
data = df, family = "gaussian")`
361361

362-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
362+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
363363
set.seed(42)
364364
365365
# Simulate a seasonal time series
@@ -433,7 +433,7 @@ plot_mvgam_smooth(fit_GS)
433433

434434
## Including dynamic components: coastal upwelling
435435

436-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
436+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
437437
cuti <- read.csv("CUTI_monthly.csv")
438438
439439
library(dplyr)
@@ -456,14 +456,18 @@ cuti_long <- cuti_subset %>%
456456
cuti_long$latitude <- as.factor(cuti_long$latitude)
457457
cuti_long$time <- sort(rep(1:447, 9))
458458
cuti_long$series <- cuti_long$latitude
459+
```
460+
461+
```{r echo=TRUE, message = TRUE, results = 'asis'}
459462
head(cuti_long)
460463
```
461464

465+
462466
## Including dynamic components: coastal upwelling
463467

464468
* for speed, we'll filter the data to be 2023--2024 (train) and use 2025 for testing
465469

466-
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=4, fig.width=7}
470+
```{r echo=FALSE, warning=FALSE, message=FALSE, fig.height=3, fig.width=7}
467471
cuti2324 <- dplyr::filter(cuti_long, year %in% 2023:2024)
468472
ggplot(cuti2324, aes(year+month/12, cuti)) +
469473
geom_line()+
@@ -505,7 +509,7 @@ plot_mvgam_smooth(fit_trend, smooth = "time")
505509

506510
* Predictions (last 3 points are forecasts)
507511

508-
```{r echo=FALSE, fig.height=4, fig.width=7}
512+
```{r echo=FALSE, fig.height=3, fig.width=7}
509513
newdata <- dplyr::filter(cuti_long, year > 2022)
510514
newdata$pred <- predict(fit_trend, newdata = newdata)[,"Estimate"]
511515
@@ -518,7 +522,7 @@ ggplot(newdata, aes(time, pred)) +
518522
## Including dynamic components: coastal upwelling
519523

520524
* Possibly some problematic residuals (autocorrelated)
521-
```{r echo=TRUE, fig.height=4, fig.width=7}
525+
```{r echo=TRUE, fig.height=3, fig.width=7}
522526
plot(fit_trend)
523527
```
524528

@@ -544,29 +548,29 @@ fit_trend_rw <- mvgam(cuti ~ s(time, bs="cr") +
544548
## Including dynamic components: coastal upwelling
545549

546550
* ACF now looks a little better (and QQ plot too)
547-
```{r echo=TRUE, fig.height=4, fig.width=7}
551+
```{r echo=TRUE, fig.height=3, fig.width=7}
548552
plot(fit_trend_rw)
549553
```
550554

551555
## Including dynamic components: coastal upwelling
552556

553-
```{r echo=TRUE, fig.height=4, fig.width=7}
557+
```{r echo=TRUE, fig.height=3, fig.width=7}
554558
plot_mvgam_smooth(fit_trend_rw)
555559
```
556560

557561
## Including dynamic components: coastal upwelling
558562

559563
* Each of the trends (here for Latitude 31) is ~ flat
560564

561-
```{r echo=TRUE, fig.height=4, fig.width=7}
565+
```{r echo=TRUE, fig.height=3, fig.width=7}
562566
plot_mvgam_trend(fit_trend_rw)
563567
```
564568

565569
## Including dynamic components: coastal upwelling
566570

567571
* Some fits to forecasts are better -- but training data actually seems worse
568572

569-
```{r echo=FALSE, fig.height=4, fig.width=7}
573+
```{r echo=FALSE, fig.height=3, fig.width=7}
570574
newdata <- dplyr::filter(cuti_long, year > 2022)
571575
newdata$pred <- predict(fit_trend_rw, newdata = newdata)[,"Estimate"]
572576
@@ -582,7 +586,7 @@ ggplot(newdata, aes(time, pred)) +
582586

583587
* Let's drop out the main smooth now and let that be explained by the trend
584588

585-
```{r echo=TRUE, fig.height=4, fig.width=7}
589+
```{r echo=TRUE, fig.height=3, fig.width=7}
586590
fit_trend_ar1 <- mvgam(cuti ~ s(time, bs="cr") +
587591
s(latitude, k=9, bs="re") +
588592
s(month, k = 12, bs="cc"),
@@ -595,7 +599,7 @@ fit_trend_ar1 <- mvgam(cuti ~ s(time, bs="cr") +
595599
## Including dynamic components: coastal upwelling
596600

597601
* Trends more flexible than the random walk appear to fit data a little better
598-
```{r echo=FALSE, fig.height=4, fig.width=7}
602+
```{r echo=FALSE, fig.height=3, fig.width=7}
599603
newdata <- dplyr::filter(cuti_long, year > 2022)
600604
newdata$pred <- predict(fit_trend_ar1, newdata = newdata)[,"Estimate"]
601605
@@ -611,7 +615,7 @@ ggplot(newdata, aes(time, pred)) +
611615

612616
* Instead, model trends as dynamic factors (n = 2)
613617

614-
```{r echo=TRUE, fig.height=4, fig.width=7}
618+
```{r echo=TRUE, fig.height=3, fig.width=6.5}
615619
fit_trend_fac <- mvgam(cuti ~ s(time, bs="cr") +
616620
s(latitude, k=9, bs="re") +
617621
s(month, k = 12, bs="cc"),
@@ -631,7 +635,7 @@ fit_trend_fac <- mvgam(cuti ~ s(time, bs="cr") +
631635
* By default mvdlm models these as AR1
632636

633637
* DFA model with AR1 trend + seasonal smooth
634-
```{r echo=TRUE, fig.height=4, fig.width=7}
638+
```{r echo=TRUE, fig.height=3, fig.width=6.5}
635639
plot(fit_trend_fac, type = 'factors')
636640
```
637641

@@ -641,7 +645,7 @@ plot(fit_trend_fac, type = 'factors')
641645

642646
* Factors for this application probably not flexible enough
643647

644-
```{r echo=FALSE, fig.height=4, fig.width=7}
648+
```{r echo=FALSE, fig.height=3, fig.width=6.5}
645649
newdata <- dplyr::filter(cuti_long, year > 2022)
646650
newdata$pred <- predict(fit_trend_fac, newdata = newdata)[,"Estimate"]
647651

‎docs/Lectures/Week 6/lec_12_multivariate_GAMs.html‎

Lines changed: 20 additions & 14 deletions
Large diffs are not rendered by default.

‎docs/index.html‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -349,7 +349,7 @@ <h2>
349349
</center>
350350
<p><br></p>
351351
<center>
352-
<em>This site was last updated at 06:57 on 08 May 2025</em>
352+
<em>This site was last updated at 09:00 on 08 May 2025</em>
353353
</center>
354354

355355

0 commit comments

Comments
 (0)