truequantile( )というRの関数を書いてみた

院生から相談があって、通常連続量とみなして分析する整数のスケールを分位数で区切ってカテゴリ化する場合に（そういう手口の論文は多数ある）、サンプルサイズがある程度大きいと、ちょうど分位数で同順位になる人が複数いて、その分位数を下と上のどちらの範囲カテゴリに入れるかで結果が大きく変わってしまうので、どちらに入れたら良いのか？　という問題を考えてみた。このカテゴリ化の意味は、順位によってほぼ同数のサンプルサイズのサブグループに分けるということだから、分位数ちょうどのサンプルをどちらの範囲カテゴリに含めてもサイズのアンバランスは不可避であって、ベストな解はないと思う（そもそもカテゴリ化によって情報量は落ちるが、分布によっては悪くない）。先行研究があるなら、それに従うべきであろう。中央値より上か下かで２つのサブグループにするだけなら、fmsbパッケージのtruemedian()関数（リンク先pdfのpp.39-41）を使えば、たぶん同順位があっても、同順位の測定値が測定単位より高い精度で測ったら等間隔に分布すると仮定した場合の真の中央値が出るので、下の区間にも上の区間にもその値を含まないことにすれば、一応問題は解決する。例えば、8という値をもつ人が多数いて、それが中央値だった場合、真の中央値が7.5以上8未満であれば8は上の区間に含まれることになり、8以上8.5未満であれば下の区間に含まれることになる（真の中央値もちょうど8になってしまった場合にどうするかという問題は残るが、タイが多いのにそうなってしまう可能性は無視できるくらい小さいと思うので、とりあえず気にしないことにする）。もちろん、等間隔に分布すると仮定すれば、同順位となる分位数ちょうどのデータの一部は下の区間、残りは上の区間に振り分けられるべきだが、実際に同順位となった値のどの個体が下の区間、どの個体が上の区間に振り分けられるべきかを判定する方法が原理的に存在しないので、それは諦めるしかない。

しかしそれさえ我慢すれば、かつて奥村さんが議論されていたように、中央値だけではなく、任意の分位数について、同じ手口が使えるはずである。fmsbパッケージのtruemedian()関数は、等間隔の分布を仮定した真の値を生成してからmedian()をとっているので、returnされるmedian(YY)をquantile(YY, probs=c(0, 1/3, 2/3, 1))に変えるだけで真の三分位数を返すことができる。というわけで、truequantile()関数を定義してみた。例えば、 truequantile(sample(1:10, 100, rep=TRUE), probs=0:3/3)のように使うことができる。

truequantile <- function(X, h=1, probs=0:4/4, ...) { 
   X <- subset(X, !is.na(X))
   YY <- rep(0, length(X))
   XX <- table(X)
   q <- length(XX)
   k <- 0
   for (i in 1:q) {
      L <- as.numeric(names(XX)[i]) - h/2
      for (j in 1:XX[[i]]) {
         k <- k+1
         YY[k] <- L + h*(2*j - 1)/(2*XX[[i]])
      }
   }
   return (quantile(YY, probs=probs, ...))
}