Showing posts with label R script. Show all posts
Showing posts with label R script. Show all posts

Saturday, July 18, 2015

After running GARLI, I wanted to compare
OR137:
OR213: lots of disagreements; redo with 100
OR4: lots of disagreements; redo with 100
OR589: lot
OR6: all agree after 8 runs
OR10: all agree after 8 runs
OR11: tree2 disagrees, 7 others agree, use this as main gene tree
OR51: lots of disagreements; redo with 100
OR52: lots of disagreements; redo with 100

100 is a lot of tree to compare, so I wrote a small script doing so:

setwd("/Volumes/Yango/Hayden_batORs/garli/OR51/trees")

#tell where to open up all of your downloaded files from CIPRES
my.files<-list.files("/Volumes/Yango/Hayden_batORs/garli/OR51/trees/")
#read in all of these trees to a list using read.nexus function
my.trees<-lapply(my.files, read.nexus)

#create a matrix to store all of the tree comparisons
tree.matrix<-matrix(nrow = length(my.trees), ncol = length(my.trees))

for(i in 1:length(my.trees)){ #for every tree in the file
  this.tree<-as.phylo(my.trees[[i]]) #look at that tree
  #compare the topology to all the other trees
  tree.matrix[i,]<-sapply(my.trees, function(my.trees) dist.topo(this.tree,my.trees))
}
tree.matrix

#print out the ones that agree with each other
for(i in 1:length(tree.matrix[1,])){
  this.list<-as.numeric(which(tree.matrix[i,]==0))
  print(this.list)
}
write.csv(tree.matrix, "tree.matrix.csv")

For OR51, are things looking better?
[1] 1
[1] 2
[1] 3
[1]  4 84
[1] 5
[1] 6
[1] 7
[1]  8 15 41 79
[1]  9 14
[1] 10 54
[1] 11
[1] 12
[1] 13 73
[1]  9 14
[1]  8 15 41 79
[1] 16
[1] 17 20 75
[1] 18
[1] 19
[1] 17 20 75
[1] 21 51 68
[1] 22 37
[1] 23 50 71
[1] 24
[1] 25
[1] 26
[1] 27
[1] 28 69 88
[1] 29
[1] 30
[1] 31
[1] 32 46
[1] 33
[1] 34 70
[1] 35 53
[1] 36
[1] 22 37
[1] 38 59
[1] 39 47 78
[1] 40
[1]  8 15 41 79
[1] 42
[1] 43
[1] 44
[1] 45
[1] 32 46
[1] 39 47 78
[1] 48
[1] 49
[1] 23 50 71
[1] 21 51 68
[1] 52
[1] 35 53
[1] 10 54
[1] 55
[1] 56
[1] 57
[1] 58
[1] 38 59
[1] 60
[1] 61
[1] 62 86 87
[1] 63
[1] 64
[1] 65
[1] 66 95
[1] 67
[1] 21 51 68
[1] 28 69 88
[1] 34 70
[1] 23 50 71
[1] 72
[1] 13 73
[1] 74
[1] 17 20 75
[1] 76
[1] 77
[1] 39 47 78
[1]  8 15 41 79
[1] 80
[1] 81
[1] 82
[1] 83
[1]  4 84
[1] 85
[1] 62 86 87
[1] 62 86 87
[1] 28 69 88
[1] 89
[1] 90
[1] 91
[1] 92
[1] 93
[1] 94
[1] 66 95
[1] 96
[1] 97
[1] 98

Womp, womp...not really--only a maximum of 4 trees out of 100 are agreeing.

OR4s looking much better:
[1]  1 22 28 42 75
[1]  2 99
[1]  3 14 18 83 89
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1] 5
[1] 6
[1]  7 15 98
[1] 8
[1]  9 19 38 80
[1] 10 50
[1] 11 67
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1] 13 86
[1]  3 14 18 83 89
[1]  7 15 98
[1] 16
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1]  3 14 18 83 89
[1]  9 19 38 80
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1] 21 31 44 62 84
[1]  1 22 28 42 75
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1] 25
[1] 26
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1]  1 22 28 42 75
[1] 29
[1] 30 33
[1] 21 31 44 62 84
[1] 32 37 61 68 81
[1] 30 33
[1]  34 100
[1] 35
[1] 36
[1] 32 37 61 68 81
[1]  9 19 38 80
[1] 39 71
[1] 40 49
[1] 41
[1]  1 22 28 42 75
[1] 43
[1] 21 31 44 62 84
[1] 45
[1] 46
[1] 47
[1] 48
[1] 40 49
[1] 10 50
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1] 52
[1] 53
[1] 54
[1] 55
[1] 56 90
[1] 57
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1] 60
[1] 32 37 61 68 81
[1] 21 31 44 62 84
[1] 63
[1] 64
[1] 65
[1] 66
[1] 11 67
[1] 32 37 61 68 81
[1] 69
[1] 70
[1] 39 71
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1] 73
[1] 74
[1]  1 22 28 42 75
[1] 76
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1] 78
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1]  9 19 38 80
[1] 32 37 61 68 81
[1] 82
[1]  3 14 18 83 89
[1] 21 31 44 62 84
[1] 85
[1] 13 86
[1] 87
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1]  3 14 18 83 89
[1] 56 90
[1] 91
[1] 92
[1] 93
[1] 94
[1] 95
 [1]  4 12 17 20 23 24 27 51 58 59 72 77 79 88 96
[1] 97
[1]  7 15 98
[1]  2 99
[1]  34 100

Friday, June 14, 2013

Finally have some meaningful(?) output. This is going to be lengthy but it will be really helpful in organizing my thought. Because I had been including morphological data for very "distantly-related babblers", all of the nodes for the rate shifts occurred, of course on the root node. I removed some more taxa to only include measurements for birds within the babbler family. Therefore, I reran the PCA. Still similar trends:

Count the number of species in my analysis and the number of specimens per species.
> babblers.species.count<-ddply(.data=new.babblers, .(Genus, Species), species.count=length(Wing.Length), summarize)
> babblers.species.count
             Genus         Species species.count
1       Actinodura        egertoni             6
2       Actinodura    morrisoniana             3
3       Actinodura      nipalensis             3
4       Actinodura         waldeni             5
5          Alcippe    brunneicauda             4
6          Alcippe          grotei            11
7          Alcippe      morrisonia             5
8          Alcippe      peracensis             2
9          Alcippe    poioicephala            11
10           Babax     lanceolatus             4
11           Cutia      nipalensis             7
12         Dumetia      hyperythra             6
13  Gampsorhynchus         rufulus             8
14        Garrulax         affinis            11
15        Garrulax     albogularis             5
16        Garrulax         canorus            14
17        Garrulax       chinensis            40
18        Garrulax      cineraceus             4
19        Garrulax erythrocephalus             1
20        Garrulax        formosus             4
21        Garrulax     leucolophus            41
22        Garrulax        lineatus             6
23        Garrulax           maesi             3
24        Garrulax       merulinus             1
25        Garrulax         milleti             1
26        Garrulax          milnei             4
27        Garrulax        mitratus             6
28        Garrulax       monileger             8
29        Garrulax        nuchalis             3
30        Garrulax       ocellatus             7
31        Garrulax       palliatus             4
32        Garrulax poecilorhynchus             4
33        Garrulax          sannio            12
34        Garrulax       squamatus             3
35        Garrulax        striatus             6
36        Garrulax     subunicolor             6
37        Garrulax      variegatum             3
38        Garrulax        virgatus             2
39    Heterophasia     auricularis             3
40    Heterophasia      capistrata             6
41    Heterophasia        gracilis             4
42    Heterophasia    maelanoleuca             9
43    Heterophasia       picaoides             7
44    Heterophasia       pulchella             3
45     Illadopsis       albipectus             4
46     Illadopsis         cleaveri             4
47     Illadopsis       fulvescens             8
48     Illadopsis           puveli             2
49     Illadopsis      pyrrhoptera             2
50     Illadopsis        rufescens             3
51     Illadopsis       rufipennis             6
52         Kenopia         striata             2
53       Leiothrix     argentauris            15
54       Leiothrix           lutea             4
55       Liocichla       phoenicea             1
56       Liocichla          steeri             3
57        Macronus         gularis            25
58        Macronus        ptilosus             3
59        Macronus     striaticeps            12
60    Malacocincla         abbotti             4
61    Malacocincla     cinereiceps             2
62    Malacocincla     malaccensis             6
63    Malacocincla        sepiaria             3
64    Malacopteron          affine             3
65    Malacopteron      albogulare             4
66    Malacopteron        cinereum             1
67    Malacopteron     magnirostre             7
68    Malacopteron          magnum             5
69    Malacopteron    palawanenese             3
70           Minla      ignotincta             6
71           Minla        strigula             7
72       Napothera    brevicaudata             4
73       Napothera          crassa             2
74       Napothera     crispifrons             7
75       Napothera     epilepidota            10
76      Pellorneum      albiventre             9
77      Pellorneum    capristratum            15
78      Pellorneum      pyrrogenys             6
79      Pellorneum        ruficeps             5
80      Pellorneum        tickelli             1
81     Phyllanthus      atripennis             6
82    Pomatorhinus    erythrogenys            18
83    Pomatorhinus    ferruginosus             4
84    Pomatorhinus      horsfieldi             6
85    Pomatorhinus      hypoleucos             4
86    Pomatorhinus        montanus             4
87    Pomatorhinus      ruficollis             9
88    Pomatorhinus        swinhoei             2
89     Ptilocichla         falcata             2
90     Ptilocichla   leucogrammica             3
91     Ptilocichla     mindanensis             3
92       Ptyrticus        turdinus             6
93    Rhopocichla         atriceps             4
94         Rimator    malacoptilus             3
95    Schoeniparus         brunnea            10
96    Schoeniparus     castaneceps             1
97    Schoeniparus         cinerea             3
98    Schoeniparus           dubia             6
99    Schoeniparus     rufogularis             9
100    Spelaeornis    chocolatinus             5
101    Spelaeornis   troglodytoide             3
102   Sphenocichla           humei             3
103      Stachyris         ambigua             4
104      Stachyris        chrysaea            11
105      Stachyris    erythroptera             8
106      Stachyris      grammiceps             2
107      Stachyris        leucotis             3
108      Stachyris        maculata             3
109      Stachyris    melanothorax             6
110      Stachyris       nigriceps             7
111      Stachyris     nigricollis             4
112      Stachyris           oglei             3
113      Stachyris    poliocephala             4
114      Stachyris        pyrrhops             3
115      Stachyris        ruficeps            10
116      Stachyris       rufifrons             3
117      Stachyris       striolata             5
118      Stachyris       thoracica             4
119        Timalia         pileata             8
120    Trichastoma         bicolor             4
121    Trichastoma       celebense            11
122    Trichastoma       rostratum             5
123      Turdoides         bicolor             3
124      Turdoides         gularis             4
125      Turdoides       jardineii            10
126      Turdoides        plebejus            16
127      Turdoides      reinwardii             3
128  Xiphirhynchus   superciliaris             3

I had to match up the names of my species with the names on the tips (huge pain) and I ended up doing it by hand--I could have probably scripted it but I think it would have taken me just as long. 

I collapsed the PCA values into means for each species but I believe I already posted the script for that. So now I have a .csv that contains all the mean PCA that corresponds to each species with the species name that matches the tip. YAY.

Now I finally have evol.rate.mcmc() running and I ran the output for the first four principle components:
#read in the PCA vector:
temp<-read.csv(file.choose())
x<-as.vector(temp$Tarsus) #vector of tarsus measurement traits

names(x)<-as.vector(temp$Species) #names the measurements
#i have to date the trees
undated.tree<-read.tree(file.choose())
#the tree was based on two calibration points (see Moyle, et al. 2012) and the root

tree<-chronopl(undated.tree, lambda=.5, age.min=c(8.4, 13, 27.1), age.max=c(11.4,17, 37.3), node = c(336, 297, 293)) 

#now time for the juicy part--using default parameters to start with as suggested by Revell, et al. 2011

result<-evol.rate.mcmc(tree,x,ngen=100000,control=list(sd2=2.0,sdlnr=3.0))
min_split<-minSplit(tree,result$mcmc[21:101,c("node","bp")]) 
mcmc<-posterior.evolrate(tree,min_split,result$mcmc[21:101,],result$tips[201:1001])
mcmc.data.1<-as.data.frame(mcmc)
#nicely summarize the mcmc for Tarsus
output.summary.1<-ddply(.data=mcmc.data.1, .(node), node.count=length(bp), mean.bp=mean(bp),mean.sig1=mean(sig1), mean.sig2=mean(sig2), summarize)

#let's visualize the nodes! first trim the tree to visualize only nodes that were used in analysis
trim.tree<-drop.tip(tree,tree$tip.label[-(match(names(x),tree$tip.label))]) 

Nodes 3 and 4 are probably root nodes and have high probability. Let's look at 186, 185, and 188.

#now extract the clades we want to look at

this.extract<-extract.clade(tree.trim, tree.trim$edge[186])
this.extract<-extract.clade(tree.trim, tree.trim$edge[58])
#perhaps the most interesting because it contains nearly the whole subfamily
this.extract<-extract.clade(tree.trim, tree.trim$edge[4])
Now time for the bill length:
#read in the PCA vector:
temp<-read.csv(file.choose())
x<-as.vector(temp$Bill.Length) #vector of bill length measurement traits

#you have to reset the names
names(x)<-as.vector(temp$Species) #names the measurements

#use the same dated tree as before
result<-evol.rate.mcmc(tree,x,ngen=100000,control=list(sd2=2.0,sdlnr=3.0))
min_split<-minSplit(tree,result$mcmc[21:101,c("node","bp")]) 
mcmc<-posterior.evolrate(tree,min_split,result$mcmc[21:101,],result$tips[201:1001])
mcmc.data.2<-as.data.frame(mcmc)
#nicely summarize the mcmc for Bill Length
output.summary.2<-ddply(.data=mcmc.data.2, .(node), node.count=length(bp), mean.bp=mean(bp),mean.sig1=mean(sig1), mean.sig2=mean(sig2), summarize)

Hmm...very strange. 250 and 251 are only tips of two sister taxa.
> this.extract<-extract.clade(tree.trim, tree.trim$edge[122])
> plot(this.extract)

Now for the Hallux:
x<-as.vector(temp$Hallux) #vector of hallux measurement traits

#you have to reset the names
names(x)<-as.vector(temp$Species) #names the measurements

#use the same dated tree as before
result<-evol.rate.mcmc(tree,x,ngen=100000,control=list(sd2=2.0,sdlnr=3.0))
min_split<-minSplit(tree,result$mcmc[21:101,c("node","bp")]) 
mcmc<-posterior.evolrate(tree,min_split,result$mcmc[21:101,],result$tips[201:1001])
mcmc.data.2<-as.data.frame(mcmc)
#nicely summarize the mcmc for hallux
output.summary.2<-ddply(.data=mcmc.data.3, .(node), node.count=length(bp), mean.bp=mean(bp),mean.sig1=mean(sig1), mean.sig2=mean(sig2), summarize)

Node 9 has high probability (low frequency):
> this.extract<-extract.clade(tree.trim, tree.trim$edge[9])
> plot(this.extract)

> this.extract<-extract.clade(tree.trim, tree.trim$edge[29])
> plot(this.extract)

> this.extract<-extract.clade(tree.trim, tree.trim$edge[235])
> plot(this.extract)

Enough of that, let's finally look at Tail Length and see if there is anything interesting:
x<-as.vector(temp$Tail.Length) #vector of tail measurement traits

#you have to reset the names
names(x)<-as.vector(temp$Species) #names the measurements

#use the same dated tree as before
result<-evol.rate.mcmc(tree,x,ngen=100000,control=list(sd2=2.0,sdlnr=3.0))
min_split<-minSplit(tree,result$mcmc[21:101,c("node","bp")]) 
mcmc<-posterior.evolrate(tree,min_split,result$mcmc[21:101,],result$tips[201:1001])
mcmc.data.2<-as.data.frame(mcmc)
#nicely summarize the mcmc for tail length
output.summary.4<-ddply(.data=mcmc.data.4, .(node), node.count=length(bp), mean.bp=mean(bp),mean.sig1=mean(sig1), mean.sig2=mean(sig2), summarize)



> this.extract<-extract.clade(tree.trim, tree.trim$edge[8])
> plot(this.extract)

Okay I am kind of petering out on uploading all of these. Gives me something to compare to though.

Important questions:
Next step is to start making my slides for my presentation and date random samples of trees and get feedback from adviser to make sure I am on the right track. 

Thursday, June 6, 2013

Had a minor panic attack because I realized my missing tarsus values were being imputed with mean of the bill length which of course would lead to extreme variance in the trait. Anyways, fixed it right up and still got somewhat similar results. Whew.
New PCA result:


Fixed script:
#impute missing values by taking mean of all other values of of other species ofthat character 

i<-1 #set the index
species.removed<-c()
#sift through all babbler species
while (i <= length(babblers$UniqueID_1)){
  #missing wing lengths; if the Wing Length is NA for that value
  if(is.na(babblers$Wing.Length[i]) == TRUE){
    #select the wing lengths of all birds which match this species with missing wing length
    wing.species.vector<-babblers$Wing.Length[which(babblers$Species[i] == babblers$Species)]
    #make a new vector without the missing value
    wing.species.filler<-wing.species.vector[which(is.na(wing.species.vector) == FALSE)]
    #calculate the mean
    wing.impute<-mean(wing.species.filler)
    #fill in missing value with the mean for that bird
    babblers$Wing.Length[i]<-wing.impute
  }
  
  #repeat the above method for all characters
  if(is.na(babblers$Tail.Length[i]) == TRUE){
    tail.species.vector<-babblers$Tail.Length[which(babblers$Species[i] == babblers$Species)]
    tail.species.filler<-tail.species.vector[which(is.na(tail.species.vector) == FALSE)]
    tail.impute<-mean(tail.species.filler)
    babblers$Tail.Length[i]<-tail.impute
  }
  
  
  
  if(is.na(babblers$Bill.Length[i]) == TRUE){
    BL.species.vector<-babblers$Bill.Length[which(babblers$Species[i] == babblers$Species)]
    BL.species.filler<-BL.species.vector[which(is.na(BL.species.vector) == FALSE)]
    BL.impute<-mean(BL.species.filler)
    babblers$Bill.Length[i]<-BL.impute
  }
  
  if(is.na(babblers$Bill.Width[i]) == TRUE){
    BW.species.vector<-babblers$Bill.Width[which(babblers$Species[i] == babblers$Species)]
    BW.species.filler<-BW.species.vector[which(is.na(BW.species.vector) == FALSE)]
    BW.impute<-mean(BW.species.filler)
    babblers$Bill.Width[i]<-BW.impute
  }
  
  if(is.na(babblers$Bill.Depth[i]) == TRUE){
    BD.species.vector<-babblers$Bill.Depth[which(babblers$Species[i] == babblers$Species)]
    BD.species.filler<-BD.species.vector[which(is.na(BD.species.vector) == FALSE)]
    BD.impute<-mean(BD.species.filler)
    babblers$Bill.Depth[i]<-BD.impute
  }
  
  if(is.na(babblers$Tarsus.Length[i]) == TRUE){
    TL.species.vector<-babblers$Tarsus.Length[which(babblers$Species[i] == babblers$Species)]
    TL.species.filler<-TL.species.vector[which(is.na(TL.species.vector) == FALSE)]
    TL.impute<-mean(TL.species.filler)
    babblers$Tarsus.Length[i]<-TL.impute
  }
  
  if(is.na(babblers$Hallux[i]) == TRUE){
    hallux.species.vector<-babblers$Hallux[which(babblers$Species[i] == babblers$Species)]
    hallux.species.filler<-hallux.species.vector[which(is.na(hallux.species.vector) == FALSE)]
    hallux.impute<-mean(hallux.species.filler)
    babblers$Hallux[i]<-hallux.impute
  }

  i<-i+1
}


PCAs and reweights of the values for each species is finally done and in a format that can be easily used to run the evol.rate.mcmc()


Script to impute data:
#impute missing values by taking mean of all other values of of other species ofthat character 

i<-1 #set the index
species.removed<-c()
#sift through all babbler species
while (i <= length(babblers$UniqueID_1)){
  #missing wing lengths; if the Wing Length is NA for that value
  if(is.na(babblers$Wing.Length[i]) == TRUE){
    #select the wing lengths of all birds which match this species with missing wing length
    wing.species.vector<-babblers$Wing.Length[which(babblers$Species[i] == babblers$Species)]
    #make a new vector without the missing value
    wing.species.filler<-wing.species.vector[which(is.na(wing.species.vector) == FALSE)]
    #calculate the mean
    wing.impute<-mean(wing.species.filler)
    #fill in missing value with the mean for that bird
    babblers$Wing.Length[i]<-wing.impute
  }
  
  #repeat the above method for all characters
  if(is.na(babblers$Tail.Length[i]) == TRUE){
    tail.species.vector<-babblers$Tail.Length[which(babblers$Species[i] == babblers$Species)]
    tail.species.filler<-tail.species.vector[which(is.na(tail.species.vector) == FALSE)]
    tail.impute<-mean(tail.species.filler)
    babblers$Tail.Length[i]<-tail.impute
  }
  
  
  
  if(is.na(babblers$Bill.Length[i]) == TRUE){
    BL.species.vector<-babblers$Bill.Length[which(babblers$Species[i] == babblers$Species)]
    BL.species.filler<-BL.species.vector[which(is.na(BL.species.vector) == FALSE)]
    BL.impute<-mean(BL.species.filler)
    babblers$Bill.Length[i]<-BL.impute
  }
  
  if(is.na(babblers$Bill.Width[i]) == TRUE){
    BW.species.vector<-babblers$Bill.Width[which(babblers$Species[i] == babblers$Species)]
    BW.species.filler<-BW.species.vector[which(is.na(BW.species.vector) == FALSE)]
    BW.impute<-mean(BW.species.filler)
    babblers$Bill.Width[i]<-BW.impute
  }
  
  if(is.na(babblers$Bill.Depth[i]) == TRUE){
    BD.species.vector<-babblers$Bill.Depth[which(babblers$Species[i] == babblers$Species)]
    BD.species.filler<-BD.species.vector[which(is.na(BD.species.vector) == FALSE)]
    BD.impute<-mean(BD.species.filler)
    babblers$Bill.Depth[i]<-BD.impute
  }
  
  if(is.na(babblers$Tarsus.Length[i]) == TRUE){
    TL.species.vector<-babblers$Tarsus.Length[which(babblers$Species[i] == babblers$Species)]
    TL.species.filler<-TL.species.vector[which(is.na(TL.species.vector) == FALSE)]
    TL.impute<-mean(BD.species.filler)
    babblers$Tarsus.Length[i]<-TL.impute
  }
  
  if(is.na(babblers$Hallux[i]) == TRUE){
    hallux.species.vector<-babblers$Hallux[which(babblers$Species[i] == babblers$Species)]
    hallux.species.filler<-hallux.species.vector[which(is.na(hallux.species.vector) == FALSE)]
    hallux.impute<-mean(hallux.species.filler)
    babblers$Hallux[i]<-hallux.impute
  }

  i<-i+1
}

Script to scale by geometric mean:
#make empty columns for geometric mean scale
babblers$gm.wing<-0
babblers$gm.tail<-0
babblers$gm.BL<-0
babblers$gm.BW<-0
babblers$gm.BD<-0
babblers$gm.tarsus<-0
babblers$gm.hallux<-0

#Size-adjust the raw data by dividing each individual by the geometric mean 

j<-1
could.not.scale<-c()
while(j<=length(babblers$UniqueID_1)){
  #grab all the data of that species and put in vector
  temp<-c(babblers$Wing.Length[j], babblers$Tail.Length[j], babblers$Bill.Length[j],
          babblers$Bill.Width[j], babblers$Bill.Depth[j], babblers$Tarsus.Length[j], 
          babblers$Hallux[j])
  temp<-na.omit(temp)
  #calculate the geometric mean of meausrements for that species
  geo.mean<-exp(mean(log(temp)))
  
  if(is.nan(geo.mean)==FALSE){
  #divide all values for that species by the geometric mean
  geo.wing<-(babblers$Wing.Length[j])/geo.mean
  babblers$gm.wing[j]<-geo.wing
  
  geo.tail<-babblers$Tail.Length[j]/geo.mean
  babblers$gm.tail[j]<-geo.tail
  
  geo.BL<-babblers$Bill.Length[j]/geo.mean
  babblers$gm.BL[j]<-geo.BL
  
  geo.BW<-babblers$Bill.Width[j]/geo.mean
  babblers$gm.BW[j]<-geo.BW
  
  geo.BD<-babblers$Bill.Depth[j]/geo.mean
  babblers$gm.BD[j]<-geo.BD
  
  geo.tarsus<-babblers$Tarsus.Length[j]/geo.mean
  babblers$gm.tarsus[j]<-geo.tarsus
  
  geo.hallux<-babblers$Hallux[j]/geo.mean
  babblers$gm.hallux[j]<-geo.hallux
 }
 else{
    append(could.not.scale, babblers$UniqueID_1[j])
  }
  #reset the temp
  temp<-c()
  
  #increment to next species
  j<-j+1
}

cat("The following species could not be scaled:", could.not.scale, "\n")

Script to do PCA with geometric mean values and put rescaled scores into new data frame:
#PCA for all measurements
pca.result <- prcomp(~babblers$gm.wing + babblers$gm.tail + babblers$gm.BL + babblers$gm.BW + babblers$gm.BD + 
       babblers$gm.tarsus + babblers$gm.hallux, data = babblers, retx = T)

#standard deviation of PCA
sd <- pca.result$sdev
#make vector of weights
loadings <- pca.result$rotation
#value of the rotated data (scores)
scores <- pca.result$x

#make a new data frame
new.babblers<-as.data.frame(babblers)

#make empty columns in the data frame
new.babblers$pca1<-0
new.babblers$pca2<-0
new.babblers$pca3<-0
new.babblers$pca4<-0
new.babblers$pca5<-0
new.babblers$pca6<-0
new.babblers$pca7<-0

#put PCA Scores into the data frame
new.babblers$pca1<-scores[,1]
new.babblers$pca2<-scores[,2]
new.babblers$pca3<-scores[,3]
new.babblers$pca4<-scores[,4]
new.babblers$pca5<-scores[,5]
new.babblers$pca6<-scores[,6]
new.babblers$pca7<-scores[,7]

new.babblers$pca# are the traits that will be used in the Bayesian analysis.
This site was really useful for understanding PCA:
http://strata.uga.edu/software/pdf/pcaTutorial.pdf