[R] plotting issues with PCA

Gavin Simpson gavin.simpson at ucl.ac.uk
Mon Oct 17 13:18:19 CEST 2011


On Mon, 2011-10-17 at 21:02 +1000, Andrew Halford wrote:
> Hi Listers,
> 
> This has a simple answer but it has been eluding me nonetheless.

ordiellipse() doesn't work the way you think it does - it can only take
a single colour at a time. Therefore you need to do /n/ calls to
ordiellipse() to draw /n/ ellipses.

Here is a *reproducible* example, taken from ?ordiellipse:

require(vegan)
data(dune)
data(dune.env)
mod <- cca(dune ~ Management, dune.env)
plot(mod, type="n", display = "sites")
text(mod, display="sites", labels = as.character(Management))

## vector of colours
cols <- c("blue","red","darkgreen","grey70")

## add ellipses
with(dune.env, ordiellipse(mod, Management, kind="se", conf=0.95,
                           lwd=2, col=cols[1], show.groups = "BF"))
with(dune.env, ordiellipse(mod, Management, kind="se", conf=0.95,
                           lwd=2, col=cols[2], show.groups = "HF"))
with(dune.env, ordiellipse(mod, Management, kind="se", conf=0.95,
                           lwd=2, col=cols[3], show.groups = "NM"))
with(dune.env, ordiellipse(mod, Management, kind="se", conf=0.95,
                           lwd=2, col=cols[4], show.groups = "SF"))

We could automate this a bit:

## set up plotting region
plot(mod, type="n", display = "sites")
text(mod, display="sites", labels = as.character(Management))

## get the levels of the factor for plotting groups
lev <- with(dune.env, levels(Management))
## vector of colours
cols <- c("blue","red","darkgreen","grey70")

## loop to draw each group
for (i in seq_along(lev)) {
with(dune.env, ordiellipse(mod, Management, kind="se", conf=0.95,
                           lwd=2,
                           col=cols[i], ## ith colour
                           show.groups = lev[i])) ## for ith group
}

This works with rda() too.

HTH

G

> I have been building a PCA plot from scratch with the ability to plot
> predefined groups in different colors. This has worked fine but when I try
> to get a polygon drawn around each of the groups it is not recognising my
> colour file correctly and is only printing the first colour in the
> file....code is below
> 
>  site.codings <-
> c(1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,2,2,5,4,4,4,4,4,3,3,3,6,6,6,6,6,6,6,6,5,5,5,5)
> > names(site.codings) <- c( "WM1", "WM2", "WM3", "NM1", "NM2", "NM3", "KH1",
> "KH2", "KH3", "LM1" ,"LM2" ,"LM3", "DB1" ,"DB2" ,"DB3", "DM1" , "DM2" ,
> "DM3" , "FI1",  "FI2",  "BKI1", "BKI2", "BKO1", "BKO2", "BKO3",
> "SUR1","MI1","MI2","MI3","BHE1","BHE2","BHE3","BHW1","BHW2","BHW3","HAL1","HAL2","HAL3","HAL4","HAL5","HAL6","HAL7","DOH1","DOH2","DOH3","DOH4","DOH5")
> > fish.pca <-rda(fish.sqrt.h)
> > fish.site <- scores(fish.pca,display="sites",scaling=3)
> > fish.spp <-
> scores(fish.pca,display="species",scaling=3)[omanfish.mrt.indval$pval<=0.05,]
> > graph <- plot(fish.pca,display=c("sites","species"),type="n",scaling=3)
> > plotcolor <-
> c("red","green","blue","aquamarine","magenta","yellow")[site.codings]
> > points(fish.site,pch=21,bg=plotcolor,cex=1.2)
> 
> #up to this point all works well but when I try to draw the polygons I cant
> get the lines to colour code the same way as the points did
> 
> > ordiellipse(graph,site.codings,kind="sd",conf=0.90,draw="polygon")
> 
> I see there is a command called show.groups but I cant work out how to use
> it to access the plotcolor file.
> 
> Any help appreciated.
> 
> 

-- 
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%
 Dr. Gavin Simpson             [t] +44 (0)20 7679 0522
 ECRC, UCL Geography,          [f] +44 (0)20 7679 0565
 Pearson Building,             [e] gavin.simpsonATNOSPAMucl.ac.uk
 Gower Street, London          [w] http://www.ucl.ac.uk/~ucfagls/
 UK. WC1E 6BT.                 [w] http://www.freshwaters.org.uk
%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%~%



More information about the R-help mailing list