국가재무제표를 처음으로 작성했더니 나랏빚이 300조 증가했다는 기사가 나왔습니다. 하지만 기사를 읽어 보아도 경제학에 무지한 저는 무슨 소린지 도통 모르겠습니다. 그래서 제가 구할 수 있는 국가통계포털 KOSIS 자료 ¹를 좀 들여다보았습니다.

2012년 1/4 분기를 기준으로 우리나라는 외국에서 받을 돈(대외채권)이 외국에 갚을 돈(대외채무)보다 많은 순채권국입니다. 그런데 자료를 좀 더 살펴보만면 통화 당국(중앙은행)을 제외한 나머지 일반정부, 예금취급기관, 기타 부분은 모두 외국에 갚은 돈이 더 많습니다. 이중에서 정부의 순대외채권의 추이를 그려보면 다음과 같습니다. ²

정부의 재정이 흑자였던 때도 있긴 했군요. 현 정권의 추이는 그렇다 하더라도 2007년의 낙하 폭에 눈이 가는군요. 무슨 일이 있었던 건지…. 또 한 가지 흥미로운 것은 정권 이양 시에는 정부 빚에 그렇게 큰 차이가 없다는 점입니다. 앞으로 약 10개월, 더도 말고 앞의 두 정권의 수준만큼이라도 회복되길 바라는 것은 허무맹랑한 욕심일지도 모르겠네요.

아무튼, 신문의 기사는 여기에 300조 원을 더해서 생각하면 된다는 이야기라고 이해하고 있습니다만 아니라면 지적해 주세요. ^^

이글과 관련된 글

• 2012/05/24 -- 우리나라 공기업의 순대외채권 추이, 1994~2012
• 2012/05/30 -- 2012 프로야구 각 팀의 표정을 보니…
• 2012/05/23 -- 당신의 삶 나아지셨나요?
• 2011/08/29 -- Visualize This
• 2011/05/29 -- OECD에서 평가한 우리나라의 웰빙지수는?

올해 프로야구가 개막되고 팀별로 약 40게임을 소화했습니다. 어찌 보면 야구만큼 통계와 밀접한 관계가 있는 스포츠도 드문 것 같습니다.
타율, 방어율, 출루율, 승률 등 이름에서도 알 수 있듯이 선수 개개인의 성적은 확률로 나타내는 경우가 많죠.

오늘은 각 팀의 성적, 기록을 다변량 데이터의 시각화 방법의 한가지인 체르노프의 얼굴 그림(Chernoff face)를 이용해 그려보았습니다. 체르노프의 얼굴 그림은 각 변수의 크기를 얼굴 각 부분의 길이와 넓이로 표시하는 방법입니다. 이 방법은 변수의 순서에 따라 그림이 바뀌게 되므로 같은 데이터를 이용한다 하더라도 다양한 결과를 얻을 수 있는 특징이 있습니다.

위 그림은 승률이 미소를 표현하게 해놓았으므로 승률이 높으면 웃는 표정이 되겠죠.
어디까지나 재미로 그려본 것이니 좋아하는 팀 표정이 이상하다고 항의하지는 마세요 ^^.

이글과 관련된 글

• 2012/05/31 -- 정부의 순대외채권 추이, 1994~2012
• 2012/05/24 -- 우리나라 공기업의 순대외채권 추이, 1994~2012
• 2012/05/23 -- 당신의 삶 나아지셨나요?
• 2011/08/29 -- Visualize This
• 2011/05/29 -- OECD에서 평가한 우리나라의 웰빙지수는?
• 2011/05/23 -- 1980년대와 2000년대 프로야구 투수들의 방어율 어느 쪽이 낮을까?

잃어버린 10년… 돌려받고 싶습니다.

이글과 관련된 글

• 2012/05/31 -- 정부의 순대외채권 추이, 1994~2012
• 2012/05/30 -- 2012 프로야구 각 팀의 표정을 보니…
• 2012/05/23 -- 당신의 삶 나아지셨나요?
• 2011/08/29 -- Visualize This
• 2011/05/29 -- OECD에서 평가한 우리나라의 웰빙지수는?

경제협력개발기구(OECD)가 2012년 판 ‘Your Better Life Index‘를 발표했습니다. 이 지수는 OECD 회원국의 각종 통계 자료들을 모아 나라별 웰빙지수를 주거, 소득, 일자리, 공동체, 교육, 환경, 시민참여, 건강, 삶의 만족도, 안전성, 일과 삶의 균형의 총 11개 분야로 측정하고 종합한 지수입니다. 각 분야에는 세부 측정항목이 있는데 올해 특히 주목할 점은 각 항목별로 성별, 빈부격차에 따른 지수의 차이가 추가되었다는 것입니다.

OECD 회원국 안에서 우리나라는 어떤 위치에 있을까요? 조사 결과에 따르면

• 우리나라는 한 사람이 한 해에 OECD 평균 22,387달러에 못 미치는 16,570달러를 벌고,
• 상위 20%의 수입이 하위 20% 수입의 5배가 넘는 심각한 빈부격차에 시달리고 있으며,
• 15~64세 국민의 63%가 직업을 가지고 있으며 (남: 74%, 여: 55%),
• OECD 연평균 1,749시간을 훨씬 웃도는 2,193 시간을 일하며,
• 25~64세의 80%가 고등학교를 졸업했고,
• OECD의 평가방법에 따른 수학과 과학의 점수가 평균보다 높고,
• 평균수명은 평균보다 한 살 많은 81세(남: 80세, 여: 84세)
• 대기 환경은 OECD 평균보다 안 좋고, 82%의 국민이 수질에 만족하고 있으며,
• 81%의 국민이 유사시 의지할 사람이 있다고 믿고 있으며,
• 소득 상위 20%의 투표율이 91%인데 반해 하위 20%의 투표율은 59%에 불과하며,
• 62%의 사람들이 매일 부정적 경험보다 긍정적 경험을 겪는다고 답했습니다.

Better Life Index는 위와 같은 여러 평가항목을 바탕으로 처음 소개한 11가지 분야에 대해 0~10점의 점수를 부과하여 우리의 삶을 평가하고 있습니다. 예를 들어 주거에 대한 평가가 완벽하다면 10점, 형편없다면 0점 이런 식이죠. 따라서 각 항목의 기대평균은 5점이라 할 수 있습니다.

작년 이맘때쯤에도 이 지수를 이용해 포스팅을 한적이 있습니다만, 이번 글에서는 작년과 올해의 평가 점수가 어떻게 달라졌는지 그래프로 만들어 보았습니다.

많은 항목에서 평가가 좋아졌지만, 아직 개선해야 할 부분들이 많은 것 같습니다. 이 글에서는 구체적으로 살펴보지 않았지만, 거의 모든 평가항목에서 성별에 따른 지수의 차이 혹은 빈부에 따른 지수의 차이가 OECD 회원국 중 최하위를 달리고 있습니다. 남녀가 평등한 사회, 모두가 잘 살 수 있는 그런 우리나라가 되길 바라지만 요즘의 정치, 경제, 사회가 돌아가는 모습을 보면 한숨만 나옵니다.

내년에는 좀 달라질까요?

이글과 관련된 글

• 2011/05/29 -- OECD에서 평가한 우리나라의 웰빙지수는?
• 2012/05/31 -- 정부의 순대외채권 추이, 1994~2012
• 2012/05/30 -- 2012 프로야구 각 팀의 표정을 보니…
• 2012/05/24 -- 우리나라 공기업의 순대외채권 추이, 1994~2012
• 2011/08/29 -- Visualize This

1년 전에 역대 정권별 소비자물가 상승률 추이라는 글을 올린 적이 있습니다. 이후 2011년 11월에 소비자물가지수 지표가 변경된 것은 많은 분이 아시리라 생각합니다. 통계청의 보도 자료를 보면 변경 이전에는 2005년의 물가를 100으로 보았을 때의 상대 값이었는데 변경 이후에는 2010년의 물가를 100으로 놓았을 때의 상대적인 값을 물가지수로 사용하게 됩니다. 그 외에도

• 조사대상 품목 수가 489개에서 481개로 변경
• 소비행태의 변화에 따른 품목 추가: 스마트폰 이용료, 떡볶이, 외식용 막걸리, 캠핑용품 등 43개 종목 추가
• 탈락 및 변경 품목: 금반지, 한복, 정수기, 캠코더, 전자사전 등 21개 종목 탈락
• 2010년 가계동향조사의 소비지출액 구성비에 따른 가중치의 재조정
• 경제협력개발기구(OECD) 방식을 도입해 품목별 가중치에 적용

등을 주요 내용으로 삼고 있습니다.

새 지수를 적용한 결과 이전 방식대로 계산할 때 4.4%이던 2011년 물가상승률이 4.0%로 떨어지게 됩니다. 작년 글에서도 언급했듯이 정부의 물가상승 억제선은 보통 4%가 됩니다. 이를 기준으로 성적을 매기게 되겠죠. 이번 변경 덕분에 참 이상하리만치 우연일지도 모르겠습니다만 작년 물가상승률은 더도 덜도 아닌 딱 4%로 역사에 남게 되었습니다. 이에 대해 언론에서도 통계청이 새 지수 적용을 11월로 앞당겨 적용해 ‘물가 꼼수’를 부리려는 것 아니냐는 의혹을 제기하기도 했습니다.

[물가통계 조사방식 개편] 새 지수 앞당겨 적용 ‘꼼수’ 논란 – 서울신문

그럼 소비자물가지수 지표 변경 전과 이후의 물가지수를 좀 더 자세하게 비교해 보면 어떻게 될까요?

지표 변경 발표 후 국가 통계 포털 및 e-나라지표의 소비자물가지수는 2010년을 기준으로 한 값으로 변환되었습니다. 그런데 다행(?)히도 제가 작년에 글을 쓰면서 사용했던 2005년 기준을 적용한 물가지수 데이터가 남아있어서 두 값을 한 번 비교해 보았습니다. 다만 아직 변환작업이 모두 안 끝난 것인지 아니면 변환 후에도 지수의 변동이 없었는지 확인할 수 없지만 2005년 이전과 2007~2009년의 소비자 물가지수(전년 같은 달 대비)의 값들은 변화가 없더군요.

위 그래프는 2003년 3월부터 2012년 3월까지 소비자물가지수(전년 같은 달 대비)의 추이를 그려본 그래프입니다. 변경 전의 지수를 실선으로 변경 후의 지수를 표현했습니다. 그래프에서도 확인할 수 있듯이 지표변경 후 지수가 높아진 때도 낮아진 때도 있습니다. 그런데 MB 정권에 들어서는 지표변경 후 소비자물가지수가 낮아진 달이 더 많아진 듯이 보입니다.

좀 더 디테일하게 살펴보죠. 다음 표는 지표변경 후 소비자물가지수가 상승한 달과 하락한 달을 세어본 것입니다. 단, 괄호 안의 %는 두 대통령의 전체 재임 기간 중 소비자물가지수에 변화가 있었던 개월 수의 비율을 의미합니다.

위 표에서 보는 바와 같이 지표변경 후 MB 정부에서 지수가 상승한 개월 수는 참여정부의 개월 수보다 그 비율이 낮고 반대로 하락한 개월 수는 참여정부보다 MB정부의 개월 수가 10% 이상 많음을 알 수 있습니다. MB정부 전체의 1/4이나 하락했네요. 두 대통령의 재임 기간 중 월별 소비자물가지수의 평균 변화를 보면 참여정부는 2.94에서 2.93으로 MB 정부는 3.65에서 3.56으로 변경되어 평균값으로 보아도 MB정부의 하락폭이 더 크다는 걸 알 수 있습니다.

또한, 소비자물가지수 상승률 4%를 기준으로 했을 때 현 정부의 4% 이상인 월 수와 미만인 월 수의 변화를 세어보면 다음 표와 같습니다.

이렇게 전후 관계가 있는 데이터에 대한 통계 검정방법으로 McNemar의 검정법이 있습니다. 위 데이터를 이용해 검정을 해보면 유의확률(p-value) 0.0143으로 유의수준 5%에서 지표변경 전과 변경 후에 4% 이상인 월 수와 미만인 월 수의 변화에는 통계적으로 유의한 차이가 있음을 확인할 수 있었습니다. 약간 과장해서 말하자면 소비자물가지수 상승률 4%를 기준으로 보았을 때 지표변경 전의 지수와 변경 후의 지수는 성격이 다르다고 할 수 있겠죠.

물론 개편으로 인한 소비자물가 하락 효과는 1991년 -0.3%, 1996년 -0.1%, 2001년 -0.3%, 2006년 -0.2% 등으로 이전의 개편에서도 있었다고 통계청의 보도자료는 밝히고 있습니다만 이번 개편은 유독 MB 정권 하에서의 물가지수가 혜택을 받고 있다는 느낌을 지울 수가 없습니다.

게다가 2012년 3월의 소비자물가지수 상승률이 2.6%로 발표된 데 대해 실제로는 3.2%라는 민간 연구기관의 발표가 있었고, 이는 통계청의 공식 트윗 계정에서도 인정했습니다.

통계청이 발표한 3월 소비자 물가 상승이 2.6%가 아니라 3.2%라는 민간연구기관의 보고서 내용은 사실, 통계청(bit.ly/HteoL0)과 기획재정부(bit.ly/HqvFiB)에서 이미 설명했던 내용입니다.

— 통계청(StatisticsKorea) (@KOSTATIN) 4월 9, 2012

현실을 제대로 반영 못 하는 통계지표. 이런 지표를 이용해 목표달성 했다고 자화자찬하는 정부와 보수언론을 보면 그냥 찌질하게 보이는 건 저뿐일까요?

이상 McNemar 검정의 예제를 생각하다가 우발적으로 생각나 적어본 글이었습니다.

이글과 관련된 글

• 2011/05/06 -- 역대 정권별 소비자물가 상승률 추이

6. 맺음말

게놈정보는 유전계승법칙에 따라 다음 세대에 전달되는 안정된 형태의 데이터이기 때문에 통계적으로 유리하다고 말할 수 있습니다. 특히 게놈정보 중 다형성의 하나인 SNP에 주목해 모든 염색체상의 SNP좌위의 유전자형을 총체적으로 관측하여 형질과의 관련성을 검토하는 GWAS에 대해 앞선 8개의 글에서 소개하였습니다.

하지만 이번에 소개한 방법은 게놈 연구의 일부에 지나지 않습니다. GWAS는 형질과의 관련성을 총체적으로 분석하는 방법이지만 변이는 모두 공통이라는 전제(Common disease-common variant)를 필요로 하죠. 따라서 개인별로 변이가 서로 다른 형질의 원인을 GWAS로 찾아내는 것은 어렵습니다. 예를 들어 혈연관계가 있는 집단을 연구대상으로 삼을 때에는 연쇄분석(linkage analysis)을 하는 것이 일반적입니다. 연쇄분석도 총체적으로 유전적 리스크를 탐색하는 방법이긴 하지만 공통으로 가지고 있는 변이가 아닌 개인별로 서로 다른 변이도 검출할 수 있기 때문에 GWAS와는 구별되게 됩니다. 또한, 최근에는 차세대 시퀀싱(Next Generation Sequencer, NGS)라 부르는 생물의 전 염기배열 혹은 전 엑솜(exom) 영역의 염기배열을 관측하는 장치가 등장하여 생물의 모든 게놈정보를 지금까지와는 비교되지 않을 정도의 속도와 비용으로 얻을 수 있게 되었습니다. NGS 데이터의 분석방법에 대해서는 이후 기회가 있을 때 소개하도록 하겠습니다.

GWAS는 원인 유전자좌 영역의 추출 및 위치를 탐색하는 것이 목적이므로 그 결과에 대한 의학적, 생물학적 기능에 대한 해석이 불가능하다면 분석결과를 신약개발, 개인의료 등에 응용하기 어렵게 됩니다.

게다가 GWAS에서는 형질과 관련성이 있는 진짜 SNP를 탐색하는 것이 아니라 매우 가까운 위치에 있는 SNP를 탐지하게 되는 경우도 많다고 앞서 설명했습니다. 이 때문에 진짜 관련 유전자 좌를 찾아내기 위해 GWAS로 얻어낸 후보 SNP좌위의 근처에 있는 유전자 좌를 더 자세히 조사해 보아야할 필요가 있습니다. 이때 각 배우자(정자 및 난자)의 염기배열(haplotype)을 추정하여 형질과의 관계성을 조사하는 방법도 많이 사용됩니다.

또한, 사람뿐 아니라 농작물이나 가축의 게놈연구도 활발히 이루어지고 있으며 벼, 소 등의 몇몇 종에 대해서는 GWAS 데이터를 얻기 위한 DAN chip도 시판되고 있습니다. 그러나 농작물 및 가축의 게놈연구는 그리 쉽지가 않은데 이들 개체는 무작위 교배가 아닌 경우가 많기 때문입니다. 예를 들어 식용 소의 대부분은 인공수정을 하므로 같은 정자로부터 태어나는 송아지가 다수 존재하게 되고 농작물은 접목 등을 통해 클론을 손쉽게 만들 수 있습니다. 이 때문에 HWE 법칙이 성립하지 않을 때가 많아서 데이터의 품질평가가 불가능하여 분석결과의 해석이 어렵게 됩니다. 그리고 농작물은 사람과 달리 2배체(diploid) 염색체가 아닌 경우도 많으므로 사람의 데이터를 분석하기 위해 개발된 소프트웨어 등을 그대로 이용할 수 없다는 단점도 있습니다.

마지막으로 게놈데이터 분석을 위해서 유닉스(특히 리눅스) 환경은 이제 필수라고 해도 과언이 아닙니다. GWAS만 해도 파일 사이즈가 보통 수 GB ~ 수십GB의 데이터를 핸들링하게 되는데 스팩이 좋은 PC라고 해도 GUI 환경에서 이를 다루기는 무리가 있습니다.

분석에 이용하는 소프트웨어도 대부분 CUI 환경을 기반으로 하고 있고 윈도즈용은 아예 없을 때도 많습니다. 게다가 요즘 NGS 데이터 분석이 점점 주류가 되고 있는데 이런 데이터를 분석하기 위해서는 고성능 병렬 계산기와 대용량의 저장공간이 필요합니다. 즉, 요즘 유행하는 빅데이터 분석을 위한 도구가 필요하죠. 이를 윈도즈 환경에서 구축하기란 불가능하지 않을까요?

이것으로 GWAS로 배우는 유전통계학 연재를 마치도록 하겠습니다.
부족한 글 보아 주셔서 감사합니다. ^^;

이글과 관련된 글

• 2012/04/12 -- GWAS로 배우는 유전통계학 – 5 분석결과의 시각화
• 2012/04/11 -- GWAS로 배우는 유전통계학 – 4 다중비교 문제
• 2012/04/07 -- GWAS로 배우는 유전통계학 – 3.2 질적, 양적형질에 대한 관련분석
• 2012/04/05 -- GWAS로 배우는 유전통계학 – 3.1 코호트 연구와 실험-대조군 연구
• 2012/04/03 -- GWAS로 배우는 유전통계학 – 2.3 게놈 데이터의 품질평가
• 2012/03/31 -- GWAS로 배우는 유전통계학 – 2.1 가계모순과 하디-베인베르크 평형 법칙
• 2012/03/30 -- GWAS로 배우는 유전통계학 – 1. 시작하며

5. 분석결과의 시각화

GWAS로부터의 검정결과는 분석에 사용하는 DNA chip에 따라 차이가 있지만 보통 50만~150만의 p-값이 계산되므로 그 결과를 하나하나 확인하는 것은 사실상 불가능합니다. 따라서 먼저 시각적으로 분석결과를 확인하고 관련성이 있다고 판단된 SNP좌위의 정보를 확인하는 것이 일반적입니다. 분석결과의 시각화방법으로는 qq-plot(quantile-quantile plot)과 Manhattan plot이 많이 사용됩니다.

Quantile-Quantile plot

만약 분석에 사용된 모든 SNP에 대해 형질과의 관련성에 대해 검정을 할 때 관련성이 없다는 귀무가설이 바르다고 하면 모든 p-값은 0과 1 사이의 균일분포(uniform distribution)를 따르게 될 것입니다. 만약 관련성이 있다는 대립가설이 바르다고 한다면 그때의 p-값은 균일분포로 부터 벗어나게 됩니다.

관련성이 있는 SNP가 있다면 그때의 p-값은 이론값보다 작아짐

이 사실을 이용하여 이론적인 분포(균일분포)에서의 p-값과 실제 계산된 p-값을 그래프로 작성한 것이 qq-plot입니다. 즉, 관련성 검정의 대상이 되는 SNP 수를 n, 검정결과 i-번째로 작은 p-값을 \(p_{(i)}\)라고 하면
\[
\left(
-log_{10}\frac{i}{n},-log_{10}p_{(i)} \right) ,~i=1,\cdot,n
\]
을 그래프로 그리게 됩니다.

p-값이 작은 부분에서는 형질과 SNP사이에 관련성이 있다는 대립가설에 따른다고 예상되므로 아래 그림과 같이 붉은 원 안의 p-값에 대응하는 SNP가 관련성을 시사하는 SNP가 됩니다.

Manhattan plot

Manhattan plot은 관련분석결과 p-값을 염색체번호, 물리적 거리순으로 늘어놓은 플롯을 말합니다.
여기서 주의해야 할 점은 연쇄불평형(linkage disequilibrium)관계에 있는 SNP 사이의 관련분석 결과는 p-값이 비슷하게 됨에 주목하여 형질과 SNP와의 관련성에 대해 평가해야 합니다. 다시 말해 GWAS에 사용되는 DNA chip은 SNP의 밀도가 높으므로 비슷한 위치에 있는 SNP 사이에는 연쇄불평형 관계가 존재합니다. 이런 때에는 검정의 p-값이 비슷하여지므로 아래 그림과 같이 형질과 관련성이 있는 유전자 영역에 있는 SNP들의 p-값은 고층빌딩과 같이 불쑥 솟아오른 모양이 됩니다. ¹

참고문헌

1. Balding D.J. (2006), Nature reviews Genetics, 7, 10, 781-791.
2. 鎌谷直之 (2007) 遺伝統計学入門, 岩波書店 (카마타니 나오유키 (2007), 유전통계학 입문, 이와나미서점 )

이글과 관련된 글

• 2012/04/07 -- GWAS로 배우는 유전통계학 – 3.2 질적, 양적형질에 대한 관련분석
• 2012/04/05 -- GWAS로 배우는 유전통계학 – 3.1 코호트 연구와 실험-대조군 연구
• 2012/03/31 -- GWAS로 배우는 유전통계학 – 2.1 가계모순과 하디-베인베르크 평형 법칙
• 2012/03/30 -- GWAS로 배우는 유전통계학 – 1. 시작하며
• 2012/04/18 -- GWAS로 배우는 유전통계학 – 6. 맺음말
• 2012/04/11 -- GWAS로 배우는 유전통계학 – 4 다중비교 문제
• 2012/04/03 -- GWAS로 배우는 유전통계학 – 2.3 게놈 데이터의 품질평가

4. 다중비교(multiple comparison) 문제

GWAS에서는 보통 50만~250만 SNP를 이용해 관련분석을 하게 되므로 반드시 다중비교의 문제가 발생합니다. 하나의 SNP를 이용한 검정의 유의수준을 \(\alpha\)라고 한다면 한 번의 검정에서 \(\alpha \times 100 \)%의 확률로 잘못된 결론을 내리게 됩니다. 만약 50만 SNP좌위를 이용해 검정을 했을 때 단 한 번이라도 잘못된 결론을 내리게 될 확률, 즉 거짓 양성(false positive)은
\[
1-(1-\alpha)^{500K} \approx 1
\]
이 되어 100% 오류를 포함하게 되는 거죠. 이러한 문제를 개선하기 위해 매우 다양한 방법이 고안, 발표되고 있습니다. 이번 포스팅에서는 가장 간단한 방법인 Bonferroni의 보정방법과 FDR(false discovery rate)를 이용한 방법에 대해 알아보도록 하겠습니다.

Bonferroni 방법

Bonferroni의 부등식에 기초한 이 방법은 분석 전체의 유의수준을 모든 검정의 수로 나누어 이것을 1회의 검정에서의 유의수준으로 삼는 방법을 말합니다. 즉, 분석에 사용되는 SNP좌위 수를 \(k\)라고 하면 1회의 검정에서 사용하는 유의수준 \(\alpha^’\)은
\[
\alpha^' = \frac{\alpha}{k}
\]
이 되고 산출된 \(\alpha^’\)보다 작은 유의확률(p-value)이 관측된 SNP좌위에 대해 유의성을 인정하는 방법입니다. 예를 들어 100만 SNP좌위의 DNA chip을 이용하여 GWAS를 수행할 때 1회의 검정에 이용되는 유의수준은 \(5 \times 10^{-8} \)이 되게 됩니다.

이 방법은 계산이 간단하다는 장점이 있지만 모든 SNP가 독립이라는 다시 말해 연쇄 평형(linkage equilibrium)이라는 가정하에서 이루어지므로 1회의 검정에서 사용하는 유의수준이 너무나 작아져 검정력이 떨어지게 된다는 단점이 있습니다. 즉, 실제로 형질과 SNP 사이에 관련성이 있어도 관련이 없다고 판단하는 오류를 범할 확률이 커지게 됩니다.

FDR을 이용한 방법

Bonferroni 방법의 문제점을 개선하기 위해 여러 가지 대안이 있는데 그중에서 형질과 관련성이 있다고 판단된 SNP 중에서 잘못된 판단의 비율을 일정 비율 이하로 억제하는 FDR을 이용한 접근법이 있습니다.

n개의 SNP좌위를 이용해 실험-대조군 연구를 한다고 할 때 관련분석의 결과 유의성이 있는 SNP의 수를 R이라 하겠습니다. 모든 SNP가 형질과 관련이 없다(귀무가설)고 하면 R의 분포는 귀무가설 하에서
\[
R \sim B(n,\alpha)
\]
의 이항분포를 따르게 됩니다. 여기서 R개의 SNP 중에 진짜로 관련성이 있는 SNP 수를 S, 관련성이 없음에도 불구하고 관련이 있다고 판단된 SNP의 수를 V라고 해보죠. 또한, n개의 SNP중에서 진짜로 관련성이 있는 SNP의 비율을 \( \pi \)라고 한다면 관련이 없음에도 관련이 있다고 잘못 판단될 비율은 \( (1-\pi)\alpha \)가 되고, V의 기댓값은 \( n(1-\pi)\alpha \)가 됩니다. 이를 표로 그려보면 다음과 같이 됩니다.

FDR이라는 것은 위 표의 V/R의 기댓값을 말하게 됩니다. 따라서
\[
Q = \frac{n \alpha (1-\pi) }{R}
\]
이라 하면 이 값(Q-value)을 일정 값 이하로 제어하는 것이 이 방법의 목적입니다. 그러나 \( \pi \)의 값은 미지의 값(신의 영역이란 건 이럴 때 쓰는 거죠!!)이기 때문에 \(Q \leq n \alpha/R \)의 관계를 이용하여 Q-값을 컨트롤 하게 됩니다.

개인적으로는 FDR을 컨트롤하는 방법 중에 Benjamini-Hochberg법(BH법)을 즐겨 사용합니다. 이 방법은 각 검정으로부터의 유의확률을 크기순으로 늘어놓고 j번째로 작은 유의확률 \( p_{(j)} \)보다 작은 p 값을 가진 SNP는 관련성이 있다고 판단 할 때, 분석 전체의 유의수준을 \( \alpha \)라고 하면
\[
\frac{n(1-\pi)p_{(j)}}{j} \leq \frac{n p_{(j)}}{j}
\]
이므로
\[
\frac{n p_{(j)}}{j} \leq \alpha
\]
\[p_{(j)} \leq \alpha \times \frac{j}{n}
\]
의 관계를 이용해 p 값이 큰 것으로부터 평가하여 최초로 부등식이 성립되게 될 때, 이보다 작은 p 값을 가지는 SNP는 모두 유의성(관련성)이 있다고 판단하는 방법입니다.

이번 포스팅에서는 다중비교 문제에 대해 비교적 계산이 간단한 두 가지 보정방법에 대해 알아보았습니다만 가장 정확하다고 할 보정방법은 permutation test를 이용한 방법이라 할 수 있겠습니다. 그러나 현실적으로는 사용하기 불가능할 정도로 계산량이 많아서 이를 해결하기 위한 방법 또한 활발히 연구가 진행되고 있습니다.
요즘 제가 즐겨 쓰는 방법은 SLIDE(a Sliding-window approach for Locally Inter-correlated markers with asymptotic Distribution Errors corrected)라는 방법이 있는데 자세한 사항은 이곳을 참조하시길 바랍니다. 저자가 한국분이신 것 같네요.

참고문헌

1. Benjamini Y, and Hochberg Y. (1995) J. Roy. Stat. Soc. B., 57, 289-300
2. 鎌谷直之 (2007) 遺伝統計学入門, 岩波書店 (카마타니 나오유키 (2007), 유전통계학 입문, 이와나미서점 )

이글과 관련된 글

• 2012/03/31 -- GWAS로 배우는 유전통계학 – 2.1 가계모순과 하디-베인베르크 평형 법칙
• 2012/03/30 -- GWAS로 배우는 유전통계학 – 1. 시작하며
• 2012/04/18 -- GWAS로 배우는 유전통계학 – 6. 맺음말
• 2012/04/12 -- GWAS로 배우는 유전통계학 – 5 분석결과의 시각화
• 2012/04/07 -- GWAS로 배우는 유전통계학 – 3.2 질적, 양적형질에 대한 관련분석
• 2012/04/05 -- GWAS로 배우는 유전통계학 – 3.1 코호트 연구와 실험-대조군 연구
• 2012/04/03 -- GWAS로 배우는 유전통계학 – 2.3 게놈 데이터의 품질평가

3.2 질적 형질에 대한 관련분석

질적 형질에 대한 관련분석을 분할표를 이용한 Pearson의 카이제곱 검정이나 Fisher의 정확 검정법을 주로 이용합니다. 어떤 SNP 좌위에 대해 가장 기본적인 관측 데이터는 질적 형질의 표현형에 따른 유전자형의 도수겠죠. 많은 경우 질적 형질은 두 개의 카테고리를 가지므로 개체의 표현형을 D(disease)와 N(non-disease)라 하고 SNP의 allele를 A, a라고 한다면 표 1과 같은 분할표를 작성할 수 있습니다.

<표 1> 유전자형에 따른 돗수의 분할표

여기서 유전계승양식의 지식을 이용하면 보다 유전학에 따른 분석을 할 수 있게 됩니다. 예를 들어 allele A에 대해 우성양식을 가정한다면 표 1은 표 2와 같이 재구성할 수 있습니다.

<표 2> allele A에 대한 우성양식의 분할표

만약 allele A에 대해 열성양식을 가정한다면 표 3과 같은 분할표를 만들 수 있습니다.

<표 3> allele A에 대한 열성양식의 분할표

또한 allele의 도수를 두 군에 대해 비교하는 방법도 가능합니다(표 4).

<표 4> 표현형과 allele 돗수의 분할표

하지만 병에 걸리는 것은 allele가 아니고 개개인이 되므로 allele 빈도를 이용한 관련분석은 개체를 기초로 하는 분석이 아님에 주의해야 합니다. Allele 빈도를 이용한 관련분석은 검정에 사용되는 표본크기가 유전자형을 이용한 관련분석의 2배가 되므로 검정력(power of test)이 높아지게 됩니다.

우성, 열성, 유전자형, allele 빈도 이외에도 주목하는 allele의 수와 관측 유전자형 돗수사이의 경향성을 이용하여 관련성을 평가하는 방법도 있습니다. 즉, 개체가 보유하고 있는 관심 allele의 수가 질병의 리스크를 높이는가(혹은 낮추는가)에 주목하고 Armitage 검정을 이용하여 경향성의 유무를 평가는 방법입니다.

만약 분석 대상이 되는 표현형의 유전계승양식이 알려져지지 않았다면 GWAS에서는 위에서 설명한 5가지 양식 각각에 대한 관련성 평가를 하게 됩니다.

3.3 양적 형질에 대한 관련분석

질적 형질에 대한 관련분석에서는 유전계승양식을 가정하고 각 유전자형 빈도의 차이에 대해 검토하게 됩니다. 또한, 개체의 배경정보(환경정보)가 형질에 미치는 영향을 고려하기 위해 로지스틱 회귀모형 등의 통계모형을 도입해 유전적 요인의 탐색을 하기도 합니다. 실제로 분석 목적에 따라 나이, 성별, 체중, 키 등의 배경정보를 선택합니다.

양적 형질에 대해서도 마찬가지로 형질을 반응변수(목적변수)로 하고 배경정보와 게놈 정보를 설명변수로 하는 회귀모형을 이용한 관련분석을 하는 것이 일반적입니다.

from Balding (2006), Nat. Rev. Genet

예를 들어, 개체의 배경정보를 나이(age), 성별(gender)이라고 하면 양적 형질에 대한 선형모형은
\[
y= \beta_0 + \beta_1 Age + \beta_2 Gender + \beta_3 SNP + \epsilon
\]
과 같이 표현할 수 있습니다. 여기서 실제로 사용하는 SNP의 값은 유전계승양식에 따라 숫자로 코딩한 값을 입력합니다. 그리고, SNP에 대한 회귀계수에 주목하여 최소제곱법에 의해 추정된 \(\beta_3\) 값에 대한 평가를 합니다. 회귀계수에 대한 검정은 “추정된 회귀계수가 0 인가? (귀무가설)”, “0 이 아닌가? (대립가설)”에 대한 평가가 됩니다.

앞에서 통계모형을 도입해서 분석할 때, 유전자형을 유전계승양식에 따라 숫자로 코딩한 값을 사용한다 했는데 유전자형이 AA, Aa, aa이고 minor allele를 a라 하면

• 우성: minor allele에 대해 우성양식을 가정하면 AA, Aa, aa를 각각 0, 1, 1로 변환한 값을 모형에 사용합니다. 이때 추정되는 회귀계수는 하나이며 AA의 형질에 대한 영향을 0으로 가정했을 때 minor allele를 하나라도 보유하고 있는 개체(Aa or aa)의 형질에 대한 영향을 추정합니다.
• 열성: minor allele에 대해 열성양식을 가정하면 AA, Aa, aa를 각각 0, 0, 1로 변환한 값을 모형에 사용합니다. 이때 추정되는 회귀계수는 하나이며 되며 AA 혹은 Aa의 형질에 대한 영향을 0으로 가정했을 때 aa의 형질에 대한 영향을 추정합니다.
• 유전자형: 3 가지의 유전자형 AA, Aa, aa를 각각 (0,0), (1,0), (0,1)로 변환한 값을 입력합니다. 이는 세 유전자형에 대해 자유도가 2가 되므로 2차원 값으로 변환할 필요가 있기 때문입니다. 이런 변환방법을 처리대비라고 하는데 이때는 추정되는 회귀계수가 2개가 됩니다. 다음 식과 같이 첫 번째 계수는 유전자형 Aa에 대한 영향, 두 번째 계수는 aa에 대한 영향을 추정합니다.
• 경향성: 3 가지의 유전자형 AA, Aa, aa를 minor allele의 갯수 0, 1, 2로 변환한 값을 사용합니다. 추정되는 회귀계수는 1개로 minor allele가 하나 증가함에 따른 영향을 평가하게 됩니다. 즉, aa가 형질에 미치는 영향은 Aa의 2배라고 가정하는 것과 같습니다.

제가 여기서 든 예는 통계모형을 이용한 가장 간단한 예의 하나에 불과합니다. 양적 형질에 대한 관련분석은 질적 형질과 달리 통계모형을 만드는 데 사용하는 배경정보에 따라 분석 결과에 차이가 있을 수 있고 배경정보와 유전자형, 배경정보와 배경정보 사이에 상호작용(interaction)이 있을 수 있기 때문에 결과의 해석이 복잡해 질수 있습니다. 양적 형질에 대한 관련분석은 아직 개선해야 할 문제들이 산적해 있기 때문에 활발한 연구가 이루어지고 있습니다.

통계모형을 이용한 분석이 왜 어려운지는 다음 인용으로 대신하도록 하겠습니다. ^^;

Essentially, all models are wrong, but some are useful.
- George E. P. Box. Empirical Model-Building and Response Surfaces (1987)

참고문헌

1. Balding D.J. (2006), Nature reviews Genetics, 7, 10, 781-791.
2. 鎌谷直之 (2007) 遺伝統計学入門, 岩波書店 (카마타니 나오유키 (2007), 유전통계학 입문, 이와나미서점 )

이글과 관련된 글

• 2012/04/12 -- GWAS로 배우는 유전통계학 – 5 분석결과의 시각화
• 2012/04/05 -- GWAS로 배우는 유전통계학 – 3.1 코호트 연구와 실험-대조군 연구
• 2012/03/31 -- GWAS로 배우는 유전통계학 – 2.1 가계모순과 하디-베인베르크 평형 법칙
• 2012/03/30 -- GWAS로 배우는 유전통계학 – 1. 시작하며
• 2012/04/18 -- GWAS로 배우는 유전통계학 – 6. 맺음말
• 2012/04/11 -- GWAS로 배우는 유전통계학 – 4 다중비교 문제
• 2012/04/03 -- GWAS로 배우는 유전통계학 – 2.3 게놈 데이터의 품질평가

3. Genome-wide association study

관련분석은 유전적 변이와 형질과의 관련성을 검출하는 것이 목적입니다. 이때 관측된 SNP좌위가 형질의 표현형(phenotype)에 직접적인 영향을 미친다는 것을 검출할 수 있다면 가장 바람직스러운 결과일 것 입니다(direct association). 그러나 실제로는 관련성을 시사하고 있다고 한다 해도 관측된 SNP좌위가 표현형과 직접 관련이 있다고는 보장할 수 없습니다. 진짜 원인이 되는 유전자 좌와 연쇄불평형(linkage disequilibrium; LD) 상태에 있는 유전자 좌도 표현형과 간접적인 관련이 있을 때가 많기 때문입니다(indirect association).

direct & indirect association from Kruglyak (2008), Nat. Rev. Genet

관련분석의 대상이 되는 형질에는 앞서 언급한 바와 같이 질적 형질과 양적 형질이 있습니다. 따라서 형질의 형태에 따라 분석 방법이 달라집니다. 또한, 수집된 데이터의 연구디자인에 따라 관련성의 지표 및 결과의 해석방법이 달라짐에도 주의해야 합니다. 이번 포스팅에서는 코호트 연구(cohort study)와 실험-대조군 연구(case-control study)에 대해 알아보고 그 후에 질적 형질에 대한 관련분석법, 양적 형질에 대한 관련분석법에 대해 설명하도록 하겠습니다.

3.1 코호트 연구와 실험-대조군 연구

코호트 연구는 연구의 대상이 되는 집단을 일정 기간에 걸쳐 추적조사를 하는 연구법을 말하여 어떤 인자를 가지고 있는 개체와 가지고 있지 않은 개체가 미래에 어떤 표현형이 되는가에 대해 연구하는 방법을 말합니다. 반면 실험-대조군 연구는 표현형에 따라 실험군과 대조군으로 분류하고 각 군에 대해 특정 인자를 포함하고 있는가를 분석하는 방법입니다. 즉, 모든 사건(event)이 이미 일어난 과거의 일을 분석하게 됩니다. 이 때문에 코호트 연구는 연구의 방향이 전향적(prospective)이고, 실험-대조군 연구는 후향적(retrospective)으로 진행됩니다.

관련성의 척도로써 코호트 연구는 상대위험도(relative risk; RR)를 실험-대조군 연구는 오즈비(odds ratio; OR)를 사용합니다. 관측된 데이터에 대해 관련분석을 할 때 일반적으로 다음과 같은 분할표를 이용합니다.

<표 1> 관련분석에서의 분할표

개체의 표현형을 D(disease)와 N(non-disease)라 한다면 Type 1의 개체가 질환에 걸릴 확률과 Type 2의 개체가 질환에 걸릴 확률의 비로 정의되는 상대위험도는
\[
RR=\frac{\frac{a}{a+c}}{\frac{b}{b+d}} = \frac{a(b+d)}{b(a+c)}
\]
로 계산할 수 있습니다.

한편 어떤 사건이 일어나지 않은 확률에 대한 사건이 일어난 확률의 비율로 오즈(odds)를 정의한다면 오즈비는 대조군의 오즈에 대한 실험군의 오즈 비율로 정의됩니다.
\[
OR = \frac{a/b}{c/d} = \frac{ad}{bc}
\]

상대위험도가 좀 더 알기 쉬운 개념이기는 하지만 실험-대조군 연구에서는 실험군 혹은 대조군의 표본 수를 연구자가 결정하게 되므로 상대위험도를 구할 수 없습니다.

코호트 연구는 원인이 되는 개체간 유전자 다형성의 차이가 처음부터 고정되고 결과가 되는 표현형을 관측하게 되므로 자연의 인과관계와 일치하게 됩니다. 또한, 추적관찰을 하게 되므로 사건의 발생순서를 알 수 있다는 점, 측정의 바이어스가 작다는 점, 복수의 결과인자를 동시에 관찰할 수 있다는 점, 표현형이 발현하는 비율로 정의되는 발병율(침투율)을 추정할 수 있다는 장점이 있습니다. 그러나 실험-대조군 실험에 비해 비용과 시간이 걸린다는 점, 발병율(침투율)이 낮은 표현형의 연구에는 표본크기가 크지 않으면 통계분석을 하기 어렵다는 점등의 문제가 있습니다.

참고문헌

1. Balding D.J. (2006), Nature reviews Genetics, 7, 10, 781-791.
2. 이재원, 박미라, 유한나 (2005) 생명과학연구를 위한 통계적 방법. 자유아카데미
3. 鎌谷直之 (2007) 遺伝統計学入門, 岩波書店 (카마타니 나오유키 (2007), 유전통계학 입문, 이와나미서점 )

이글과 관련된 글

• 2012/04/12 -- GWAS로 배우는 유전통계학 – 5 분석결과의 시각화
• 2012/04/07 -- GWAS로 배우는 유전통계학 – 3.2 질적, 양적형질에 대한 관련분석
• 2012/03/31 -- GWAS로 배우는 유전통계학 – 2.1 가계모순과 하디-베인베르크 평형 법칙
• 2012/03/30 -- GWAS로 배우는 유전통계학 – 1. 시작하며
• 2012/04/18 -- GWAS로 배우는 유전통계학 – 6. 맺음말
• 2012/04/11 -- GWAS로 배우는 유전통계학 – 4 다중비교 문제
• 2012/04/03 -- GWAS로 배우는 유전통계학 – 2.3 게놈 데이터의 품질평가