1. Introduction
안녕하세요 박세훈입니다. 지난번 포스팅한 MLE는 여러가지 용도로 쓰입니다. 그 중 하나는 GLM 회귀분석에서 모수 추정치를 찾을 때 MLE를 사용한 Hessian 행렬이 매우 중요한 경우가 있죠.
Hessian 행렬은 다변량 편미분 값을 행렬로 나타내거나, 최적값 근처의 Log-Likelihood 표면의 국소적인 형태를 나타내기도 합니다. Hessian을 이용하여 모수의 공분산 행렬을 추정할 수도 있고, 모수의 표준오차 추정치를 얻기도 하죠. 분석가들은 최적의 솔루션에서 Hessian Matrix를 얻는 방법을 묻곤 합니다.
이번 장에서는 MLE Solution에서 Hessian Matrix를 얻는 3가지 방법을 포스팅 해보겠습니다.
2. Method
- PROC PLM의 Show Hessian 구문을 사용하여 Hessian 행렬 표현
- MODEL 구문에서 COVB 옵션을 사용하고, 공분산행렬의 역을 이용하여 Hessian 표현
- 일반 회귀 모형을 해결할 수 있는 PROC NLMIXED을 사용하고, HESS옵션으로 표현
3. data
Data는 지난 PROC PLM포스팅에서 사용했던 Neuralgia 로, PainModel 이라는 이름으로 Store된 데이터셋을 사용하겠습니다. 자세한 사항은 <code1>을 참고하시면 됩니다.
<code1>
Data Neuralgia;
input Treatment $ Sex $ Age Duration Pain $ @@;
PainY = (Pain^='Yes');
datalines;
P F 68 1 No B M 74 16 No P F 67 30 No P M 66 26 Yes B F 67 28 No B F 77 16 No
A F 71 12 No B F 72 50 No B F 76 9 Yes A M 71 17 Yes A F 63 27 No A F 69 18 Yes
B F 66 12 No A M 62 42 No P F 64 1 Yes A F 64 17 No P M 74 4 No A F 72 25 No
P M 70 1 Yes B M 66 19 No B M 59 29 No A F 64 30 No A M 70 28 No A M 69 1 No
B F 78 1 No P M 83 1 Yes B F 69 42 No B M 75 30 Yes P M 77 29 Yes P F 79 20 Yes
A M 70 12 No A F 69 12 No B F 65 14 No B M 70 1 No B M 67 23 No A M 76 25 Yes
P M 78 12 Yes B M 77 1 Yes B F 69 24 No P M 66 4 Yes P F 65 29 No P M 60 26 Yes
A M 78 15 Yes B M 75 21 Yes A F 67 11 No P F 72 27 No P F 70 13 Yes A M 75 6 Yes
B F 65 7 No P F 68 27 Yes P M 68 11 Yes P M 67 17 Yes B M 70 22 No A M 65 15 No
P F 67 1 Yes A M 67 10 No P F 72 11 Yes A F 74 1 No B M 80 21 Yes A F 69 3 No
;
/* Store문을 사용하여 모델 저장하기*/
title 'Logistic Model on Neuralgia';
proc logistic data=Neuralgia;
class Sex Treatment;
model Pain(Event='Yes')= Sex Age Duration Treatment;
store PainModel / label='Neuralgia Study';
ods select ParameterEstimates;
run;
4. Creating Hessain Matrix
다양한 방법을 사용하여 Hessian Matrix를 추출해보겠습니다. 참고로 공분산 행렬 추정치와 Hessian Matrix는 역의 관계에 있습니다.
4.1 PROC PLM 방법
전통적 방식인 PROC MEANS를 사용해보겠습니다.
<code2>
proc plm restore=PainModel; /* Restore문을 사용하여 모델 꺼내오기*/
show Hessian CovB;
ods output Cov=CovB;
run;
공교롭게도 Store문을 사용했을 때 모든 모델에 Hessian Matrix가 연산 되는 것은 아닙니다. 불가능한 PROC PLM에 Hessian을 요청할 경우 오류가 발생하게 되죠.
<TABLE_1> - PROC PLM 문에 의해 생성된 “CovB” 데이터 셋
4.2 COVB option 사용
대부분의 SAS Regression PROC은 COVB 옵션을 제공합니다. COVB를 통해 공분산 행렬을 구한 후 INV 함수를 이용하여 Hessian 행렬을 계산합니다. 데이터는 위 4.1에서 추출한 Covb 행렬을 사용하겠습니다.
<code3>
proc iml;
use CovB nobs p; /* 데이터 읽기*/
cols = "Col1"😞"Col"+strip(char(p)));/* 변수 이름 지정*/
read all var cols into Cov; /* <table_1>의COVB matrix 읽기 */
read all var "Parameter"; /* 파라미터 이름 불러오기*/
close;
Hessian = inv(Cov); /* 역행렬 계산*/
print Hessian[r=Parameter c=Cols F=BestD8.4];
quit;
Cov에 inv 함수를 적용하여 Hessian을 계산했습니다. 4.1의 결과와 동일합니다.
공분산 행렬은 분류변수의 매개변수에 따라 달라지는데, logistic PROC은 일반적으로 Effect를 사용합니다. Reference등의 매개변수를 사용하면 다른 결과를 얻을 수도 있습니다.
4.3 나만의 Log-Likelihood 함수 정의하기.
<함수 정의 절차>
1. Outdesign 옵션을 사용하여 설계 매트릭스를 생성
2. NLMIXED에 적용하기 위해 숫자형 변수로 인코딩
3. PROC NLMIXED 적용을 위해 로그 유사성 Binary 함수 관점에서 로지스틱 회귀모형 정의
<code4>
/* 매개변수화 및 아웃풋 설계(design)*/
proc logistic data=Neuralgia outdesign=Design outdesignonly;
class Pain Sex Treatment;
model Pain(Event='Yes')= Sex Age Duration Treatment;
run;
data Design;
set Design;
PainY = (Pain='Yes'); /* 문자형 변수를 더미변수화 */
run;
ods exclude IterHistory;
proc nlmixed data=Design COV HESS;
parms b0 -18 bSexF bAge bDuration bTreatmentA bTreatmentB 0;
eta = b0 + bSexF*SexF + bAge*Age + bDuration*Duration +
bTreatmentA*TreatmentA + bTreatmentB*TreatmentB;
p = logistic(eta); /* p 대신 1-p 를 사용해도 됩니다. */
model PainY ~ binary(p);
run;
결과를 위의 4.1 , 4.2와 비교하였을 때, 작은 수의 단위를 제외하고는 상당히 근사한 것을 확인 할 수가 있습니다. 저희가 직접 설계하고 추출한 행렬이 PROC LOGISTIC으로 부터 계산된 행렬과 상당히 유사함을 보입니다.
4. Conclusion
우리는 MLE 추정치에 대해 최적으로 Hessian을 구하는 3가지 방법을 실습 해보았습니다. 고급 통계분석이나 고급 프로그래밍을 진행할 때, 이처럼 세부적인 파라미터를 직접 계산하고 추출해야 하는 과정을 거쳐야할 때가 있습니다. 단순히 모델을 이해하는 것 보다 모델이 작동되는 원리를 파악하고 세밀한 핸들링을 할 수 있다면, 더 성장한 분석가가 될 수 있지않나 하는 생각이 듭니다.
Join us for SAS Innovate 2025, our biggest and most exciting global event of the year, in Orlando, FL, from May 6-9. Sign up by March 14 for just $795.