Pages

Saturday, September 15, 2007

Neighbor Joining Tree with Ape

Today I used R to create a neighbor joining phenogram, using the ape library. The phenogram will be used to visualize a species genetic differences between six different locations in Costa Rica. Genetic distances were calculated from microsatellite data using GenAlEx (GenAlEx), because I haven't figured out how to calculate Nei's genetic dissimilarities with R. This is the procedure I used:

First I imported the genetic distance matrix from Exc..l. The matrix used was symmetric, meaning that the same distances were mirrored above and below the diagonal. I copied the matrix into the clipboard from Exc..l using Ctrl-C. The distance matrix included row and coloumn names, the six different locations. The matrix was imported in R by assigning it to an object named m:

> m <- as.matrix(read.table("clipboard", head=T, row.names=1))

The previous command imports the data stored in the clipboard and transforms the data frame into a matrix format. The options of the read.table() command: head=T and row.names=1, tell R that the first line of the matrix is a header and that row names are stored in the first column, respectively. After checking the matrix was imported correctly, I proceeded to load the ape library.

> library(ape)

The ape (analysis in phylogenetics and evolution) contains various methods for the analysis of genetic and evolutionary data. Ape also provides commands designed to calculate distances in DNA sequences. I will tackle those in another post. Now back to my tree. After importing the distance matrix, and including the ape library, I proceeded to create the phylogenetic tree with the nj() command:

> arbol <- nj(as.dist(m), "unrooted") > plot(arbol)

The "unrooted" parameter provided to the nj() command produces, as expected, an unrooted tree. The resulting tree topography was stored in the 'arbol' object. The plot command produced:

We can observe how 'Puriscal' is located between two well demarcated groups, one including 'Guapiles' and 'Guatuso' and the other three populations in the second group.


5 comments:

  1. Thanks, this post was useful to me.

    hannes

    ReplyDelete
  2. adegenet package has a function to estimate Neis Genetic distance.

    dist.genpop(...,method=1)

    Gives Nei's distance, metod=3 Cavalli Sforza....

    Nice post thanks!

    ReplyDelete
  3. when i define arbol, it comes like this :

    Error in nj(as.dist(M), "unrooted") : unused argument ("unrooted")

    What should i do?

    ReplyDelete
    Replies
    1. This comment has been removed by the author.

      Delete
    2. Sorry, I made a mistake. The unrooted option is for the plot command:

      > plot(arbol, type="unrooted")

      Delete