# pcoa questions

**URL:** https://forum.mothur.org/t/pcoa-questions/2441
**Category:** Commands in mothur
**Created:** [July 7, 2015, 3:47pm UTC](https://forum.mothur.org/t/pcoa-questions/2441 "2015-07-07T15:47:50Z")
**Posts on this page:** 3
**Page:** 1

<div class="post-metadata">

### Author: ![suttee2011](https://avatars.discourse-cdn.com/v4/letter/s/898d66/32.png) [@suttee2011](https://forum.mothur.org/u/suttee2011)
#### Post date: [July 7, 2015, 3:47pm UTC](https://forum.mothur.org/t/pcoa-questions/2441/1 "2015-07-07T15:47:50Z")

</div>

🙂 Hi everyone，  
I have some questions about “pcoa” command.

1. The output ofthis command are two files ending in \*pcoa, and \*pcoa.loadings,but I can’t find the file ending in \*pcoa. However, there is a file ending in \*axes.
2. Most PCoA are based on unweighted unifrac or weighted unifrac. Is this command also based on unifrac? If not, how can I do PCoA based on unifrac?
3. Who can give me a R script to plot PCoA results  
Thanks for your help! 🙂

---

<div class="post-metadata">

### Author: ![dwaite](https://avatars.discourse-cdn.com/v4/letter/d/ecccb3/32.png) [@dwaite](https://forum.mothur.org/u/dwaite)
#### Post date: [July 7, 2015, 8:21pm UTC](https://forum.mothur.org/t/pcoa-questions/2441/2 "2015-07-07T20:21:22Z")

</div>

1. The axes file contains the plotting coordinates for graphing the results so this is what you need. Maybe the \*.pcoa is a typo in the logfile or something.

2. PCoA is just run on any distance matrix. It’s up to you what you want to run it on, you just need a phylip-formatted matrix to proceed. As you note, you can get these from the unifrac.unweighted/unifrac.weighted commands, but also from the dist.shared command (Bray-Curtis, Yue-Clayton theta, Jaccard etc).

3. If you want anything fancy you’ll need to have a play around a bit, but for a very basic plot:

```nohighlight
pcoa_data = read.table(file="your.pcoa.file.axes", header=T, sep='\t')
plot(pcoa_data$Axis2~pcoa_data$Axis1,xlab='Axis 1',ylab='Axis 2',main='Blah blah PCoA...')

```

If you have a design file that matches the groups processed in your distance matrix you can use the following to get group-based colouring or symbols:

```nohighlight
design_file = read.table(file="your.design.file", sep='\t',stringsAsFactors=FALSE)

uniq_groups = unique(design_file[,2])
num_uniq_groups = length(uniq_groups)

avail_cols = rainbow(num_uniq_groups)
avail_symbols = 1:num_uniq_groups

colour_map = list()
symbol_map = list()

for(i in 1:num_uniq_groups) {
 current_group = uniq_groups[i]
 colour_map[[current_group]] = avail_cols[i]
 symbol_map[[current_group]] = avail_symbols[i]
}

col_vector = rep(NA,nrow(design_file))
symbol_vector = rep(NA,nrow(design_file))

for(i in 1:nrow(design_file)) {
 group = design_file[i,2]
 col_vector[i] = colour_map[[group]]
 symbol_vector [i] = symbol_map [[group]] 
}

```

Then you could plot either by colours:

```nohighlight
plot(pcoa_data$Axis2~pcoa_data$Axis1,xlab='Axis 1',ylab='Axis 2',main='Blah blah PCoA...',pch=21,col='black',bg=col_vector)
legend(X,Y,uniq_groups,avail_cols)

```

Or symbols:

```nohighlight
plot(pcoa_data$Axis2~pcoa_data$Axis1,xlab='Axis 1',ylab='Axis 2',main='Blah blah PCoA...',pch=symbol_vector)

```

You’ll just need to manually work out the X/Y values for placing your legend (hopefully there will be a bit of an empty spot on your graph where you can stick it).

---

<div class="post-metadata">

### Author: ![suttee2011](https://avatars.discourse-cdn.com/v4/letter/s/898d66/32.png) [@suttee2011](https://forum.mothur.org/u/suttee2011)
#### Post date: [July 8, 2015, 2:12am UTC](https://forum.mothur.org/t/pcoa-questions/2441/3 "2015-07-08T02:12:47Z")

</div>

Hi dwaite,  
Thank you for the quick response. Your reply give me invaluable help! 🙂
