在R中線性回歸分析的函數是lm()。
(1)一元線性回歸
我們可以根據以上數據來分析合金的強度是否與碳含量有關系。
首用以下命令把數據讀取到R中:
x <- c(seq(0.10,0.18,by = 0.01),0.20,0.21,0.23)
y <- c(42.0,43.5,45.0,45.5,45.0,47.5,49.0,53.0,50.0,55.0,55.0,60.0)
plot(x,y)
通過畫圖得到想x,y兩個變量之間存在某種線性關系
所以,就可以用lm()函數來擬合直線,通過回歸函數lm()得到如下結果:
lm.sol <- lm(y~1+x) ##lm()函數返回擬合結果的對象,可以用summary()函數查看其內容。
summary(lm.sol)
回歸結果:
Call:
lm(formula = y ~ 1 + x)
Residuals:
Min 1Q Median 3Q Max
-2.0431 -0.7056 0.1694 0.6633 2.2653
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 28.493 1.580 18.04 5.88e-09 ***
x 130.835 9.683 13.51 9.50e-08 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 1.319 on 10 degrees of freedom
Multiple R-squared: 0.9481, Adjusted R-squared: 0.9429
F-statistic: 182.6 on 1 and 10 DF, p-value: 9.505e-08
結果分析:
其中,兩個回歸系數分別是28.493和130.835,結果里面還有t值,以及兩個P值,P值越小,回歸效果越顯著,並且后面星級越高。
倒數第二行R-squared數字越接近於1,回歸效果越好。
所以,本例回歸分析效果顯著,回歸直線為:y=28.493+130.835x。
做完回歸分析之后還可以進行做預測,也就是說給定一個x值,可以求出Y值的概率為0.95的相應區間。在R中可以用predict()函數實現:
> new <- data.frame(x = 0.16) ##注意,一個值也要寫出數據框的形式。
> lm.pred <- predict(lm.sol,new,interval = "prediction",level = 0.95) ##加上interval = "prediction",表示同時給出相應的預測區間。
> lm.pred
fit lwr upr
1 49.42639 46.36621 52.48657
(2)多元線性回歸
首先把數據讀取到R中:
> x1 <- c(76.0,91.5,85.5,82.5,79.0,80.5,74.5,79.0,85.0,76.5,82.0,95.0,92.5)
> x2 <- c(50,20,20,30,30,50,60,50,40,55,40,40,20)
> y <- c(120,141,124,126,117,125,123,125,132,123,132,155,147)
> mydata <- data.frame(x1,x2,y)
做線性回歸:
> lm.sol <- lm(y ~ x1 + x2,data = mydata)
> summary(lm.sol)
結果:
Call:
lm(formula = y ~ x1 + x2, data = mydata)
Residuals:
Min 1Q Median 3Q Max
-4.0404 -1.0183 0.4640 0.6908 4.3274
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -62.96336 16.99976 -3.704 0.004083 **
x1 2.13656 0.17534 12.185 2.53e-07 ***
x2 0.40022 0.08321 4.810 0.000713 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 2.854 on 10 degrees of freedom
Multiple R-squared: 0.9461, Adjusted R-squared: 0.9354
F-statistic: 87.84 on 2 and 10 DF, p-value: 4.531e-07
從以上結果來看,回歸系數和回歸方程的檢驗都是顯著的,所以回歸方程為:y=-62.96366+2.13656x1+0.40022x2
與一元回歸一樣,也可以用predict()函數做預測:
下面來預測以為體重80KG,年齡40歲的男子的血壓:
> new <- data.frame(x1 = 80,x2 = 60)
> lm.pred <- predict(lm.sol,new ,interval = "prediction", level = 0.95)
> lm.pred
fit lwr upr
1 131.9743 124.6389 139.3096