アールメカブ

_ 2007 08 30

John Verzani p.336 の構造モデルから判断すると，すべてベース (Intercept) との切片の差ということになる．

# 二元配置の分散分析

frogs3<- read.csv(
"http://150.59.18.68/frogs3.csv", 
header = FALSE)

# header = FALSE で，列名はファイルに未設定と指示

# なお列名が未定義の場合，自動的に V1, V2, V3 などの名前が付加される

# 二つの要因がある場合，それらをチルダ記号の右側に ＋ 記号で指定する

frogs3.aov<- aov(V1 ~ V2 + V3, data = frogs3)
summary(frogs3.aov)
summary.lm(frogs3.aov)

# Intercept は V2 = 12H かつ V3 = 100ug の場合．繰り返し数 3

# この標準偏差は sqrt(7.51/6)．これは V2 V3 の自由度の積か

# 2行目の V224H は sqrt(2 *7.51/9)．9 は V2 の繰り返し数か

# Intercept は V2 が 12H で V3 が 100 ug の場合

# 2行目 V224H は V2 が 24H の場合の Intercept(V2=12Hかつ V3=100ugの場合) との差

同じく，p.332 によれば共分散分析では，連続量はスロープを表す．

regrowth <- read.table(
"http://www.bio.ic.ac.uk/research/mjcraw/
therbook/data/ipomopsis.txt", 
header = T)
ancova1<- lm(Fruit ~ Grazing * Root)
summary(ancova1)
anova(ancova1)

共分散分析での各パラメータの標準誤差の計算は Crawlye The R Book p. 492 - 498 また共分散分析での「修正平均」の求めかたは Crawley S-PLUS p.292 中 を参照

# Faraway (2006) よりデータを借用

babyfood <- read.table(file =  
http://150.59.18.68/babyfood.txt")
babyfood

# データから要因別に罹患比率を求めて分割表にする．xtabs() 関数を利用

xtabs(disease/(disease+nondisease) ~ sex + food,  
babyfood)

# ロジスティック回帰分析を実行する

# 目的変数を 2 項分布とした一般化線形モデル glm() による

model1 <- glm(cbind(disease, nondisease) ~  sex + food,  
family = binomial, data = babyfood)

# glm は一般化線形モデルを実行する関数．family は分布を指定する

summary(model1)              # 要約を見る
drop1(model1, test = "Chi")  # 各項は有意か
exp(-.669)                # 母乳の効果を確認する
model.matrix(model1)

# Intercept は Boy で Bottle

# sexGirl は Girl の場合の Intercept との差

# foodBreast は Intercept (Boy Bottle) の場合に比べての差

# foodSuppl は Intercept (Boy Bottle) の場合に比べての差