고자리꽃파리는 양파 및 마늘 등 수선화과(Amaryllidaceae) 부추속(Allium)에 속하는 농작물에 중요한 해충으로(Finch and Eckenrode, 1985;Collier et al., 2020), 전 세계적으로 온대 지역에서 경제적 해충으로 취급하고 있다. 방제약제를 처리하지 않는 경우 62.9%에서 거의 100%까지 달할 수 있다(Nault et al., 2006a;Walters and Eekenrode, 1996;Moretti and Nault, 2020).
특히 고자리꽃파리의 피해는 1세대 발생기엔 마늘이나 양파의 재식 초기에 높은 것으로 보고되었다(Werling, 2004). 이것은 월동 세대의 발생 시기 예측이 초기 방제를 위하여 매우 중요하다는 것을 나타낸다. 따라서 1세대를 겨냥한 약제살포 시기를 결정하기 위해 봄 세대 발생 예측 모형이 개발되었다(미국: Eckenrode et al., 1975, 캐나다: Liu et al., 1982, 노르웨이: Thöming, et al., 2019). 우리나라에서는 주로 마늘이나 양파를 가을철 파종하기 때문에 상황이 조금 다를 수 있으나, 봄철 고자리꽃파리 성충발생을 예측하여 적기에 방제하는 것도 연간 방제체계의 중요한 요소라고 할 수 있다.
Masaki (1980)는 휴면여부와 번식기를 기준으로 곤충의 생활사 유형을 A~E형으로 구분하면서 동면(월동휴면)과 하면(여름휴면)의 행동양식(syndrome) 동시에 가지고 있는 꽃파리과(Anthomyiidae) 해충을 D-형으로 구분하고, 가을철에만 번식하는 종은 Da-형, 봄철 번식 종은 Db-형, 그리고 봄과 가을에 모두 번식하는 종은 Dc-형으로 세분화하였다. 이 중 고자리꽃파리는 Dc-형에 해당되며, 초여름 번데기 형성 초기에 24℃ 이상의 온도에 자극을 받아서 하면에 들어가는 것으로 알려져 있다 (Ishikawa et al., 2000). 국내에서 Park et al. (1990)에 따르면 진주지역에서 고자리파리는 5월 상순 나타난 번데기는 43.8%가 하면으로 들어가고 6월 상순에 생성된 번데기는 100% 하면에 들어가는 것으로 나타났다. 하면에 들어간 번데기는 22.1°C 이하 온도에서 일정기간 이상의 자극을 받아야 하면에서 깨어나서 성충으로 우화한다(Ishikawa et al., 2000;Shin et al., 2023). 따라서 우리나라와 같은 전형적인 온대지역에서 고자리꽃파리는 봄에서 초여름까지 1~2세대를 경과하여 1차 발생을 보이고(전기세대), 여름철 발생이 중단되었다가 가을철 9월 하순에서 10월 2차 발생(3세대, 즉 후기세대)하는 양상을 보인다(Choi, 1981;Park et al., 1988;Park, 1989). 하지만 제주지역과 같이 여름철부터 늦가을까지 고온이 지속되는 경우 하면번데기의 하면타파가 늦어져서 우화시기가 크게 지연되고, 소수만이 우화하기 때문에 후기세대의 발생이 미약한 것으로 나타났다(Shin et al., 2023). 이러한 현상은 강한 고온효과로 설명되었으며, 고온효과가 약한 미국이나 독일 북부지역의 경우 전기세대와 후기세대가 크게 이격되지 않고 발생하고, 우리나라 진주와 같은 지역은 중도효과 지역으로써 전기와 후기세대가 뚜렷하게 분리되어 피크를 형성하는 특징이 있었다.
이렇게 따뜻한 아열대의 제주지역에서는 후기발생이 미약하여 자연히 월동밀도 형성이 낮아짐으로 다음해 봄 월동세대의 최성기가 명확하지 않을 수 있다. 신뢰할 수 있는 월동번데기 우 화모형을 활용할 수 있다면, 고자리꽃파리의 월동후 방제적기의 결정뿐만 아니라 우화예측시기를 기준으로 봄철부터 여름철까지 발생양상을 해석하는 데 유용하게 활용할 수 있을 것이다.
고자리꽃파리는 휴면상태의 번데기로 월동하며 온도(저온) 및 광주기에 의하여 조절되고, 50% 휴면유도 온도는 광주기 16L:8D에서 14.0°C, 12L:12D에서는 18.5°C로 같은 온도의 경우 짧은 낮 길이에서 휴면율이 높아진다(Ishikawa et al., 1987). 휴면상태의 번데기가 휴면을 종료하려면 4~15°C에서 130일이 필요한 것으로 알려져 있다(Nomura, unpublished data; Reviewed in Ishikawa et al., 2000). 국내(진주)에서 고자리꽃파리 월동휴면 개시는 10월 중하순부터 시작되어 해에 따라 11월 하순에서 12월 하순 최성기에 보였다(Park et al., 1987). 또한 50% 휴면종료(타파) 시기는 10월 상순 형성된 번데기는 12월 중순이었고, 11월 하순에서 12월 하순 사이 형성된 번데기는 다음해 1월 상순 경이었다. 진주지역에서 100% 휴면타파 완료 시기는 휴면개시와는 상관없이 1월 중순이었다.
따라서 본 연구는 고자리꽃파리 월동 번데기는 1월 상순에 휴면이 대부분 종료되고 발육이 가능한 온도 조건이 조성되면 언제든지 휴면 후 발육(즉 번데기의 성충으로 발육)을 개시할 수 있다고 전제하고, 고자리꽃파리 월동 번데기 온도발육 자료를 이용하여 온도와 발육과의 관계를 구명하고 우화시기를 예측할 수 있는 모형을 개발하고자 수행하였다.
재료 및 방법
월동 번데기 발육자료
모형개발에 사용한 자료는 1986년 1월 17일 진주(N35°12’) 지역에서 채집하여 가온 사육한 것이었다(Park et al., 1987 농시논문집). 채집한 월동번데기를 축축한 모래가 담긴 페트리디쉬(9 × 3 cm)에 1 cm 깊이로 묻고 항온기에 처리하여 매일 성충 우화여부를 조사하였다. 이렇게 하여 얻은 가온 사육일부터 성충 우화까지 기간(암수 평균발육기간)은 10°C 39.6일, 15°C 22.0일, 20°C 13.5일, 25°C 10.2일, 30°C 9.5일 이었다(n = 60).
온도발육모형
선형모형
번데기 발육기간에 역수를 취하여 발육률(1/평균발육기간, 즉 발육속도)로 간주하였으며, 이하 모든 매개변수의 추정은 TableCurve2D 프로그램(Jandel Scientific, 2002)을 사용하였다.
선형모형을 적용하여 온도(x)와 고자리꽃파리 월동용(번데기)의 발육률(y)과의 관계를 추정하였다. 이 선형모형의 매개변수를 이용하여 Campbell et al. (1974)이 제시한 방법에 따라 발육영점온도는 -b/a, 발육 완료에 필요한 적산온도(Degree days, 즉 온량상수 thermal constant)는 1/a의 계산식으로 추정하였다.
여기서 y는 발육률, x는 온도를 나타내며, a는 기울기이고 b는 절편이다. 선형모형 추정시 선형구간에서 벗어난 30°C 자료른 분석에서 제외하였다. 발육영점온도와 온량상수의 표준오차는 Campbell et al. (1974)의 방법에 따라 계산하였다.
비선형 모형
온도에 따른 발육률의 비선형적 관계는 Lactin et al. (1995) 의 방정식을 변형하여 적용하였다. 이들은 기본 3개의 매개변수 모형과 곡선이 온도 축(즉 x-축)을 교차하도록 매개변수 하나를 추가한 4개의 매개변수 모형을 제시하였다. 상기 발육자료를 이용하여 추정한 3개 매개변수 모형의 곡선이 수식의 속성상 저온에서 발육률이 과대평가되었기 때문에 4번째 매개변수(λ)를 선형모형의 절편인 b값으로 고정시키고 나머지 매개변수를 추정하여 비교하였다.
여기서 r(T)는 대기온도 T °C에서 발육률이며, Tmax = 발육 상한온도(상한임계온도)을 의미하는 매개변수로써 더 이상 생명 현상을 지속할 수 없는 일종의 치사온도, ΔT = 온도상한 부근에서 급격한 발육붕괴가 일어나는 온도 범위(발육붕괴 온도범위)(the temperature range over which thermal breakdown became the overriding influence), ρ = 곤충발육에 관여하는 주요 효소가 촉매하는 생화학적 반응에 대한 합성치(the composite value for critical enzyme-catalyzed biochemical reactions in insect development)로 정의할 수 있다(Logan et al., 1976;Damos and Savopoulou-Soultani, 2008). 마지막 4번째 매개변수 b값은 앞에서 정의했듯이 선형모형의 절편을 나타낸다.
발육기간 분포모형
고자리꽃파리 번데기가 발육을 완료하는 확률분포를 작성하기 위해서는 발육기간별 발육을 완료하는 개체수(빈도수) 자료가 필요하다. 그러나 본 연구에서는 발육기간의 원자료를 이용할 수 없었기 때문에 평균과 표준편차를 이용하여 분포모형의 매개변수를 추정하는 방법(Kim and Kim, 2018;Kim and Lee, 2022)을 적용하였다. 이 방법은 발육기간의 빈도분포가 정규분포를 따른다고 가정하고 평균과 표준편차를 이용하여 표준누적 확률분포도를 작성하는 방법이다. Park et al. (1987)이 보고한 자료의 월동번데기 평균발육기간과 표준오차 자료를 이용하다. 암컷과 수컷의 자료를 분리하여 표준오차를 표본수에 따라 표준편차로 환산하고 표준연령을 계산하였으며(Supp. Sheet 1, Online available), 최종 두 자료를 종합하여 Weibull 함수에 적용시켜 매개변수를 추정하였다.
여기서 f(px)는 표준연령(또는 생리적 연령) px까지 누적발육 완료확률이며, α, β는 추정해야 할 매개변수이다.
고자리꽃파리 월동 번데기의 성충우화모형
적산온도 모형
선형모형에서 얻은 발육 영점온도와 온량상수를 이용하여 50% 성충 우화시기를 추정하였다. 즉 야외 온도자료를 이용하여 적산 온도가 온량상수에 도달하는 날을 50% 우화시기로 간주하였다.
발육률 적산모형
발육률 적산모형은 비선형 발육모형과 발육기간 분포모형을 융합한 방법으로 발육단계전이모형(Kim et al., 2001)과 동일한 구조를 갖는다. 즉 번데기의 생리적 연령 i와 i+1 사이 동안 성충으로 전이되는 개체수는 발육기간 분포모형을 이용하여 다음과 같이 계산하였다.
여기서 번데기의 생리적 연령은 로써 번데기 온도발 육모형의 출력값을 누적하여 계산하였다.
모형 시뮬레이션 및 포장적합
포장 적합자료의 수집
과거 우리나라 진지주역에서 조사된 자료를 수집하여 모형 적합에 사용하였다. 월동 번데기 집단의 우화 자료는 Park et al. (1987)이 보고한 1984년 및 1987년 자료를 이용하였다. 1984년 자료는 진주에서 3월 6일 채집하여 야외 포장 토양 5~15 cm 깊이에 매몰하고 조사한 자료이었다(n=99). 1987년 자료는 진주에서 1월 17일 채집하여 야외 포장 토양 5 cm 깊이에 매몰하고 조사한 자료이었다(n=210, 처리 개체수 기준).
또한 진주지역에서 1985년과 1986년 야외 채집한 월동번데기를 5cm 깊이 토양에 매몰하고 조사하여 얻은 50% 우화시기를 포장 적합에 이용하였다(Park et al., 1987). 1985년은 4월 5일, 1986년은 각각 4월 6일이었다.
고자리꽃파리 연간 발생양상에서 월동 번데기 우화시기의 위치를 이해하기 위하여 2021년(3월 하순~10월 중순)과 2022년(3월 상순~12월 하순) 제주 대정지역에서 수집된 자료(Shin et al., 2023)와 1987년 진주자료(Park et al., 1988)를 이용하여 비교하였다. 제주 자료는 백색 점착판 트랩(100 × 250 mm, 영국 Russell IPM)을 원통형(높이 100 mm, 직경 76 mm)으로 말아서 철제 고추 지지대를 사용하여 지표면에서 1 m 높이에 수직으로 설치한 트랩을 사용하였다(2021년 2개 지점, 2022년 3개 지점, 각 1주일 간격조사). 진주 자료는 원추형 트랩(즉, 콘형 트랩 밑면직경 28 cm, 높이 31.5 cm)을 이용하여 1987년 3월 하순부터 11월 중순까지 10일 간격으로 조사한 것이었다(트랩수 = 10 개). 이 자료를 이용하여 월동번데기 우화 예측시기 피크에서 공통발육영점온도(3.9°C) 기반의 알-알세대기간 온량상수를 적용하여(Shin et al., 2023) 이후의 세대 발달을 추적하였다. 고자리꽃파리의 각 발육단계 완료를 위한 DD(즉, 온량상수)는 Shin et al. (2023)이 제시한 바와 같이 알 53.6DD, 1령 유충 79.3DD, 2령 유충 73.6DD, 3령 유충 159.2DD, 번데기 217.0DD를 적용 하였다. 산란전기간 완료에 필요한 적산온도는 20°C에서 6.5일 (양파), 7.56일(부추), 8.53일(마늘)(Ning et al., 2017), 22°C에서 6.4일(Spencer and Miller, 2002) 값으로 재계산하여 119.9 DD를 얻었고, 결과적으로 고자리꽃파리의 알에서 알까지의 세대기간으로 702.6DD를 적용하였다.
모형 시뮬레이션 및 포장적합
선형 온도발육모형을 바탕으로 한 적산온도 모형의 경우 적산온도는 일일 평균온도를 이용한 단순 계산법과 일일 최고 최저온도를 이용하는 사인곡선법(수평 제거법)을 적용하여 계산하였다. 사인곡선법은 시간별 온도를 반영하는 성격을 가지고 있다. 즉 하루 최고, 최저온도를 이용하여 사인곡선을 만들고 발육영점온도와 상한온도 범위의 온도만을 누적하는 방법이다(Higley et al., 1986). 상한온도로 월동번데기 비선형발육모형의 매개변수 Tm을 이용하였다. 계산은 Baskerville and Emin (1969), Zalom et al. (1983) 등이 제시한 수식을 바탕으로 제작한 엑셀시트(DegDay.xls, Version 1.01: Snyder, 2002)를 이용 하였다(Supp. Sheet 2, Online available). 실제 관측자료(실측치)와 모형 예측자료(예측치) 모두 누적분포도를 작성하고 50% 우화시기를 비교하였다.
발육률 적산모형은 선형 또는 비선형 발육모형을 이용하여 포장온도를 발육율로 변환시키고 이를 누적하여 발육기간 분포 모형에 입력하여 우화시기를 예측하였다. 온도자료는 일일 평균온도와 시간별 온도를 추정하여 이용하여 비교하였다. 즉 시간별 온도자료는 하루 최고와 최저온도를 이용하여 사인곡선으로 시간별 온도를 추정하는 방법(Supp. Sheet 3, Online available)을 적용하였다(Jung et al., 2023). 이 사인곡선식에서 최고와 최저 온도는각각 15:00, 06:00에 나타난다고 가정하였다.
모형의 정확도는 실측치 및 예측치의 10, 25, 50, 75, 90% 누 적발생기에 해당되는 날짜(Julian date)를 추출하여 차이를 검정 하였다(5분위 평균값, five quantile mean). 즉 차이가 3일 또는 5일과 차이가 없다는 귀무가설 하에 단일평균 t-검정(one sample t-test)을 실시하였다. 포장 관측일이 각 퍼센트 발생기와 일치하지 않는 경우 두 관측시기간 내삽을 통하여 결정하였다.
유의성이 없는 경우는 정해진 3일 또는 5일과 차이가 없다고 평가하였으며, 유의성이 있는 경우 평균편차가 3일 또는 5일보다 작으면 이들 값보다 편차가 통계적으로 유의하게 작다고 판단하였다. 반면, 유의성이 있는 상태에서 평균편차가 3일 또는 5일보다 크면 차이가 이보다 크다고 평가하였다.
더불어 최종 검증된 월동번데기 우화시기 예측모형의 결과를 제주(2021, 2022년) 및 진주지역(1987년)에서 조사된 포장 발생시기와 비교하여 초기 페놀로지(계절발생)를 해석하는데 활용하였다.
기상자료
기상자료(일별 최고 및 최저온도)는 기상청(the Korea Meteorological Administration)의 종관관측 자료를 이용하였으며, 트랩조사 위치에서 가장 가까운 진주 기상대 자료, 제주는 고산 기상대 자료를 사용하였다. 기상청에서 제공하는 토양온도는 농장의 온도를 대표하기 어렵기 때문에 대기온도를 이용하여 실제 농경지에서 경험적으로 추정한 보정식으로 계산하였다. 즉 대기온도(Tair)와 토양온도(Tsoil, 5 cm 깊이)는 다음과 같은 관계가 있었다(Choi et al., 2019).
여기서 Taυe.air 는 대기 일평균온도를 의미한다.
결과
모형 매개변수 추정
선형모형
추정 된 회귀모형은 통계적 유의성이 있었으며(F = 1151.51; df = 1,3; P < 0.0001; R2 = 0.99), 온도와 발육률 사이에 선형회귀 선(Fig. 1A)을 사용하여 월동번데기의 발육영점온도는 5.2°C, 온량상수는 202.4DD (degree-days, DD)로 추정되었다(Table 1).
비선형 모형
고자리꽃파리 번데기 발육률과 온도 사이의 관계는 Lactin et al. (1995) 3-매개변수 모형에 잘 적합되었으며, 결정계수(R2)는 0.99로 변이에 대한 설명력이 높았다(Fig. 1B). 또한 회귀모형은 통계적으로 유의하였다(F = 298.722; df = 2,2; P = 0.0033, R2 = 0.99). 추정된 곡선은 최적 온도까지 온도 증가에 따라 증가하였고 그 후 높은 온도 영역에서 감소하였다. 고온 임계온도(매개 변수 Tm)는 34.8°C로 추정되었다(Table 1). 선형모형의 절편값으로 매개변수 하나를 고정시킨 4-매개변수 모형 또한 통계적으로 유의하였다(F =1242.31; df 2,2; P = 0.0008, R2 = 0.99).
저온 영역인 10°C 이하 구간에서 4-매개변수 발육률 곡선은 3-매개변수 곡선보다 선형모형에 더 근접하였고, 25°C 이상에서는 더 빨리 하강하였다.
발육기간 분포모형
매개변수 2개의 Weibull 식을 이용하여 번데기 발육기간 분포모형의 매개변수를 추정하였으며(Table 1), 추정된 회귀선은 통계적으로 유의하였다(F = 8852.0; df = 1,88; P < 0.0001, R2 = 0.99). 분포모형의 기울기를 나타내는 매개변수(β)는 11.3689이었고, 10, 25, 50, 75, 90% 누적우화율에서 생리적 연령은 각 각 0.8525, 0.9312, 1.0061, 1.0694, 1.1182에 대응되었다.
모형의 포장 적합
고자리꽃파리 월동번데기의 성충우화시기 50% 예측에서 각 예측모형의 정확도는 적용한 발육모형의 종류와 입력온도의 형태에 따라서 차이가 있었다(Table 2). 전체적으로 일일 평균온도를 입력한 경우 시간별 온도를 사용한 것과 비교하여 실제 관측일과 모형 예측일의 차이가 증가하였다. 온도와 발육률 간 선형 발육모형(선형모형)에서 얻은 발육영점온도와 온량상수(적산온도)를 이용한 적산온도 모형은 오차(차이)가 평균 9.8일로 3일 또는 5일보다 통계적으로 유의하게 큰 차이가 발생하였다. 반면, 일일 평균온도를 이용한 선형모형 발육률 적산모형의 예측 정확도는 3일 또는 5일과 차이가 없었으나, 실제 차이의 평균은 7.9일로 연도간 변이가 매우 심하였다. 일일 평균온도와 비선 형모형을 이용한 발육률 적산모형의 경우 3-매개변수 모형의 예측일 편차는 실측일과 평균 15.5일로 매우 차이가 있었다. 4- 매개변수 모형의 경우는 평균 5.9일 편차로 3일 또는 5일과 차이가 없었다.
시간별 온도를 입력값으로 한 경우 3-매개변수 발육모형을 이용한 경우 모형 예측치와 실측치의 편차가 20.4일로 큰 차이가 있었다. 그러나 적산온도 모형(사인곡선법 적용), 선형모형 기반의 발육률 적산모형, 4-매개변수 비선형 모형 기반의 발육률 적산모형은 모두 실측일과 예측일의 평균편차가 3일과 차이가 없었다. 즉 각각 1.5일, 1.8일, 및 2.5일의 작은 차이가 있었다.
1984년 및 1987년 고자리꽃파리 월동번데기 우화집단 자료를 이용하여 누적우화율 0~100%까지 모형의 누적우화곡선과 비교한 결과 Fig. 2와 같았고, 그 예측정확도 통계량은 Table 3 에 제시하였다(Table 2의 50% 예측에서 정확도가 높았던 시간별 온도를 이용한 선형모형과 4-매개변수 비선형 모형 기반의 발육률 적산모형만 제시되었음). 선형 발육률 적산모형의 경우 10, 25, 50, 75, 90% 누적우화율에서 평균편차(즉 관측치-예측치)는 1984년 3.1일, 1987년 2.5일로 3일과 차이가 없었다. 비선형모형(4-매개변수) 발육률 적산모형은 1984년은 평균편차가 0.8일로 3일보다 유의하게 작았고, 1987년에는 비록 평균편차가 6.5일로 다소 증가하였으나 3일과 통계적 차이가 없었다.
제주와 진주지역에서 연간 고자리꽃파리 성충 발생자료와 월동번데기 우화예측시기를 비교한 결과 Fig. 3과 같았다(월동 번데기 n = 1,000개 가정). 제주지역에서는 예측된 월동번데기 우화시기에 트랩의 성충 유살량이 매우 적었다(2022년 3월 18일 3마리 성충 포획). 반면, 진주에서는 월동번데기 우화시기와 첫 번째 피크 앞 부분이 일치하였고, 큰 발생피크를 형성하였다. 또한 고온기가 시작하기 전까지 진주지역에서는 부분적 2세대까지 발생이 가능하였으며, 제주지역에서는 알-알세대기간으로 추정했을 때 월동세대를 포함하여 3~4세대의 출현이 가능하였다.
고찰
기존 고자리꽃파리 번데기의 발육영점온도로 비휴면 상태의 경우 4.0°C를 적용한 바 있었다(Otto and Hommes, 2000). 외국에서 월동번데기의 성충우화 시기 추정 시 적용했던 발육영점 온도는 미국 뉴욕 4.44°C (Eckenrode et al., 1975), 미시건 4.0°C (Carruthers, 1979), 그리고 캐나다 4.38°C (reviewed in Liu et al., 1982) 그리고 노르웨이 3.9°C (Thöming et al., 2019)로 약간의 변이를 보였다. 본 연구에서는 5.2°C로 기존 보고의 값들과 비교하여 가장 높게 추정되었다. 결과에서 보여준 바와 같이 이 발육영점온도를 기준으로 적용했을 때 우리나라 환경에서는 좋은 예측력을 보였다.
본 연구에서는 선형모형과 비선형모형, 그리고 일평균온도 및 시간별 온도자료(사인곡선 추정) 등 다양한 조합으로 월동번데기 우화시기를 예측하여 비교하였다. 우선 발육모형의 적용에 있어서 선형모형을 이용한 경우 예측치와 관측치간 정확도가 가장 양호했다(발육률 적산모형 기준). 락틴 3-매개변수 모형은 우화시기를 매우 과대평가하였다. 그것은 선형모형으로 추정한 발육영점온도 이하에서도 발육률이 높게 유지되기 때문으로 해석할 수 있다(Fig. 1A 참조). 저온영역에서 발육률 곡선이 0으로 떨어지지 않는 문제를 개선하기 위하여 제시된 락틴 4- 매개변수 모형은 저온에서 발육의 과대평가 부분을 상당히 개선하였다 하더라도 여전히 선형모형과 비교하여 높은 발육률을 보이고 있으며, 이것이 월동번데기 발육을 과대평가하는 원인이 되는 것으로 보인다. 즉 선형모형의 절편을 락틴 모형의 4번 째 매개변수로 대체하였으나 저온영역에서 선형모형보다 높은 발육률을 나타내고 있고 이것이 과대평가의 원인으로 되었다. 이렇게 선형모형이 번데기 발육완료 예측에 정확도가 가장 높았다는 것은 온도에 따른 고자리꽃파리의 실제 월동번데기 발육이 선형적인 관계에 있었다는 것을 말한다. 보통 곤충의 발육은 저온영역에서 급격히 저온저해가 일어나거나 또는 느리게 0에 수렴하는 형태를 보이는데(Wagner et al., 1985), 고자리꽃파리 월동번데기는 선형적인 관계를 보이는 것으로 판단된다. 또한 발육기간분포 모형의 경우 정규분포를 가정하고 평균과 표준편차를 이용하여 추정하였기 때문에 실제 변이를 과소평가할 수도 있으나(상대적으로 큰 매개변수 β 값), 모형의 정확도에는 큰 영향을 미치지 않은 것으로 보인다.
모형의 정확도에 큰 영향을 미치는 요인은 온도자료를 어떻게 처리하는가 하는 것이었다. 시간별 온도자료를 이용했을 때와 일평균온도 자료를 사용했을 때 예측 정확도의 차이가 매우 컸다. 예를 들어 50% 우화기를 예측하는 선형모형 기반의 발육률 적산모형에서 전자는 평균편차가 1.8일이었으나 후자는 7.9 일로 편차가 매우 컸다. 이를 통해서 볼 때 주야간 외부온도가 발육영점온도 상에서 변동하고 하는 경우 시간별 온도자료를 이용해야 할 필요가 있다. 주야간 온도 변화가 심하고 만일 평균 온도가 발육영점온도 이하가 되면 그날의 발육률이 0가 되기 때 문에 실제 번데기의 발육을 과소평가하게 된다. 그러나 시간별 온도자료를 이용하게 되면 하루 중 발육영점온도 이상의 시간 동안 진행된 발육을 반영할 수 있다. 봄철 비교적 저온 시기 유효온도 또는 발육률을 누적할 때는 대상 해충의 발육을 과소평가할 수 있으므로(Jeon et al., 1996) 발육모형의 형태뿐만 아니라 시간별 온도자료를 이용한 것이 필요해 보인다.
모형의 적합에서 토양온도는 Choi et al. (2019)이 제안한 방법에 따라 대기 온도를 보정하여 사용하였다. 이 보정 방법은 실제 토양온도 측정자료(토양깊이 5 cm)를 이용하여 경험적으로 추정한 것으로 토양 해충 발생 예측에서 온도보정에 유용하게 활용될 수 있을 것으로 보인다(Shin et al., 2023). 이와 유사하게 검거세미밤나방 발생예측에서 단순히 대기온도에 0.8°C를 더하여 토양온도를 보정한 경우도 활용성이 있었다(Choi et al., 2023). 이러한 토양온도의 보정을 비롯하여 현장온도가 아닌 기상청의 온도자료의 이용(최고 및 최저온도)과 시간별 온도의 추정, 포장 관측일의 내삽 등은 모형 예측치와 실측치 간 편차를 발생시킬 수 있는 요인이 될 수 있다. 하지만, 적절한 발육모형의 선택과 시간별 온도자료의 적용을 통하여 예측정확도가 높은 결과를 얻을 수 있었다.
본 연구에서 제시한 선형발육모형과 시간별 온도자료를 적용한 고자리꽃파리 월동번데기 우화시기 예측모형은 실측치와 비교하여 평균편차가 3일과 차이가 없을 정도로 정확하였다. 이러한 정확한 예측 능력은 실제 포장상태에서 월동번데기의 성충으로 우화시기를 예측하여 방제체계를 수립하는데 유용할 것 이다. 캐나다 Keswick Marsh에서 적산온도를 이용한 고자리꽃파리 예측방제는 관행방제 연 7~12회 살포와 비교하여 연 2회로 방제가 가능할 정도로 예측모형은 효과적이었다(Liu et al., 1982). 또한 봄철 고자리꽃파리 첫 발생시기에 일반적으로 양파와 마늘의 경제적 피해가 가장 심하기 때문에 월동세대 우화예 측은 매우 필요하다(Nault et al., 2006a, b;Werling, 2004). 미국 뉴욕에서 Eckenrode et al. (1975)은 712DD 되는 시기를 첫 번째 성충 우화기의 50% 우화기로 예측하여 방제에 응용하였다. 즉 고자리꽃파리 발생세대에서 초기 밀도가 증가하기 시작하여 보통 7일 후 피해가 나타나므로 예측모형의 결과는 예방조치를 강구하는데 활용할 수 있다. 다만 Eckenrode et al. (1975) 의 모형은 50% 발생기만을 예측하므로 한계가 있지만, 본 연구의 예측모형은 0~100%까지 우화시기를 예측할 수 있기 때문에 더욱 폭넓게 응용이 가능하다. 이뿐만 아니라 고자리꽃파리의 봄철 페놀로지를 해석하는데 기준점으로 사용할 수 있다.
진주지역의 경우(Fig. 3C) 월동번데기 50% 우화 예측일은 4월 4일 이었고, 실측 포장 발생시기의 앞부분에 배치되었다. 포장 발생 최성기는 4월 19일로 모형보다 약 2주 후에 나타나고 있는데, 고자리꽃파리 성충 수명이 길어서 개체군밀도가 계속 누적되기 때문에 뒤쪽으로 지연되었을 것이라는 점으로 이해할 수 있다. 이렇게 월동세대(1세대) 피크가 앞쪽에 크게 나타나고, 반면 후속세대(2세대)는 고온기를 맞아 대부분 하면에 들어가게 되므로 알-알세대기간 도달지점에서 부분적으로 나타났다고 해석할 수 있다.
포장적합 통계량에서 보듯이 월동번데기 우화모형이 정확하다는 전제하에 제주지역에서 우화시기는 진주지역보다 훨씬 앞쪽에 위치하고 있다(Fig. 3A and 3B). 이렇게 포장에서 4월 중순부터 증가하기 시작한 개체군은 월동개체군(즉 월동번데기에서 우화)이 아니라 번식 개체군으로 해석할 수 있다. 예측된 월동번데기 우화시기에 실제 포장에서는 뚜렷한 발생피크가 형성되지 않고 있는 특징이 있다. 제주지역에서 또 하나의 특징은 5월~6월에 형성된 발생피크가 크게 나타나고 있다는 점이다. 이러한 현상은 기주식물인 마늘과 양파의 재배작형과 고자리꽃 파리의 하면기작을 감안하면 해석하기 난해한 부분이다. 즉 마늘과 양파는 빠르면 5월 상순부터 수확을 시작하여 6월에는 기주식물이 존재하지 않는다. 다만, 마늘 수확후 방치된 잔재물은 먹이원이 될 수 있으므로 부분적으로 6월 피크를 설명할 수 있다. 실제 포장에 버려진 비상품 마늘 잔재물 내지 친환경 재배지 같은 번식처가 존재하고 있는 것이 확인된다(Personal observation on 26 May, 2024). 이런 먹이원을 바탕으로 하면임계온도(24°C; Ishikawa et al., 2000) 도달 전에 번데기로 된 개체들이 존재하고 후에 우화하는 집단이 있는 것으로 추정된다. 또 하나의 가능성은 비래집단의 존재일 수 있다. 고자리꽃파리의 비래성에 대해서는 아직 보고가 없으나, 근연종인 씨고자리꽃파 리는 시리아 등 중동지역에서 유럽으로 비래하여 이주하는 꽃파리과(Anthomyiidae) 곤충의 대표종이라고 보고되었다(Hawkes et al., 2022). 향후 고자리꽃파리의 장거리 비래 가능성에 대하여 검토가 필요할 것이다.
제주지역에서 월동개체군의 크기가 매우 작은 것은 서론에서 언급했듯이 여름철 고온으로 하면타파가 지연되어 후기발생이 미약하다는 것으로 설명될 수 있다. 이렇게 기후적으로 따뜻 한 제주지역에서 고자리꽃파리 발생 양상은 마늘재배 작형과 관련하여 방제체계 수립에 유용하게 적용할 수 있다. 즉 제주에서는 가을 마늘 파종 후 유묘기 피해가 상대적으로 적게 되며, 난지형 마늘의 경우 4월 중하순 수확이 완료되기 때문에 고자리꽃파리 피해를 상당히 회피할 수 있다(일반적으로 5월 중순 수확). 그럼에도 불구하고 마늘재배에서 고자리꽃파리의 피해를 방지하기 위하여 예방적으로 토양살충제를 과도하게 살포하는 경향이 있다(Lee and Choi, 1999;Nault et al., 2006a, b;Moretti, 2020;Moretti and Nault, 2020). 앞으로 포장실험을 통하여 검증이 필요하지만, 장기간 토양 생태계에 악영향을 미칠 수 있는 토양살충제 대신 월동번데기에 우화하는 시기 성충을 표적으로 하는 선택적인 살충제 지상살포로 효과적인 방제를 달성할 수 있을 것이다.