Last updated: -October 2009. Kajiyama                    [ 目次に戻る ]

線形単回帰分析 A simple linear regression analysis


  1. 線形単回帰分析の関数 lm( )     lm=Linear Modelの略

    単回帰分析は,2変数の散布の直線モデル Y = a + bX の係数 a とb を推定する.

    Y に目的変数の測定値(データ)を与え,X に説明変数の測定値(データ) を与えたときに,

    目的変数の測定値とその予測値の差を最小とする最小二乗法で係数 a と b を推定する.

      Y:目的変数,被説明変数,従属変数,外生変数と呼ぶ.

      X:説明変数,      独立変数,内生変数と呼ぶ.

    Y = f( X ) 説明変数が1つの場合を単回帰分析という.2つ以上は線形重回帰分析を行う.

    単回帰分析の関数は, lm ( Y ~ X )  直線式が当てはまらない場合,非線形(=曲線)回帰分析が用意されている.関数 nls ( )

    質問 「なぜ回帰分析と名づけたのか」調べておくこと.


  2. 「iris」データの「setosa」品種のデータを抜き取り分析に使用する

    (1)「iris」データの確認

    ・「コンソール」 画面に,

    str (iris) #  と記述すると,「iris」のデータ構造が表示され,5変数の名前と,型が表示される.

    ・Sepal.Length: num 数値型 ・Sepal.Width : num 数値型 ・Petal.Length: num 数値型 ・Petal.Width : num 数値型
    ・Species  : Factor 因子型

      あやめの大きな3枚のはなびらは,「Sepal がく片」で,小さな3枚のはなびらが,「Petal 花びら」である.

    Sepal「がく片」の長さと幅・Petal「花びら」の長さと幅

    Species「あやめ3品種 [setosa・versicolor・virginica]」

    (2)「iris」データの「setosa」品種のデータを抜き取り分析に使用

    ・「コンソール」 画面に,

    A <- iris[iris$Species=="setosa", ] #   「setosa」品種を抽出してデータ作成

    と記述すると,「iris」データの第5変数「Species」の「"setosa"」品種,50サンプルが選ばれて,これに新しい名前を付ける

    str( A ) #   データ構造を確認する

    plot (A$Sepal.Width, A$Sepal.Length) #   変数Sepal.Widthと変数Sepal.Lengthの散布図を描く

    > A <- iris[iris$Species=="setosa", ] 
    
    > str( A ) 
    'data.frame':   50 obs. of  5 variables:
     $ Sepal.Length: num  5.1 4.9 4.7 4.6 5 5.4 4.6 5 4.4 4.9 ...
     $ Sepal.Width : num  3.5 3 3.2 3.1 3.6 3.9 3.4 3.4 2.9 3.1 ...
     $ Petal.Length: num  1.4 1.4 1.3 1.5 1.4 1.7 1.4 1.5 1.4 1.5 ...
     $ Petal.Width : num  0.2 0.2 0.2 0.2 0.2 0.4 0.3 0.2 0.2 0.1 ...
     $ Species     : Factor w/ 3 levels "setosa","versicolor",..: 1 1 1 1 1 1 1 1 1 1 ...
    
    > plot (A$Sepal.Width, A$Sepal.Length) 
    
     
  3. 「iris」の「setosa」品種のデータを相関分析

    ・「コンソール」 画面に,

    cor.test( A$Sepal.Width, A$Sepal.Length) #  変数Sepal.Widthと変数Sepal.Lengthの相関分析

    > cor.test( A$Sepal.Width, A$Sepal.Length, method="pearson" ) 
    
            Pearson's product-moment correlation
    
    data:  A$Sepal.Width and A$Sepal.Length 
    t = 7.6807, df = 48, p-value = 6.71e-10           # t値,自由度,p値
    alternative hypothesis: true correlation is not equal to 0   # 対立仮説 ρ(ロー)≠0
    95 percent confidence interval:
     0.5851391 0.8460314            # 95%信頼区間 [0.585, 0.846] 0.5851391 ≦r≦ 0.8460314
    sample estimates:
          cor                       # 相関係数
    0.7425467 
    


    ・相関係数はr = 0.743であり,2変数間に強い相関があることを示している.

    ・相関係数rを用いて母集団相関係数(母相関)ρに関する検定を行う.

     これはt分布を用いた相関係数の検定となる.結果はt値= 7.6807, df・自由度 = 48, p-value・p値 = 6.71e-10 < 0.001

     p値が有意水準0.05(5%)より大きい小さいので母集団の2つの変数「がく片の長さ」と「がく片の幅」の間には有意な相関があるといえる.


  4. 関数 lm( ) による線形単回帰分析目的変数「がく片の長さ」Sepal.Length 説明変数「がく片の幅」Sepal.Width

    (1) 単回帰分析

    「がく片の長さ」を「がく片の幅」で予測することを考える.回帰分析を使用する.

     目的変数「がく片の長さ」Sepal.Length

     説明変数「がく片の幅」Sepal.Width

    ・Rの「コンソール」画面に,

    B <- lm( Sepal.Length ~ Sepal.Width, data=A ) #  と記述し,線形単回帰分析

    B #  と記述し結果を表示

    summary ( B ) #  と記述し,より詳細な結果を表示

    (2) 散布図と回帰直線の作図

    plot ( A$Sepal.Width, A$Sepal.Length, main="Y=0.69X+2.64" ) #  と記述し,目的変数説明変数の散布図を描く

    abline( B, col = "blue", lwd=2 ) #  と記述し,回帰直線を挿入する

    (3) 残差分析

    plot ( 測定値Sepal.Width, residuals( B )) #  と記述し,残差-説明変数の散布図を作成

    > B <- lm( Sepal.Length ~ Sepal.Width, data=A )   #  回帰分析
    
    > B  
    Call:モデル式
    lm(formula = Sepal.Length ~ Sepal.Width, data = A) 
    
    Coefficients:
      切片      傾き
    (Intercept)  Sepal.Width  
         2.6390       0.6905 
    
    > summary ( B )                                 #  回帰分析の要約
    Call:モデル式
    lm(formula = Sepal.Length ~ Sepal.Width, data = A)
    
    Residuals:残差の四分位数 残差は測定値と予測値の差
         Min       1Q   Median       3Q      Max 
    -0.52476 -0.16286  0.02166  0.13833  0.44428 
    
    Coefficients:
           回帰係数 標準誤差    t値       P値
                Estimate  Std. Error  t value   Pr(>|t|)    
    (Intercept)   2.6390     0.3100     8.513     3.74e-11 *** 切片の推定値とその検定結果 帰無仮説は切片が0
    Sepal.Width   0.6905     0.0899     7.681     6.71e-10 *** 傾きの推定値とその検定結果 帰無仮説は傾きが0
    ---
    Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
    
    Residual standard error: 0.2385 on 48 degrees of freedom
    Multiple R-squared(決定係数R2): 0.5514,   Adjusted R-squared(自由度調整済みの決定係数R2): 0.542 
    F-statistic: 58.99 on 1 and 48 DF,  p-value: 6.71e-10 
    
    > plot ( A$Sepal.Width, A$Sepal.Length, main="Y=0.69X+2.63" )      # 測定した2変数の分布図 X軸は説明変数
     
    > abline( B, col = "blue", lwd=2 )              # 推定の回帰直線を描く.lwd=2は回帰直線の幅
    
    > plot ( 測定値Sepal.Width, residuals( B ))     # 残差の分布図 X軸は説明変数
    


    (4) 結果と考察

    Y切片= 2.6390,Xの傾き= 0.6905である.切片と傾きのP値はいずれも0.1%より小さく有意である.

    Y = 0.69X + 2.64 Y 目的変数「がく片の長さ」Sepal.Length X 説明変数「がく片の幅」Sepal.Width

    決定係数: 0.5514, 自由度調整済みの決定係数: 0.542 である.決定係数が1に近いと直線傾向がより強いと判断できる.

    切片を除いたモデル全体の当てはまりの良さの検定結果はF値=58.99であり P値が0.1%より小さく有意である.

    残差は,その散布図によるとほぼ均一なデータの散らばりであり,二次曲線や系統的な分布は観察できない.


  5. 測定値と予測値の確認

    ・Rの「コンソール」画面に,

    測定値Sepal.Length <- A$Sepal.Length

    測定値Sepal.Width <- A$Sepal.Width

    予測値Sepal.Length <- predict( B )

    残差 <- residuals( B )

    C <- data.frame ( 測定値Sepal.Width, 測定値Sepal.Length, 予測値Sepal.Length, 残差 )

    C #  と記述し,データの内容を表示

    > 測定値Sepal.Length <- A$Sepal.Length    
    
    > 測定値Sepal.Width  <- A$Sepal.Width     
    
    > 予測値Sepal.Length  <- predict( B )     
    
    > 残差  <- residuals( B )     
    
    > C <- data.frame ( 測定値Sepal.Width, 測定値Sepal.Length, 予測値Sepal.Length, 残差 )   
    
    > C     # と記述し,データの内容を表示
    
       測定値Sepal.Width	測定値Sepal.Length	予測値Sepal.Length	残差
    1                3.5                5.1					5.055715			 0.04428474
    2                3.0                4.9					4.710470			 0.18952960
    3                3.2                4.7					4.848568			-0.14856834
    4                3.1                4.6					4.779519			-0.17951937
    5                3.6                5.0					5.124764			-0.12476423
    6                3.9                5.4					5.331911			 0.06808885
    7                3.4                4.6					4.986666			-0.38666629
    8                3.4                5.0					4.986666			 0.01333371
           ・・・・・・・・・・・・・・・・・・・・・・
    43               3.2                4.4					4.848568			-0.44856834
    44               3.5                5.0					5.055715			-0.05571526
    45               3.8                5.1					5.262862			-0.16286217
    46               3.0                4.8					4.710470			 0.08952960
    47               3.8                5.1					5.262862			-0.16286217
    48               3.2                4.6					4.848568			-0.24856834
    49               3.7                5.3					5.193813			 0.10618680
    50               3.3                5.0					4.917617			 0.08238268
    

線形単回帰分析                          [ 目次に戻る ]