します つまり 操 作 をなぞる 場 合 は 赤 字 部 分 をコピー&ペーストします 入 力 コマンドの 各 行 の# 以 降 右 はコメントで なくても 操 作 に 影 響 しません library(qtl) # qtl ライブラリーを 起 動 する testmap <- read.cross(

Transcription

1 R/qtl による 連 鎖 地 図 の 作 成 清 水 顕 史 R の 利 用 R の qtl ライブラリーは 連 鎖 地 図 も 作 成 できます 1 そのための 遺 伝 子 型 および 表 現 型 デ ータは 次 のフォーマットで 準 備 します データは 表 現 型 (QTL 解 析 用 ) 遺 伝 子 型 の 順 に 上 から 並 べます 1 列 目 は 表 現 型 名 またはマーカー 名 を 記 入 します 表 現 型 データは 最 低 1 つ 記 入 しておく 必 要 があるのでデータがない 場 合 も 何 かデタラメな 値 を 入 力 します 2 列 目 はマーカーの 座 乗 する 染 色 体 番 号 または 連 鎖 群 番 号 ( 任 意 の 数 字 )を 記 入 します( 不 明 の 場 合 も 何 か 記 入 し ておく) 3 列 目 は 必 要 な 場 合 に cm 単 位 で 連 鎖 地 図 位 置 を 記 入 します( 連 鎖 地 図 情 報 が 記 入 されていない 場 合 マーカーは 連 鎖 地 図 順 に 並 べられていると 見 なしてデータが 読 み 込 まれます マーカーは 物 理 地 図 位 置 の 順 に 入 力 します) 表 現 型 データの 2 3 列 目 は 空 白 にしておきます 遺 伝 子 型 データは 母 親 ホモ 接 合 体 を A 父 親 ホモ 接 合 体 を B ヘテロ 接 合 体 を H とし 優 性 マーカーの 場 合 は A と C(B or H)か B と D(A or H)と 表 記 します ( 注 意 1) 表 現 型 データや 遺 伝 子 型 データは 測 定 の 際 にコンタミやエラーを 起 こさないよ うに 注 意 を 払 うのは 当 然 ですが 記 録 したものをパソコンに 入 力 する 際 にもタイプミスが 生 じます 入 力 したデータは 印 刷 し 必 ず 見 直 しましょう( 複 数 人 で 確 認 し 合 うとよい) ( 注 意 2) 解 析 (マッピング)は 一 つの 個 体 に 由 来 する 分 離 集 団 ごとに 計 算 します たとえ 同 じ 交 配 組 合 せの F2 集 団 でも 由 来 する F1 個 体 が 異 なる 場 合 は 別 の 集 団 として 解 析 します ( 減 数 分 裂 は 環 境 要 因 の 影 響 を 受 けることがあり 分 離 比 が 歪 むかもしれない) 2-1)データ 入 力 例 題 用 データ(test.csv)を 例 えば C ドライブの Rdata フォルダに 置 きます R 起 動 後 に[フ ァイル]メニューの[ディレクトリの 変 更 ]からデータのあるフォルダ(ここでは C: Rdata) を 指 定 します 以 下 では R への 入 力 コマンドは 赤 字 で 記 し 出 力 結 果 は 青 字 で 示 すことに 1 研 究 室 のパソコンに 導 入 済 みです R を 自 力 でインストールする 場 合 は 参 考 にしましょう その 後 qtl ライブライ-もインストールします R の 使 い 方 は 研 究 室 に 参 考 図 書 があります 1

2 します つまり 操 作 をなぞる 場 合 は 赤 字 部 分 をコピー&ペーストします 入 力 コマンドの 各 行 の# 以 降 右 はコメントで なくても 操 作 に 影 響 しません library(qtl) # qtl ライブラリーを 起 動 する testmap <- read.cross("csvr", file="test.csv", estimate.map=false) summary(testmap) # 交 配 タイプ 個 体 数 マーカー 数 表 現 型 数 を 確 認 する #データ 読 み 込 み データが 正 常 に 読 み 込 めているかどうか 交 配 タイプや 個 体 数 マーカー 数 表 現 型 数 を 確 認 しましょう 染 色 体 数 や 染 色 体 毎 のマーカー 数 もみておきましょう 2-2) 欠 損 データの 確 認 plot.missing(testmap) # 欠 損 データをプロットする par(mfrow=c(1,2), las=1) plot(ntyped(testmap), ylab="no. typed markers", main="no. genotypes by individual") plot(ntyped(testmap, "mar"), ylab="no. typed individuals", main="no. genotypes by marker") 2

3 欠 損 データの 多 い 2 マーカーや 個 体 は とりあえず 解 析 から 外 します 必 要 があれば ジ ェノタイピングし 直 すか 近 傍 マーカーを 作 り 直 します データの 削 除 は 以 下 の 通 りです 3 例 えばデータ 数 が 50 以 上 の 個 体 のみを 採 用 する 場 合 testmap <- subset(testmap, ind=(ntyped(testmap)>50)) 例 えば 個 体 数 が 70 未 満 のマーカーを 解 析 から 除 外 する 場 合 nt.bymar <- ntyped(testmap, "mar") todrop <- names(nt.bymar[nt.bymar < 70]) testmap <- drop.markers(testmap, todrop) 2-3) 分 離 比 の 歪 んだマーカーを 検 出 する 各 マーカー 遺 伝 子 座 の 分 離 比 が 期 待 値 から 大 きくずれる 場 合 4 遺 伝 子 型 の 判 別 方 法 に 問 題 5 がある 可 能 性 が 高 い 期 待 比 からのズレは geno.table 関 数 を 使 って 調 べることができ ボン フェローニ 補 正 6 した 5% 水 準 以 下 のマーカーは 以 下 のように 検 出 できる gt <- geno.table(testmap) gt[gt$p.value < 0.05/totmar(testmap),] この 例 の 場 合 染 色 体 3 のこの 領 域 の 隣 接 マーカー(RM7 や RM251)の p 値 も 小 さいた め 歪 みの 原 因 はマーカー 遺 伝 子 型 の 問 題 ではなく 遺 伝 的 な 要 因 であるとみなせます こ の 種 のマーカーの 歪 みは 連 鎖 地 図 の 作 成 や 形 質 マッピングの 妨 げにはならないでしょう 隣 接 マーカーから 突 出 した 歪 みを 示 すマーカーは 削 除 するのが 無 難 です もう 一 度 ジェ ノタイピングしなおすか 代 替 マーカーを 探 すのがよいですが さしあたっては 解 析 から 外 します 例 えばボンフェローニ 補 正 前 の p 値 が 以 下 の 歪 みの 激 しいマーカーを 削 除 する 場 合 には 以 下 のようにします todrop <- rownames(gt[gt$p.value < 1e-10,]) testmap <- drop.markers(testmap, todrop) 2 明 確 な 基 準 はありませんが データの 20% 以 上 が 欠 損 する 場 合 は 怪 しいと 考 えましょう 3 Excel で 入 力 ファイルを 編 集 してもよいでしょう 4 F2 分 離 集 団 では 共 優 性 マーカーの 遺 伝 子 型 の 分 離 が A:H:B=1:2:1 に 優 性 マーカーは A:C または B:D=1:3 に 分 離 することが 期 待 でされる 5 マーカーが 単 一 遺 伝 子 座 由 来 でない 場 合 や ホモ ヘテロの 判 別 が 難 しく 誤 った 遺 伝 子 型 を 含 んでいる 場 合 などが 考 えられる 6 ボンフェローニ 補 正 は χ 2 適 合 度 検 定 による p 値 に 調 査 した 全 マーカー 数 を 掛 ける 3

4 2-4) 連 鎖 地 図 の 検 定 formlinkagegroups() 関 数 で 連 鎖 群 を 推 定 できます 指 定 したmax.rf( 組 換 価 の 上 限 )と min.lod(lodスコアの 下 限 )を 満 たす 任 意 の2つのマーカーは 同 じ 連 鎖 群 に 割 り 当 てます 7 lg <- formlinkagegroups(testmap, max.rf=0.45, min.lod=5) table(lg[,2]) 連 鎖 群 別 に 割 り 当 てられたマーカー 数 が 降 順 に 示 されます 例 データの 場 合 連 鎖 群 数 はちょうどイネの12の 染 色 体 数 に 等 しくなりました その 連 鎖 群 順 に 8 マーカー 対 の 組 換 価 とLODスコアをプロットすると 連 鎖 群 ごとのマーカー 順 序 の 精 度 を 図 示 できます testmap <- formlinkagegroups(testmap, max.rf=0.45, min.lod=5, reorgmarkers=true) plot.rf(testmap, alternate.chrid=true) 図 の 左 上 半 分 は 組 換 価 のプロットで 右 下 半 分 は LOD スコアのプロットです マーカー 間 で 連 鎖 群 作 成 の 閾 値 よりも 組 換 価 が 小 さい もしくは LOD スコアが 大 きい 場 合 を 赤 色 で 示 し( 連 鎖 する 可 能 性 が 非 常 に 高 い) 青 は 逆 に 組 換 価 が 大 きいまたは LOD スコアが 小 さいことを 示 します( 独 立 である 可 能 性 が 非 常 に 高 い) 上 例 では 対 角 線 上 に 赤 色 がプロ 7 置 換 系 統 を 用 いたマッピングの 場 合 は より 厳 しい 条 件 でも 連 鎖 地 図 は 作 成 できるでしょう 8 組 換 価 と LOD スコアのプロットは 染 色 体 番 号 順 ではないことに 注 意 4

5 ットされ 連 鎖 群 の 異 なるマーカー 間 は 全 て 青 色 となり 連 鎖 群 の 割 り 当 ておよびマーカ ーの 並 びに 矛 盾 がないきれいな 連 鎖 地 図 であることを 示 しています ちなみに 例 データ で マーカー 順 序 をわざとシャッフルした 間 違 った 連 鎖 地 図 では 以 下 のようなプロットに なります 連 鎖 地 図 の 検 証 これまでの 操 作 で 欠 損 データや 分 離 比 の 歪 みがみられたら 入 力 データからおかしなマー カーを 削 除 します 必 要 に 応 じてジェノタイピング データのエラーチェック 再 ジェノ タイピング マーカーの 差 し 替 えを 行 います 連 鎖 群 内 のマーカー 順 序 の 検 証 ( 修 正 )データの 再 読 み 込 み testmap <- read.cross("csvr", file="test.csv", estimate.map=false) #データ 再 読 込 例 データでは 染 色 体 毎 に 連 鎖 群 が 定 まりマーカー 順 序 も 物 理 地 図 通 りになりました(p6) そこで 以 下 では 染 色 体 ごとにマーカー 順 序 を 検 定 します そのさい 染 色 体 番 号 ( 予 め 入 力 用 データの 2 列 目 に 記 入 しておく test.csv を 参 照 )を chr オプションで 指 定 します 連 鎖 群 内 のマーカー 数 が 30 以 下 程 度 の 場 合 のマーカー 順 序 の 探 索 ordermarkers() 関 数 では 任 意 のマーカー 対 を 取 り 上 げ 無 作 為 に 選 んだマーカーを 最 良 の 場 所 に1つずつ 追 加 します 5

6 testmap <- ordermarkers(testmap, chr=12) pull.map(testmap, chr=12) ripple() 関 数 を 使 うと 指 定 したwindowsサイズの 数 だけマーカー 順 序 の 並 べ 替 えを 検 討 し てくれます ripple() 関 数 の"likelihood"メソッドを 用 いると 尤 度 が 一 番 大 きいマーカー 順 序 を 最 善 の 並 びとして 検 索 してくれます その 際 error.probを 指 定 することで 遺 伝 子 型 デー タに 含 まれるエラー 確 率 ( 推 定 値 )を 考 慮 できます windowsサイズは 大 きい 程 ( 最 大 値 は 連 鎖 群 内 のマーカー 数 )より 緻 密 な 検 索 ができます ただし windowsサイズを7 以 上 9 にす ると ( 用 いるPCの 性 能 によっては) 計 算 に 数 時 間 以 上 要 する 場 合 があります 組 換 価 と LODスコアのプロットで 連 鎖 の 順 序 に 不 具 合 が 見 られない 場 合 は windowsサイズを4~7 程 度 で 並 びを 探 索 するとよいでしょう rip12lik <- ripple(testmap, chr=12, window=7, method="likelihood", error.prob=0.005) summary(rip12lik) 例 の 場 合 2 と 1 を 入 れ 替 えた 並 びが 2 番 目 に LOD が 高 く 染 色 体 長 の 推 定 値 が 少 し 長 くなることが 分 かります マーカー 順 序 は LOD が 高 いほどよく 染 色 体 長 はなるべく 短 い ものを 選 びます 選 ばれるマーカー 順 序 は error.prob で 指 定 した 遺 伝 子 型 データのエラー 率 に 影 響 を 受 けま す エラー 率 を 変 えた 場 合 に 最 適 マーカー 順 序 が 変 わるか 否 かを 検 討 する 場 合 compareorder 関 数 を 用 いて 新 たに 指 定 したマーカー 順 序 (new)と 初 期 入 力 順 序 (orig)を 比 較 します 染 色 体 12 で マーカー2 と 1 を 入 れ 替 えた 順 序 との 比 較 は エラー 率 0.01 で 以 下 のように 計 算 します compareorder(testmap, chr=12, c(2,1,3:10), error.prob=0.01) compareorder(testmap, chr=12, c(2,1,3:10), error.prob=0.001) 9 オーバーナイトで 計 算 させる 場 合 があるとしても せいぜい windows サイズは 10 程 度 までにする 6

7 compareorder(testmap, chr=12, c(2,1,3:10), error.prob=0) 結 果 からマーカー 順 序 は 初 期 順 序 の LOD が 最 大 で 染 色 体 の 長 さも 短 いことがわかりま す もし 初 期 順 序 よりも LOD スコアの 高 いマーカーがみつかった 場 合 は switch.order() 関 数 でマーカー 順 序 を 入 れ 替 えます この 関 数 は est.map()と 同 じ 引 数 を 使 用 して マーカー 順 序 を 交 換 したマーカー 順 序 で 染 色 体 地 図 を 推 定 します testmap <- switch.order(testmap, chr=x, c(2,1,3:10), error.prob=0.005) 例 データの 場 合 全 ての 染 色 体 で 最 良 になるマーカー 順 序 を 初 期 値 としています 連 鎖 地 図 の 検 証 ギャップの 検 証 summary.map() 関 数 を 用 いると 各 染 色 体 の 最 大 ギャップがわかります summary.map(testmap) プロットplot.map() 関 数 で 引 数 にshow.marker.names=TRUEを 用 いると マーカー 名 付 きの 連 鎖 地 図 でギャップ( 連 鎖 地 図 のマーカー 密 度 が 低 い 箇 所 )を 表 示 できます plot.map(testmap, show.marker.names=true) 7

8 例 データの 一 番 大 きなギャップは 34.0cM でした 94 個 体 程 度 以 上 のデータ 数 がある 場 合 50cM 未 満 のマーカー インターバルであればマーカー 間 の 連 鎖 は 有 意 といえるでしょ う ただし 形 質 マッピングで 40cM 以 上 のマーカー 間 に QTL が 検 出 された 場 合 は マ ーカー 間 (ギャップの 中 央 付 近 がよい)に 座 乗 する 新 たなマーカーを 追 加 しマーカー 密 度 を 高 めた 方 がよいでしょう 1 マーカー 除 去 ジェノタイピング エラーを 起 こしやすいマーカー( 例 えばヘテロが 区 別 しにくい)を 検 出 する 方 法 に droponemarker() 関 数 を 用 いた 1 マーカー 除 去 があります 一 度 に 一 つのマー カーだけ 除 去 した 場 合 の LOD スコアと 染 色 体 長 の 減 少 程 度 を 観 察 することで 怪 しいマー カーを 検 出 することができます dropone <- droponemarker(testmap, error.prob=0.005) par(mfrow=c(2,1)) plot(dropone, lod=1, ylim=c(-100,0)) plot(dropone, lod=2, ylab="change in chromosome length") 8

9 結 果 は 横 軸 に 染 色 体 1~12 が 連 鎖 地 図 順 に 並 び 地 図 上 のその 位 置 のマーカーを 取 り 除 い たときの LOD スコアおよび 染 色 体 長 の 減 少 程 度 がプロットされています 取 り 除 くと LOD スコアが 上 昇 するマーカーや 取 り 除 くと 染 色 体 長 が 大 きく 減 少 するマーカーが 怪 しいマー カーです 各 染 色 体 で 一 番 染 色 体 長 を 減 少 させるマーカーは 次 のように 検 出 できます summary(dropone, lod.column=2) 染 色 体 の 両 端 に 位 置 するマーカーは 取 り 除 くと 染 色 体 長 が 減 少 するのは 当 たり 前 なので これは 割 り 引 いて 考 える 必 要 があります 例 の 場 合 染 色 体 6 と 11 で 検 出 されたマーカー のみ 染 色 体 の 内 側 に 位 置 するマーカーになりました 試 しにこれらのマーカーを 取 り 除 いて みる 場 合 は 以 下 のようにします badmar <- c(rownames(summary(dropone, lod.column=2))[6],rownames(summary(dropone, lod.column=2))[11]) testmap <- drop.markers(testmap, badmar) 9

10 newmap <- est.map(testmap, error.prob=0.005) testmap <- replace.map(testmap, newmap) summary.map(testmap) # 連 鎖 地 図 の 再 推 定 取 り 除 くと 遺 伝 距 離 が 10cM 以 上 減 少 する 染 色 体 の 内 側 のマーカーについては 遺 伝 子 型 の 再 チェックか 必 要 に 応 じて 代 わりになるマーカーの 遺 伝 子 型 を 調 べます 問 題 個 体 の 探 索 countxo() 関 数 を 用 いると 各 個 体 の 染 色 体 の 乗 換 え 数 ( 推 定 値 )を 表 示 できます plot(countxo(testmap), ylab="number of crossovers") 集 団 の 中 で 突 出 して 染 色 体 乗 換 え 数 が 多 い 個 体 は 怪 しいので 削 除 します 例 えば 染 色 体 乗 換 え 数 が 80 以 上 のマーカーを 除 いて 連 鎖 地 図 をつくるには testmap <- subset(testmap, ind=(countxo(testmap) < 80)) とします 削 除 後 のマーカー 順 序 は 連 鎖 群 毎 の ripple 関 数 を 用 いた 検 証 (p7)を 再 度 行 い ます ジェノタイピング エラー 率 を 見 積 もる 実 際 に 調 査 したマーカーの 遺 伝 子 型 データにはどうしてもジェノタイピングのエラーが 含 まれます 以 下 のコードはエラー 率 0.001~ の 範 囲 で 連 鎖 地 図 の LOD を 計 算 し プ ロットしてくれます LOD が 最 大 になるエラー 率 を 採 用 します loglik <- err <- c(0.001, , 0.005, , 0.01, , 0.015, , 0.02) for(i in seq(along=err)) { cat(i, "of", length(err), "\n") tempmap <- est.map(testmap, error.prob=err[i]) loglik[i] <- sum(sapply(tempmap, attr, "loglik")) } 10

11 lod <- (loglik - max(loglik))/log(10) plot(err, lod, xlab="genotyping error rate", xlim=c(0,0.02), ylab=expression(paste(log[10], " likelihood"))) 例 データの 場 合 はエラー 率 で LOD スコアが 最 大 ( 付 近 )になりました データに よっては エラー 率 の 検 索 範 囲 を 広 げる 必 要 があるでしょう ジェノタイピング エラーの 可 能 性 がある 個 体 の 探 索 calc.errorlod() 関 数 を 用 いると ジェノタイピング エラーの 可 能 性 がある 個 体 を 検 出 する ことができます top.errorlod() 関 数 で エラーを 持 つ 可 能 性 が 高 いマーカーおよび 個 体 番 号 をリスト 表 示 できます LODスコアが4~5 以 上 の 場 合 エラーの 可 能 性 が 高 いと 考 えるのが 一 般 的 です 例 えばLOD5 以 上 のリスト 表 示 は 以 下 のようにします(ここでは マーカー 順 序 をデタラメに 並 べたerror.csvを 使 用 します) testmap <- read.cross("csvr", file="error.csv", estimate.map=false) #データ 読 み 込 み testmap <- calc.errorlod(testmap, error.prob=0.015) print(toperr <- top.errorlod(testmap, cutoff=4)) plot.geno() 関 数 を 用 いると 検 出 された 怪 しい 遺 伝 子 型 を 表 示 できます 染 色 体 1についてみ てみると 以 下 のようになります 結 果 の 図 では 左 縦 軸 の 個 体 番 号 (5,9,11) 中 に 見 出 され た 怪 しい 遺 伝 子 型 が 赤 色 四 角 で 示 されています 横 軸 は 染 色 体 の 連 鎖 地 図 位 置 (cm)です plot.geno(testmap, chr=1, ind=toperr$id[toperr$chr==1], cutoff=5, include.xo=false) 11

12 ジェノタイピング エラーの 検 出 では 隣 接 マーカーが 二 重 組 換 え 型 のものを 検 出 していま す(Lincoln and Lander 1992 Genomics 14: ) 検 出 されたエラー 遺 伝 子 型 については 再 調 査 を 行 うべきです error.csv の 場 合 全 ジェノタイピング 数 (94 個 体 138 マーカー- 51 欠 損 データ)のうち LOD 5 の 閾 値 で 検 出 されるジェノタイピング エラーは 4 つのみ (0.03%)であり 連 鎖 地 図 と 形 質 マッピングには 影 響 を 及 ぼしません 10 検 出 されたジェノタイピング エラーを 削 除 する 場 合 は 次 のコードが 使 えます testmap.clean <- testmap for(i in 1:nrow(toperr)) { chr <- toperr$chr[i] id <- toperr$id[i] mar <- toperr$marker[i] testmap.clean$geno[[chr]]$data[testmap$pheno$id==id, mar] <- NA } # 分 離 比 の 歪 みを 再 チェック 最 後 に もう 一 度 分 離 比 の 歪 みをチェックします testmap <- read.cross("csvr", file="test.csv", estimate.map=false) #データ 再 読 み 込 み gt <- geno.table(testmap, scanone.output=true) par(mfrow=c(2,1)) plot(gt, ylab=expression(paste(-log[10], "P-value"))) plot(gt, lod=3:5, ylab="genotype frequency") abline(h=c(0.25, 0.5), lty=2, col="gray") 10 マーカー 二 重 組 換 個 体 は その 領 域 の 准 同 質 遺 伝 子 系 統 作 成 に 非 常 に 有 益 で 遺 伝 子 単 離 に 使 えます 12

13 単 一 マーカーによる 分 離 比 の 歪 みはみられません( 上 段 ) 遺 伝 子 型 の 比 ( 下 段 の 黒 青 赤 線 は それぞれ AA, AB, BB に 対 応 )は 染 色 体 7 に1つ 染 色 体 8 に1つ 染 色 体 9 に 2 つある 優 性 マーカーで 歪 んで 見 えますが 分 離 比 自 体 は 歪 みません 最 終 的 な 地 図 を 描 きます plot.map(testmap, show.marker.names=true) 13

14 最 終 的 に 得 られた 連 鎖 地 図 のマーカー 順 序 や 染 色 体 群 および 地 図 位 置 は 入 力 データに 反 映 させましょう(その 際 マーカー 名 と 遺 伝 子 型 データが 入 れ 違 いを 起 こさないように 注 意 すること) 形 質 マッピングに 必 要 な 連 鎖 地 図 情 報 とは 染 色 体 全 域 に 渡 り 20cM 程 度 の 遺 伝 距 離 で 連 鎖 する DNA マーカーの 遺 伝 子 型 情 報 です 当 然 ですが 染 色 体 断 片 置 換 系 統 と 反 復 親 との 交 雑 後 代 の 場 合 は 置 換 領 域 のみの 遺 伝 子 型 情 報 が 必 要 になります 14