散点图的非参数分位数回归曲线

Posted 2023-03-24

【中文标题】散点图的非参数分位数回归曲线

【英文标题】：Nonparametric quantile regression curves to scatterplot

【发布时间】：2013-02-07 13:29:03

【问题描述】：

我用IV=time、DV=concentration 创建了一个散点图（多组 GRP）。我想将分位数回归曲线(0.025,0.05,0.5,0.95,0.975) 添加到我的绘图中。

顺便说一句，这就是我创建散点图的方法：

attach(E)  ## E is the name I gave to my data
## Change Group to factor so that may work with levels in the legend
Group<-as.character(Group)
Group<-as.factor(Group)

## Make the colored scatter-plot
mycolors = c('red','orange','green','cornflowerblue')
plot(Time,Concentration,main="Template",xlab="Time",ylab="Concentration",pch=18,col=mycolors[Group])

## This also works identically
## with(E,plot(Time,Concentration,col=mycolors[Group],main="Template",xlab="Time",ylab="Concentration",pch=18))

## Use identify to identify each point by group number (to check)
## identify(Time,Concentration,col=mycolors[Group],labels=Group)
## Press Esc or press Stop to stop identify function

## Create legend
## Use locator(n=1,type="o") to find the point to align top left of legend box
legend('topright',legend=levels(Group),col=mycolors,pch=18,title='Group')

由于我在此处创建的数据是我的较大数据的一小部分，因此它可能看起来可以近似为矩形双曲线。但我还不想称我的自变量和因变量之间存在数学关系。

我认为quantreg包中的nlrq可能是答案，但是当我不知道我的变量之间的关系时，我不明白如何使用该函数。

我从一篇科学文章中找到了这张图，我想做完全相同的图：

再次感谢您的帮助！

更新

Test.csv 有人指出我的样本数据不可重现。这是我的数据示例。

library(evd)
qcbvnonpar(p=c(0.025,0.05,0.5,0.95,0.975),cbind(TAD,DV),epmar=T,plot=F,add=T)

我也试过qcbvnonpar::evd，但曲线似乎不是很平滑。

【问题讨论】：

如果您无法提供自己的数据，请尝试创建一个随机数数据集并演示您的问题。向我们展示您的尝试。它为我们提供了一些可以合作的东西，同时也是我们诚信的标志。

哦。对不起——我会编一些数字。它可能相当大。

这可以帮助您生成数据。 ***.com/questions/5963269/…

几位cmets：您显示的数字似乎来自works.bepress.com/phil_reiss/16？那篇论文似乎有一个相关的 R 包（还没有看过......）要获得完全非参数的平滑分位数将相当困难（我认为）。两种可能的解决方案是（1）拟合广义加法模型（例如library(splines); rq(y~s(x,5),tau=0.9)）； (2) 使用分位数的运行估计。

也许 quantreg 的 rqss() 也适合你？

【参考方案1】：

也许看看 quantreg:::rqss 来平滑样条曲线和分位数回归。 抱歉，示例数据不太好：

set.seed(1234)
period <- 100
x <- 1:100
y <- sin(2*pi*x/period) + runif(length(x),-1,1)


require(quantreg)
mod <- rqss(y ~ qss(x))
mod2 <- rqss(y ~ qss(x), tau=0.75)
mod3 <- rqss(y ~ qss(x), tau=0.25)
plot(x, y)
lines(x[-1], mod$coef[1] + mod$coef[-1], col = 'red')
lines(x[-1], mod2$coef[1] + mod2$coef[-1], col = 'green')
lines(x[-1], mod3$coef[1] + mod3$coef[-1], col = 'green')

【讨论】：

+1 for rqss() - 如果非参数是必需的，并且示例图像表明它是基于样条的拟合，那么 rqss() 肯定是我首先要看的地方。跨度>

它与您的示例配合得非常好，但我不知道为什么，即使我检查我的 x 和 y 具有相同的 n，我仍然会收到 Error in xy.coords(x, y) : 'x' and 'y' lengths differ 警告我的数据集。仍在进行错误调试：P

你能多使用一点你的数据吗？您上面的示例数据显然不合适。

【参考方案2】：

我过去经常与rqss 斗争，而我的问题几乎总是与分数的顺序有关。

您在不同的时间点进行了多次测量，这就是您获得不同长度的原因。这对我有用：

dat <- read.csv("~/Downloads/Test.csv")

library(quantreg)
dat <- plyr::arrange(dat,Time)
fit<-rqss(Concentration~qss(Time,constraint="N"),tau=0.5,data = dat)
with(dat,plot(Time,Concentration))
lines(unique(dat$Time)[-1],fit$coef[1] + fit$coef[-1])

似乎有必要在拟合模型之前对数据框进行排序。

【讨论】：

这行得通！太感谢了。我不知道订购可能是个问题。

【参考方案3】：

如果你想要ggplot2 图形...

我的这个例子是基于@EDi 的。我增加了x 和y，这样分位数线就不会那么摆动了。由于这种增加，我需要在某些调用中使用unique(x) 代替x。

这是修改后的设置：

set.seed(1234)
period <- 100
x <- rep(1:100,each=100)
y <- 1*sin(2*pi*x/period) + runif(length(x),-1,1)


require(quantreg)
mod <- rqss(y ~ qss(x))
mod2 <- rqss(y ~ qss(x), tau=0.75)
mod3 <- rqss(y ~ qss(x), tau=0.25)

这是两个图：

# @EDi's base graphics example
plot(x, y)
lines(unique(x)[-1], mod$coef[1] + mod$coef[-1], col = 'red')
lines(unique(x)[-1], mod2$coef[1] + mod2$coef[-1], col = 'green')
lines(unique(x)[-1], mod3$coef[1] + mod3$coef[-1], col = 'green')

# @swihart's ggplot2 example:
## get into dataset so that ggplot2 can have some fun:
qrdf <- data.table(x       = unique(x)[-1],
                   median =  mod$coef[1] +  mod$coef[-1],
                   qupp   = mod2$coef[1] + mod2$coef[-1],
                   qlow   = mod3$coef[1] + mod3$coef[-1]
)

line_size = 2
ggplot() +
  geom_point(aes(x=x, y=y),
             color="black", alpha=0.5) +
  ## quantiles:
  geom_line(data=qrdf,aes(x=x, y=median),
            color="red", alpha=0.7, size=line_size) +
  geom_line(data=qrdf,aes(x=x, y=qupp),
            color="blue", alpha=0.7, size=line_size, lty=1) +
  geom_line(data=qrdf,aes(x=x, y=qlow),
            color="blue", alpha=0.7, size=line_size, lty=1)