pm25_korea <- read.csv(
"./data/첨부1_2025_PM2.5 시간자료.txt",
fileEncoding = "CP949", sep = "\t"
)
airkorea_station <- read_xlsx("./data/첨부2_측정소_인근_대기관측소(AWS).xlsx")[, c(1:3, 5)]② 전국 단기 환경기준 달성 측정소 추출
532개 도시대기 측정소 중 일평균 99백분위 ≤ 35 ㎍/㎥ 달성 지역
1 개요
전국 도시대기 측정소 전체를 순회하며, 2025년 일평균 PM2.5의 99백분위수가 단기 환경기준(35 ㎍/㎥) 이하인 측정소를 추출한다.
노트전체 532개 도시대기 측정소 중 몇 %가 될까?
내가 사는 동작구가 단기 환경기준 미달인 것을 보면, 환경부에서 정한 환경기준이 상당히 엄격한 기준인 것으로 판단된다. 이 기준을 달성 가능한 도시대기 측정소가 많지는 않을 것이라 생각되었다.
2 데이터 로드
데이터는 동작구 분석 시의 데이터와 동일하고, 필요한 컬럼만 뽑아서 사용했다.
3 전국 측정소 단기기준 달성 스캔
측정소별로 일평균을 구한 뒤 오름차순 99백분위 값이 35 ㎍/㎥ 이하이면 ’달성’으로 판정한다. 단, 환경기준 판정에는 아래 자료 유효성 기준을 함께 적용한다.
- 일평균 유효성: 하루 24시간 중 18시간(75%) 이상 측정된 날만 유효 일평균으로 인정
- 연간 유효성: 연간 유효일이 274일(75%) 이상인 측정소만 판정 대상으로 인정하고, 미달 측정소는 ’자료부족’으로 판정에서 제외
자료가 부족한 측정소를 그대로 판정에 넣으면 관측이 빠진 기간이 통계를 왜곡할 수 있어, 유효자료 요건을 갖춘 측정소끼리만 비교한다.
dir.create("outputs", showWarnings = FALSE)
airkorea_stn <- unique(pm25_korea$AREA_INDEX)
pass_stn <- c() # 단기 기준 달성 측정소
insufficient_stn <- c() # 연간 유효자료 부족(<274일) → 판정 제외
for (i in airkorea_stn) {
pm25 <- pm25_korea %>% filter(AREA_INDEX == i)
pm25$date <- substr(pm25$YYYYMMDDHH, 1, 8)
pm25_daily <- daily_valid_mean(pm25) # ① 하루 18시간 이상 측정된 유효 일평균만
# ② 연간 유효일이 274일(75%) 미만인 측정소는 통계적으로 유효하지 않아 판정 제외
if (nrow(pm25_daily) < ANNUAL_MIN_DAYS) {
insufficient_stn <- c(insufficient_stn, i)
next
}
n_valid <- nrow(pm25_daily)
ascending_order <- order(pm25_daily$daily_mean)
critical_value <- pm25_daily[ascending_order[p99_index(n_valid)], ]
if (critical_value$daily_mean > PM25_DAILY_STD) {
# 99백분위 값이 35 초과 → 달성 아님
} else {
pass_stn <- c(pass_stn, i)
g1 <- pm25_daily %>%
mutate(DateTime = ymd(as.character(date))) %>%
ggplot() +
geom_line(aes(x = DateTime, y = daily_mean), lwd = 1.0) +
labs(title = paste0("(", i, ") 측정소 일평균 PM2.5 농도 변화 추이"),
y = "PM2.5 농도(㎍/㎥)", x = "날짜") +
geom_hline(yintercept = PM25_DAILY_STD, col = "red")
ggsave(file.path("outputs", paste0(i, "_측정소_일평균_pm25_농도_시계열_추이.png")),
g1, width = 12, height = 7, units = "in", dpi = 150)
g2 <- pm25_daily %>%
mutate(DateTime = ymd(as.character(date))) %>%
ggplot() +
geom_bar(aes(x = reorder(DateTime, daily_mean), y = daily_mean),
lwd = 1.0, stat = "identity") +
scale_y_continuous(limits = c(0, 100), breaks = seq(0, 100, 20)) +
labs(title = paste0("(", i, ") 측정소 일평균 PM2.5 (오름차순)"),
y = "PM2.5 농도(㎍/㎥)", x = "날짜") +
theme(axis.text.x = element_blank()) +
geom_hline(yintercept = PM25_DAILY_STD, col = "blue") +
geom_vline(xintercept = p99_index(n_valid), col = "red")
ggsave(file.path("outputs", paste0(i, "_측정소_일평균_pm25_농도_오름차순_추이.png")),
g2, width = 12, height = 7, units = "in", dpi = 150)
}
}
length(airkorea_stn) # 전체 측정소[1] 532length(insufficient_stn) # 자료부족으로 판정 제외[1] 15length(airkorea_stn) - length(insufficient_stn) # 실제 판정 대상[1] 517length(pass_stn) # 단기 기준 달성[1] 32pass_stn [1] 238183 238203 238212 336111 336112 336127 336133 336135 336221 336357
[11] 336452 336471 336481 336491 336501 336511 336541 336551 336561 336571
[21] 336572 339122 339124 422212 437113 437421 437431 437542 632131 632151
[31] 632441 632451전체 532개 도시대기 측정소 중 연간 유효자료 기준을 충족한 517개를 대상으로 판정한 결과(15개는 자료부족으로 제외), 32개 측정소가 단기 환경기준을 달성한 것으로 나온다. 그런데 에어코리아 측정소 번호로는 어딘지 감이 안 온다.
4 달성 측정소 목록 (지역명 매핑)
그래서 간략하게 에어코리아 측정소가 위치한 지역을 찾아보았다.
stn_name <- data.frame()
for (j in pass_stn) {
stn <- airkorea_station %>% filter(station_code == j)
stn_name <- rbind(stn_name, stn)
}
stn_name# A tibble: 32 × 4
location station_name station_code mang_name
<chr> <chr> <chr> <chr>
1 경남 장유동 238183 도시대기
2 경남 고현동 238203 도시대기
3 경남 향촌동 238212 도시대기
4 전남 용당동 336111 도시대기
5 전남 부흥동 336112 도시대기
6 전남 여천동(여수) 336127 도시대기
7 전남 순천만 336133 도시대기
8 전남 신대 336135 도시대기
9 전남 화양면 336221 도시대기
10 전남 봉강면 336357 도시대기
# ℹ 22 more rows위 표는 앞의 10개 측정소만 보여주고 나머지 22개는 생략되므로, 32개 달성 측정소가 어느 지역에 주로 위치하는지 파악할 수 있도록 시·도별로 묶어 전체 목록을 정리했다.
stn_name %>%
group_by(location) %>%
summarise(n = n(),
지역명 = paste(station_name, collapse = ", "),
.groups = "drop") %>%
arrange(desc(n), location) %>%
mutate(line = sprintf("- **%s** (%d개소): %s", location, n, 지역명)) %>%
pull(line) %>%
cat(sep = "\n")- 전남 (18개소): 용당동, 부흥동, 여천동(여수), 순천만, 신대, 화양면, 봉강면, 화순읍, 장흥읍, 진도읍, 신지면, 함평읍, 고흥읍, 강진읍, 곡성읍, 구례읍, 벌교읍, 보성읍
- 강원 (4개소): 주문진읍, 천곡동, 상리, 양양읍
- 경북 (4개소): 장량동, 울진군, 석포면, 울릉읍
- 경남 (3개소): 장유동, 고현동, 향촌동
- 제주 (2개소): 성산읍, 강정동
- 대구 (1개소): 군위읍
32개 중 18개소가 전남 지역에 위치하고, 나머지도 강원·경북·경남·제주 등 비수도권 청정 지역이다. 아쉽게도 수도권(서울·인천·경기) 측정소는 한 곳도 포함되지 않았으며, 광역시 중에는 대구 군위읍이 유일하게 단기환경기준을 달성했다.
5 사례: 대구 군위읍 (422212)
이전의 동작구 분석 시에 만든 시각화 기법을 활용해서, 달성 측정소 중 유일하게 광역시로 구분되는 측정소인 대구 군위읍의 초미세먼지 농도 추이와 단기 환경기준 달성 여부 판별 시각화를 동시에 수행해 보았다.
동작구와 비교해 보면, 일별 환경기준(35 ㎍/㎥)을 초과하는 날 자체가 적은 것을 알 수 있다.
# 군위읍(422212) 유효 일평균 — 스캔과 동일한 유효성 기준 적용
gunwi_daily <- pm25_korea %>%
filter(AREA_INDEX == 422212) %>%
mutate(date = substr(YYYYMMDDHH, 1, 8)) %>%
daily_valid_mean() %>%
mutate(DateTime = ymd(date))
gunwi_n <- nrow(gunwi_daily)
plot3_1 <- gunwi_daily %>%
ggplot() +
geom_line(aes(x = DateTime, y = daily_mean), lwd = 1.0) +
geom_vline(xintercept = ymd("20250313"), col = "red") +
labs(title = "대구 군위읍 일평균 PM2.5 농도 변화 추이", y = "PM2.5 농도(㎍/㎥)", x = "날짜") +
geom_hline(yintercept = PM25_DAILY_STD, col = "blue") +
annotate("text", x = ymd("20250313"), y = 60,
label = "고농도 황사유입일(3/13)", col = "red", hjust = 0)
plot3_2 <- gunwi_daily %>%
ggplot() +
geom_bar(aes(x = reorder(DateTime, daily_mean), y = daily_mean),
lwd = 1.0, stat = "identity") +
scale_y_continuous(limits = c(0, 100), breaks = seq(0, 100, 20)) +
labs(title = "대구 군위읍 일평균 PM2.5 (오름차순 정렬)", y = "PM2.5 농도(㎍/㎥)", x = "날짜") +
theme(axis.text.x = element_blank()) +
geom_hline(yintercept = PM25_DAILY_STD, col = "blue") +
geom_vline(xintercept = p99_index(gunwi_n), col = "red")
plot3_1 / plot3_2
6 분석 결과
- 전국 도시대기 측정소 중 단기 환경기준(99백분위 일평균 ≤ 35) 달성 측정소는 모두 32개로, 연간 유효자료 기준을 충족한 517개 측정소의 약 6.2%(32/517) 수준으로 나타났다. (자료부족으로 판정에서 제외된 측정소는 15개)
- 당연한 이야기겠지만, 도시 지역(예: 동작구) 대비, 청정 지역 측정소(예: 강원 양양읍, 대구 군위읍 등)에서 단기 환경기준을 달성하는 경향이 보였다.