한국형 말라리아의 수학 모형과 응용프로그램 개발

 Mathematical modelling on P. vivax malaria transmission and development of its application program

 질병관리본부 전염병대응센터 역학조사팀   
경북대학교 수학과  
 

Ⅰ. 들어가는 말
  우리나라의 말라리아는 1993년 파주지역 DMZ(Demilitarized Zone)에 복무 중이던 군인에서 재출현한 이후 매년 지속적으로 발생하고 있다. 발생 양상을 보면 2000년 이후 감소 추세에서 2005년 다시 증가하기 시작하였으나, 2008년 1,017건으로 전년도 동기간(‘07년 2,192건) 대비 53.6% 감소하였다.
  전염병 전파의 역동성을 조사하기 위한 수학 모형은 1760년 Bernoulli가 천연두(smallpox)의 수학 모형을 개발한 이후로 사용되어 왔다[1]. 모형은 현상들 사이의 관계를 표현하고 분석할 수 있도록 도와주고 특정한 사례들을 통합하고 확장하는 원리를 발견할 수 있도록 도와준다. 또한 눈으로 확인할 수 없는 기전을 추론할 수 있게 한다. 수학은 정신적인 것, 질적인 것을 대신할 수 없지만 구성 요소와 그들의 중간 연결, 가정과 목표가 정확한 형식적이고 양적인 범위로 확장할 수 있다. 그러나 수학 모형이 문제의 모든 특징을 다 포함할 수는 없다. 말라리아 모형의 유용성은 수학 모형의 수학적 결과가 데이터를 얼마나 정확하게 설명하고 예상하는 지에 의존한다. 모델링은 예측을 조절하고 현실성을 강화하여 애매함을 제거함으로써 의사결정, 의사의 요약과 전달, 실행을 관찰하고 그에 알맞게 전략을 수립할 수 있게 한다. 또한 모델링은 최적의 투자를 위한 수단을 제공하고, 실제적인 것과 이론적인 것 사이에서 접근하는 가능성들에 풍부한 범위를 제공할 수 있다. 그러나 말라리아 모형이 관찰이나 경험을 대신하지 않기 때문에 중재를 대신하지는 않는다.
  말라리아 발생은 환자 요인뿐만 아니라 기후, 모기 밀도 등과도 연관이 있으므로 통계학적 또는 수학적 모형을 이용한 발생 예측은 관리 사업의 효율적 수행을 위하여 관심도가 높아지고 있다. 말라리아 수학 모형을 통해 우리는 효과적인 중재 수단을 찾고, 관리 사업을 위한 정책 결정에 유용한 도구로 이용할 수도 있지만, 한편 판단에 대한 오류를 범하게 할 수도 있다. 그 이유는 현실에서는 모형에서 필요한 신뢰할 만한 모수를 추정하는 것이 매우 어렵기 때문이다[2]. 본 연구에서는 최근 외국의 말라리아에 관한 수리적 모델 연구를 기반으로 국내 현실에 맞는 수학 모형 개발과 시뮬레이션으로 예상 발생률을 분석하여 말라리아 관리를 위하여 예측할 수 있는 사용자 중심(end user)의 프로그램을 개발하고자 하였다.

Ⅱ. 몸 말

1. 한국형 말라리아의 수학 모형

삼일열 말라리아가 열대열 말라리아와 다른 특징은 hypno-zoite에 의해 유발되는 장·단기 잠복기이다. 감염 모기가 사람을 물면 혈관을 통해 인체감염형인 포자소체가 간으로 들어가서 분열소체로 발육한 다음 cryptomerozoite가 터져 나오고 이들 대부분은 적혈구로 침입하지만 간세포로 들어가기도 한다. 간세포로 들어간 포자소체가 분열증식 없이 장기간 잠복하면서 그대로 경과하는 수도 있는데, 이때의 충체를 hypnozoite라 하며 재발의 원인이 되기도 한다. 열대열 말라리아 모형으로는 hypno-zoite에 의한 이러한 장·단기 잠복기를 설명할 수 없다. 이번에 개발한 수학 모형은 delay term을 적용하여 장·단기 잠복기를 설명하고자 하였다. 또한 1979년 말라리아 박멸 선언 이후 1993년 현역 군인 1명으로 재출현한 한국의 상황을 설명하지 못한다. 한국에서의 발병에 영향을 미치는 주요 인자로는 무증상(asymptomatic) 환자와 모기의 sporozoite rate가 있다. 무증상 환자는 모기에 의해 감염은 되었지만 증상이 나타나지 않는 환자들로 이들도 모기를 감염 시킬 수 있기 때문에 우리 모형에서는 이러한 집단을 Infectious한 부류에서 고려하였다. 또한 모기의 sporozoite rate이라는 것은 모기의 침샘에서 sporozoite가 발견되는 암컷 모기의 수를 의미하는 것으로 모기의 수명과 감염된 사람을 흡혈하는 횟수와 관계가 있다. 즉, 매일 감염된 사람을 흡혈하는 횟수가 증가 할수록 sporozoite rate는 증가하고 하루 동안의 모기 생존율이 높을수록 증가한다[3].
  연구자들은 2006년부터 2008년까지 질병관리본부로 보고된 말라리아 사례조사서 DB 분석과 문헌고찰을 통하여 수학 모형에 필요한 모수를 추정하였다. 그리고 한국형 말라리아 수학 모형 개발을 위하여 기존의 열대열 말라리아 수학 모형과 삼일열 말라리아 수학 모형을 조사하고 지연형 미분방정식 해법을 검토 및 비교 분석하여 수학 모형을 개발하였다. 이 수학 모형을 공학용 소프트웨어인 Matlab으로 시뮬레이션 하였으며, 한국형 말라리아의 특성을 기반으로 말라리아 예상환자 수를 추정할 목적으로 C++6.0을 기반으로 MFC(Microsoft Foundation Class Library)를 이용하여 사용자 중심의 프로그램(P.Vivaxsim)을 개발하였다.
  이 연구의 수학 모형에서는 사람의 전체 집단을 감수성 집단(susceptible humans), 잠복기 집단(exposed humans), 전염성 집단(infectious humans), 치료 집단(treated humans)으로 나누고, 모기 전체 집단을 3개의 부류(susceptible mos-quitoes, exposed mosquitoes, infectious mosquitoes)로 나누었다(Figure 1). 전염성 집단은 일반적으로는 증상이 있는 집단을 의미하지만 이 모형에서는 모기가 흡혈을 하였을 경우 감염시킬 수 있는 상태의 인구 집단을 의미하는 것으로 유증상자와 무증상 감염자, 회복 후 면역을 가진 사람까지 포함하였다. 치료 집단은 전염성 집단에서 진단을 받고 치료를 시작한 사람을 의미한다. 다른 수학 모형들처럼 회복된 집단(recovered humans)을 고려하지 않고 치료 집단(treated humans)을 고려한 이유는 단순히 병에 걸리고 회복된 집단들 간의 이동만을 보는 것이 아니라 감염된 사람들을 빨리 진단하고 치료하는 것이 말라리아를 예방 관리하는데 얼마나 효과적인지 알아보기 위해서이다.

  감수성 집단의 사람이 감염된 모기에 물렸을 경우 잠복기 집단으로 들어가게 되는데 사람이 모기에 물렸을 경우 바로 모기를 감염시킬 수 있는 것이 아니기 때문이다. 잠복기 집단에 있던 사람은 증상이 나타나거나 증상이 나타나지 않더라도 분열된 원충이 혈액 속에 존재해서 모기를 감염시킬 수 있는 상태가 되면 전염성 집단으로 들어가게 된다. 또한 병이 회복되었더라도 재발하는 경우가 생기는데 이 경우는 원충이 hypnozoite 상태로 간 속에 숨어 있다가 다시 혈액으로 나오게 되는 경우로 전염성 집단에서 잠복기 집단에 속해 있다가 증상이 나타나게 되면 다시 전염성 집단에 포함된다. 전염성 집단에서 몇몇의 사람들은 자연적으로 치유되어 감수성 집단으로 들어가고 몇몇의 사람들은 병을 진단 받아 치료 집단으로 들어가게 된다. 치료 집단에 속한 사람들은 병원에서 치료를 받고 온전히 회복되어 감수성 집단에 들어가거나 hypnozoite가 간 속으로 숨어들어 잠복기 집단으로 들어가게 된다(Figure 2).

  Figure 2에서 사람 집단들 간의 이동량을 보면, 먼저 Λh만큼 유입되는 인구와 전체 인구 Nh(=Sh+Eh+Ih+T)에서 Ψh의 비율로 사람이 태어나면 ΨhNh의 수 만큼 감수성 집단으로 들어가게 되고, λh의 비율로 감염되면 λhSh만큼 감수성 집단에서 잠복기 집단으로 빠져 나가게 된다(a+b+c=1). 전염성 집단을 보게 되면 a,b,c의 비율로 나뉘게 된다. 그중 aIh의 인구에서 β의 비율로 재발하게 되면 aβIh만큼 다시 잠복기 집단으로 들어가게 되고 bIh의 인구에서 α의 비율로 감염된 사람이 병을 진단받게 되면 bαIh만큼 치료 집단으로 들어가게 된다. 여기서 재발률이란 자연 치료된 사람 중 원충이 완전히 사라지지 않고 간 속으로 hypnozoite가 숨어 들어가게 되는 비율을 말한다. 또 cIh의 인구에서 ρh의 비율로 자연 회복된다면 cρhIh만큼 다시 감수성 집단으로 들어가게 된다. 치료 집단은 e,f의 비율로 나뉘는데(e+f=1) 먼저 eT의 인구 중 γ의 비율로 치료가 되면 eγT의 인구는 감수성 집단으로 들어가고 fT의 인구 중 δ의 비율로 재발이 되면 fδT의 인구는 치료 집단으로 들어가게 된다. 여기서 재발률인 δ는 치료 후 재발률이므로 앞서 말한 β와는 다르다. φh는 사망률로 각각의 집단에서 φh의 비율만큼 빠져나간다.

모기 집단 간 개체수의 흐름은 Figure 3에서 보여진다. 유입되는 모기의 수는 생각하지 않고 모기의 산란율인 Ψv는 모기가 자라서 성충이 되는 비율로 ΨvNv의 수 만큼 감수성 집단으로 들어가게 된다. 모기에서는 수직감염이 없고 감염된 사람에게서 흡혈할 경우 λv의 비율로 감염되어 잠복기 집단으로 들어가게 되다. 그 이유는 모기가 흡혈할 때 사람의 몸속에 있던 원충이 모기의 몸속으로 들어가서 생식분열을 한 뒤 다시 사람을 감염 시킬 수 있는 sporozoite가 침샘까지 오는데 시간이 걸리기 때문이다. 잠복기 집단의 모기는 Vv의 비율로 전염성 집단으로 들어가는데 Vv는 잠복기의 역수와 같다. φv는 모기의 사멸률로 각 집단에서 모기가 φv의 비율로 빠져나간다. 모형에서 각 집단 별로 들어오는 양을 (+)로 하고 나가는 양을 (-)로 하여 시간에 따른 변화량을 식으로 만들면 다음과 같은 연립미분방정식을 얻을 수 있다.

  Table 1에서 처음 4개의 식은 사람에 관한 식이고 다음 3개는 모기에 관한 식이다. λh는 모기로부터 삼일열 말라리아 원충이 사람에게 전염되는 비율이며, βhv는 모기가 사람을 감염시킬 확률이고 σvh는 모기의 biting rate이다. 반대로 λv는 사람으로부터 모기에게 원충이 전염되는 비율로 모기가 감염될 확률은 βvh이다. p와 q는 각각 단기 잠복기, 장기 잠복기의 비율이고(p+q=1),  τ1,τ2는 각각 단기 잠복기, 장기 잠복기이다.

2. 모수 추정

 1) 진단율(a)=1/8.216
 말라리아 환자 자료상의 증상 발생일과 진단일의 시간차를 계산한 결과 평균 8.216일이었다.

 2) 잠복기-모기(Sporogenic development(extrinsic rate))(Vv)
 삼일열 말라리아 모기의 경우 잠복기 x와 온도 y와의 관계식을 나타내면 아래와 같다. 수학 모형에서 Vv= 1/x 이므로, 온도 y에 관한 식으로 표현이 가능하다. 따라서 수학 모형을 컴퓨터 시뮬레이션 할 때, 그 지역의 실제 기온을 이용할 수 있다.
                    

 3) 자연회복률(ρh)=0.0016
   감염된 환자가 자연 치료되는데 걸리는 기간은 2년 반에서 3년 정도로 알려져 있다[4].

 4) 사람에서 모기로, 모기에서 사람으로 감염시킬 확률
               

 모기가 사람을 물었을 때, 사람에서 모기로, 모기에서 사람으로 감염시킬 확률은 0.5 정도인 것으로 알려져 있다〔4〕. 이 확률에 대한 우리나라의 실험적 연구는 부족한 실정이다.

 
5) 시간당 모기가 사람을 무는 율
                
  모기가 사람을 무는 율은 온도에 의해 다소 영향을 받으나 일반적으로 0.3인 것으로 알려져 있다[5]. 이에 대해서는 우리나라의 간헐적인 실험 자료가 있으나, 수학 모형을 시뮬레이션 하기에는 유용하지 못하다.

 6) 모기의 산란율 및 사멸률
 이 연구에서는 기온과 강수량 중에서 기온에 따른 모기의 산란율과 사멸률을 구하는 것에 중점을 두었다. 만약 기온의 절대적 수치로써 산란율과 사멸률을 구한다면 5, 6월과 9월은 평균기온이 비슷하기 때문에 두 시기 모두 모기 수가 증가하거나 감소하여야 한다. 하지만 실제적으로 5, 6월은 증가하는 반면 9월은 감소한다. 이런 현상을 설명하기 위해서 기온의 변화량에 따른 모기의 산란율과 사멸률을 구하였다.

 여기서 4를 곱한 것은 철원과 파주의 모기 자료를 이용하여 모기의 대략적인 증가 추이를 살펴본 결과 5월의 모기 수와 7, 8월의 모기 수가 대략 100-1,000배 증가 하는데 4를 곱하면 대략적으로 그 정도의 증가 추이를 보이기 때문이다. 위 식에서 산란율과 사멸률의 정확한 식을 구하려면 산란율과 사멸률 중에서 하나의 값에 대한 정보가 더 필요하다. 성충 모기가 보통 10-30일 정도 살기 때문에 모기의 생존 여건이 점점 좋아지는 시기(5-7월 초)의 사멸률은 0.03-0.05, 장마 시기에는 0.2, 기온이 하강하는 시기(8-10월)에는 0.1로 고정된 값을 주었다. 그리고 산란율은 위 식에서 사멸률을 더하면 다음과 같다.

 7) 잠복기-사람
 기본적인 방법은 우리나라 말라리아 발생 데이터를 이용하여 잠복기를 계산한 논문의 방법을 따르고 몇 가지 부분을 보완했다. 모델 식에 들어가는 단기 잠복기과 장기 잠복기 환자의 비율(식의 p와 q)은 사례 조사를 통해 계산하지 않고 실험 데이터를 이용했다. 북한 strain P. vivax 실험에 따르면, 77명 중에 58명(75.3%)이 장기 잠복기를 가진 것으로 되어 있다[4].

3. 수학 모형을 이용한 말라리아 유행 지역의 시뮬레이션 결과
  각각의 모수들을 이용하여 우리나라 말라리아 위험지역 중 일부 지역을 대상으로 시뮬레이션 하였다. 각 지역별 기상, 모기 발생현황, 잠복기의 환자 수가 다르기 때문에 지역별 모수가 전체 모수를 대변할 수 없어 전국의 시뮬레이션은 의미가 없다. 모든 시뮬레이션은 0인 시점이 5월 1일이며 일 단위이다.

 1) 인천광역시 강화군의 2008년 예상 T값과 환자 수
 2007년 자료를 이용한 2008년 5월 1일부터 다음해 4월 30일까지의 실제 T와 예상 T를 matlab 프로그램을 사용하여 하루 단위로 시뮬레이션 하였다. 그 시뮬레이션에 대한 변수 값은 아래와 같다.

 여기서 실제 T의 값을 전날 환자 수와 당일의 환자 수의 합으로 구하였다. 이를 이용하면 예상한 T를 이용하여 역으로 예상 환자수를 구할 수 있는데 실질적으로 T 그래프의 반값이 되는 연속적인 그래프가 나온다. 하지만 실제 환자 수는 불연속적으로 되어 있기 때문에 실제 T값과 환자 수가 같을 수도 있다. 그래서 예상 T를 이용한 예상 환자 수에 비해 실제 환자 수는 돌출되는 부분이 나올 수 있다. 아래 그래프는 실제로 우리가 궁극적으로 알고 싶은 것이 바로 미래의 환자 수이므로 실제 환자 수와 시뮬레이션 한 예상 환자 수를 그래프로 나타낸 것이다. 아래의 그림은 2007년도 자료를 이용하여 2008년 강화군의 예상 T값과, 실제 T값 일별 실제 환자 수, 예상 환자 수를 나타낸 것이다.

 2) 응용프로그램 P.vivaxsim을 이용한 말라리아 발생 예측 결과
  예상하고자 하는 지역의 전년도의 5월 1일부터 그해 4월 30일
까지의 감수성자 수, 잠복기 환자, 전염성 환자, 치료중 환자, 감수성 모기, 잠복기 모기, 전염성 모기 수에 대한 자료를 데이터 파일(history 파일)로 저장해 놓고 기온 자료도 데이터 파일로 저장하여 사용한다. 아래의 그래프는 P.vivaxsim을 이용한 2009년 지역별(강화군, 파주시, 김포시, 철원군) 말라리아 발생 예측 결과이다. 지역별 환자수는 강화군 182명, 파주시 48명, 김포시 89명, 철원군 67명, 월별 환자수가 가장 많은 달은 7월로 예측하고 있다.

Ⅲ. 맺는 말

  기생충 질환의 모델링은 1934년 Kostitzin에 의해 시작된 이후 MacDonald의 연구에 의해 활성화되었다. 말라리아 연구에 획기적인 도약을 이룩한 사람은 로널드 로스(Ronald Ross)로 말라리아의 수학적 모델은 1957년 MacDonald의 책에 소개되어진 Ross의 모델과 그 모델의 확장으로 시작되었다. 기존의 말라리아 모형은 대부분 Ross-MacDonald의 모형을 기반으로 한 열대열 말라리아에 관한 것이어서 삼일열 말라리아의 주 특징인 장·단기 잠복기를 잘 표현하지 못하고 있다. 그 외 2-3개의 삼일열 말라리아 전파 수학 모형이 알려져 있지만 아주 초기단계로 수학적 분석 단계에 지나지 않고 있다. 이렇듯 말라리아 모델링의 역사는 길지만, 불규칙해서 아직 표준이 되는 모형은 만들어지지 않았다. 기존의 삼일열 말라리아의 수학 모형로는 P. Pongsumpun and I. Tang의 모델[6], H. Isikawa 모델[7], Zoysa 모델[8]등이 보고된 바 있다.
  삼일열 말라리아 전파에 가장 적합한 수학 모형은 장·단 잠복기와 재발을 잘 설명할 수 있는 지연형 미분방정식이다. 지연형 미분방정식은 현재 과학의 모든 영역에서 중요한 위치를 차지하고 있지만 생물학 영역의 질병 역학분야에서는 특히 더 중요한 위치를 차지하고 있다. 그런 지연형 미분방정식은 일반적인 상미분 방정식보다 좀 더 현실적인 예측을 위해서 사용한다. 일반적인 지연형 미분방정식은 다른 상미분방정식보다 더 많은 정보가 필요하다. 왜냐하면 지연형 미분방정식 안에서 미분계수는 이전시간의 해에 의존하므로 이전 시간(t-τj)의 해의 값을 나타내는 history함수가 있어야 하기 때문이다. 대부분의 지연형 미분방정식에서는 이 history함수가 고정된 값으로 되어 있지만 고정되지 않은 history함수도 종종 볼 수 있다. 여러 가지 프로그램들은 지연형 미분방정식과 상미분방정식 시스템을 통합하는 것으로 정확히 Runge-Kutta method를 사용을 한다. 이 연구에서 사용하는 C++프로그램 또한 Runge-Kutta method를 기본적으로 사용하고 있다.
  이 연구에서는 수학 모형에서 환자 발생에 가장 중요한 영향을 미치는 모기와 관련된 데이터 부족으로 말라리아 유행 5개 지역만을 선정하여 시뮬레이션 하였다. 그리고 이 모형의 모수 중 자연회복률, 사람에서 모기로, 모기에서 사람으로 감염시킬 확률, 시간당 모기가 사람을 무는 횟수는 국내에 축적된 데이터가 없기 때문에 외국의 자료를 사용하였다. 이러한 현실에서 얻어내기 어려운 모수들이 국내 말라리아 발생 예측에 한계점으로 작용할 수밖에 없다. 이번 시뮬레이션 연구 결과의 활용에 대한 주의점은 말라리아 수학 모형에 의해서 결과로 얻어지는 환자 수는 실제 발생 환자수와는 다소 차이가 날 수 있다는 것이다. 그 이유는 말라리아 매개모기 관련 자료의 부족으로 좀 더 엄밀한 식의 적용에 어려움이 있었을 뿐 아니라, 수학 모형은 말라리아 환자 발생의 당해연도의 추이를 알아보고자 하는 것이기 때문이다. 즉, 수학 모형은 특정 날짜에 환자 수가 얼마인가보다는 그 해에 언제 환자 수가 가장 많이 발생할 것이며, 다음 해에는 그 시기가 앞당겨질 것인지 늦어질 것인지와 그 해에 비해 다음 해에 환자 수가 줄어들 것인지 늘어날 것인지 등의 전체 경향을 보여주는 것이다. 2009년에도 2008년도와 유사한 기온 분포를 보인다면, 이 모형에 의하면 말라리아 환자 수는 2008년과 비슷한 수준이거나 더 적게 나올 것이라 예상된다. 3월 중 기상청에서 2009년 여름 예상 기온 자료가 발표되면, 2009년 온도 자료를 P.Vivaxsim에 사용하여 2009년도 지역별 환자 추이를 좀 더 정확히 예측해 볼 수 있다.
  개발된 P.Vivaxsim의 정확한 예측 시뮬레이션과 효과적인 활용을 위해서는 모수의 특성별 심층 연구, 예를 들어 예방약을 지속적으로 군인에게 투여하였을 때 모기에서 사람으로의 감염률의 감소 정도 연구, 기피제를 사용한 경우 모기가 사람을 무는 확률(biting rate)의 변화 연구 등이 필요하다. 말라리아 수학 모형에서 가장 큰 영향을 미치는 것은 모기임에도 불구하고 현재 수집되어 이 연구에 사용한 모기에 대한 자료는 주별 자료이고 채집지역도 전국을 대표할 수 없으므로 일별 시뮬레이션에 적용하여 예측하는 데에는 한계가 있다. 따라서 적어도 6월부터 9월 초순까지 일별로 동일지역에서 모기 밀도뿐만 아니라 및 감염 모기 수의 조사가 꾸준히 이루어져야 할 것이다.

Ⅳ. 참고문헌
 1. Richard J. Maude, The role of mathematical modelling in malaria elimination and eradicatioin, Trans R soc Trop Med Hyg, 2009
 2. David G. Regan, Modelling sexually trandmitted infections, Trans R soc Trop Med Hyg, 102 : 207-208, 2008
 3. George MacDonald, C.M.G., M.D., F.R.C.P., The epidemiology and Control of Malaria, London Oxford University Press, 1957
 4. NA Tiburskaja, OS Vrublevskaja : The course of infection caused by the north Korean strain of P. vivax, 1977
 5. Shigui Ruan, Dongmei Xiao, John C. Beier : On the Delayed Ross–Macdonald Model for Malaria Transmission, Mathematical Biology, vol. 70, No.4, 2008
 6. Puntani Pongsumpun and I-Ming Tang: Mathematical Model for the transmission of Plasmodium Vivax Malaria, International Journal of Mathematical
    Models and Methods in Applied Sciences, vol 1, 2007
 7. Hirofumi Ishikawa, Akira Ishii, Nobuhiko Nagai, Hiroshi Ohmae, Masakazu Harada, Setsuo Suguri, Judson Leagasia, A mathematical model for the
     transmission of Plasmodium vivax malaria, Parasitology International 52, 81~93, 2003
 8. De Zoysa AP, Mendis C, Gamage-Mendis AC, Weerasinghe S, Herath PR, Mendis KN : A mathematical model for Plasmodium vivax malaria
     transmission : estimation of the impact of transmission-blocking immunity in an endemic area, Bull World Health Organ, 69(6) : 725~734, 1991