R-如何在haploNet haplotyp网络中绘制正确的饼图{pegas} {ape} {adegenet} [英] R- How to plot correct pie charts in haploNet haplotyp Networks {pegas} {ape} {adegenet}
问题描述
使用haploNet软件包在单倍型网络上绘制某些图时， 我使用互联网上可用的脚本来执行此操作.但是我认为这是有问题的.该脚本以woodmouse示例的形式提供.我使用的代码是:
When using the haploNet package to make some plots on a haplotype network, I used a script available on the internet to do so. However I think there is something wrong. The script is available in form of the woodmouse example. The code I used is:
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时，您会注意到某些行不在正确的位置.您可以在这里看到它:
However when plotting ind.hap, you can notice that some rows are not in the right place. You can see this here:
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行进行绘制以绘制IX的饼图，该饼图是VIII的数据.结果如下: (由于我的信誉低于10，所以我无法插入图片...，无论如何您都可以通过执行整个文件来获取图片)
You can see that row IX is not on its right place. This would not be too much of a problem, but the program takes row 9 to make the pie plot for IX, which is the data of VIII. This is the result: (I could not insert the image since my reputation is below 10..., you get the image by executing the whole file anyway)
您可以看到，直到IX为止，V都不是应有的状态(这些是交换的行).例如:IX中只有1个单倍型，但是有一个饼图用于2个单倍型(它们都占图表的50％)，这是使用VIII数据生成的.由于行是按字母顺序而不是升序排序的，但这是包固有的，所以我不知道该怎么做. 我离R的掌握还很远，所以不要太抽象，而要提供代码.
You can see that for V until IX it's not as it should be (these are the swapped rows). For instance: IX has only 1 haplotype in it, but there's a pie chart for 2 haplotypes (which both have 50% of the chart), which is generated using the VIII data. Since the rows are sorted alphabetically instead of ascending, but this is inherent to the package, I don't know what to do. I'm far from a master in R, so try not to be too abstract, but provide code instead.
如果有人对这个软件包非常了解，请解释一下为什么真实图表后面有这些奇怪的多余线(上面有数字)，因为它们在木鼠标示例中是不可见的(也许是因为还有什么问题吗?)
If there is someone who knows this package very well, please explain also why there are these weird extra lines behind the real charts (these with the numbers on them), as they were not visible in the woodmouse example (maybe that's because of what's wrong too?)
先感谢
推荐答案
我一直在努力解决同一问题，但相信我提出了一个解决方案.
I've struggled with the same issue, but believe I came up with a solution.
问题在于，对每个人口"的单倍型计数table的步骤按字母顺序对单倍型进行排序.因此，例如，单倍型"IX"在"V"之前.另一方面，函数haplotype()按单倍型的数字"顺序对其进行排序.这就是在绘图时产生差异的原因.
The problem is that the step making the table of haplotype counts per "population" orders the haplotypes alphabetically. So, for example, haplotype "IX" comes before "V". On the other hand, the function haplotype() sorts the haplotypes by their "numerical" order. And this is what creates a discrepancy when plotting.
可以通过按标签"对单倍型对象进行排序来解决此问题，如?haplotype帮助中所述.
This can be solved by sorting the haplotype object by "label", as explained in ?haplotype help.
我将使用woodmouse示例数据进行举例说明:
I'll use the woodmouse example data to exemplify:
# Sample 9 distinct haplotypes
library(pegas)
data(woodmouse)
x <- woodmouse[sample(9, 100, replace = T), ]
To simplify, I create a function to create the count table of haplotypes (based on this post):
countHap <- function(hap = h, dna = x){
with(
stack(setNames(attr(hap, "index"), rownames(hap))),
table(hap = ind, pop = attr(dna, "dimnames")[[1]][values])
)
}
现在，让我们在不对单元型进行排序的情况下查看结果:
Now, let's see the result without sorting haplotypes:
h <- haplotype(x) # create haplotype object
net <- haploNet(h) # create haploNet object
plot(net, pie = countHap(), size = attr(net, "freq"), legend = T)
现在，让我们看一下计数表，检查以下结果:
Now, let's look at our count table, to check these results:
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"的标签).
Things do not match: for example, haplotype "V" should occur in individual "No0908S", but instead is coloured as individual "No0913S" (which should be the label for haplotype "IX").
现在，让我们对单倍型进行排序:
Now, let's sort haplotypes:
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并没有要求这样做，但我认为如果其他人都感兴趣，可以将其留在这里. 有时，我发现按频率标记单倍型很方便.可以通过将单倍型标签更改为等于其频率来实现:
Although this is not requested by the OP, I thought of leaving it here if it is of interest for anyone else. Sometimes, I find it convenient to label haplotypes by their frequency. This can be done by changing the haplotype labels to be equal to their frequencies:
attr(h, "labels") <- attr(h, "freq")
plot(net, pie = countHap(), size = attr(net, "freq"), legend = T)
这篇关于R-如何在haploNet haplotyp网络中绘制正确的饼图{pegas} {ape} {adegenet}的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
