plinkでペアワイズ連鎖不平衡解析｜非モデル生物の場合
はじめに
RAD-seqやMIG-seqなどで取得できるSNP情報は、もろもろの解析に使用する前に連鎖不平衡を調べて除去してやる必要がある。
ここではplinkを使用してSNP間の相関係数をペアワイズで計算し、閾値（こちらで設定できる）よりも相関係数が大きな組み合わせを連鎖不平衡のあるSNPとして検出する。
また、今回は非モデル生物でRAD-seqやMIG-seqをおこなってde novoでアライメントした場合を想定し、各locusの染色体番号や位置情報を使わずに連鎖解析をおこなう。
インプットファイル
.pedファイルと.mapファイルが必要
.pedファイルは.strなどからPGDspiderで変換可能
ただし、同じ座位に3つ以上アリルがあるとはじかれるので、事前に取り除いておく必要がある
.mapファイルはエクセルで自分で作れる
左から順に
染色体番号|座位の番号|遺伝的距離|ポジション
を半角スペース区切りで記述していくが、de novoの場合、座位の番号以外は0と記入すれば良い
0 loc1 0 0
0 loc2 0 0
0 loc3 0 0
0 loc4 0 0といった感じ
.pedファイルに記述されている座位と.mapファイルに記述されている座位は一致する必要があるかも
実行
.pedファイルと.mapファイルを同じファイル名で（拡張子は異なる）同じディレクトリに格納する
例えば、
mkdir ~/Desktop/plinkで、できたフォルダに
○○.ped
○○.mapを移すといった感じ
その後、
plink --file ○○ --r2 --ld-window 1000000 --ld-window-r2 0※〇〇：拡張子をつけてはいけない
と打ち込めば、同ディレクトリに
plink.ldが作成される
コマンドのオプションを解説すると、
-- file ○○ インプットファイルの指定
-- r2 相関係数の計算
-- ld-window 1000000 何個隣のSNPまで相関係数を計算するか …※
--ld-window-r2 0 相関係数が指定した数値以上のものを出力する（全て出力するときは0）
※入力したデータが実際の配列順で並んでいる場合、あまり遠くの座位と相関係数を出しても意味が無い。そのため、デフォルトでは9個隣の座位までしか相関係数を計算しないようになっている。de novoの場合、実際の配列順に並んでいるわけではないので、全てのSNPを総当りで調べる必要がある。よって、ここには手持ちのSNP数以上の数字を入れる必要がある。SNP数以上であれば数字はなんでもよい。
アウトプットファイル
plink.ldをmiでもエクセルでもなんでもいいので開いてみると、
CHR_A BP_A SNP_A CHR_B BP_B SNP_B R2
0 0 loc1 0 0 loc2 0.0152985
0 0 loc1 0 0 loc3 0.0515134
0 0 loc1 0 0 loc4 1こんな感じになっているはず
見たままだが、loc1とloc2の相関係数が0.0152985であるという感じに読む
ちなみに、たとえば
--ld-window-r2 0
の0の数字を0.2にすれば、相関係数が0.2以上のものだけをリストアップしてくれる
