画一个sloan图~

 1#随便找了一个OTU表
 2otu = read.table(file = "otu.txt",sep="\t",row.names=1,header=T)
 3
 4library(MicEco)
 5library(ggplot2)
 6res = neutral.fit(t(otu))
 7
 8##########画图
 9m = res[[1]][1]
10r2 = res[[1]][3]
11out = res[[2]]

1p1 = ggplot() +
2    geom_line(data = out,aes(x=log(p),y=freq.pred),size = 1.2,linetype = 1)+
3    geom_line(data = out,aes(x=log(p),y=Lower),size = 1.2,linetype = 2)+
4    geom_line(data = out,aes(x=log(p),y=Upper),size = 1.2,linetype = 2)+
5    xlab("log10(mean relative abundance)")+ylab("Occurrence frequency")
6p1

 1#判断点的位置，分成三组，红-绿-黄
 2require(plyr)
 3out2 = mutate(out, group=NA)
 4out2$group[out[,2]<out[,4]]="yellow" ##低于下界
 5out2$group[out[,2]>out[,5]]="red" ##高于上界
 6out2$group[(out[,2]>=out[,4])&(out[,2]<=out[,5])]="green"##中间
 7
 8mycols<-c("green","red1","orange")
 9p2 = p1 + geom_point(data = out2,aes(x=log(p),y=freq,color = group),size = 2)+
10         scale_colour_manual(values = mycols)+
11      annotate("text",x=-12.5,y=0.25,label=paste("m = ",round(m,3),sep=''),size=7)+
12      annotate("text",x=-12.5,y=0.3,label=paste("R2 = ",round(r2,3),sep=''),size=7)
13p2

 1plot_theme = theme(panel.background=element_blank(),
 2                   panel.grid=element_blank(),
 3                   axis.line.x=element_line(size=.5, colour="black"),
 4                   axis.line.y=element_line(size=.5, colour="black"),
 5                   axis.ticks=element_line(color="black"),
 6                   axis.text=element_text(color="black", size=24),
 7                   legend.position="none",   ##right
 8                   legend.background=element_blank(),
 9                   legend.key=element_blank(),
10                   legend.text= element_text(size=24),
11                   text=element_text(family="sans", size=24)
12)
13p3 = p2 + plot_theme;p3

 1low = nrow(out2[out2[,6]== "yellow",])
 2med = nrow(out2[out2[,6]== "green",])
 3high = nrow(out2[out2[,6]== "red",])
 4
 5type <- c('med','high','low')
 6nums <- c(med,high,low)
 7df <- data.frame(type = type, nums = nums)
 8label_value <- paste('', round(df$nums/sum(df$nums) * 100, 1), '%', sep = '')
 9label_value
10label <- paste(df$type, label_value, sep = ' ')
11label
12p4 <- ggplot(data = df, aes(x = 1, y = nums, fill = type)) + 
13  geom_bar(stat = 'identity', position = 'stack', width = 0.5)+
14  scale_fill_manual(name='', 
15                    labels=c(label[2], label[3], label[1]),
16                    values=c("red1","orange","green"))+
17  coord_polar(theta = 'y')+
18  labs(x = '', y = '', title = '')+
19  theme(axis.text = element_blank())+
20  theme(axis.ticks = element_blank())+
21  theme(legend.position = "right")
22 
24p4
25p4 = p4+theme(panel.background=element_blank(),
26              panel.grid=element_blank(),
27              legend.background=element_blank(),
28              legend.key=element_blank(),
29              legend.text= element_text(size=24))
30
31p4

1require(ggimage)
2g <- p3 + geom_subview(subview = p4 + theme_void(), x=-11.5, y=0.75, w=5, h=5)
3g
4#dev.off()