R 関数を作ってみる (1) 正規分布 塗りつぶし

graph.norm <- function(x1,x2,m,s,oneway=FALSE,rev=FALSE) {
  sig <- sqrt(s)
  xmin <- m-4*sig
  xmax <- m+4*sig
  dx <- (xmax-xmin)/1000
  x <- seq(xmin,xmax,by=dx)
  y <- dnorm(x,mean=m,sd=sig)
  xl <- c(xmin,xmax)
  yl <- c(0,max(y)*1.05)

  title = paste("Normal Distribution N(",m,",",s,")",sep="")
  
  plot(x,y,type="l",xlim=xl,ylim=yl,xlab="x",ylab="y",main=title)
  if (oneway==TRUE) {
    if (rev==FALSE) {
      pp1 <- 1-pnorm(x1,mean=m,sd=sig)
      px1 <- round(pp1,digits=4)
      y1 <- 0; y2 <- 0; x3 <- xmax; x4 <- x1
      y3 <- dnorm(x3,mean=m,sd=sig)
      y4 <- dnorm(x4,mean=m,sd=sig)
      xx <- rev(seq(x4,x3,by=(x3-x4)/1000))
      yy <- dnorm(xx,mean=m,sd=sig)
      xx <- c(x1,x1,x3,xx,x4); yy <- c(y1,y2,y3,yy,y4)
      polygon(xx,yy,col=gray(0.5))
      x0 <- (x1+xmax)/2; y0 <- max(yy)/2; text(x0,y0,px1)      
    } else {
      pp1 <- pnorm(x1,mean=m,sd=sig)
      px1 <- round(pp1,digits=4)
      y1 <- 0; y2 <- 0; x3 <- x1; x4 <- xmin
      y3 <- dnorm(x3,mean=m,sd=sig)
      y4 <- dnorm(x4,mean=m,sd=sig)
      xx <- rev(seq(x4,x3,by=(x3-x4)/1000))
      yy <- dnorm(xx,mean=m,sd=sig)
      xx <- c(x4,x3,x3,xx,x4); yy <- c(y1,y2,y3,yy,y4)
      polygon(xx,yy,col=gray(0.5))
      x0 <- (x1+xmin)/2; y0 <- max(yy)/2; text(x0,y0,px1)
    }

    } else {
      pp2 <- pnorm(x2,mean=m,sd=sig)-pnorm(x1,mean=m,sd=sig)
      px2 <- round(pp2,digits=4)
      y1 <- 0; y2 <- 0; x3 <- x2; x4 <- x1
      y3 <- dnorm(x3,mean=m,sd=sig)
      y4 <- dnorm(x4,mean=m,sd=sig)
      xx <- rev(seq(x4,x3,by=(x3-x4)/1000))
      yy <- dnorm(xx,mean=m,sd=sig)
      xx <- c(x1,x2,x3,xx,x4); yy <- c(y1,y2,y3,yy,y4)
      polygon(xx,yy,col=gray(0.5))
      x0 <- (x1+x2)/2; y0 <- max(yy)/2; text(x0,y0,px2)
    }
  }

#使ってみる
par(mfrow=c(3,3))
graph.norm(7,,10,4,oneway=TRUE,rev=TRUE)　     #左上
graph.norm(7,,10,4,oneway=TRUE,rev=FALSE)　　　#中央上
graph.norm(7,11,10,4,oneway=FALSE,rev=FALSE)　#中央右

graph.norm(9,,10,4,oneway=TRUE,rev=TRUE)　　　#左中央
graph.norm(9,,10,4,oneway=TRUE,rev=FALSE)　　 #中央
graph.norm(9,12,10,4,oneway=FALSE,rev=FALSE)　#中央右

graph.norm(11,,10,4,oneway=TRUE,rev=TRUE)　　　#左下
graph.norm(11,,10,4,oneway=TRUE,rev=FALSE)　　 #中央下
graph.norm(11,13,10,4,oneway=FALSE,rev=FALSE)　#右下

正規分布のグラフを描き、区間を指定して塗った上で、領域の確率を返すことができる