패널모형을 활용한 논벼 수량성과 온도조건 간의 상관성 분석: 3개 농업기후지대의 지역성 비교를 중심으로

1. 서론

2019년 한반도의 평균기온은 역대 두 번째로 높았으며, 최근 10년간(2010년 ~ 2019년) 한반도의 평균기온은 평년(1981년 ~ 2010년)에 비하여 0.5℃ 높았다. 이와 같은 기후변화는 폭염일수의 150% 상승과 기후변동성의 증가를 유발하고 있다(KMA et al., 2020).

한국의 주식작물인 논벼의 수량성에 있어서, 국립식량과학원 온도구배온실에서의 실험은 이산화탄소 농도 상승으로 인한 기후변화의 긍정적 효과가 기온의 상승이 유발하는 등숙율 저하, 낟알 무게 감소 등의 부정적인 영향에 의해 상쇄됨을 확인하였다(Lee et al., 2013). 다만, DSSAT 4.0-CERES-Rice 작물모형에 기반을 둔 분석은 이와 같은 상쇄의 정도가 지역에 따라 서로 상이할 것으로 예상하였다(Kim et al., 2013). 또한, APEX-Paddy 모델을 활용한 이천과 김제 지역에서의 분석 역시 지역과 RCP 시나리오에 따라 변화의 방향과 크기가 상이할 것으로 예상하였다(Choi et al., 2017).

작물생육모형에 기반을 둔 분석은 그 특성상 높은 농학적 타당성과 정밀성을 갖지만, 계산과 검증에 풍속, 일사량 등 다양한 데이터의 활용이 요구된다는 점에서 광범위한 지역에 걸친 통계실증분석에는 적용되기 어렵다는 한계를 갖는다. 또한, 작물생육모형은 작물의 생리학적 특성을 정확히 예측한다는 장점을 갖지만, 실제 재배에서는 작물이 자연적인 조건에만 노출되는 것이 아니라, 농가 차원의 기상조건 적응시도(예컨대, 논벼 농가에서 냉해피해 우려 시 물 깊이대기 실시)도 발생한다는 점에서 실제 수량성 통계와는 차이를 갖는다는 단점을 갖는다. 다시 말해, 단순한 작물생육모형에 기반을 둔 기후변화 예측은 ‘어리석은 농가 시나리오(dumb-farmer scenario)’를 가정한다는 점에서, 수량성 통계 또는 농가단위 생산액에 기반을 둔 계량경제학적 분석은 농가 수준의 적응을 반영할 수 있다는 장점을 갖는다(Mendelsohn et al., 1994)

한편, 생육온도일수(growing degree days, GDD)는 작물생육기간의 예측이나 품종의 조만성, 특정지역에서의 신품종 재배가능 여부 등을 판단하는데 주로 사용되었다(Jong et al., 1986; MAF 2006; Han et al., 2011). GDD를 계산하는 방법은 매우 다양한데, 이 중 일 최고온도와 최저온도의 평균에서 생육 최저한계온도인 기준온도를 뺀 나머지를 누적하여 계산하는 잉여온도형이 많이 이용되고 있다(Hyun et al., 1994; Lee 1983; MAF 2006). 하지만 이는 고정된 기준온도(논벼의 경우 통상적으로 10℃)를 활용하기 때문에, 정밀한 온도구간별 소요일수의 계산에는 적합하지 않다.

최근의 계량경제학적 연구는 고해상도 기온 및 강수데이터와 일중 온도분포에 대한 정밀한 분석의 결합이 기후인자와 식량작물 수량성 간의 연관성을 높은 설명력으로 유추할 수 있음을 제시한다(Schlenker and Roberts 2009; Blanc and Schlenker 2017). 구체적으로, 해당 접근은 농림기상학적 이론(Snyder 1985)에 기반을 두고 일중 최저ㆍ최고온도를 3℃ (또는 5℃) 단위로 분절한다. 예컨대 일중 최저온도가 4℃이고, 최고온도가 12℃인 날의 경우, 일중 온도의 종형(sinusoidal) 분포를 상정하여 4℃ ~ 7℃의 온도 구간에 0.5467일이 포함되었음을 추정할 수 있다.¹⁾ 이러한 접근은 일중 최저·최고온도를 활용한다는 점에서 GDD와 유사하지만, 온도구간의 분절을 통해 보다 유연한 분석이 가능하다는 차별성을 갖는다. 이와 같은 계산은 모든 생육주기에 걸쳐서 반복되며, 이를 통하여 특정 생육주기에서 해당 작물이 주어진 온도구간을 경험한 소요일수를 합산할 수 있다. 해당 계산을 활용한 계량경제학적 분석은 옥수수, 대두 등의 식량작물과 기후조건 간의 상관성을 분석하는 데에 상당한 타당성을 가짐이 입증되었다(Liu et al., 2016).

본 연구의 목표는 이와 같은 접근방법을 국내 논벼 데이터에도 적용, 기초지자체 단위에서도 기후변화 취약성 분석을 위해 활용 가능한 계량경제학적 모형을 구현하는 데에 있다. 이 과정에서 선행연구의 결론을 참조, 상이한 지리적 특성의 시군구를 선정해 지역 간 차이를 확인하였으며, 논벼의 생육주기별로 온도와 수량성 간의 상관성이 갖는 차이에 대해서도 분석을 진행하였다.

또한, 기상변수를 설정하는 데에 있어서, 본 연구는 개별 관측소의 데이터가 아닌 1 km 해상도의 정밀 자료를 통해 유효기상데이터를 추산하였다. 이와 같이 정밀한 기상데이터는 한국에서의 식량작물에 대한 계량분석에 적용된 바 없으며, 따라서 향후 기후변화 조건 하 작물의 취약성을 분석하는 데에 유의미한 시사점을 제공할 것으로 예상한다.

2. 연구방법

2.1 연구대상 지역 및 기간

본 연구는 전국 19개 농업기후지대 중 7개 이상의 시군구를 갖고, 지리적으로 이격되어있는 지대를 선정하였다. 이를 통하여 총 세 개의 상이한 농업기후지대(중북부 내륙, 소백서부 내륙, 영남 내륙, Choi and Yun, 1989)에 걸친 27개 시군구가 연구대상 지역으로 설정되었으며, 18년(2000년 ~ 2017년)에 걸친 시군구별 논벼 수량성 및 기상자료를 수집하였다(Fig. 1).

Fig. 1.

A map of Korean mainland sigungus, with the three surveyed agroclimactic zones (as defined by Choi and Yun. 1985) colored in green, brown and orange

이앙이후 논벼 전체 생육주기(5월 말 ~ 9월 말, Table 2 참조)에 걸쳐서 중북부내륙지대(이하 ‘북부지역’)는 평균 27.68℃와 17.35℃의 일중 최고･최저온도를 보인 반면, 소백서부내륙지대(이하 ‘중부지역’)는 28.43℃와 18.35℃의 평균 최고·최저 온도를, 영남 내륙지대(이하 ‘남부지역’)는 평균 28.95℃와 19.34℃의 최고·최저 온도를 보였다. 생육주기 중 연평균 누적강수량은 각각 1061 mm, 901 mm, 893 mm로 상이하였다. 본 연구에서 종속변수로 활용한 수량성의 18년 평균값은 북부지역 627 kg/10 a(최대 696 kg/10 a, 최소 594 kg/10 a), 중부지역 평균 688 kg/10 a(최대 710 kg/10 a, 최소 662 kg/10 a), 남부지역 평균 653 kg/10 a(최대 666 kg/10 a, 최소 636 kg/10 a)로 각 지역 내에서 비교적 균일하였다.

2.2 분석자료

본 연구는 통계청의 『농작물생산조사』를 참조하여 18년(2000년 ~ 2017년)에 걸친 분석대상 27개 시군구의 논벼 수량성 데이터를 수집하였다(단위: kg/10 a). 해당 데이터는 시군별 표본농가에 대한 현장실측을 통하여 기록되며, 표본추출틀은 매년 갱신된다(Statistics Korea 2019).

본 연구는 기상청(KMA 2020)의 1 km 단위 격자 기상데이터(gridded weather data)를 활용하였다. 해당 데이터는 실측자료에 해양도, 고도, 경사방향 등의 미기상적인(micrometeorological) 요소를 반영 후 재생산된 정밀기상데이터로, 이는 농림축산식품부 제공 팜맵(AGIS 2020)과 결합되어 시군구별 논벼 생산지 기상자료를 형성하는 데에 활용되었다. 즉, GIS를 활용하여 각 격자에 포함된 논벼 재배면적을 전체 시군구 논벼 재배면적 대비 비율로 분할한 뒤, 이를 가중치로 활용하여 최저온도, 최고온도 및 강수량의 가중평균을 도출하였다. 이와 같이 시군구 논벼 재배지의 지리적 분포를 감안한 유효기온 및 강수량의 도출은 선행연구(Schlenker and Roberts 2009)의 핵심적인 특성 중 하나로, 단일 관측소에서의 관측치를 활용하는 경우와 달리 기상변수의 공간적 분포를 반영한다는 장점을 갖는다.

예컨대 어떤 시군구의 논벼 재배지역이 5개 격자에 걸쳐 분포하고, 팜맵 조회 결과 각 격자 내 논벼 면적이 0.2 km², 0.3 km², 0.5 km², 0.4 km². 그리고 0.6 km²였다면 해당 시군구의 총 논벼 재배면적은 2 km²였으므로, 각각의 격자 기상데이터에 0.1, 0.15, 0.25, 0.2, 그리고 0.3의 가중치를 곱하여 시군구 내 논벼재배지의 평균 기상 데이터 값을 구하는 것이다. 이와 같이 추출된 기온변수는 기상청의 실측자료와 비교하여 평균 0.996(최저 0.985, 최고 0.998)의 Pearson 상관계수를 보였다.

본 연구에서 활용한 기상변수는 작물생육기간 중 기상변수와 극한기상 변수로 분류된다. 전자에는 선행연구(Schlenker and Roberts 2009)를 참조하여 온도구간별 누적소요시간과 강수량이 포함되었으며, 후자에는 기상청이 기후변화 시나리오 하에서 제공 중인 극한기후지수 중 강수량과 관련된 3종(5일 최다강수량, 최대 무강수지속기간, 호우일수)이 포함되었다.

극한기상의 수량성에 대한 영향력을 평가하기 위하여 각 격자별로 집중호우(일 강수량 80 mm 이상) 발생 여부를 확인한 뒤, 이 역시 앞선 방식과 동일하게 가중평균하여 해당 시군구 논 재배면적 중 집중호우 피해가 발생한 면적비중을 계산하였다.²⁾ 이외의 5일간 최대 누적강우량의 최대치와 연속 무강수일 계산은 시군구 단위로 종합된 강수량 데이터를 활용하였다(Table 1).

Table 1.

List of data surveyed

2.3 분석모형

본 연구는 논벼의 생육주기별로 변화하는 기온 민감성을 반응하기 위하여 논벼가 노지 조건에 노출되는 시기³⁾를 분얼활착기(이하 ‘1기’), 유수형성 이후 출수개화 전(이하 ‘2기’), 그리고 등숙기(이하 ‘3기’)로 구분하였다. 또한, 지역별(북부, 중부, 남부)로 별도의 분석을 시도하여 기온조건에 대한 지역별 차이 역시 확인할 수 있도록 하였다.

논벼의 생육주기는 품종의 조만성(早晩性), 이앙시기 및 작기의 기후특성에 따라 달라진다. 본 연구는 분석의 일관성을 위하여 국립종자원의 정부 보급종 중 각 지역에 적합한 것으로 판단된 품종(북부지역 ‘오대’, 중부지역 ‘삼광’, 남부지역 ‘새일미’, KSVS, 2020)을 선택, 각 품종의 품종보고서(NICS, 2004; 2006; 2012)를 참고하여 생육주기를 선정하였다(Table 2).

Table 2.

Dates designated to three growth stages for each region

전술된 바와 같이, 본 연구는 생육주기 내 온도구간별 소요일수를 독립변수로 활용하여 논벼 수량성을 추정한다. 온도구간을 설정하기 위해서는 구간설정의 시작온도와 종결온도를 정하는 것이 필요한데, 선행연구에서는 결빙점 이하 온도를 배제하는 방식을 선택하였다(Schlenker and Roberts, 2009)⁴⁾. 본 연구에서는 각 지역의 온도분포에 맞춰진 온도구간을 선정하기 위해 세 지역별로 각 생육주기의 일 최저온도와 일 최고온도를 집계하였다. 이어서 시작온도와 종결온도 간의 구간이 전체 집계일의 99.9%를 포함할 수 있도록 구간을 선정하였다(Table 3).

Table 3.

Range of temperature bins designated to each region (based on quantiles)

즉, 북부지역에서는 생육주기별로 각각 10개, 8개, 10개의 3도 너비(width) 온도구간을 설정하고, 중부지역에서는 9개, 8개, 11개의 온도구간을, 남부지역에서는 10개, 8개, 11개의 온도구간을 선정하되, 각 구간은 지역적 온도 분포를 반영하도록 했다.

이외의 기상변수(강수량 및 극한기상 3종)는 생육주기별로 서로 다른 계수를 가질 수 있도록 설정하였다. 끝으로, 시간에 따른 논벼 재배기술의 변화를 포함하기 위하여 추세항과 추세제곱항을 추가하였다.

계량경제학적 모형으로는 선행연구를 참조하여 미관측 효과 패널모형(unobserved effects panel model)을 활용하였다. 해당 모형의 추정은 대상 지역의 시군구별로 기상데이터와 추세항으로는 확인되지 않는 시간불변 특성(예컨대 토성(土城))이 존재하여 논벼 수량성에 영향을 미친다는 가정 하에 실시된다. 이 미관측 특성이 독립변수(기상, 극한기상 및 추세항)와 상관성을 갖는지 여부에 따라 해당 모형은 고정효과 모형 또는 임의효과 모형의 형태로 구체화된다(Wooldridge, 2010). 시군구별 미관측 특성의 존부를 확인하기 위해 F-test를 우선 실시하였으며, 해당 특성과 독립변수간의 상관성을 검정하기 위해 Hausman 검정을 실시하였다. 이를 종합한 모형의 형태는 아래와 같다.

$\[ \mathrm{d} ln\left(y_{i,t}\right) =\begin{array}{l}{c}_i+\sum _{j=1}^{10}{\beta }_j{h}_{s=1,j,i,t}+\sum _{ j=11}^{18}{\beta }_j{h}_{s=2,j,i,t}\\ +\sum _{j=19}^{28}{\beta }_j{h}_{s=3,j, i,t}+\sum _{s=1}^3{\gamma }_s{p}_{s,i,t}\\ +\sum _{s=1}^3{\delta }_s^1{HR}_{s,i,t} +\sum _{s=1}^3{\delta }_s^2{cumprcp}_{s,i,t}\\ +\sum _{s=1}^3{\delta }_s^3{dconsec}_{s,i,t}+{\zeta }_1 T_t+{\zeta }_2{T}_t^2+{\epsilon }_{i,t}\end{array} \]$

y_it는 i 시군구 t년도에서의 수량성을, c_i는 i 시군구의 불변특성으로 인한 효과⁵⁾를, h_s,j,i,t는 i 시군구 t년도 s 생육주기 j번째 온도구간에서의 누적소요일수를 의미한다. p_s,i,t, HR_s,i,t, cumprcp_s,i,t, dconsec_s,i,t는 각각 s기 i 시군구 t년도에서의 강수량, 집중호우 피해면적 비율, 최대 5일 누적강수량, 최대 무강수일을 의미한다. T_t과 T_t²는 각각 추세항과 추세제곱항을 뜻하고, ε_i,t는 오차항을 지칭하며, E(ϵ_i,t|x_i, c_i)가정한다. 해당 수식에서 각 생육주기별 온도구간의 수(10개, 8개, 10개)는 북부지역의 사례이며, 중부지역과 남부지역의 생육주기별 온도구간 수는 상이하다.

이와 같이 계산된 모형의 오차항은 시간적(serial) 및 횡단면적(cross-sectional) 상관성을 가질 수 있다. 다시 말해, 특정 시군구에서 t년도의 수량성 충격은 t′년도(t≠t′)의 수량성 충격과 영향을 가질 수 있으며(시간적 상관성), 특정 연도에서 i 시군구의 수량성 충격은 i′ 시군구(i≠i′)의 수량성 충격과 상관성을 가질 수 있다(횡단면적 상관성). 각 상관성의 검정에는 Breusch-Godfrey 검정(Breusch 1978; Godfrey, 1978)과 Pesaran의 CD검정(Pesaran, 2004)을 활용하였으며, 귀무가설 기각 기준은 p-value 0.05였다. 오차항의 상관성이 존재할 경우 엄밀 표준편차(robust standard deviation)의 추정이 필요한데, 이를 위해서는 White와 Driscoll and Kraay 기법을 활용하였다(White, 1984; Driscoll and Kraay, 1998).

3. 연구 결과

3.2 지역별 기온-수량성 관계

추정 결과는 Table 4와 같다. 각 온도구간별 소요시간의 형태로 표현된 기온변수는 세 지역 모두에서 강수와 극한기상에 비하여 큰 중요성을 갖는 변수로 선정되었다.

Table 4.

Bin-based linear regression results

세 지역은 각 생육주기 기온조건의 역할에 대해 판이한 해석을 제시하였다. 북부지역에서는 기온변수가 뚜렷한 영향력을 갖는 시기가 부재하였던 반면, 중부지역에서는 3기가, 남부지역에서는 2기가 수량성에 긍정적으로 기여하는 시기인 것으로 드러났다.

이와 같은 수량성 반응에서의 차이를 감안하면, 저온피해･고온피해 취약성 여부는 각 지역 수준에서 판단하는 것이 적합하다. 다시 말해, 세 지역이 각 온도구간에 대해 갖는 계수의 범위가 상이하므로, 온도장해 발생은 장해온도 구간에서의 계수 추정치가 동일 생육주기 내 다른 온도구간의 그것과 비교하였을 때 차이를 보이는지의 여부를 통하여 확인되어야 한다.

예컨대, 남부지역의 경우(Fig. 2)에서는 각 생육주기의 온도구간이 다른 온도구간에 비하여 어떠한 계수를 갖는가를 통하여 저온피해와 고온피해의 존부가 판단되어야 하는데, Fig. 2를 통하여 남부지역에서는 1기, 2기, 3기 모두 최말온도구간에서 고온피해가 관측되고, 저온피해는 3기에서만 관측되었음을 확인할 수 있다. 같은 맥락에서, 북부지역에서는 1기 고온피해, 3기 저온피해가 뚜렷하며, 중부지역에서는 1기 저온피해, 2기 저온피해가 드러난다.

Fig. 2.

An illustration of the yield-temperature relationship, in the case of Southern region

각각 생육주기를 살펴보면, 북부지역과 남부지역은 1기 최말온도구간(bin 10, 각각 32℃ ~ 35℃, 35℃ ~ 38℃)에서 고온피해의 규모가 증가하는 것으로 나타났는데, 이는 33℃ 이상의 온도에서 논벼의 분얼발생이 억제된다(Yoshida 1973)는 선행연구와 일치한다. 또한, 냉해피해가 중부지역에서 발견된 구간(8℃ ~ 11℃) 역시 분얼활착기 논벼의 냉해피해 임계온도로 알려진 13℃(RDA 2014) 이하로, 농학분야 연구와 상당부분 일치한다.

2기에서는 냉해피해온도(15℃ ~ 17℃, RDA 2014) 이하 온도구간에서 중부지역에서 냉해피해가 드러났다. 이는 중부지역 이는 삼광벼 재배지역에서의 유수형성기 ~ 출수기 간 냉해피해가 임실률에 미치는 영향에 대한 농학분야 선행연구(Kim et al., 2009)와 일치한다.

3기에서의 고온피해온도는 남부에서 드러났는데, 이는 등숙기 하의 고온조건이 수량성에 부정적인 영향을 가진다는 선행연구(Lee et al., 2013)와 일치한다. 또한 냉해피해가 북부와 남부에서 관측된 온도구간(6℃ ~ 9℃, 5℃ ~ 8℃) 역시 선행연구의 등숙기 냉해피해 임계온도인 14℃(RDA, 2014) 이하로, 농학분야 연구와 일관적이다.

4. 고찰

본 연구를 통해 드러난 지역적 특성에 따른 논벼 수량성-기온변수 간의 관계에서의 차이점은 지역성 반영의 필요성을 입증한다. 또한, 대부분의 지역에서 생육주기별로 나타난 상이한 수량성-기온 관계는 공간적 이질성과 함께 시간적 이질성 또한 감안해야 함을 보여준다. 특히 남부지역의 경우 북·중부지역과 달리 고온피해의 발생이 실제 통계를 통하여 확인되고 있어, 향후 기후변화 조건 하 적응을 위한 정책의 도입이 필요할 것으로 보인다.

선행연구(Schlenker and Roberts, 2009)는 기후인자와 옥수수, 대두, 목화 수량성 간의 상관성 분석에 고정효과 패널모형을 활용하였으나, 본 연구에서는 Hausman 검정결과 시군구별 미관측 개별특성과 기상데이터 간의 상관성 및 독립성이 타당하지 않은 것으로 드러나 임의효과 모형을 활용하였다. 이와 같은 차이는 본 연구가 선행연구에 비하여 짧은 시간(56년 대 18년)에 걸쳐 진행되었다는 점에서 설명가능하다. 예컨대 개별 시군구의 논벼 수량성에 기여하는 미관측 특성으로 지형과 도시 인접성이 기여한다면, 후자는 18년에 걸친 분석에서만 불변 요소로 판단될 수 있다. 그런데 PRISM 데이터는 두 요소 중 지형적 요소만을 직접적으로 반영하므로, 이와 같은 사례에서는 미관측 불변요소가 장기에서는 독립변수와의 종속성을, 단기에는 독립성을 보일 수 있다.⁶⁾

기상 데이터의 제한성은 고온피해의 강도를 추정하는 데 한계를 유발한다. 냉해피해의 경우 현재 남한 논벼 재배지, 특히 중·북부지역의 논벼 재배지에서 장해온도를 경험할 확률이 높기에 비교적 뚜렷하게 확인되는 반면, 고온피해는 남부지역에서만 명확하게 확인되었다. 수원지역의 화성벼와 다산벼에 대해 대기온도 대비 3.0℃의 온도 상승 조건에서야 논벼 수량에 대한 유의한 영향이 나타난 연구결과(Lee et al., 2015)를 고려하면, 온도구배 하우스나 작물생육모형이 아닌 실제 집계된 논벼 수량성 통계 및 실측기상자료에서 고온피해의 강도를 정확히 예측하기 위해서는 향후 추가적인 기상데이터, 특히 극한고온조건 하에서의 기상데이터 축적이 필요함을 알 수 있다. 2018년은 전국 각지에서 역대 최고 수준의 온도가 기록되었던 해인만큼, 향후 연구에서 해당 해의 데이터 제공 및 활용이 요구된다.

기상 데이터 활용의 측면에서, 본 연구는 광범위한 연구대상 지역으로 인하여 1 km 격자 기상 데이터를 활용하였으나, 현재 농업재해 조기경보시스템은 30 m 해상도 격자 기상데이터를 활용하여 기상재해를 예보 중에 있다(Ahn et al., 2014). 향후 연구에서는 지형적 특성과 미기상적 특성을 고려한 분석을 위해 필지단위 격자 데이터를 활용하는 것 또한 유의미할 것으로 예상된다.

또한, 본 연구는 시군별 논벼의 생육주기를 고정하였으나, 실제 논벼의 생육주기는 이앙시기와 각 연도별 기후조건에 따라 신축적으로 변화한다. 해당 연구방법이 향후 작물생육모형과 결합된다면 분석의 정밀도를 향상시킬 수 있을 것으로 예상한다.

기후변화 정책의 관점에서, 본 연구의 추세 및 추세제곱항은 주어진 기상조건(기온 및 강우)과는 무관하게 기술적인 요소가 수량성에 미치는 영향을 분석하므로, 넓은 의미에서 기후변화 조건에 대한 기술적 대처역량을 보인다고 할 수 있다. 이와 같은 추세선이 중부지역에서는 2008년 이후 상승하는 형태를, 남부지역에서는 2013년 이후 하락하는 형태를 보였음은 향후 연구의 분석대상이 되어야 할 것으로 보인다. 실제 2000년 ~ 2017년 사이 「농업생산비조사」의 ‘시도별 논벼 농가 생산비 조사’ 항목 중 ‘농구비’ 항목이 물가조정 이후에도 두 지역에서 상당히 다른 흐름을 보인 점(Fig. 3)은 이를 설명할 수 있는 단서가 될 수 있을 것으로 보인다. 다만, 이와 같은 분석이 단순한 상관관계를 넘어 인과관계를 규명하는 분석이 되기 위해서는 추가적인 정교한 계량경제학적 접근이 요구된다.

Fig. 3.

Cost of farming equipment per unit area in the Southern and the Central region, showing distinctly different trends

끝으로, 본 연구는 단순한 기온 및 강수량 뿐만이 아닌, 집중호우 피해면적 및 강도 또한 독립변수로 포함하였다. 기상데이터를 이용한 패널분석에 있어서는 시공간적 상관성을 반영하는 한에서 종속변수와 연관성을 가질 수 있는 데이터를 최대한 활용하는 것이 필요하므로(Auffhammer et al., 2013) 향후 연구에서는 통계적 보정을 거친 1 km 격자단위의 풍속 데이터(Kim et al., 2011), 일사량 데이터 등을 활용할 경우 정밀도를 개선할 수 있을 것으로 예상된다. 이와 같이 정교화된 모형과 기후변화 시나리오 하 격자 단위 기상자료의 결합은 향후 시군구 단위 논벼 생산성을 추정하는 데에 활용될 수 있을 것이다.

5. 요약 및 결론

본 연구는 기초지자체 단위에서 농업분야의 기후변화 취약성을 평가하기 위한 수단으로서 미관측 효과 패널모형의 타당성을 확인하였다. 지역성 고려의 필요성을 확인하기 위하여 분석대상 지역을 3개 농업기후지대의 27개 시군구로 설정하였으며, 지역별로 온도구간에 차이를 두고 각 온도구간에서의 누적시간이 수량성과 형성하는 관계에 대한 분석을 실시하였다. 각 모형 내에서 시군구 특성의 임의성과 오차항의 시·공간적 자기상관성을 검정하였다.

연구결과 모든 기후지대에서 임의효과 모형이 선정되었으며, 기온조건에 대하여 지역 간 상당한 수준의 반응 차이가 관측되었다. 특히 남부지역은 모든 생육주기에 걸쳐 고온에 대한 취약성을 보였으며, 저온피해 취약성은 중부와 북부 분얼 활착기에서 주로 발견되었다. 극한기상은 전반적으로 모형의 설명력을 개선하였으며, 집중호우 피해면적과 집중호우 강도가 유의하였다.

향후 연구에 있어서는 포함된 변수와 데이터의 시공간적 개선을 활용한 분석이 요구된다. 고온조건 하에서의 고해상도 기상데이터의 축적은 향후 기후변화 조건 하에서 각 시군구 농업부문의 취약성을 분석하는 데에 유의미한 기여를 제공할 것으로 기대된다.

Acknowledgments

본 연구는 환경부의 재원으로 환경산업기술원의 기후변화대응사업 지원을 받아 수행되었습니다(과제번호: 201800131002).

References

• Agricultural Geographical Information System. [accessed 2020 Sep 1]. https://agis.epis.or.kr
• Ahn JB, Hur J, Lim AY. 2014. Estimation of Fine-Scale Daily Temperature with 30m-Resolution Using PRISM. Atmos Kor Meteorol Soc 24(1): 101-110. [https://doi.org/10.14191/Atmos.2014.24.1.101]
• Auffhammer M, Hsiang SM, Schlenker W, Sobel A. 2013. Using Weather Data and Climate Model Output in Economic Analyses of Climate Change. Rev Environ Econ Policy 7(2): 181-198. [https://doi.org/10.1093/reep/ret016]
• Blanc E, Schlenker W. 2017. The Use of Panel Models in Assessments of Climate Impacts on Agriculture. Rev Environ Econ Policy 11(2):258-279. [https://doi.org/10.1093/reep/rex016]
• Breusch TS. 1978. Testing for autocorrelation in dynamic linear models. Australian Econ Pap 17(31): 334-355. [https://doi.org/10.1111/j.1467-8454.1978.tb00635.x]
• Choi DH, Yun SH. 1989. Agroclimatic Zone and Characters of the Area Subject to Climatic Disaster in Korea(in Korean with English abstract). Kor J Crop Sci 34(2): 13-33.
• Choi SK, Kim MK, Jeong J, Choi D, Hur SO. 2017. Estimation of Crop Yield and Evapotranspiration in Paddy Rice with Climate Change using APEX-Paddy Model(in Korean with English abstract). J Kor Soc Agric Eng 59(4): 27-42.
• D’Agostino AL, Schlenker W. 2016. Recent weather fluctuations and agricultural yields: implications for climate change. Agric Econ 47 supplement: 159-171. [https://doi.org/10.1111/agec.12315]
• Driscoll JC, Kraay AC. 1998. Consistent covariance matrix estimation with spatially dependent panel data. Rev Econ Stat 80(4): 549-560. [https://doi.org/10.1162/003465398557825]
• Godfrey LG. 1978. Testing against general autoregressive and moving average error models when the regressors include lagged dependent variables. Econometrica 46(6): 1293-1301. [https://doi.org/10.2307/1913829]
• Han JH, IC Son, IM Choi, SH Kim, JG Cho, SK Yun, HC Kim, TC Kim. 2011. Predicting harvest date of ‘Niitaka’ pear by using full bloom date and growing season weather. Kor J Hort Sci Technol 29(6):549-554.
• Hyun SW, YM Park, MS Ko, YK Kang. 1994. Effect of seeding date on growth and yield in oats. Kor J Crop Sci 39:359-365.
• Jong SK, SS Lee, KY Park. 1986. Methods of estimating growing degree days to predict growth duration in maize. Kor J Crop Sci 31:186-194.
• Kim DS, Song J, Lee JI, Chun A, Jeong EG, Kim JT, Hur OS, Kim SL, Suh SJ, 2009, Rice Quality Characterization According to Damaged Low Temperature in Rice Plant, Kor J Crop Sci 54(4): 452-457.
• Kim HJ, Kim HS, Choi YJ, Lee SW, Seo BK. 2011. A Statistical Tuning Method to Improve the Accuracy of 1km × 1km Resolution-Wind Data of South Korea Generated from a Numerical Meteorological Model (in Korean with English abstract). Kor J Appl Stat 24(6): 1225-1235. [https://doi.org/10.5351/KJAS.2011.24.6.1225]
• Kim HY, Ko J, Kang S, Tenhunen J. 2013. Impacts of climate change on paddy rice yield in a temperate climate. Glob Change Biol 19:548-562. [https://doi.org/10.1111/gcb.12047]
• Kim MK, Han MS, Jang DH, Baek SG. 2012. Production Technique of Observation Grid Data of 1km Resolution (in Korean with English abstract). J Clim Res 7(1): 55-68.
• Korea Meteorological Administration, 2020. Climate Data Portal. [accessed 2020 June 23]. http://www.climate.go.kr/
• Korea Seed & Variety Service (KSVS). 2020. 2020 Government-provided cultivar: names and characteristics. Gimcheon, Korea: KSVS.
• Korea Meteorological Administration, APEC Climate Center, Ministry of Agriculture, Food and Rural Affairs, Rural Development Administration and 18 others. 2020. 2019 Extreme Weather Report. Seoul, Korea: KMA.
• Lee SS. 1983. Utilization of growing degree days as an index of growth duration of rice varieties. Kor J Crop Sci 28:173-183.
• Lee CK, Son JY, Gu B, Mo YJ. 2013. Investigation of rice yield and quality change and development of adjustment technique to climate change (in Korean with English abstract). Jeonju, Korea: Rural Development Administration. Research Report.
• Lee KJ. Nguyen DN, Choi DH, Ban HY, Lee BW. 2015. Effects of Elevated Air Temperature on Yield and Yield Components of Rice. Kor J Agric For Meteorol 17(2):156-164. [https://doi.org/10.5532/KJAFM.2015.17.2.156]
• Liu B, Asseng S, Müler C, Ewert F, Elliott J, Lobell DB and 56 others. 2016. Similar estimates of temperature impacts on global wheat yield by three independent methods. Nat Clim Change 6:1130-1136. [https://doi.org/10.1038/nclimate3115]
• Mendelsohn R, Nordhaus WD, Shaw D. 1994. The Impact of Global Warming on Agriculture: A Ricardian Analysis. Am. Econ. Rev. 84(4): 753-771.
• Ministry for Agriculture and Forestry(MAF). 2006. Prediction of harvesting time of super sweet and waxy corn hybrids with growing degree days. Gwacheon, Korea: MAF. p. 49-62.
• National Institute of Crop Science(NICS). 2004. Samgwangbyeo cultivar report(in Korean). Jeonju, Korea: NICS.
• National Institute of Crop Science(NICS). 2006. Odae1-ho cultivar report(in Korean). Jeonju, Korea: NICS.
• National Institute of Crop Science(NICS). 2012. Saeilmi cultivar report(in Korean). Jeonju, Korea: NICS.
• Pesaran M. 2004. General Diagnostic Tests for Cross Section Dependence in Panels. CESifo Work Pap, 1229.
• Rural Development Administration. 2014. Essential Rice Production Techniques: Questions and Answers (in Korean). Jeonju, Korea: Rural Development Administration.
• Schlenker W, Roberts MJ. 2009. Nonlinear Temperature Effects Indicate Severe Damage to U.S. Crop Yields under Climate Change. PNAS 106(37): 15594-15598. [https://doi.org/10.1073/pnas.0906865106]
• Snyder RL. 1985. Hand Calculating Degree Days. Agric For Meteorol 35: 353-358. [https://doi.org/10.1016/0168-1923(85)90095-4]
• Statistics Korea, 2019. Crop production survey report(in Korean). Daejeon, Korea: Statistics Korea.
• White H. 1984. Asymptotic Theory for Econometricians. New York: Academic Press.
• Wooldridge JM. 2010. Econometric Analysis of Cross Section and Panel Data. Cambridge, MA: MIT Press.
• Yoshida S. 1973. Effects of temperature on growth of rice plant(Oryza sativa L.) in a controlled environment. Soil Sci Plant Nutr 19(4):299-310. [https://doi.org/10.1080/00380768.1973.10432599]