嗨,我试着使用R中的分段包来拟合分段线性回归模型来估计我的数据中的断点。我使用了下面的代码来获得这个图形。
library(segmented)
set.seed(5)
x <- c(1:10, 13:22)
y <- numeric(20)
## Create first segment
y[1:10] <- 20:11 + rnorm(10, 0, 1.5)
## Create second segment
y[11:20] <- seq(11, 15, len=10) + rnorm(10, 0, 1.5)
## fitting a linear model
lin.mod <- lm(y~x)
segmented.mod <- segmented(lin.mod, seg.Z = ~x, psi=15)
summary(segmented.mod)
plot(x,y, pch=".",cex=4,xlab="x",ylab="y")
plot(segmented.mod, add=T, lwd = 3,col = "red")

理论计算表明,断点附近两条线的斜率应该相等,但在符号上相反。我刚开始学习lm和glms。我希望能有一种方法来估计受此关系约束的斜坡的断点,slope2 1=-slope2 2。
发布于 2016-02-13 20:21:25
分段包中不支持这一点。
nls2与"plinear-brute"算法可结合使用。在输出中,.lin1和.lin2分别是常数项和斜率。这将尝试x范围内的每个值,作为一个可能的bp来拟合每个值的线性回归。
library(nls2)
st <- data.frame(bp = seq(min(x), max(x)))
nls2(y ~ cbind(1, abs(x - bp)), start = st, alg = "plinear-brute")给予:
Nonlinear regression model
model: y ~ cbind(1, abs(x - bp))
data: parent.frame()
bp .lin1 .lin2
14.000000 9.500457 0.709624
residual sum-of-squares: 45.84213
Number of iterations to convergence: 22
Achieved convergence tolerance: NA下面是另一个例子,它可以澄清这一点,因为它是从相同的模型中生成数据的:
library(nls2)
set.seed(123)
n <- 100
bp <- 25
x <- 1:n
y <- rnorm(n, 10 + 2 * abs(x - bp))
st <- data.frame(bp = seq(min(x), max(x)))
fm <- nls2(y ~ cbind(1, abs(x - bp)), start = st, alg = "plinear-brute")给予:
> fm
Nonlinear regression model
model: y ~ cbind(1, abs(x - bp))
data: parent.frame()
bp .lin1 .lin2
25.000 9.935 2.005
residual sum-of-squares: 81.29
Number of iterations to convergence: 100
Achieved convergence tolerance: NA注:在上面的中,我们假设bp是x范围内的一个整数,但是如果不需要这样的条件,我们可以通过使用nls2的结果作为nls优化的起始值,即nls(y ~ cbind(1, abs(x - bp)), start = coef(fm)[1], alg = "plinear")来放松。
https://stackoverflow.com/questions/35383554
复制相似问题