線形モデルの前提の再考

　これまで様々なパターンの線形モデルを取り扱ってきた。今ここで、検定の時と同じように、線形モデルを使える前提について考えてみよう。線形モデルは、（1）残差の分布が正規分布に従う、（2）残差の分散は説明変数の値によらず一定、（3）説明変数と被説明変数が直線関係（線形）であることを仮定している。これらの仮定が重要であることは、線形回帰にて最小二乗法を使って係数の推定を行っていれば、おのずと理解できる。ゆえに本来は、これらのいずれかの仮定が満たされないときに単純な線形モデルに当てはめて解析するのは問題がある。（3）に関しては以前、非線形最小二乗法でも取り扱った。今回は、さらに（1）や（2）も絡んでくる状況を考えよう。

例1：等分散性が満たされない

　例えば、下記のようなデータを解析することを考える。

------------------------------------------------------

library(plotn)

x1 <- c(6.9, 29.4, 24.6, 25.7, 4.2, 11.3, 26.2, 19.1, 

        15, 18.5, 20.9, 19.2, 13.3, 28.3, 1.3, 12.9, 

        10.7, 18.7, 18.4, 26.4, 5.3, 3.2, 9.1, 26.2, 

        8.5, 26.9, 9.4, 10.5, 0.9, 0.7, 10.3, 15.6, 5,

        9, 23.5, 5.8, 2.4, 27, 16.6, 6.9, 25.8, 1.3,

        13.3, 13.5, 9.9, 18, 21.2, 20.5, 27.4, 24.2, 

        7.1, 18.3, 26.4, 5.7, 14.1, 27.8, 20.3, 17.1, 

        13.6, 29.2, 18.7, 19.9, 15.6, 7.3, 25.2, 6, 9.9, 

        20.8, 20.6, 12.2, 26.7, 18.7, 23.4, 9.1, 10, 

        23.8, 28.2, 13.9, 24.1, 22.6, 5.5, 1.1, 17.1,

        20.4, 20.1, 24.8, 4, 13.5, 10.1, 29.6, 5.4, 12,

        1.5, 23.5, 13.2, 8, 15.1, 6.8, 19.6, 21.8)

y <- c(3.7, 15.1, 5.7, 10.5, 2.4, 4.8, 8, 6.5, 3.5, 12.5,

       11.1, 10.5, 5.9, 16.8, 0, 9.4, 6.6, 13.8, 7.1, 15.6, 

       3.1, 0.8, 5.4, 17, 5.8, 13.2, 6, 7.3, 1.1, 0, 3.4,

       6.2, 2.9, 2.4, 11.9, 1.4, 1.5, 7.1, 9.9, 2.1, 8.3, 

       0.4, 4.5, 5.5, 5.9, 12.7, 6.7, 12.6, 10.1, 12.2, 

       3.5, 12.7, 17.8, 2, 7.1, 13.4, 14.8, 5.7, 8, 18.2,

       18.4, 19.4, 9.1, 1.2, 10.7, 3.9, 10.2, 10.1, 6, 7.3, 

       8.4, 9.1, 7.7, 2.9, 3.8, 10.6, 13.8, 5.3, 16.5, 11,

       3.6, 0.4, 8.6, 7.7, 6.5, 10.4, 1.3, 18.3, 2.4, 16.9, 

       2.9, 3.6, 0, 13.1, 5.9, 3.2, 9, 1.8, 13.7, 15.1)

d <- data.frame(x1 = x1, y = y)

head(d)

##     x1    y

## 1  6.9  3.7

## 2 29.4 15.1

## 3 24.6  5.7

## 4 25.7 10.5

## 5  4.2  2.4

## 6 11.3  4.8

plotn(y ~ x1, data = d) #図1の描画

------------------------------------------------------

明らかにこのデータはxが大きくなるにつれて、分散が大きくなっている。もちろん、このようなデータであっても形式的には線形回帰を行うことが可能である。解析後、併せて回帰曲線を描くことにする。

------------------------------------------------------

fit1 <- lm(y ~ x1, data = d) #線形回帰

fit2 <- summary(fit1)

fit2

## 

## Call:

## lm(formula = y ~ x1, data = d)

## 

## Residuals:

##     Min      1Q  Median      3Q     Max 

## -6.6832 -1.7653 -0.3118  1.2104 11.4369 

## 

## Coefficients:

##             Estimate Std. Error t value Pr(>|t|)    

## (Intercept)  0.14957    0.65350   0.229    0.819    

## x1           0.49730    0.03707  13.414   <2e-16 ***

## ---

## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

## 

## Residual standard error: 3.054 on 98 degrees of freedom

## Multiple R-squared:  0.6474, Adjusted R-squared:  0.6438 

## F-statistic: 179.9 on 1 and 98 DF,  p-value: < 2.2e-16

coefs <- coef(fit2)

xx <- seq(min(d$x1), max(d$x1), length = 200)

yy <- coefs[1] + coefs[2]*xx

plotn(y ~ x1, data = d)

overdraw(points(xx, yy, col = "red", type = "l")) #図2の描画

------------------------------------------------------

最小二乗法に従っているため、データの中央を通るような回帰曲線になっている。実はこのデータ、下記のようにガンマ分布に従うように生成したデータである。詳しい性質は、当該ページを見てもらうとよいが、この分布に従って生成された値は0よりも大きく、平均が大きくなるほど分散も大きくなる。今回は、線形予測子の切片を0とし、xの係数を0.5とした。

------------------------------------------------------

b1 <- 0

b2 <- 0.5

n <- 100

x1 <- runif(n, 0, 30)

y <- rgamma(n = n, shape = b1 + b2 * x1, scale = 1)

x1

##   [1] 10.1047424  9.4317910  2.3436613  5.2843793 17.3159607 15.1307419

##   [7] 12.6825164  5.1162131 28.9944106 18.4355306 19.1081254 28.4372772

##  [13] 18.4272964 16.8662314 23.4053531  2.6083536 13.9791594  6.3959398

##  [19]  8.7019284 19.9495560 10.0785974 26.2398214 22.2744329  0.4358120

##  [25] 20.6069826 10.9650592 15.6620012 17.4292989 14.9723749 12.8403838

##  [31] 19.7399289 24.6924959 12.9843380 20.4821414 28.0146524 13.6872696

##  [37] 10.3006881  0.5593573  7.3278870 23.5961861 17.4800555  0.1462855

##  [43]  7.7375858 29.8039046  1.2969873 27.8417893  9.5622219  7.1520469

##  [49] 24.7085548 22.1263876  0.3407775  4.5830740 24.0108868 11.3252800

##  [55] 24.9733557 17.7208741 21.2518413 14.6085185  4.0727465 17.7910232

##  [61] 23.3867083 21.7513056  3.3739654  5.9572262  2.3404273 24.7383164

##  [67] 16.7341782 17.6337523 29.0084312  8.4603074 14.6759334 25.3527017

##  [73] 15.5824682  5.8881744  7.6075867 20.3449419 17.3244139 18.8521406

##  [79] 19.1027849 15.8654722 29.9780254 19.4873486 14.7773346  0.7174303

##  [85] 27.6160480  9.7934722  0.2793857 25.5542374 25.3807799 16.1705580

##  [91]  9.9202025  6.0556224 26.7562156 20.9985778 17.5914787  8.6946109

##  [97]  7.2574454 18.7914068 21.8030127  6.8573973

y

##   [1] 9.921449e+00 5.770857e+00 1.661283e+00 3.680733e+00 9.838678e+00

##   [6] 5.638701e+00 9.812146e+00 2.558509e+00 1.214186e+01 5.909798e+00

##  [11] 7.626572e+00 9.685944e+00 6.868889e+00 8.688919e+00 1.076761e+01

##  [16] 1.520142e+00 8.409799e+00 3.492446e+00 5.738852e+00 6.230552e+00

##  [21] 7.152797e+00 1.165255e+01 9.101552e+00 1.974275e-01 1.085297e+01

##  [26] 4.726012e+00 7.189375e+00 1.466522e+01 8.552982e+00 5.132835e+00

##  [31] 9.494570e+00 1.002141e+01 5.564657e+00 8.433057e+00 1.858311e+01

##  [36] 4.868214e+00 3.962124e+00 8.695255e-03 1.038696e+00 1.386070e+01

##  [41] 1.048838e+01 8.411577e-02 5.315866e+00 1.922044e+01 9.718462e-01

##  [46] 1.401155e+01 5.021364e+00 5.447531e+00 1.225042e+01 1.266820e+01

##  [51] 5.728057e-05 2.129705e+00 9.759901e+00 6.495971e+00 1.175501e+01

##  [56] 7.679863e+00 1.423744e+01 8.126174e+00 1.357999e+00 6.309421e+00

##  [61] 1.366505e+01 9.347699e+00 3.985169e+00 1.618706e+00 2.805452e+00

##  [66] 8.295418e+00 5.430346e+00 7.873904e+00 1.664956e+01 2.544609e+00

##  [71] 8.966028e+00 1.888579e+01 9.867566e+00 2.521715e+00 4.442281e+00

##  [76] 9.669529e+00 5.709712e+00 7.038061e+00 9.632330e+00 8.086906e+00

##  [81] 1.514704e+01 7.010325e+00 7.069339e+00 2.775649e-01 1.423114e+01

##  [86] 3.401713e+00 1.054347e-01 1.168921e+01 1.008653e+01 3.858754e+00

##  [91] 4.033299e+00 3.893945e+00 1.299766e+01 6.315076e+00 7.530335e+00

##  [96] 4.390782e+00 6.397352e+00 1.133896e+01 7.794367e+00 2.874802e+00

------------------------------------------------------

最小二乗法に従っているため、データの中央を通るような回帰曲線になっている。実はこのデータ、下記のようにガンマ分布に従うように生成したデータである。詳しい性質は、当該ページを見てもらうとよいが、この分布に従って生成された値は0よりも大きく、平均が大きくなるほど分散も大きくなる。今回は、線形予測子の切片を0とし、xの係数を0.5とした。