生命情報科学人材養成コンソーシアム バイオインフォマティクス実習コース 遺伝医学 2010年10月29日（金曜日） 京都大学大学院医学研究科 附属ゲノム医学センター 統計遺伝学分野 教授 山田　亮 © 2010 山田　亮

本日の構成 「遺伝統計学」 遺伝統計学概観 遺伝学基礎 遺伝統計学でのデータの扱い 分割表検定 連鎖解析 たくさんの検定

目標 データ解析をする立場で、研究対象を捉えることにイメージを持つ データを眺めるスタイルを理解する 既存の手法や用語に対して、自分の理解の役に立つのかどうかという視点で客観的(批判的)に対する

受講の動機は？ 出身・背景は？ 内容は、受講者の希望に応じて変更しながら進みます せっかくなので知りたいことを知って帰りましょう 個別研究？/研究全般のアシスト？/ ”遺伝統計学”研究?/シリーズ受講なのでついでに? 出身・背景は？ 生物・遺伝系？/統計・計算機系？

今日の教材の元ネタ 遺伝統計学の基礎 Rによる遺伝因子解析・遺伝子機能解析 ISBN 978-4-274-06822-5 Download Sources of R http://www.genome.med.kyoto-u.ac.jp/wiki_tokyo/index.php/StatGenetOhm

※　使用するRコードの入手とその使用法 -"http://www.genome.med.kyoto-u.ac.jp/wiki_tokyo/index.php/StatGenetOhm"から、 -Rsrc.zipファイルをダウンロードして解凍 -Rの作業ディレクトリを解凍してできたRsrcディレクトリに指定する -"Rコードのソースを読み込み"で、"StatGenetDemo.R"を読み込むと、デモ開始 -必要なパッケージのインストールが行われ、Rsrcディレクトリ内のRソースを順に実行する -ソースファイルを順送りするときと、表示図を確認するためのプロンプトが現れるので、"Enter"キー等でデモを進める

遺伝統計学概観

似ている 似ていない Alike, not alike

10

離散的 ４カテゴリ：Ａ、Ｔ、Ｇ、Ｃ 長さ　30億塩基対

acagccagag ggacgagcta gcccgacgat ggcccagggg acattgatcc gtgtgacccc agagcagccc acccatgccg tgtgtgtgct gggcaccttg actcagcttg acatctgcag ctctgcccct gaggactgca cgtccttcag catcaacgcc tccccagggg tggtcgtgga tattgcccac ggccctccag ccaagaagaa atccacaggt tcctccacat ggcccctgga ccctggggta gaggtgaccc tgacgatgaa agtggccagt ggtagcacag gcgaccagaa ggttcagatt tcatactacg gacccaagac tccaccagtc aaagctctac tctacctcac cggggtggaa atctccttgt gcgcagacat cacccgcacc ggcaaagtga agccaaccag agctgtgaaa gatcagagga cctggacctg gggcccttgt ggacagggtg ccatcctgct ggtgaactgt gacagagaca atctcgaatc ttctgccatg gactgcgagg atgatgaagt gcttgacagc gaagacctgc aggacatgtc gctgatgacc ctgagcacga agacccccaa ggacttcttc acaaaccata cactggtgct ccacgtggcc aggtctgaga tggacaaagt gagggtgttt caggccacac ggggcaaact gtcctccaag tgcagcgtag tcttgggtcc caagtggccc tctcactacc tgatggtccc cggtggaaag cacaacatgg acttctacgt ggaggccctc gctttcccgg acaccgactt cccggggctc attaccctca ccatctccct gctggacacg tccaacctgg agctccccga ggctgtggtg ttccaagaca gcgtggtctt ccgcgtggcg ccctggatca tgacccccaa cacccagccc ccgcaggagg tgtacgcgtg cagtattttt gaaaatgagg acttcctgaa gtcagtgact actctggcca tgaaagccaa gtgcaagctg accatctgcc ctgaggagga gaacatggat gaccagtgga tgcaggatga aatggagatc ggctacatcc aagccccaca caaaacgctg cccgtggtct tcgactctcc aaggaacaga ggcctgaagg agtttcccat caaacgcgtg atgggtccag attttggcta tgtaactcga gggccccaaa cagggggtat cagtggactg gactcctttg ggaacctgga agtgagcccc ccagtcacag tcaggggcaa ggaatacccg ctgggcagga ttctcttcgg ggacagctgt tatcccagca atgacagccg gcagatgcac caggccctgc aggacttcct cagtgcccag caggtgcagg cccctgtgaa gctctattct gactggctgt ccgtgggcca cgtggacgag ttcctgagct ttgtgccagc acccgacagg aagggcttcc ggctgctcct ggccagcccc aggtcctgct acaaactgtt ccaggagcag cagaatgagg gccacgggga ggccctgctg ttcgaaggga tcaagaaaaa aaaacagcag aaaataaaga acattctgtc aaacaagaca ttgagagaac ataattcatt tgtggagaga tgcatcgact ggaaccgcga gctgctgaag cgggagctgg gcctggccga gagtgacatc attgacatcc cgcagctctt caagctcaaa gagttctcta aggcggaagc ttttttcccc aacatggtga acatgctggt gctagggaag cacctgggca tccccaagcc cttcgggccc gtcatcaacg gccgctgctg cctggaggag aaggtgtgtt ccctgctgga gccactgggc ctccagtgca ccttcatcaa cgacttcttc acctaccaca tcaggcatgg ggaggtgcac tgcggcacca acgtgcgcag aaagcccttc tccttcaagt ggtggaacat ggtgccctga gcccatcttc cctggcgtcc tctccctcct ggccagatgt cgctgggtcc tctgcagtgt ggcaagcaag agctcttgtg aatattgtgg ctccctgggg gcggccagcc ctcccagcag tggcttgctt tcttctcctg tgatgtccca gtttcccact ctgaagatcc caacatggtc ctagcactgc acactcagtt ctgctctaag aagctgcaat aaagtttttt taagtcactt tgtac

染色体伝達図 系統樹

遺伝統計学の位置づけ

The Lady Tasting Tea: How Statistics Revolutionized Science in the Twentieth Century [ペーパーバック] David Salsburg  ISBN-13: 978-0805071344 統計学を拓いた異才たち―経験則から科学へ進展した一世紀 [単行本] デイヴィッド サルツブルグ (著), David S. Salsburg (原著), 竹内 惠行 (翻訳), 熊谷 悦生 (翻訳)  ISBN-13: 978-4532351946

遺伝学基礎

血縁関係　→　似る 似る程度にばらつき 似る特徴もあれば似ない特徴もある 似るか似ないかには「理由」がある

遺伝的に伝わることと多様であること Heredity and Variation 何が伝わり、何が多様か？ “もの”と”こと” 遺伝子型と表現型

ジェノタイプとフェノタイプ

フェノタイプ分類のこつ

ジェノタイプ

染色体

遺伝的多様性の３要素 変異 交叉・組換え 遺伝的浮動

遺伝子多型・変異

遺伝子多型 構造分類

遺伝子多型 サイズ分類

遺伝子座 アレル ハプロタイプ ディプロタイプ

交叉

組み換え 起源の一致(ＩＢＳ)は、２本の染色体ができるまでの交叉の回数の偶奇で決まる

交叉回数 ポアッソン分布 R2-6.R

#R2-6.R # 可能箇所すべてで交叉がおきるかどうかを試す方法 RecombSim<-function(L=10000,r=0.001,Niter=1000){ # Lは配列長,rは箇所あたりの交叉確率,Niterはシミュレーション試行回数 # 行数Niter、列数L-1箇所の行列にする m<-matrix(rbinom((L-1)*Niter,1,r),nrow=Niter) apply(m,1,sum) } # ポアッソン分布を使う方法 RecombPois<-function(L=10000,r=0.001,Niter=1000){ rpois(Niter,(L-1)*r) # rpois() 関数については付録A.5 確率分布関数・疑似乱数列の発生を参照

#R2-6.R(続き) L<-10000;r<-0.0001;Niter<-1000 NumSim<-RecombSim(L=L,r=r,Niter=Niter) NumPois<-RecombPois(L=L,r=r,Niter=Niter) ylim<-c(0,max(NumSim,NumPois)) plot(ppoints(Niter,a=0),sort(NumSim),ylim=ylim,col=gray(6/8)) par(new=T) plot(ppoints(Niter,a=0),sort(NumPois),type="l",ylim=ylim)

交叉間距離 指数分布 R2-7.R

#R2-7.R Niter<-1000 # シミュレーション回数 L<-1000000 #染色体の長さ r<-0.0001 #塩基間あたりの交叉確率 # 交叉箇所数をポアッソン分布からの乱数で指定し、交叉箇所をsample関数で指定する crosses<-sort(sample(1:(L-1),rpois(1,(L-1)*r),replace=FALSE)) # 交叉間距離のベクトルを作る A<-c(0,crosses) # 染色体の始点と交叉箇所のベクトル B<-c(crosses,L) # 交叉箇所と染色体の終点のベクトル C<-B-A #交叉間距離のベクトル # 平均がmean(C)の指数分布からの乱数をlength(C)の数だけ発生させてプロット rexps<-rexp(length(C),1/mean(C)) # rexp() 関数については付録A.5 確率分布関数、疑似乱数列の発生を参照 # 交叉間距離をソートしてプロット ylim<-c(0,max(C,rexps)) plot(sort(C),ylim=ylim,cex=0.5,pch=15) #交叉間距離の昇順プロット par(new=T) plot(sort(rexps),col="red",ylim=ylim,type="l") # 指数分布からの乱数の昇順プロット

遺伝的浮動

#”a”と”b”の集団を作る Na<-1;Nb<-11;k<-4 Ns<-Na+Nb A<-c(rep("a",Na),rep("b",Nb)) A #k倍する B<-rep(A,k) #Ns個抜き取る sample(B,Ns) Niter<-1000 #Niter回、繰り返して Numa<-rep(0,Niter) for(i in 1:Niter){ S<-sample(B,Ns) Numa[i]<-length(which(S=="a")) } #"a"の抜き取られ数の分布をみる hist(Numa) #"a"が２個抜き取られた確率は？ length(which(Numa==2))/Niter #n!を使って計算してみる gamma(12+1)/(gamma(2+1)*gamma(10+1))*gamma(36+1)/(gamma(2+1)*gamma(34+1))*gamma(4+1)*gamma(44+1)/gamma(48+1) 多型性を失う ドリフトアウト

#"a"の抜き取り数は、0,1,2,3,4のいずれか #その確率を計算して足し合わせてみる P0<-gamma(12+1)/(gamma(0+1)*gamma(12+1))*gamma(36+1)/(gamma(4+1)*gamma(32+1))*gamma(4+1)*gamma(44+1)/gamma(48+1) P1<-gamma(12+1)/(gamma(1+1)*gamma(11+1))*gamma(36+1)/(gamma(3+1)*gamma(33+1))*gamma(4+1)*gamma(44+1)/gamma(48+1) P2<-gamma(12+1)/(gamma(2+1)*gamma(10+1))*gamma(36+1)/(gamma(2+1)*gamma(34+1))*gamma(4+1)*gamma(44+1)/gamma(48+1) P3<-gamma(12+1)/(gamma(3+1)*gamma(9+1))*gamma(36+1)/(gamma(1+1)*gamma(35+1))*gamma(4+1)*gamma(44+1)/gamma(48+1) P4<-gamma(12+1)/(gamma(4+1)*gamma(8+1))*gamma(36+1)/(gamma(0+1)*gamma(36+1))*gamma(4+1)*gamma(44+1)/gamma(48+1) P0+P1+P2+P3+P4

酔歩 Random walk R16-sup1.R

#R16-sup1.R nstep<-100 rwalk<-matrix(0,nstep,2) rtheta<-rnorm(nstep-1) stepx<-cos(rtheta) stepy<-sin(rtheta) for(i in 1:(nstep-1)){ rwalk[i+1,1]<-rwalk[i,1]+stepx[i] rwalk[i+1,2]<-rwalk[i,2]+stepy[i] } plot(rwalk,type="l")

差が大きいときには急速に 差が小さいときにはゆっくりと 差と「時間微分」が比例 拡散方程式 拡散 R6-4.R

差が大きいときには急速に 差が小さいときにはゆっくりと 差と「時間微分」が比例 拡散方程式

# R6-4.R # pa,pb:２集団の人口,d:単位時間あたりの移住人数,t:世代 pa<-9000 pb<-1000 d<-100 t<-0:100 fa<-(pa+pb*exp(-2*d*(pa+pb)/(pa*pb)*t))/(pa+pb) fb<-(pa*(1-exp(-2*d*(pa+pb)/(pa*pb)*t)))/(pa+pb) plot(t,fa,ylim=c(0,1),type="l",xlab="t",ylab="frequency") par(new=T) plot(t,fb,ylim=c(0,1),type="l",xlab="t",ylab="frequency")

アレルの行く末 どちらかのアレルのみになる アレル頻度 R2-8.R 時間 片方のアレルに固定 低頻度ながら、こちら側にも固定

遺伝統計学でのデータの扱い

カテゴリカルデータ 全体を覆う 重なりがない

カテゴリと順序 量的形質とカテゴリ カテゴリと次元

カテゴリカルデータの表現 分割表

自由度

自由度　記述するための変数の数 ２ｘ２ table

いくつの数値で表を説明するか？

テーブルを説明するのに 必要な変数の数 (自由度)

自由度　１

有意性を判断する 説明をするのに、変数を増やすか増やさないか、それが問題

珍しさで順序をつける δが大きさで順序をつける δが大きいほど「珍しい」

独立を仮定したとき 順列(パーミュテーション) 項目１ 項目２ ○１ ■１ ○２ ■２ ○３ □３ ●４ □４ ●５ □５ ●６ ■６ 項目１ 項目２ ○１ □５ ○２ □４ ○３ ■１ ●４ □３ ●５ ■２ ●６ ■６ 項目１と項目２のラベル付けは独立だと仮定すれば 並べ替えない 並べ替える 順列(パーミュテーション) N個の並べ替え：N!

独立を仮定したとき 順列(パーミュテーション) 項目１ 項目２ ○１ ■１ ○２ ■２ ○３ □３ ●４ □４ ●５ □５ ●６ ■６ 項目１ 項目２ ○１ □５ ○２ □４ ○３ ■１ ●４ □３ ●５ ■２ ●６ ■６ 項目１と項目２のラベル付けは独立だと仮定すれば 順列(パーミュテーション) N個の並べ替え：N! 並べ替えごとに、分割表ができる δの分布がとれる 珍しいδとありきたりなδがわかる 観測データのδのありきたりな程度を0-1の値で表わす

完全パーミュテーション と モンテカルロパーミュテーション 全順列(N!)は膨大 しらみつぶしにするか 1部を使うか(乱数を使う：モンテカルロ)

パーミュテーションの網羅と サンプリング library(gtools) n<-3 permutations(n,n) permutations(n,n,repeats=TRUE) sample(1:n,n)

パーミュテーション(順列)の計算 パーミュテーション(順列)の足し合わせ N!は膨大 モンテカルロをするにしても膨大 モンテカルロの回数とパーミュテーションテストのｐ値 便法はないか？

便法 ないこともある

便法 ないこともある 分割表にはある

分割表の正確生起確率 珍しさの計算 G1,G2をn..人に割り付ける場合はn..!/(n1.! n2.!)通り A,aをn..人に割り付ける場合はn..!/(n.1! n.2!) “G1-A”,”G1-a”,”G2-A”,”G2-a”をn..人に割り付ける場合はn..!/(n11! n12!n21!n22!)

正確確率と完全パーミュテーションの結果は同じ ラベルの枚数を考慮するかしないか しない パーミュテーション する 正確確率計算

δと生起確率 生起確率 δ R13-2.R

正確確率の計算 正確確率の足し合わせ 計算自体が面倒臭い 足し合わせるためには、「膨大な」観測可能テーブルのすべてについて計算が必要 自由度が上がると非現実的 便法はないか？

正確確率のプロットを近似できる関数があれば・・・ 生起確率 δ R13-2.R

２つの便法 確率 ピアソンの独立性検定のカイ自乗値 尤度 尤度比検定

正確生起確率とはなんだったか 「ある仮説(帰無仮説)」の下で 「あるデータ」が観察される 確率

正確生起確率はなんだったか 「ある仮説(帰無仮説)」の下で 「あるデータ」 「別のデータ」 … 「すべてのデータ」が観察される 確率 →　確率分布に照らす(カイ自乗検定)

正確生起確率はなんだったか 「ある仮説(帰無仮説)」 「別のある仮説(対立仮説)」 … 「すべての仮説(対立仮説と帰無仮説)」の下で 「あるデータ」が観察される 確率 →　仮説に基づく確率(尤度)の比較

確率と尤度 Probability and likelihood 仮説を固定、観察を動かす 確率：G1,G2に差がないときにn11=x (x=0,1,2,…)という観察をする確率 尤度：G1ではAの割合がp1でG2ではAの割合がp2であるという仮定のもとでn11=n11という観察をする確率(p1=0～1,p2=0～1) 仮説を動かす、観察を固定

observation hypothesis

ピアソンの独立性検定のカイ自乗値 独立仮説でのデータ観測確率を知る

尤度による方法 確率pの事象がn回続けて起きる確率 確率pの事象がn回続けて起きて、そのあとm回続けて起きない確率 p^n x (1-p)^m n+m回のうちn回起きてm回起きない確率 C x p^n (1-p)^m C: (n+m)からnを取り出す場合の数

２ｘ２表の尤度 A a All G1 N1A N1a N1 G2 N2A N2a N2 Gall NA Na T A a All G1 P(1A) P(1a) 1 G2 P(2A) P(2a) Gall P(A) P(a)

ピアソンのカイ自乗検定の自由度 尤度比検定の自由度 表の中身を変更するときに、自由に動かせるセルの数 尤度比検定 比較する２つの仮説で使う変数の差 ２ｘ２表の場合はなぜ、自由度が１？

検定３種 Three types of tests 正確確率検定 Exact tests パーミュテーションテスト Permutation-based テーブルの正確生起確率による Exact Probability based on table ピアソンの独立性検定 Pearson's independence test 尤度比検定 Likelihood ratio test

だいたい同じ 少し違う Similar each other but a bit different

漸近近似 Nが大きいとほとんど一緒

カテゴリカルデータの検定 検定の手法 Ｐ値に持ち込む形式の違いでわける 正確確率に基づく検定 期待値表からの距離に基づく検定 尤度に基づく検定

カテゴリと順序 量的形質とカテゴリ カテゴリと次元

ディプロタイプというカテゴリ AaとaA AA, Aa(=aA),aaの３つ Aを父からaを母から Aを母からaを父から アレルの並び順を区別するかしないか、と、エピジェネティック効果 AA, Aa(=aA),aaの３つ

ディプロタイプというカテゴリ ＳＮＰの場合 ３タイプが対等 ３タイプに順序(アレルについて評価) アレルの本数の力 相乗的 相加的 それ以外 優性　劣性　ヘテロが特別

多アレルの場合 順序あり リピート数など（タンデムリピート、ＣＮＶ) 順序なし ハプロタイプ

順序あり・多アレルの ディプロタイプの場合 1 2 (0,0) (0,1) (0,2) (1,0) (1,1) (1,2) (2,0) (2,1) (2,2) 順序をつける 足し算・・・ (0,0)=0, (0,1)=(1,0)=1,…,(2,2)=2

順序なし(かもしれない)多アレル ハプロタイプというカテゴリ ＳＮＰのアレルの組み合わせ 組み合わせは多次元 多次元のまま扱う 順序を入れられる？ ディプロタイプにすると？ Na+Na(Na-1)/2 カテゴリ

普通の表 サンプル Ｇ Ｐ S1 G1 P1 S2 P2 S3 G2 S4 S5 G3 … Sn P1 P2 G1 N11 N12 G2

普通の表 2x2, 2x3, …, NxM 表 どれも同じ 自由度は(N-1)(M-1) 自由度が大きくなると正確確率検定は大変すぎる

正確確率検定の重さ #行列を作る m<-matrix(c(10,20,30,40,50,60),nrow=2,nrow=3) m #作る行列の値の納め方を変える m<-matrix(c(10,20,30,40,50,60),nrow=2,nrow=3,byrow=TRUE) fisher.test(m) #大きな表、大きな自由度 nrow=4;ncol=4 m<-matrix(round(runif(nrow*ncol)*10,0),nrow=nrow,ncol=ncol)

普通の表の検定 検定の手法 Ｐ値に持ち込む形式の違いでわける 正確検定・ピアソン・尤度比検定

普通の表の検定 検定の手法 Ｐ値に持ち込む形式の違いでわける 正確検定・ピアソン・尤度比検定 表の形とカテゴリの順序関係の違いでわける

普通の表の検定 検定の手法 確率・尤度 生起確率を元に、それをまねる確率分布(カイ自乗のような)を使う方法 傾向性検定、ＫＷ、ＪＴ 比較するべき仮説(帰無と対立)を変数でモデル化する方法(尤度比検定のような) 線形回帰・ロジスティック・それ以外のもっといろいろな回帰(一般化線形回帰とか)

普通の表 サンプル Ｇ Ｐ S1 G1 P1 S2 P2 S3 G2 S4 S5 G3 … Sn P1 P2 G1 N11 N12 G2

普通ではない表 ハーディ・ワインバーグ平衡検定を例に

ハーディワインバーグ平衡(HWE)とは AA Aa aa ALL N(AA) N(Aa) N(aa) N(ALL)

ランダムメイティング

HWEの期待値表 P(AA)=P(A) ^2 P(Aa)=2P(A)p(a) P(aa)=P(a)^2 AA Aa aa ALL N(AA) N(ALL) E(AA)=P(A)^2 N(ALL) E(Aa)=2P(A)P(a) N(ALL) E(aa)=P(a)^2 N(ALL) P(AA)=P(A) ^2 P(Aa)=2P(A)p(a) P(aa)=P(a)^2

ＨＷＥと近交係数 N(AA)/N(ALL) = E(AA) + f p(A)p(a) f:近交係数 (HWEの指数) f=0 : HWE f=1 : 全員ホモ接合体

HWEの観察表からカイ二乗検定 アレル本数 M(A)=(2N(AA)+N(Aa))/2 M(a)=(N(Aa)+2(aa))/2 ALL N(AA) N(Aa) N(aa) N(ALL) アレル本数 M(A)=(2N(AA)+N(Aa))/2 M(a)=(N(Aa)+2(aa))/2 アレル頻度 (の推定値) P(A)=M(A)/(2 N(ALL)) P(a)=M(a)/(2 N(ALL))

HWEの期待値からのずれ カイ自乗値 自由度:1 (自由変数:f) N(AA)/N(ALL) = E(AA)+ f p(A)p(a) E(AA)=P(A)^2 N(ALL) E(Aa)=2P(A)P(a) N(ALL) E(aa)=P(a)^2 N(ALL) カイ自乗値 自由度:1 (自由変数:f) N(AA)/N(ALL) = E(AA)+ f p(A)p(a) N(Aa) /N(ALL)= E(Aa) – 2f p(A)p(a) N(aa)/N(ALL) = E(aa) + f p(A)p(A)

ＨＷＥ検定カイ自乗値と近交係数 N(AA)/N(ALL) = E(AA) + f p(A)p(a) f:近交係数 (HWEの指数) f=0 : HWE f=1 : 全員ホモ接合体 f^2 x N(ALL)=カイ自乗値 Chi^2 = N f^2

普通ではない表　HWEの表 AA Aa aa ALL N(AA) N(Aa) N(aa) N(ALL)

ＨＷＥ正確確率検定 正確率検定はできる データの生起確率を計算できて すべての場合を網羅できればよい

ＨＷＥ正確確率検定 R4-1.R 正確率検定はできる データの生起確率を計算できて

ＨＷＥ正確確率検定 正確率検定はできる すべての場合を網羅できればよい ヘテロの人数は奇数か偶数かのどちらか R4-1.R

ハーディ・ワインバーグ平衡 ↓ 連鎖平衡

アレル関連 連鎖不平衡 連鎖平衡 ハプロタイプ

アレル関連 数を確認 アレル数が２の多型 k個 ハプロタイプ種類数 2^k ディプロタイプ種類数 3^k

アレルの関係が独立とは 独立は直交 内積がゼロ～相関がゼロ 独立は掛け算

ハプロイド：２ローカスの独立 ハプロイド 同一染色体上 同一配偶子(精子・卵子)内

LE状態 2SNP haplotype ((1)AB,(2)Ab,(3)aB,(4)ab) H1 = pA pB H2 = pA pb

２ローカスが独立でなくなるとき 同一分子上にあるから、自由でない 交差してよければ・・・ 遠ければ遠いほど 完全な自由：無限遠：異なる染色体

２ローカスが独立でなくなるとき 「無限」に遠くても不自由 メイティングが不自由(HWD,集団の構造化)

LE状態からのずれを表す 2SNP haplotype ((1)AB,(2)Ab,(3)aB,(4)ab) H1 = pA pB +d d=r SQRT(pA pa pB pb) r^2 : LD 指標

HWD と LD H1 = pA pB +d MM : p(M)^2 +f p(M)p(m) H2 = pA pb -d Mm : 2p(M)p(m) – fp(M)p(m) mm : p(m)^2 +f p(M)p(m) MM : p(M)p(M) + d Mm : p(M)p(m) – d mM : p(m)p(M) – d mm : p(m)p(m) + d d= f sqrt(p(M)p(m)p(M)p(m)) H1 = pA pB + d H2 = pA pb – d H3 = pa pB – d H4 = pa pb + d d=r sqrt(pA pa pB pb)

HWE LE

指数と統計量 MM : p(M)p(M) +d Mm : p(M)p(m) – d mM : p(m)p(M) – d d= f sqrt(p(M)p(m)p(M)p(m)) H1 = pA pB +d H2 = pA pb -d H3 = pa pB –d H4 = pa pb + d d=r sqrt(pA pa pB pb) r^2 : LD index N : No. samples Chi^2 = N r^2 r : 相関係数 Chi^2 = N f^2

LDは染色体上に広がる

染色体上に広がるLE/LDの具合 K個のマーカー ペアワイズな比較 K^2 な分布 K(K-1)/2 ペア Mx とMy (x , yは異なる) との関係(K(K-1)/2 ペア) Mx とMx (自身) との関係 分散・共分散行列　　　　　　　　　 (K^2の要素)

多数のＳＮＰ ペアについてのr^2から絵を描く R5-2.R R5-3.R m<-matrix(rbinom(120,1,0.5),20,6) heatmap(m) cormatrix<-cor(m);rsqmatrix<-cormatrix^2 image(1:nrow(rsqmatrix),1:ncol(rsqmatrix),rsqmatrix,col=gray((100:0)/100))

ペアワイズLDプロット ペアワイズ情報を用いる マーカーの位置を変えない パターンから染色体に沿った交叉・組み換えの様子を読み取る

LDプロットができるまで M1 M2 M3 M4 … N1 1 N2 O N3 N4 N5 N6 NxM MxM

NxMのデータの処理 ペアワイズな関係

多数の要素のペアワイズ関係の図示 要素の順番を並びかえれば、それは階層的クラスタリング R5-1.R

階層的クラスタリングは 樹形図の１方法

クラスタリング 階層的クラスタリング 非階層的クラスタリング クラスタ化する方法の違い クラスタ化して残す情報の違い

『階層的』クラスタリング 樹形図化 ペアワイズ比較情報 ５つの要素(葉)が ３つの結合点(節)を介して ７本の線(辺)で結びつけられている 葉への要素の割り当てパターンと辺の長短がクラスタを定めている

『非階層的』クラスタリング 要素の数だけ点がある 点が空間に配置されている 位置の情報が、帰属クラスタを決めている R6-2.R

類別してプロットする 主成分分析 PCA

分散の分解 Decomposition of variance into pieces R7-5.R

適切な軸 Appropriate axes

固有値分解・主成分分析 (PCA) R7-5.R 正規直交基底 どうして「直交」　 分散が基底成分の分散に分解できるから

エッセンス　連鎖解析 木の解析 尤度の解析

家系図 Pedigree

家系の観察 マーカーＸに想定される木 マーカーＹに想定される木 表現型に想定される木

木の解析 尤度の解析 マーカーのジェノタイプの伝達木 個体のフェノタイプの伝達木 異なる伝達木を一致させるためのねじれ 交叉 交叉の起きやすさ(尤度)は距離の関数 面倒くさいのは、木の場合分けが多いから

教科書　第10章 伝達の木 たくさんあるけれど、数えられる マーカーＸの木(Tx)、表現型の木(Tp)、マーカーＹの木(Ty)が異なるときには、間で『組換え』が起きたと考える 距離が遠ければ、組換えが起きやすい 木のパターンによって、Tx,Tp,Tyの遠近関係として「ありそうな位置」が推定できる

マルチプルテスティング たくさんの検定 教科書　第17章

均一分布 N<-1000 X<-runif(N) plot(X) plot(sort(X)) plot(ppoints(N,a=0),sort(X))

均一分布からのk個の乱数の最小値はどれくらい小さいか

N<-10000 k<-5 Xs<-matrix(runif(N*k),N,k) Mins<-apply(Xs,1,min) mean(Mins) 1/(k+1) kの値を振ってやってみる

１番小さいｐ値の期待値 k個の独立なテスト 1/(k+1) i番目に小さいｐ値の期待値?

N<-1000 k<-10 Xs<-matrix(runif(N*k),N,k) Sorted<-matrix(apply(Xs,1,sort),N,k,byrow=TRUE) plot(apply(Sorted,2,mean)) kの値を振ってやってみる

期待値は1/(k+1) 分布は？ N<-1000 k<-100 Xs<-matrix(runif(N*k),N,k) Mins<-apply(Xs,1,min) mean(Mins) plot(density(Mins)) abline(v=mean(Mins))

このグラフをどう読む

plot(sort(Mins)) abline(h=mean(Mins))

i番目に大きいｐ値の期待値はi/(k+1) 分布は？ doubleSorted<-apply(Sorted,2,sort) matplot(doubleSorted,type="l")

独立なp値

Sidak’s family-wise error rate Bonferroni’s method すべての検定を「排他的」に考える もっとも保守的 p_corrected= kp Sidak’s family-wise error rate すべての検定を「独立」と考える p_corrected=1-(1-p)^k k<-100000 plot(1:k,1/k*(1:k),ylim=c(0,1),type="l") par(new=TRUE) plot(1:k,1-(1-1/k)^(1:k),ylim=c(0,1),col="red",type="l")

実際には、すべての検定は排他的でもなく、独立でもない パーミュテーションで、「全通り」を調べる モンテカルロ・パーミュテーションで「全通り」に近い分布を調べる

対立仮説が正しいとき 非心カイ自乗分布 カイ自乗分布に１変数追加

強弱いろいろな対立仮説が成り立つとき たとえば、集団構造化があるときのＧＷＡＳ 「集団構造化」を１追加変数でモデル化すれば、観察データから、その変数の値を推定することができて、その推定値に基づいて補正ができる

たくさんの対立仮説が正しいとき 1変数でモデル化するほど単純ではない 大雑把に「正しい対立仮説の数」が分かっている(と仮定する) False Discovery Rate（FDR)