■google

■最近のコメント
■最近のトラックバック
■最近の記事
■月別アーカイブ
■ブログランキング
■ブログ検索

■ブロとも申請フォーム
■リンク
■RSSフィード
スポンサーサイト
上記の広告は1ヶ月以上更新のないブログに表示されています。
新しい記事を書く事で広告が消せます。


スポンサー広告 | --:--:--
ステップワイズ法による変数選択
重回帰分析を行うときには
変数の選択が大きな問題になります。

そのときには、
意味を良く考えて
自分の考えを反映させて、
変数を設定するべきです。

しかし、ここでは
変数の選び方に一つの示唆を
与えるものとして、
ステップワイズ法
stepAIC
を用いた変数選択を紹介しましょう。

library(MASS)で
使用できます。


まず、変数選択の方法としては
変数増加法
変数減少法
変数増減法があります。


変数増加法は
まず、定数のみの回帰式を設定し、
そこにAICを最も改善させる
一つ変数を追加させて行きます。
AICが改善しなくなったらそこで
ストップです。

result<-lm(y~1)
result2<-stepAIC(result,direction="forward",scope=list(upper=~x1+x2+x3+x4))

定数項は1で表します。
そこに加える変数を
scope=listで指定、
変数増加法はdirection="forward"と表します。
ここで、x1+x2+x3+x4でなく、
x1*x2*x3*x4
とすれば、交互作用項もOKです。
(ただし、あんまりたくさんの変数で行うと重いです)



変数減少法は
まず、変数を全て取り込んだ回帰式を設定し、
そこにAICを最も改善させる
一つ変数を取り除いて行きます。
AICが改善しなくなったらそこで
ストップです。

direction="backward"

変数増減法は
まず、定数のみの回帰式を設定し、
そこにAICを最も改善させる
一つ変数を加えます。
次に変数を一つ加えるか、
一つ取り除くかをAICが
改善するように追加あるいは削除します。
AICが改善しなくなったらそこで
ストップです。

direction="both"


オススメ文献




スポンサーサイト
変数選択 | 22:48:31 | Trackback(0) | Comments(0)
単回帰分析
まずは、
やっぱり単回帰分析から。
使用するデータは
以下のコマンドで得られる車のデータ。

>library(rpart)
>data(car.test.frame)


目的変数に燃費、説明変数に重量として単回
帰分析を例に挙げる。

>result<-lm(Mileage~Weight)
>summary(result)


lmが回帰のコマンドです。
lm(目的変数~説明変数)

もし、重回帰なら、
説明変数のところを
説明変数1+説明変数2+・・・
というように
"+"でつないでください。

summaryは回帰の結果の
要約です。

その結果は、
summary

(クリックすると大きくなります)

用語の日本語訳は、
Residuals:残差
Estimate:係数
Std.error:標準誤差
t value:t値
Pr:有意確率
Multiple R-Squared:決定係数
Ajusted R-Squared:修正済み決定係数

これを解釈して・・・
決定係数 0.72
回帰式
Y=48.3(24.43)-0.0082(-12.89 )X +ε
ただし()内はt値、
εは誤差項(独立性,等分散性,正規性を仮定)

という感じ。

決定係数72%なので説明力のあるモデル
と言えます。

回帰式を描画するには、
>plot(Weight,Mileage)
>abline(result)


回帰


さて、ここまで
誤差項に独立性,等分散性,正規性
という仮定をおいてきました。

回帰式が出来た後に
必ず確認しましょう。

視覚的に確認するには、
残差プロットと
正規QQプロットが良く使われます。

残差プロット
plot(resid(result))
残差プロット
(クリックで大きくなります)

これは均一になっているかを見ます。
どうでしょうか。
少し2次関数に見えるかもしれませんが、
均一と言えそうです。

正規Q-Qプロット
qqnorm(resid(result))
qqline(resid(result))

QQ

(クリックで大きくなります)
直線に近いほど、
正規分布に従っています。
今回はそう言えそうですね。


残差プロットと
QQプロットで
回帰診断を行いましたが、
このとき、残差の仮定に
沿わないときは、
層別の必要や
何か異常値を含んでいます。
そういうときは、
元データに戻って
じっくり見てみましょう。

例えば、車のデータでは
プレミアムの付いた高級車が
あるかもしれません。

Rではユーザに優しく、
そういったデータは
プロットしたときに
データの番号が表示されています。
(今回のQQプロットのように)




回帰分析 | 22:46:22 | Trackback(0) | Comments(0)
ボックスプロット
ヒストグラムで用いたデータを
引き続き使う。

車のタイプ(Small, Sporty, Compact, Medium)別に
価格の分布を比べる。

まず、データを分ける作業を行う。
subset()
は部分集合を抽出する命令。

subset(データ名,条件)

>Small<-subset(x,Type=="Small")[,1]
>Sporty<-subset(x,Type=="Sporty")[,1]
>Compact<-subset(x,Type=="Compact")[,1]
>Medium<-subset(x,Type=="Medium")[,1]


ボックスプロットを
描くコマンドは
boxplot()

>boxplot(Small,Sporty,Compact,Medium,names=c("Small",
"Sporty","Compact","Medium"),ylab="価格($)")


namesでラベルの指定を行っている。
(デフォルトでは1,2,3・・・とでてくる)
ylabはy軸の名前。


boxplot

(クリックすると大きくなります)


参考:
ボックスプロットの読み方や
意味合いについて書いてある.


データの視覚化 | 15:05:17 | Trackback(0) | Comments(0)
Rの性質
Rを使う上で
最低限理解しておくべきことがあります。

それは・・・

Rでは保存されている関数やデータをオブジェクトと呼ぶ。

Rは基礎部分とパッケージ部分によって構成されている。

関数およびデータを分野別に分類してパッケージとしてまとめられている。

実行した結果をオブジェクトとして保存するには“<-”を使う。
>オブジェクト<-関数(引数)





Rの性質 | 10:03:51 | Trackback(0) | Comments(0)
ヒストグラム
ヒストグラムはhist(x)を
使って描きます。

ここでは、
rpartパッケージに保存されている、
自動車のデータを使って紹介します。

>library(rpart)
>data(car.test.frame)

を打つとデータを得ることができます。

このデータの中から、
priceに相当するものを抜き出します。

>Price<-car.test.frame$Price

これでPriceを得ました。


それではヒストグラムを書いてみましょう。

>hist(price,breaks=c(0,5000,10000,15000,20000,25000),   
 xlab="価格",ylab="度数",col="gray")


breaksは階級の幅を決めるものです。
区切り位置を指定しましょう。

xlabはx軸の名前、
ylabはy軸の名前です。
colは柱の色です。

結果です。

hist.jpg



データの視覚化 | 09:53:33 | Trackback(0) | Comments(0)
代表値
Rでは以下の関数を使って基本統計量を求めることができる。

代表値

(クリックすると大きくなります)


データの要約 | 17:25:13 | Trackback(0) | Comments(0)

FC2Ad

上記広告は1ヶ月以上更新のないブログに表示されています。新しい記事を書くことで広告を消せます。