标签: quantreg

当我在 R 中运行分位数回归时，使用quantreg包，然后运行summary(quantregObject)，我收到此错误消息：

base::backsolve(r, x, k = k, upper.tri = upper.tri, transpose = transpose, : 'backsolve' 中的奇异矩阵中的错误。对角线 [1] 中的第一个零

任何建议我怎么能解决这个问题？

我想获得分位数回归的回归系数的95%置信区间.您可以使用R中包的rq功能quantreg(与OLS模型相比)计算分位数回归:

library(quantreg)
LM<-lm(mpg~disp, data = mtcars)
QR<-rq(mpg~disp, data = mtcars, tau=0.5)

我可以使用confint函数获得线性模型的95%置信区间:

confint(LM)

当我使用分位数回归时,我理解以下代码会产生引导标准错误:

summary.rq(QR,se="boot")

但实际上我想要95%的置信区间.也就是说,有些东西可以解释为:"概率为95%,间隔[...]包括真实系数".当我使用summary.lm()计算标准误差时,我只需乘以SE*1.96并得到与confint()类似的结果.但是使用自举标准错误是不可能的.所以我的问题是如何获得分位数回归系数的95%置信区间？

我使用quantregR 中的包来计算分位数回归模型。在模型中，因变量（Y）为NAS_DELAY，自变量（Xs）为SEANSON1TO4、SEANSON2TO4、SEANSON3TO4。
模型是：

NAS_DELAY=aSEANSON1TO4+bSEANSON2TO4+cSEANSON3TO4+d

的SEANSON1TO4，SEANSON2TO4，SEANSON3TO4是虚拟变量，0或1。我用R计算表明，截距和其它回归系数，但结果

“rq.fit.br(x,y,tau=tau,....) 单一设计矩阵中的错误；此外：警告消息 1:in summary.rq(xi,....):278951 non-positivr菲斯”。

我不明白为什么。

"fit2<-summary(rq(NAS_DELAY ~SEASON1TO4+SEASON2TO4+SEASON3TO4,tau=c(0.1,0.2,0.3,0.4,0.5),data=fddata))
Error in base::backsolve(r, x, k = k, upper.tri = upper.tri, transpose = transpose,  :   singular matrix in 'backsolve'. First zero in diagonal [1]"
In addition: Warning messages:
1: In rq.fit.br(x, y, tau = tau, ...) : Solution may be nonunique
2: In rq.fit.br(x, y, tau = tau, ...) : Solution may …

使用 R 和包quantreg，我正在对我的数据执行分位数回归分析。

我可以使用汇总函数中的 se（标准误差）估计器来访问 p 值，如下所示，但是我只能得到 5 位小数，并且想要更多。

model <- rq(outcome ~ predictor)
summary(model, se="ker")

Call: rq(formula = outcome ~ predictor)

tau: [1] 0.5

Coefficients:
            Value    Std. Error t value  Pr(>|t|)
(Intercept) 78.68182  2.89984   27.13312  0.00000
predictor    0.22727  0.03885    5.84943  0.00000

如何获取 p 值的更多小数位？

更新

好的，所以我可以通过选择包含数值结果矩阵的子对象来获得更多的小数位；

> summary(model, se="ker")[[3]]
                 Value Std. Error   t value     Pr(>|t|)
(Intercept) 78.6818182 3.13897835 25.066059 0.000000e+00
predictor    0.2272727 0.04105681  5.535567 4.397638e-08

然而，当该值 <1e-12 时，P 值仍四舍五入为 0（以上输出是简化的示例模型）。通过应用 @seancarmody 的建议，我可以得到更多；

format(summary(model, se="ker")[[3]], digits=22)

但如果 P …

我正在使用该quantreg包在R中运行以下分位数回归:

bank <-rq(gekX~laggekVIXclose+laggekliquidityspread+lagdiffthreeMTBILL+
lagdiffslopeyieldcurve+lagdiffcreditspread+laggekSPret, tau=0.99)

并通过提取系数和摘要统计

bank$coefficients
summary(bank)

我得到的结果是

Call: rq(formula = gekX ~ laggekVIXclose + laggekliquidityspread + 
lagdiffthreeMTBILL + lagdiffslopeyieldcurve + lagdiffcreditspread + 
laggekSPret, tau = 0.99)

tau: [1] 0.99

Coefficients:
                       Value    Std. Error t value  Pr(>|t|)
(Intercept)            -0.03005  0.01018   -2.95124  0.00319
laggekVIXclose          0.00471  0.00069    6.81515  0.00000
laggekliquidityspread  -0.01295  0.01619   -0.79976  0.42392
lagdiffthreeMTBILL     -0.12273  0.12016   -1.02136  0.30717
lagdiffslopeyieldcurve -0.13100  0.06457   -2.02876  0.04258
lagdiffcreditspread    -0.21198  0.15659   -1.35377  0.17592
laggekSPret            -0.01205  0.46559   -0.02588  0.97936

但是,我想知道R ^ 2 /调整后的R ^ …

我想使用该quantreg包从分位数回归中提取系数以及上限和下限。这是帮助文件中的示例。

data(engel)
attach(engel)
taus <- c(.05,.1,.25,.75,.9,.95)
f <- rq((foodexp)~(income),tau=taus)
sf <- summary(f)
sf[1]
#[[1]]

#Call: rq(formula = (foodexp) ~ (income), tau = taus)

#tau: [1] 0.05

#Coefficients:
#            coefficients lower bd  upper bd 
#(Intercept) 124.88004     98.30212 130.51695
#income        0.34336      0.34333   0.38975

我知道我可以用来coefficients()获取系数。

cf <- t(data.frame(coefficients(f)))    # transpose for better arrangement
cf
#              (Intercept)    income
#tau..0.05   124.88004 0.3433611
#tau..0.10   110.14157 0.4017658
#tau..0.25    95.48354 0.4741032
#tau..0.75    62.39659 0.6440141
#tau..0.90    67.35087 0.6862995
#tau..0.95    64.10396 0.7090685

但我不知道如何获得出现在summary(). 我看了 …

这是我可重复的数据:

library("ggplot2")
library("ggplot2movies")
library("quantreg")    
set.seed(2154)
msamp <- movies[sample(nrow(movies), 1000), ]

我正在尝试熟悉stat_quantile,但文档中的示例提出了几个问题.

mggp <- ggplot(data=msamp, mapping=aes(x=year, y=rating)) + 
    geom_point() + 
    stat_quantile(formula=y~x, quantiles=c(0, 0.25, 0.50, 0.75, 1)) + 
    theme_classic(base_size = 12) + 
    ylim(c(0,10))
mggp

1. 根据我的理解,分位数将数据分成小于定义的截止值的部分,对吗？如果我像下面的代码中那样定义分位数,我得到五行.为什么？它们代表什么？

2. 似乎分位数是基于y轴上的因变量(评级)计算的.有可能扭转这种局面吗？我的意思是根据"年份"中的分位数分割数据？

我想将中值样条和相应的置信区间带添加到ggplot2散点图中.我正在使用'quantreg'包,更具体地说是rqss函数(Additive Quantile Regression Smoothing).

在ggplot2我能够添加中值样条曲线,但不是置信区间带:

fig = ggplot(dd, aes(y = MeanEst, x = N, colour = factor(polarization)))
fig + stat_quantile(quantiles=0.5, formula = y ~ qss(x), method = "rqss") +
  geom_point()

该quantreg-package本身自带的绘图功能; plot.rqss.我可以在哪里添加置信区间(bands=TRUE):

plot(1, type="n", xlab="", ylab="", xlim=c(2, 12), ylim=c(-3, 0)) # empty plot
plotfigs = function(df) {
  rqss_model = rqss(df$MeanEst ~ qss(df$N))
  plot(rqss_model, bands=TRUE, add=TRUE, rug=FALSE, jit=FALSE)
  return(NULL)
}
figures = lapply(split(dd, as.factor(dd$polarization)), plotfigs)

然而,quantreg-package …

我有一个大约250个年度最大降雨量测量的数据系列,最大值[,]并且想要一次性对所有系列应用分位数回归并获得R中每个回归模型的显着性.

library(quantreg)


qmag <- array(NA, c(250,4))
taus <- c(0.05, 0.1, 0.95, 0.975)

for(igau in 1:250){
qure <- rq(maxima[,igau+1]~maxima[,1], tau=taus)
qmag[igau,] <- coef(qure)[2,]

}

我试过了

summary(qure, se="boot")$p.value
ci(qure)

和其他类似的变化但获得NULL值.实际上是否可以自动从quantreg中将p值提取到表中,而不是仅在每个模型的summary()中单独查看它们？

我有一个分位数回归模型,我有兴趣估计.25,.5和.875分位数的效果.我的模型中的系数彼此不同,其方式符合我的模型的实质性实质理论.

下一步是测试一个分位数的特定解释变量的系数是否与另一个分位数的估计系数显着不同.我该如何测试？此外,我还想测试给定分位数的该变量的系数是否与OLS模型中的estimnate显着不同.我怎么做？

我对任何答案感兴趣,虽然我更喜欢一个涉及R的答案.这里有一些测试代码:(注意:这不是我的实际模型或数据,但是一个简单的例子,因为数据在R安装中可用)

data(airquality)
library(quantreg)
summary(rq(Ozone ~ Solar.R + Wind + Temp, tau = c(.25, .5, .75), data = airquality, method = "br"), se = "nid")

tau: [1] 0.25

Coefficients:
            Value     Std. Error t value   Pr(>|t|) 
(Intercept) -69.92874  12.18362   -5.73957   0.00000
Solar.R       0.06220   0.00917    6.77995   0.00000
Wind         -2.63528   0.59364   -4.43918   0.00002
Temp          1.43521   0.14363    9.99260   0.00000

Call: rq(formula = Ozone ~ Solar.R + Wind + Temp, tau = c(0.25, 0.5, 
    0.75), data = airquality, method = "br")

tau: [1] 0.5 …