R中查找曲线的峰值

目的：计算密度曲线对应的峰值并作图。

library(ggplot2)
library(reshape2)

mytheme2 <- theme_bw() + theme(legend.title=element_blank(), panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
                               axis.text.x=element_text(size=6,angle=0), legend.position="top")

argv<-commandArgs(TRUE)
data <- read.table(argv[1], header = TRUE, sep="\t")

densFindPeak <- function(x){
 td <- density(x)
 maxDens <- which.max(td$y)
 list(x=td$x[maxDens],y=td$y[maxDens])
}


peak =  densFindPeak(data$depth)

peak

p1 <- ggplot(data,aes(depth, alpha=0.8)) +
      geom_density() + 
      geom_rug() +
      geom_vline(xintercept = peak[[1]],colour="red",linetype="dashed") +
      geom_text(aes(x=peak[[1]], y= peak[[2]]+0.3,label=paste("peak =",round(peak[[1]],4))), color="black") +
      ylab("Density") + 
      mytheme2

ggsave(file=argv[2], plot=p1, width = 5, height = 4)