knitr::opts_chunk$set(echo = TRUE)
library(readxl)
## Warning: package 'readxl' was built under R version 4.5.3
library(sf)
## Linking to GEOS 3.13.1, GDAL 3.11.0, PROJ 9.6.0; sf_use_s2() is TRUE
library(dplyr)
## Warning: package 'dplyr' was built under R version 4.5.3
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.5.3
library(viridis)
## Loading required package: viridisLite
library(geodata)
## Warning: package 'geodata' was built under R version 4.5.3
## Loading required package: terra
## Warning: package 'terra' was built under R version 4.5.3
## terra 1.9.34
# =========================================================
# USDA 2022 CATTLE PRODUCTION MAP
# =========================================================
# Load USDA livestock data
livestock <- read_excel(
"data.xlsx",
sheet = "Livestock and Animals"
)
# Load county-name lookup
county_names <- read_excel(
"data.xlsx",
sheet = "County Names"
)
# ---------------------------------------------------------
# California county boundaries
# ---------------------------------------------------------
us_counties <- geodata::gadm(
country = "USA",
level = 2,
path = "data"
)
california <- st_as_sf(us_counties)
california <- california[
california$NAME_1 == "California",
]
# ---------------------------------------------------------
# Find the county-name and FIPS columns
# ---------------------------------------------------------
fips_column <- names(county_names)[
grepl("fips", names(county_names), ignore.case = TRUE)
][1]
name_column <- names(county_names)[
grepl("county|name", names(county_names), ignore.case = TRUE)
][1]
# Create clean lookup table
county_lookup <- county_names %>%
select(
FIPS = all_of(fips_column),
County = all_of(name_column)
) %>%
mutate(
FIPS = as.character(FIPS),
FIPS = gsub("\\D", "", FIPS),
FIPS = sprintf("%05d", as.numeric(FIPS)),
County = trimws(County),
County = gsub(
"\\s+County$",
"",
County,
ignore.case = TRUE
)
) %>%
distinct(FIPS, .keep_all = TRUE)
# ---------------------------------------------------------
# Extract 2022 cow/heifer inventory
# ---------------------------------------------------------
cattle_2022 <- livestock %>%
mutate(
FIPS = as.character(FIPSTEXT),
FIPS = gsub("\\D", "", FIPS),
FIPS = sprintf("%05d", as.numeric(FIPS))
) %>%
select(
FIPS,
cattle = y22_M099_valueNumeric
) %>%
left_join(
county_lookup,
by = "FIPS"
) %>%
filter(
!is.na(County),
substr(FIPS, 1, 2) == "06"
) %>%
select(
County,
cattle
)
# ---------------------------------------------------------
# Join cattle data to California counties
# ---------------------------------------------------------
cattle_map <- california %>%
left_join(
cattle_2022,
by = c("NAME_2" = "County")
)
# ---------------------------------------------------------
# Check data match
# ---------------------------------------------------------
print(
cattle_map %>%
st_drop_geometry() %>%
summarise(
total_counties = n(),
counties_with_data = sum(!is.na(cattle)),
counties_without_data = sum(is.na(cattle))
)
)
## total_counties counties_with_data counties_without_data
## 1 58 55 3
# ---------------------------------------------------------
# Create dark California cattle map
# ---------------------------------------------------------
ggplot(cattle_map) +
geom_sf(
aes(fill = cattle),
color = "#34383A",
linewidth = 0.25
) +
scale_fill_viridis_c(
option = "magma",
trans = "sqrt",
na.value = "#151819",
name = "Cows & heifers\ninventory"
) +
labs(
title = "Where Cattle Production Is Concentrated",
subtitle = "Cows and heifers that had calved — California, 2022",
caption = "Source: USDA National Agricultural Statistics Service, 2022 Census of Agriculture."
) +
theme_void() +
theme(
plot.background = element_rect(
fill = "#0B0D0E",
color = NA
),
panel.background = element_rect(
fill = "#0B0D0E",
color = NA
),
plot.title = element_text(
color = "#F2EFE9",
size = 21,
face = "bold"
),
plot.subtitle = element_text(
color = "#8E9294",
size = 11
),
plot.caption = element_text(
color = "#6F7477",
size = 8,
hjust = 0
),
legend.background = element_rect(
fill = "#0B0D0E",
color = NA
),
legend.text = element_text(
color = "#E8E5DF"
),
legend.title = element_text(
color = "#E8E5DF",
face = "bold"
)
)
library(readxl)
farms <- read_excel(
"data.xlsx",
sheet = "Farms"
)
economics <- read_excel(
"data.xlsx",
sheet = "Economics"
)
producers <- read_excel(
"data.xlsx",
sheet = "Producers"
)
cat("FARMS:\n")
## FARMS:
print(names(farms))
## [1] "FIPS" "FIPSTEXT" "y22_M001_valueNumeric"
## [4] "y22_M001_valueText" "y22_M001_classRange" "y22_M002_valueNumeric"
## [7] "y22_M002_valueText" "y22_M002_classRange" "y22_M050_valueNumeric"
## [10] "y22_M050_valueText" "y22_M050_classRange" "y22_M051_valueNumeric"
## [13] "y22_M051_valueText" "y22_M051_classRange" "y22_M052_valueNumeric"
## [16] "y22_M052_valueText" "y22_M052_classRange" "y22_M053_valueNumeric"
## [19] "y22_M053_valueText" "y22_M053_classRange" "y22_M054_valueNumeric"
## [22] "y22_M054_valueText" "y22_M054_classRange" "y22_M055_valueNumeric"
## [25] "y22_M055_valueText" "y22_M055_classRange" "y22_M056_valueNumeric"
## [28] "y22_M056_valueText" "y22_M056_classRange" "y22_M057_valueNumeric"
## [31] "y22_M057_valueText" "y22_M057_classRange" "y22_M058_valueNumeric"
## [34] "y22_M058_valueText" "y22_M058_classRange" "y22_M059_valueNumeric"
## [37] "y22_M059_valueText" "y22_M059_classRange" "y22_M064_valueNumeric"
## [40] "y22_M064_valueText" "y22_M064_classRange" "y22_M065_valueNumeric"
## [43] "y22_M065_valueText" "y22_M065_classRange" "y22_M066_valueNumeric"
## [46] "y22_M066_valueText" "y22_M066_classRange" "y22_M067_valueNumeric"
## [49] "y22_M067_valueText" "y22_M067_classRange" "y22_M068_valueNumeric"
## [52] "y22_M068_valueText" "y22_M068_classRange" "y22_M069_valueNumeric"
## [55] "y22_M069_valueText" "y22_M069_classRange" "y22_M070_valueNumeric"
## [58] "y22_M070_valueText" "y22_M070_classRange" "y22_M071_valueNumeric"
## [61] "y22_M071_valueText" "y22_M071_classRange" "y22_M079_valueNumeric"
## [64] "y22_M079_valueText" "y22_M079_classRange" "y22_M080_valueNumeric"
## [67] "y22_M080_valueText" "y22_M080_classRange" "y22_M081_valueNumeric"
## [70] "y22_M081_valueText" "y22_M081_classRange" "y22_M174_valueNumeric"
## [73] "y22_M174_valueText" "y22_M174_classRange"
cat("\nECONOMICS:\n")
##
## ECONOMICS:
print(names(economics))
## [1] "FIPS" "FIPSTEXT" "y22_M003_valueNumeric"
## [4] "y22_M003_valueText" "y22_M003_classRange" "y22_M004_valueNumeric"
## [7] "y22_M004_valueText" "y22_M004_classRange" "y22_M005_valueNumeric"
## [10] "y22_M005_valueText" "y22_M005_classRange" "y22_M006_valueNumeric"
## [13] "y22_M006_valueText" "y22_M006_classRange" "y22_M007_valueNumeric"
## [16] "y22_M007_valueText" "y22_M007_classRange" "y22_M008_valueNumeric"
## [19] "y22_M008_valueText" "y22_M008_classRange" "y22_M009_valueNumeric"
## [22] "y22_M009_valueText" "y22_M009_classRange" "y22_M010_valueNumeric"
## [25] "y22_M010_valueText" "y22_M010_classRange" "y22_M011_valueNumeric"
## [28] "y22_M011_valueText" "y22_M011_classRange" "y22_M012_valueNumeric"
## [31] "y22_M012_valueText" "y22_M012_classRange" "y22_M013_valueNumeric"
## [34] "y22_M013_valueText" "y22_M013_classRange" "y22_M014_valueNumeric"
## [37] "y22_M014_valueText" "y22_M014_classRange" "y22_M015_valueNumeric"
## [40] "y22_M015_valueText" "y22_M015_classRange" "y22_M016_valueNumeric"
## [43] "y22_M016_valueText" "y22_M016_classRange" "y22_M017_valueNumeric"
## [46] "y22_M017_valueText" "y22_M017_classRange" "y22_M018_valueNumeric"
## [49] "y22_M018_valueText" "y22_M018_classRange" "y22_M019_valueNumeric"
## [52] "y22_M019_valueText" "y22_M019_classRange" "y22_M020_valueNumeric"
## [55] "y22_M020_valueText" "y22_M020_classRange" "y22_M021_valueNumeric"
## [58] "y22_M021_valueText" "y22_M021_classRange" "y22_M022_valueNumeric"
## [61] "y22_M022_valueText" "y22_M022_classRange" "y22_M023_valueNumeric"
## [64] "y22_M023_valueText" "y22_M023_classRange" "y22_M024_valueNumeric"
## [67] "y22_M024_valueText" "y22_M024_classRange" "y22_M025_valueNumeric"
## [70] "y22_M025_valueText" "y22_M025_classRange" "y22_M026_valueNumeric"
## [73] "y22_M026_valueText" "y22_M026_classRange" "y22_M027_valueNumeric"
## [76] "y22_M027_valueText" "y22_M027_classRange" "y22_M028_valueNumeric"
## [79] "y22_M028_valueText" "y22_M028_classRange" "y22_M029_valueNumeric"
## [82] "y22_M029_valueText" "y22_M029_classRange" "y22_M030_valueNumeric"
## [85] "y22_M030_valueText" "y22_M030_classRange" "y22_M031_valueNumeric"
## [88] "y22_M031_valueText" "y22_M031_classRange" "y22_M032_valueNumeric"
## [91] "y22_M032_valueText" "y22_M032_classRange" "y22_M033_valueNumeric"
## [94] "y22_M033_valueText" "y22_M033_classRange" "y22_M034_valueNumeric"
## [97] "y22_M034_valueText" "y22_M034_classRange" "y22_M035_valueNumeric"
## [100] "y22_M035_valueText" "y22_M035_classRange" "y22_M036_valueNumeric"
## [103] "y22_M036_valueText" "y22_M036_classRange" "y22_M037_valueNumeric"
## [106] "y22_M037_valueText" "y22_M037_classRange" "y22_M038_valueNumeric"
## [109] "y22_M038_valueText" "y22_M038_classRange" "y22_M039_valueNumeric"
## [112] "y22_M039_valueText" "y22_M039_classRange" "y22_M040_valueNumeric"
## [115] "y22_M040_valueText" "y22_M040_classRange" "y22_M041_valueNumeric"
## [118] "y22_M041_valueText" "y22_M041_classRange" "y22_M042_valueNumeric"
## [121] "y22_M042_valueText" "y22_M042_classRange" "y22_M043_valueNumeric"
## [124] "y22_M043_valueText" "y22_M043_classRange" "y22_M044_valueNumeric"
## [127] "y22_M044_valueText" "y22_M044_classRange" "y22_M045_valueNumeric"
## [130] "y22_M045_valueText" "y22_M045_classRange" "y22_M046_valueNumeric"
## [133] "y22_M046_valueText" "y22_M046_classRange" "y22_M047_valueNumeric"
## [136] "y22_M047_valueText" "y22_M047_classRange" "y22_M060_valueNumeric"
## [139] "y22_M060_valueText" "y22_M060_classRange" "y22_M061_valueNumeric"
## [142] "y22_M061_valueText" "y22_M061_classRange" "y22_M062_valueNumeric"
## [145] "y22_M062_valueText" "y22_M062_classRange" "y22_M063_valueNumeric"
## [148] "y22_M063_valueText" "y22_M063_classRange" "y22_M175_valueNumeric"
## [151] "y22_M175_valueText" "y22_M175_classRange"
cat("\nPRODUCERS:\n")
##
## PRODUCERS:
print(names(producers))
## [1] "FIPS" "FIPSTEXT" "y22_M048_valueNumeric"
## [4] "y22_M048_valueText" "y22_M048_classRange" "y22_M049_valueNumeric"
## [7] "y22_M049_valueText" "y22_M049_classRange" "y22_M072_valueNumeric"
## [10] "y22_M072_valueText" "y22_M072_classRange" "y22_M073_valueNumeric"
## [13] "y22_M073_valueText" "y22_M073_classRange" "y22_M074_valueNumeric"
## [16] "y22_M074_valueText" "y22_M074_classRange" "y22_M075_valueNumeric"
## [19] "y22_M075_valueText" "y22_M075_classRange" "y22_M076_valueNumeric"
## [22] "y22_M076_valueText" "y22_M076_classRange" "y22_M077_valueNumeric"
## [25] "y22_M077_valueText" "y22_M077_classRange" "y22_M078_valueNumeric"
## [28] "y22_M078_valueText" "y22_M078_classRange" "y22_M082_valueNumeric"
## [31] "y22_M082_valueText" "y22_M082_classRange" "y22_M083_valueNumeric"
## [34] "y22_M083_valueText" "y22_M083_classRange" "y22_M084_valueNumeric"
## [37] "y22_M084_valueText" "y22_M084_classRange" "y22_M085_valueNumeric"
## [40] "y22_M085_valueText" "y22_M085_classRange" "y22_M086_valueNumeric"
## [43] "y22_M086_valueText" "y22_M086_classRange" "y22_M087_valueNumeric"
## [46] "y22_M087_valueText" "y22_M087_classRange" "y22_M088_valueNumeric"
## [49] "y22_M088_valueText" "y22_M088_classRange" "y22_M089_valueNumeric"
## [52] "y22_M089_valueText" "y22_M089_classRange" "y22_M090_valueNumeric"
## [55] "y22_M090_valueText" "y22_M090_classRange" "y22_M091_valueNumeric"
## [58] "y22_M091_valueText" "y22_M091_classRange" "y22_M092_valueNumeric"
## [61] "y22_M092_valueText" "y22_M092_classRange" "y22_M093_valueNumeric"
## [64] "y22_M093_valueText" "y22_M093_classRange" "y22_M094_valueNumeric"
## [67] "y22_M094_valueText" "y22_M094_classRange" "y22_M095_valueNumeric"
## [70] "y22_M095_valueText" "y22_M095_classRange" "y22_M096_valueNumeric"
## [73] "y22_M096_valueText" "y22_M096_classRange" "y22_M097_valueNumeric"
## [76] "y22_M097_valueText" "y22_M097_classRange"
library(readxl)
library(dplyr)
library(plotly)
##
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
##
## last_plot
## The following object is masked from 'package:stats':
##
## filter
## The following object is masked from 'package:graphics':
##
## layout
# =========================================================
# CALIFORNIA'S LIVESTOCK LANDSCAPE
# USDA 2022 CENSUS OF AGRICULTURE
# =========================================================
# Load livestock data
livestock <- read_excel(
"data.xlsx",
sheet = "Livestock and Animals"
)
# Load county-name lookup
county_names <- read_excel(
"data.xlsx",
sheet = "County Names"
)
# ---------------------------------------------------------
# Create county lookup
# ---------------------------------------------------------
fips_col <- names(county_names)[
grepl("fips", names(county_names), ignore.case = TRUE)
][1]
county_col <- names(county_names)[
grepl("county|name", names(county_names), ignore.case = TRUE)
][1]
county_lookup <- county_names %>%
select(
FIPS = all_of(fips_col),
County = all_of(county_col)
) %>%
mutate(
FIPS = gsub("\\D", "", as.character(FIPS)),
FIPS = sprintf("%05d", as.numeric(FIPS)),
County = trimws(County),
County = gsub(
"\\s+County$",
"",
County,
ignore.case = TRUE
)
) %>%
distinct(FIPS, .keep_all = TRUE)
# ---------------------------------------------------------
# Prepare livestock variables
# ---------------------------------------------------------
livestock_3d <- livestock %>%
mutate(
FIPS = gsub("\\D", "", as.character(FIPSTEXT)),
FIPS = sprintf("%05d", as.numeric(FIPS))
) %>%
left_join(
county_lookup,
by = "FIPS"
) %>%
filter(
substr(FIPS, 1, 2) == "06",
!is.na(County)
) %>%
transmute(
County = County,
cattle_density =
y22_M098_valueNumeric,
cattle_inventory =
y22_M099_valueNumeric,
cattle_sold =
y22_M101_valueNumeric,
dairy_share =
y22_M100_valueNumeric
) %>%
filter(
!is.na(cattle_density),
!is.na(cattle_inventory),
!is.na(cattle_sold)
)
# ---------------------------------------------------------
# Create hover information
# ---------------------------------------------------------
livestock_3d <- livestock_3d %>%
mutate(
hover_text = paste0(
"<b>", County, "</b><br>",
"Cattle per 100 acres: ",
format(round(cattle_density, 1), big.mark = ","),
"<br>",
"Cows & heifers: ",
format(round(cattle_inventory), big.mark = ","),
"<br>",
"Cattle sold: ",
format(round(cattle_sold), big.mark = ","),
"<br>",
"Cows/heifers as % of cattle: ",
format(round(dairy_share, 1)),
"%"
)
)
# ---------------------------------------------------------
# Create interactive 3D visualization
# ---------------------------------------------------------
livestock_3d_plot <- plot_ly(
data = livestock_3d,
x = ~cattle_density,
y = ~cattle_inventory,
z = ~cattle_sold,
type = "scatter3d",
mode = "markers",
marker = list(
size = 7,
color = ~dairy_share,
colorscale = "Magma",
opacity = 0.9,
showscale = TRUE,
colorbar = list(
title = "Cows & heifers<br>% of cattle",
titlefont = list(
color = "#E8E5DF"
),
tickfont = list(
color = "#E8E5DF"
)
)
),
text = ~hover_text,
hoverinfo = "text"
) %>%
layout(
# -----------------------------------------------------
# TITLE
# -----------------------------------------------------
title = list(
text = "California's Livestock Landscape",
font = list(
color = "#F2EFE9",
size = 24
),
x = 0.03,
xanchor = "left",
y = 0.98,
yanchor = "top"
),
# -----------------------------------------------------
# EXTRA SPACE ABOVE THE GRAPH
# -----------------------------------------------------
margin = list(
t = 110,
r = 40,
b = 40,
l = 40
),
paper_bgcolor = "#0B0D0E",
# -----------------------------------------------------
# 3D SCENE
# -----------------------------------------------------
scene = list(
bgcolor = "#0B0D0E",
xaxis = list(
title = "Cattle per 100 acres",
color = "#E8E5DF",
gridcolor = "#303438",
zerolinecolor = "#555555",
titlefont = list(
color = "#E8E5DF"
)
),
yaxis = list(
title = "Cows & heifers inventory",
color = "#E8E5DF",
gridcolor = "#303438",
zerolinecolor = "#555555",
titlefont = list(
color = "#E8E5DF"
)
),
zaxis = list(
title = "Cattle sold",
color = "#E8E5DF",
gridcolor = "#303438",
zerolinecolor = "#555555",
titlefont = list(
color = "#E8E5DF"
)
),
camera = list(
eye = list(
x = 1.5,
y = 1.5,
z = 1.2
)
)
)
)
# Display the visualization
livestock_3d_plot