![]() |
ノート/R/回帰分析https://pepper.is.sci.toho-u.ac.jp:443/pepper/index.php?%A5%CE%A1%BC%A5%C8%2FR%2F%B2%F3%B5%A2%CA%AC%C0%CF |
![]() |
ノート/ノート
ノート/R
訪問者数 2125 最終更新 2016-01-11 (月) 13:42:50
線形での依存性を見つけたい。Excelだと説明変数が16個までで具合が悪い。
手順を忘れるので、メモしておく。
(あとで)
まず、作業ディレクトリの変更が必要なら、作業ディレクトリの変更
> setwd("ディレクトリ指定")
CSVファイルの読み込みは、
> x0 <- read.csv("全体検討.csv")
で読み込む。ファイルからデータを読み込む、 Rメモ:read.csv/write.csvで読み込んだり書き出したりするときの注意点
この時の元データは、ExcelからCSV形式で書き出したもの。
デフォルトで1行目をタイトル行(列タイトルの指定)と見なすが、漢字もOK。 但し、漢字を含む場合は、Linux上のRを使うのであれば、CSVファイル上でUTF-8に 変換しておくのがよい。(でないと後で引っかかる)
行番号(行タイトル)は、今回は指定しないで上から順に打たせる。行タイトルを 含む欄を指定したい場合は、row.names="列名"かrow.names=列番号を指定する。Rメモ:read.csv/write.csvで読み込んだり書き出したりするときの注意点
列名の変更は、行列に対してもデータフレームに対しても、
> colnames(dat) <- c("A","B","C","D")
のような形で良いらしい。行列やデータフレームの列名・変数名の変更
> x1 <- as.data.frame(x0)
で変換する。データ型とデータ構造
念のため、x1の内容を確認しておくこと。
要は、「各列が変量となっているデータ行列 x の各変量の単位が異なる場合に,各変量を平均が 0 ,分散が 1 になるように変換すること.規準化または標準化という.」ということ。
で、やり方は、
> x2 <- scale(x1)
あとは、回帰分析の関数lmに食わせればよいだけ。
> x3 <- lm(目的変数~説明変数, data=データフレーム)
重回帰分析(複数の説明変数の線形結合で説明したい)の場合は、説明変数の指定を
> x3 <- lm(目的変数~説明変数1+説明変数2+... , data=データフレーム)
のように書けばよい。R による重回帰分析
まず、womenのデータ自体は、Rを起動すると変数womenとして既にあるので、そのまま使ってしまおう。
> women height weight 1 58 115 2 59 117 3 60 120 4 61 123 5 62 126 6 63 129 7 64 132 8 65 135 9 66 139 10 67 142 11 68 146 12 69 150 13 70 154 14 71 159 15 72 164
次に標準化する
> w0 <- scale(women) > w0 height weight [1,] -1.5652476 -1.4022687 [2,] -1.3416408 -1.2732255 [3,] -1.1180340 -1.0796608 [4,] -0.8944272 -0.8860962 [5,] -0.6708204 -0.6925315 [6,] -0.4472136 -0.4989668 [7,] -0.2236068 -0.3054021 [8,] 0.0000000 -0.1118374 [9,] 0.2236068 0.1462489 [10,] 0.4472136 0.3398136 [11,] 0.6708204 0.5978998 [12,] 0.8944272 0.8559861 [13,] 1.1180340 1.1140723 [14,] 1.3416408 1.4366802 [15,] 1.5652476 1.7592880 attr(,"scaled:center") height weight 65.0000 136.7333 attr(,"scaled:scale") height weight 4.472136 15.498694
このままlmに食わせると、データフレームでない、というエラーが出るので、データフレームに変換する。
> w1 <- scale(w0) > w1 height weight 1 -1.5652476 -1.4022687 2 -1.3416408 -1.2732255 3 -1.1180340 -1.0796608 4 -0.8944272 -0.8860962 5 -0.6708204 -0.6925315 6 -0.4472136 -0.4989668 7 -0.2236068 -0.3054021 8 0.0000000 -0.1118374 9 0.2236068 0.1462489 10 0.4472136 0.3398136 11 0.6708204 0.5978998 12 0.8944272 0.8559861 13 1.1180340 1.1140723 14 1.3416408 1.4366802 15 1.5652476 1.7592880
これで、lmに食わせると、
> w2 <- lm(weight~height, data=w1) > w2 Call: lm(formula = weight ~ height, data = w1) Coefficients: (Intercept) height 8.597e-16 9.955e-01
更に、summaryを出力すると
> summary(w2) Call: lm(formula = weight ~ height, data = w1) Residuals: Min 1Q Median 3Q Max -0.11184 -0.07312 -0.02473 0.04785 0.20109 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 8.597e-16 2.541e-02 0.00 1 height 9.955e-01 2.630e-02 37.85 1.09e-14 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 0.0984 on 13 degrees of freedom Multiple R-squared: 0.991, Adjusted R-squared: 0.9903 F-statistic: 1433 on 1 and 13 DF, p-value: 1.091e-14
次に、漢字のカラム名が動作するかどうかを確認しておきたい。
> colnames(w0) <- c("身長","体重")
でカラム名を書き換えられるのだが(うまく書き換えられる)、どうも裏に英語のカラム名 が潜んでいそうなので、全面的にデータを作り直してみる。まずこの書き変えたデータ(行列)をcsvファイルに書き出し、それを再度読み込んでみる。
> write.csv(w1, "women-j.csv", quote=TRUE, row.names=TRUE) > w11 <- read.csv("women-j.csv",header=TRUE, row.names=1) > w11 身長 体重 1 -1.5652476 -1.4022687 2 -1.3416408 -1.2732255 3 -1.1180340 -1.0796608 4 -0.8944272 -0.8860962 5 -0.6708204 -0.6925315 6 -0.4472136 -0.4989668 7 -0.2236068 -0.3054021 8 0.0000000 -0.1118374 9 0.2236068 0.1462489 10 0.4472136 0.3398136 11 0.6708204 0.5978998 12 0.8944272 0.8559861 13 1.1180340 1.1140723 14 1.3416408 1.4366802 15 1.5652476 1.7592880
これでよさそうなので、データフレームに変換
> w21 <- as.data.frame(w11) > w21 身長 体重 1 -1.5652476 -1.4022687 2 -1.3416408 -1.2732255 3 -1.1180340 -1.0796608 4 -0.8944272 -0.8860962 5 -0.6708204 -0.6925315 6 -0.4472136 -0.4989668 7 -0.2236068 -0.3054021 8 0.0000000 -0.1118374 9 0.2236068 0.1462489 10 0.4472136 0.3398136 11 0.6708204 0.5978998 12 0.8944272 0.8559861 13 1.1180340 1.1140723 14 1.3416408 1.4366802 15 1.5652476 1.7592880
これでlmに食わせる
> w31 <- lm(体重~身長,data=w21) > w31 Call: lm(formula = 体重 ~ 身長, data = w21) Coefficients: (Intercept) 身長 1.663e-15 9.955e-01 > summary(w31) Call: lm(formula = 体重 ~ 身長, data = w21) Residuals: Min 1Q Median 3Q Max -0.11184 -0.07312 -0.02473 0.04785 0.20109 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 1.663e-15 2.541e-02 0.00 1 身長 9.955e-01 2.630e-02 37.85 1.09e-14 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 0.0984 on 13 degrees of freedom Multiple R-squared: 0.991, Adjusted R-squared: 0.9903 F-statistic: 1433 on 1 and 13 DF, p-value: 1.091e-14
これで、lmの引数に、漢字のコラム名を直接書いてOKらしいことが分かる。
いよいよ、目的とする「社会的強み」の各項目と成績(単位数)との重回帰分析をしてみる。
データの読み込みは、csvファイルからとし、UTF-8に変換して置くのを忘れないこと。
> x0 <- read.csv("全体検討.csv")
あとで標準化(scale)するときに、学生番号は具合が悪い(平均が取れない?)ので、覗いてしまう。
x0[,2:(ncol(x0))]
更に、データは部分的に穴があいている(欠落)と具合が悪いらしいので、na.omitで取り除く。その始末はRに任せてしまった。要確認。
na.omit(x0[,2:(ncol(x0))])
その上で、標準化(scale)
> x1 <- scale(na.omit(x0[,2:(ncol(x0))]))
データフレームに変換
> x2 <- as.data.frame(x1)
これで、回帰分析lmを適用する。
> x3 <- lm(単位数 ~ 社会的強み平均(偏差値)+ ... ,data=x2)
ここでエラー。
エラー: 想定外の入力です in "x3 <- lm(単位数 ~ 社会的強み平均・"
どうも、欄名 「社会的強み平均(偏差値)」 のカッコ (偏差値) がいけないらしい。 面倒なので(どうせこの平均との相関はあまりないことが分かっているので)説明変数から外してしまおう、ということで、やってみると今度はOK。
> x3 <- lm(単位数 ~ 意欲_偏差値+ ... ,data=x2)
で、結果が得られた。この先は結果を参照