1

我正在尝试计算相当分散的数据的百分位数或分位数。使用 Loess 函数可以很好地呈现平均值，但是，我无法从此函数中获得百分位数/分位数。

我试图将 quantreg 与黄土结合起来。此图显示线性曲线而不是黄土平滑曲线。

我想得到类似这样的结果： 在此处输入图像描述

data(cars)
plot(cars)
lmodel <- loess(cars$dist~cars$speed,span = 0.3, degree = 1)
lpred<-predict(lmodel, newdata= 5:25,se=TRUE)
lines(5:25, lpred$fit,col='#000066',lwd=4)
lines(5:25, lpred$fit - qt(0.975, lpred$df)*lpred$se, lty=2)
lines(5:25, lpred$fit + qt(0.975, lpred$df)*lpred$se, lty=2)


#### combination of quantreg with loess

plot(cars$speed,cars$dist)
xx <- seq(min(cars$speed),max(cars$speed),1)
f <- coef(rq(loess(cars$dist~cars$speed,span = 0.3, degree = 1), tau=c(0.1,0.25,0.5,0.75,0.9)) )
yy <- cbind(1,xx)%*%f
for(i in 1:length(taus)){
  lines(xx,yy[,i],col = "gray")
}

我也尝试了建议的代码，但是，我无法更改平滑的设置。线条显示波浪形路径。

library(quantreg)
data(cars)
taus <- c(0.1, 0.25, 0.5, 0.75, 0.9)
lmodel <- loess(dist ~ speed, data = cars, span = 0.9, degree = 1)
rqmodel <- rq(lmodel, tau = taus, data = cars)
f <- coef(rqmodel)
xx <- seq(min(cars$speed), max(cars$speed), length.out = nrow(cars))
yy <- predict(rqmodel)
plot(cars)
matlines(xx, yy, col = "grey",lwd=3)

在此处输入图像描述

Loess 函数不像 rg 那样提供分位数数据。

然而，黄土函数允许获得没有锯齿形的曲线。请查看代码片段。与 Loess 函数相比，使用 rg 函数产生相同结果的 tau=0.5 的设置是什么。

data(cars)
lmodel <- loess(dist ~ speed, data = cars, span = 0.9 )
plot(cars)
lines( x=4:25 , y=predict(lmodel, newdata= data.frame(speed=4:25)) ,col="Blue")

在此处输入图像描述

4

2 回答 2

1

我相信问题中的代码混合了黄土和分位数回归，因为它们是不同的方法，而后者不需要前者。

我将尝试同时拟合两者并绘制各自的结果。在下面的代码中，我将使用matlines，而不是for循环。

这些代码行很常见。

library(quantreg)

data(cars)

xx <- seq(min(cars$speed), max(cars$speed), length.out = nrow(cars))

首先是黄土模型。 

lmodel <- loess(dist ~ speed, data = cars, span = 0.5, degree = 1)
ls_yy <- predict(lmodel, se = TRUE)
ls_yy <- cbind(ls_yy$fit, 
               ls_yy$fit - 2*ls_yy$se.fit, 
               ls_yy$fit + 2*ls_yy$se.fit)

plot(cars)
matlines(xx, ls_yy, col = "darkgrey")

在此处输入图像描述

现在分位数回归。 

taus <- c(0.1, 0.25, 0.5, 0.75, 0.9)
rqmodel <- rq(dist ~ speed, tau = taus, data = cars)

rq_yy <- predict(rqmodel)

plot(cars)
matlines(xx, rq_yy, col = "darkgrey")

在此处输入图像描述

于 2019-06-01T17:12:08.110 回答
1

下面的代码（取自“答案”）不正确，不应包含在正确的解决方案中。这将为拟合提供 95% 的置信区间，以及区间落在真实趋势线上的概率。它不对应于根据该移动平均线范围内的数据计算的分位数。推荐的基于正态的近似值需要将 ls_yy$se.fit 乘以 sqrt(ni)，其中 ni 是特定跨度中的观察数。不幸的是 loess 不返回 ni，所以这不是一个站得住脚的解决方案，除非跨度覆盖整个数据集并且 ni 可以设置为 n 并且不存在异方差。

data(cars)
plot(cars)

lmodel <- loess(dist ~ speed, data = cars, span = 0.5, degree = 1)
ls_yy <- predict(lmodel, se = TRUE)

#wrong - this does not denote quantiles for the input data:
ls_yy <- cbind(ls_yy$fit, 
               ls_yy$fit - 2*ls_yy$se.fit, 
               ls_yy$fit + 2*ls_yy$se.fit)
plot(cars)
matlines(xx, ls_yy, col = "darkgrey")

我们可以使用具有更多观察结果的样本数据集使这一点更加明显。样本 1 和 2 是相同的，除了样本大小（500 和 1500 个观测值）之外，因此它们应该具有非常相似的分位数。

set.seed(1)
x1 = runif(500,0,10)
y1 = x1 + rnorm(length(x1))

x2 = runif(1500,0,10)
y2 = x1 + rnorm(length(x2))

dfpd = data.frame(x=1:9)

lmodel1 <- loess(y ~ x, data = data.frame(x=x1,y=y1), span = 0.5, degree = 1)
ls_yy1 <- predict(lmodel1, newdata=dfpd, se = TRUE)

lmodel2 <- loess(y ~ x, data = data.frame(x=x2,y=y2), span = 0.5, degree = 1)
ls_yy2 <- predict(lmodel2, newdata=dfpd, se = TRUE)

#the only difference between lmodel1 and lmodel2 is the number of observations
#the quantiles should be very similar, but their se values are a function of sample
#size and are thus quite different
ls_yy1$se
ls_yy2$se


ls_yy1$se / ls_yy2$se

我们可以看到 se 值的比率约为 60%，这证实了它们不能按原样用于分位数计算

se值之比

于 2021-03-24T23:52:44.130 回答