Rプログラミングの小ネタ

tapplyを使えば、データフレームのカテゴリの列でグループ分けしつつ、値の列に対して関数を適用するということができます。 ・・・と言っても分かりにくいと思いますので、具体的なサンプルをあげてみますと、 > # サンプルデータを作る > 名前 <- c("A", "B",  "C", "D",  "E", "F") > 性別 <- c("男", "男", "女", "男", "女", "女" ) > 身長 <- c(175,  165, 165,  170, 160, 155  ) > df <- data.frame(名前, 性別, 身長) > df   名前 性別 身長 1    A   男  175 2    B   男  165 3    C   女  165 4    D   男  170 5    E   女  160 6    F   女  155 > > # ここから本題 > tapply(df$身長, df$性別, mean)    女  男 160 170 性別でグループ分けして、身長の平均を算出する（mean関数を適用）という感じです。 第3引数にhist関数を指定してやれば、グループ分けを適用したあとに、複数の度数分布表を算出したり、複数の度数分布図（ヒストグラム）を書いたりすることもできます。 もう少し大きいサンプルデータを使いましょう。↓この本に載っていたサンプルを使わせていただきます。 ↓こちらからダウンロードできる、"年収.csv"を使います。 データマイニング入門－Rで学ぶ最新データ解析 - 東京図書 ↓こんな感じのサンプルデータです。 > d <- read.csv("年収.csv") > head(d)   学位 地位 性別 現職年数 学位年数  年収 1 博士 教授 男性       25       3

こんな例を考えます。 10人の元々のテストの平均点は50点でした。点数をあげる施策を行った後、再度テストをしてみると、 after <- c(44, 47, 48, 51, 54, 56, 59, 62, 63, 65) という点数が得られました。施策は効果があったかどうかを、有意水準5%で検定しましょう。 【Aさんの検定】 Aさんは施策に懐疑的でした。場合によっては施策の悪影響が出て、 点数が下がる可能性もある と考えていますので、 両側検定 を行います。 > t.test(after, mu=50,  alternative="two.sided")         One Sample t-test data:  after t = 2.1198, df = 9, p-value = 0.06306 alternative hypothesis: true mean is not equal to 50 95 percent confidence interval:  49.67088 60.12912 sample estimates: mean of x      54.9 p値が0.05よりも大きいので、 帰無仮説は棄却できず、「施策は効果があるとは言えない」 と結論付けました。 【Bさんの検定】 Bさんは施策の効果を確信していました。 点数が下がることはありえない と考えていますので、 片側検定 を行います。 > t.test(after, mu=50,  alternative="greater")         One Sample t-test data:  after t = 2.1198, df = 9, p-value = 0.03153 alternative hypothesis: true mean is greater than 50 95 percent confidence interval:  50.66264      Inf sample estimates: mean of x      54.9 p値が0.05よりも小さいので、 帰無仮説は棄却でき、「施策は効果があると言える」 と結論付けました。 う～ん、気持ち悪い。 何が気持ち悪いって、効果に懐

例えば↓こんなデータがあったとします。40人のクラスがあって、1度目のテスト実施後に、なんらかの施策をして、その後、2度目のテストを実施したと。で、施策の効果は点数に表れているのか？ みたいなのを想像していただけると理解しやすいかと。行名の 1, 2, ... , 40 はその学生の出席番号みたいなイメージで。 > df    first second 1     56     75 2     50     59 3     53     61 4     60     62 5     45     59 6     51     60 7     48     70 8     55     79 9     50     67 10    44     50 11    50     64 12    56     49 13    42     53 14    52     53 15    57     41 16    39     63 17    50     54 18    45     54 19    54     59 20    48     74 21    48     58 22    47     59 23    47     57 24    43     54 25    56     76 26    54     39 27    48     59 28    45     43 29    50     59 30    51     52 31    50     58 32    48     70 33    52     61 34    51     50 35    43     62 36    44     39 37    49     64 38    47     45 39    52     68 40    46     64 で、いろいろな可視化をしてみようと思います。 1つ目は直線の傾きで見る感じです。matplot関数をつかえば直線が右肩上がりか、右肩下がりかで効果を把握できます。加えて、色分けすることで分かりやすくしてみました。 # 変化に応じて線の色を決める dif <- df $second - df $first color <

例えば↓こんなサンプルデータがあったとしましょう。属するグループと観測値の組みたいな感じですね。 > x1    group value 1      A  2.66 2      A  2.24 3      A  2.41 4      A  1.12 5      A  2.48 6      A  2.33 7      A  2.90 8      A  4.35 9      A  3.82 10     A  3.45 11     B  6.65 12     B  4.61 13     B  5.26 14     B  5.24 15     B  7.26 16     B  5.97 17     B  2.90 18     B  2.37 19     B  4.94 20     B  4.08 21     C  5.63 22     C  5.75 23     C  4.87 24     C  8.14 25     C  1.01 26     C  1.58 27     C  3.76 28     C  8.12 29     C  5.95 30     C  4.48 で、グループごとの平均やばらつきにどんな傾向があるか見るために、エラーバー付きでプロットしたいような状況だと考えてください。 こういうときは、gplotsパッケージのplotmeans関数が便利です。 install.packages ( "gplots" ) library ( gplots ) plotmeans ( value~ group , data =x1 , connect=F , ylim= c ( 0 , 10 ) ) plotmeansによるエラーバー付きプロット ところで、エラーバーの長さは何を表しているのでしょうか？ > ?plotmeans でマニュアルの引数のところを見ると、 p     confidence level for error bars. Defaults to 0.95. と書いてあります。p= で指定することによって、エラーバーの信頼水準を指定することができ、そのデフォルト値は0.95であると書いてあるので、つ