to be continued...

Rや心理統計学の備忘録的な

There are three kinds of lies: lies, damned lies, and statistics.

- 嘘には三種類ある。嘘、大嘘、そして統計だ -
Benjamn Disraeli(19世紀のイギリス首相)

Rで分位分け関数

Rで分位分け関数

どうも、最近は「女子高生の無駄遣い」見ています。

「輪! moon! dass! cry!」が頭から離れない。

今回は、ショート記事。

Rで分位分け関数を作ります

(分位分けという表現あってるかわからない)

というのも、機械学習や統計分析をしてるとデータの10のランクに分けたりして、データの上位40%はどうなのかを見たり、中央についてみたりします。

ググったけどなかったので、個人的なメモとして残しておきます。

(余裕でdplyrにありました)


使う関数は

  • ceiling():x未満ではない最小の整数を返す

    例:0→0, 1→1, 1.001→2, 2.999→3

  • rank():データの大きさの順番を返す関数(順序尺度変換)

  • length():ベクトルの要素数を返す関数

そして関数がこれ

sep_tile <- function(x,tile = 10){
  res <- ceiling(rank(x)/length(x)*tile)
  return(res)
}

ほんまに使えるの?


データ生成

> set.seed(42)
> umr <- runif(100,0,1)
> head(umr)
[1] 0.9148060 0.9370754 0.2861395 0.8304476 0.6417455 0.5190959

シードを固定して、[0,1]の一様乱数を100個生成したものをumrに格納

仮想データに自作関数を適用してみると

> sort(sep_tile(umr))
  [1]  1  1  1  1  1  1  1  1  1  1  2  2  2  2  2  2  2  2  2  2  3  3  3  3  3  3  3  3  3  3  4  4  4  4  4  4  4  4  4
 [40]  4  5  5  5  5  5  5  5  5  5  5  6  6  6  6  6  6  6  6  6  6  7  7  7  7  7  7  7  7  7  7  8  8  8  8  8  8  8  8
 [79]  8  8  9  9  9  9  9  9  9  9  9  9 10 10 10 10 10 10 10 10 10 10

わー、どうやら出来てるっぽい。


関数の説明

tileの引数を10として、10分位に分けます

  • rank()でデータを順序尺度に変換します

    1~100のランクが与えられます

rank(umr)
  [1]  88  92  28  80  59  50  71  13  61  67  42  68  91  25  43  93  98  12  44  51  85  14 100  94   9  46  37  86  41
 [30]  82  72  78  36  65   3  81   4  20  87  55  35  40   6  97  39  95  84  58  96  56  31  32  38  77   7  74  64  17
 [59]  26  48  63  99  75  53  83  18  27  79  66  24   8  15  22  45  19  69   5  34  47   2  54  16  33  60  76  52  23
 [88]  11  10  29  62   1  21  90  89  70  30  49  73  57
  • length()で割ります

    0~1の値に変換されます

rank(umr)/length(umr)
  [1] 0.88 0.92 0.28 0.80 0.59 0.50 0.71 0.13 0.61 0.67 0.42 0.68 0.91 0.25 0.43 0.93 0.98 0.12 0.44 0.51 0.85 0.14 1.00
 [24] 0.94 0.09 0.46 0.37 0.86 0.41 0.82 0.72 0.78 0.36 0.65 0.03 0.81 0.04 0.20 0.87 0.55 0.35 0.40 0.06 0.97 0.39 0.95
 [47] 0.84 0.58 0.96 0.56 0.31 0.32 0.38 0.77 0.07 0.74 0.64 0.17 0.26 0.48 0.63 0.99 0.75 0.53 0.83 0.18 0.27 0.79 0.66
 [70] 0.24 0.08 0.15 0.22 0.45 0.19 0.69 0.05 0.34 0.47 0.02 0.54 0.16 0.33 0.60 0.76 0.52 0.23 0.11 0.10 0.29 0.62 0.01
 [93] 0.21 0.90 0.89 0.70 0.30 0.49 0.73 0.57
  • 10を掛けます

    ここで掛ける値が分けたい分位数になります(0~10の値に変換)

    今回は10分位に分けました

rank(umr)/length(umr)*10
  [1]  8.8  9.2  2.8  8.0  5.9  5.0  7.1  1.3  6.1  6.7  4.2  6.8  9.1  2.5  4.3  9.3  9.8  1.2  4.4  5.1  8.5  1.4 10.0
 [24]  9.4  0.9  4.6  3.7  8.6  4.1  8.2  7.2  7.8  3.6  6.5  0.3  8.1  0.4  2.0  8.7  5.5  3.5  4.0  0.6  9.7  3.9  9.5
 [47]  8.4  5.8  9.6  5.6  3.1  3.2  3.8  7.7  0.7  7.4  6.4  1.7  2.6  4.8  6.3  9.9  7.5  5.3  8.3  1.8  2.7  7.9  6.6
 [70]  2.4  0.8  1.5  2.2  4.5  1.9  6.9  0.5  3.4  4.7  0.2  5.4  1.6  3.3  6.0  7.6  5.2  2.3  1.1  1.0  2.9  6.2  0.1
 [93]  2.1  9.0  8.9  7.0  3.0  4.9  7.3  5.7
  • ceiling()でx未満ではない最小の整数を返します

    こんな感じ

> ceiling(0)
[1] 0
> ceiling(0.1)
[1] 1
> ceiling(1.03)
[1] 2

これで完成

ceiling(rank(umr)/length(umr)*10)
  [1]  9 10  3  8  6  5  8  2  7  7  5  7 10  3  5 10 10  2  5  6  9  2 10 10  1  5  4  9  5  9  8  8  4  7  1  9  1  2  9
 [40]  6  4  4  1 10  4 10  9  6 10  6  4  4  4  8  1  8  7  2  3  5  7 10  8  6  9  2  3  8  7  3  1  2  3  5  2  7  1  4
 [79]  5  1  6  2  4  6  8  6  3  2  1  3  7  1  3  9  9  7  3  5  8  6

ほんとに分位分け出来てるの?

> ceiling(rank(umr)/length(umr)*10)
  [1]  9 10  3  8  6  5  8  2  7  7  5  7 10  3  5 10 10  2  5  6  9  2 10 10  1  5  4  9  5  9  8  8  4  7  1  9  1  2  9
 [40]  6  4  4  1 10  4 10  9  6 10  6  4  4  4  8  1  8  7  2  3  5  7 10  8  6  9  2  3  8  7  3  1  2  3  5  2  7  1  4
 [79]  5  1  6  2  4  6  8  6  3  2  1  3  7  1  3  9  9  7  3  5  8  6
> sort(sep_tile(umr,tile = 5))
  [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
 [59] 3 3 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5
> sort(sep_tile(umr,tile = 3))
  [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
 [59] 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
> sort(sep_tile(umr,tile = 100))
  [1]   1   2   3   4   5   6   8   8   9  10  11  12  13  15  15  16  17  18  19  20  21  22  23  24  25  26  27  29  29
 [30]  30  31  32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54  56  57  57  58
 [59]  59  60  61  62  63  64  65  66  67  68  69  70  71  72  73  74  75  76  77  78  79  80  81  82  83  84  85  86  87
 [88]  88  89  90  91  92  93  94  95  96  97  98  99 100
> sort(sep_tile(umr,tile = 20))
  [1]  1  1  1  1  1  2  2  2  2  2  3  3  3  3  3  4  4  4  4  4  5  5  5  5  5  6  6  6  6  6  7  7  7  7  7  8  8  8  8
 [40]  8  9  9  9  9  9 10 10 10 10 10 11 11 11 11 11 12 12 12 12 12 13 13 13 13 13 14 14 14 14 14 15 15 15 15 15 16 16 16
 [79] 16 16 17 17 17 17 17 18 18 18 18 18 19 19 19 19 19 20 20 20 20 20

こげな感じです。

今回のコードまとめ

#シードを固定
set.seed(42)
#[0,1]の一様乱数を100個→umrに格納
umr <- runif(100,0,1)
#10分位に分ける
ceiling(rank(umr)/length(umr)*10)
#分位分け自作関数
sep_tile <- function(x,tile = 10){
  res <- ceiling(rank(x)/length(x)*tile)
  return(res)
}
#関数使ってみた
sep_tile(umr)

わーい。


{dplyr}に搭載されていた...

> library(tidyverse)
> ntile(umr,10)
  [1]  9 10  3  8  6  5  8  2  7  7  5  7 10  3  5 10 10  2  5  6  9  2 10 10  1  5  4  9  5  9  8  8  4  7  1  9  1  2  9  6  4  4  1 10  4 10  9  6
 [49] 10  6  4  4  4  8  1  8  7  2  3  5  7 10  8  6  9  2  3  8  7  3  1  2  3  5  2  7  1  4  5  1  6  2  4  6  8  6  3  2  1  3  7  1  3  9  9  7
 [97]  3  5  8  
to be continued...