ポケモンGO バンギラスの厳選はどこで妥協するか。Rで計算してみる
#注意：2018.10.10 近い将来CPなどのデータが調整されるとのことです。
#CP値などのデータは更新前の2018.10.09現在の値なのでご注意ください。
ポケモンgoで岩祭り開催中。
ヨーギラスがよく出るので、厳選したいが、どの程度で妥協すればよいだろうか。
バンギラスの各攻撃、防御、HP(スタミナ）毎での、最大CPやレア度を知りたい。
Rで計算してみる。
#順番をCP強い方が先頭に来るように改定。CP順位も追加。同値の扱いを変更。
#野生を何匹以上捕まえるか、あるいは何個以上卵を孵化させれば80%あるいは
90%以上の確率でそのCP以上のヨーギラスが得られるかも追加。
固体値などは
http://www.pokemontrash.com/pokemon-go/niveau-pokemon-poussieres/
などを参考にした。
とりあえず、全個体値の分と、卵に限った場合で計算した。
CP3580超えを卵で孵化させるには
AT DF HP kotaichi_per Bangiras_max_cp rank joui_per N80per N90per 25 14 12 15 91.11111 3580 25 11.574074 14 19
14体ヨーギラスが孵化すれば80％以上の確率で
19体ヨーギラスが孵化すれば90％以上の確率で
3580超えの固体値のヨーギラスが得られる。
これまで3体しか生まれていないのでまだ先が長いが
そのうちなんとかなるか？
CP3600超えとなると卵産でも
AT DF HP kotaichi_per Bangiras_max_cp rank joui_per N80per N90per 7 15 15 13 95.55556 3600 7 3.240741 49 70
上位3.2％とかなり厳しいですな。8割の確率でゲットするには卵49個も孵化せねばならぬ。
野生ならCP3600超は
AT DF HP kotaichi_per Bangiras_max_cp rank joui_per N80per N90per 7 15 15 13 95.55556 3600 7 0.17089844 941 1347
上位0.17%で無理ゲー。8割の確率でゲットするにはヨーギラス941匹もゲットせねばならぬ。
以下バンギラスの各攻撃、防御、HP(スタミナ）毎の計算値を記す。
ヨーギラスの時点から固体値は固定なので、ヨーギラス厳選の参考にしてください。
＃2017.5.21現在のネット上の情報で作成。間違っていても責任は持ちません。
それにしても、この表を見ると攻撃力大事なのね。
＃全個体値（卵に限らず）
攻撃、防御、HP(スタミナ）、固体値（％）、バンギラスLv39でのmaxCP、maxCP順位、maxCP上位何％でゲットできるか、何匹野生のヨーギラスを捕獲すれば80%以上の確率でそのmaxCP以上のヨーギラスが得られるか、同９０％以上
AT DF HP kotaichi_per Bangiras_max_cp rank joui_per N80per N90per 1 15 15 15 100.00000 3617 1 0.02441406 6592 9431 2 15 15 14 97.77778 3609 2 0.07324219 2197 3143 3 15 14 15 97.77778 3609 2 0.07324219 2197 3143 4 14 15 15 97.77778 3604 4 0.09765625 1648 2357 5 15 14 14 95.55556 3601 5 0.14648438 1098 1571 6 15 13 15 95.55556 3601 5 0.14648438 1098 1571 7 15 15 13 95.55556 3600 7 0.17089844 941 1347 8 14 14 15 95.55556 3596 8 0.19531250 824 1178 9 14 15 14 95.55556 3595 9 0.21972656 732 1047 10 15 14 13 93.33333 3593 10 0.29296875 549 785 11 15 13 14 93.33333 3593 10 0.29296875 549 785 12 15 12 15 93.33333 3593 10 0.29296875 549 785 13 15 15 12 93.33333 3592 13 0.31738281 507 725 14 13 15 15 95.55556 3590 14 0.34179688 471 673 15 14 13 15 93.33333 3588 15 0.36621094 439 628 16 14 15 13 93.33333 3587 16 0.41503906 387 554 17 14 14 14 93.33333 3587 16 0.41503906 387 554 18 15 13 13 91.11111 3585 18 0.48828125 329 471 19 15 12 14 91.11111 3585 18 0.48828125 329 471 20 15 11 15 91.11111 3585 18 0.48828125 329 471 21 15 15 11 91.11111 3584 21 0.53710938 299 428 22 15 14 12 91.11111 3584 21 0.53710938 299 428 23 13 15 14 93.33333 3582 23 0.58593750 274 392 24 13 14 15 93.33333 3582 23 0.58593750 274 392 25 14 12 15 91.11111 3580 25 0.61035156 263 377 26 14 15 12 91.11111 3579 26 0.68359375 235 336 27 14 14 13 91.11111 3579 26 0.68359375 235 336 28 14 13 14 91.11111 3579 26 0.68359375 235 336 29 15 12 13 88.88889 3577 29 0.78125000 206 294 30 15 11 14 88.88889 3577 29 0.78125000 206 294 31 15 10 15 88.88889 3577 29 0.78125000 206 294 32 12 15 15 93.33333 3577 29 0.78125000 206 294 33 15 14 11 88.88889 3576 33 0.83007812 194 277 34 15 13 12 88.88889 3576 33 0.83007812 194 277 35 15 15 10 88.88889 3575 35 0.85449219 188 269 36 13 14 14 91.11111 3574 36 0.90332031 178 254 37 13 13 15 91.11111 3574 36 0.90332031 178 254 38 13 15 13 91.11111 3573 38 0.92773438 173 248 39 14 12 14 88.88889 3572 39 0.97656250 165 235 40 14 11 15 88.88889 3572 39 0.97656250 165 235 41 14 14 12 88.88889 3571 41 1.02539062 157 224 42 14 13 13 88.88889 3571 41 1.02539062 157 224 43 14 15 11 88.88889 3570 43 1.04980469 153 219 44 15 11 13 86.66667 3569 44 1.14746094 140 200 45 15 10 14 86.66667 3569 44 1.14746094 140 200 46 15 9 15 86.66667 3569 44 1.14746094 140 200 47 12 14 15 91.11111 3569 44 1.14746094 140 200 48 15 13 11 86.66667 3568 48 1.22070312 132 188 49 15 12 12 86.66667 3568 48 1.22070312 132 188 50 12 15 14 91.11111 3568 48 1.22070312 132 188 51 15 15 9 86.66667 3567 51 1.26953125 126 181 52 15 14 10 86.66667 3567 51 1.26953125 126 181 53 13 14 13 88.88889 3566 53 1.34277344 120 171 54 13 13 14 88.88889 3566 53 1.34277344 120 171 55 13 12 15 88.88889 3566 53 1.34277344 120 171 56 13 15 12 88.88889 3565 56 1.36718750 117 168 57 14 11 14 86.66667 3564 57 1.41601562 113 162 58 14 10 15 86.66667 3564 57 1.41601562 113 162 59 14 13 12 86.66667 3563 59 1.48925781 108 154 60 14 12 13 86.66667 3563 59 1.48925781 108 154 61 11 15 15 91.11111 3563 59 1.48925781 108 154 62 14 15 10 86.66667 3562 62 1.53808594 104 149 63 14 14 11 86.66667 3562 62 1.53808594 104 149 64 15 10 13 84.44444 3561 64 1.63574219 98 140 65 15 9 14 84.44444 3561 64 1.63574219 98 140 66 15 8 15 84.44444 3561 64 1.63574219 98 140 67 12 13 15 88.88889 3561 64 1.63574219 98 140 68 15 12 11 84.44444 3560 68 1.73339844 93 132 69 15 11 12 84.44444 3560 68 1.73339844 93 132 70 12 15 13 88.88889 3560 68 1.73339844 93 132 71 12 14 14 88.88889 3560 68 1.73339844 93 132 72 15 14 9 84.44444 3559 72 1.78222656 90 129 73 15 13 10 84.44444 3559 72 1.78222656 90 129 74 15 15 8 84.44444 3558 74 1.87988281 85 122 75 13 13 13 86.66667 3558 74 1.87988281 85 122 76 13 12 14 86.66667 3558 74 1.87988281 85 122 77 13 11 15 86.66667 3558 74 1.87988281 85 122 78 13 15 11 86.66667 3557 78 1.92871094 83 119 79 13 14 12 86.66667 3557 78 1.92871094 83 119 80 14 10 14 84.44444 3556 80 1.97753906 81 116 81 14 9 15 84.44444 3556 80 1.97753906 81 116 82 14 12 12 84.44444 3555 82 2.07519531 77 110 83 14 11 13 84.44444 3555 82 2.07519531 77 110 84 11 15 14 88.88889 3555 82 2.07519531 77 110 85 11 14 15 88.88889 3555 82 2.07519531 77 110 86 14 14 10 84.44444 3554 86 2.12402344 75 108 87 14 13 11 84.44444 3554 86 2.12402344 75 108 88 14 15 9 84.44444 3553 88 2.24609375 71 102 89 15 9 13 82.22222 3553 88 2.24609375 71 102 90 15 8 14 82.22222 3553 88 2.24609375 71 102 91 15 7 15 82.22222 3553 88 2.24609375 71 102 92 12 12 15 86.66667 3553 88 2.24609375 71 102 93 15 11 11 82.22222 3552 93 2.36816406 68 97 94 15 10 12 82.22222 3552 93 2.36816406 68 97 95 12 15 12 86.66667 3552 93 2.36816406 68 97 96 12 14 13 86.66667 3552 93 2.36816406 68 97 97 12 13 14 86.66667 3552 93 2.36816406 68 97 98 15 13 9 82.22222 3551 98 2.41699219 66 95 99 15 12 10 82.22222 3551 98 2.41699219 66 95 100 15 14 8 82.22222 3550 100 2.51464844 64 91 101 13 12 13 84.44444 3550 100 2.51464844 64 91 102 13 11 14 84.44444 3550 100 2.51464844 64 91 103 13 10 15 84.44444 3550 100 2.51464844 64 91 104 15 15 7 82.22222 3549 104 2.61230469 61 87 105 13 14 11 84.44444 3549 104 2.61230469 61 87 106 13 13 12 84.44444 3549 104 2.61230469 61 87 107 10 15 15 88.88889 3549 104 2.61230469 61 87 108 13 15 10 84.44444 3548 108 2.68554688 60 85 109 14 9 14 82.22222 3548 108 2.68554688 60 85 110 14 8 15 82.22222 3548 108 2.68554688 60 85 111 14 11 12 82.22222 3547 111 2.78320312 58 82 112 14 10 13 82.22222 3547 111 2.78320312 58 82 113 11 14 14 86.66667 3547 111 2.78320312 58 82 114 11 13 15 86.66667 3547 111 2.78320312 58 82 115 14 13 10 82.22222 3546 115 2.85644531 56 80 116 14 12 11 82.22222 3546 115 2.85644531 56 80 117 11 15 13 86.66667 3546 115 2.85644531 56 80 118 14 15 8 82.22222 3545 118 3.02734375 53 75 119 14 14 9 82.22222 3545 118 3.02734375 53 75 120 15 8 13 80.00000 3545 118 3.02734375 53 75 121 15 7 14 80.00000 3545 118 3.02734375 53 75 122 12 12 14 84.44444 3545 118 3.02734375 53 75 123 15 6 15 80.00000 3545 118 3.02734375 53 75 124 12 11 15 84.44444 3545 118 3.02734375 53 75 125 15 10 11 80.00000 3544 125 3.12500000 51 73 126 15 9 12 80.00000 3544 125 3.12500000 51 73 127 12 14 12 84.44444 3544 125 3.12500000 51 73 128 12 13 13 84.44444 3544 125 3.12500000 51 73 129 15 12 9 80.00000 3543 129 3.19824219 50 71 130 15 11 10 80.00000 3543 129 3.19824219 50 71 131 12 15 11 84.44444 3543 129 3.19824219 50 71 132 15 14 7 80.00000 3542 132 3.34472656 48 68 133 15 13 8 80.00000 3542 132 3.34472656 48 68 134 13 11 13 82.22222 3542 132 3.34472656 48 68 135 13 10 14 82.22222 3542 132 3.34472656 48 68 136 13 9 15 82.22222 3542 132 3.34472656 48 68 137 10 14 15 86.66667 3542 132 3.34472656 48 68 138 15 15 6 80.00000 3541 138 3.44238281 46 66 139 13 13 11 82.22222 3541 138 3.44238281 46 66 140 13 12 12 82.22222 3541 138 3.44238281 46 66 141 10 15 14 86.66667 3541 138 3.44238281 46 66 142 13 15 9 82.22222 3540 142 3.51562500 45 65 143 13 14 10 82.22222 3540 142 3.51562500 45 65 144 14 7 15 80.00000 3540 142 3.51562500 45 65 145 14 10 12 80.00000 3539 145 3.66210938 44 62 146 14 9 13 80.00000 3539 145 3.66210938 44 62 147 11 14 13 84.44444 3539 145 3.66210938 44 62 148 14 8 14 80.00000 3539 145 3.66210938 44 62 149 11 13 14 84.44444 3539 145 3.66210938 44 62 150 11 12 15 84.44444 3539 145 3.66210938 44 62 151 14 12 10 80.00000 3538 151 3.73535156 43 61 152 14 11 11 80.00000 3538 151 3.73535156 43 61 153 11 15 12 84.44444 3538 151 3.73535156 43 61 154 14 14 8 80.00000 3537 154 3.88183594 41 59 155 14 13 9 80.00000 3537 154 3.88183594 41 59 156 15 6 14 77.77778 3537 154 3.88183594 41 59 157 12 11 14 82.22222 3537 154 3.88183594 41 59 158 15 5 15 77.77778 3537 154 3.88183594 41 59 159 12 10 15 82.22222 3537 154 3.88183594 41 59 160 14 15 7 80.00000 3536 160 4.05273438 39 56 161 15 9 11 77.77778 3536 160 4.05273438 39 56 162 15 8 12 77.77778 3536 160 4.05273438 39 56 163 12 13 12 82.22222 3536 160 4.05273438 39 56 164 15 7 13 77.77778 3536 160 4.05273438 39 56 165 12 12 13 82.22222 3536 160 4.05273438 39 56 166 9 15 15 86.66667 3536 160 4.05273438 39 56 167 15 11 9 77.77778 3535 167 4.15039062 38 55 168 15 10 10 77.77778 3535 167 4.15039062 38 55 169 12 15 10 82.22222 3535 167 4.15039062 38 55 170 12 14 11 82.22222 3535 167 4.15039062 38 55 171 15 13 7 77.77778 3534 171 4.29687500 37 53 172 15 12 8 77.77778 3534 171 4.29687500 37 53 173 13 10 13 80.00000 3534 171 4.29687500 37 53 174 13 9 14 80.00000 3534 171 4.29687500 37 53 175 13 8 15 80.00000 3534 171 4.29687500 37 53 176 10 13 15 84.44444 3534 171 4.29687500 37 53 177 15 14 6 77.77778 3533 177 4.41894531 36 51 178 13 12 11 80.00000 3533 177 4.41894531 36 51 179 13 11 12 80.00000 3533 177 4.41894531 36 51 180 10 15 13 84.44444 3533 177 4.41894531 36 51 181 10 14 14 84.44444 3533 177 4.41894531 36 51 182 15 15 5 77.77778 3532 182 4.51660156 35 50 183 13 14 9 80.00000 3532 182 4.51660156 35 50 184 13 13 10 80.00000 3532 182 4.51660156 35 50 185 14 6 15 77.77778 3532 182 4.51660156 35 50 186 13 15 8 80.00000 3531 186 4.71191406 34 48 187 14 10 11 77.77778 3531 186 4.71191406 34 48 188 14 9 12 77.77778 3531 186 4.71191406 34 48 189 14 8 13 77.77778 3531 186 4.71191406 34 48 190 11 13 13 82.22222 3531 186 4.71191406 34 48 191 14 7 14 77.77778 3531 186 4.71191406 34 48 192 11 12 14 82.22222 3531 186 4.71191406 34 48 193 11 11 15 82.22222 3531 186 4.71191406 34 48 194 14 12 9 77.77778 3530 194 4.80957031 33 47 195 14 11 10 77.77778 3530 194 4.80957031 33 47 196 11 15 11 82.22222 3530 194 4.80957031 33 47 197 11 14 12 82.22222 3530 194 4.80957031 33 47 198 14 13 8 77.77778 3529 198 4.93164062 32 46 199 15 5 14 75.55556 3529 198 4.93164062 32 46 200 12 10 14 80.00000 3529 198 4.93164062 32 46
＃卵産のみ
攻撃、防御、HP(スタミナ）、固体値（％）、バンギラスLv39でのmaxCP、maxCP順位、maxCP上位何％でゲットできるか、何個卵を孵化させれば80%以上の確率でそのmaxCP以上のヨーギラスが得られるか、同９０％以上
AT DF HP kotaichi_per Bangiras_max_cp rank joui_per N80per N90per 1 15 15 15 100.00000 3617 1 0.462963 347 497 2 15 15 14 97.77778 3609 2 1.388889 116 165 3 15 14 15 97.77778 3609 2 1.388889 116 165 4 14 15 15 97.77778 3604 4 1.851852 87 124 5 15 14 14 95.55556 3601 5 2.777778 58 82 6 15 13 15 95.55556 3601 5 2.777778 58 82 7 15 15 13 95.55556 3600 7 3.240741 49 70 8 14 14 15 95.55556 3596 8 3.703704 43 62 9 14 15 14 95.55556 3595 9 4.166667 38 55 10 15 14 13 93.33333 3593 10 5.555556 29 41 11 15 13 14 93.33333 3593 10 5.555556 29 41 12 15 12 15 93.33333 3593 10 5.555556 29 41 13 15 15 12 93.33333 3592 13 6.018519 26 38 14 13 15 15 95.55556 3590 14 6.481481 25 35 15 14 13 15 93.33333 3588 15 6.944444 23 32 16 14 15 13 93.33333 3587 16 7.870370 20 29 17 14 14 14 93.33333 3587 16 7.870370 20 29 18 15 13 13 91.11111 3585 18 9.259259 17 24 19 15 12 14 91.11111 3585 18 9.259259 17 24 20 15 11 15 91.11111 3585 18 9.259259 17 24 21 15 15 11 91.11111 3584 21 10.185185 15 22 22 15 14 12 91.11111 3584 21 10.185185 15 22 23 13 15 14 93.33333 3582 23 11.111111 14 20 24 13 14 15 93.33333 3582 23 11.111111 14 20 25 14 12 15 91.11111 3580 25 11.574074 14 19 26 14 15 12 91.11111 3579 26 12.962963 12 17 27 14 14 13 91.11111 3579 26 12.962963 12 17 28 14 13 14 91.11111 3579 26 12.962963 12 17 29 15 12 13 88.88889 3577 29 14.814815 11 15 30 15 11 14 88.88889 3577 29 14.814815 11 15 31 15 10 15 88.88889 3577 29 14.814815 11 15 32 12 15 15 93.33333 3577 29 14.814815 11 15 33 15 14 11 88.88889 3576 33 15.740741 10 14 34 15 13 12 88.88889 3576 33 15.740741 10 14 35 15 15 10 88.88889 3575 35 16.203704 10 14 36 13 14 14 91.11111 3574 36 17.129630 9 13 37 13 13 15 91.11111 3574 36 17.129630 9 13 38 13 15 13 91.11111 3573 38 17.592593 9 12 39 14 12 14 88.88889 3572 39 18.518519 8 12 40 14 11 15 88.88889 3572 39 18.518519 8 12 41 14 14 12 88.88889 3571 41 19.444444 8 11 42 14 13 13 88.88889 3571 41 19.444444 8 11 43 14 15 11 88.88889 3570 43 19.907407 8 11 44 15 11 13 86.66667 3569 44 21.296296 7 10 45 15 10 14 86.66667 3569 44 21.296296 7 10 46 12 14 15 91.11111 3569 44 21.296296 7 10 47 15 13 11 86.66667 3568 47 22.685185 7 9 48 15 12 12 86.66667 3568 47 22.685185 7 9 49 12 15 14 91.11111 3568 47 22.685185 7 9 50 15 14 10 86.66667 3567 50 23.148148 7 9 51 13 14 13 88.88889 3566 51 24.537037 6 9 52 13 13 14 88.88889 3566 51 24.537037 6 9 53 13 12 15 88.88889 3566 51 24.537037 6 9 54 13 15 12 88.88889 3565 54 25.000000 6 9 55 14 11 14 86.66667 3564 55 25.925926 6 8 56 14 10 15 86.66667 3564 55 25.925926 6 8 57 14 13 12 86.66667 3563 57 27.314815 6 8 58 14 12 13 86.66667 3563 57 27.314815 6 8 59 11 15 15 91.11111 3563 57 27.314815 6 8 60 14 15 10 86.66667 3562 60 28.240741 5 7 61 14 14 11 86.66667 3562 60 28.240741 5 7 62 15 10 13 84.44444 3561 62 29.166667 5 7 63 12 13 15 88.88889 3561 62 29.166667 5 7 64 15 12 11 84.44444 3560 64 31.018519 5 7 65 15 11 12 84.44444 3560 64 31.018519 5 7 66 12 15 13 88.88889 3560 64 31.018519 5 7 67 12 14 14 88.88889 3560 64 31.018519 5 7 68 15 13 10 84.44444 3559 68 31.481481 5 7 69 13 13 13 86.66667 3558 69 32.870370 5 6 70 13 12 14 86.66667 3558 69 32.870370 5 6 71 13 11 15 86.66667 3558 69 32.870370 5 6 72 13 15 11 86.66667 3557 72 33.796296 4 6 73 13 14 12 86.66667 3557 72 33.796296 4 6 74 14 10 14 84.44444 3556 74 34.259259 4 6 75 14 12 12 84.44444 3555 75 36.111111 4 6 76 14 11 13 84.44444 3555 75 36.111111 4 6 77 11 15 14 88.88889 3555 75 36.111111 4 6 78 11 14 15 88.88889 3555 75 36.111111 4 6 79 14 14 10 84.44444 3554 79 37.037037 4 5 80 14 13 11 84.44444 3554 79 37.037037 4 5 81 12 12 15 86.66667 3553 81 37.500000 4 5 82 15 11 11 82.22222 3552 82 39.814815 4 5 83 15 10 12 82.22222 3552 82 39.814815 4 5 84 12 15 12 86.66667 3552 82 39.814815 4 5 85 12 14 13 86.66667 3552 82 39.814815 4 5 86 12 13 14 86.66667 3552 82 39.814815 4 5 87 15 12 10 82.22222 3551 87 40.277778 4 5 88 13 12 13 84.44444 3550 88 41.666667 3 5 89 13 11 14 84.44444 3550 88 41.666667 3 5 90 13 10 15 84.44444 3550 88 41.666667 3 5 91 13 14 11 84.44444 3549 91 43.055556 3 5 92 13 13 12 84.44444 3549 91 43.055556 3 5 93 10 15 15 88.88889 3549 91 43.055556 3 5 94 13 15 10 84.44444 3548 94 43.518519 3 5 95 14 11 12 82.22222 3547 95 45.370370 3 4 96 14 10 13 82.22222 3547 95 45.370370 3 4 97 11 14 14 86.66667 3547 95 45.370370 3 4 98 11 13 15 86.66667 3547 95 45.370370 3 4 99 14 13 10 82.22222 3546 99 46.759259 3 4 100 14 12 11 82.22222 3546 99 46.759259 3 4 101 11 15 13 86.66667 3546 99 46.759259 3 4 102 12 12 14 84.44444 3545 102 47.685185 3 4 103 12 11 15 84.44444 3545 102 47.685185 3 4 104 15 10 11 80.00000 3544 104 49.074074 3 4 105 12 14 12 84.44444 3544 104 49.074074 3 4 106 12 13 13 84.44444 3544 104 49.074074 3 4 107 15 11 10 80.00000 3543 107 50.000000 3 4 108 12 15 11 84.44444 3543 107 50.000000 3 4 109 13 11 13 82.22222 3542 109 51.388889 3 4 110 13 10 14 82.22222 3542 109 51.388889 3 4 111 10 14 15 86.66667 3542 109 51.388889 3 4 112 13 13 11 82.22222 3541 112 52.777778 3 4 113 13 12 12 82.22222 3541 112 52.777778 3 4 114 10 15 14 86.66667 3541 112 52.777778 3 4 115 13 14 10 82.22222 3540 115 53.240741 3 4 116 14 10 12 80.00000 3539 116 55.092593 3 3 117 11 14 13 84.44444 3539 116 55.092593 3 3 118 11 13 14 84.44444 3539 116 55.092593 3 3 119 11 12 15 84.44444 3539 116 55.092593 3 3 120 14 12 10 80.00000 3538 120 56.481481 2 3 121 14 11 11 80.00000 3538 120 56.481481 2 3 122 11 15 12 84.44444 3538 120 56.481481 2 3 123 12 11 14 82.22222 3537 123 57.407407 2 3 124 12 10 15 82.22222 3537 123 57.407407 2 3 125 12 13 12 82.22222 3536 125 58.333333 2 3 126 12 12 13 82.22222 3536 125 58.333333 2 3 127 15 10 10 77.77778 3535 127 59.722222 2 3 128 12 15 10 82.22222 3535 127 59.722222 2 3 129 12 14 11 82.22222 3535 127 59.722222 2 3 130 13 10 13 80.00000 3534 130 60.648148 2 3 131 10 13 15 84.44444 3534 130 60.648148 2 3 132 13 12 11 80.00000 3533 132 62.500000 2 3 133 13 11 12 80.00000 3533 132 62.500000 2 3 134 10 15 13 84.44444 3533 132 62.500000 2 3 135 10 14 14 84.44444 3533 132 62.500000 2 3 136 13 13 10 80.00000 3532 136 62.962963 2 3 137 14 10 11 77.77778 3531 137 64.814815 2 3 138 11 13 13 82.22222 3531 137 64.814815 2 3 139 11 12 14 82.22222 3531 137 64.814815 2 3 140 11 11 15 82.22222 3531 137 64.814815 2 3 141 14 11 10 77.77778 3530 141 66.203704 2 3 142 11 15 11 82.22222 3530 141 66.203704 2 3 143 11 14 12 82.22222 3530 141 66.203704 2 3 144 12 10 14 80.00000 3529 144 66.666667 2 3 145 12 12 12 80.00000 3528 145 67.592593 2 3 146 12 11 13 80.00000 3528 145 67.592593 2 3 147 12 14 10 80.00000 3527 147 68.518519 2 2 148 12 13 11 80.00000 3527 147 68.518519 2 2 149 13 10 12 77.77778 3526 149 69.444444 2 2 150 10 12 15 82.22222 3526 149 69.444444 2 2 151 13 12 10 77.77778 3525 151 71.296296 2 2 152 13 11 11 77.77778 3525 151 71.296296 2 2 153 10 14 13 82.22222 3525 151 71.296296 2 2 154 10 13 14 82.22222 3525 151 71.296296 2 2 155 10 15 12 82.22222 3524 155 71.759259 2 2 156 11 12 13 80.00000 3523 156 73.148148 2 2 157 11 11 14 80.00000 3523 156 73.148148 2 2 158 11 10 15 80.00000 3523 156 73.148148 2 2 159 14 10 10 75.55556 3522 159 74.537037 2 2 160 11 14 11 80.00000 3522 159 74.537037 2 2 161 11 13 12 80.00000 3522 159 74.537037 2 2 162 11 15 10 80.00000 3521 162 75.000000 2 2 163 12 12 11 77.77778 3520 163 76.388889 2 2 164 12 11 12 77.77778 3520 163 76.388889 2 2 165 12 10 13 77.77778 3520 163 76.388889 2 2 166 12 13 10 77.77778 3519 166 76.851852 2 2 167 10 12 14 80.00000 3518 167 77.777778 2 2 168 10 11 15 80.00000 3518 167 77.777778 2 2 169 13 11 10 75.55556 3517 169 79.629630 2 2 170 13 10 11 75.55556 3517 169 79.629630 2 2 171 10 14 12 80.00000 3517 169 79.629630 2 2 172 10 13 13 80.00000 3517 169 79.629630 2 2 173 10 15 11 80.00000 3516 173 80.092593 1 2 174 11 12 12 77.77778 3515 174 81.481481 1 2 175 11 11 13 77.77778 3515 174 81.481481 1 2 176 11 10 14 77.77778 3515 174 81.481481 1 2 177 11 14 10 77.77778 3514 177 82.407407 1 2 178 11 13 11 77.77778 3514 177 82.407407 1 2 179 12 11 11 75.55556 3512 179 83.333333 1 2 180 12 10 12 75.55556 3512 179 83.333333 1 2 181 12 12 10 75.55556 3511 181 83.796296 1 2 182 10 11 14 77.77778 3510 182 84.722222 1 2 183 10 10 15 77.77778 3510 182 84.722222 1 2 184 13 10 10 73.33333 3509 184 86.111111 1 2 185 10 13 12 77.77778 3509 184 86.111111 1 2 186 10 12 13 77.77778 3509 184 86.111111 1 2 187 10 15 10 77.77778 3508 187 87.037037 1 2 188 10 14 11 77.77778 3508 187 87.037037 1 2 189 11 11 12 75.55556 3507 189 87.962963 1 2 190 11 10 13 75.55556 3507 189 87.962963 1 2 191 11 13 10 75.55556 3506 191 88.888889 1 2 192 11 12 11 75.55556 3506 191 88.888889 1 2 193 12 10 11 73.33333 3504 193 89.351852 1 2 194 12 11 10 73.33333 3503 194 89.814815 1 2 195 10 11 13 75.55556 3502 195 90.740741 1 1 196 10 10 14 75.55556 3502 195 90.740741 1 1 197 10 13 11 75.55556 3501 197 91.666667 1 1 198 10 12 12 75.55556 3501 197 91.666667 1 1 199 10 14 10 75.55556 3500 199 92.129630 1 1 200 11 10 12 73.33333 3499 200 92.592593 1 1
使用したRプログラムはこちら
#pokemon go options(max.print = 99999) AT <- DF <- HP <- 0:15 #0:15, 卵なら10:15 x <- expand.grid(AT, DF, HP) kotaichi <- rowSums(x) N <- nrow(x) K <- kotaichi /45 * 100 hist(K, breaks = seq(0, 100, by = 1)) dat <- data.frame(AT = x[,1], DF = x[,2], HP = x[,3], kotaichi_per = K) cplv.list <- structure(c(1, 0.094, 1.5, 0.135137432009048, 2, 0.166397869998387, 2.5, 0.192650918995991, 3, 0.215732469994667, 3.5, 0.236572660994461, 4, 0.255720049996085, 4.5, 0.273530380999991, 5, 0.290249899999983, 5.5, 0.306057377006338, 6, 0.321087600040861, 6.5, 0.335445036034221, 7, 0.34921268003897, 7.5, 0.362457751055209, 8, 0.375235589996471, 8.5, 0.387592406014359, 9, 0.399567279941689, 9.5, 0.411193551019468, 10, 0.4225, 10.5, 0.432926419036769, 11, 0.443107550037234, 11.5, 0.453059957842226, 12, 0.462798390014486, 12.5, 0.472336082995996, 13, 0.481684950045151, 13.5, 0.490855800006478, 14, 0.499858439960755, 14.5, 0.508701764986126, 15, 0.517393950003283, 15.5, 0.52594251102188, 16, 0.534354330009592, 16.5, 0.542635766974496, 17, 0.550792690038639, 17.5, 0.558830576024612, 18, 0.566754519964332, 18.5, 0.574569153018155, 19, 0.582278909973562, 19.5, 0.589887916980845, 20, 0.597400009959826, 20.5, 0.604818814026151, 21, 0.612157289999882, 21.5, 0.619399365030349, 22, 0.626567130002843, 22.5, 0.633644532999378, 23, 0.640652949965892, 23.5, 0.647576425991558, 24, 0.654435629989688, 24.5, 0.66121480601995, 25, 0.667934000032937, 25.5, 0.67457753698148, 26, 0.681164919971662, 26.5, 0.687680647975497, 27, 0.694143649974557, 27.5, 0.700538673022411, 28, 0.706884209966526, 28.5, 0.713164995986202, 29, 0.719399090004985, 29.5, 0.72557155201951, 30, 0.7317, 30.5, 0.734741008995687, 31, 0.737769479986805, 31.5, 0.740785573968608, 32, 0.74378943001363, 32.5, 0.746781210998241, 33, 0.749761039998745, 33.5, 0.752729086989469, 34, 0.755685509984147, 34.5, 0.758630377983903, 35, 0.761563840002925, 35.5, 0.764486065013614, 36, 0.767397169984357, 36.5, 0.770297265995408, 37, 0.773186500011478, 37.5, 0.776064961971612, 38, 0.778932749985517, 38.5, 0.781790055001981, 39, 0.784636970005875, 39.5, 0.787473577969445, 40, 0.790300009996204, 40.5, 0.7931164), .Dim = c(2L, 80L)) cplv <- function(LV){ cplv.list[2, which(cplv.list[1,]==LV)] } cp.cal <- function(AT, DF, HP, SHUAT, SHUDF, SHUHP, LV){ cp <- (AT + SHUAT)*sqrt(SHUDF + DF)*sqrt(SHUHP + HP)*(cplv(LV)^2) / 10 floor(cp) } dat$Bangiras_max_cp <- apply(dat, 1, function(x) cp.cal(x[1], x[2], x[3], 251, 212, 200, 39)) dat <- dat[order(-dat$Bangiras_max_cp),] rownames(dat) <- 1:N dat$rank <- rank(-dat$Bangiras_max_cp, ties.method = c("min")) dat$joui_per <- rank(-dat$Bangiras_max_cp, ties.method = c("max")) / N * 100 calc.egg.N <- function(joui_per, p){ ceiling(log(1 - p) / log(1 - joui_per / 100)) } dat$N80per <- calc.egg.N(dat$joui_per, 0.8) dat$N90per <- calc.egg.N(dat$joui_per, 0.9) dat[1:200,] hist(dat$Bangiras_max_cp, freq = FALSE, breaks = seq(3180, 3620, by = 10), main = "", xlab = "バンギラスMAX CP")
R3.4.0にアップデート。Rが高速に。
R3.4.0が出たと聞いた。
JITバイトコンパイラになり、loopが高速になるらしい。
適当にフィボナッチ数を計算する関数を作って試してみる。
fibonacci <- function(n){ if(n == 1) 1 else if(n == 2) 2 else (fibonacci(n-2) + fibonacci(n-1)) } fibonacci2 <- function(n){ f1 <- 1 f2 <- 2 f3 <- 1 if(n >= 3){ for(i in 3:n){ f3 <- f1 + f2 f1 <- f2 f2 <- f3 } return(f3) } return(n) } library(rbenchmark) benchmark(fibonacci(20),fibonacci2(20), replications = 100) benchmark(fibonacci2(200000), replications = 100)
#R3.3.0 > benchmark(fibonacci2(200000), replications = 100) test replications elapsed relative user.self sys.self user.child sys.child 1 fibonacci2(2e+05) 100 17.03 1 16.94 0 NA NA > benchmark(fibonacci(20),fibonacci2(20), replications = 100) test replications elapsed relative user.self sys.self user.child sys.child 1 fibonacci(20) 100 1.81 NA 1.8 0 NA NA 2 fibonacci2(20) 100 0.00 NA 0.0 0 NA NA
#R3.4.0 > benchmark(fibonacci2(200000), replications = 100) test replications elapsed relative user.self sys.self user.child sys.child 1 fibonacci2(2e+05) 100 2.25 1 2.18 0 NA NA > benchmark(fibonacci(20),fibonacci2(20), replications = 100) test replications elapsed relative user.self sys.self user.child sys.child 1 fibonacci(20) 100 0.81 NA 0.8 0 NA NA 2 fibonacci2(20) 100 0.00 NA 0.0 0 NA NA
fibonacciを見ると再帰関数も、2倍以上高速になっている。
fibonacci2をみると、for文での関数は7倍以上高速になっている。
単純に足し算でも比較
f.loop <- function(n){ x <- 0 for(i in seq_len(n)){ x <- x + 1 } return(x) } f.loop2 <- function(n){ x <- rep(1, n) y <- 0 for(i in seq_len(n)){ y <- y + x[i] } return(y) } f.vec <- function(n){ x <- rep(1, n) sum(x) } library(rbenchmark) benchmark(f.vec(100000),f.loop(100000), f.loop2(100000), replications = 1000)
#R3.3.0 > benchmark(f.vec(100000),f.loop(100000), f.loop2(100000), replications = 1000) test replications elapsed relative user.self sys.self user.child sys.child 2 f.loop(1e+05) 1000 45.19 72.887 43.75 0.14 NA NA 3 f.loop2(1e+05) 1000 67.94 109.581 65.74 0.25 NA NA 1 f.vec(1e+05) 1000 0.62 1.000 0.42 0.19 NA NA
#R.3.4.0 > benchmark(f.vec(100000),f.loop(100000), f.loop2(100000), replications = 1000) test replications elapsed relative user.self sys.self user.child sys.child 2 f.loop(1e+05) 1000 4.24 6.424 4.07 0.00 NA NA 3 f.loop2(1e+05) 1000 6.69 10.136 6.26 0.17 NA NA 1 f.vec(1e+05) 1000 0.66 1.000 0.42 0.22 NA NA
loopは10倍程度高速ですね。
ここまで高速になると良いですねぇ。他言語との絶望的な差が縮まった感じ。。
アップデート後、実際今までのコードが動かなくなったりしているので面倒だけど
（パッケージのバージョンの問題もある）、今後のことを考えてR3.4系に移行するとしよう。
探偵ナイトスクープ 得点した選手の背番号の合計は偶数が多い？
ナイトスクープでサッカーの試合で得点した人の背番号の合計が
偶数有利か奇数有利かということを検証していた。
偶数背番号の選手も奇数背番号の選手もおなじ得点率なら同じになるだろう。
しかし、実際は異なる。
偶数の選手が得点することがたまたま多かったから、偶数が勝つことが多かった
との説明であったが、もやもやする。得点を量産するストライカーに奇数背番号の
人しかいなくても、2点入れば偶数になってしまう。一方、得点を量産するストラ
イカーに偶数背番号の人しかいない場合、何点取っても永遠に偶数だ。
これでは偶数有利じゃね？
Rで分析してみる
日本代表の得点は
JFAhttp://www.jfa.jp/samuraiblue/schedule_result/
から、2008-2016の日本の得点を目視で集計した。間違いがあるかも。
tokuten <- c(45, 20, 24, 10, 6, 4, 1, 1) #0点の28試合は除く1-8点の試合数 tokuten.rate <- tokuten / sum(tokuten) #各試合の1-8点になる確率(0点は除く) guu <- function(p, n){1/2 + (p - 1/2) * ( (2 * p - 1)^(n -1) )} #得点n点で偶数になる確率 nmax <- length(tokuten) p <- seq(0, 1, 0.05) gu.win <- numeric(length(p)) ct <- 1 for(i in seq(0, 1, 0.05) ){ gu.win[ct] <- sum(guu(i, 1:nmax) * tokuten.rate) # plot(guu(i, 1:nmax), type = "l", xlim = c(1,nmax), ylim = c(0, 1)) ct <- ct + 1 } plot(p, gu.win, ylim = c(0,1), type = "l", xlab = "偶数背番号が得点する率", ylab = "偶数にかけて勝つ確率") abline(h = 0.5, v = 0.5, lty = 2)
偶数背番号が得点する率は、せいぜい0.3-0.7くらいだとすると、
それほど差はなさそうだ。
ちなみに、試合結果の総得点が１点から８点まで同じ確率で起こるとすると
tokuten.rate <- rep(1/8, 8) gu.win <- numeric(length(p)) ct <- 1 for(i in seq(0, 1, 0.05)){ gu.win[ct] <- sum(guu(i, 1:nmax) * tokuten.rate) # plot(guu(i, 1:nmax), type = "l", xlim = c(1,nmax), ylim = c(0, 1)) ct <- ct + 1 } plot(p, gu.win, ylim = c(0,1), type = "l", xlab = "偶数背番号が得点する率", ylab = "偶数にかけて勝つ確率") abline(h = 0.5, v = 0.5, lty = 2)
奇数背番号が得点する率が高くても、それほど0.5からズレない。
一方、偶数背番号の得点率が高いと偶数に賭けると勝率半端ない。
圧倒的に偶数に賭けたほうが有利になる。
tokuten.rate <- rep(1/8, 8) i <- 0.2 #偶数背番号が得点する率0.2 plot(guu(i, 1:nmax), type = "l", xlim = c(1,nmax), ylim = c(0, 1), xlab = "得点", ylab = "偶数が勝つ確率") abline(h = 0.5, lty = 2)
奇数背番号が得点する率が80%の場合(つまり偶数背番号の選手が20%)
1点では偶数にかけた人の勝率は20%、しかし２点では逆転。
8点とか入るとほぼ五分五分だ。奇数の優位性が減っちゃう。
したがって２番めのグラフのようになる。
実際は１－３点が多いので、一番上のグラフのようになる。
得点分布などに大きく影響されるので、スポーツごとでも違うことになる。
結論としては、得点する人の背番号が偶数が多いか奇数が多いかを分析して
賭けるのが正解。その情報がない場合には、偶数にかけとけば良いんじゃね？（適当
受験番号と合格率の関係
受験番号と合格率の関係
rs_fan_jp
2016年4月3日
2016年4月3日20:20念のため追記。言うまでもなく受験番号から合格への直接的因果関係はないからね。統計上有意ですが、単なる見かけ上の関係です。間違っても受験番号１ゲットなどという不毛なことを目指さないように。追記終了
すでに、旬は終わっているが、かつての上司のお子さんの合否が気になり 大学のサイトをみたものの、当たり前だが合格した人の受験番号しか書かれてなくて 合否は分からず。 その時に合格した受験番号だけダウンロードしておいた。
これを使って遊んでみる。
受験番号が若いほど合格率が高いのではないか？と気になったのである。 受付開始と同時に願書を出す人と、、ぎりぎりまで悩んだ人とでは 多分学力が違いそうという仮説である。
まずは元データの再現
hist.goukaku <- function(x, bw = 10){
check.ok <- function(y){any(x == y)}
goukaku <- sapply(range(x)[1]:range(x)[2], check.ok)
hist(x, breaks = seq(range(x)[1], (2 + range(x)[2] %/% bw) * bw, by= bw))
}
###
sougouningen_bunkei <- c(1, 2, 6, 8, 12, 13, 14, 15, 16, 21, 25, 29, 30, 32, 33, 35,
37, 44, 49, 55, 56, 61, 63, 65, 69, 70, 72, 74, 75, 76, 83, 89,
90, 92, 94, 106, 107, 116, 118, 119, 125, 126, 131, 141, 146,
147, 148, 149, 150, 154, 156, 159, 163, 165, 169, 174, 178, 182,
185, 187, 195, 207, 209, 211, 217)
sougouningen_rikei <- c(3007, 3011, 3012, 3013, 3015, 3016, 3020, 3021, 3024, 3025,
3032, 3034, 3037, 3039, 3041, 3042, 3043, 3044, 3049, 3052, 3057,
3061, 3063, 3073, 3076, 3078, 3083, 3091, 3092, 3093, 3095, 3097,
3101, 3103, 3104, 3106, 3109, 3110, 3112, 3113, 3128, 3131, 3132,
3139, 3148, 3151, 3155, 3158, 3159, 3160, 3164, 3169, 3186, 3187
)
igaku <- c(1, 3, 4, 5, 7, 8, 9, 15, 16, 24, 25, 27, 33, 40, 47, 48, 54,
58, 59, 61, 62, 67, 69, 71, 74, 75, 78, 79, 83, 84, 85, 86, 88,
89, 90, 94, 97, 99, 100, 102, 103, 104, 105, 111, 113, 114, 115,
116, 118, 125, 131, 132, 133, 134, 141, 143, 145, 146, 147, 151,
152, 153, 155, 156, 158, 162, 165, 168, 169, 170, 171, 172, 179,
185, 187, 188, 190, 192, 195, 197, 198, 201, 202, 203, 206, 207,
209, 210, 211, 214, 219, 223, 224, 225, 226, 227, 230, 233, 235,
237, 242, 245, 248, 256, 258, 270, 272, 276, 278, 294, 304, 314
)
bungaku <- c(1, 6, 8, 10, 12, 13, 14, 16, 18, 22, 23, 25, 26, 29, 31, 33,
34, 35, 37, 41, 42, 43, 45, 51, 52, 54, 55, 56, 57, 58, 61, 62,
63, 65, 66, 69, 72, 75, 77, 85, 86, 91, 94, 98, 101, 103, 104,
107, 110, 112, 114, 126, 132, 135, 136, 137, 138, 139, 141, 148,
155, 156, 159, 161, 165, 167, 168, 169, 170, 171, 172, 179, 180,
181, 185, 186, 194, 195, 196, 197, 200, 204, 205, 206, 208, 209,
211, 214, 215, 216, 222, 223, 226, 228, 229, 230, 234, 237, 238,
239, 246, 251, 253, 254, 255, 256, 258, 261, 264, 265, 269, 273,
281, 286, 290, 300, 302, 303, 305, 306, 307, 308, 310, 311, 315,
317, 322, 323, 329, 331, 332, 336, 337, 338, 339, 340, 345, 352,
353, 357, 359, 361, 363, 364, 366, 367, 372, 375, 376, 378, 382,
383, 393, 397, 399, 405, 406, 407, 411, 423, 425, 426, 428, 432,
434, 436, 438, 446, 448, 449, 450, 454, 456, 457, 460, 461, 463,
464, 470, 477, 481, 484, 487, 491, 493, 494, 495, 502, 505, 507,
508, 512, 513, 515, 516, 519, 522, 524, 528, 532, 535, 537, 538,
540, 551, 557, 562, 563, 568, 578, 580, 595, 598, 612, 613, 614
)
kyoiku_bunkei <- c(1, 2, 3, 4, 8, 11, 14, 15, 18, 19, 20, 21, 24, 25, 26, 28,
35, 36, 37, 38, 39, 41, 42, 44, 48, 50, 54, 60, 63, 64, 68, 69,
72, 80, 82, 83, 86, 95, 96, 99, 101, 104, 110, 114, 124, 128,
136)
kyoiku_rikei <- c(3001, 3003, 3007, 3008, 3015, 3016, 3018, 3022, 3028, 3038)
hougaku <- c(1, 3, 4, 5, 9, 11, 12, 18, 20, 23, 24, 27, 29, 36, 37, 38,
39, 40, 41, 43, 44, 46, 48, 50, 51, 52, 53, 55, 57, 59, 68, 70,
71, 76, 78, 81, 82, 84, 85, 86, 87, 88, 89, 90, 93, 94, 95, 96,
97, 98, 99, 100, 102, 103, 104, 105, 107, 108, 109, 112, 116,
117, 123, 124, 125, 128, 133, 134, 136, 137, 139, 141, 145, 146,
149, 150, 152, 153, 155, 156, 161, 162, 163, 167, 168, 172, 175,
177, 180, 182, 183, 184, 185, 189, 191, 192, 193, 195, 196, 197,
202, 206, 210, 212, 213, 216, 219, 222, 224, 225, 228, 229, 232,
234, 238, 240, 244, 246, 248, 249, 250, 251, 253, 254, 255, 256,
257, 258, 259, 260, 262, 263, 264, 265, 266, 267, 269, 270, 272,
275, 277, 279, 283, 284, 289, 290, 295, 296, 297, 300, 302, 303,
304, 305, 306, 311, 312, 315, 316, 318, 319, 321, 323, 324, 325,
326, 332, 334, 336, 338, 339, 342, 343, 345, 350, 351, 353, 355,
356, 359, 360, 362, 363, 364, 365, 367, 368, 370, 371, 372, 373,
374, 376, 379, 382, 383, 386, 387, 389, 390, 392, 393, 399, 400,
401, 403, 404, 407, 410, 411, 412, 421, 424, 428, 433, 435, 440,
441, 442, 444, 446, 448, 450, 451, 452, 457, 462, 463, 464, 467,
479, 480, 484, 486, 487, 488, 507, 517, 519, 528, 531, 537, 538,
540, 541, 543, 548, 550, 562, 563, 577, 580, 582, 595, 598, 599,
600, 602, 603, 605, 610, 625, 627, 629, 632, 636, 638, 639, 641,
645, 649, 652, 654, 655, 656, 663, 666, 673, 674, 682, 686, 693,
694, 695, 697, 698, 705, 709, 712, 718, 721, 722, 726, 728, 729,
737, 742, 746, 748, 752, 753, 755, 760, 761, 764, 766, 776, 787,
789, 801, 804, 814)
keizai_bunkei <- c(2, 3, 6, 7, 9, 11, 17, 19, 20, 22, 24, 28, 30, 31, 36, 37,
40, 41, 44, 45, 46, 47, 48, 50, 51, 52, 58, 59, 60, 62, 65, 66,
67, 68, 69, 70, 72, 73, 74, 75, 77, 79, 80, 81, 82, 86, 89, 90,
91, 92, 93, 94, 96, 97, 102, 103, 104, 105, 106, 107, 108, 109,
110, 111, 112, 113, 114, 116, 118, 123, 124, 125, 128, 130, 131,
132, 134, 135, 137, 142, 144, 145, 148, 149, 160, 163, 165, 168,
169, 170, 171, 181, 183, 187, 188, 191, 193, 194, 197, 198, 200,
201, 203, 204, 206, 215, 218, 219, 223, 224, 225, 229, 232, 236,
239, 240, 244, 245, 246, 247, 252, 255, 256, 259, 261, 264, 273,
274, 278, 279, 280, 281, 282, 284, 286, 289, 290, 296, 299, 303,
305, 309, 310, 312, 319, 320, 321, 322, 324, 325, 326, 327, 332,
336, 338, 340, 346, 349, 357, 360, 363, 370, 374, 378, 380, 384,
388, 391, 393, 394, 400, 401, 402, 405, 406, 410, 414, 415, 416,
419, 421, 429, 433, 436, 442, 447, 450, 452, 462, 468)
keizai_rikei <- c(4002, 4003, 4005, 4006, 4024, 4031, 4033, 4040, 4042, 4045,
4046, 4051, 4054, 4055, 4058, 4059, 4071, 4072, 4081, 4084, 4097,
4102, 4110, 4112, 4131)
rigaku <- c(1, 3, 9, 10, 11, 12, 14, 15, 17, 22, 24, 26, 32, 35, 36, 37,
39, 41, 43, 44, 45, 46, 50, 51, 55, 56, 58, 59, 61, 62, 63, 67,
68, 72, 75, 76, 77, 82, 83, 85, 86, 93, 97, 101, 102, 103, 104,
105, 110, 111, 113, 117, 119, 120, 123, 124, 126, 127, 128, 129,
130, 132, 133, 134, 135, 139, 140, 145, 147, 148, 150, 151, 153,
154, 156, 157, 160, 162, 164, 165, 171, 172, 174, 178, 181, 182,
183, 184, 185, 187, 189, 190, 192, 193, 196, 197, 199, 200, 202,
203, 204, 206, 207, 218, 219, 222, 229, 230, 233, 234, 236, 238,
239, 240, 241, 242, 254, 256, 258, 263, 264, 265, 267, 268, 271,
277, 281, 283, 286, 287, 289, 290, 292, 298, 299, 304, 307, 310,
312, 314, 316, 317, 318, 321, 322, 324, 326, 327, 328, 332, 334,
336, 340, 341, 342, 343, 345, 351, 352, 354, 357, 358, 360, 361,
363, 365, 366, 372, 375, 380, 381, 384, 385, 386, 388, 389, 391,
392, 393, 394, 395, 400, 402, 408, 409, 413, 414, 415, 420, 423,
429, 432, 436, 440, 449, 451, 453, 454, 458, 463, 469, 473, 474,
475, 476, 478, 481, 482, 484, 485, 493, 500, 509, 517, 522, 524,
526, 533, 536, 540, 541, 542, 551, 555, 557, 559, 564, 565, 569,
573, 574, 579, 580, 582, 590, 592, 596, 604, 605, 607, 608, 610,
613, 616, 618, 620, 624, 627, 633, 634, 635, 636, 637, 644, 651,
653, 655, 656, 657, 661, 669, 670, 672, 673, 676, 679, 680, 683,
684, 687, 689, 690, 692, 693, 696, 698, 700, 701, 706, 707, 711,
721, 722, 728, 732, 734, 735, 739, 741, 742, 745, 754, 758, 760,
766, 769, 772, 774, 781, 788, 791, 798, 799, 804, 809, 810, 822,
824, 827, 829)
kango <- c(2001, 2002, 2003, 2010, 2013, 2015, 2016, 2017, 2018, 2020,
2021, 2023, 2025, 2027, 2031, 2033, 2034, 2035, 2036, 2037, 2038,
2039, 2040, 2042, 2043, 2044, 2046, 2047, 2048, 2049, 2050, 2052,
2053, 2055, 2056, 2057, 2061, 2062, 2063, 2064, 2065, 2066, 2067,
2071, 2075, 2078, 2079, 2080, 2088, 2091, 2094, 2095, 2096, 2097,
2098, 2100, 2101, 2102, 2105, 2106, 2107, 2110, 2113, 2114, 2115,
2117, 2118, 2120, 2123, 2124, 2125, 2129, 2131)
kensagijutu <- c(4001, 4002, 4005, 4006, 4008, 4012, 4013, 4014, 4016, 4017,
4019, 4023, 4025, 4027, 4029, 4031, 4034, 4036, 4037, 4041, 4044,
4048, 4049, 4050, 4053, 4054, 4055, 4057, 4058, 4066, 4068, 4070,
4071, 4072, 4076, 4077, 4081, 4084, 4085)
rigakuryouhou <- c(5002, 5003, 5005, 5006, 5007, 5008, 5009, 5012, 5014, 5015,
5019, 5020, 5025, 5029, 5032, 5033, 5034, 5036)
sagyouryouhou <- c(6001, 6002, 6003, 6004, 6005, 6006, 6007, 6009, 6010, 6011,
6012, 6015, 6018, 6019, 6021, 6022, 6024, 6025)
yakka <- c(1, 3, 7, 8, 9, 10, 11, 12, 15, 17, 19, 20, 25, 26, 27, 29,
31, 32, 34, 37, 39, 46, 49, 50, 51, 52, 55, 56, 57, 59, 60, 62,
64, 67, 69, 72, 73, 78, 83, 89, 91, 92, 93, 97, 102, 103, 107,
110, 113, 116, 117, 119, 121)
yakugaku <- c(3001, 3003, 3004, 3005, 3007, 3008, 3010, 3013, 3014, 3016,
3018, 3020, 3021, 3023, 3025, 3027, 3028, 3031, 3032, 3036, 3040,
3043, 3044, 3047, 3048, 3050, 3054, 3057, 3058, 3059, 3067)
kougaku <- c(5, 6, 10, 16, 21, 22, 23, 26, 28, 30, 32, 36, 44, 45, 46, 48,
50, 51, 52, 53, 54, 55, 56, 57, 59, 60, 61, 62, 63, 66, 71, 73,
80, 85, 86, 87, 88, 89, 94, 98, 99, 101, 106, 108, 109, 112,
113, 114, 115, 116, 118, 119, 120, 122, 124, 129, 130, 132, 133,
137, 138, 142, 143, 147, 148, 156, 157, 160, 161, 164, 167, 171,
172, 179, 185, 191, 192, 193, 196, 198, 200, 203, 207, 209, 211,
212, 213, 221, 227, 231, 232, 234, 239, 243, 247, 249, 250, 253,
254, 255, 257, 259, 262, 264, 265, 270, 274, 275, 277, 278, 282,
285, 292, 293, 298, 303, 304, 305, 306, 311, 313, 316, 317, 318,
319, 320, 321, 323, 329, 331, 332, 333, 341, 344, 345, 346, 347,
348, 349, 350, 353, 354, 357, 358, 359, 360, 361, 362, 363, 368,
369, 372, 373, 374, 377, 381, 383, 385, 387, 388, 390, 391, 397,
400, 401, 403, 404, 406, 407, 409, 412, 414, 418, 421, 423, 425,
426, 430, 432, 433, 436, 437, 438, 439, 450, 453, 454, 458, 459,
462, 463, 465, 466, 467, 471, 478, 485, 487, 490, 491, 496, 497,
499, 500, 501, 503, 504, 510, 511, 514, 515, 517, 518, 521, 527,
530, 532, 533, 536, 542, 545, 549, 550, 555, 558, 560, 562, 563,
565, 566, 568, 571, 575, 579, 582, 585, 586, 587, 588, 589, 593,
594, 607, 613, 616, 622, 624, 628, 629, 630, 631, 635, 636, 637,
638, 643, 645, 647, 651, 652, 654, 655, 657, 659, 661, 662, 665,
669, 671, 672, 673, 678, 680, 681, 696, 698, 700, 704, 707, 708,
711, 712, 714, 718, 719, 720, 722, 723, 728, 729, 731, 732, 734,
737, 738, 740, 741, 742, 743, 745, 747, 748, 750, 753, 756, 760,
761, 762, 766, 767, 769, 770, 773, 774, 777, 778, 781, 787, 791,
795, 797, 801, 802, 813, 815, 817, 818, 822, 825, 828, 830, 834,
838, 840, 841, 842, 843, 845, 846, 848, 852, 854, 855, 857, 858,
859, 866, 867, 872, 873, 874, 875, 876, 877, 878, 881, 884, 885,
888, 889, 890, 892, 895, 897, 899, 901, 902, 908, 909, 910, 913,
914, 915, 916, 917, 918, 919, 922, 923, 927, 929, 932, 936, 938,
939, 945, 946, 949, 950, 952, 955, 958, 968, 970, 971, 973, 976,
979, 984, 985, 988, 989, 992, 993, 994, 998, 999, 1001, 1002,
1003, 1007, 1013, 1014, 1016, 1017, 1018, 1020, 1021, 1024, 1026,
1028, 1029, 1034, 1035, 1040, 1041, 1043, 1045, 1048, 1051, 1052,
1053, 1058, 1062, 1063, 1076, 1078, 1079, 1083, 1084, 1088, 1090,
1091, 1092, 1101, 1102, 1106, 1110, 1112, 1116, 1117, 1118, 1122,
1123, 1124, 1126, 1127, 1128, 1129, 1130, 1131, 1133, 1136, 1138,
1139, 1140, 1142, 1144, 1145, 1149, 1151, 1152, 1156, 1166, 1168,
1169, 1170, 1174, 1176, 1185, 1189, 1191, 1194, 1195, 1197, 1199,
1203, 1206, 1208, 1214, 1220, 1221, 1223, 1224, 1225, 1226, 1230,
1231, 1234, 1235, 1240, 1252, 1253, 1255, 1256, 1257, 1258, 1261,
1262, 1265, 1266, 1267, 1268, 1273, 1274, 1275, 1277, 1278, 1279,
1283, 1293, 1294, 1296, 1302, 1303, 1306, 1309, 1311, 1312, 1313,
1314, 1317, 1318, 1319, 1320, 1321, 1324, 1326, 1329, 1330, 1331,
1332, 1336, 1338, 1339, 1340, 1341, 1342, 1347, 1349, 1351, 1352,
1353, 1354, 1355, 1356, 1360, 1362, 1366, 1367, 1370, 1373, 1374,
1375, 1378, 1385, 1386, 1387, 1390, 1391, 1393, 1399, 1400, 1403,
1404, 1408, 1409, 1411, 1417, 1418, 1421, 1425, 1428, 1431, 1432,
1433, 1437, 1440, 1441, 1444, 1452, 1456, 1465, 1466, 1467, 1468,
1478, 1479, 1482, 1484, 1485, 1489, 1496, 1497, 1499, 1500, 1505,
1512, 1513, 1515, 1516, 1520, 1522, 1525, 1527, 1531, 1536, 1539,
1544, 1546, 1551, 1558, 1560, 1562, 1565, 1566, 1568, 1569, 1577,
1579, 1580, 1583, 1588, 1600, 1606, 1609, 1616, 1617, 1618, 1622,
1624, 1628, 1634, 1639, 1644, 1645, 1646, 1652, 1654, 1655, 1660,
1669, 1671, 1676, 1677, 1685, 1686, 1689, 1690, 1691, 1692, 1693,
1704, 1706, 1710, 1714, 1715, 1720, 1721, 1724, 1725, 1727, 1733,
1735, 1737, 1738, 1741, 1745, 1746, 1750, 1751, 1753, 1754, 1755,
1756, 1761, 1763, 1764, 1766, 1767, 1769, 1772, 1774, 1783, 1796,
1797, 1798, 1804, 1805, 1806, 1807, 1811, 1813, 1814, 1820, 1822,
1827, 1828, 1830, 1831, 1838, 1840, 1841, 1842, 1846, 1851, 1855,
1856, 1857, 1859, 1860, 1864, 1865, 1875, 1876, 1877, 1882, 1883,
1886, 1888, 1889, 1896, 1897, 1903, 1915, 1916, 1917, 1920, 1921,
1922, 1928, 1931, 1932, 1935, 1946, 1952, 1962, 1963, 1964, 1965,
1967, 1969, 1972, 1973, 1974, 1975, 1977, 1981, 1983, 1985, 1987,
1990, 1991, 1992, 1994, 1996, 1997, 1999, 2017, 2023, 2026, 2027,
2028, 2036, 2038, 2046, 2048, 2050, 2055, 2058, 2059, 2061, 2063,
2068, 2069, 2073, 2074, 2075, 2076, 2077, 2083, 2084, 2086, 2090,
2092, 2093, 2100, 2103, 2104, 2105, 2107, 2109, 2110, 2112, 2115,
2119, 2120, 2124, 2127, 2129, 2131, 2132, 2134, 2141, 2144, 2148,
2154, 2156, 2165, 2170, 2171, 2172, 2177, 2186, 2190, 2192, 2193,
2197, 2198, 2200, 2204, 2205, 2208, 2210, 2221, 2225, 2231, 2232,
2237, 2241, 2242, 2243, 2251, 2261, 2265, 2266, 2267, 2269, 2274,
2287, 2291, 2292, 2295, 2298, 2300, 2304, 2305, 2312, 2316, 2319,
2320, 2323, 2324, 2326, 2329, 2330, 2333, 2334, 2336, 2340, 2357,
2358, 2359, 2360, 2362, 2364, 2365, 2367, 2372, 2375, 2381, 2383,
2384, 2386, 2390, 2394, 2397, 2406, 2414, 2415, 2420, 2421, 2429,
2438, 2439, 2441, 2443, 2446, 2449, 2452, 2455, 2459, 2466, 2474,
2475, 2486, 2487, 2489, 2492, 2493, 2495, 2515, 2516, 2522, 2524,
2529, 2537, 2540, 2545, 2553, 2558, 2559, 2571, 2581, 2585, 2591,
2596, 2597, 2607, 2613, 2616, 2621, 2623, 2630, 2643, 2647, 2651,
2652, 2667, 2691, 2705, 2716, 2723, 2725, 2727)
kougaku_gakka <- structure(c(3L, 3L, 4L, 6L, 6L, 6L, 6L, 1L, 1L, 1L, 1L, 2L, 1L,
3L, 6L, 3L, 3L, 3L, 3L, 3L, 3L, 6L, 3L, 1L, 3L, 4L, 1L, 3L, 1L,
4L, 5L, 6L, 1L, 2L, 2L, 2L, 2L, 3L, 1L, 2L, 2L, 2L, 1L, 3L, 1L,
3L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 3L, 3L, 1L, 6L, 3L, 1L, 4L,
6L, 6L, 6L, 5L, 4L, 4L, 6L, 3L, 3L, 3L, 6L, 6L, 6L, 5L, 6L, 6L,
6L, 6L, 1L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 6L, 6L, 6L, 5L, 5L,
4L, 1L, 2L, 2L, 4L, 4L, 4L, 6L, 4L, 6L, 3L, 1L, 1L, 1L, 1L, 1L,
2L, 2L, 2L, 6L, 3L, 3L, 2L, 2L, 1L, 1L, 2L, 2L, 3L, 3L, 3L, 3L,
3L, 3L, 3L, 3L, 4L, 4L, 4L, 1L, 6L, 6L, 5L, 3L, 3L, 3L, 1L, 4L,
4L, 5L, 5L, 4L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 5L, 4L, 6L, 6L, 3L,
6L, 3L, 1L, 3L, 3L, 6L, 3L, 3L, 3L, 3L, 3L, 3L, 5L, 3L, 3L, 5L,
1L, 4L, 3L, 3L, 5L, 5L, 4L, 3L, 3L, 1L, 4L, 6L, 3L, 1L, 3L, 3L,
5L, 5L, 5L, 5L, 5L, 6L, 2L, 6L, 4L, 4L, 4L, 4L, 4L, 3L, 3L, 3L,
3L, 1L, 4L, 4L, 4L, 6L, 6L, 6L, 1L, 5L, 5L, 1L, 1L, 2L, 6L, 6L,
1L, 6L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 6L, 3L, 4L, 6L, 6L, 6L, 6L,
6L, 5L, 1L, 5L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 1L, 1L, 6L, 3L, 6L, 1L, 3L, 3L, 3L, 3L, 6L, 6L, 5L, 6L,
5L, 1L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 3L, 6L,
3L, 3L, 3L, 3L, 6L, 1L, 4L, 3L, 3L, 6L, 6L, 3L, 3L, 6L, 1L, 3L,
3L, 4L, 4L, 4L, 4L, 6L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 1L, 3L, 3L,
5L, 1L, 6L, 6L, 6L, 6L, 6L, 1L, 3L, 3L, 1L, 1L, 6L, 3L, 1L, 3L,
3L, 3L, 1L, 4L, 4L, 4L, 4L, 6L, 6L, 6L, 6L, 6L, 6L, 1L, 3L, 3L,
3L, 3L, 3L, 3L, 1L, 3L, 3L, 3L, 3L, 6L, 3L, 3L, 3L, 4L, 4L, 6L,
6L, 6L, 6L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 5L, 4L, 4L, 4L, 4L, 6L,
6L, 6L, 5L, 5L, 1L, 6L, 3L, 6L, 3L, 3L, 6L, 6L, 6L, 6L, 3L, 4L,
3L, 6L, 6L, 1L, 3L, 3L, 6L, 1L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 2L,
3L, 3L, 3L, 3L, 3L, 6L, 3L, 4L, 4L, 6L, 6L, 6L, 4L, 4L, 4L, 4L,
5L, 5L, 5L, 6L, 6L, 6L, 2L, 2L, 2L, 4L, 4L, 4L, 6L, 6L, 6L, 3L,
3L, 1L, 3L, 1L, 3L, 1L, 4L, 3L, 3L, 6L, 3L, 6L, 3L, 6L, 3L, 6L,
6L, 6L, 1L, 6L, 3L, 3L, 3L, 6L, 6L, 6L, 6L, 1L, 3L, 4L, 3L, 6L,
2L, 2L, 5L, 3L, 3L, 3L, 1L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 3L,
3L, 3L, 3L, 3L, 3L, 6L, 6L, 1L, 6L, 1L, 1L, 3L, 2L, 1L, 2L, 4L,
6L, 6L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 6L, 5L, 5L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 3L, 6L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 6L, 1L, 4L,
5L, 6L, 6L, 6L, 6L, 6L, 6L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L,
3L, 3L, 6L, 3L, 3L, 3L, 3L, 6L, 6L, 5L, 1L, 4L, 4L, 4L, 4L, 4L,
3L, 3L, 3L, 4L, 1L, 1L, 1L, 1L, 3L, 3L, 1L, 2L, 2L, 2L, 2L, 6L,
2L, 1L, 1L, 6L, 6L, 6L, 6L, 6L, 5L, 3L, 3L, 6L, 4L, 5L, 5L, 6L,
6L, 6L, 6L, 6L, 1L, 3L, 3L, 1L, 4L, 4L, 4L, 5L, 1L, 1L, 2L, 1L,
1L, 1L, 1L, 3L, 3L, 4L, 6L, 1L, 1L, 1L, 1L, 2L, 4L, 5L, 5L, 1L,
1L, 6L, 6L, 6L, 6L, 6L, 3L, 3L, 1L, 1L, 6L, 1L, 6L, 6L, 5L, 5L,
5L, 6L, 1L, 3L, 1L, 1L, 1L, 5L, 5L, 1L, 5L, 3L, 4L, 3L, 6L, 1L,
6L, 4L, 6L, 3L, 1L, 1L, 1L, 5L, 5L, 5L, 4L, 4L, 3L, 3L, 3L, 3L,
3L, 3L, 1L, 3L, 6L, 3L, 3L, 6L, 1L, 6L, 6L, 6L, 5L, 5L, 5L, 4L,
4L, 4L, 4L, 6L, 1L, 6L, 6L, 3L, 3L, 3L, 6L, 6L, 6L, 3L, 6L, 2L,
2L, 3L, 1L, 3L, 3L, 3L, 4L, 4L, 4L, 2L, 2L, 2L, 1L, 2L, 1L, 3L,
3L, 1L, 2L, 6L, 4L, 6L, 3L, 3L, 4L, 3L, 3L, 2L, 4L, 5L, 5L, 4L,
4L, 6L, 4L, 3L, 6L, 3L, 4L, 6L, 4L, 5L, 5L, 4L, 3L, 3L, 6L, 3L,
1L, 3L, 2L, 1L, 6L, 6L, 6L, 6L, 5L, 4L, 4L, 4L, 4L, 1L, 1L, 1L,
1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 3L, 3L, 4L,
4L, 4L, 4L, 3L, 6L, 1L, 3L, 1L, 6L, 6L, 6L, 3L, 3L, 3L, 3L, 6L,
2L, 2L, 2L, 2L, 2L, 6L, 5L, 5L, 4L, 6L, 4L, 6L, 4L, 3L, 6L, 6L,
6L, 1L, 3L, 1L, 1L, 6L, 6L, 5L, 5L, 2L, 1L, 1L, 1L, 4L, 5L, 5L,
5L, 5L, 6L, 6L, 6L, 3L, 3L, 1L, 4L, 1L, 4L, 4L, 5L, 3L, 3L, 3L,
5L, 4L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 4L, 1L, 6L, 6L, 6L, 6L, 6L,
6L, 6L, 3L, 4L, 3L, 6L, 6L, 6L, 1L, 3L, 6L, 3L, 3L, 3L, 4L, 2L,
6L, 5L, 5L, 3L, 6L, 4L, 1L, 1L, 6L, 4L, 1L, 1L, 6L, 6L, 6L, 1L,
5L, 1L, 1L, 5L, 5L, 4L, 6L, 6L, 1L, 1L, 1L, 3L, 6L, 4L, 5L, 3L,
4L, 6L, 4L, 5L, 6L, 6L, 3L, 6L, 3L, 1L, 1L, 1L, 3L, 1L, 3L, 4L,
1L, 6L, 1L, 6L), .Label = c("A", "B", "C", "D", "E", "F"),
class = "factor")
nougaku <- c(1, 3, 9, 11, 13, 14, 15, 16, 18, 20, 21, 22, 23, 24, 25, 28,
30, 35, 36, 38, 43, 52, 54, 55, 56, 60, 61, 63, 67, 69, 70, 71,
72, 76, 77, 78, 82, 83, 93, 95, 96, 97, 98, 100, 103, 106, 107,
109, 116, 119, 120, 123, 124, 128, 129, 130, 134, 142, 143, 145,
147, 149, 150, 151, 152, 154, 155, 158, 159, 163, 164, 166, 168,
172, 173, 175, 177, 178, 179, 182, 185, 186, 188, 189, 190, 194,
197, 198, 202, 205, 206, 207, 211, 212, 218, 219, 224, 225, 228,
231, 233, 236, 239, 244, 248, 249, 251, 252, 253, 256, 260, 261,
262, 266, 268, 271, 272, 274, 286, 291, 292, 296, 311, 316, 318,
319, 320, 325, 326, 329, 333, 334, 336, 337, 339, 347, 348, 350,
352, 353, 359, 366, 369, 371, 374, 376, 377, 378, 379, 381, 385,
386, 387, 391, 394, 395, 399, 400, 406, 407, 427, 428, 429, 432,
434, 435, 438, 441, 443, 447, 451, 452, 457, 460, 461, 463, 466,
468, 475, 478, 480, 482, 485, 486, 492, 495, 498, 501, 504, 505,
506, 510, 512, 521, 523, 525, 527, 531, 532, 533, 536, 538, 539,
542, 545, 546, 552, 553, 554, 558, 567, 569, 583, 584, 587, 589,
590, 591, 593, 595, 597, 598, 602, 603, 605, 609, 612, 613, 614,
618, 621, 624, 625, 626, 627, 634, 640, 642, 647, 649, 651, 652,
653, 654, 655, 660, 662, 663, 664, 665, 668, 669, 671, 672, 673,
676, 677, 680, 681, 683, 685, 686, 688, 692, 693, 695, 700, 701,
706, 707, 708, 713, 717, 720, 725, 727, 728, 731, 735, 736, 738,
745, 746, 748, 752, 754, 759, 760, 761, 765, 766, 770, 772, 773,
783, 785, 786, 791, 793, 805, 810, 819, 824, 826, 831, 833, 835,
842, 844, 845)
nougaku_gakka <- c(1, 3, 9, 11, 13, 14, 15, 16, 18, 20, 21, 22, 23, 24, 25, 28,
30, 35, 36, 38, 43, 52, 54, 55, 56, 60, 61, 63, 67, 69, 70, 71,
72, 76, 77, 78, 82, 83, 93, 95, 96, 97, 98, 100, 103, 106, 107,
109, 116, 119, 120, 123, 124, 128, 129, 130, 134, 142, 143, 145,
147, 149, 150, 151, 152, 154, 155, 158, 159, 163, 164, 166, 168,
172, 173, 175, 177, 178, 179, 182, 185, 186, 188, 189, 190, 194,
197, 198, 202, 205, 206, 207, 211, 212, 218, 219, 224, 225, 228,
231, 233, 236, 239, 244, 248, 249, 251, 252, 253, 256, 260, 261,
262, 266, 268, 271, 272, 274, 286, 291, 292, 296, 311, 316, 318,
319, 320, 325, 326, 329, 333, 334, 336, 337, 339, 347, 348, 350,
352, 353, 359, 366, 369, 371, 374, 376, 377, 378, 379, 381, 385,
386, 387, 391, 394, 395, 399, 400, 406, 407, 427, 428, 429, 432,
434, 435, 438, 441, 443, 447, 451, 452, 457, 460, 461, 463, 466,
468, 475, 478, 480, 482, 485, 486, 492, 495, 498, 501, 504, 505,
506, 510, 512, 521, 523, 525, 527, 531, 532, 533, 536, 538, 539,
542, 545, 546, 552, 553, 554, 558, 567, 569, 583, 584, 587, 589,
590, 591, 593, 595, 597, 598, 602, 603, 605, 609, 612, 613, 614,
618, 621, 624, 625, 626, 627, 634, 640, 642, 647, 649, 651, 652,
653, 654, 655, 660, 662, 663, 664, 665, 668, 669, 671, 672, 673,
676, 677, 680, 681, 683, 685, 686, 688, 692, 693, 695, 700, 701,
706, 707, 708, 713, 717, 720, 725, 727, 728, 731, 735, 736, 738,
745, 746, 748, 752, 754, 759, 760, 761, 765, 766, 770, 772, 773,
783, 785, 786, 791, 793, 805, 810, 819, 824, 826, 831, 833, 835,
842, 844, 845)
nougaku_gakka <- structure(c(5L, 1L, 3L, 5L, 1L, 2L, 2L, 3L, 3L, 2L, 6L, 2L, 4L,
5L, 1L, 4L, 1L, 5L, 3L, 3L, 4L, 5L, 6L, 1L, 5L, 3L, 1L, 5L, 6L,
1L, 1L, 6L, 1L, 1L, 2L, 2L, 1L, 5L, 5L, 3L, 6L, 6L, 1L, 4L, 2L,
3L, 4L, 5L, 5L, 5L, 6L, 1L, 2L, 2L, 4L, 6L, 2L, 5L, 1L, 5L, 2L,
6L, 1L, 1L, 1L, 2L, 6L, 5L, 4L, 1L, 1L, 2L, 1L, 1L, 1L, 4L, 4L,
5L, 1L, 2L, 4L, 6L, 2L, 1L, 2L, 6L, 1L, 4L, 2L, 3L, 2L, 3L, 5L,
6L, 5L, 2L, 5L, 3L, 4L, 6L, 3L, 1L, 1L, 1L, 3L, 5L, 2L, 1L, 2L,
5L, 2L, 1L, 1L, 4L, 1L, 1L, 1L, 5L, 1L, 4L, 1L, 5L, 4L, 1L, 5L,
6L, 1L, 1L, 1L, 5L, 3L, 2L, 2L, 1L, 1L, 3L, 1L, 5L, 2L, 1L, 3L,
6L, 5L, 3L, 1L, 4L, 1L, 6L, 4L, 1L, 5L, 4L, 4L, 4L, 1L, 4L, 1L,
1L, 5L, 6L, 5L, 1L, 1L, 1L, 5L, 2L, 3L, 2L, 3L, 6L, 3L, 5L, 1L,
3L, 3L, 6L, 5L, 1L, 2L, 5L, 6L, 5L, 2L, 5L, 1L, 2L, 6L, 1L, 3L,
3L, 2L, 4L, 2L, 3L, 1L, 1L, 1L, 3L, 5L, 1L, 6L, 3L, 4L, 2L, 6L,
3L, 1L, 1L, 1L, 5L, 2L, 5L, 6L, 1L, 1L, 2L, 4L, 1L, 2L, 5L, 2L,
5L, 1L, 1L, 1L, 1L, 5L, 5L, 3L, 1L, 5L, 1L, 3L, 1L, 5L, 1L, 4L,
1L, 3L, 4L, 2L, 2L, 5L, 1L, 1L, 5L, 3L, 4L, 3L, 5L, 5L, 1L, 1L,
1L, 3L, 5L, 6L, 6L, 6L, 1L, 3L, 5L, 2L, 6L, 1L, 6L, 3L, 1L, 2L,
1L, 1L, 1L, 3L, 5L, 1L, 2L, 1L, 4L, 1L, 1L, 3L, 5L, 2L, 6L, 2L,
2L, 5L, 2L, 6L, 1L, 4L, 6L, 4L, 1L, 2L, 2L, 6L, 5L, 5L, 5L, 3L,
1L, 2L, 4L, 3L, 5L, 5L, 6L, 1L, 1L),
.Label = c("A", "B", "C", "D", "E", "F"), class = "factor")これをデータフレーム化
name1 <- c("bungaku", "hougaku", "igaku", "kango", "keizai_bunkei",
"keizai_rikei", "kensagijutu", "kougaku", "kyoiku_bunkei", "kyoiku_rikei",
"nougaku", "rigaku", "rigakuryouhou", "sagyouryouhou", "sougouningen_bunkei",
"sougouningen_rikei", "yakka", "yakugaku")
hist.goukaku2 <- function(x, bw = 10){
eval(parse(text = paste0("z <- ", x, sep="", collapse = "")))
check.ok <- function(y){any(z == y)}
goukaku <- sapply(range(z)[1]:range(z)[2], check.ok)
hist(z, breaks = seq(range(z)[1], (2 + range(z)[2] %/% bw) * bw, by= bw))
}
data1 <- data.frame() #data frame化
for(i in name1){
eval(parse(text = paste0("tmp_data <-", i)))
data_one <- tmp_data - (tmp_data[1] %/% 10) *10 #はじめの10人が連続不合格した学科はないと仮定
#最後の受験番号は不明。最後と最後から二番目合格者の間の人数/2（切り上げ）を追加
temp_length <- length(data_one)
temp_delta <- data_one[temp_length] - data_one[temp_length - 1]
ceiling(temp_delta / 2)
num_one <- seq_len(max(data_one) + ceiling(temp_delta / 2))
goukaku <- sapply(num_one, function(x)any(x == data_one))
data_rate <- num_one / max(num_one)
tmp <- data.frame(gakka = i, num_one = num_one, num_rate = data_rate, goukaku = goukaku)
data1 <- rbind(data1, tmp)
}まずは普通のロジスティック回帰
res_1 <- glm(goukaku ~ num_rate, data = data1, family = binomial(logit))
summary(res_1)##
## Call:
## glm(formula = goukaku ~ num_rate, family = binomial(logit), data = data1)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.1564 -0.9665 -0.8186 1.3172 1.6730
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.04933 0.04627 -1.066 0.286
## num_rate -1.06729 0.08327 -12.817 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 10247 on 7833 degrees of freedom
## Residual deviance: 10079 on 7832 degrees of freedom
## AIC: 10083
##
## Number of Fisher Scoring iterations: 4res_2 <- glm(goukaku ~ num_one, data = data1, family = binomial(logit))
summary(res_2)##
## Call:
## glm(formula = goukaku ~ num_one, family = binomial(logit), data = data1)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.0145 -0.9829 -0.8923 1.3673 1.6837
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -3.959e-01 3.188e-02 -12.42 < 2e-16 ***
## num_one -2.728e-04 3.462e-05 -7.88 3.28e-15 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 10247 on 7833 degrees of freedom
## Residual deviance: 10182 on 7832 degrees of freedom
## AIC: 10186
##
## Number of Fisher Scoring iterations: 4res_3 <- glm(goukaku ~ num_rate + num_one, data = data1, family = binomial(logit))
summary(res_3)##
## Call:
## glm(formula = goukaku ~ num_rate + num_one, family = binomial(logit),
## data = data1)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.1565 -0.9678 -0.8189 1.3159 1.7075
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -4.908e-02 4.627e-02 -1.061 0.289
## num_rate -9.965e-01 9.868e-02 -10.098 <2e-16 ***
## num_one -5.413e-05 4.092e-05 -1.323 0.186
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 10247 on 7833 degrees of freedom
## Residual deviance: 10077 on 7831 degrees of freedom
## AIC: 10083
##
## Number of Fisher Scoring iterations: 4res_4 <- glm(goukaku ~ num_rate + gakka, data = data1, family = binomial(logit))
summary(res_4)##
## Call:
## glm(formula = goukaku ~ num_rate + gakka, family = binomial(logit),
## data = data1)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.6469 -0.9622 -0.8043 1.3038 2.0736
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.0866914 0.0941919 -0.920 0.357380
## num_rate -1.0810523 0.0837775 -12.904 < 2e-16 ***
## gakkahougaku 0.1308349 0.1122017 1.166 0.243586
## gakkaigaku 0.0002669 0.1461376 0.002 0.998543
## gakkakango 0.8493735 0.1966870 4.318 1.57e-05 ***
## gakkakeizai_bunkei 0.2276737 0.1277555 1.782 0.074732 .
## gakkakeizai_rikei -0.9348748 0.2379133 -3.929 8.51e-05 ***
## gakkakensagijutu 0.4424097 0.2352474 1.881 0.060024 .
## gakkakougaku -0.0223312 0.0945930 -0.236 0.813373
## gakkakyoiku_bunkei -0.0672868 0.2000343 -0.336 0.736587
## gakkakyoiku_rikei -0.5800056 0.3739524 -1.551 0.120898
## gakkanougaku 0.0673537 0.1118000 0.602 0.546876
## gakkarigaku 0.0983553 0.1120666 0.878 0.380134
## gakkarigakuryouhou 0.5864504 0.3436677 1.706 0.087925 .
## gakkasagyouryouhou 1.4775357 0.4377294 3.375 0.000737 ***
## gakkasougouningen_bunkei -0.2591714 0.1719577 -1.507 0.131764
## gakkasougouningen_rikei -0.2993656 0.1838186 -1.629 0.103399
## gakkayakka 0.3614791 0.2035938 1.775 0.075817 .
## gakkayakugaku 0.3737973 0.2567470 1.456 0.145421
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 10246.7 on 7833 degrees of freedom
## Residual deviance: 9990.3 on 7815 degrees of freedom
## AIC: 10028
##
## Number of Fisher Scoring iterations: 4とりあえず、受験番号そのものよりも、各学科での“個人の受験番号/各学科最大の受験番号” num_rateの方がより良い予測因子のようだ。 係数は負であり、各学科での受験番号が大きいほど、合格率は低いようだ。
次に混合モデル。学科ごとにランダムな合格率があると仮定。 本当は正しい合格者数があるはずだが、調べるのが面倒なのでこれで行く。
library(lme4)## Loading required package: Matrixlibrary(lmerTest)## Warning: package 'lmerTest' was built under R version 3.2.4## Warning: replacing previous import by 'ggplot2::unit' when loading 'Hmisc'## Warning: replacing previous import by 'ggplot2::arrow' when loading 'Hmisc'## Warning: replacing previous import by 'scales::alpha' when loading 'Hmisc'##
## Attaching package: 'lmerTest'
##
## 以下のオブジェクトは 'package:lme4' からマスクされています:
##
## lmer
##
## 以下のオブジェクトは 'package:stats' からマスクされています:
##
## stepres_mix1 <- glmer(goukaku ~ num_rate + (1 + num_rate|gakka), data = data1, family = binomial(logit))
summary(res_mix1)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula: goukaku ~ num_rate + (1 + num_rate | gakka)
## Data: data1
##
## AIC BIC logLik deviance df.resid
## 10057.5 10092.3 -5023.7 10047.5 7829
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.2820 -0.7676 -0.6236 1.1626 2.5547
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## gakka (Intercept) 0.12325 0.3511
## num_rate 0.07179 0.2679 -0.01
## Number of obs: 7834, groups: gakka, 18
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.02136 0.10618 0.201 0.841
## num_rate -1.10296 0.12397 -8.897 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## num_rate -0.411同じく、各学科での受験番号が大きいほど、合格率は低いようだ。
“個人の受験番号/各学科最大の受験番号”ごとの合格者のヒストグラム。 大まかに合格率と比例。
data1_goukaku <- subset(data1, goukaku == TRUE)
hist(data1_goukaku$num_rate, breaks = seq(0, 1, by = 0.05))
なんだか、“個人の受験番号/各学科最大の受験番号”が0.5くらいまでは 合格率は良好。その後下がり、最後の駆け込み１０％くらいは激しく 悪い感じ。
線形のモデルは相性が悪いと考えて、最後にgbmでの解析を追加する
library(gbm)## Loading required package: survival
## Loading required package: lattice
## Loading required package: splines
## Loading required package: parallel
## Loaded gbm 2.1.1res_gbm <- gbm(goukaku ~ num_rate + gakka, data = data1, distribution = "bernoulli",
n.trees = 10000,
interaction.depth = 3,
n.minobsinnode = 10,
shrinkage = 0.001,
bag.fraction = 0.5,
train.fraction = 1.0,
cv.folds=5,
keep.data = TRUE,
verbose = "CV",
class.stratify.cv=NULL,
n.cores = 4)
best.iter <- gbm.perf(res_gbm,method="cv")
print(best.iter)## [1] 4212summary(res_gbm, n.trees=best.iter) 
## var rel.inf
## num_rate num_rate 59.60418
## gakka gakka 40.39582print(pretty.gbm.tree(res_gbm,res_gbm$n.trees))## SplitVar SplitCodePred LeftNode RightNode MissingNode ErrorReduction
## 0 1 1.257000e+04 1 8 9 0.6199666
## 1 0 1.581571e-02 2 3 7 0.5989936
## 2 -1 6.670651e-04 -1 -1 -1 0.0000000
## 3 0 6.823725e-02 4 5 6 0.6161752
## 4 -1 -5.788456e-04 -1 -1 -1 0.0000000
## 5 -1 -9.316886e-05 -1 -1 -1 0.0000000
## 6 -1 -1.167603e-04 -1 -1 -1 0.0000000
## 7 -1 -1.021398e-04 -1 -1 -1 0.0000000
## 8 -1 2.800834e-05 -1 -1 -1 0.0000000
## 9 -1 -4.055205e-06 -1 -1 -1 0.0000000
## Weight Prediction
## 0 3917 -4.055205e-06
## 1 965 -1.021398e-04
## 2 18 6.670651e-04
## 3 947 -1.167603e-04
## 4 46 -5.788456e-04
## 5 901 -9.316886e-05
## 6 947 -1.167603e-04
## 7 965 -1.021398e-04
## 8 2952 2.800834e-05
## 9 3917 -4.055205e-06plot(res_gbm,1:2,best.iter, type = "response")
クロスバリデーションありのgbmではこんな感じになった。 やはり受験番号が大きいほど合格率は悪い感じ。
最後に受験番号ごとの合格率（最小値～＋９９まで）
bangoubetu <- function(x){
temp <- subset(data1, data1$num_one == x)
sum(temp$goukaku) / nrow(temp)
}
barplot(sapply(1:100, bangoubetu), names.arg = 1:100)
abline(h = sum(data1$goukaku) / nrow(data1))
受験番号1とか3001とかその学科最小の受験番号の合格率は ７割超えであった。 有意だが、100個の検定を繰り返したことを調整すると、有意 ではなくなる程度のもの。グラフを見る限り、たまたまの印象。
レアなタクシーに遭遇する確率
レアなタクシーに遭遇する確率
京都で学会があり、レアなタクシーに遭遇する機会があった。
ヤサカタクシーは通常三つ葉なのだが、四つ葉と双葉
が存在する。それぞれ1200台中4台と2台しか無いらしい。
今回同じ日に四つ葉と双葉を見た。なんだか良い気分。
多分２０－３０台のタクシーを見たので、その場合に
四つ葉と双葉を両方見るのがどの位レアな事象なのか
気になった。
Rで計算してみる
四つ葉1台と双葉１台自身とそれ以上にレアな組合せ
すべての確率の合計を計算する。
同じタクシーに複数回遭遇することもあるという前提
で各遭遇は独立とする。
library(gmp) f0 <- function(m=30, x=1, y=1){ n = 1200 #n 総台数 #m = m #m 出会った台数 a = 2 #a ふたば b = 4 #b よつば f1 <- function(x,y){ #x 出会ったふたば台数 #y 出会った四つ葉台数 library(gmp) n = as.bigz(n) p <- factorialZ(m)*(a^x)*(b^y)*((n-a-b)^(m-x-y)) / (factorialZ(x)*factorialZ(y)*factorialZ(m-x-y)*n^m) return(p) } list1 <- matrix(0, nrow = (m+1)*(m+2)/2, ncol = 2) ct <- 1 for(j in 0:m){ for(k in 0:m){ if(m - j - k >= 0){ list1[ct,] <- c(j,k) ct <- ct + 1 } } } pxy <- matrix(0, nrow = (m+1)*(m+2)/2, ncol = 3) bigxy <- as.bigq(numeric((m+1)*(m+2)/2)) for(i in 1:((m+1)*(m+2)/2)){ pxy[i,] <- c(list1[i,], gmp::asNumeric(f1(list1[i,1],list1[i,2]))) bigxy[i] <- f1(list1[i,1],list1[i,2]) } f2 <- function(x=1, y=1){ sum(bigxy[which(pxy[,3] <= pxy[which(pxy[,1]==x & pxy[,2]==y),3])]) } return(gmp::asNumeric(f2(x,y))) } f0(20,1,1) f0(30,1,1) L2 <- 2:50 res <- numeric(length(L2)) ct1 <- 1 for(m1 in L2){ res[ct1] <- f0(m=m1, x=1, y=1) ct1 <- ct1 + 1 } plot(res)
> f0(20,1,1)
[1] 0.004473894
> f0(30,1,1)
[1] 0.009909649
20台で0.45%しか起きない珍しい事象
30台でも0.99%しか起きない珍しい事象であった。
でも30台見れば1%位の確率で起こるのね。想像より
レアじゃかなったが、まあレアだと言えよう。
ババ抜きで１３枚以上のスタートになる確率
ババ抜きで13枚スタートになる確率
異動が決まり、週末は長距離移動かあるいはフルに仕事。
サンデープログラマーなので最近記事がかけていない。
秋からは精力的に更新したいものだ。
某アニメでババ抜きのシーンが有った
４人でゲームを行い、なんと13枚スタートであった。
これは、相当にレアなのではないか？
Rで計算＆シミュレートしてみる
ある個人が１３枚配られた時に１３枚スタートになる確率を計算する。
勝率をコントロールするため、自分で配って有利な枚数のものを取る
輩がいなければ、この結果を3/4倍すれば求める確率になる。
１３枚以上スタートになる確率だと１４枚スタートになる確率を
足してやる必要がある。
まずは１３枚配られた時に１３枚スタートになる確率
#小数で library(Rmpfr) fact1 <- function(x){ one <- mpfr(1,1000) fact.x <- Reduce("*", c(one, 1:x)) } four <- mpfr(4,1000) a13 <- (four ^ 13 + 13 * four ^ 12) * fact1(13) * fact1(40) / fact1(53) a13 #分数で library(gmp) fact2 <- function(x){ one <- as.bigz(1) fact.x <- Reduce("*", c(one, 1:x)) } four2 <- as.bigz(4) a13.2 <- (four2 ^ 13 + 13 * four2 ^ 12) * fact2(13) * fact2(40) / fact2(53) a13.2
小数
0.0003389767722882 以下略
分数
8388608/24746851955
2950.05464017391204833984375回に１回起こる。
約３０００回に一回なので相当レアだ。
シミュレーションで確認する
tra <- c(rep(1:13, 4), 14) sim = 1000000 res1 <- numeric(sim) for(i in 1:sim){ temp <- sample(tra, 13, replace = FALSE) temp2 <- table(temp) res1[i] <- length(temp2[temp2 == 1 | temp2 == 3]) } sum(res1 > 12) / sim table(res1) hist(res1, freq = FALSE, xlim = c(0,14))
> sum(res1 > 12) / sim
[1] 0.000332
0.0003389767722882 と極めて近いのでOKだろう。
> binom.test(sum(res1 > 12) , sim)
95 percent confidence interval:
0.0002972492 0.0003696965
0.0003389767722882を含んでいる。
次に１４枚スタートになる確率
a14 <- four ^ 13 * fact1(13) * fact1(40) / fact1(53) a14 1 / a14
> ans14
[1] 7.97592405以下略e-5
> 1 / ans14
[1] 12537.73222073以下略
１３枚以上スタートになる確率は
ans <- a13 * 3 / 4 + a14 * 1 / 4 ans 1 /ans
> ans <- a13 * 3 / 4 + a14 * 1 / 4
> ans
[1] 0.00027417238935以下略
> 1 /ans
[1] 3647.34028以下略
ある個人が１３枚以上スタートになるのは、およそ3650回に一回のレアな現象であった。
４人の誰か１人以上に起こる、ということだと
ans2 <- 1 - (1 - a13)^3 * (1 - a14)
ans2
1/ans2
> ans2
[1] 0.00109626以下略
> 1/ans2
[1] 912.189202以下略
900回に一回程度であった。
RMST生存曲線下面積２
#RMST生存曲線下面積２
前回の記事で(surv2sampleComp)パッケージを使ってみた
http://r-statistics-fan.hatenablog.com/entry/2014/08/04/225135
今回survRM2パッケージが出たことを知った。
使ってみる
前回と同じく、イレッサのシミュレーションデータを使用する。
#RMST2 library(survival) set.seed(1) N <- 1000 #全症例 NAK <- N%/%4 #count N NBK <- N%/%4 NCK <- N%/%4 NDK <- N - NAK * 3 AK <- 0.55 #survival rate at AKT AKT <- 8 #x Months survive BK <- 0.2 BKT <- 4 CK <- 0.5 CKT <- 6 DK <- 0.4 DKT <- 6 tinib <- c(rep(1, NAK + NBK), rep(0, NCK + NDK)) mut <- c(rep(1, NAK), rep(0, NBK), rep(1, NCK), rep(0, NDK)) data <- data.frame(tinib = tinib, mut = mut) data$event <- rep(1, N) data$cov <- rep(1, N) AKk <- -log(AK)/AKT BKk <- -log(BK)/BKT CKk <- -log(CK)/CKT DKk <- -log(DK)/DKT AK <- rexp(NAK, rate = AKk) #指数分布で BK <- rexp(NBK, rate = BKk) CK <- rexp(NCK, rate = CKk) DK <- rexp(NDK, rate = DKk) data$sur <- c(AK, BK, CK, DK) ## 分子標的薬ありなし。遺伝子変異考慮なし all <- survfit(Surv(sur, event) ~ 1, data = data, conf.type = "none") tinib_1 <- survfit(Surv(sur, event) ~ 1, subset(data, data$tinib == 1), conf.type = "none") tinib_0 <- survfit(Surv(sur, event) ~ 1, subset(data, data$tinib == 0), conf.type = "none") plot(all, xlim = c(0, 50), ylim = c(0, 1), col = "black") par(new = TRUE) plot(tinib_0, xlim = c(0, 50), ylim = c(0, 1), conf.int = TRUE, col = "red") par(new = TRUE) plot(tinib_1, xlim = c(0, 50), ylim = c(0, 1), conf.int = TRUE, col = "purple") par(new = TRUE) title(main = "遺伝子変異考慮せず。黒線は全員の曲線") legend(25, 1, legend = c("All", "分子標的薬なし", "分子標的薬あり"), col = c("black", "red", "purple"), lty = 1) #######RMST2########## library(survRM2) res.RMST2 <- rmst2(time = data$sur, status = data$event, arm = data$tinib, tau=12) print(res.RMST2) ############ library(surv2sampleComp) # 調整なしRMST res.RMST <- rmstaug(y = data$sur, x = data$cov, arm = data$tinib, delta = data$event, tau = 12, type = "difference") res.RMST
結果
> res.RMST2 <- rmst2(time = data$sur, status = data$event, arm = data$tinib, tau=12) > print(res.RMST2) ##survRM2 result The truncation time: tau = 12 was specified. Restricted Mean Survival Time (RMST) by arm Est. se lower .95 upper .95 RMST (arm=1) 5.346 0.200 4.955 5.738 RMST (arm=0) 6.450 0.189 6.079 6.821 Restricted Mean Time Lost (RMTL) by arm Est. se lower .95 upper .95 RMTL (arm=1) 6.654 0.200 6.262 7.045 RMTL (arm=0) 5.550 0.189 5.179 5.921 Between-group contrast Est. lower .95 upper .95 p RMST (arm=1)-(arm=0) -1.103 -1.643 -0.564 0 RMST (arm=1)/(arm=0) 0.829 0.755 0.910 0 RMTL (arm=1)/(arm=0) 1.199 1.097 1.310 0 ##surv2sampleComp result > res.RMST $result.ini coef se(coef) z p lower .95 upper .95 intercept 6.449593 0.1892441 34.080817 0.000000e+00 6.078674 6.8205112 arm -1.103421 0.2751565 -4.010159 6.067775e-05 -1.642728 -0.5641147 $result.aug coef se(coef) z p lower .95 upper .95 intercept 6.449593 0.1893388 34.063772 0.000000e+00 6.078489 6.8206968 arm -1.103421 0.2752942 -4.008154 6.119524e-05 -1.642998 -0.5638448
Est. lower .95 upper .95 p
RMST (arm=1)-(arm=0) -1.103 -1.643 -0.564 0
coef se(coef) z p lower .95 upper .95
arm -1.103421 0.2752942 -4.008154 6.119524e-05 -1.642998 -0.5638448
結果は同じ感じになった。
