한반도 겨울철 사망자에 미치는 한파의 영향 분석

2.1. 자료

본 논문에서는 한반도에서 한파–사망자 영향 연구를 위해 기상청과 국가데이터처의 자료를 통합하여 활용하였다(Table 1). 먼저, 한파와 관련된 기상 요인을 파악하기 위해 기상청 종관관측시스템(Automated Synoptic Observing System, ASOS)의 일자료를 수집하였다. 사용된 변수는 기온(평균·최저·최고), 풍향·풍속, 강수량, 일사량 등이며, 분석 기간은 국가데이터처 자료 확보가 가능한 2000년부터 2023년까지 최근 23년간의 겨울철(12월 ~ 2월)이다. 연구 대상지역은 행정구역 기준 전국 17개 시·도이다. 전국 ASOS 관측소 중 분석 기간 동안 연속 관측이 가능했던 97개 관측소를 우선 선정하였으며(Table 2), 각 시·도의 대표 관측소는 인구 규모(국가데이터처 자료 기준)가 가장 큰 지역의 ASOS 지점을 기준으로 설정하여 위치를 시각적으로 제시하였다(Fig. 1). 한편, 사망자료가 17개 시·도 행정구역 단위로 제공되는 점을 고려하여, 실제 분석에서는 97개 ASOS 관측소 자료를 시·도 단위로 구분한 후 권역별 평균값을 산출하여 기상 변수로 사용하였다. 다만 서울, 부산, 인천 등 일부 대도시는 단일 대표 관측소를 통해 사망자료와 기상자료 간의 1:1 대응이 가능하였으며(Table 2), 그 외 지역은 권역 평균값을 적용하였다. 한편, 세종 관측소는 2019년 6월부터 자료가 제공되므로 본 연구에서는 제한적으로 사용하였다.

Table 1.

Summary of meteorological and statistical data used in this study

Fig. 1.

Map of the study area. The blue circles and red squares show locations of automated synoptic observing system (ASOS) station and representative station, respectively

기상자료와 함께, 본 연구에서는 사망자 자료를 국가데이터처에서 집계·제공하는 마이크로데이터인 사망원인통계 A형 원시자료를 수집하여 활용하였다. 국가데이터처 사망자 자료는 실제 사망시각을 포함하고 있어 질병관리청 자료에 비해 신뢰도가 높으며, 장기간의 자료 확보가 가능하다. 사망원인통계 A형에는 사망일자, 시·도 구분 코드, 사망원인 코드 등이 일별로 제공된다(Table 1). 사망원인 코드 분류는 국제사인분류(International Classification of Disease, ICD) 기준을 적용하였다. 특히, 최근 선행 연구들(Achebak et al., 2024; Lv et al., 2025; Pehlivan et al., 2024)에서 기온과 사망률 간의 연관성이 보고된 심혈관계질환(Cardiovascular Disease, CVD)인 I00–I99 코드와 호흡기계질환(Respiratory Disease, Resp)인 J00–J99 코드에 해당하는 사망자를 포함하였다. 아울러 한파로 인한 건강 영향은 뇌졸중, 급성 심근경색, 천식 등 기저질환자의 간접적 사망뿐만 아니라 저체온증(T68), 동창(T69)과 같은 직접적인 영향으로도 나타날 수 있으므로, 본 연구에서는 이러한 원인도 분석에 추가적으로 반영하였다. 사망원인별 통계는 1997년부터 제공되지만, 사망시각 자료는 2000년 이후부터 구축되어 있다. 이에 따라 기상청 자료와의 연구 기간을 일치시키기 위하여 2000 ~ 2023년을 분석 기간으로 설정하였다.

2.2. 연구방법

국가데이터처 사망원인통계 자료는 ICD 코드를 포함하여 행정구역(전국 17개 시·도) 단위로 사망자 주소지가 입력되어 있어, 기상자료와 직접적인 연계에는 한계가 존재한다. Table 2와 같이 서울, 부산, 인천 등 대규모 도시는 사망자 자료와 기상자료 간의 1:1 연계가 가능하였으나, 경기도, 강원도 등 기타 지역의 경우 권역 내 ASOS 관측소들의 평균값을 대체하여 사용하였다. 한편 풍향은 원형(circular) 변수의 특성상 단순 일평균값으로 사용할 경우 왜곡이 발생할 수 있으며(Rutgersson et al., 2015; Satari et al., 2015), 특히 한파 연구에서는 특정 방향의 기류 유입이 중요한 의미를 갖는다. 이러한 점을 고려하여 본 연구에서는 각 시·도의 풍향에 대해서는 일최다풍향 자료를 활용하였다.

Table 2.

Classification of 17 administrative districts and 97 ASOS stations used

또한 지역 간 인구 규모 차이를 보정하기 위하여 사망자수를 정규화하였다. 구체적으로, 각 시·도의 일별 사망자수는 해당 지역 총인구로 나누어 인구 100만 명당 사망자수로 환산하였으며, 이는 Garner et al. (2020)이 기상재해 취약성 분석에서 제시한 방법론을 준용한 것이다. 정규화된 사망자수는 식 (1)과 같이 정의된다.

여기서, D_(j)는 시·도 j의 일별 사망자수, P_(j)는 해당 시·도의 총인구, M_norm(j)는 인구 100만 명당 사망자수를 의미한다. 이러한 정규화 과정을 통해 지역 간 인구 규모 불균형으로 인한 왜곡을 최소화하고, 한파에 따른 사망 영향의 상대적 크기를 비교할 수 있도록 하였다. 다음으로 한파의 시간지연 효과(time-lagged effect)를 탐색적으로 파악하기 위하여, 특정 시·도에서 한파일수로 정의된(아침 최저기온 ≤ –12℃) 날짜를 기준일(t_i)로 설정하였다. 이후 기준일로부터 최대 7일까지(t_i+1, ..., t_i+7) 정규화된 사망자수를 추적하였으며, 지연일별 평균 사망자수(M_lag=k,j)를 식 (2)와 같이 산출하였다.

여기서, D_j,ti+k는 시·도 j에서 한파 발생일 t_i로부터 k일 뒤의 정규화된 사망자수, N_j(k)는 해당 lag k일 이후 값이 존재하는 한파 발생일의 수를 의미한다. 한파일수가 여러 날 연속해서 발생하면 각 lag 간 정보가 겹치므로, 연속 한파일자 중 첫 날만 선택하였다. 이를 통해 지연일별 평균 사망자수의 변화를 비교하고, 어느 시차에서 사망 증가의 연관성이 가장 크게 나타나는지를 확인하였다. 여기서 사용한 절대 기준(아침 최저기온 ≤ -12℃)은 인체가 경험하는 절대적인 한랭 노출 수준을 반영하는 지표로, 급격한 기온 저하에 따른 단기적 건강 반응과 시간지연 효과를 평가하는 데 초점을 두었다.

또한 한파일의 강도와 지속일수(duration)에 따른 초과 사망 특성을 종합적으로 평가하기 위해, 최저기온 백분위를 기반으로 한 강도 분석과 연속 한파일 지속일수에 기반한 분석의 두 가지 접근법을 병행하였다(Kim et al., 2022; Meng et al., 2021). 한파일 강도는 각 관측소별 겨울철 최저기온 분포의 하위 백분위수(percentile)를 기준으로 정의하였으며, 하위 1% 이하 구간(P1)과 하위 1 ~ 5% 구간(P5)을 각각의 분석 구간으로 설정하였다. 이후 한파 발생일을 기준으로 0 ~ 7 일의 초과 사망자(ΔM_k,c)를 식 (3)과 같이 산출하여, 한파일 강도별 시간지연 효과를 분석하였다.

여기서, $$ {\stackrel{-}{M}}_{k,c}$$는 사망원인군 c의 lag k일 평균 정규화 사망자수, $$ {\stackrel{-}{M}}_{noncold,c}$$는 해당 원인군의 비한파일 평균 사망자수를 의미한다. 백분위수 기반 정의는 각 지역의 기후적 평년 상태를 고려한 상대적인 한파 강도 수준을 반영하는 지표로, 동일한 절대 기온 조건에서도 한파의 영향이 어떻게 달라지는지를 비교·평가하는 데 활용된다. 다음으로, 한파일의 지속성에 따른 누적 건강영향을 평가하기 위해 한파일이 연속적으로 지속된 일수를 기준으로 구간화하였다. 각 한파일 구간에 대해 사망자 평균을 산출하고, 비한파일 평균을 기준값(baseline)으로 한 초과 사망자를 계산하였다. Baseline 정의는 강도 분석과 동일하게 지점별(station) 특성을 반영한 독립적 평균값을 적용하였으며, 두 접근법 모두 겨울 간 교차 영향을 방지하기 위해 경계월/연도 처리를 수행하였다. 즉, 12월부터 다음 해 2월까지의 기간을 하나의 겨울 시즌으로 그룹화하여, 서로 다른 겨울 시즌 간 한파일이 연속되거나 중복되지 않도록 하였다.

마지막으로, 한파와 관련된 사망자 영향에 대한 다요인적 기상요소의 기여도를 평가하기 위하여 앙상블 학습 기법 중 하나인 Gradient Boosting Machine (GBM) 모형을 활용하였다. 앙상블 기법은 크게 Bagging과 Boosting으로 구분되는데, Bagging 기법의 대표적인 알고리즘인 Random Forest는 여러 개의 독립적인 결정트리를 병렬적으로 학습시켜 분산을 줄이고 과적합을 방지하는 데 효과적이다. 그러나 Random Forest는 단순 평균을 통해 결과를 도출하기 때문에 편향(bias)을 줄이는 데 한계가 있으며, 복잡한 비선형 관계를 충분히 반영하지 못할 수 있다. 반면 GBM은 Boosting 계열 알고리즘으로, 약한 학습기(weak learner)를 순차적으로 결합하면서 이전 단계에서 발생한 오차를 보정하는 방식으로 학습한다. 이러한 순차적 보정 과정은 모형의 분산뿐만 아니라 편향을 동시에 줄일 수 있어, 기상 변수와 같이 상호작용이 복잡하고 비선형성이 강한 자료에서 더욱 높은 예측 성능을 기대할 수 있다. 또한 GBM은 학습률(learning rate) 조정, 정규화(regularization) 기법, 표본 추출(sampling) 등을 통해 과적합 문제를 완화할 수 있어 안정성과 설명력을 동시에 확보할 수 있다. 본 연구에서는 R 프로그램을 이용하여 GBM 모형을 구현하였다.

3.1. 한파일 및 사망자 발생 추세 분석

한파일은 시기와 지역에 따라 발생 특성이 달라지며, 누적 규모는 사망자 분석에 의미 있는 참고자료가 될 수 있다. 따라서 장기적인 한파일 추세를 먼저 파악할 필요가 있다. Fig. 2(a)는 2000년부터 2023년까지 전국 97개 ASOS 관측소를 기준으로 연도별 겨울철 한파일수의 누적 합계를 보여준다. 2010년, 2012년, 2017년 겨울철은 분석 기간 중 한파일수가 가장 높은 수준을 기록하였다. 반면, 2019년 겨울철은 시베리아 고기압이 상대적으로 약화되었고(Zhang et al., 2021), 이에 따라 분석 기간 중 가장 적은 약 60일의 한파일수가 관측되었다. Locally Estimated Scatterplot Smoothing (LOESS) 추세선은 2011년을 전후로 한파일수가 점차 감소하는 경향을 보여주었다.

Fig. 2.

(a) Annual total number of cold-wave days aggregated across all ASOS stations in Korea from 2000 to 2023. Light blue bars indicate yearly totals, and the red line represents the smoothed trend estimated using LOESS regression. (b) Province-level distribution of cold-wave days across the 17 administrative regions (sido) during the same period. Bars show the mean number of cold-wave days, and error bars indicate standard errors (SE)

Fig. 2(b)는 동일한 기간 동안 시·도별 평균 한파일수와 표준오차를 시각화한 결과로, 한파일수의 전국적 분포 경향을 파악하기 위해 제시하였다. 전국적으로 뚜렷한 지역적 차이가 나타났으며, 특히 경기, 강원, 충북 지역은 각각 평균 250일 이상으로 가장 높은 한파일수를 기록하였다. 이들 지역은 중부 내륙에 위치하거나 지형고도가 높아 겨울철 기온 하강이 빈번한 기후적 특성을 지니며, 이에 따라 한파일 기준을 충족하는 날이 다른 지역에 비해 더 자주 나타난 것으로 해석된다. 반면, 부산, 광주, 울산, 전남 등 남부 해안권 지역은 평균 한파일수가 5일 미만으로 낮게 나타났으며, 특히 제주도에서는 분석 기간 동안 한파일이 발생하지 않은 것으로 확인되었다. 이러한 결과는 절대 임계값 기반 한파 정의의 영향이 반영된 것으로, 상대적으로 온화한 남부 및 해안 지역에서는 한파 발생이 과소평가되었을 가능성이 있다. 인천은 평균값에 비해 표준오차가 상대적으로 크게 나타났으며, 이는 해당 지역에서의 연도별 한파일 발생 빈도의 변동성이 컸음을 시사한다. 한편, 월별 전체 한파일수는 1월(5,961건)이 가장 많았고, 12월(2,977건), 2월(2,051건) 순으로 나타나 겨울철 중 1월에 한파일이 가장 빈번하게 발생하는 경향이 확인되었다.

전체 사망자는 ICD 코드를 기준으로 심혈관계질환(CVD) 및 호흡기계질환(Resp)에 의한 사망자와, 저체온증(T68) 및 동창(T69)과 같은 한파 관련 사망자로 분류한 후, 성별, 연령군, 하루 중 시간대별로 구분하여 분석하였으며, 그 결과는 Table 3에 제시하였다. CVD+Resp 사망자수는 남성과 여성 모두 약 27만명 수준으로 유사하게 나타났으며, 이 중 65세 이상 고령층이 남성은 전체의 약 77%, 여성은 약 91%를 차지하였다. 특히 여성의 경우, CVD+Resp 사망자의 대부분이 고령층에서 발생한 것으로 나타났다. T68+T69 사망자에서는 고령층 비중이 여성(약 73%)에서 높게 나타난 반면, 남성은 약 34%로 상대적으로 낮았다. 이는 2023년 기준 여성 고령자 인구 비중이 더 높음(성비 78.7%)에도 불구하고, 전체 T68+T69 사망자수는 남성이 여성보다 2배 이상 많았던 점에서, 남성은 비고령층에서도 외부 활동, 산업재해, 주거환경 등 한파로 인한 위험에 더 많이 노출되었을 가능성을 시사한다. 시간대별 분포에서는 CVD+Resp 및 T68+T69 사망자 모두 오전 시간대(06:00 ~ 11:59 KST)에 가장 많이 집중되었으며, 이는 이른 아침 급격한 기온 변화나 생리적 스트레스 반응 등이 건강에 영향을 미쳤을 가능성과 관련이 있는 것으로 판단된다.

Table 3.

Number of deaths related to cardiovascular and respiratory diseases (CVD+Resp) and cold-related illnesses (T68+T69), stratified by gender, age group, and time period

3.2. 최저기온과 사망자 간의 연관성 분석

최저기온은 한파 환경을 나타내는 대표적인 지표 중 하나로, 저온 조건에서의 기온 특성을 반영한다. 이에 따라 본 연구에서는 사망자 발생일의 일 최저기온 분포를 질병 유형별로 비교하여, 기온 조건에 따른 사망 발생 특성을 탐색적으로 살펴보고자 하였다.

Fig. 3은 CVD+Resp 사망자와 T68+T69 사망자에 대해 사망 발생일의 일 최저기온 분포를 Boxplot으로 비교한 결과이다. 두 질병군 모두 전반적으로 음의 기온 구간에서 사망자가 많이 분포하는 경향을 보였으나, T68+T69 사망자는 보다 낮은 기온 범위에 집중되어 있었고, 평균(녹색 삼각형) 및 중앙값 또한 –5℃ 이하에서 형성되었다. 반면, CVD+Resp 사망자는 상대적으로 더 높은 기온 범위에서도 분포가 나타났으며, 75백분위수(Q3) 기준으로는 0℃ 이상의 기온에서도 사망이 관측되었다. 다만 본 결과는 사망 발생 시점의 기온 분포를 단순 비교한 것으로, 실제 사망 위험과의 직접적인 인과관계를 의미하기보다는 기온 조건에 따른 분포 특성을 보여주는 기술적 결과로 해석될 필요가 있다.

Fig. 3.

Boxplots of daily minimum air temperatures (T_MIN) by cause of death. The green triangles denote the mean value

Fig. 4는 한반도 16개 시·도에서 겨울철 T68+T69 사망자와 CVD+Resp 사망자가 발생한 날짜의 기온 특성을 비교한 결과이다. 연구 기간 동안 세종(시·도코드: 29)에서는 T68+T69 사망자가 보고되지 않아 분석에서 제외하였다. 시·도별 T68+T69 사망일의 최저기온(T_MIN) 평균을 기준으로 오름차순 정렬한 결과, 경기·강원·충북 등 중부 내륙 지역에서 가장 낮은 기온 수준(–7.8 ~ –8.8℃)이 나타났고, 부산, 제주 등 남부 및 해안 지역에서는 상대적으로 높은 기온 조건에서 사망자가 발생하였다. 이는 Fig. 2(b)에서 확인된 바와 같이, 경기·강원·충북 등이 평균 한파일수가 가장 많은 지역이라는 점과 일치하며, 한파 발생 빈도가 높은 지역일수록 사망 발생일의 절대적 저온 수준 또한 더욱 낮아짐을 시사한다. CVD+Resp 사망자군의 경우에도 지역별 기온 패턴은 유사하였으나, 최고기온(T_MAX), 평균기온(T_AVE), 최저기온(T_MIN) 모두 T68+T69 사망자군보다 일관되게 높아 상대적으로 덜 낮은 기온 조건에서도 사망자가 발생함을 보였다. 한편, 인천은 사망 발생일의 T_MAX와 T_AVE가 낮고 T_MIN의 하강 폭이 작아 일교차가 작은 특징을 보였다. 이는 서해에 인접한 지리적 특성으로 인해, 낮 동안 기온 상승이 억제되고 밤 동안 급격한 냉각이 제한되는 기후 특성이 반영된 결과로, 인천에서의 한파일은 내륙 지역에서 흔히 나타나는 극저온 기반 한파와는 다르게 완만한 저온이 지속되는 형태의 한랭 환경으로 형성되었을 가능성을 시사한다.

Fig. 4.

Mean daily maximum (T_MAX), average (T_AVE), and minimum (T_MIN) air temperatures on days with cold-related deaths (T68+T69; solid lines) and cardiovascular-respiratory deaths (CVD+Resp; dotted lines) across the 16 provinces of Korea

Fig. 5는 겨울철 기간을 세부 구간으로 구분하여, 두 질병군(T68+T69, CVD+Resp)의 사망자수와 사망 발생일의 평균 T_MIN 변화를 비교한 결과이다. 먼저 CVD+Resp 사망자수는 12월 초순(Early Dec)부터 2월 중순(Mid Feb)까지 약 4만2천 ~ 4만7천 명 수준으로 유지하여 겨울철 전반에 걸쳐 큰 변동 없이 높은 수준을 보였다. 다만 2월 하순(Late Feb)에는 사망자수가 약 3만 중반 수준까지 감소하여 계절 후반부에 완만한 하락 패턴을 나타냈다. 반면 T68+T69 사망자수는 전체 규모는 작지만 뚜렷한 계절적 변동성을 보였으며, 12월 초순 이후 증가하여 12월 하순에 최대 사망자수(477명)에 도달한 뒤, 1월 하순까지 비교적 높은 수준을 유지하였다가 2월 초순 이후 감소하는 흐름을 나타냈다. 두 사망군의 평균 T_MIN은 겨울철 경과에 따라 유사한 U자형 계절 패턴을 보였으나, 절대적 기온 수준은 모든 시점에서 T68+T69가 CVD+Resp보다 일관되게 낮았다. 특히 기온이 최저에 도달하는 1월 초·중순에는 두 사망군 간 T_MIN의 차이가 약 2.5℃까지 확대되어, 직접 한랭 사망(T68+T69)이 상대적으로 더 극저온 조건에서 발생함을 보여주었다. 또한 CVD+Resp 사망자수는 T_MIN 변화와 뚜렷한 연관성이 거의 없었던 반면, T68+T69 사망자수는 기온 하강 구간에서 증가하고 기온 상승 구간에서 감소하는 반대 양상(out of phase)을 나타냈다. 다만, 최저기온이 겨울철 가장 낮아지는 1월 초·중순에 T68+T69 사망자수가 추가적으로 증가하지 않고 일정 수준에서 유지된 점은 기온 요인만으로 설명하기 어려우며, 초기 한랭 노출에 대한 취약성 발현, 개인 및 사회적 난방 대비 강화, 한랭 스트레스에 대한 단기 적응 등 여러 비기온적 요인이 복합적으로 작용했을 가능성이 있다.

Fig. 5.

Winter sub-period variations in death counts (bars; left y-axis) and mean daily minimum air temperatures (dashed lines; right y-axis) for T68+T69 and CVD+Resp groups. The number of T68+T69 deaths was multiplied by 100

3.3. 한파에 의한 초과 사망자 분석

앞선 3.2절에서는 겨울철 기온 특성과 사망 발생 간의 연관성을 지역별·시기별로 비교함으로써, 한파가 사망 위험을 증가시키는 주요 기후적 조건임을 확인하였다. 그러나 한파의 영향은 반드시 즉시 나타나는 것이 아니라, 노출 이후 수일에 걸쳐 지연적으로 발생하거나 한파일이 일정 기간 지속될 때 누적적으로 확대될 수 있다. 이에 3.3절에서는 한파가 사망에 미치는 영향을 정량적으로 평가하기 위해, 한파 발생 이후 나타나는 시간지연 효과(time-lagged effect)와 한파일이 연속적으로 지속될　때 나타나는 누적 영향(duration effect)에 초점을 두어 초과 사망자수를 산출하였다. 초과 사망자수는 인구 규모 차이를 보정하기 위해 인구 100만 명당 정규화 사망자수를 기반으로 하였으며, 연속된 시계열 정보가 필요한 점을 고려하여 분석 대상은 전체(All-cause) 사망자와 심혈관계·호흡기계(CVD+Resp) 사망자로 한정하였다. 이를 통해 한파의 단기 영향뿐만 아니라 한랭 노출의 지속성과 누적성이 사망 위험 증가에 어떠한 기여를 하는지 규명하고자 하였다.

Fig. 6은 한파 발생 이후 0 ~ 7일 동안의 초과 사망자수를 한파 강도(P1, P5)에 따라 구분하여 제시한 결과이다. 전체(All-cause) 사망자의 경우, 강한 한파일(P1)에서는 초과 사망자수가 모든 지연일에서 양(+)의 값을 보였으며, 특히 lag 1일에서 인구 100만 명당 약 2.5명으로 가장 크게 나타나 강한 한파 직후 초과 사망이 뚜렷하게 증가함을 보여주었다. 이후 lag 2 ~ 5일까지도 약 1.7 ~ 2.0명 수준의 증가가 지속되어, 강한 한파의 영향이 노출 당일 직후뿐 아니라 최대 5일까지 이어지는 것으로 나타났다. 반면 P5 구간에서는 초과 사망자수가 P1보다 훨씬 작은 규모로 나타났다. lag 1일과 lag 6일에서 약 0.6명 이상의 증가가 나타나는 bi-modal 형태가 관찰되었지만 대체로 0.1 ~ 0.3명 수준에서 유지되었다. 이는 중간 강도의 한파에서 초과 사망 변화가 뚜렷하게 지속되지 않고 일부 지점에서만 일시적 증가가 나타난 결과로 해석된다. CVD+Resp 사망자군의 경우에도 전반적인 시간지연 구조는 All-cause와 유사하였으나, 초과 사망자수의 절대 규모는 더 작고 변동폭도 완만하게 나타났다. P1 구간에서는 lag 1 ~ 4일 사이에 약 0.3 ~ 0.5명 수준의 초과 사망이 나타났으며, P5 구간에서는 대부분의 지연일에서 ±0.2명 내외의 미미한 변화만 나타났다. 즉, 한파 강도가 클수록(P1 > P5) 시간지연 초과사망 영향이 커지는 경향은 동일하게 나타났으나, 그 영향의 크기와 지속성은 All-cause에서 더욱 뚜렷하게 나타났다. 이는 All-cause가 더 큰 모집단을 포함하여 한파 노출 시 영향을 받을 수 있는 잠재적 위험군이 더 넓기 때문에 CVD-Resp에 비해 초과 사망자가 크게 나타난 것으로 보인다.

Fig. 6.

Lagged excess mortality (ΔMk) following cold-wave days for All-cause (solid lines) and CVD+Resp (dashed lines) deaths, stratified by cold intensity―P1 (lowest 1%; orange lines) and P5 (1-5% range; green lines)

Fig. 7은 한파 지속시간(duration)에 따른 초과 사망자수를 지속기간별 구간(1일, 2 ~ 3일, 4 ~ 5일)으로 구분하여 제시한 결과이다. 전체(All-cause) 사망자의 경우, 1일 지속 시 초과 사망자수는 인구 100만 명당 약 0.6명, 2 ~ 3일 지속 시에는 약 0.55명으로 나타나 두 구간 사이의 차이는 크지 않았다. 반면 4 ~ 5일 지속 구간에서는 초과 사망자가 약 2.1명으로 급격히 증가하여, 한파가 장기간(4일 이상) 이어질 때 초과사망이 비선형적으로 확대되는 특징이 확인되었다. CVD+Resp 사망자군에서는 1일 지속 구간에서 약 -0.1명으로 초과 사망이 나타나지 않았으나, 2 ~ 3일 지속 구간에서는 약 0.2명으로 음(–)에서 양(+)으로 전환되는 증가가 관찰되었다. 이어 4 ~ 5일 지속 구간에서는 약 0.4명으로 추가 상승하여, 전체 사망자군에 비해 절대 규모는 작지만 지속기간 증가에 따른 초과 사망의 점진적·선형적 증가 경향이 뚜렷하게 나타났다. 즉, 심혈관계·호흡기계 사망자군에서는 한파 발생 직후의 단기 영향보다 한파가 며칠간 지속되는지가 사망 위험 증가를 설명하는 데 보다 중요한 요인으로 작용한 것으로 판단된다.

Fig. 7.

Excess mortality (ΔM) by cold-wave duration category, showing mean ± standard error (SE) relative to station- and cause-specific baselines, with blue and orange bars indicating All-cause and CVD+Resp deaths, respectively

3.4. GBM 모형 기반 기상요소 영향 평가

한파와 관련한 기상요소 간의 복합적 상호작용을 고려하기 위해, 본 절에서는 GBM 모형을 활용하여 다요인적 기상요소가 겨울철 사망자 변동에 기여하는 상대적 중요도를 평가하였다. 주요 기상요소는 모두 당일 관측값으로 구성하였으며, 사망자수는 당일(lag 0)과 1일 지연(lag 1)을 대상으로 분석하였다. 앞서 3.3절의 시간지연 효과 분석에서 한파 발생 이후 lag 1 시점에서 초과 사망자수가 가장 뚜렷하게 나타났고, lag 0과 lag 1 분석 결과 간 비선형 반응과 상호작용 패턴 또한 전반적으로 유사하게 확인되었다. 이에 따라 본 논문에서는 한파의 단기 지연 효과를 보다 명확히 반영하는 lag 1 결과를 중심으로 제시하였다. 한편, 종속변수로는 한파에 취약한 주요 사망원인군으로 알려진 CVD+Resp 사망자를 설정하여 모형을 구축하였다. 직접 한랭 관련 사망원인(T68+T69)은 사례수가 극히 제한적이므로 본 분석에서는 제외하였다. 사망자 예측모형의 신뢰성을 확인하기 위해 전체 자료를 80%의 학습자료와 20%의 검증자료로 분리하여 교차검증을 수행하였으며, GBM 모형의 성능 평가는 Table 4에 제시하였다. 모형 학습 과정에서는 반복 횟수(n.trees)를 200으로 설정하고, 100회의 Bootstrap 반복을 적용하여 변수 기여도의 안정성을 확보하고자 하였다. lag 0 및 lag 1에 대한 성능 평가 결과, 전체 자료의 약 69 ~ 70%가 factor-of-two (FAC2) 범위 내에 분포하여 관측값의 분포 특성을 비교적 잘 재현하는 것으로 나타났다. RMSE와 상관계수(r)는 학습자료와 검증자료에서 거의 동일한 수준을 보였으며, 이는 모형의 과적합(overfitting)이 발생하지 않았음을 의미한다. 편향 지표(NMB = 0.00 ~ 0.01)는 매우 낮은 수준이었고, 기상-사망 변동 구조에 대한 학습 정도를 나타내는 IOA 역시 0.57 ~ 0.58 범위로 나타나 모형이 비선형적 관계를 안정적으로 학습하고 있음을 확인하였다(Table 4). 이러한 결과를 종합할 때, 본 연구에서 구축한 GBM 모형은 기상요소와 겨울철 사망자 간의 복합적인 관계를 분석하기에 합리적인 수준의 신뢰성을 확보한 것으로 판단하였다.

Table 4.

Performance statistics of the GBM model

Fig. 8은 GBM 모형을 통해 산출된 주요 기상 변수들의 부분의존도(Partial Dependence Plots, PDP)를 제시한 것으로, 1일 지연(lag 1) 사망자수에 대한 각 변수의 상대적 기여도와 비선형 반응 특성을 보여준다. 각 PDP 하단의 러그(rug)는 관측치 분포를 나타내며, 관측 빈도가 낮은 구간에서는 불확실성이 상대적으로 크게 나타난다. 전체적으로 일사량(26.1%)과 최다풍향(25.9%)이 겨울철 사망자 변동에 가장 높은 기여도를 보였으며, 강수량은 기여도(5.2%)가 가장 낮고 불확실성이 크게 나타났다. 일사량은 값이 증가할수록 사망자수가 감소하는 경향을 보였으며, 풍향은 북서풍 및 동풍 계열에서 사망자수가 상대적으로 높게 나타나는 비선형적 반응이 확인되었다. 또한 일최저기온(T_MIN)은 기온이 낮아질수록 사망자수가 증가하는 경향을 보였으며, 극저온 구간에서는 자료 수의 제한으로 불확실성이 확대되는 특징을 보였다. 이슬점온도(T_DEW) 역시 낮은 값에서 사망자수가 증가하는 경향을 보였고, 일교차(Diurnal Temperature Range, DTR)는 증가할수록 사망자수가 증가하는 양상을 나타냈다. 이러한 결과는 겨울철 사망자 변동이 단일 기온 변수뿐 아니라 다양한 기상요소들의 복합적 작용에 의해 영향을 받음을 보여준다.

Fig. 8.

Partial dependence plots of key meteorological variables for lag 1 winter death. Shaded areas indicate uncertainty, and relative variable influence (%) is shown in each panel

실제로 극한 저온 조건에서는 기온과 함께 일사량, 풍향 등 다른 기상 요인들이 동시에 작용하며, 이들 간의 상호작용에 따라 사망 위험이 증폭되거나 완화될 가능성이 존재한다. 이러한 점을 보다 구체적으로 확인하기 위해 Fig. 9에서는 이변량 상호작용 분석(plot2Way)을 통해 일최저기온과 주요 기상 변수 간의 결합 효과를 살펴보았다. 먼저 T_MIN과 일사량의 상호작용에서는 T_MIN이 0℃ 미만, 특히 -10℃ 이하이면서 일사량이 2 MJ m^-2 이하인 조건에서 lag 1 사망자수가 약 25 ~ 30명 수준으로 가장 높게 나타났으며, 동일한 저온 조건에서도 일사량이 10 MJ m^-2 이상으로 증가할 경우 사망자수는 약 20명 내외로 감소하였다(Fig. 9(a)). 풍향과의 상호작용에서는 저온 조건에서 북서풍 계열(약 280 ~ 330°)과 동풍 계열에서 사망자수가 상대적으로 높게 나타나는 경향이 확인되었다(Fig. 9(b)). 더욱이, 이러한 풍향별 차이는 Fig. 10에 제시된 최다풍향·풍속에 따른 사망자수 분포에서도 보다 직관적으로 관찰되었다. 4 m s^-1 이하의 북서풍 조건에서 사망자수(20명 이상)가 집중적으로 분포하는 패턴이 확인되었으며, 2 m s^-1 이하의 동풍 계열에서도 일정 수준 이상의 사망자 증가가 나타났다. 즉, 겨울철 시베리아 고기압 확장과 연계된 찬 공기 유입이 한반도에 직접적으로 영향을 미치는 상황에서 한파의 건강 부담이 가장 크게 증가하였고, 영동지역을 포함한 동해안 일부 지역에서는 북동 기류에 의한 차가운 공기가 유입되는 상황이 한파 사망 반응에 일부 기여했던 것으로 추정된다. T_MIN과 T_DEW의 상호작용에서는 T_MIN이 –15℃ 이하이면서 T_DEW가 –25℃ 이하인 조건에서 사망자수는 약 25 ~ 30명 수준으로 가장 높게 분포하였다(Fig. 9(c)). 또한 T_MIN과 DTR의 상호작용에서는 T_MIN이 –15℃ 이하이면서 DTR이 15℃ 이상인 경우 사망자수가 약 25명 내외로 높게 나타났으며, T_MIN이 0℃ 이상인 조건에서는 DTR 변화에 따른 사망자수의 민감도가 크게 감소하였다(Fig. 9(d)). 이는 겨울철 급격한 기온 변동성이 한파 상황에서 사망 위험을 추가적으로 증폭시키는 요인으로 작용했던 것으로 사료된다.

Fig. 9.

Two-way interaction plots of lag 1 winter death between T_MIN　(units: ℃) and (a) insolation (units:　MJ·m-2), (b) mode wind direction (units: °), (c) T_DEW (units: ℃), and (d) DTR (units: ℃) from the GBM model

Fig. 10.

Polar plot of lag 1 winter death according to mode wind direction (units: °) and wind speed (m·s-1) during the study period