一般化線形混合モデル のRによる実施例です。
重回帰分析 や ロジスティック回帰分析 は、単独で紹介されることが一般的ですが、 一般化線形モデル の一種です。
下記は、glmを使っているところは同じですが、familyとlinkの部分が違います。 下記の中には、デフォルトのlinkなのでlinkを省略しても大丈夫なものもありますが、ここでは省略しないで書いています。
任意の位置について、予測値を計算したい時は、Data2.csvというファイルを別に作り、そこに任意の位置を記入します。 Data2.csvのファイルに含まれる変数名は、Data.csvと同じにする必要があります。 Yの変数はあってもなくても良く、あっても使われません。
library(MASS)
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data1 <- read.csv("Data.csv", header=T) # データを読み込み
gm <- step(glm(Y~., data=Data1, family= gaussian(link = "identity"))) # 一般化線形モデル
summary(gm) # 結果の出力
#ここから予測の手順
library(ggplot2)
Data2 <- read.csv("Data2.csv", header=T) # データを読み込み
s2 <- predict(gm,Data2)# 作業用ディレクトリを変更
Data2s2 <- cbind(Data2,s2)# データを読み込み
ggplot(Data2s2, aes(x=X01, y=X02)) + geom_point(aes(colour=s2)) + scale_color_viridis_c(option = "D")
ポアソン回帰分析は、目的変数が度数（頻度、カウントデータ）の時によく使われます。
library(MASS)
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data1 <- read.csv("Data.csv", header=T) # データを読み込み
gm <- step(glm(Y~., data=Data1, family= poisson(link = "log"))) # 一般化線形モデル
summary(gm) # 結果の出力
#ここから予測の手順
library(ggplot2)
Data2 <- read.csv("Data2.csv", header=T) # データを読み込み
s2 <- predict(gm,Data2,type="response")# 作業用ディレクトリを変更
Data2s2 <- cbind(Data2,s2)# データを読み込み
ggplot(Data2s2, aes(x=X01, y=X02)) + geom_point(aes(colour=s2)) + scale_color_viridis_c(option = "D")
対数線形分析 は、ポアソン回帰分析の一種ですが、一般的に二次以上の項も見ます。 Rによる対数線形分析 に実施例があります。
ロジスティック回帰分析では、Yは「0]と「1」の２種類の数字になっている必要があります。 Yが「OK」と「NG」のような質的変数の場合、２種類しかなくてもエラーになります。 この場合、それぞれを0と1に変換しておきます。
library(MASS)
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data1 <- read.csv("Data.csv", header=T) # データを読み込み
gm <- step(glm(Y~., data=Data1, family= binomial(link = "logit"))) # 一般化線形モデル
summary(gm) # 結果の出力
#ここから予測の手順
library(ggplot2)
Data2 <- read.csv("Data2.csv", header=T) # データを読み込み
s2 <- predict(gm,Data2,type="response")# 作業用ディレクトリを変更
Data2s2 <- cbind(Data2,s2)# データを読み込み
ggplot(Data2s2, aes(x=X01, y=X02)) + geom_point(aes(colour=s2)) + scale_color_viridis_c(option = "D")
「type="response"」を入れると、予測が0から1までの値になります。 Yが「1」のサンプルの発生確率と解釈することもできますし、0と1の期待値と解釈することもできます。
線形混合モデル の実施例です。
事前に「lme4」というパッケージをインストールして、読み込んでおきます。
下記は、コピーペーストで、そのまま使えます。 この例では、Cドライブの「Rtest」というフォルダに、 「Data.csv」という名前でデータが入っている事を想定しています。
データは、１列目は「Y」という列名で目的変数、２列目は「X1」という列名で説明変数、３列目は「C1」という列名でカテゴリが入っていることを想定しています。
setwd("C:/Rtest") # 作業用ディレクトリを変更
Data <- read.csv("Data.csv", header=T) # データを読み込み
library(lme4)
lmer <- lmer(Y ~ X1 + (1 + X1|C1) + (1|C1), data=Data) # 線形混合モデル（傾きと切片に変量効果あり）
summary(lmer) # モデルの作成結果の出力
切片や、傾きの片方だけでに、変量効果がある場合も作れます。
lmer <- lmer(Y ~ X1 + (1|C1), data=Data) # 線形混合モデル（切片に変量効果あり）
lmer <- lmer(Y ~ X1 + (1 + X1|C1), data=Data) # 線形混合モデル（傾きに変量効果あり）
R Tips
http://cse.naro.affrc.go.jp/takezawa/r-tips/r/72.html
上記以外の、familytとlinkについても紹介しています。
同志社大学 金明哲先生のページ
一般化線形モデル
https://www1.doshisha.ac.jp/~mjin/R/Chap_16/16.html
Tweet