시뮬레이션 함수 추정과 그 예
송 한 식*
|
I. 서론 V. 시뮬레이션 응용 2 II. GMDH 방법 VI. 결론 III. GMDH를 이용한 시뮬레이션 함수추정 VII. 참고문헌 IV. 시뮬레이션 응용 1 |
I. 서론
본 연구는 시뮬레이션 과정을 하나의 함수로 추정하여 표현하고 그것을 활용하여 시스템의 모수(통제변수)를 정하는 것에 관한 내용이다. 주어지는 조건에 따라서 투입을 각기 다른 산출 결과로 변환하는 시뮬레이션 과정은 하나의 함수로 볼 수 있다. 문제는 그 함수를 분석적인 함수 형태로 구체적으로 표현하는 것이다. 본 연구에서는 주어지는 조건과 투입에 따른 산출 결과를 GMDH (group method of data handling) 방법에 적용하여 명시적인 함수(explicit form)로 표현하는 것을 다루었다.
일반적으로 시스템의 작동은 "투입(I) - 변환처리(T) - 산출(O)"의 과정을 따른다. 이러한 변환 처리 과정을 설계하고 그 시스템을 운영하는 데 있어서 정책적인 요소를 시스템 모수(systems parameters)라고 부른다. 예를 들면, 대기행렬 시스템에서 서버의 서비스율이나 서버의 수는 시스템 모수이다. 시스템의 설계와 운영을 결정하는 것은 바로 그 시스템의 모수를 정하는 것이 된다. 이러한 뜻에서 시스템 모수를 우리는 통제변수(controllable variable)라고 부른다.
통제변수에 대한 적정한 값을 정해주는 것은 시스템의 운영을 책임진 사람(의사결정자)에게는 중요한 과제가 된다. 의사결정자는 시스템을 수리모형으로 표현하거나 시뮬레이션 모형으로 표현한 다음, 시스템의 운영 효율을 최적화할 수 있도록 모수를 찾게된다.
_____________________________________________________________
* 동아대학교 경영학부 교수
따라서 시스템을 모형화할 때에는 이러한 통제변수를 포함하여야 하고, 통제변수가 시스템의 작동에 미치는 영향을 효과적으로 표현하고 있어야 한다. 시스템을 시뮬레이션 모형으로 표현한 경우에 그것에서 최적의 통제변수를 찾는 것을 우리는 시뮬레이션 최적화(simulation optimization)이라고 부른다. 시뮬레이션 최적화는 다음 3 가지 방법으로 접근할 수 있다.
(1) 발견적 방법(heuristic method)
가능한 통제변수 값을 임의로 정한 다음, 난수 흐름에 의한 투입에 대해서 시스템의 산출을 시뮬레이션을 통해서 얻는다. 이러한 과정을 반복하면서 그중 가장 나은 산출결과를 주는 통제변수 값을 최선의 것으로 정한다. 그러나 이때에는 그것이 최선의 것인지 확신할 수 없다는 문제와 통제변수의 여러 가지 조합마다 시스템 시뮬레이션을 계속 반복해 보아야 한다는 불편함이 있다. 만일 시스템 시뮬레이션에 시간이 많이 걸린다면 이것은 더욱 큰 문제가 된다. 발견적 방법보다는 더 체계적인 탐색을 하는 방법으로서 심플렉스 탐색법(simplex search method)이 연구되고 있다(D. A.Humphrey & J. Wilson(2000); R. R. Barton & J. S. Ivey, Jr.(1996)).
(2) 뉴런 넷에 의한 시뮬레이션 흉내내기
다른 통제변수마다 시뮬레이션을 반복해야 하는 경우 시간을 절약하기 위해서는 시뮬레이션을 뉴런 넷으로 제2차 모형화를 하는 방법이 있다. 뉴런 넷으로 모형화할 때에는 시스템의 투입과 통제변수 값을 투입노드(input nodes)로 표현하고서, 시뮬레이션의 산출을 뉴런 넷의 산출노드(output node)로 한다. 뉴런 넷으로 모형화하기 위해서는 시뮬레이션 결과를 많이 얻어서 그것을 샘플 포인트로 하여 뉴런 넷을 훈련시킨다. 뉴런 넷으로 모형화 한 다음에는 가능한 통제변수 값을 뉴런 넷에 투입하면 금방 시뮬레이션의 산출을 비교해 볼 수 있다. 이때에도 발견적인 방법말고는 적절한 통제변수 값을 찾을 수 있는 방법이 없다. 그럼에도 최근에 뉴런 넷의 활용과 응용이 다양한 분야에 시도되면서 시뮬레이션 흉내내기를 이용한 연구가 나타나고 있다(Gupta, et al., 2000).
(3) 유전자 알고리즘을 이용한 시뮬레이션 최적화
최근에는 시뮬레이션 최적화에 유전자 알고리즘(genetic algorithm)을 이용하고 있다. RISKOpimizerⓡ(미국 Palisade사 제품)는 시스템의 운영이 확률적인 과정을 거치는 경우, 시스템의 산출결과에 대한 여러 가지 통계치(평균, 분산, 최소값, 최대값 등)를 목적치로 삼아서 그것을 최적화하는 통제 변수 값을 유전자 알고리즘으로 찾도록 설계되어 있다. 이 방법은 최적해에 가까운 해를 아주 잘 찾아주는 것으로 받아들여지고 있다. 실제로 많은 시험을 해보면 상당히 효과가 있다. 그러나 두 가지 단점이 있다. 최적의 통제변수 값을 찾는 것 이외 추가적인 분석(예: 민감도 분석)을 해 볼 수 없다. 또, 통제변수에 대응한 시스템의 반응(response)을 전체적으로 이해하는데 큰 도움을 주지 못한다. 두 번째 단점은 유전자 알고리즘으로 최적해를 찾는데는 시간이 많이 걸린다는 것이다. 통제변수가 6개 또는 10개 정도만 되어도, 그리고 통제변수의 적절한 범위를 모르는 경우, 최적해를 찾는데는 아주 많은 시간이 걸린다. 또, 그것을 찾아도 과연 그것이 최적해에 얼마나 가까운 것인지 확신하지 못한다는 단점이 있다.
이상의 문제점을 해결하는 다른 한가지 방법은 시뮬레이션을 하나의 함수로 보고서 그 함수를 추정해서 얻는 방법이다. 시뮬레이션을 "투입 -> 함수변환 -> 결과산출"로 보고 (투입 -> 결과)의 대응을 하나의 함수 대응관계(functional correspondence 또는, sample point)로 보는 것이다. 이때, 투입은 난수흐름을 이용한 시스템의 투입과 시스템 운영에 관한 정책적인 모수를 포함한다. 시뮬레이션을 통해서 이러한 함수 대응관계를 나타내는 데이터 포인트를 많이 얻은 다음 그 함수관계를 가장 적절하게 표현하는 함수를 추정하는 것이다. 만일 함수 추정이 가능하다면 그 함수에다 기존의 최적화 기법(주로 비선형기법)을 적용하는 것이 가능하다. 이 방법은 아직 연구되고 있지 않다. 그 이유는 함수를 추정하는 적절한 방법이 없었기 때문이다. 본 연구에서는 GMDH 방법으로 함수추정을 살펴보고 그 예를 살펴보았다.
II. GMDH 방법
GMDH 방법은 1966년 러시아 수학자 A. G. Ivakhnenko에 의해서 처음 소개되었다(A. G. Ivakhnenko, 1966, 1971. 태풍의 진로를 예측하거나, 어군의 이동을 예측하거나, 경제모형에서 이자율이나 주가를 예측하는데 활용되었다(S. J. Farlow, 1984). 그러나 그 이후에는 별로 주목을 받지 못하고 있다. 최근에는 주어진 데이터에서 유용한 정보를 찾아내는 Knowledge Mining에 활용하고 있다. (상용으로 KnowledgeMinerⓡ라는 소프트웨어가 선보이고 있다. (참조, Julian Miller, http://www.knowledgeminer.net)
GMDH 방법은 주어진 독립변수들로서 종속변수를 가장 잘 표현하는 함수를 찾는 방법이다. 독립변수를 '현세대'의 각각의 개체(하나의 생명개체)로 간주한다. 개체 가운데 임의의 두 개의 개체(독립변수)를 선택하여 그 독립변수의 값으로 종속변수를 표현하는 식을 만든다. 그 식으로 표현된 것을 '새로운 세대'의 개체로 삼는다. 이것은 생명체의 암수가 결합하여 후손을 낳는 것에 비유할 수 있다. 만들어진 개체(후손 식)에다 독립변수 값을 투입하여 얻은 '개체의 값'을 새로운 데이터 포인트로 삼는다. 이것을 독립변수 가운데 가능한 두 개의 조합마다 반복한다. 이러한 과정은 세대교체를 이루는 것과 같다. 다음에는 생성된 개체들 가운데서 임의의 두 개의 개체를 골라서 종속변수를 표현하는 새로운 식을 만든다. 이러한 세대교체 과정을 반복해 나가면서 후손 세대에서 최적의 개체를 찾는 것이다. GMDH 방법은 세대교체 과정에서 부적합한 개체는 도태시키고 적자만 살려서 후손을 생성시켜나가는 적자생존의 유전적 과정을 따른다.
GMDH의 모형화 방법은 다음 4가지 단계로 나누어 설명할 수 있다. 이중 두 번째와 세 번째 단계를 반복하여 세대를 교체하면서 중단시점이 될 때까지 계속한다(S. J. Farlow(1984); H. R. Madala & A. G. Ivakhnenko(1994) 참조)
Srep 1) 준비단계
m개의 독립변수 (x1,x2, ... xm)와 종속변수 y가 있는 것으로 가정한다. m개 독립변수 데이터 값과 종속변수 값의 짝으로 된 데이터 포인트 (xi1, xi2, ... xim, yi)를 N개 다음 그림 1과 같이 준비한다.
N개 데이터 포인트를 훈련세트(training set: Nt개)와 체크세트(checking set: Nc개)로 구분한다(N = Nt + Nc). X = (x1,x2, ... xm)을 '현 세대'로 삼고, ( xi )는 그 세대의 한 '개체'가 된다. xi의 데이터 값은 '개체값'이 된다. 준비단계에서는 주어지는 독립변수를 현세대의 개체로 삼는다. 이후 반복에서는 새로 생성된 세대(G+1세대)를 현세대(G세대)로 교체하여 Step 2)와 Step 3)를 반복한다.
Step 2) 세대생성단계
G번째 세대의 개체들 XG = (x1,x2, ... )로서 다음 G+1번째 세대를 생성한다. G세대의 훈련세트의 데이터 포인트(Nt)를 이용하여 임의의 두 개 변수 (xg, xh)에 대해서 y에 대한 다음의 추정식 fk(xg, xh)을 정한다.
fk(xg, xh) = a0 + a1xg + a2xh + a3xg2 + a4xh2 + a5xg*xh
이것은 비선형회귀식 추정에 해당한다. 이 추정식에 (xg, xh)의 N개 데이터 값을 각각 투입하여 얻은 값을 새로운 개체의 값으로 삼고, 그 개체를 zk = (z1k, z2k. ... zNk)으로 나타낸다. 새로운 개체들로 된 Z = (z1, z2, ... zk, ...)을 G+1번째의 '새로운 세대'로 삼는다. G세대에서 개체가 m개 인 경우, G+1세대에서는 개체 zk가 최대 mC2 = m(m-1)/2 개수만큼 생긴다.
Step 3) 평가 및 도태단계
G+1세대에서 체크세트의 데이터 포인터(Nc)를 이용하여 zk에 대해서 다음 평가값을 구한다.
이렇게 구한 평가값이 사전에 정한 기준 R 보다 큰 경우에는 zk을 도태시킨다.
이렇게 구한 평가값 중 최소인 을 선택하고 그것을 그 세대의 '최소평가값'으로 삼고, zb 를 그 세대의 '최적개체'로 부른다. 가 직전 세대의 '최소평가값' 보다 크면, 직전 세대의 최적개체를 시스템 전체의 '최적개체'으로 선정하고 Step 4)의 '중단단계'로 간다.
최소 평가값이 그전 세대의 최소 평가값 보다 작거나 같으면 생성된 세대의 개체 ZG+1 = (z1, z2, ... zk, .. )를 XG로 대체하여 Step 2)의 '세대생성단계'로 간다.
Step 4) 중단단계
최적개체로 선정된 개체를 생성하는 데 투입된 이전 세대의 개체를 역으로 추적하여 최적개체 추정식을 만든다. 이것을 'GMDH 모형식'으로 삼는다. GMDH모형식은 원래의 독립변수에 대한 다항식(polynomial)으로 표현된다.
이 과정에서, fk(xg, xh)의 계수는 다음과 같이 Gauss normal equation으로 구할 수 있다. 우선 데이터 포인터 (xi1, xi2, ... xim, yi) 에 대해서 다음의 식을 적는다.
y1 = a0 + a1x1g + a2x1h + a3x1g2 + a4x1h2 + a5x1g*x1h
：
yNt = a0 + a1xNt,g + a2xNt,h + a3xNt,g2 + a4xNt,h2 + a5xNt,g*xNt,h
이것을 간단히 행렬로 나타내면
Y = X A
가 된다. Y는 Nt x 1 , X 는 Nt x 6, A 는 6 x 1 행렬이다. 행렬 X의 첫째 열은 (1, 1, 1 ... 1)T 으로 된 열벡터이다. 양변에 전치행렬 XT를 곱하면
XT Y = (XTX)A
이다. (XTX)는 6 x 6 행렬이고, X의 각 열이 1차 독립이라면 역행렬이 존재하므로, 그 역행렬을 양변에 곱하면,
A = (XTX)-1XTY
으로서 A = (a0, a1, a2, a3, a4, a5)을 구할 수 있다.
한편 '평가 및 도태단계'에서 각 세대의 '최소평가값'은 세대를 거듭할수록 감소한다. 그러다가 어느 순간 증가하게된다. 이는 그 세대에서 개체의 식이 '훈련세트'의 데이터 포인트에 지나치게 근접한 나머지 '체크세트'에 대해서는 추정값의 오차 가 커졌음을 뜻한다. 이것은 과다모형화(overfitting, overspecification)가 이루어짐을 의미한다. GMDH 방법은 데이터 포인터를 훈련세트와 체크세트로 나눈 때문에 이러한 과다모형화를 방지할 수 있다(H. R. Madala & A. G. Ivakhnenko, p. 11). 경험적으로 훈련세트와 체크세트는 5:5 또는 6:4 정도의 비율로 나눈다.
본 논문에서는 GMDH 방법에서 한가지 변형을 추가하였다. 원래 GMDH는 복잡한 시스템(complex system)을 대상으로 하며 독립변수가 4 - 5개 이상인 경우를 대상으로 한다. 독립변수가 두 개인 경우는 제2세대에서 더 이상 세대교체를 할 수 없기 때문이다. ZG -> Z(G+1)세대로 교체할 때, 원래의 GMDH방법은 ZG세대 내의 결합 (ZG x ZG)만 이용하였다. 그래서 본 논문에서는 변형을 추가하여서, ZG -> Z(G+1)세대로 교체할 때, ZG 세대 뿐만 아니라, 원래 주어진 X1의 개체를 재차 사용하였다. 따라서 ( X1, ZG) -> Z(G+1)으로 세대교체를 하였다. 두 번째 이후 세대에서는, (X1 x ZG )와 ( ZG x ZG )의 결합을 고려하였다( (X1 x ZG-1 )의 결합은 고려하지 않음).
본 논문에서는 GMDH 방법을 엑셀에서 사용할 수 있도록 Visual Basic 프로그램으로 만들어서 사용하였으며, 앞에서 설명한 변형을 프로그램에 추가하였다.
III. GMDH를 이용한 시뮬레이션 함수 추정
투입(I)을 받아들여서 산출(O)로 변환하는 시스템의 변환처리과정은 Fx(I) = O의 함수로 볼 수 있다. 여기서 x는 시스템 모수(통제변수)이고, 시스템 운영이나 설계상의 정책을 반영하게 된다. 이 시스템을 시뮬레이션으로 분석하기 위해서는 시스템을 시뮬레이션 모형으로 표현한 후 시뮬레이션 반복에 의한 결과를 이용한다. 즉, 난수(random number)를 발생시켜서 여러 가지 가능한 투입 I를 생성시키고, 그것에 의하여 얻어지는 산출 O에 대한 통계치 를 시뮬레이션으로 얻게 된다. 보통 는 O의 기대값 E[O]을 시뮬레이션 결과로 산출하지만, 때로는 최대값이나 최소값, 또는 분산이나 표준편차를 시뮬레이션 산출로 얻어서 이용한다.
이 때, 시뮬레이션 산출로 얻어지는 통계치 는 통제변수 X에 대한 함수가 된다. 통제변수 X를 독립변수로 간주한다면, 시뮬레이션은 독립변수 X에 대한 종속변수 의 값을 얻어내는 과정으로 볼 수 있다. 종속변수 를 으로 나타내어서 시뮬레이션 함수(simulation function)로 표기하자. 그리고, i번째 시뮬레이션을 통하여 관찰된 값 를 으로 표기하자.
독립변수 X가 (x1,x2, ... xm)로 여러 개인 경우를 생각해 보자. i번 째 시뮬레이션에서는 통제변수 값을 ( xi1, xi2, ... xim )으로 고정하고, 난수발생을 통하여 여러 가지 가능한 투입 ( I )를 만들어서 반복투입(iteration)하고, 그 결과 얻어지는 통계치 ( = )를 얻게 된다. 통제변수 값을 여러 가지로 바꾸면서 이러한 시뮬레이션을 N회 반복하면, N개 데이터 포인트 ( xi1, xi2, ... xim, ), i = 1, 2, ..., N 를 얻을 수 있다. 그러면 앞에서 설명한 GMDH 방법을 이용하여 함수: X -> 을 추정할 수 있게된다. 단, 이러한 시뮬레이션 함수를 추정할 때에는 시뮬레이션 함수 가 X의 적절한 구간 안에서 연속(continuous)이라는 것을 전제로 한다. 그리고, GMDH로 추정된 다항식 함수(polynomial function)이 미분가능하기 때문에, 원래 함수 = 가 미분가능하지 않더라도 미분가능한 함수로 추정한다는 것을 전제로 하고 있다.
이상의 GMDH에 의한 시뮬레이션 함수 추정을 실제 시뮬레이션에 적용해 보자.
IV. 시뮬레이션 응용 1
변수가 2개인 경우의 "가상항공회사"의 over booking 문제에 적용해 보았다. 이 예제는 Cliff T. Ragsdale(1998)에 있는 예제를 RISKOptimizer guide book에서 발췌한 것을 인용한 것이다.
예제: 이 항공사는 좌석이 19개인데 정규 요금은 $195이고 할인요금은 $85이다. 단, 예약 후에 나타나지 않을(불현) 경우, 정규요금으로 산 항공권은 100% 환불해 주지만, 할인 요금으로 구입한 항공권은 $50을 공제하고 환불해 준다. 만일 예약 후, 좌석이 없을 경우(최과 부킹)은 확률 (0,1, 0.3, 0.4, 0.2)으로서 각각 $150, $200, $250, $300 씩 보상하게 된다.
이 항공회사의 운영은 다음과 같은 확률적인 과정을 따른다.
. 정규요금 예약자의 불현율은 평균이 0.20이고 표준편차가 0.03인 정규분포를 따른다.
. 정규요금 탑승 수요자는 (3,7,5)인 삼각형 분포(triangular)를 따른다.
. 할인 요금예약자 중 불현율은 평균이 0.10이고 표준편차가 0.01인 정규분포를 따른다.
. 할인요금 탑승 수요자는 (12, 20, 35, 10, 90)인 삼각형 분포(triangular)를 따른다. 하위 10%는 12명, 가장 가능성이 높게는 20명, 90%까지는 35명이 된다.
. 예약 판매 예정인 좌석수는 19 -35사이에 조정한다(통제변수).
. 정규요금 판매 비율은 0% ~100%사이에 조정한다(통제변수).
. 판매 순이익의 분포에 대해서 표준편차는 400이하이어야 한다.
이러한 조건하에서 변상비용을 제외한 순이익을 최대화하려면, 예약판매 좌석수와 정규요금 판매비중을 어떻게 정하면 되겠는가하는 것이 문제이다. 그림 2는 이것을 RISKOptimizer 패키지로 시뮬레이션하여 최적해를 찾는 엑셀 모형이다.
이 시뮬레이션 모형에서 통제변수는 예약좌석수(x1)와 정규요금 판매비율(x2)이고, 시뮬레이션 산출은 순이익의 기대값(y)이다. 이것을 이용하여 다음의 절차로 시뮬레이션 함수식을 추정하였다.
① 한 쌍의 통제변수 값에 대하여 500회 난수발생을 하여 예약-탑승-변상의 과정을 거쳐서 순이익을 500번 발생시키고 그것의 평균값을 순이익 기대값으로 구하였다.
② RIDKOptimizer는 이러한 시뮬레이션을 150회 반복하면서 유전자알고리즘으로 기대이익을 최대로 해주는 예약좌석수(x1)과 판매비율(x2)을 스스로 찾는다. 이과정에서 모두 122개의 데이터 포인터를 관찰하였다. (나머지는 순이익 분포의 표준편차가 400을 넘어서 무시하였다.) 시뮬레이션 결과는 예약판매좌석을 33석으로 하고, 정규요금 비율을 58.29% 일 때 얻어진 순이익 평균이 $2,325이 최대로 생각되는 것이었다.
③ 이러한 과정에서 얻어진 데이터 포인터(이것을 로그파일이라고 함)를 GMDH 방법에 투입하여서 시뮬레이션 함수식을 추정하였다. 그림 3의 상단에 있는 B,C,D열의 4 ~ 125번 행은 122개의 데이터포인터 (예약좌석수, 정규요금판매비율, 순이익 기대값)를 나타낸다. 122개 데이터 포인터 중, 70개는 훈련세트, 52개는 체크세트로 구분하여 활용하였다.
④ 이것을 GMDH 프로그램에 투입하여 순이익(y)을 예약좌석수(x1)와 정규요금 판매비율(x2)에 대한 함수로 표현해 보았다. 그림 3의 하단에는 GMDH에 의한 추정식의 계수를 세대별로 순차적으로 적은 것이다. 아래는 그 식이다. 이 식에서는 독립변수 x1, x2에서 유도된 식(개체) x3, x4, x5, x6, x7, x21을 차례로 투입하여 최적개체 x61의 식을 유도할 수 있다. (여기서는 그 식의 표현을 생략하였다.)
x3 =1041.339 +35.69506x1+ 2680.233x2 -1.03653x1^2 -3723.8x2^2 +40.12604x1*x2
x4 =2358.084 -22.3416x1 -1.40111x3 -0.63484x1^2 +0.000438x3^2+0.031879x1*x3
x5 =2437.338 -1278.67x2 -1.91431x3 +1062.977x2^2+ 0.00089x3^2 +0.100827x2*x3
x6 =672.0297 -6.1961x1 +0.325359x4 +0.209038x1^2 + 0.000214x4^2-0.00338x1*x4
x7 =765.3232 -35.9736x1+ 0.683575x5 -0.09464x1^2 -9.9*10^(-05)*x5^2
+0.023821x1*x5
x21 =7.061734 +1.091683x6 -0.08509x7 -0.00681x6^2 -0.00602x7^2 +0.012828x6*x7
x61 = 48.5697 -46.9216x2 +0.959096x21 +32.1401x2^2+1.32*10^(-05)*x21^2
-0.0059x2*x21
그림 3의 E열은 x61의 식에 해당 x1, x2 값을 투입하여 얻은 순이익 추정값( )을 적은 것이다. 이 값은 시뮬레이션 산출값인 D열의 값과 상당히 유사하다.
추정된 식을 그래프로 그려 보았다. 이그래프는 이 시뮬레이션의 결과를 투입값에 따라 그린 시뮬레이션 반응 표면(simulation response surface)이다.
그림 4는 Mathematicaⓡ 에 위의 식을 투입하여 변수 x1, x2에 대한 순이익 x61의 분포모양을 그린 것이다. 그림을 보면, 순 이익은 예약 판매 좌석 수 30 ~ 35 사이와 정규요금 비중을 60% 부근으로 할 때 순이익 평균이 대체로 높음을 알 수 있다. 실제로 예약 좌석을 35석으로 하고, 정규요금 판매비중을 60%으로 하고 500회 반복 투입하여 시뮬레이션을 해 보아도 같은 순이익 평균값 $2325을 얻을 수 있었다. 이 그림을 통해서 항공사의 over booking 문제를 더욱 잘 이해할 수 있다. (Airyield_논문용.xls, Airyield_논문용_로그파일.xls, Airyield_gmdh함수추정.xls)
V. 시뮬레이션 최적화 응용 2
- 변수가 여러 개인 경우
GMDH에 의한 함수 추정을 통제변수가 여러 개인 경우에 적용해 보았다.
예제문제: 신년도 칼렌더를 4종류 구입해서 판매하여 이익을 남기려고 한다. 칼렌더 1,2,3,4의 구입 비용은 각각 $3.00, $4.00, $5.00, $6.00 이고, 판매가격은 $6.00, $7.00, $7.50, $8.50 이다. 그러나 수요는 불확실하며, 각각, 차례로,(1800, 200, 2500), (1300, 1700, 2000), (900, 1600, 1800), (1400, 1600, 1900) 개인 삼각형 분포(triangular distribution)를 따른다. 팔다 남은 칼렌더는 각각 $1.00, $1.20, $1.40, $1.60씩 되팔 수 있다. 단 칼렌더 구입예산은 $2500을 넘을 수 없다. 각각 몇 개씩 구입하면 이익이 최대가 되겠는가?
수요가 확률적인 분포를 띄므로 선형계획으로 풀지 못한다. RISKOptimizer를 이용하면 유전자 알고리즘으로 최적해를 찾는데, 찾은 최적해는 각각1987개, 1435개, 939개, 1434개일 때 기대되는 최대이익은 $16,113.66 이었다. (이상은 Wayne Winston의 "Decision Making under Uncertainty with RISKOptimizer" 의 pp.23-26에서) 그림 5는 이것을 엑셀으로 표현한 시뮬레이션 모형이다.
(multinews-논문용.xls)
이 문제에 대해서 다음과 같이 GMDH로 시뮬레이션 함수를 추정하였다.
① RISKOptimizer로 시뮬레이션을 약 400회 수행하였다. 1회 시뮬레이션은 이익값의 평균이 수렴할 때까지 반복(iteration)을 계속하였다. 그 중 340개를 샘플 데이터 포인터로 뽑았다. 이 샘플에는 4개 칼렌더의 구입개수(x1, x2 x3, x4)와, 그때의 기대되는 최대 이익 평균값(y)이 포함되었다. 여기에는 전체 예산 조건을 만족하지 못하는 것(invalid)도 포함하였다. 이 데이터 포인터에는 RIDKOptimizer가 최적해 부근에서 찾은 데이터 포인터(1987, 1435, 939, 1434, 이익평균값)들이 많이 포함된다.
② 샘플된 340개 데이터 포인터를 가지고 200개를 훈련세트로, 140개를 체크세트로 하여 GMDH로 함수를 추정하였다(그림 6).
③ 추정된 함수를 가지고, 엑셀의 해찾기(Solver)로서 최적해를 찾았다. GMDH 추정식의 계수를 차례로 대입한 추정식은 그림7 에서 차례로 읽을 수 있다. 찾은 최적해는 2035, 1645, 1060, 1168개이다. 그 때 해찾기로 구한 최대이익평균은 $16,388 이었다(그림 7).
④ ③에서 찾아진 해(칼렌더 구입수량)를 가지고 다시 시뮬레이션 모형에 투입하여 기대 이익의 평균치를 을 구한 결과 평균은 $16295 ~ $16326 사이에 분포하였다(그림5 참조). 시뮬레이션 결과는 GMDH 함수 추정에 의한 예측치와 일치하지 않다. 따라서 GMDH 함수 추정이 시뮬레이션과 정확히 일치하지는 않음을 알 수 있다. 그러나, 새로 투입한 해(2035개, 1645개, 1060개, 1168개)는 RISKOptimizer로 구한 해보다 더 개선된 해임을 알 수 있다. 기대 이익의 평균치가 $16,113 보다 높기 때문이다.
(참고: multinews_논문용.xls, multinews_논문용_gmdh.xls, multinews_논문용_최적해.xls)
VI. 결론
본 연구에서는 시뮬레이션으로 얻어진 데이터 포인터를 가지고 그 시뮬레이션을 함수로 추정하는 것을 시도하였다. 이것을 구체적인 시뮬레이션에 적용하여 살펴보았다. 시뮬레이션의 통제변수가 2개인 경우의 예는 GMDH로 추정된 함수를 그래프로 그려서 시뮬레이션 반응을 시각적으로 볼 수 있었다. 그 결과, 두 개 통제변수에 대해서 시뮬레이션 결과치(순이익 평균값)가 최적해 부근에서 민감도가 비교적 낮았으며, 비교적 최적해를 포함하는 주변의 많은 영역에서 결과치를 비슷하게 얻을 수 있음을 알 수 있었다. 이는 시뮬레이션 결과를 가지고 시뮬레이션 함수로 추정한 것이 시뮬레이션 반응표면(simulation response surface)을 시각적으로 이해하는 데 활용할 수 있는 예이다.
통제변수가 여럿인 경우에도 시뮬레이션 결과로 얻어진 데이터포인터를 GMDH에 투입하여 시뮬레이션 함수를 추정하였다. 얻어진 함수를 가지고 다시 최적해를 비선형계획법으로 푼 결과(이것은 엑셀의 해찾기(Solver)기능으로 하였음)를 재차 시뮬레이션에 투입하였더니, 그 앞에 시뮬레이션으로 얻은 해보다 더 나은 결과를 주었다. 이는 GMDH에 의한 시뮬레이션 함수 추정을 시뮬레이션 최적화에 재차 투입하여 현재 사용하고 있는 유전자 알고리즘이 발견하지 못하는 최적해를 재차 발견할 수 있음을 보여주는 예이다.
이상, 두 가지의 경우의 예에서 알 수 있듯이, GMDH 방법을 이용하여 시뮬레이션을 하나의 함수로 추정한 것을 유용하게 활용할 수 있음을 알 수 있다. 물론, 시뮬레이션 함수로 추정하는 데에는 두 가지 전제조건이 있다. 하나는 시뮬레이션 자체가 통제변수에 대해서 연속 함수(continuous function)이어야 한다는 것과, 시뮬레이션이 미분불가능인 함수인 경우에 그것을 미분가능한 함수로 추정하면서 약간 스무딩을 해준다는 것이다.
본 연구에서는 시뮬레이션 함수 추정을 GMDH 방법으로 하고, 그 결과를 활용하여 시스템 설계 운영에 활용하고자 하는 것이다. 이러한 연구는, 현재 활용되고 있는 RISKOptimizerⓡ 같은 소프트웨어의 성능을 더욱 개선시킬 수 있다. 또, 위험 관리를 위한 시뮬레이션에서 최적의 모수를 찾는 데에도 활용할 수 있다. 뉴런 넷 응용과 병행하거나 대체 방법으로도 활용할 수 있다. 최근, 뉴런넷 방법이나 RISKOptimizerⓡ과 같이 엑셀 기반하의 응용 소프트웨어가 개발되어 공급되고 있다. 이러한 소프트웨어에 이러한 함수추정 기능이 추가된다면 분석이 더욱 확대될 것으로 기대할 수 있다.
VII. 참고문헌
R. R. Barton & J. S. Ivey, Jr, "Nelder-Mead Simplex Modification for Simulation Optimization," Management Science, Vol 42, No. 7, July 1996
S. J. Farlow(edited), Self-Organizing Method in Modeling: GMDH Type Algorithm, Marcell Dekker, Inc., 1984:
J. N. D. Gupta, et al., "Selecting scheduling Heuristics Using Neural networks'" INFORMS Journal on Computing, Vol. 17, No.2, Spring 2000)
D. A. Humphrey & J. Wilson의 "A Revised Simplex Search Procedure for Stochastic Simulation Response Surface Optimization, INFORMS Journal on Computing, Vol.12, No.4, Fall, 2000;
A. G. Ivakhnenko, "Group method of data handling - a rival of of the method of stochastic approximation," Soviet Automatic Control, 13, 43 -71
A. G. Ivakhnenko, "Polynomia Theory of Complex System," IEEE Transaction on Systems, Man, and Cybernetics, Vol. SMC-1, No. 4, October, 1971
H. R. Madala & A. G. Ivakhnenko, Inductive Learning Algorithm for complex Systems Modeling, CRC Press, 1994
Palisade, Guide to RISKOptimizerⓡ: Simulation Optimizzation for Microsoft Excel, Palisade, May, 2000
Cliff T. Ragsdale, Spreadsheet Modeling and Decision Analysis, South-Western Publishing, 1998, 2001
W. Winston, Decision Making Under Uncertainty with RISKOptimizer, Palisade, 2000
Julian Miller, "KnowledgeMiner," http://www.knowledgeminer.net)
본 논문에 언급된 엑셀 파일은 http://home.donga.ac.kr/~hssong 에 있음.