R言語を用いて回帰分析を行なってみます。回帰分析は連続的な目的変数を説明変数から求めることで一つの説明変数に対する予測を単回帰、複数の説明変数を用いて行う方法を重回帰といいます。今回は回帰問題でよく使われるボストンの住宅価格のデータを使って行なってみます。

1. ボストン住宅価格データの確認

まず今回使用するボストン住宅価格データの中身を見てみます。

#データを取得するためのライブラリ

library(MASS)



#データの取得表示

data <- Boston

head(data)

dim(data)

説明変数の意味を下に記載します。

この一番下の「mev:住宅価格」を目的変数として使います。

また、各変数の情報はsummary()関数を用いて確認ができます。

library(MASS)

data <- Boston

summary(data)

2. 相関係数と線形回帰について

次に回帰を実装するために変数を選びます。変数を選ぶには変数同士の相関関係がなくてはなりません。まず、相関係数と線形単回帰について軽く説明します。

相関係数は-1~1の値を取ります。上の図のように相関係数が正の場合一方の変数が増えるともう一方が増える性質を持ち、負の場合一方が増えるともう一方が減っていくことを示します。また、相関係数の絶対値が1に近づくほどその関係が強くなっていき、逆に０に近いほど一方の変数に対するもう一方の関係がないということになります。つまり相関係数の絶対値が大きいほど上図のような直線での近似がしやすくなり、線形回帰分析が行えるということになります。

では、住宅の価格であるmevと相関関係がないと線形回帰分析が行えないため相関係数を調べます。

#相関係数行列の計算

library(MASS)

data <- Boston

correlation <- cor(data)



#相関係数行列の可視化

library(corrplot)

corrplot(correlation, method="shade", shade.col=NA, tl.col="black", tl.srt=45, 

               addCoef.col="black", order="AOE")

cor()関数を用いて相関係数行列を作成し、corrplotライブラリのcorrplot()利用して相関係数行列をヒートマップとして表示しました。これをみるとmedv:住宅価格と最も相関係数が高いのがlstat:低所得者人口の割合であり、最も直線で近似できそうであることがわかります。

3. Rを使って線形単回帰の実装

3.1 線形単回帰について

簡単に線形単回帰について説明すると一つの目的変数と一つの説明変数の関係を直線の式で近似することで以下の式で近似します。

この a と b を求めることが線形単回帰です。説明変数と目的変数の相関係数の絶対値が大きいほど近似がうまくいきます。線形単回帰については以下の記事を参考にしていただくとより理解が深まると思います。

3.2 Rでの線形単回帰を実装

では、上の相関関係の分析結果からlstat:低所得者人口の割合がmedv:住宅価格と相関関係が強いことがわかったためこれを用いて線形回帰を行いたいと思います。R言語では線形回帰は非常に簡単に行えます。

result <- lm(medv~lstat, data=data)

このように書きます。lm()の第一引数に目的変数と説明変数を”~”を使って指定します。”~”の左側に目的変数、右側に説明変数を書きます。重回帰の場合は右側の説明変数を”+”を使って追加します。今回は単回帰なので一つの説明変数lstatを指定しています。

まず、データセットをtrainデータとtestデータにわけます。その後trainデータを元に線形回帰分析を行い、結果を用いてtestデータの予測を行います。

#データの取得

library(MASS)

data <- Boston



#回帰分析

result <- lm(medv~lstat, data=data)



#結果の表示

result

結果を表示してみました。

1つ目の(Intercept) が切片で の に対応し、2つ目のlstatが傾きに対応しており、 の に対応しています。

さらに詳しくみたい場合は以下のようにsummary()を使うといいです。

summary(result)

3.3 近似した直線とデータのプロット

線形単回帰を行ったので結果をグラフにプロットしてみます。

#データの取得

library(MASS)

data <- Boston



####データの分割#####

#　データの順番をシャッフル

set.seed(42)

shuffle <- sample(1:506)

data <- data[shuffle,]


# trainデータとtestデータの分類

split_data <- 1:10

train <- data[-split_data, ]
test <- data[split_data, ]


#####################



#回帰分析

result <- lm(medv~lstat, data=train)



#グラフの描画

lmpredict <- predict(result)

plot(train$medv, train$lstat, xlim=c(-10,60), xlab="lstat", ylab="medv")

par(new=T)

plot(lmpredict, train$lstat,type="l", col="red",xlim=c(-10,60),xlab="lstat", ylab="medv")

(グラフの描画ではRのデフォルトの関数を使用して描画しました。par(new=T)で直線のグラフとデータプロットを重ねています。Rでの詳しい可視化については以下の記事で紹介しているのでよかったら参考にしてみてください)

可視化してみるとある程度線形回帰で予測できることがわかりました。

3.4 テストデータで予測

最後にテストデータに対して予測を行い結果を確認してみます。

#データの取得

library(MASS)

data <- Boston



####データの分割#####

#　データの順番をシャッフル

set.seed(42)

shuffle <- sample(1:506)

data <- data[shuffle,]


# trainデータとtestデータの分類

split_data <- 1:10

train <- data[-split_data, ]
test <- data[split_data, ]


#####################



#回帰分析

result <- lm(medv~lstat, data=train)

predictdata <- predict(result, test)



#予測データの表示

cbind(predictdata,data.frame(test$medv))

一番左側はデータ番号で真ん中のデータが予測データ(predictdata)、右側のデータが正しいテストデータ(test.medv)です。値が結構異なる場所もあり、完全とは言えませんがある程度予測はできています。相関係数が±1に近いですが単回帰で完全に予測しようとするのは難しそうです。

4. まとめ

Rには統計処理やデータ分析を行うツールが揃っており、このようなちょっとした回帰分析はすぐに実装できるようになっています。今回は回帰分析のRを使ったやり方を紹介しました。簡単なので遊んでみるといいですね。