fine

3dplot

al<-matrix(0,36,2) al[,1]<-c(rep(2,9),rep(1,12),rep(0,15)) al[,2]<-c(2,2,1,1,1,0,0,0,0,2,2,2,1,1,1,1,0,0,0,0,0,2,2,2,2,1,1,1,1,1,0,0,0,0,0,0) sr<-al[,1]*1+al[,2]*1 P<-sr+rnorm(36,0,0.1) rgl.clear() rgl.bg(col="white") x<-cbind(P,al) d<-dis…

■

遺伝率0.5で、リスク０とリスクが±1くらい（正規分布でちょっとランダムに揺らす）にしたらこうなった。リスク0ではパフォーマンスがよくなくてほっとした。次は遺伝率を変えてみよう

■

リスクのあたえかたを正規分布にしてみた。右は平均0で、左は平均が1のものと-1のものを組み合わせたものになる。RVにリスクがある時は検出されることがわかったけれど、RVにリスクがない時検出されないかが分からない。今遺伝率１なのでリスク0にするとフェ…

■

200人、バリアント80個、ただし200人について観測されるRV(mafフェノタイプの値を、検定の時に、-0.5~0.5くらいの間になるようにした。最大値と最小値の絶対値の和で割る。

2011-10-05

■

フェノタイプの単位とジェノタイプの単位は違って、フェノタイプ=2とジェノタイプ=2の意味は違うのだけど、MSTで処理するときに、同じものとしてあつかわれているので、ちょっと微妙になっているかもしれない。困った。リスクの与え方を工夫しないといけない…

2011-10-05

■

source("MST.R") par(ask=F) Np<-100 HN<-500 df<-2000 N<-280 #バリアントの数 i ka<-1000 coa<-c(0,1000,2000,3000,4000,5000) lcoa<-length(coa) syu<-matrix(200,Np,lcoa) syu2<-matrix(0,Np,lcoa) colnames(syu)<-coa nrvcv<-matrix(200,Np,lcoa) afb<-…

2011-10-04

■

1000人、N=200、観測されるRV150個くらい、krを変えたとき 500人、n=300、観測されるRV150個くらい、krを変えたとき観測されるRVの数が揺れていて微妙なので、何人はこのRVを持ってる、と決めてしまって、誰が持っているかはランダムに決めるのがいいと思う…

2011-09-30

■

十分な人数、十分なバリアントでリスクの与え方を変えたいのだが、人数を2倍にすると4倍くらいの時間がかかる気がする。どのくらいにするか難しいところ。1000人、RV150個で、大体72秒。あるkrについて100回試行して2時間。

2011-09-30

距離のはかり方について

距離のはかり方をなんとなくmanhattanからeuclidにしてみたら上手くいきそうな感じ。。。。manhattanでCVをちょっと増やすとかより、こっちの違いについて考えるのが重要な感じがする。euclidで色々な条件で試してみて、それで上手くいくならそれでいいし、…

2011-09-30

■

サンプルを200人とってきて、こういったRVセットについて考える、としたとき、どんなRVを選んで調べて見てもいいのだから、アレル頻度のところは実際は融通がききそうな気がしてきた。といっても、みつかるRVのアレル頻度は指数分布に従うだろうから、10-50%…

2011-09-29

■

サンプル数がある程度少ないとき、世界全体でのアレル頻度分布と、少ないサンプル数でのアレル頻度分布は変わるかもしれない。実際のプログラミングについて考えると、最初100個RVを設定していても、多分100種類も観察されていなくて、重要なのはどちらかと…

2011-09-28

■

http://www.hsph.harvard.edu/~xlin/SKAT-paper-and-supplement.pdf のアレル頻度モデルのもととなった実際のデータ。RV93個。アレル頻度かなり低い。これだけ低いと、RVをもってて1,2個だから、上手く結果がでないかもしれない。大体ka=300くらい。kaとkr…

■

レアバリアントが高リスクのとき、フェノタイプと関係性があるとでて、レアバリアントが低リスクのとき、フェノタイプと関係性がないという結果が出したい。この高リスク、低リスクは、CVに比べて、というものなので、CVについてもう少し考えて組み込んで、…

メモ

CVはリスク0で、アレル頻度20~50%。４個。 RVは高リスクで、アレル頻度0~10%。40個。CVはリスク0だから、フェノタイプにまったく寄与せず、遺伝要因に占める割合は0。RVは、持っているとフェノタイプに大きく寄与する。というモデルにおいて、MST検定を行う…

■

CV4個で固定してRVの数を変えてみるとこうなった。こうなる理由がつかめない。モデルとしては、4個、と固定するより、RVの数によって変わる方が説明しやすそう。たとえば、グラフの変化が小さい順にka=20,40,60,80 上のモデルのような確率分布に従う時、検…

■

#バリアントの数が60のときの、CNを変える Date<-"09261" source("MST.R") Np<-100 #MST何回か HN<-50000 df<-200 cnnu<-c(20,30,40,50,60,70) Rnu<-length(cnnu) syu<-matrix(0,Np,Rnu) #boxplot用 colnames(syu)<-cnnu for (nnp in 1:Rnu){ CN<-cnnu[nnp] …

リスク０のCVを入れてみる　RV40,60個のとき

RV40個 CVの数を変えてみる。 0のときを加えて、CVの数が小さいときを調べてみる。CVについて、今アレル頻度を一様分布で与えているが、数が少ないから適当に設定すべきだ。RV60個

N=60のとき、アレル頻度の刻み幅を変える

レアバリアントの数を60以下にしたとき検出力がよくなかったのはなぜかについて。段の高さの違いは、afの与え方のパラメタka（af=rexp(N,ka)）で変わるので、段の幅のとりかたによるのではないかと考えた。RVの数が少ないとき、アレル頻度の設定の仕方を、1/…

euclid,manhattan

オレンジと緑のところで、木のでき方が違うことになる。木のでき方がちがっても、距離は常にマンハッタン距離の方が長かった。差がプラスかマイナスかしか調べていないが。

■

p3を円上で動かしたときの、p1,p2,p3の全域木の長さを調べる。 [:] eはユークリッド距離、mはマンハタン距離 123だったら、p1-p2-p3と結んだときの長さ最小全域木の長さあまり特徴がない。今度は楕円でやってみる。これを円でやると、p3-p1-p2と結んだと…

バリアントの数を変えたとき２

100個で調べたときに少し平均と分散があがっているのは、試行回数が少ないからだと考えられる。120個あったら十分で、80個でも使えるかもしれない。

■

バリアントの数が少ないとき、リスク0のコモンバリアントを考えるとMSTがうまくいくのではないか、という予想。