2016년 4월 30일 토요일

아빠가 들려 주는 [통계] 생존분석에서 샘플 수 계산


오늘 보라매 병원에서 생존분석 강의를 했습니다.
지난번 인천에서도 생존분석 강의를 했는데,
그 때는 주로 다변수 분석인 Cox regression 강의를 했었구요.

오늘은 RCT 그러니까 log-rank test를 주로 사용하는 것에 대해서 강의했구요.
통계적 방법은 사실 오히려 쉽습니다.

생존분석에서 샘플 수 계산하는 것이 조금 어렵고 이해가 잘 되는 면이 있어서
이부분 강의를 다시 녹화해 보았습니다.


다음.

그 다음..

동영상 마지막 부분에 두 집단의 비율이 2가 아니라, 0.5로 바뀌어야 되는 것같은데..
어쨌든.. 그것은 나중에 수정하기로 하고.
전체적 방법이 어떤 식인지는 이해 되셨을 것으로 생각됩니다.

2016년 4월 27일 수요일

2016년 4월 26일 화요일

mulptiple test experiment :

multiple test experiment

아래는 '다중성 문제'가 어떤 것인지 실험적으로 알려주고자 만든 자료입니다.

아래의 두개의 열은 평균 0 표준편차 1인 하나의 집단에서 무작위 추출된 것입니다.
각 집단은 100개씩의 샘플로 되어 있고요.

당연히 한 집단에서 추출된 것이니 두 집단을 t-test 하면 차이가 없어야죠.
즉 p>0.05 이어야 겠죠.

그런데 p의 특성상 여러번 반복하다 보면
우연에 의해 차이가 커지는 경우가 생기고 p값은 작아집니다.

대략 100번정도 해보면 5번 정도가 그러하죠.
그것이 p값의 특징입니다.

우선 한번 실험해 보세요.
1을 2로 바꾸었다가, 다시 1로 바꾸었다가 다시 바꾸어 보세요.
그 때마다 p값이 계속 바뀝니다.

그러다가 어쩌다 한두번.. p값이 정말 작아질 때가 있을 겁니다.


multiple test experiment

우리가 여러 가지 변수를 두고 검정하게 되면,
아무런 차이가 없을지라도 우연에 의해 p값이 작아지는
경우가 발생합니다.

A수술과 B수술을 하고서,
이것저것 기능점수, 각도차이, 크기차이, 통증차이 등등을 조사합니다.
그 때 어떤 하나가 p=0.008 이 나왔다고 합시다.
그런 경우에 어떻게 해석해야 할까요?

잘 생각해 보면,
p값이 무엇인지?
1종오류, 2종오류가 무엇인지.

다중성이 왜 문제가 되는지 아실 수 있으실 것입니다.

=============================
이렇게 한번 p=0.008이 나왔다고 해도
다른 사람이 실험하면 다시 p 값이 크게 나올 가능성이 훨씬 많아집니다.

그것이 p값의 성격이고,
다중 검정을 했기 때문에 발생하는 문제점입니다.

이것이 재현가능하지 않은 연구의 주된 원인 중의 하나라고 저는 생각합니다.



2016년 4월 22일 금요일

2X5 Bang's Blinding Index

2X5 Bang's Blinding Index

무작위 대조연구에서 중요한 것 중의 하나가 얼마나 가림 blinding이 잘 되었는가
하는 이슈이죠.
Blinding Index 중에서 간단한 2X3으로 하는 것은 지난 번에 올린 적이 있고요,
이번에는 2X5 로 된 것 올리도록 하겠습니다.

역시 동일하게 노란 칸에만 입력하는 것입니다.
계산이 좀 복잡했습니다. 엑셀 아래쪽으로 내려보면 어떻게 계산했는지..
아실 수 있으실 것입니다.



 혹시, R로 계산하시고 싶으신 분은..
아래의 것을 복사해서 사용하셔도 될 것같습니다.
저도 이거 받아서 이걸로 엑셀로 전환한 것입니다.
========================
 #2by3;
bifun<-function(n11,n12,n13){
n1<-n11+n12+n13
bi<-(2*(n11/(n11+n12))-1)*(n11+n12)/n1
p1.1<-n11/n1
p2.1<-n12/n1
var.bi<-(p1.1*(1-p1.1)+p2.1*(1-p2.1)+2*p1.1*p2.1)/n1
se.bi<-sqrt(var.bi)
ub<-bi+1.96*se.bi
lb<-bi-1.96*se.bi
#return(bi,lb,ub,se.bi)
print(bi)
print(lb)
print(ub)
#print(se.bi)
}
bifun(n11=82,n12=25,n13=170)
bifun(n11=29,n12=27,n13=83)




#2by5+2x2;
bifun<-function(n11,n12,n13,n14,n15,nt,np,w11,w12,w13,w14,wt,wp){
n1<-n11+n12+n13+n14+n15
bi<-(w11*n11+w12*n12+wt*nt+wp*np+w13*n13+w14*n14)/n1
p1.1<-n11/n1
p2.1<-n12/n1
p3.1<-n13/n1
p4.1<-n14/n1
pt.1<-nt/n1
pp.1<-np/n1
var.bi<-(w11**2*p1.1*(1-p1.1)+w12**2*p2.1*(1-p2.1)+w13**2*p3.1*(1-p3.1)+w14**2*p4.1*(1-p4.1)+wt**2*pt.1*(1-pt.1)+wp**2*pp.1*(1-pp.1)
+2*w11*w12*p1.1*p2.1+2*w11*w13*p1.1*p3.1+2*w11*w14*p1.1*p4.1+2*w11*wt*p1.1*pt.1+2*w11*wp*p1.1*pp.1
                    +2*w12*w13*p2.1*p3.1+2*w12*w14*p2.1*p4.1+2*w12*wt*p2.1*pt.1+2*w12*wp*p2.1*pp.1
                                        +2*w13*w14*p3.1*p4.1+2*w13*wt*p3.1*pt.1+2*w13*wp*p3.1*pp.1
                                                            +2*w14*wt*p4.1*pt.1+2*w14*wp*p4.1*pp.1
                                                                               +2*wt*wp*pt.1*pp.1
)
var.bi<-var.bi/n1
se.bi<-sqrt(var.bi)
ub<-bi+1.96*se.bi
lb<-bi-1.96*se.bi
print(bi)
print(lb)
print(ub)
#print(se.bi)
}
bifun(n11=38,n12=44,n13=21,n14=4,n15=170,nt=79,np=86,w11=1,w12=0.5,w13=-0.5,w14=-1,wt=0.25,wp=-0.25)
bifun(n11=8,n12=21,n13=16,n14=11,n15=83,nt=45,np=36,w11=1,w12=0.5,w13=-0.5,w14=-1,wt=0.25,wp=-0.25)



2016년 4월 13일 수요일

아빠가 들려 주는 [통계] 데이터가 많아지면, 설명력이 올라갈까?

 
이런 질문을 한번 해 볼께요, 너무 당연한 결론이지만, 가끔 착각하는 경우가 있기도 해서요.
 
기울기 3 절편 3
잔차가 정규분포로 분포하는 모집단을 만들자

 
무작위로 10개의 샘플을 뽑았더니,
P=0.035087로 기울기가 0이 아닐 것으로 생각되고,
기울기의 값은 1.7정도
. 원래값인 3과는 꽤 거리가 있군요.
그런데, 기울기의 95%신뢰구간이 0.15에서 3.33까지이므로,
어쨌든 그 범위에는 포함되는 군요.
이렇게 큰 범위로 추정한다면 추정하나 마나 한 것아닌가 하는 생각이 듭니다.
이렇게 샘플 수가 작으면 95%신뢰구간도 넓어지죠.
설명력은 44%정도, 예상보다는 꽤 나쁘지만, 사실 의학이나, 경제학이나 실제에서는 이보다 훨씬 적게 나오기 쉽죠.
 
이제 20개의 샘플을 뽑았습니다. P=0.0000000136 무지 작아졌습니다.
기울기도 많이 근접했군요. 설명력도 84%나 됩니다.

 
이제는 30개를 뽑았습니다.
P값은 이제 세기도 어려울 정도로 작아졌고요.
기울기도 3.5, 20개로 추정했던 3.02보다 참값에 더 멀어지긴 했지만.
95%신뢰구간 안에는 있습니다.
, 샘플 수가 늘어나면서 95%신뢰구간이 점점 줄어 들지만,
더 많은 샘플 수가 더 정확한 점추정을 하지는 않군요.
그렇지만, 더 많은 샘플 수를 가기게 되면 이 구간이 점점 더 줄어서 추정하기는 좋겠군요.
설명력은 85%입니다. 꽤 높군요.
 
이제 샘플 수를 확 늘여서 121개로 잡았습니다.
P값은 훨씬 적어지고, 95% 구간도 줄어서 3개 가까워 지고 있습니다.
설명력은 86%..
. 이렇게 점점 늘여 가다보면,
기울기는 점차 참값에 가깝게 추정할 수 있겠군요.


 
이제 484개의 샘플을 추출해서 실험해 보았더니, p값은 상상할 수 없을 정도로 작아졌습니다.
95% 신뢰구간은 어느 정도 줄어 들었고요,
그런데, 설명력은 83%로 거의 늘고 있지 않습니다.
왜 그럴까요?
이 회귀모형이 원래 설명할 수 있는 것이 83% 정도밖에 안되기 때문에 아무리 샘플을 모은다 하더라도 100%가 될 수 없습니다.


 
공식을 굳이 말하고 싶진 않지만, 최소한의 것이라
정도는 언급해야 할 것같군요.
설명력의 공식은 이와 같습니다.
수학식만 나오면 거부감이 있는 분을 위해서
이 공식은 이어서 설명하겠습니다.


 
여러 개의 관측된 점들이 있다고 치고, 그리고, 빨간 점선의 회귀직선.
파란 직선의 Y값들의 평균입니다.


 
평균이 무엇일까요,
특별한 이유가 없다면, 우리는 모두 키가 같고, 몸무게가 같을 것입니다. 차이가 나야할 이유가 없으니까요.(이론상)
그런데, 유전적인 요인, 후천적인 요인, 등등 온갖 요인에 의해서
키가 달라집니다.
그 달라짐의 정도가 큰 사각형이다.
달라짐을 제곱한 것이지요.
작은 사각형은 무엇인가요.
우리가 구한 회귀식에 의해 달라진 정도입니다.
우리가 구한 회귀식은
몇 가지 측정된 변수들에 의해 설명되는 것입니다.
즉 설명이 가능한 정도이다.


 
이렇게 모든 점들에게서 큰 사각형과 작은 사각형을 그려 볼 수 있을 것입니다.
큰 사각형은 점들이 변화 정도,
작은 사각형은 회귀선에 의해서 설명되는 정도.


 
큰 사각형들의 합을 분모로 놓고,
작은 사각형의 합을 분자로 놓으면
전체 변화량에 대해서 회귀식에 의한 변화량을 구하면,
(‘무엇에 대해서이렇게 하면 무엇을 분모에 놓는 거죠)
이것을 설명력이라고 말합니다.
그래서 전체의 사각형을 분포로, 회귀식의 사각형을 분자로 놓은 것이 설명력이 되고,
그것을 수학적으로 표현한 것이 위의 공식이 됩니다.

만일 모든 점들이 회귀식 위에 올라간다면, 큰 사각과 작은 사각이 크기가 같아질 것이고,
1이 될 것입니다.
회귀직선이 수평선이 되면 빨간 사각형들이 모두 0 되어서 설명력은 0 됩니다.
 
위키피디아에 나오는 그림인데, 찬찬히 보면 모두 설명이 가능합니다.
1행의 것들이 왜 직선에서 멀어질 때 점점 설명력이 0에 가까워 지는지 이해되시죠?
그리고,
2행은 모두 직선이기 때문에 모두 1입니다. 기울기와 상관없이..
단 가운데 있는 것은 분모가 0이 되어서 뭐라고 말할 수가 없습니다.
3행의 것들은 모두 회귀식이 수평이라 모두 0이 됩니다.
 
다시 원래의 예로 돌아가 봅시다.
샘플 수가 많아지면, 우리는 P값이 점점 작아진다는 것을 확인했습니다.
그렇지만, 샘플 수가 많아진다고 해서,
점들이 점점 직선에 가까워 지지 않습니다.
그렇기 때문에 설명력은 점점 1에 수렴하는 것이 아닙니다.
원래의 퍼진 정도에 있을 뿐이죠.
지금은 단순히 X가 하나뿐이지만, X가 여러 개로 늘어난다면,
즉 설명 변수가 더 많아진다면, 직선위로 모여질까요?
그렇지 않습니다.
 
설명력을 더 높이기 위해서는 변수를 더 찾아 내는 것도 하나의 방법이지만, 잔차 residual 또는 오차 error를 줄이는 것이 중요한 방법입니다.
이 잔차라는 개념을 생각한 사람은 정말 천재. 아마도 피어슨이 만든 것같은데……
가급적 잔차를 줄이기 위해서 조금이라도 더 정확한 값을 찾는 것이 설명력을 높이는 방법이 될 것입니다.
예를 들면 나이를 측정할 때 그냥 5년 단위로나 1년 단위로 할 수도 있고, 개월로도 할 수 있겠지만,
조금이라도 더 정확한 값을 찾아 보는 것이죠.
수입이라면 대충의 값보다는 더 정확한 값을,
체중도 마찬가지…….
그렇게 하더라도 실제로 측정하지 못한 값도 많기 때문에 설명력에는 한계가 있습니다.
 
그런데, 만일 년도를 날짜 단위로 측정한다면 정말 설명력이 더 올라갈까?
해보진 않았지만, 조금 아주 아주 조금 올라갈 것같습니다.
나이 요인이 전체에서 얼마나 될까?
마치 600화소 카메라과 1200화소 카메라가 확대해 보지 않으면
우리 눈에 차이를 주지 못하는 것과도 같이 더 정밀하게 조사해도 차이가 미비할 수 있습니다.
(300DPI 이상은 실제적인 차이가 없다고들 말한다)
그렇지만, 연구하는 입장에서 가능한 정밀하고 정확하게 조사하는 것이
잔차를 줄이고, 설명력을 높이는 방법입니다.
또는 숫자를 엄청 많이 늘인다면,
아마도 p값이 엄청 작아지면서,
무의미했던 변수들이 유의미하게(p값이 작아지는) 되는 효과는 있을지언정
설명력의 변화는 크지 않을 가능성이 있겠지요.
숫자를 아주 많이 하면 그동안 발견되지 않은 유의미한 변수를 찾아 낼 것같은 느낌이 들지만.
정말 중요한 변수들은 숫자가 적을 때도 발견이 되지요.
즉 통계적으로는 유의미한(p가 작은) 변수이지만,
임상적으로는 무의미한(설명력이 거의 없는) 변수일 가능성이 많습니다.

Mata 분석을 좀 싫어 하는 분은 이런 말을 하기도 합니다.
그 비판 중에 하나가, 숫자를 아주 많이 모아서 유의미한 결과를 만들어 낸다는 것이 과연 진짜 얼마나 중요한 것일까.. 하는 것이죠.
(참고로 저는 meta 분석을 싫어하진 않습니다. 그런 면도 고려해야 한다는 정도…)
빅데이터에 대해서 비판적인 근거 중의 하나도 이런 면이 있습니다.
천명 정도의 데이터에서 무의미하게 나온 변수가 1억명의 빅데이터로 유의미하게 되었다면,
과연 임상적으로 의미가 있을까 하는 회의론
그런데 한편, 작은 차이라도 판매에 영향을 미치는 것을 찾고자 한다면 그것도 의미는 있을 수 있지요.
참 오늘의 주제는 메타분석이나 빅데이터 이야기가 아니고, 설명력에 관한 것이었습니다.
설명력이 무엇인가와 어떻게 올릴 수 있는가 하는 이야기

2016년 4월 7일 목요일

아빠가 들려 주는 [통계] Multiple P values : Bonferroni, Holm, Hochberg

 
여러 번의 통계적 방법을 사용한 경우에는 다중성 문제를 해결하기 위해 전체의 1종 오류를 5%로 통제하기 위해서 p값의 경계를 조절하게 됩니다.
만일 조절하지 않게 되면 1종 오류가 매우 커지게 되지요.
그 방법 중에 가장 많이 알려진 방법이 Bonferroni의 방법이며 쉽습니다.
그런데 너무 엄격하고 보수적이라서
Holm-Bonferroni method의 방법과 Hochberg-Bonferroni method의 방법이 사용되기도 합니다.
그것에 대해서 알아 보겠습니다

To solve the multiple test problem, we should control the margin of p values.
If we do not, the type 1 error( alpha error) will increase dramatically.
The most common method would be the Bonferroni method, which is considered too much conservative.
So there are two modifications from the Bonferroni method.
Those are ‘Holm-Bonferroni method’ and ‘Hochberg-Bonferroni method’.
I will explain what they are.

처음 만든 것은 201496일이었는데, 엑셀 파일로 보급했었지요.
오늘 다시 웹 버전으로 만들면서 조금 더 편하게 볼 수 있도록 수정했습니다

This excel sheet was made and distributed first at September 6th, 2014
And now I made some modification and distribute again at web-version.

 
만약 여러분이 20개의 p값을 가지고 있다면, 그것을 작은 순서부터, 큰 순서로 정렬하여
노란 칸에 입력합니다.
엑셀에 붙여넣기를 한다면 더 좋습니다.

If you have 20 P values, at first arrange them from the smallest to the largest.
And copy and paste them into the yellow cells.
I recommend to paste “values only”

그렇게 하면 오른쪽에 보이는 그래프를 얻게 될 것입니다.
제일 아래쪽에 보이는 보라색 선이 Bonferroni의 경계선으로 이것보다 아래에 있는 것만
유의하다고 판단합니다.
Bonferroni의 방법에 의하면 말이지요

You will get the chart as above.
The purple line is the line for the Bonferroni method.
If the p vales are under the purple line, you can consider they are significant according to the Bonferroni method.

 
연두색 선은 변형된 것입니다. 이 연두색과 p값들과 비교하게 되는데,
왼쪽에서 읽어서 즉 가장 작은 것에서부터 읽어서 일단 한번이라도 연두색 선 위로 올라오게 되면 그 오른쪽은 모두 not significant 한 것으로 간주하는 방법이 
Holm-Bonferroni method입니다

Now we will see the green line. This is the modified cut-line.
From the left side(that is from the smallest ones) once p value rise above the green line, that and the larger( right sided ) are considered as not significant.
That is the Holm-Bonferroni method.

반대로 오른쪽에서부터 읽어서 즉 가장 큰 값부터 읽어서 한번이라도 연두색 선 아래로 내려오면
그보다 왼쪽에 있는 것은 모두 significant 한 것으로 간주하는 방법이
Hochberg-Bonferroni method입니다

From the right side(that is from the largest ones) once p value merge below the green line, that and the smaller( left sided ) are considered as significant.
That is the Hochberg-Bonferroni method.

그래서 Holm-Bonferroni methodHochberg-Bonferroni method보다 더 보수적인 방법이 됩니다.
As a result, the Holm-Bonferroni method is conservative than the Hochberg-Bonferroni method.



그리고 여러분은 이 화면을 보게 될 것입니다.
제일 오른쪽에 차트가 있습니다.
가운데 3개의 빨간색 문자열을 보게 될 것입니다.
각각 Bonferroni, Holm, Hochberg 방법의 결과들입니다.
어떤 것이 가장 보수적인지 금방 할 수 있을 것입니다

You will see this sheet in reality.
The chart will be right.
There are three red columns.
Relatively Bonferroni, Holm, Hochberg’s.
You will notice which is the most conservative.