Worksheet 06 — Real Climate Data in R

Downloading, summarizing, and modeling Duluth weather station data

r-basics

project-setup

Summarizing data and using an API and more R

Author

Bill Perry

Published

September 10, 2026

Downloading and Summarizing Climate Data

Recap from Worksheet 05

Loaded the paper calibration data and made a scatter plot
Fit a linear regression with lm(area_cm2 ~ mass_g, data = paper_df)
Read every line of summary() — slope, intercept, R², p-value
Checked assumptions — residual vs. fitted plot, QQ plot, Shapiro-Wilk on residuals
Predicted leaf area from tracing mass using predict() with confidence and prediction intervals

Today’s Objectives

Download real NOAA weather data directly into R using GSODR
Explore the raw daily data and see why it’s hard to read
Summarize the data by month and by year using group_by() + summarize()
Plot the raw data and each summary level side by side
Fit a regression (lm()) to estimate the rate of temperature change over time
Use case_when() to build a summer vs. winter season variable
Fit separate models for each season and compare warming rates
Interpret a real long-term climate trend in your own words

How to use this worksheet

Work through each part in order. Type the code into a new R script in Positron and run it line by line.

Code blocks marked ▶ Run this should be executed as written.

Blocks marked ✏️ Your turn ask you to write, modify, or interpret.

The Going further section is optional.

🔮 Predict before you run — and type, don’t paste

Before each ▶ Run this block, cover the output and predict what R will print. Then type the code yourself.

Committing to a guess — how many rows survive a group_by(), which sign the slope will have, which season wins — gives you something to be right or wrong about, and being wrong is what you’ll actually remember. Typing the code instead of pasting it also slows you down just enough to catch the errors this dataset is full of: a lowercase temp where NOAA wants TEMP, an == where you needed %in%.

🧩 Chunk 1 — Download & explore raw data (after lecture Chunk 1)

Parts 1–3: download the Duluth record, trim the columns, and plot the raw daily data.

Part 1 · Load libraries and download the data

Libraries

▶ Run this at the top of your script:

# Load packages at the top — always ----------------------
library(tidyverse)   # data manipulation + ggplot2
library(GSODR)        # download NOAA weather station data

Download Duluth’s weather record

We will use the Duluth International Airport station (727450-14913), which has daily records back to 1948.

▶ Run this (this may take a minute — it’s downloading real data):

# Download daily Duluth weather data, 1948-2025 ----------
duluth_df <- get_GSOD(years = 1948:2025,
                       station = "727450-14913")

✏️ Your turn: Run dim(duluth_df) and glimpse(duluth_df). How many rows and columns does the raw download have?

# Write your code here:

Rows:
Columns:
What does one row represent (a year? a month? a day?):

Part 2 · Trim to the columns we need

▶ Run this:

# Keep only the columns we actually need -------------------
dlh_temp_df <- duluth_df %>%
  select(
    STNID, NAME, YEARMODA,
    YEAR, MONTH, DAY,
    TEMP, MXSPD, PRCP, I_SNOW_ICE
  )

head(dlh_temp_df)

✏️ Your turn: Why do we use select() here instead of just working with the full duluth_df? (Hint: how many columns did get_GSOD() actually return?)

Your answer:

Part 3 · Plot the raw daily data

🔮 Predict first: this plots ~28,000 daily points. Before you run it — will a warming trend be visible, or will something else dominate? Write your guess.

▶ Run this:

# Plot every single daily temperature -----------------------
raw_temp_plot <- dlh_temp_df %>%
  ggplot(aes(x = YEARMODA, y = TEMP)) +
  geom_point(alpha = 0.2, size = 0.5) +
  geom_line(alpha = 0.3) +
  labs(
    title = "Duluth Daily Mean Temperature, 1948–2025",
    x     = "Date",
    y     = "Temperature (°C)"
  ) +
  theme_minimal()

raw_temp_plot

✏️ Your turn: Describe what you see. Can you tell whether Duluth is warming over time just by looking at this plot? Why or why not?

What dominates the plot (seasonal cycle / long-term trend / both)?
Can you see a clear warming trend?  Y / N
Why or why not:

💡 Key idea: The seasonal swing between a Duluth summer and a Duluth winter is roughly 30°C. A warming trend of a few tenths of a degree per decade is real, but at this scale it’s buried under swings a hundred times larger — that’s why we have to summarize before we can see it.

🧩 Chunk 2 — Summarize to reveal the trend (after lecture Chunk 2)

Parts 4–5: collapse the daily data to monthly and yearly means.

Part 4 · Summarize by month

▶ Run this:

# Mean temperature per year-month ----------------------------
dlh_month_df <- dlh_temp_df %>%
  group_by(YEAR, MONTH) %>%
  summarize(
    YEARMODA = first(YEARMODA),
    TEMP     = mean(TEMP, na.rm = TRUE)
  )

head(dlh_month_df)

# Plot the monthly means --------------------------------------
month_temp_plot <- dlh_month_df %>%
  ggplot(aes(x = YEARMODA, y = TEMP)) +
  geom_point(alpha = 0.4, size = 0.8) +
  geom_line(alpha = 0.4) +
  geom_smooth(method = "lm", color = "darkblue") +
  labs(
    title = "Duluth Monthly Mean Temperature",
    x = "Date", y = "Temperature (°C)"
  ) +
  theme_minimal()

month_temp_plot

✏️ Your turn: How many rows does dlh_month_df have compared to dlh_temp_df? What pattern is now visible that wasn’t visible in the raw daily plot?

Rows in dlh_temp_df (daily):
Rows in dlh_month_df (monthly):
New pattern now visible:

Part 5 · Summarize by year

▶ Run this:

# Mean temperature per year ------------------------------------
dlh_year_df <- dlh_temp_df %>%
  group_by(YEAR) %>%
  summarize(
    YEARMODA = first(YEARMODA),
    TEMP     = mean(TEMP, na.rm = TRUE)
  )

head(dlh_year_df)

yearly_temp_plot <- dlh_year_df %>%
  ggplot(aes(x = YEAR, y = TEMP)) +
  geom_point(size = 2) +
  geom_line() +
  geom_smooth(method = "lm", color = "firebrick") +
  labs(
    title = "Duluth Annual Mean Temperature, 1948–2025",
    x = "Year", y = "Mean Temp (°C)"
  ) +
  theme_minimal()

yearly_temp_plot

✏️ Your turn: Fill in the table comparing all three levels of summary:

Level      Rows     What it shows clearly      What it hides
--------   ------   -------------------------  -------------------------
Daily                
Monthly              
Yearly

✏️ Your turn: Which summary level would you use if you wanted to know about a single extreme heat wave? Which would you use to study climate change? Explain why they’re different.

Best level for a single heat wave:
Best level for climate change:
Why they're different:

🧩 Chunk 3 — Model the rate of change (after lecture Chunk 3)

Part 6: fit lm(TEMP ~ YEAR) and read the warming slope.

Part 6 · Model the rate of change with `lm()`

🔮 Predict first: will the slope of TEMP ~ YEAR be positive or negative? Roughly how many °C per decade? Guess before you run summary().

▶ Run this:

# Fit a regression: TEMP predicted by YEAR ----------------------
yearly_temp_model <- lm(TEMP ~ YEAR, data = dlh_year_df)

summary(yearly_temp_model)

✏️ Your turn: This is the exact same lm() syntax from Worksheet 05. Fill in what X and Y are in this new context:

Worksheet 05:  X = mass_g       Y = area_cm2
Worksheet 06:  X = _____        Y = _____

▶ Run this:

# Pull the slope and convert to °C per decade --------------------
yearly_intercept <- coef(yearly_temp_model)[1]
yearly_slope     <- coef(yearly_temp_model)[2]

cat("Slope (b)      =", round(yearly_slope, 4), "°C/year\n")
cat("Warming rate   =", round(yearly_slope * 10, 3), "°C per decade\n")

✏️ Your turn: Report the results from your summary() output:

Slope (b) =
p-value for the slope =
R-squared =
Warming rate in °C per decade =
Is the warming trend statistically significant at α = 0.05?  Y / N

📖 Whitlock & Schluter §17.3: the slope’s p-value tests H₀: β = 0. A significant slope means temperature is changing systematically with year — not just by chance.

🧩 Chunk 4 — Seasons: is winter warming faster? (after lecture Chunk 4)

Parts 7–11: build the season variable, model each season, and write the results paragraph.

Part 7 · Build a season variable with `case_when()`

🔮 Predict first: which season category (summer / winter / shoulder) will have the most rows in the count()? Predict before running.

We now ask a sharper question: is winter warming at the same rate as summer?

▶ Run this:

# Label each day as summer, winter, or shoulder season -------------
dlh_season_df <- dlh_temp_df %>%
  mutate(
    season = case_when(
      MONTH %in% c(6, 7, 8)   ~ "summer",
      MONTH %in% c(12, 1, 2)  ~ "winter",
      TRUE                     ~ "shoulder"
    )
  )

# Quick check - did it work? -----------------------------------------
dlh_season_df %>%
  count(season)

✏️ Your turn: How many rows fall into each season category? Does the count roughly make sense given there are 12 months and 3 are assigned to summer, 3 to winter?

n(summer):
n(winter):
n(shoulder):
Does this make sense?  Y / N

✏️ Your turn: In your own words, explain what each line of the case_when() is doing. What does TRUE ~ "shoulder" mean?

What MONTH %in% c(6, 7, 8) ~ "summer" does:
What TRUE ~ "shoulder" does:

Part 8 · Summarize each season by year

▶ Run this:

# Keep only summer and winter, then summarize by year ----------------
season_year_df <- dlh_season_df %>%
  filter(season %in% c("summer", "winter")) %>%
  group_by(YEAR, season) %>%
  summarize(
    TEMP = mean(TEMP, na.rm = TRUE)
  )

head(season_year_df)

✏️ Your turn: What does group_by(YEAR, season) do differently from group_by(YEAR) alone? Why do we need both?

Your answer:

Part 9 · Plot both seasons together

▶ Run this:

# Plot both seasons with separate trend lines --------------------------
season_temp_plot <- season_year_df %>%
  ggplot(aes(x = YEAR, y = TEMP, color = season)) +
  geom_point(size = 1.8, alpha = 0.7) +
  geom_smooth(method = "lm", se = TRUE) +
  scale_color_manual(values = c(
    "summer" = "firebrick",
    "winter" = "darkblue"
  )) +
  labs(
    title = "Duluth Summer vs. Winter Temperature Trends",
    x = "Year", y = "Mean Temp (°C)",
    color = "Season"
  ) +
  theme_minimal()

season_temp_plot

✏️ Your turn: Just from looking at the plot, which trend line looks steeper — summer or winter? Make a prediction before you run the models in the next part.

Visual prediction — steeper line:

Part 10 · Fit a separate model for each season

🔮 Predict first: commit now — which season warms faster, summer or winter? Write your guess before the two slopes print.

▶ Run this:

# Summer-only model -------------------------------------------------
summer_df <- season_year_df %>% filter(season == "summer")
summer_model <- lm(TEMP ~ YEAR, data = summer_df)
summer_slope <- coef(summer_model)[2]
summary(summer_model)

# Winter-only model -------------------------------------------------
winter_df <- season_year_df %>% filter(season == "winter")
winter_model <- lm(TEMP ~ YEAR, data = winter_df)
winter_slope <- coef(winter_model)[2]
summary(winter_model)

# Compare the two slopes ----------------------------------------------
cat("Summer warming rate:", round(summer_slope * 10, 3), "°C/decade\n")
cat("Winter warming rate:", round(winter_slope * 10, 3), "°C/decade\n")

✏️ Your turn: Fill in the comparison table:

Season    Slope (°C/year)   Rate (°C/decade)   p-value   Significant? (Y/N)
-------   ----------------  ------------------  --------  -------------------
Summer
Winter

✏️ Your turn: Was your visual prediction from Part 9 correct? Which season is warming faster in Duluth’s data? Does this match what you might expect from the news or your own experience with winters here?

Faster-warming season:
Did this match your prediction?  Y / N
Does this match your real-world experience?  Y / N
Why might winter and summer warm at different rates? (Hint: think about ice cover, snow reflectivity, cloud cover)

Part 11 · Write a short results paragraph

✏️ Your turn: Write a short paragraph (4–6 sentences) summarizing what you found. Include:

The overall annual warming rate (°C/decade) and whether it was significant
The summer warming rate and the winter warming rate, compared
One sentence connecting this to something you’ve personally noticed about Duluth winters (optional, but encouraged)

Write your results paragraph here:

Part 12 · Review and checkpoint

At this point you should be able to:

Download real station data into R with get_GSOD()
Explain why raw daily data can hide long-term trends
Use group_by() %>% summarize() to collapse data to monthly or yearly means
Fit lm(Y ~ YEAR, data) and extract the slope as a rate of change
Convert a yearly slope into a more interpretable “per decade” rate
Use case_when() to build a new multi-condition categorical variable
Filter to specific categories with filter(x %in% c(...))
Fit and compare two group-specific regression models

✏️ Your turn — before you move on: Run your entire script with Ctrl/Cmd + Shift + Enter. Does it run from top to bottom without errors?

Ran cleanly?  Y / N
If not, what error appeared:

Part 13 · Going further

This section is optional — work through it if you finish early or want to push deeper.

Try a different season split

▶ Try this: Instead of strict 3-month summer/winter, try defining seasons using meteorological vs. astronomical boundaries, or add a “spring” and “fall” category to the case_when().

# Add all four seasons instead of just summer/winter -----------------
dlh_four_season_df <- dlh_temp_df %>%
  mutate(
    season = case_when(
      MONTH %in% c(3, 4, 5)   ~ "spring",
      MONTH %in% c(6, 7, 8)   ~ "summer",
      MONTH %in% c(9, 10, 11) ~ "fall",
      MONTH %in% c(12, 1, 2)  ~ "winter"
    )
  )

dlh_four_season_df %>% count(season)

✏️ Your turn: Does adding spring and fall change your interpretation at all?

Your answer:

Look at snow instead of temperature

▶ Try this:

# Explore snow/ice flag data -------------------------------------------
dlh_snow_df <- dlh_temp_df %>%
  mutate(
    snow_mm = case_when(
      I_SNOW_ICE == 1 ~ PRCP,
      TRUE             ~ NA
    )
  ) %>%
  mutate(snow_mm = ifelse(snow_mm == 0, NA, snow_mm))

yearly_snow_df <- dlh_snow_df %>%
  group_by(YEAR) %>%
  summarize(sum_snow = sum(snow_mm, na.rm = TRUE))

yearly_snow_df %>%
  ggplot(aes(x = YEAR, y = sum_snow)) +
  geom_point() +
  geom_line() +
  geom_smooth(method = "lm") +
  labs(x = "Year", y = "Total Yearly Snow (mm)") +
  theme_minimal()

✏️ Your turn: Has total yearly snowfall changed over time the same way temperature has? Why might these two trends differ?

Your answer:

What your `figures/` folder should contain after this worksheet

figures/
├── raw_temp_plot.png        ← from Part 3
├── month_temp_plot.png      ← from Part 4
├── yearly_temp_plot.png     ← from Part 5
├── season_temp_plot.png     ← from Part 9

Extension — out of class (~30–40 min)

Turn this in with your worksheet. You measured Duluth’s warming rate. Now run the same workflow for a city you pick and compare.

Important

E1 and E3 must be handwritten on paper, photographed, and embedded (![caption](my_photo.jpg)). Typed answers get at most half credit, even if correct — E1 is your prediction before you have the numbers, E3 is your own reasoning after.

E1 · Predict, then check — ✍️ by hand (3 pts) — do this FIRST

Before downloading anything, on paper: name your city and roughly where it sits relative to Duluth (latitude, coastal vs inland). Predict whether its annual warming rate will be faster, slower, or about the same as Duluth’s, and why (one sentence of real reasoning). Photograph the page now.

E2 · Your city’s warming rate (4 pts)

Using the same GSODR workflow as the in-class part: pick a station with a long record, download it, summarize to one mean temperature per year, and fit lm(TEMP ~ YEAR).

# same as class, new station id
# get_GSOD(years = ..., station = "your-station") -> group_by(YEAR) ->
#   summarize(mean_temp = mean(TEMP, na.rm = TRUE)) -> lm(mean_temp ~ YEAR)

Report: your city + station id, years of record, the slope in °C per decade (slope × 10), its p-value, R², and whether your city is warming faster or slower than Duluth.

E3 · Explain it — ✍️ by hand, with YOUR numbers (3 pts)

Get your E1 page back out. Using your actual numbers:

Was your prediction right? Quote your real warming rate against your guess.
A single station’s significant trend does not prove why a city is warming. Name one confound that could inflate one station’s trend with nothing to do with global climate (station moved, instrument change, urban heat island).
Using your R², explain what fraction of year-to-year temperature variation is explained by year alone — and why a low R² can still accompany a real, significant warming trend.

Getting unstuck

get_GSOD() taking forever / failing: check your internet connection — this function reaches out to NOAA’s servers. If it keeps failing, ask for the saved .csv backup.
select() error — object not found: run names(duluth_df) first to confirm exact column names; capitalization matters (TEMP, not temp).
case_when() returns all NA: check that your conditions use %in% for multiple values (MONTH %in% c(6,7,8)), not == with a vector.
group_by() summary looks wrong: always check with head() immediately after summarizing — did you group by the right combination of variables?
Plot legend missing or wrong colors: confirm color = season is inside aes(), and that scale_color_manual() values exactly match your case_when() labels ("summer", "winter").
Cheat sheets — https://posit.co/resources/cheatsheets/

💡 Key idea: Notice you never opened Excel for any of this — download, trim, plot, summarize, and model all happened inside the same R session. That’s the habit worth keeping for your final project, whatever data source you end up pulling from.

End of the Weather Data worksheet.

Other Formats

Worksheet 06 — Real Climate Data in R

Downloading and Summarizing Climate Data

Recap from Worksheet 05

Today’s Objectives

How to use this worksheet

🧩 Chunk 1 — Download & explore raw data (after lecture Chunk 1)

Part 1 · Load libraries and download the data

Libraries

Download Duluth’s weather record

Part 2 · Trim to the columns we need

Part 3 · Plot the raw daily data

🧩 Chunk 2 — Summarize to reveal the trend (after lecture Chunk 2)

Part 4 · Summarize by month

Part 5 · Summarize by year

🧩 Chunk 3 — Model the rate of change (after lecture Chunk 3)

Part 6 · Model the rate of change with `lm()`

🧩 Chunk 4 — Seasons: is winter warming faster? (after lecture Chunk 4)

Part 7 · Build a season variable with `case_when()`

Part 8 · Summarize each season by year

Part 9 · Plot both seasons together

Part 10 · Fit a separate model for each season

Part 11 · Write a short results paragraph

Part 12 · Review and checkpoint

Part 13 · Going further

Try a different season split

Look at snow instead of temperature

What your `figures/` folder should contain after this worksheet

Extension — out of class (~30–40 min)

E1 · Predict, then check — ✍️ by hand (3 pts) — do this FIRST

E2 · Your city’s warming rate (4 pts)

E3 · Explain it — ✍️ by hand, with YOUR numbers (3 pts)

Getting unstuck

Other Formats

Downloading and Summarizing Climate Data

Recap from Worksheet 05

Today’s Objectives

How to use this worksheet

🧩 Chunk 1 — Download & explore raw data (after lecture Chunk 1)

Part 1 · Load libraries and download the data

Libraries

Download Duluth’s weather record

Part 2 · Trim to the columns we need

Part 3 · Plot the raw daily data

🧩 Chunk 2 — Summarize to reveal the trend (after lecture Chunk 2)

Part 4 · Summarize by month

Part 5 · Summarize by year

🧩 Chunk 3 — Model the rate of change (after lecture Chunk 3)

Part 6 · Model the rate of change with lm()

🧩 Chunk 4 — Seasons: is winter warming faster? (after lecture Chunk 4)

Part 7 · Build a season variable with case_when()

Part 8 · Summarize each season by year

Part 9 · Plot both seasons together

Part 10 · Fit a separate model for each season

Part 11 · Write a short results paragraph

Part 12 · Review and checkpoint

Part 13 · Going further

Try a different season split

Look at snow instead of temperature

What your figures/ folder should contain after this worksheet

Extension — out of class (~30–40 min)

E1 · Predict, then check — ✍️ by hand (3 pts) — do this FIRST

E2 · Your city’s warming rate (4 pts)

E3 · Explain it — ✍️ by hand, with YOUR numbers (3 pts)

Getting unstuck

Part 6 · Model the rate of change with `lm()`

Part 7 · Build a season variable with `case_when()`

What your `figures/` folder should contain after this worksheet