r 语言关于 ses 命令结果显示问题，关键就是要pd的sd

各位大虾，本人在用r语言做论文，想求通过ses命令求得sd等值，但运行下列命令后，不知道如何看到具体结果，求助啊？？？？？？？？
命令如下
`ses.pd` <-
function (samp, tree, null.model = c("taxa.labels", "richness", "frequency", "sample.pool",
    "phylogeny.pool", "independentswap", "trialswap"), runs = 999, iterations = 1000, ...)
{
    pd.obs <- as.vector(pd(samp, tree, ...)$PD)
    null.model <- match.arg(null.model)
    pd.rand <- switch(null.model,
            taxa.labels = t(replicate(runs, as.vector(pd(samp, tipShuffle(tree), ...)$PD))),
            richness = t(replicate(runs, as.vector(pd(randomizeMatrix(samp,null.model="richness"), tree, ...)$PD))),
            frequency = t(replicate(runs, as.vector(pd(randomizeMatrix(samp,null.model="frequency"), tree, ...)$PD))),           
            sample.pool = t(replicate(runs, as.vector(pd(randomizeMatrix(samp,null.model="richness"), tree, ...)$PD))),
            phylogeny.pool = t(replicate(runs, as.vector(pd(randomizeMatrix(samp,null.model="richness"),
                    tipShuffle(tree), ...)$PD))),
            independentswap = t(replicate(runs, as.vector(pd(randomizeMatrix(samp,null.model="independentswap", iterations), tree,

...)$PD))),
            trialswap = t(replicate(runs, as.vector(pd(randomizeMatrix(samp,null.model="trialswap", iterations), tree, ...)$PD)))
    )
    pd.rand.mean <- apply(X = pd.rand, MARGIN = 2, FUN = mean, na.rm=TRUE)
    pd.rand.sd <- apply(X = pd.rand, MARGIN = 2, FUN = sd, na.rm=TRUE)
    pd.obs.z <- (pd.obs - pd.rand.mean)/pd.rand.sd
    pd.obs.rank <- apply(X = rbind(pd.obs, pd.rand), MARGIN = 2,
        FUN = rank)[1, ]
    pd.obs.rank <- ifelse(is.na(pd.rand.mean),NA,pd.obs.rank)   
    data.frame(ntaxa=specnumber(samp),pd.obs, pd.rand.mean, pd.rand.sd, pd.obs.rank,
        pd.obs.z, pd.obs.p=pd.obs.rank/(runs+1),runs=runs, row.names = row.names(samp))
}




前期还有一部分命令，主要是想得到pd的sd，大虾也可以参考一下
samp<-read.csv("d:\\127comm110.csv",row.names=1)
phy<-read.tree("d:\\127tree.new")
traits<-read.table("d:\\127traits.csv",header=TRUE,row.names=1,sep=",")
plot(phy)
plot(x=phy,type="phylogram",show.node.label=TRUE, cex=.75)
prunedphy<-prune.sample(samp,phy)
samp<-samp[,prunedphy$tip.label]
traits<-traits[prunedphy$tip.label,]
par(mfrow=c(2,3),cex=.75)   #side by side plots
for(i in row.names(samp))
{
        plot(prunedphy,show.tip.label=TRUE,main=i)
        tiplabels(tip=which(samp[i, ]>0), pch=19,col="blue",cex=.7)
        }
par(mfrow=(c(1,1)))
par(mfrow=c(1,1),cex=1)   #just one plot
plot(prunedphy, show.tip.label=TRUE, show.node.label=TRUE, cex=.7)
tiplabels(tip=which(traits$habit=="tree"),pch=19,cex=1,col="blue")
tiplabels(tip=which(traits$habit=="shrub"),pch=19,cex=1,col="yellow")
tiplabels(tip=which(traits$habit=="vine"),pch=19,cex=1,col="red")
plot(prunedphy, show.tip.label=TRUE,
tip.color=as.numeric(traits$habit),edge.width=3)
phydist<-as.data.frame(cophenetic(prunedphy))


species<-row.names(phydist)
p<-cbind(traits,species,phydist)
pd.result<-pd(samp,prunedphy,include.root=TRUE)
#ses.mpd(samp,phydist,null.model="taxa.labels",abundance.weighted=FALSE,runs=999)

#ses.mntd(samp,phydist,null.model="taxa.labels",abundance.weighted=FALSE,runs=999)
psv.result<-psv(samp,prunedphy,compute.var=TRUE)
psv.result
pd.result
cdnt<-comdistnt(samp,phydist)
cdnt
library(cluster)
cclust<-hclust(cdnt)
plot(cclust)

pd.result<-pd(samp,prunedphy,include.root=TRUE,compute.var=TRUE)
psv.result<-psv(samp,prunedphy,compute.var=TRUE)

library(picante)
   community.txt (264 Bytes)

2012-12-3 11:45:06 上传

comm<-readsample("community.txt")
#create a simple phylogeny
cat("primates(((Saimiri1:2, Saimiri2:2):1, (Cebus1:2,Cebus2:2):1):4,
   ((Alouatta1:2, Alouatta2:2):1,(Ateles1:2, Ateles2:2):1):4);",
   file ="ex.tre", sep = "\n")
tree <- read.tree("ex.tre")
str(tree)
plot(tree)
tree<-prune.sample(comm, tree)

comm<-comm[,tree$tip.label]

par(mfrow = c(2, 2))
plot(tree)
for (i in rownames(comm))
{plot(tree, show.tip.label = FALSE, main = i)
tiplabels(tip = which(comm[i, ] > 0), pch = 19, cex = 2, col ="red")
legend("topleft" , i, bty = "n")}

#Standardized effect size of PD
ses.pd(comm,tree, null.model="taxa.labels", runs=99)

           ntaxa pd.obs pd.rand.mean pd.rand.sd pd.obs.rank   pd.obs.z pd.obs.p runs
community1     4     14     18.19192   2.302070         8.0 -1.8209348    0.080   99
community2     3     17     14.19192   3.835061        81.5  0.7322129    0.815   99
community3     4     19     18.54545   1.785844        57.0  0.2545269    0.570   99

呵呵，终于做出来了，实在实在太感谢了！！！