Triad sou.

やっつけ jitter 関数 (箱ひげ図用)

ある本を読んでいたら、箱ひげ図に jitter を付加した散布図を重ねて描いていたので、再現してみようと思い適当に作りました。
散布図のプロットの際に、

  • 近くに点があり重なりそうな場合には、jitter を強めに付加する
  • 特に tie になる場合はちゃんと重ならないように処理する
  • 近くに点がない場合は、jitter があまり付加されない (軸に近づく) ようにする

という方針で作っています。

myBoxplotJitter <- function(x, range, overlap=1, lambda=5) {
  index <- seq(1, length(x))
  jitter <- rep(0, length(x))
  tie <- rep(FALSE, length(x))
  df <- data.frame(x, jitter, index, tie)
  df <- df[order(df$x),]
  for (i in 1:length(x)) {
    region <- c(seq(i - overlap, i), seq(i + 1, i + overlap))
    region <- region[region>0 & region<=length(x)]
    df$jitter[i] <- abs(
      sqrt(sum((df$x[region] - df$x[i])^2)) / (length(region) - 1)
    )
    df$tie[i] <- (length(table(df$x[region])) < length(region))
  }
  attenuation <- exp(-lambda*df$jitter) / exp(-lambda*min(df$jitter))
  attenuation <- ifelse(df$tie, 1, attenuation)
  df$jitter <- runif(length(x), min=-abs(range), max=abs(range)) * attenuation
  df <- df[order(df$index),]
  df$jitter
}

set.seed(1192)

x <- c(
  runif(95), c(2,3,4,5,6),
  runif(95), c(-1,-1,-1,-1,-2)
)
# as.numeric(factor(c(rep(1, 100), rep(2, 100))))
g <- c(rep(1, 100), rep(2, 100))

jitter <- unlist(by(x, g, function(x) myBoxplotJitter(x, 0.3, lambda=50)))

xlim <- c(1-0.5, 2+0.5)
ylim <- c(min(x)-0.1, max(x)+0.1)


postscript(
  file="jitter.eps", width = 18, height = 18, pointsize = 40,
  horizontal = FALSE, onefile = FALSE, paper = "special",
  family="Times"
)

par(mar=c(3.2, 3.2, 1.5, 0.5), mgp=c(2.0, 0.7, 0))

boxplot(x ~ g, xlim=xlim, ylim=ylim, outline=FALSE)
par(new=T)
plot(
  g + jitter, x,
  xlab="group + jitter", ylab=expression(italic(x)), xaxt = "n", yaxt = "n",
  pch = 20, cex = 1.25,
  xlim=xlim, ylim=ylim,
  main="myJitter"
)

dev.off()


使い方
myBoxplotJitter(x, range, overlap=1, lambda=5)

引数
x 	 データ (numeric 型のベクトル)
range    jitter を付加する範囲、±range の幅の間の値が返ります (numeric)
overlap  前後何点を考慮するか (integer)
lambda   1 以上の値を取る重み係数、値が大きいと軸に近づきやすい (numeric)