首先,当我们要创建一串虚拟随机数的时候,我们必须总要记得写set.seed(n)。这样做,随机数生成器总能产生同等数目的数据。
set.seed(20)
预测变量q:使用seq来快速产生等间距的序列:
q <- seq(from=0, to=20, by=0.1)
预测y值:
y <- 500 + 0.4 * (q-10)^3
我们现在产生一些噪音并把它添加到模型中:
noise <- rnorm(length(q), mean=10, sd=80)
noisy.y <- y + noise
对噪声数据进行画图:
plot(q,noisy.y,col='deepskyblue4',xlab='q',main='Observed data')
lines(q,y,col='firebrick1',lwd=3)
下面的这个图根据观测数据进行模拟。其中,模拟的图的散点是蓝色的,而红色线则是信号(信号是一种术语,它通常用于表示我们感兴趣的东西的通常变化趋势)。
我们得出的模型应当是 y = aq + bq2 + c*q3 + cost。
现在,我们用R对此进行模拟。要拟合一个多项式模型,你也可以这样用:
model <- lm(noisy.y ~ poly(q,3))
或者:
model <- lm(noisy.y ~ x + I(X^2) + I(X^3))
然而,我们要知道q,I(q^2),I(q^3)存在相关的关系,而这些相关变量很有可能引起某些问题的产生。这时,使用poly()可以避免这个问题,因为它是创建一个垂直的多项式。因此,我喜欢第一种方法:
summary(model)
Call:
lm(formula = noisy.y ~ poly(q, 3))
Residuals:
Min 1Q Median 3Q Max
-212.326 -51.186 4.276 61.485 165.960
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 513.615 5.602 91.69 <2e-16 ***
poly(q, 3)1 2075.899 79.422 26.14 <2e-16 ***
poly(q, 3)2 -108.004 79.422 -1.36 0.175
poly(q, 3)3 864.025 79.422 10.88 <2e-16 ***
Signif. codes: 0 ‘’ 0.001 ‘’ 0.01 ‘’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 79.42 on 197 degrees of freedom
Multiple R-squared: 0.8031, Adjusted R-squared: 0.8001
F-statistic: 267.8 on 3 and 197 DF, p-value: 0
我们可以使用confint()来获得一个模型的参数的置信区间。
一下是模型参数的置信区间:
confint(model, level=0.95)
2.5 % 97.5 %
(Intercept) 502.5676 524.66261
poly(q, 3)1 1919.2739 2232.52494
poly(q, 3)2 -264.6292 48.62188
poly(q, 3)3 707.3999 1020.65097
现在,我们要作一个拟合VS残差图。如果这是一个拟合效果比较不错的模型,我们应该看不到任何一种模型的模式特征:
plot(fitted(model),residuals(model))
整体来说,这个模型的拟合效果还是不错的,毕竟残差为0.8。第一和第三个订单序列的系数,在统计学当中,是相当这样的,这样在我们的意料之中。现在,我们可以使用predict()函数来获得拟合数据以及置信区间,这样,我们可以不按照数据来作图。
下面是预测值和预测置信区间:
predicted.intervals <- predict(model,data.frame(x=q),interval='confidence',level=0.99)
在已有的图像中添加拟合线:
lines(q,predicted.intervals[,1],col='green',lwd=3)
lines(q,predicted.intervals[,2],col='black',lwd=1)
lines(q,predicted.intervals[,3],col='black',lwd=1)
添加图例:
legend("bottomright",c("Observ.","Signal","Predicted"),
col=c("deepskyblue4","red","green"), lwd=3)
下面是它的拟合图像:
我们可以看到我们的模型在数据的拟合方面做的不错,我们也因此感到非常满意。
注意:多项式回归是一种更能强大的工具。可是,我们也可能得到事与愿违的结果:在这个例子中,我们知道我们的信号是使用三次多项式而产生的,然而,当我们在分析实际数据的时候,我们通常对此不知情,因此,正因为多项式次数n大于4的时候会产生过度拟合的情况,我们要在这里注意一下。但你的模型取了噪音而不是信号的时候会产生过拟合的情况;甚至,当你在现有的数据进行模型优化的时候,当你要尝试预测新的数据的时候就不好了,它会导致缺失值的产生。
版权声明:本文内容由阿里云实名注册用户自发贡献,版权归原作者所有,阿里云开发者社区不拥有其著作权,亦不承担相应法律责任。具体规则请查看《阿里云开发者社区用户服务协议》和《阿里云开发者社区知识产权保护指引》。如果您发现本社区中有涉嫌抄袭的内容,填写侵权投诉表单进行举报,一经查实,本社区将立刻删除涉嫌侵权内容。