単回帰分析 - 太郎丸博のホームページ

# 例題 USArrests data を使い、都市人口比率と 3 種類の犯罪の
# 発生率の散布図を作り、OLS で回帰直線を推定しなさい。
# その結果をもとに、都市人口比ともっとも関連が大きいのは、
# どの種類の犯罪か述べなさい。
# また、もっとも関連の大きい犯罪に関して、都市人口比率が 10% 30% 50% 70% 90%
# のときの犯罪率の予測値をもとめなさい。
例題の解答例
都市人口比率と犯罪発生率の関係
太郎丸博
1973 年の米国 50 州の都市人口比率と殺人、暴力犯罪、レイプの発生率の関係を示した
のが、図 1 であり、単回帰分析を行った結果が表１である。図 1 を見ると殺人と暴力犯罪
は都市人口比率とは関係ないように見えるが、レイプに関しては都市人口比率の多い州で
は発生率が高い州もある程度あることが分かる。表 1 を見ると、殺人と暴力犯罪では回帰
直線の傾きは有意ではないが、レイプに関しては１%水準で有意であり、都市人口比率が
１％上がると、レイプの発生率が 0.27 上がるという推測結果である。モデルから予測され
るレイプの発生率は、Y = 3.79 + 0.27 × X で計算できる。その結果を示したのが表 2 で
250
50
150
Assault
10
5
Murder
15
ある。
30 40 50 60 70 80 90
30 40 50 60 70 80 90
UrbanPop
10 20 30 40
Rape
UrbanPop
30 40 50 60 70 80 90
UrbanPop
図 1 1973 年の米国 50 州の都市居住者人口比率と殺人、暴力犯罪、レイプの発生率（人口
10 万人当たりの逮捕件数）の散布図
1
表 1 都市居住者人口比率と殺人、暴力犯罪、レイプの単回帰分析の結果
================================================
lm.Murder lm.Assault lm.Rape
-----------------------------------------------(Intercept)
6.42*
73.08
3.79
(2.91)
(53.85)
(5.71)
UrbanPop
0.02
1.49
0.27**
(0.04)
(0.80)
(0.09)
-----------------------------------------------N
50
50
50
================================================
（）内は標準誤差. ** p < 0.01, * p < 0.05,
表 2 モデルから予測される都市居住者人口比ごとのレイプ発生率
都市居住者人口比
10%
30%
50%
70%
90%
予測値
6.4
11.8
17.1
22.4
27.7
以下は分析のための R スクリプト
USArrests # データの中身を確認
?USArrests # データの解説を見る
attach(USArrests) # スクリプトを簡単に書けるように USArrests 内の変数にパスを通す
par(mfrow=c(2,2)) # 一度に複数のプロットを表示するためにプロット画面を2 行2 列に分割
lm.Murder <- lm(Murder~UrbanPop) # 回帰分析の結果に名前を付ける
plot(UrbanPop, Murder) # 散布図を作成
abline(lm.Murder) # 回帰直線を追加
summary(lm.Murder) # 回帰分析の結果をチェック
lm.Assault <- lm(Assault~UrbanPop)
plot(UrbanPop, Assault)
abline(lm.Assault)
summary(lm.Assault)
lm.Rape <- lm(Rape~UrbanPop)
plot(UrbanPop, Rape)
abline(lm.Rape)
summary(lm.Rape)
# エラーなどのせいで図の表示がずれた時はpar からやりなおすと、リセットされる
library(memisc) #memisc パッケージを読み込む。PC によってインストールされていない場合もある
# mtable 関数で回帰分析の結果を奇麗にまとめる。memisc パッケージを読みこまないと使えない
# 結果をまとめたい回帰分析の結果を指定し、そのあと "digits=2" で小数点以下2 桁まで表示するように指定
mtable(lm.Murder, lm.Assault, lm.Rape, digits=2)
# 予測したい UrbanPop の値を含む data frame を作る
d0 <- data.frame(UrbanPop = 10 + 0:4 * 20)
pred.rape <- round(predict(lm.Rape, d0), 1) # 回帰分析の結果と説明変数の値を含むデータを指定して予測
bind(d0$UrbanPop, pred.rape)
# あとはできた図や表をコピーして、レポートを完成させる。
2

Download Report