区間高次元化回帰分析 のRによる実施例です。
この実施例は、
区間高次元化回帰分析
のページと、ほぼ同じ手順をRで実行しています。
Rによる一般化線形混合モデル にあるような一般化線形混合モデルや、交互作用モデルを使うと、コードはシンプルにする方法もあります。 しかし、この方法を使うと、xの元データの項がモデルに入って来て、係数の評価がわかりにくくなります。
library(MASS)
library(dummies)
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data <- read.csv("Data.csv", header=T) # データを読み込み
Data1 <- Data
DataY <- Data
DataY$X <- NULL
Data1$Y <- NULL
DataX <- Data1
Data1[,1] <- droplevels(cut(Data1[,1], breaks = 3,include.lowest = TRUE))# 1列目の量的変数を１次元クラスタリングで３つカテゴリに分ける
Data2 <- dummy.data.frame(Data1)# ダミー変換
Data3 <-Data2*DataX[,1]# 交互作用項の作成
colnames(Data3)<-paste0("X:",colnames(Data3))# 変数名の修正
Data4 <-cbind(DataY,Data2,Data3)
gm <- step(glm(Y~., data=Data4, family= gaussian(link = "identity"))) # 重回帰分析
summary(gm) # 結果の出力
#ここから予測の手順
library(ggplot2)
s2 <- predict(gm,Data4)# 作業用ディレクトリを変更
Data4s2 <- cbind(Data4,s2)# データを読み込み
ggplot(Data4s2, aes(x=Y, y=s2)) + geom_point() + labs(x="Y",y="predicted Y")
この実施例は、 クラスター高次元化回帰分析 のページのデータを使っています。 クラスター分析の部分がEXCELではできないところですが、Rなら簡単にできます。
ダミー変数と、元の説明変数との交互作用項を作るのに、 線形混合モデル を使っています。 こうすると、コードがシンプルになります。
ここでは、クラスター分析としてk-means法を使っています。
他の方法は、
Rによるクラスター分析
のページが参考になります。
library(MASS)
library(dummies)
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data <- read.csv("Data.csv", header=T) # データを読み込み
Data10 <- Data
Data10$Y <- NULL# Yの列を削除
Data11 <- Data10 # 出力先の行列を作る
for (i in 1:ncol(Data10)) {
Data11[,i] <- (Data10[,i] - min(Data10[,i]))/(max(Data10[,i]) - min(Data10[,i]))# 正規化
}
km <- kmeans(Data11,3) # k-means法で分類。これは3個のグループ分けの場合
cluster <- km$cluster # 分類結果の抽出
cluster <- as.character(cluster)# 文字列に変換
cluster <- as.data.frame(cluster)# データフレームに変換
cluster <- dummy.data.frame(cluster)# ダミー変換
Data4 <-cbind(Data,cluster)# データセットの作成
gm <- step(glm(Y~.^2, data=Data4, family= gaussian(link = "identity"))) # 線形混合モデル
summary(gm) # 結果の出力
このコードだと、モデルに元の説明変数が入るため、結果の解釈がちょっと面倒です。
結果にcluster3の場合が入っていないですが、cluster3の単独の係数はX01やX02という所の係数に表れています。
cluster1とcluster2が少々面倒で、cluster1のX01の係数というのは、X01とX01:cluster1というのの和です。 この場合だと、前者が約-2で、後者が約4なので、「約2 （-2+4)」というのが、知りたい値になります。
#ここから予測の手順
library(ggplot2)
s2 <- predict(gm,Data4)# 作業用ディレクトリを変更
Data4s2 <- cbind(Data4,s2)# データを読み込み
ggplot(Data4s2, aes(x=Y, y=s2)) + geom_point() + labs(x="Y",y="predicted Y")
ほぼ一直線に並ぶので、非常に精度が高く予測できていることが、わかります。
Tweet
