多元線性回歸分析 (Multiple Linear Regression)

套路43: 多元線性回歸分析 (Multiple Linear Regression)

1. 使用時機: 以多個獨立自變項預測一個應變項。
2. 分析類型: 母數分析(parametric analysis)。
3. 範例資料: 咪路測量菌菌的生長條件與菌數，資料如下表。咪路希望能得到一回歸方程式可用來預測菌菌生長。
item
Temp
pH
hour
ml
CFU(106)
1
6
5.7
1.6
2.12
9.9
2
1
6.4
3
3.39
9.3
3
-2
5.7
3.4
3.61
9.4
4
11
6.1
3.4
1.72
9.1
5
-1
6
3
1.8
6.9
6
2
5.7
4.4
3.21
9.3
7
5
5.9
2.2
2.59
7.9
8
1
6.2
2.2
3.25
7.4
9
1
5.5
1.9
2.86
7.3
10
3
5.2
0.2
2.32
8.8
11
11
5.7
4.2
1.57
9.8
12
9
6.1
2.4
1.5
10.5
13
5
6.4
3.4
2.69
9.1
14
-3
5.5
3
4.06
10.1
15
1
5.5
0.2
1.98
7.2
16
8
6
3.9
2.29
11.7
17
-2
5.5
2.2
3.55
8.7
18
3
6.2
4.4
3.31
7.6
19
6
5.9
0.2
1.83
8.6
20
10
5.6
2.4
1.69
10.9
21
4
5.8
2.4
2.42
7.6
22
5
5.8
4.4
2.98
7.3
23
5
5.2
1.6
1.84
9.2
24
3
6
1.9
2.48
7
25
8
5.5
1.6
2.83
7.2
26
8
6.4
4.1
2.41
7
27
6
6.2
1.9
1.78
8.8
28
6
5.4
2.2
2.22
10.1
29
3
5.4
4.1
2.72
12.1
30
5
6.2
1.6
2.36
7.7
31
1
6.8
2.4
2.81
7.8
32
8
6.2
1.9
1.64
11.5
33
10
6.4
2.2
1.82
10.4
4. 執行回歸。
第一步: 輸入建立資料。
  m <- read.table(header = T, text = "
Temp pH hour ml CFU
6 5.7 1.6 2.12 9.9
1 6.4 3 3.39 9.3
-2 5.7 3.4 3.61 9.4
11 6.1 3.4 1.72 9.1
-1 6 3 1.8 6.9
2 5.7 4.4 3.21 9.3
5 5.9 2.2 2.59 7.9
1 6.2 2.2 3.25 7.4
1 5.5 1.9 2.86 7.3
3 5.2 0.2 2.32 8.8
11 5.7 4.2 1.57 9.8
9 6.1 2.4 1.5 10.5
5 6.4 3.4 2.69 9.1
-3 5.5 3 4.06 10.1
1 5.5 0.2 1.98 7.2
8 6 3.9 2.29 11.7
-2 5.5 2.2 3.55 8.7
3 6.2 4.4 3.31 7.6
6 5.9 0.2 1.83 8.6
10 5.6 2.4 1.69 10.9
4 5.8 2.4 2.42 7.6
5 5.8 4.4 2.98 7.3
5 5.2 1.6 1.84 9.2
3 6 1.9 2.48 7
8 5.5 1.6 2.83 7.2
8 6.4 4.1 2.41 7
6 6.2 1.9 1.78 8.8
6 5.4 2.2 2.22 10.1
3 5.4 4.1 2.72 12.1
5 6.2 1.6 2.36 7.7
1 6.8 2.4 2.81 7.8
8 6.2 1.9 1.64 11.5
10 6.4 2.2 1.82 10.4")  # 資料間以空白分隔
第二步: 使用基本模組(base)的attach及names函數賦予m中的資料名稱(標題)。
  attach(m)  # 讓R能透過變數名稱搜尋資料。
  names(m)  # 賦予資料名稱(標題)。
第三步: 使用基本模組(base)的lm函數代入m資料求回歸方程式，結果儲存至變數x。
  x <- lm(CFU ~ Temp + pH + hour + ml, data = m)
  # CFU因變數，Temp、pH、hour及ml為自變數。
  # CFU ~ Temp + pH + hour + ml即回歸方程式的組成:
     CFU = 係數*Temp + 係數*pH + 係數*hour + 係數*ml
第四步: 顯示判讀結果。
  coefficients(x)
 # 以coefficients 函數顯示求得之回歸方程式係數(model coefficients)。
(Intercept)        Temp        pH        hour         ml
 13.9229254   0.1116533  -0.9924694   0.3241469  -0.2109881

  summary(x)
  # 以summary 函數顯示求得之回歸方程式係數及相關資料p值。
Call:
lm(formula = CFU ~ Temp + pH + hour + ml, data = m)
Residuals:
   Min    1Q  Median    3Q    Max
-2.2849 -1.0231  0.1202  0.9970  2.5673
Coefficients:
                Estimate Std. Error t value Pr(>|t|)  
(Intercept)  13.9229     4.1799   3.331  0.00244 **
Temp          0.1117     0.1065   1.048  0.30364  
pH           -0.9925     0.6695  -1.482  0.14939  
hour          0.3241     0.2534   1.279  0.21129  
ml           -0.2110     0.6321  -0.334  0.74102  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 1.42 on 28 degrees of freedom
Multiple R-squared:  0.2001,    Adjusted R-squared:  0.08578
F-statistic: 1.751 on 4 and 28 DF,  p-value: 0.167

  confint(x, level = 0.95)
  # 以confint函數求信賴區間(CIs for model parameters)。
                   2.5 %     97.5 %
(Intercept)    5.3606934 22.4851573
Temp        -0.1066003  0.3299069
pH          -2.3638398  0.3789010
hour         -0.1948710  0.8431647
ml          -1.5057751  1.0837989
   
  anova(x)
  # 顯示anova table
Analysis of Variance Table
Response: CFU
              Df Sum Sq Mean Sq F value  Pr(>F) 
Temp       1  7.631  7.6311  3.7822 0.06191 .
pH         1  2.924  2.9245  1.4495 0.23870 
hour       1  3.348  3.3480  1.6594 0.20823 
ml         1  0.225  0.2248  0.1114 0.74102 
Residuals 28 56.494  2.0176                 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

  vcov(x)
  # covariance matrix for model parameters
            (Intercept)         Temp           pH        hour          ml
(Intercept)  17.4719516 -0.120085075 -2.544985069  0.30277089 -1.08175295
Temp         -0.1200851  0.011352466 -0.005645347 -0.01089930  0.05271535
pH           -2.5449851 -0.005645347  0.448205285 -0.03677360  0.01105980
hour          0.3027709 -0.010899303 -0.036773602  0.06419957 -0.08129613
ml           -1.0817529  0.052715348  0.011059800 -0.08129613  0.39954353

  # 其他相關資料可以下列函數取得。
  fitted(x) # predicted values
  residuals(x) # residuals
  influence(x) # regression diagnostics

5. 挑選變數，重新求回歸方程式。
第一步: 建立初始模型CFU~1及全模型CFU~. (全模型就是用上所有的變數)。
  null <- lm(CFU~1, data = m)
  null
Call:
lm(formula = CFU ~ 1, data = m)
Coefficients:
(Intercept) 
      8.885
  full <- lm(CFU~., data = m)
  full
Call:
lm(formula = CFU ~ ., data = m)
Coefficients:
(Intercept)         Temp          pH         hour          ml 
    13.9229       0.1117      -0.9925       0.3241      -0.2110 
第二步: 使用基本模組(base)的step函數代入初始模型null及全模型full以forward addition方法挑選變數。
  step(null, scope = list(lower = null, upper = full), direction = "forward")
Start:  AIC=27.11
CFU ~ 1
       Df Sum of Sq    RSS    AIC
+ Temp  1    7.6311 62.991 25.334
<none>              70.622 27.108
+ ml    1    3.1752 67.447 27.590
+ hour  1    2.2965 68.326 28.017
+ pH    1    1.4950 69.127 28.402
Step:  AIC=25.33
CFU ~ Temp
       Df Sum of Sq    RSS    AIC
<none>              62.991 25.334
+ pH    1   2.92447 60.067 25.765
+ hour  1   1.88749 61.104 26.330
+ ml    1   0.12047 62.871 27.271
Call:
lm(formula = CFU ~ Temp, data = m)
Coefficients:
(Intercept)         Temp 
     8.3186       0.1271
第三步: 使用基本模組的step函數以backward elimination方法挑選變數。
  step(full, data = m, direction = "backward")
Start:  AIC=27.74
CFU ~ Temp + pH + hour + ml
       Df Sum of Sq    RSS    AIC
- ml    1    0.2248 56.719 25.873
- Temp  1    2.2156 58.710 27.011
- hour  1    3.3021 59.796 27.616
<none>              56.494 27.742
- pH    1    4.4341 60.928 28.235
Step:  AIC=25.87
CFU ~ Temp + pH + hour
       Df Sum of Sq    RSS    AIC
- hour  1     3.348 60.067 25.765
<none>              56.719 25.873
- pH    1     4.385 61.104 26.330
- Temp  1     8.928 65.647 28.697
Step:  AIC=25.77
CFU ~ Temp + pH
       Df Sum of Sq    RSS    AIC
- pH    1    2.9245 62.991 25.334
<none>              60.067 25.765
- Temp  1    9.0606 69.127 28.402
Step:  AIC=25.33
CFU ~ Temp
       Df Sum of Sq    RSS    AIC
<none>              62.991 25.334
- Temp  1    7.6311 70.622 27.108
Call:
lm(formula = CFU ~ Temp, data = m)
Coefficients:
(Intercept)         Temp 
     8.3186       0.1271 
第四步: 使用基本模組的step函數以stepwise regression方法挑選變數。
  step(null, scope = list(upper = full), data = m, direction = "both")
Start:  AIC=27.11
CFU ~ 1
       Df Sum of Sq    RSS    AIC
+ Temp  1    7.6311 62.991 25.334
<none>              70.622 27.108
+ ml    1    3.1752 67.447 27.590
+ hour  1    2.2965 68.326 28.017
+ pH    1    1.4950 69.127 28.402
Step:  AIC=25.33
CFU ~ Temp
       Df Sum of Sq    RSS    AIC
<none>              62.991 25.334
+ pH    1    2.9245 60.067 25.765
+ hour  1    1.8875 61.104 26.330
- Temp  1    7.6311 70.622 27.108
+ ml    1    0.1205 62.871 27.271
Call:
lm(formula = CFU ~ Temp, data = m)
Coefficients:
(Intercept)         Temp 
     8.3186       0.1271 
# 步驟二到四三種挑選模型的方法都挑出CFU ~ Temp為最佳模型。
第五步: 比較兩組回歸方程式統計上有無差別(就是省略某些變數不影響)。
  y1 <- lm(CFU ~ Temp + pH + hour + ml, data = m)
  y2 <- lm(CFU ~ Temp, data = m)
  anova(y1, y2)
第六步: 判讀結果。
Analysis of Variance Table
Model 1: CFU ~ Temp + pH + hour + ml
Model 2: CFU ~ Temp
  Res.Df    RSS Df Sum of Sq      F Pr(>F)
1     28 56.494                          
2     31 62.991 -3   -6.4973 1.0734  0.3763
  # Pr < 0.05，H0: 兩個方程式效果相當，不成立。
  # Pr > 0.05，H0: 兩個方程式效果相當，成立。

來勁了嗎? 想知道更多?? 補充資料(連結):
1. 關於Regression analysis (https://en.wikipedia.org/wiki/Regression_analysis)
2. 關於Linear regression (https://en.wikipedia.org/wiki/Linear_regression)
3. STHDA Linear Regression Essentials in R (http://www.sthda.com/english/articles/40-regression-analysis/165-linear-regression-essentials-in-r/)
4. 關於R基礎，R繪圖及統計快速入門:
   a. R Tutorial: https://www.tutorialspoint.com/r/index.htm
   b. Cookbook for R: http://www.cookbook-r.com/
   c. Quick-R: https://www.statmethods.net/
   d. Statistical tools for high-throughput data analysis (STHDA): http://www.sthda.com/english/
e. The Handbook of Biological Statistics: http://www.biostathandbook.com/
f. An R Companion for the Handbook of Biological Statistics: http://rcompanion.org/rcompanion/index.html
5. Zar, JH. 2010. Biostatistical Analysis, Fifth Edition, Pearson.

留言

張貼留言

這個網誌中的熱門文章

統計不球人 目錄 (Table of Contents)

用 R 做生醫資料科學   (Biomedical Data Science using R) 如何選擇統計方法 1 如何選擇統計方法  2 ***************************   第一段 : 起手式   ****************************** 套路 0: 安裝 R 不球人 套路 1: R 算統計不球人 ， 使用 R 函數 套路 2: 使用 ggplot2 繪圖 (Graphics using ggplot2) 套路 3: 敘述性統計 (Descriptive Statistics) 套路 4: 如何檢查資料是否接近常態分布 (Check Normality using R) 套路 5: 機率分布 (Probability Distributions) ***************************   第二段 : 比較平均值、中數   *************************** 套路 6: 單一樣本 Z 檢定 (One-Sample Z test) 套路 7: 單一樣本 t 檢定 (One-Sample t test) 套路 8: 單一樣本 Wilcoxon signed-rank 檢定 ( One-Sample Wilcoxon signed-rank test ， non-parametric ) 套路 9: 單一樣本符號檢定 (One Sample Sign test) 套路 10: 兩組獨立樣本 Z 檢定 (Two-Sample Z Test) 套路 11: 兩組變異數相同獨立樣本 t 檢定 ( t Test for Two-Samples with Equal Variances) 套路 12: 兩組變異數不同樣本 t 檢定 ( t Test for Two-Samples with Unequal Variances ， Welch Two Sample t test) 套路 13: 成對樣本 t 檢定 (Paired Sample t test) 套路 14: 兩組獨立樣本曼恩 - 惠尼 U 檢定 (...
閱讀完整內容

如何選擇統計方法 1

如何選擇統計方法 1 第一步 : 選擇統計方法可 依照統計方法用途的分類 來選擇統計方法 。 第二步 : 要了解確認統計方法的用途可看他的虛無假設 (Null hypothesis) H 0 。 以下為統計方法用途分類及本站網頁聯結 。 類別 1. 對資料進行基本分析 ， 以了解資料的大概情況 (summary) ， 例如平均值 、 中位數 、 四分位數 、 有無離群值 (outlier) 及是否為常態分布等 … 敘述性統計 (Descriptive Statistics) 如何檢查資料是否接近常態分布 (Normality Test using R) 類別 2. 比較平均值 (mean) 或中位數 (median) 2.1. 單一樣本 ( 與期望值或特定值相比 ) 單一樣本 Z 檢定 (One-Sample Z test) 單一樣本 t 檢定 (One-Sample t test) 單一樣本 Wilcoxon signed-rank 檢定 ( One-Sample Wilcoxon signed-rank test) 單一樣本符號檢定 (One Sample Sign test) 2.2. 兩組獨立樣本 ( 兩組獨立樣本相比 ) 兩組獨立樣本 Z 檢定 (Two-Sample Z Test) 兩組變異數相同樣本 t 檢定 (t Test for Two-Samples with Equal Variances) 兩組變異數不同樣本 t 檢定 (t Test for Two-Samples with Unequal Variances ， Welch Two Sample t-test) 兩組獨立樣本曼恩 - 惠尼 U 檢定 (Two-Sample Mann-Whitney U Test) 比較兩組獨立次序資料 (Two-Sample Test for Ordinal Data) 兩組獨立樣本 Mood's 中位數檢定 (Two-Sample Mood's Median Test) 2.3. 兩組成對樣本 ( 兩組成對樣本相比 ) 成對樣本 t 檢定 (Pai...
閱讀完整內容

單因子多樣本中位數差異檢定 (Kruskal-Wallis test)

圖片
套路 21: 單因子多樣本中位數差異檢定 (Kruskal-Wallis test) 變異數分析適用連續資料之差異分析，若自變項有一個，就是單因子變異數分析。 多樣本中位數差異檢定是對應 單因子變異數分析的無母數分析方法。當數據不符合單因子變異數分析的前提假設 ( 相同變異數，數據為常態分布 ) 時，可選擇 多樣本中位數差異檢定來做為替代分析方法 。 ** 注意 ** ” 比較變異數 ” (comparing variances) 和 ” 變異數分析 ” (analysis of variance) 不同。變異數分析是多組資料比較平均值或中位數。統計假設檢定檢定什麼 ? 看 H 0 。例如多組獨立樣本假設檢定 H 0 : k 組樣本中位數相同， H A : 至少有一組中位數不同 ， 是檢定多組資料的中位數是否相同。假設相等時為雙尾 (two-tailed test) 檢定。 1. 使用時機 : 用於比較多組樣本中位數 (median) 。 2. 分析類型 : 無母數分析 (non-parametric analysis) 。 直接使用資料數值算統計叫 parametric 方法，把資料排序之後用排序的名次算統計叫 non-parametric 方法。 3. 假設前提 : 無。 4. 範例資料 : 咪路調查餵食不同飼料的肉雞體重 (g) ， 資料如下： 飼料 1 飼料 2 飼料 3 飼料 4 61.8 78.8 70.5 60.3 65.1 79.5 72.6 63.8 61.7 76.0 71.7 64.1 63.3 73.4 72.0 61.4 77.3 71.1 60.9 不同飼料是否效果不同 ? H 0 : 餵食不同飼料的肉雞體重 ( 中位數 ) 相同。 H A : 餵食不同飼料的肉雞體重 ...
閱讀完整內容