glmnet part 04 - Poisson Models

Glmnet Vignette

Trevor Hastie and Junyang Qian
Stanford September 13, 2016

포아송 모델 (Poisson Models)

포아송 회귀(Poisson regression)는 포아송 에러의 추정 하에 셈 데이터의 모델을 구축하기 위해 사용한다, 또는 그렇지 않으면 평균과 분산이 비율인 경우 비-음수 데이터. 가우시안과 이항 모델 처럼, 포아송은 분포가 지수인 군의 일원이다. 우리는 로그 척도 logμ(x) = β0 + β'x 에서 그것의 양의 평균으로 일반적으로 모델을 생성한다. 관측 {xi, yi}1N의 로그-우도는 다음으로 주어진다:
이전 처럼, 우리는 패널티 로그-우도를 최적화한다:
Glmnet은 outer Newton loop를 사용한다, 그리고 이 기준을 최적화하기 위해 내부 가중치 최소-제곱 루프(inner weighted least-squares loop) (로지스틱 회구에서 처럼).
먼저, 우리는 포아송 데이터의 미리 생성된 집합을 적재한다. 
data(PoissonExample)
우리는 "poisson" 옵션을 가지는 함수 glmnet을 적용한다. 
fit = glmnet(x, y, family = "poisson")
"poisson" 군에서 glmnet의 옵션 입력 인수는 다른 것에서 그것들과 비슷한다.
offset은 포아송 모델에서 특히 유용한 인수이다.
포아송 모델에서 비율 데이터를 처리할 때, 수집된 셈(count)은 종종 다른 노출에 기반한다, 관측된 시간의 길이, 넓이(area) 그리도 년 수(years). 포아송 비율 μ(x)는 단위 노출 시간에 연관한다, 그래서 만약 관측 yi가 시간 단위 Ei에서 노출 된다면, 그러면 기대 셈(count)은 Eiμ(x)가 될 것이다, 그리고 로그 평균은 log(Ei) + log(μ(x))가 될 것이다. 이 같은 경우에, 우리는 각 관측에서 offset log(Ei)를 적용할 것이다. 그러므로 offset은 선형 예측자를 포함하는 길이 nobs의 벡터이다. 다른 군은 또한 이 옵션을 사용할 수 있다, 전형적으로 다른 이류로.
(경고: offset이 glmnet에 제공된다면, offsets은 합리적인 예측을 생성하기 위해 또한 predict에 적용되어야 한다.)
다시, 우리는 결과의 첫 느낌을 가지고 계수를 그린다. 
plot(fit)
이전 처럼, 우리는 계수를 추출하고 그리고 coef 그리고 predict 각각 사용하여 어떤 λ에서 예측을 만든다. 옵션 입력 인수는 다른 군에서와 비슷하다. 함수 predict 에서, 옵션 type, 요구하는 예측의 타입, 포아송 군에서 그 자신이 특별함을 가진다. 그것은, * "link" (기본)은 다른 것들 처럼 선형 예측자를 제공한다, *"response"는 적합된 평균을 제공, *"coefficients"는 s의 요구 값에서 계수를 계산, coef 함수 * "nonzero"가 s의 각 값에서 0이 아닌 색인의 목록을 반환하여 실제화할 수 있음.
예로, 우리는 아래를 수행할 수 있다. 
coef(fit, s = 1)

## 21 x 1 sparse Matrix of class "dgCMatrix"
##                       1
## (Intercept)  0.61123371
## V1           0.45819758
## V2          -0.77060709
## V3           1.34015128
## V4           0.04350500
## V5          -0.20325967
## V6           .
## V7           .
## V8           .
## V9           .
## V10          .
## V11          .
## V12          0.01816309
## V13          .
## V14          .
## V15          .
## V16          .
## V17          .
## V18          .
## V19          .
## V20          .


predict(fit, newx = x[1:5,], type = "response", s = c(0.1,1))

##               1          2
## [1,]  2.4944232  4.4263365
## [2,] 10.3513120 11.0586174
## [3,]  0.1179704  0.1781626
## [4,]  0.9713412  1.6828778
## [5,]  1.1133472  1.9934537
우리는 최적의 λ를 찾기 위해 교차-검증을 또한 사용할 수 있다 그리고 추론한다. 
cvfit = cv.glmnet(x, y, family = "poisson")
옵션은 type.measure을 제외하고 가우시안 군과 거의 같다, *"deviance" (기본)은 편차를 제공, *"mse" 평균 제곱 에러를 나타냄, *"mae"는 평균 절대 에러이다.
우리는 cv.glmnet를 그릴 수 있다. 
plot(cvfit)
우리는 최적 λ 그리고 해당하는 계수를 볼 수 있다 . 
opt.lam = c(cvfit$lambda.min, cvfit$lambda.1se)
coef(cvfit, s = opt.lam)

## 21 x 2 sparse Matrix of class "dgCMatrix"
##                        1            2
## (Intercept)  0.031262578  0.185696196
## V1           0.619053001  0.575373801
## V2          -0.984550161 -0.932121975
## V3           1.525234026  1.470567302
## V4           0.231590890  0.196923579
## V5          -0.336659414 -0.304694503
## V6           0.001025945  .
## V7          -0.012829626  .
## V8           .            .
## V9           .            .
## V10          0.015983078  .
## V11          .            .
## V12          0.030866986  0.025850171
## V13         -0.027970571  .
## V14          0.032750270  .
## V15         -0.005932709  .
## V16          0.017505507  .
## V17          .            .
## V18          0.004026211  .
## V19         -0.033579041  .
## V20          0.012048941  0.009929548
predict 메서드는 비슷하다 그리고 여기서 그것을 반복하지 않는다.

Cox Models

콕스 비례 위험 모델(Cox proportional hazards model)은 예측 변수와 생존 시간 사이에 관계의 연구에서 일반적으로 사용된다. 일반적인 생존 분석 프레임워크에서, 우리는 (y1, x1, δ1), ... , (yn, xn, δn) 형태의 데이터를 갖는다, yi, 관측 시간, 만약 δi가 1이면 실패 시간 또는 δi가 0이면 바른-검열(right-censoring)인 경우. 우리는 또한 t1 < t2 < ... < tm은 한 건의 실패가 발생한 증가하는 시간의 목록이다, 그리고 ji는 시간 ti에서 관측한 실패의 색인을 나타낸다.
콕스 모델은 위험의 준-모수적 형태(semi-parametric form)를 추정한다.
hi(t)는 시간 t에서 환자 i의 위험, h0(t)는 공통 기준 위험(shared baseline hazard), 그리고 β는 고정값, 길이 p 벡터인 경우. 전통적인 설정 n ≥ p인 경우에, 추론은 부분 우도(partial likelihood)로 한다:
Ri는 yj ≥ ti를 가지는 색인의 집합인 경우(그러므로 시간 ti에서 위험).
콕스 모드에서는 교차점(intercept)이 없음을 주목하라 (기본 위험에서 그것의 생성, 그리고 그것 처럼, 부분 우도에서 취소될 것이다).
우리는 표본 데이터의 미리 생성된 집합 그리고 응답을 사용한다. 사용자는 그들 자신의 데이터를 적재할 수 있다 그리고 비슷한 절차를 따른다. 이 경우에 x는 공변량(covariate) 값의 n × p 행렬이 되어야 한다 - 각 행은 환자에 해당하고 그리고 각 열은 공변량. y는 n × 2 행렬, 실패/검열(failure/censoring) 시간의 열 "“time"을 가짐, 그리고 "status" 0/1 지시자, 1이 의미하는 시간은 실패 시간, 그리고 0는 검열(censoring) 시간. 
data(CoxExample)
y[1:5,]


##            time status
## [1,] 1.76877757      1
## [2,] 0.54528404      1
## [3,] 0.04485918      0
## [4,] 0.85032298      0
## [5,] 0.61488426      1
패키지 survival에 있는 Surv 함수는 그런 행렬을 생성할 수 있다. 주의, 어쨌든, coxph 그리고 연관된 선형 모델은 구간(interval)을 처리할 수 있다 그리고 다른 검열의 형태, glmnet은 그것의 표현 형태로 오직 바른-검열(right censoring)만 처리할 수 있다.
우리는 기본 설정 하에 해결 경로를 계산하기 위해 glmnet 함수를 적용한다. 
fit = glmnet(x, y, family = "cox")
모든 표준 옵션은 alpha, weights, nlambda 그리고 standardize 처럼 사용 가능하다. 그들의 사용은 가우시안에서 처럼 비슷하다 그리고 우리는 여기서 상세를 뺐다. 사용자는 도움말 파일 help(glmnet)에서 또한 참조할 수 있다.
우리는 계수를 그릴 수 있다. 
plot(fit)
이전 처럼, 우리는 λ의 어떤 값에서 계수를 추출할 수 있다. 
coef(fit, s = 0.05)

## 30 x 1 sparse Matrix of class "dgCMatrix"
##               1
## V1   0.37693638
## V2  -0.09547797
## V3  -0.13595972
## V4   0.09814146
## V5  -0.11437545
## V6  -0.38898545
## V7   0.24291400
## V8   0.03647596
## V9   0.34739813
## V10  0.03865115
## V11  .
## V12  .
## V13  .
## V14  .
## V15  .
## V16  .
## V17  .
## V18  .
## V19  .
## V20  .
## V21  .
## V22  .
## V23  .
## V24  .
## V25  .
## V26  .
## V27  .
## V28  .
## V29  .
## V30  .
콕스 모델이 예측에서 일반적으로 사용되지 않기 때문에, 우리는 예측에 설명 예제를 제공하지 않는다. 만약 필요하면, 사용자는 help(predict.glmnet)을 입력하여 도움말을 참조할 수 있다.
또한, 함수 cv.glmnet는 콕스 모델에서 k-배 교차-검증을 계산하기 위해 사용할 수 있다. 사용법은 2개의 주요 다른 것을 제외하고 다른 군에서 그것과 비슷하다.
하나는 type.measure는 오직 "deviance"(기본)만 지원한다는 것이다, 부분-우도를 제공함.
다른 하나는 옵션 grouped에 있다. grouped = TRUE는 빼기(subtraction)로 k번째 배(fold)에서 CV 부분 우도를 얻을 수 있다; 전체 데이터 셋에서 추정된 로그 부분 우도에서 (K-1)/K 테이터셋에서 추정된 것을 뺄샘으로. grouped=FALSE 가지는 로그 부분 우도는 오로지 k번째 배(fold)에서 계산된다, 만약 각 fold가 큰 관측 수를 가진다면 유일하게 합리적이다. 
cvfit = cv.glmnet(x, y, family = "cox")
일단 적합한다, 우리는 최적의 λ 값을 볼 수 있다 그리고 교차 검증 에러로 우리의 모델을 평가하기 위해 그릴 수 있다. 
plot(cvfit)
이전 처럼, 우리의 그림에서 왼쪽 수직선은 우리에게 CV-error 곡선은 그것의 최소값과 만나는 경우를 보여준다. 오른쪽 수직선은 우리에게 최소의 1 표준 편차 안에 CV-error를 가지는 가장 정규화된 모델을 모여준다. 우리는 또한 그런 최적 λ를 추출한다. 
cvfit$lambda.min
## [1] 0.01594373

cvfit$lambda.1se
## [1] 0.04868986
우리는 모델에서 활성 공변량(active covariates)을 확인한다 그리고 그들의 계수를 본다. 
coef.min = coef(cvfit, s = "lambda.min")
active.min = which(coef.min != 0)
index.min = coef.min[active.min]

index.min

##  [1]  0.491297030 -0.174600537 -0.218648771  0.175112406 -0.186673099
##  [6] -0.490250398  0.335197365  0.091586519  0.450169329  0.115921793
## [11]  0.017594799 -0.018364649 -0.002805776 -0.001422613 -0.023429230
## [16]  0.001687662 -0.008236046

coef.min

## 30 x 1 sparse Matrix of class "dgCMatrix"
##                1
## V1   0.491297030
## V2  -0.174600537
## V3  -0.218648771
## V4   0.175112406
## V5  -0.186673099
## V6  -0.490250398
## V7   0.335197365
## V8   0.091586519
## V9   0.450169329
## V10  0.115921793
## V11  .
## V12  .
## V13  0.017594799
## V14  .
## V15  .
## V16  .
## V17 -0.018364649
## V18  .
## V19  .
## V20  .
## V21 -0.002805776
## V22 -0.001422613
## V23  .
## V24  .
## V25 -0.023429230
## V26  .
## V27  0.001687662
## V28  .
## V29  .
## V30 -0.008236046

댓글 없음:

댓글 쓰기

홈
피드 구독하기: 글 (Atom)