ka=kr=40,サンプル集団の人数を変えたとき。母集団からランダムに200人取り出して検定を行っている。母集団120人の時は200人選べないので100人選んだ。ランダムに何人取り出すかを変えないとあまり意味がない、と反省。特に、母集団が少ないと、あたえるアレ…

白いほど検出力が高い。それぞれ100回ずつ。今週やること： ・アレル頻度とリスクのプロットなど、グラフを分かりやすく使うための勉強 −軸を対数に ・ユークリッド距離とマンハッタン距離の違いをトポロジーから ・尤もらしいモデルで人数取り出す人数をか…

ka=60,kr=0 横軸ユークリッド距離、縦軸マンハッタン距離。 正の相関。マンハッタン距離の方が数値がかなり大きい。正の相関なので、どちらを使っても結果は似ると思われる。ka↑、あるいはkr↑で、距離はどちらも短くなる。kr↑で短くなる理由はなんだろう。リ…

・ユークリッド距離とマンハッタン距離について、違いをプロット・ka[0-90],kr[0-180]位でまわす。・アレル頻度、リスク、遺伝率、マーカーの数をそれなりのものに固定したとき、サンプルの数による検出力の挙動・人数固定してマーカーの数を変えてみるマー…

アレル頻度とリスクの係数を10刻みで調べてみた 横軸がアレル頻度の係数、奥行きがリスクの係数、高さが値。 左手前の四角の隅っこが(0,0,0)。平均値 分散 どちらも低いといい。分散が低いのが検出力が高いといえるかもしれない。常にある値らへんが出るわけ…

フェノタイプとジェノタイプに関係があるかを、遺伝率で調整してもいいのだけど、リスクとアレル頻度に絞って考えることにする。たとえばリスク=0では絶対にフェノタイプとジェノタイプに関係性がなくなるので。リスクを、risk=exp(-10000*af)にしてみると、…

par(ask=TRUE) N<-200 #バリアントの数 i HN<-100000 #人数 j ir<-al<-matrix(0,HN,N) #各ローカスにおけるリスク、持っているアレル wei<-wei2<-P<-sr<-rep(0,HN) #重み、フェノタイプ、リスクの総和 risk<-af<-iaf<-rep(0,N) #真のアレル頻度 Np<-50 #MST…

ユークリッド距離かマンハッタン距離かで悩むところなのだけれど、どちらも検出力がものすごく変わるわけではなさそうなので、マンハッタン距離でやる。２つやると処理にものすごく時間がかかるので。k1=60,80あたりで、遺伝率がどのくらいならどれくらい検…

マンハッタン距離で。 下のことに気づいたので途中やめしたけれど、リスク＝±１のときが一番いいなぁ。。。パーミテーションテストでランダムに作った値をソートしたとき、下から何番目の枠に入るか、で調べているけど、分布関数を予測してｎ分位点で切るほ…

横軸はrisk=±exp(-k1*af)のk1で 1,2,3,4・・・にk1=0,10,20,30・・・が対応 縦軸はパーミテーションテストでランダムに作った値をソートしたとき、下から何番目の枠に入るか、の分布軸の書き方が必要。。。パーミテーションの回数か、取るサンプル数か、試行…

risk=±exp(-k1*af),af=rexp(N,k2) ・まず適当なk2についk1を動かしてみる。 ・k1、k2をいろいろ動かしてみて、関係性があるかどうかk1について、k1作ってなかったので１だ、と思っていたけど、パラメータわりあてて40にしていた。10刻みで0〜80まで調べてみ…

遺伝率ごとの検出力。0.68まではポジコンとしていけるのかもしれない。ポジコンとしていけるということはつまり、ある条件下なら、遺伝率が0.7以上のフェノタイプについて常に解析ができるということ。しかし、ポジコンの条件を決める際は、遺伝率は0.9で固…

・ポジコンはリスク=±exp(-アレル頻度)にしてみる。 ・意味としては、頻度が低いほどリスクが高い。 ・今、アレル頻度によってリスクを決めているので、アレル頻度の分布自体は結果にあまりかかわらない。アレル頻度が高いものが増えると、リスク０のものが…

0/1× MST[0,0,7,6,8](400,rexp(N,60),abs(rnorm)) 0/1× MST[18,2,0,0,0](400,rexp(N,60),exp(-40*af)*sign(rnorm(N, mean=0, sd=1))) 0/1× MST[0,0,2,5,13] afはアレル頻度(400,rexp(N,60),rnorm(N,0,exp(-40*af))) 0/1× MST[0,0,10,6,4]

0/1×（有意な結果は得られなかった。） MST[0,0,10,4,6] 上のほうに数値が出なかったら、5段階でかく。(400,rexp(N,60),+1) 0/1× MST[20,0,0,0,0] 常に高い。 ここまで極端だと何かおかしい気がする。考えてみる。アレル頻度をrexp(N,60)としているのは、こ…

０／１のほうは、 パラメータ(バリアントの数,アレル頻度,リスク,閾値) 結果は[有意に上、間、有意に下]遺伝率は常に0.5で。それに環境要因（適当な正規分布）をランダムに加える。という感じで結果を書いておくことにする。 閾値は上位10〜2%について0.5%刻…

P2<-sr*sdfrac+rnorm(Npt,sd=sqrt(var(P))*(1-sdfrac)) というところがあって、多分これは遺伝率を考えるためのものだと思うので、加えてみる。 N<-200 #バリアントの数 i HN<-100000 #人数 j #アレル頻度 al<-matrix(0,HN,N) #持っているアレル aa<-60.0 #…

最小全域木の考え方のほうのプログラムを、いままで作ってきたほうと一緒に動くようにする。 主要ハプロタイプのところはリスク0で考えているので、とりあえず除いて考える。ＭＳＴのほうはぼんやりとしかイメージをもてていないが、なるべくそのままがんば…

リスクをアレルと独立にしてみたら、ケースとコントロールの重みの分散の差は、１０００だったり０だったりー１０００だったり、ばらばらというかランダムな数が出た。時々、分散の差が（適当に選んだ100人と100人の差）、すべてほとんど0になるときがある。…

パーミュテーションテストには、全部調べる方法と、ランダムに選んで分布を推定する方法があり、シミュレーションではランダムパーミュテーションを使おうと思う。集中検定だと、アレル頻度の逆数を重みとして、その分散の、ケースとコントロールとの差が、…

N<-100 #バリアントの数 i HN<-10000 #人数 j #アレル頻度 al<-matrix(0,HN,N) #持っているアレル aa<-5.0 #アレル頻度パラメータ af<-rexp(N,aa) #アレル頻度を与える #リスク sr<-numeric(HN) #リスクの合計点 ir<-matrix(0,HN,N) #個人の各アレルのリスク…

N<-40 p<-0.001 #変異率 T<-80 #時間 Npt<-80 #試行回数 K<-10 #最初の変異アレルの本数 jun<-array(0,c(2,N,2)) for(i in 1:2) for (j in 1:N) jun[i,j,]<-c(i,j) sme<-matrix(0,T,Npt) #時間t、nパターン目におけるアレル頻度 for (we in 1:Npt){ Gpa<-mat…

N<-100 p<-0.0001 #変異率 T<-20 #時間 Npt<-100 #試行回数 sme<-matrix(0,T,Npt) #時間t、nパターン目におけるアレル頻度 for (we in 1:Npt){ Gpa<-matrix(0,2,N) #親ゲノム Gpa[1,1]<-1 Gch<-matrix(0,2,N) #子ゲノム Gse<-matrix(0,2,N) #子がどの親から…

中立の場合について、時間がたつとどんな感じになるかシュミレーションして、あるtで切って頻度の分布を作って、それが人数によってどう代わるか調べて、expで近似してみる。アレル頻度が何なのかわからなくなってきた。。。 というかなんとなくSNPについて…

http://www.ncbi.nlm.nih.gov/pubmed/11404818 アレル頻度をもっともらしいものにしたいと思ったので調べていると、この論文の１２７ＰのWright's formulaが使えるかもしれない。ただ、２つのアレルについてしか考えていないので、３つ以上のときどうしよう…

先生からいただいたソースを読み解いていった。USBを忘れたのでメモ書き。 Ns<-100 #アレルの数 Fs<-runif(Ns)*0.05 #頻度 Rs<-runif(Ns) #リスク #Fs<-sort(Fs) #Rs<-sort(Rs,decreasing=TRUE)^1*sample(c(-1,1),Ns,replace=TRUE) plot(Fs,Rs) Npt<-100000 …

N<-100 #バリアントの数 i HN<-100 #人数 j ir<-matrix(numeric(HN*N),HN,N) #個人の各アレルのリスク（２本分足した値） al<-matrix(numeric(HN*N),HN,N) #持っているアレル（複数のtypeなら０，１，２，４を足す感じで） alt<-1 #アレルの種類 aa<-5.0 #ア…

尤度比検定、カイ自乗検定、正確確率検定などがあるが、サンプルが少ないときは正確確率検定がいいらしい。 p値が大きめに出ることを保守的という。正確確率検定は保守的。 遺伝形式が優性、劣性の場合は２×２表でいいけれど、相加的な場合は２×３なので、フ…

N<-1000 #バリアントの数 i HN<-1000 #人数 j aa<-5.0 #アレル頻度パラメータ al<-matrix(numeric(HN*N),HN,N) #アレルの有無０〜２ アレルがいろいろあるときはそれぞれ1,2,4を足すかもしれない sr<-numeric(HN) #リスクの合計点 af<-exp(-aa*runif(N,0,)) …

cc<-quantile(sr,probs=seq(0.95,0.95,1)) #上位５％点 byo<-numeric(HN) for (j in 1:HN) if(sr[j]>cc) byo[j]<-1 tt<-matrix(numeric(N*6),6,N) #各アレルについて、caseMM,,,contmm for (j in 1:HN){ for (i in 1:N){ if (byo[j]==1){ if (al[j,i]==2) tt…