Phn tch hi quy tuyn tnh a bin

Phân tích hồi quy tuyến tính đa biến sử dụng phần mềm R Tuan V. Nguyen Garvan Institute of Medical Research, Sydney, Australia

Số liệu Tác động của nhiệt độ (x 1), thời gian (x 2) lên CO 2 (y) Id y X 1 X 2 X 3 1 36. 98 5. 1 400 51. 37 2 13. 74 26. 4 400 3 10. 08 23. 8 4 8. 53 5 X 4 X 5 X 6 X 7 4. 24 1484. 83 2227. 25 2. 06 72. 33 30. 87 289. 94 434. 90 1. 33 400 71. 44 33. 01 320. 79 481. 19 0. 97 46. 4 400 79. 15 44. 61 164. 76 247. 14 0. 62 36. 42 7. 0 450 80. 47 33. 84 1097. 26 1645. 89 0. 22 6 26. 59 12. 6 450 89. 90 41. 26 605. 06 907. 59 0. 76 7 19. 07 18. 9 450 91. 48 41. 88 405. 37 608. 05 1. 71 8 5. 96 30. 2 450 98. 60 70. 79 253. 70 380. 55 3. 93 9 15. 52 53. 8 450 98. 05 66. 82 142. 27 213. 40 1. 97 … … … … … 25 11. 19 11. 5 450 77. 88 25. 20 663. 09 994. 63 1. 61 26 75. 62 5. 2 470 75. 50 8. 66 1464. 11 2196. 17 4. 78 27 36. 03 10. 6 470 83. 15 22. 39 720. 07 1080. 11 5. 88 Chú thích: y = sản lượng CO 2; X 1 = thời gian (phút); X 2 = nhiệt độ (C); X 3 = phần trăm hòa tan; X 4 = lượng dầu (g/100 g); X 5 = lượng than đá; X 6 = tổng số lượng hòa tan; X 7 = số hydrogen tiêu thụ.

Nhập số liệu vào R y <- c(36. 98, 13. 74, 10. 08, 8. 53, 36. 42, 26. 59, 19. 07, 5. 96, 15. 52, 56. 61, 26. 72, 20. 80, 6. 99, 45. 93, 43. 09, 15. 79, 21. 60, 35. 19, 26. 14, 8. 60, 11. 63, 9. 59, 4. 42, 38. 89, 11. 19, 75. 62, 36. 03) x 1 <- c(5. 1, 26. 4, 23. 8, 46. 4, 7. 0, 12. 6, 18. 9, 30. 2, 53. 8, 5. 6, 15. 1, 20. 3, 48. 4, 5. 8, 11. 2, 27. 9, 5. 1, 11. 7, 16. 7, 24. 8, 24. 9, 39. 5, 29. 0, 5. 5, 11. 5, 5. 2, 10. 6) x 2 <- c(400, 450, 450, 400, 425, 450, 450, 460, 450, 470) x 3 <- c(51. 37, 72. 33, 71. 44, 79. 15, 80. 47, 89. 90, 91. 48, 98. 60, 98. 05, 55. 69, 66. 29, 58. 94, 74. 74, 63. 71, 67. 14, 77. 65, 67. 22, 81. 48, 83. 88, 89. 38, 79. 77, 87. 93, 79. 50, 72. 73, 77. 88, 75. 50, 83. 15) x 4 <- c(4. 24, 30. 87, 33. 01, 44. 61, 33. 84, 41. 26, 41. 88, 70. 79, 66. 82, 8. 92, 17. 98, 17. 79, 33. 94, 11. 95, 14. 73, 34. 49, 14. 48, 29. 69, 26. 33, 37. 98, 25. 66, 22. 36, 31. 52, 17. 86, 25. 20, 8. 66, 22. 39) x 5 <- c(1484. 83, 289. 94, 320. 79, 164. 76, 1097. 26, 605. 06, 405. 37, 253. 70, 142. 27, 1362. 24, 507. 65, 377. 60, 158. 05, 130. 66, 682. 59, 274. 20, 1496. 51, 652. 43, 458. 42, 312. 25, 307. 08, 193. 61, 155. 96, 1392. 08, 663. 09, 1464. 11, 720. 07) x 6 <- c(2227. 25, 434. 90, 481. 19, 247. 14, 1645. 89, 907. 59, 608. 05, 380. 55, 213. 40, 2043. 36, 761. 48, 566. 40, 237. 08, 1961. 49, 1023. 89, 411. 30, 2244. 77, 978. 64, 687. 62, 468. 38, 460. 62, 290. 42, 233. 95, 2088. 12, 994. 63, 2196. 17, 1080. 11) x 7 <- c(2. 06, 1. 33, 0. 97, 0. 62, 0. 22, 0. 76, 1. 71, 3. 93, 1. 97, 5. 08, 0. 60, 0. 90, 0. 63, 2. 04, 1. 57, 2. 38, 0. 32, 0. 44, 8. 82, 0. 02, 1. 72, 1. 88, 1. 43, 1. 35, 1. 61, 4. 78, 5. 88) REGdata <- data. frame(y, x 1, x 2, x 3, x 4, x 5, x 6, x 7) Muốn lưu giữ file, cần phải dùng lệnh Write. table(REGdata, ”địa chỉ cần lưu giữ file/data. csv”)

Hoặc gọi file từ R Nhập số liệu từ bảng biểu trình bày trên vào file xcel. Sau đó save as dưới dạng file csv. Cách này có thể lưu giữ số liệu gốc. Rồi gọi file số liệu trưc tiếp vào R: REGdata <- read. csv(“c: /đường dẫn đến nơi lưu file/ data. csv, header=T)

Kết quả phân tích reg <- lm(y ~ x 1+x 2+x 3+x 4+x 5+x 6+x 7, data=REGdata) summary(reg) Residuals: Min 1 Q Median 3 Q Max -20. 035 -4. 681 -1. 144 4. 072 21. 214 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 53. 937016 57. 428952 0. 939 0. 3594 x 1 -0. 127653 0. 281498 -0. 453 0. 6553 x 2 -0. 229179 0. 232643 -0. 985 0. 3370 x 3 0. 824853 0. 765271 1. 078 0. 2946 x 4 -0. 438222 0. 358551 -1. 222 0. 2366 x 5 -0. 001937 0. 009654 -0. 201 0. 8431 x 6 0. 019886 0. 008088 2. 459 0. 0237 * x 7 1. 993486 1. 089701 1. 829 0. 0831. --Signif. codes: 0 '***' 0. 001 '**' 0. 01 '*' 0. 05 '. ' 0. 1 ' ' 1 Residual standard error: 10. 61 on 19 degrees of freedom Multiple R-Squared: 0. 728, Adjusted R-squared: 0. 6278 F-statistic: 7. 264 on 7 and 19 DF, p-value: 0. 0002674

Tương quan chéo

Hồi quy loại từng bước (stepwise) reg <- lm(y ~. , data=REGdata) step(reg, direction=”both”)

Hồi quy loại từng bước Start: AIC= 134. 07 y ~ x 1 + x 2 + x 3 + Df Sum of Sq - x 5 1 4. 54 - x 1 1 23. 17 - x 2 1 109. 34 - x 3 1 130. 90 <none> - x 4 1 168. 31 - x 7 1 377. 09 - x 6 1 681. 09 x 4 + x 5 + x 6 + x 7 RSS AIC 2145. 37 132. 13 2164. 00 132. 36 2250. 18 133. 42 2271. 74 133. 68 2140. 83 134. 07 2309. 14 134. 12 2517. 92 136. 45 2821. 92 139. 53 Step 2: AIC= 130. 42 y ~ x 2 + x 3 + x 4 + x 6 + x 7 Df Sum of Sq RSS - x 2 1 96. 8 2264. 9 - x 3 1 122. 0 2290. 0 <none> 2168. 1 - x 4 1 187. 4 2355. 5 + x 1 1 22. 7 2145. 4 + x 5 1 4. 1 2164. 0 - x 7 1 385. 0 2553. 1 - x 6 1 1526. 2 3694. 3 AIC 129. 6 129. 9 130. 4 130. 7 132. 1 132. 4 132. 8 142. 8 Step 1: AIC= 132. 13 y ~ x 1 + x 2 + x 3 + x 4 + x 6 Df Sum of Sq RSS - x 1 1 22. 7 2168. 1 - x 2 1 113. 8 2259. 1 - x 3 1 133. 5 2278. 9 <none> 2145. 4 - x 4 1 170. 8 2316. 2 + x 5 1 4. 5 2140. 8 - x 7 1 375. 7 2521. 1 - x 6 1 1058. 5 3203. 8 + x 7 AIC 130. 4 131. 5 131. 8 132. 1 132. 2 134. 1 134. 5 141. 0 Step 3: AIC= 129. 59 y ~ x 3 + x 4 + x 6 + x 7 Df Sum of Sq RSS - x 3 1 25. 4 2290. 3 - x 4 1 90. 9 2355. 8 <none> 2264. 9 + x 2 1 96. 8 2168. 1 + x 5 1 8. 3 2256. 5 + x 1 1 5. 7 2259. 1 - x 7 1 384. 9 2649. 7 - x 6 1 2015. 6 4280. 5 AIC 127. 9 128. 7 129. 6 130. 4 131. 5 131. 8 144. 8

Hồi quy loại từng bước Step 4: AIC= 127. 9 y ~ x 4 + x 6 + x 7 Df Sum of Sq - x 4 1 73. 5 <none> + x 3 1 25. 4 + x 1 1 11. 3 + x 5 1 6. 3 + x 2 1 0. 3 - x 7 1 486. 6 - x 6 1 1993. 8 RSS 2363. 8 2290. 3 2264. 9 2279. 0 2284. 0 2290. 0 2776. 9 4284. 1 AIC 126. 7 127. 9 129. 6 129. 8 129. 9 131. 1 142. 8 Call: lm(formula = y ~ x 6 + x 7, data = REGdata) Coefficients: (Intercept) x 6 x 7 2. 52646 0. 01852 2. 18575 Step 5: AIC= 126. 75 y ~ x 6 + x 7 Df Sum of Sq RSS <none> 2363. 8 + x 4 1 73. 5 2290. 3 + x 1 1 33. 4 2330. 4 + x 3 1 8. 1 2355. 8 + x 5 1 7. 7 2356. 1 + x 2 1 7. 3 2356. 6 - x 7 1 497. 3 2861. 2 - x 6 1 4477. 0 6840. 8 AIC 126. 7 127. 9 128. 4 128. 7 129. 9 153. 4

Phương pháp lựa chọn mô hình tối ưu bằng Bayesian Model Average (BMA) library(BMA) xvars <- REGdata[, -1] co 2 <- REGdata[, 1] bma <- bicreg(xvars, co 2, strict=FALSE, OR=20) summary(bma)

Bayesian Model Average (BMA) analysis 16 models were selected Best 5 models (cumulative posterior probability = 0. 6599 ): p!=0 EV SD model 1 model 2 model 3 Intercept 100. 0 5. 75672 14. 6244 2. 5264 6. 1441 8. 6120 x 1 12. 4 -0. 01807 0. 1008. . . x 2 10. 4 -0. 00075 0. 0282. . . x 3 10. 7 0. 00011 0. 0791. . . x 4 20. 2 -0. 03059 0. 1020. . -0. 1419 x 5 10. 5 -0. 00023 0. 0030. . . x 6 100. 01815 0. 0040 0. 0185 0. 0193 0. 0164 x 7 73. 7 1. 60766 1. 2821 2. 1857. 2. 1628 model 4 7. 5936 -0. 1393. . 0. 0162 2. 1233 model 5 7. 3537. . -0. 0572. . 0. 0179 2. 2382 n. Var r 2 BIC post prob 3 0. 704 -22. 9721 0. 072 3 0. 701 -22. 6801 0. 063 2 0. 700 -25. 8832 0. 311 1 0. 636 -24. 0238 0. 123 3 0. 709 -23. 4412 0. 092