使者兄，根据软件提示，应该是带宽h-tau<0所导致的。summary(fit1)命令其实是通过呼叫summary.rq，然后看下summary.rq内部codes：
else if (se == "nid") {
        h <- bandwidth.rq(tau, n, hs = hs)
        if (tau + h > 1)
            stop("tau + h > 1:  error in summary.rq")
        if (tau - h < 0)
            stop("tau - h < 0:  error in summary.rq")
        bhi <- rq.fit.fnb(x, y, tau = tau + h)$coef
        blo <- rq.fit.fnb(x, y, tau = tau - h)$coef
        dyhat <- x %*% (bhi - blo)
        if (any(dyhat <= 0))
            warning(paste(sum(dyhat <= 0), "non-positive fis"))
        f <- pmax(0, (2 * h)/(dyhat - eps))
        fxxinv <- diag(p)
        fxxinv <- backsolve(qr(sqrt(f) * x)$qr[1:p, 1:p, drop = FALSE],
            fxxinv)
        fxxinv <- fxxinv %*% t(fxxinv)
        cov <- tau * (1 - tau) * fxxinv %*% crossprod(x) %*%
            fxxinv
        scale <- mean(f)
        serr <- sqrt(diag(cov))
    }
    else if (se == "ker") {
        h <- bandwidth.rq(tau, n, hs = hs)
        if (tau + h > 1)
            stop("tau + h > 1:  error in summary.rq")
        if (tau - h < 0)
            stop("tau - h < 0:  error in summary.rq")
        uhat <- c(y - x %*% coef)
        h <- (qnorm(tau + h) - qnorm(tau - h)) * min(sqrt(var(uhat)),
            (quantile(uhat, 0.75) - quantile(uhat, 0.25))/1.34)
        f <- dnorm(uhat/h)/h
        fxxinv <- diag(p)
        fxxinv <- backsolve(qr(sqrt(f) * x)$qr[1:p, 1:p, drop = FALSE],
            fxxinv)
        fxxinv <- fxxinv %*% t(fxxinv)
        cov <- tau * (1 - tau) * fxxinv %*% crossprod(x) %*%
            fxxinv
        scale <- mean(f)
        serr <- sqrt(diag(cov))
    }
    else if (se == "boot") {
        B <- boot.rq(x, y, tau, ...)
        cov <- cov(B)
        serr <- sqrt(diag(cov))
    }


从你取的tau可以看出，你的tau+h很可能会超过1.所以解决方法有2种方法：
第一是换过一种方法产生标准误：
y=c(2167,2509,3530,4669,5868,6763,6820,2866,8248,8868,9336,10464,11040,12631,13773,14762,17255,19398,20992,23200)
income=c(2486,3009,4277,5868,7172,8159,8439,8773,10932,11718,12883,13250,14867,16683,18645,20668,23623,26675,28838,21754)
cpi=c(110.5,110.0,120.2,123.9,118.7,109.2,102.8,100.0,101.5,102.5,100.0,100.5,100.1,102.2,101.0,101.2,103.2,105.8,99.6,103.1)
library(quantreg)
fit1=rq(y~income+cpi,tau=c(0.2,0.5,0.8,0.99))
summary(fit1,se="boot")

第二种是改变tau值，我只尝试tau=0.5可以运行


y=c(2167,2509,3530,4669,5868,6763,6820,2866,8248,8868,9336,10464,11040,12631,13773,14762,17255,19398,20992,23200)
income=c(2486,3009,4277,5868,7172,8159,8439,8773,10932,11718,12883,13250,14867,16683,18645,20668,23623,26675,28838,21754)
cpi=c(110.5,110.0,120.2,123.9,118.7,109.2,102.8,100.0,101.5,102.5,100.0,100.5,100.1,102.2,101.0,101.2,103.2,105.8,99.6,103.1)
library(quantreg)
fit1=rq(y~income+cpi,tau=c(0.5))
summary(fit1,se="nid")

tau=c(0.2,0.5,0.8,0.99)
n=20
h=bandwidth.rq(tau, n, hs=TRUE, alpha=0.05)

#h
[1] 0.21062684 0.35792541 0.21062684 0.02586753
所以产生了tau + h > 1 和tau-h<0错误。

谢谢!