이원 로지스틱 회귀에 대한 개념설명은 http://statnmath.tistory.com/86 여기 참조하세요~ 

이번 포스팅은 이원 로지스틱 회귀에 대한 두번째 포스팅으로 모델에 대해 정리해보고자 합니다.

이원 로지스틱 모델은 다음과 같아요.

$Y_{i}|X_{i}$ = 1 if response is in category of interest, 0 otherwise. (Binary!)

where $Y_{i}|X_{i}$  is distributed by Bernoulli ($\pi_{i}$), where $\pi_{i}$ is a probability of success.

 

즉 일반 선형 모델 (General linear Model)과 다르게 일반화 선형 모델(Generalized Linear Model)이기 때문에 Y값이 continuous하지 않습니다. 이해가 잘 안가시면 그 전 포스팅 참조하세요. 이원 로지스틱이라서 y값이 1 혹은 0을 가지게 되죠. 이때 y는 Bernoulli distribution을 띄기 때문에 y의 평균과 분산 역시 Bernoulli distribution을 따르게 됩니다.  

Then $E[Y_{i}|X_{i}]=\pi_{i}$, and $Var[Y_{i}|X_{i}]=\pi_{i} \cdot (1-\pi_{i})$

이원 로지스틱 회귀 조건 중에 참조할 사항으로 분산이 일반 선형 모델과 다르게 일정하지 않다라고 했었죠. 

 

일반화 선형 모델은 선형 모델을 만들어주기 위한 link function이 필요합니다. 그래서 logit과 explanatory variable (X 변수) 가 linear 관계가 되는거죠. 

Logistic Regression Model : $log(\frac{\pi}{1-\pi})= \beta_{0}+\beta_{1}X_{1}+...+\beta_{p}X_{p}$

Logistic Function : $\pi=\frac{\exp(\beta_{0}+\beta_{1}X_{1}+...+\beta_{p}X_{p})}{1+\exp(\beta_{0}+\beta_{1}X_{1}+...+\beta_{p}X_{p})}$, where $\exp(\beta_{0}+\beta_{1}X_{1}+...+\beta_{p}X_{p})\in (-\infty, \infty )$

 

그렇다면 logistic regression model에서 beta값은 어떻게 구할까요? 바로 MLE를 통해 구하게 됩니다  

Data is $Y_{i}$ where 1 if response is in category of interest, 0 otherwise.

Model : $P(Y_{i}=y_{i})=\pi_{i}^{y_{i}}\cdot (1-\pi_{i})^{1-y_{i}}$ 

Joint Density: $P(Y_{i}=y_{i},...,Y_{n}=y_{n})=\prod_{i=1}^{n}\pi_{i}^{y_{i}}\cdot (1-\pi_{i})^{1-y_{i}}$

           where  $\pi=\frac{\exp(\beta_{0}+\beta_{1}X_{1}+...+\beta_{p}X_{p})}{1+\exp(\beta_{0}+\beta_{1}X_{1}+...+\beta_{p}X_{p})}$ and $1-\pi_{i}= \frac{1}{1+\exp(\beta_{0}+\beta_{1}X_{1}+...+\beta_{p}X_{p})}$

Joint density 값을 구하기 위해서 데이타 값이 independent해야해요. (서로 결과값에 영향을 미치면 안됩니다.)  

 

Joint Density를 구하면 베타에 대한 likelihood 을 구할 수 있고요~ 

Likelihood Function will be $L(\beta_{0},...,\beta_{p})= \prod_{i=1}^{n}\pi_{i}(\beta)^{y_{i}}(1-\pi_{i}(\beta))^{1-y_{i}}$

Log 를 씌운 뒤 최대값을 구할 수 있게 되지요. $(\hat{\beta_{0}},...,\hat{\beta_{p}})=\arg \max \left \{ \log L(\beta_{0},...,\beta_{p}) \right \}$

하지만 이 값이 딱 떨어지는 값이 아니게 됩니다. 그래서 Newton-Raphson algorithm 이나 혹은 Fisher Scoring Algorithm을 통해 맥시멈 값을 구하게 됩니다. 물론 컴퓨터가 계산해주죠. 단 한가지 명심해야할 것은, MLE를 통해 구하는거는 데이타 샘플이 충분히 있어야 한다는 점이예요. 그래야 이 beta값이 신뢰할 수 있다라고 합니다. 이해 안가시면 댓글 주세요~

 

짧게 정리했는데 다음엔 Wald test와 Likelihood Ratio test에 대해 정리해보도록 할게요.

http://statnmath.blogspot.ca/2015/08/the-binary-logistic-regression.html 참고하세요~

 

 

 

 

 

[4] Wald Procedure

 

How can we know whether a explanatory has effect on log-odds or not?

 

We can use Wald procedure to test whether $\beta's$ are zero or not!!

 

 

 

Hypothesis : $H_{0}:\beta_{j}=0$ which means $X_{j}$ has no effect on log-odds!!, $H_{1}:\beta_{j}\neq 0$

 

Test Statistics : $Z_{obs}=\frac{\hat{\beta_{j}}}{se(\hat{\beta_{j}})}$

 

                        where $\hat{\beta_{j}}$ is a maximum likelihood estimate.

 

Note that if there is large enough sample, MLE's are normally distributed so that under t$H_{0}$, our test statistics, $Z_{obs}$, is an observation from an approximate Normal(0,1) distribution!!

 

 

 

95% Confidence interval : $\hat{\beta_{j}}\pm 1.96 \cdot se(\hat{\beta_{j}})$

 

 

 

 

 

[5] Likelihood Ratio Test

 

In order to compare which model is appropriate, we can use likelihood ratio test.

 

Likelihood Ratio : $\frac{L_{R}}{L_{F}}$, where $L_{R}$ is reduced model, and $L_{F}$ is full model of same data.

 

 

 

Hypothesis :

 

$H_{0}: \beta_{1}=...=\beta_{k}=0$ Reduced model is appropriate so that it fits data as well as full model.

 

$H_{1}:$ at least one $\beta_{1}=...=\beta_{k}\neq 0$

이원 로지스틱 회귀에 대한 개념설명은 http://statnmath.tistory.com/86 여기 참조하세요~ 

 

Test Statistics : $G^2=-2 \log L_{R}-(-2\log L_{F})=-2 \log \frac{L_R}{L_F}$

 

Note that, under the null hypothesis, $G^2$ is an observation from a chi-square distribution with k degrees of freedom for large n! k is the number of parameter fewer in reduced model.

 

 

 

 

 

Remark!! Both Wald Test and Likelihood Ratio Test are testing whether $\beta$ parameter is zero or not! But those two procedures are different so that each following distribution is also different!! Just in case they do not agree, we use the Likelihood Ratio Test since it is more reliable.

 

반응형

이번 포스팅은 이원 로지스틱 회귀(Binary Logistic Regression)에 대해 짧게 정리해보도록 하겠습니다. 

 

이 이원 로지스틱 회귀는 Generalized Linear Model의 한 예로 볼 수 있는데요. General Linear Model (일반선형모델)과 다르게 일반화선형모델은 Y와 X간의 linear 관계가 아닙니다. 그래서 link function으로 linear 관계를 만들어주게 됩니다.  이원 로지스틱 회귀에서는 logit이라는 link function이 사용되는데 이건 나중에 정리해보도록 할게요.

 

일반 선형 모델은 y값이 continuous인데요. 이원 로지스틱 회귀에서 y값은 binary입니다. 즉 1 아니면 0 이예요. 성공 아니면 실패, 합격 아니면 불합격 이런식의 두가지로 이뤄진 데이타 입니다. y값이 이렇게 두 가지로 이뤄져있기 때문에 Bernoulli distribution을 따르게 됩니다.

 

그럼 이원 로지스틱 회귀 조건에 대해 정리해볼게요.

우선 데이타 y값끼리 independent 해야하고요. (결과값이 서로 영향을 미치지 말아야하고요)

그리고 샘플 데이타가 충분히 많이 있어야 신뢰할만한 추론을 할 수 있어요. 왜냐하면 일반선형모델과 다르게 일반화선형모델에서 수학적인 식을 만들때 MLE (Maximum Likelihood Estimate)를 사용하게 됩니다. 이때 MLE의 조건이 unbiased 하고 정규분포를 띄기 위해서 샘플사이즈가 충분히 있어야해요. 그래서 일반 선형 모델과 다르게 일반화선형모델은 충분한 샘플 사이즈가 필요하다고 하는거예요~

 

일반 선형 모델은 outlier가 있는지 살펴보았잖아요~ 일반화 선형 모델은 그럴필요가 없어요. 왜냐!! 어짜피 Y값이 1아니면 0 이렇게 둘 중 하나라서 outlier라는 자체가 없습니다.

 

그리고 일반 선형 모델은 variance가 줄어들거나 늘어나거나 하지 않고 그룹간의 비슷한 variance를 가져야한다고 했는데요. 이원 로지스틱 회귀에서는 variance가 일정하지 않습니다. 왜냐!! Bernoulli distribution에서 분산 구할때 식을 생각해보면 성공할 확률 곱하게 실패할 확률이잖아요. 여기도 마찬가지입니다. 성공할 확률x실패할 확률로 분산이 계산되기때문에 explanatory variable에 따라 y가 성공할때 확률이 달라지므로 분산값이 일정하지 않게 됩니다.

 

 

다음 포스팅에서는 이원 로지스틱 회귀에 대한 모델과 Wald Test & Likelihood Ratio Test에 대해 정리해보도록 할게요.

http://statnmath.blogspot.ca/2015/08/the-binary-logistic-regression.html 참고하세요~

반응형

http://statnmath.blogspot.ca/2015/08/case-study-pygmalion-effect-in-sas-two.html

자세한 내용은 위 링크 참조하세요. 이해 안가는 부분은 댓글 남겨주세요. 짧게 요약한거라 많은 내용이 빠져있긴 합니다.  

 

Case Study: The Pygmalion Effect

Reference : Ramsey and Schafer (1997) The Statistical Sleuth, Duxbury Press, p. 365. The Pygmalion effect.

예제내용을 짧게 설명하자면, 10개의 군부대가 있고 한개의 군부대에는 3개의 소대가 있습니다. 3 소개 중 하나의 소대 리더에게 너네 소대가 킹왕짱이라고 말을 해둡니다. 나머지 두개 소대에는 아무말을 하지 않고요. 그 후, 기본적인 무기 교육을 실시한 다음 테스트를 해서 성적을 가지고 10개의 군부대가 차이가 있는지 및 Pygmalion 효과가 있는지 살펴보는 실험입니다. 

이때 y (response variable)값은 무기 성적이 되겠고요~ x (explanatory variable)는 10개의 군부대와 treatment (한 소대에는 Pygmalion, 나머지 두 소대는 control groups) 그리고 이 군부대와 treatment간의 interaction까지 19개가 있게 됩니다. 데이타 모습은 아래와 같아요~ 

 

 

Two Way ANOVA 인 이유는 factor가 두개가 있기 때문입니다. factor 두개와 interaction까지 포함한 glm (General Linear Model)을 SAS로 돌려보도록 할게요. 그래서 모델은 다음과 같아요.

The model with interaction :

$Y_{i}=\beta_{0}+ \beta_{1}I_{pyg,i}+ \beta_{2}I_{comp1,i}+...+\beta_{10}I_{comp9,i}$

                                        $+\beta_{11}I_{pyg,i}\cdot I_{comp1,i}+...+\beta_{19}I_{pyg,i}\cdot I_{comp9,i}+e_{i}$

 

여기서 interaction term이 의미가 있는지 없는지부터 살펴봅시다. 가설은 다음과 같아요.

$H_{0}:\beta_{11}=\beta_{12}=...=\beta_{19}=0$ (all possible interactions)

$H_{1}$ : at least $\beta_{1i}$ is not equal to 0

 

* SAS Code and Result

Interaction을 만들려면, SAS code에 score = company | treatment 라고 적으면 됩니다.

혹은 score = company treatment company*treatment 라고 해도 돼요. SAS Code는 이미지사진으로 올릴게요.

 

 

이렇게 SAS를 돌리면 다음과 같은 ANOVA 값이 뜹니다. 

 

우리가 관심있는데 interaction을 보는거라서 company와 treatment간의 interaction의 P-value값을 보면 되어요. F값이 어떻게 계산되는지 적혀있고, 이해안가시면 댓글주세요. 이때 P-value값이 0.7221이라서 null hypothesis를 fail to reject하게 됩니다. 즉, interaction terms의 beta값이 모두 0으로 간주해도 된다는 말이예요. 그러니까 company와 treatment간의 interaction이 없다는걸 의미하고, 그렇다면 Y값에 영향을 주는게 company인지, treatment인지 알아봐야겠죠.

 

그래서 첫번째 질문은, treatment가 Y값에 영향을 주는지, 두번째 질문은 company가 Y값에 영향을 주는지 알아보려고 합니다. 마지막으로 이 데이타에서 assumption이 충족되었는지 살펴보도록 할게요. 그 뒤 내용은 링크 참조하시면 나와있습니다. 

 

http://statnmath.blogspot.ca/2015/08/case-study-pygmalion-effect-in-sas-two.html

읽어보시고 이해가 안가거나 궁금한 사항 있으면 댓글 남겨주세요. 감사합니다~

 

 

반응형

http://statnmath.blogspot.ca/2015/08/two-way-anova.html 자세한 내용은 링크 참조하세요.

[1] Two Way Classification or Two-Way Analysis of Variance

Two Way ANOVA는 general linear model (GLM 일반 선형 모델)에 속하고요. GLM을 잘 모르신다면 그 전 포스팅 참조하세요~ TWO WAY는 factor가 두개 있는 경우를 말합니다. 그래서 이 factor가 y값에 영향이 있는지 없는지 살펴보는게 바로 TWO WAY ANOVA예요.

factor는 categorical predictor variable인데요. 예를들어서 설명해보자면~ 임산부가 아이를 낳았을때 그 아기 몸무게를 예측해보려고 합니다. factor로는 임산부의 흡연여부, 임산부의 나이, 임산부의 몸무게 등등 여러 요인이 있겠죠. 이렇게 두개의 factor가 있다면 Two Way ANOVA를 이용할 수 있어요. 이때 신생아 아이의 몸무게는 continuous variable이고요. factor 역시 continuous variable이 가능합니다. 혹은 임산부의 흡연여부처럼 categorical variable도 가능하고요. 

 

[2] Assumptions

The samples must be independent, and selected by randomization condition. The equal variance assumption and normal error assumption should be satisfied.

Assumption 과 table 설명은 건너띄겠습니다. 이해 안가시면 댓글 주셔요 :D  

 

[3] Model and the Expected Values Example  

Consider the model for a two-way analysis of variance with two levels of each factor (a 2x2) classification. $Y_{i}=\beta_{0}+\beta_{1}I_{factor 1,i}+ \beta_{2}I_{factor 2,i}+\beta_{3}I_{factor 1,i}I_{factor 2,i}+e_{i}$ where $I_{factor 1,i}$ if the ith observation is in the first group of factor 1 and is 0 otherwise. 

질문은 TWO WAY ANOVA일때 평균 Y값을 찾는건데요. 답변은 링크 참조하세요. http://statnmath.blogspot.ca/2015/08/two-way-anova.html 설명은 건너띄지만 혹시 이해 안가시면 댓글 달아주셔요.

 

[4] The Full Model & The Reduced Model and Test Statistics

The Full Model은 explanatory variables을 다 포함한걸 Full Model이라고 합니다. 반면 Reduced Model은 우리가 관심없는 explanatory variable 몇가지를 제외한 모델을 말해요. 다른말로는 Additive Model이라고도 합니다.  

$F_{obs}= \frac{(SStrt_{full}-SStrt_{reduced}) / \#of\beta's \ being \ tested }{MSE_{full}}$

 

다음 포스팅엔 예를들면서 설명해볼게요. :D

 

 

 

반응형

http://statnmath.blogspot.ca/2015/08/case-study-spock-conspiracy-trial-in_8.html 자세한 내용은 링크 참조하세요.

 

데이타는 The Spock Conspiracy Trial 입니다. 출처는 여기 참조하세요. Reference: http://www.inside-r.org/node/159733

목적은 6 judges 들의 평균이 차이가 있는지 살펴보려고 합니다. One-way ANVOA를 이용할 예정이고 (SAS: proc glm), 6개의 그룹이 있기에 pairwise 비교가 있으므로 The Bonferroni Adjustment 명령어를 넣으려고 합니다.

 

1. 가설과 모델 

$H_{0}=\mu_A = \mu_B= ...=\mu_E = \mu_F$

$H_{1}$ = At least one judge is different from others.

Model : $Y_{i}=\beta_0 + \beta_{1} \cdot I_{A_{,i}} + \beta_{2} \cdot I_{B_{,i}} + ... + \beta_{5} \cdot I_{F_{,i}} + e_{i}$, 여기서 I는 indicator variable입니다.

 

[2] SAS Code

 

set all 에서 all은 제가 그 전에 이 데이타 불렀을때 이름이었어요. 참고로 infile로 데이타를 부르고요. 제가 그 전에 SPOCKS과 OTHERS로 분류를 해놓았습니다. 이 데이타에서 SPOCKS를 제외한 OTHER 그룹에는 A부터 F까지 judges가 있는데요. 각각 indicator variable로 분류했습니다. 

 

그리고 proc glm을 이용해서 One Way ANOVA를 돌려봅니다. A부터 F까지 pairwise로 비교를 하다보니 (즉, A와 B, B와 C, C와 D...등등) adjust=bon 명령어를 추가로 넣었습니다.

 

 

*cldiff statement reference: ttps://support.sas.com/documentation/cdl/en/statug/63347/HTML/default/viewer.htm#statug_glm_sect018.htm

*pdiff statement reference : http://support.sas.com/documentation/cdl/en/statug/63347/HTML/default/viewer.htm#statug_glm_sect016.htm

 

 

[3.1] SAS Result Overall ANOVA 

F값이 1.22 인데 계산은 65.45789/53.591280으로 계산되어서 나온 값이고요. Variance의 ratio를 말합니다. 이해가 안가시면 제가 그 전에 포스팅해놓은걸 참고하세요. http://statnmath.tistory.com/81 P-value 값이 0.3239이므로 0.05보다 크기때문에 null hypothesis가 fail to reject 됩니다.

 

The null hypothesis가 A부터 F까지 mean에서 차이가 없다! 였잖아요. 그런데 P value값이 이 가설이 아니라고 충분한 증거가 없다라는 겁니다. 즉, A부터 F까지 mean의 차이가 없다라고 결론내릴 수 있겠네요. 만약 차이가 있었으면 어디와 어디가 차이가 있는지 pairwise로 비교해서 알아볼 수 있는데 그게 아니므로 여기서 결론맺으면 됩니다.

 

 

어짜피 The Bonferroni correction 명령어를 넣었으니 정말 차이가 없는지 보도록 할게요.

 

 

[3.2] SAS Result - LSMeans and Difference Matrix with the Bonferroni Correction  

LSMeans을 보면 각 그룹의 평균값이 나오고요. Difference Matrix를 보면 pairwise로 비교한 p-value값이 나와요. 여기서 null hypothesis (H0)을 보면 i번째의 그룹과 j번째의 그룹의 평균이 같다라고 적혀있네요. 그리고 매트릭스로 P-value값이 나와있고요. 모든값이 1이네요. 즉 0.05보다 크니까 null hypothesis가 fail to reject 하게 됩니다. 즉, pairwise로 비교해도 서로의 평균이 다르다는 충분한 증거가 없다라는 의미예요.  

 

도움되셨나요? >_<

 

 

반응형

+ Recent posts