3

haploNet パッケージを使用してハプロタイプ ネットワーク上でいくつかのプロットを作成する場合、インターネットで入手できるスクリプトを使用して実行しました。しかし、私は何かが間違っていると思います。スクリプトは、woodmouse の例の形式で入手できます。私が使用したコードは次のとおりです。

x <- read.dna(file="Masto.fasta",format="fasta")
h <- haplotype(x)
net <- haploNet(h)
plot(net)

plot(net, size = attr(net, "freq"), fast = TRUE)
plot(net, size = attr(net, "freq"))
plot(net, size=attr(net, "freq"), scale.ratio = 2, cex = 0.8

table(rownames(x))

ind.hap<-with(
    stack(setNames(attr(h, "index"), rownames(h))), 
    table(hap=ind, pop=rownames(x)[values])
)
ind.hap 

plot(net, size=attr(net, "freq"), scale.ratio = 2, cex = 0.8, pie=ind.hap)
legend(50,50, colnames(ind.hap), col=rainbow(ncol(ind.hap)), pch=20)

legend(x=7,y=10,c("Baeti ero","Felege weyni","Golgole naele","Hagare selam","Ruba feleg","Ziway"),c("red","yellow","green","turquoise","blue","magenta"))

ただし、ind.hap をプロットすると、一部の行が適切な場所にないことがわかります。これはここで見ることができます：

      pop
hap    Baetiero ETH022 ETH742 Felegeweyni Golgolenaele Rubafeleg
  I           0      0      1           0            0         0
  II          0      1      0           0            0         0
  III         1      0      0           1            0         1
  IV          2      0      0           0            0         3
  IX          0      0      0           1            0         0
  V           4      0      0           0            2         0
  VI          4      0      0           1            0         4
  VII         2      0      0           1            0         0
  VIII        0      0      0           1            0         1
  X           3      0      0           0            1         0
  XI          0      0      0           0            1         1
  XII         0      0      0           1            0         0
  XIII        0      0      0           0            0         1

行 IX が適切な場所にないことがわかります。これはそれほど問題にはなりませんが、プログラムは行 9 を使用して、VIII のデータである IX の円グラフを作成します。結果は次のとおりです: (評判が 10 を下回っているため、画像を挿入できませんでした...とにかく、ファイル全体を実行することで画像を取得できます)

V から IX までは、本来あるべき状態ではないことがわかります (これらはスワップされた行です)。例: IX にはハプロタイプが 1 つしかありませんが、VIII データを使用して生成された 2 つのハプロタイプ (どちらもチャートの 50% を占めます) の円グラフがあります。行は昇順ではなくアルファベット順にソートされますが、これはパッケージ固有のものであるため、どうすればよいかわかりません。私は R の達人にはほど遠いので、抽象的になりすぎないようにして、代わりにコードを提供してください。

このパッケージをよく知っている人がいる場合は、woodmouse の例では見えなかったので、実際のチャートの後ろにこれらの奇妙な余分な線がある理由も説明してください (数字が付いています)。それも？）

事前にサンクス

4

1 に答える 1

5

私は同じ問題に苦労しましたが、解決策を思いついたと信じています。

問題は、table「人口」ごとのハプロタイプ数を作成するステップで、ハプロタイプがアルファベット順に並べられることです。したがって、たとえば、ハプロタイプ「IX」は「V」の前に来ます。一方、関数haplotype()はハプロタイプを「数値」順にソートします。そして、これがプロット時に矛盾を生み出すものです。

これは、ヘルプで説明されているように、ハプロタイプ オブジェクトを「ラベル」でソートすることで解決できます?haplotype。

サンプル データを使用してwoodmouse例を示します。

# Sample 9 distinct haplotypes
library(pegas)
data(woodmouse)
x <- woodmouse[sample(9, 100, replace = T), ]

簡単にするために、ハプロタイプのカウント テーブルを作成する関数を作成します (この投稿に基づく)。

countHap <- function(hap = h, dna = x){
    with(
        stack(setNames(attr(hap, "index"), rownames(hap))),
        table(hap = ind, pop = attr(dna, "dimnames")[[1]][values])
    )
}

それでは、ハプロタイプをソートせずに結果を見てみましょう。

h <- haplotype(x) # create haplotype object
net <- haploNet(h) # create haploNet object

plot(net, pie = countHap(), size = attr(net, "freq"), legend = T)

ここに画像の説明を入力

それでは、カウント テーブルを見て、これらの結果を確認してみましょう。

countHap(h, x)

      pop
hap    No0906S No0908S No0909S No0910S No0912S No0913S No304 No305 No306
  I          0       0       0       0       0       0     0     8     0
  II         0       0       0       0       0       0     9     0     0
  III        0       0       0       0       0       0     0     0    10
  IV        16       0       0       0       0       0     0     0     0
  IX         0       0       0       0       0       8     0     0     0
  V          0      12       0       0       0       0     0     0     0
  VI         0       0      10       0       0       0     0     0     0
  VII        0       0       0      13       0       0     0     0     0
  VIII       0       0       0       0      14       0     0     0     0

一致しないもの: たとえば、ハプロタイプ「V」は個体「No0908S」に存在するはずですが、代わりに個体「No0913S」として色付けされています (これはハプロタイプ「IX」のラベルである必要があります)。

それでは、ハプロタイプを並べ替えましょう。

h <- haplotype(x)
h <- sort(h, what = "labels") # This is the extra step!!
net <- haploNet(h)

plot(net, pie = countHap(), size = attr(net, "freq"), legend = T)

ここに画像の説明を入力

そして、すべてがうまくいっています！

おまけ：

これはOPから要求されたものではありませんが、他の誰かが興味を持っている場合はここに残すことを考えました. 時々、頻度によってハプロタイプにラベルを付けると便利だと思います。これは、ハプロタイプ ラベルをその頻度と等しくなるように変更することで実行できます。

attr(h, "labels") <- attr(h, "freq")
plot(net, pie = countHap(), size = attr(net, "freq"), legend = T)

ここに画像の説明を入力

于 2015-10-13T15:55:04.090 に答える