2017-01-10 7 views
0

データを区分線形モデルに適合させようとしています。セグメント化されたパッケージは、係数に制約を適用する方法が見つからない限り、本当にうまく機能します。線形セグメントの傾きをある範囲の値に制限する必要があります(たとえば、セグメント1では0〜0.1、セグメント2では> 0.5)。ここに制約のない例があります。ありがとう。セグメント化された区分線形フィットの係数を制限する

library(segmented) 

#generate data for fit 
c <- 1 
m <- 0.47 #slope of 2nd line 
d <- 4 
n<- 30 
sd <- 0.2 
b <- c-m*d 
x<- runif(n,0,10) 
y<- ifelse(x<=d,c+rnorm(n,0,sd),m*x+b+rnorm(n,0,sd)) #piecewise data for fit 
plot(x,y) 

lin.mod <- lm(y~x) 
segmented.mod <- segmented(lin.mod, seg.Z = ~x, psi=6) 
summary(segmented.mod) 
plot(segmented.mod, add=T) 

答えて

0

これはエレガントな方法ではありませんが、私はconstrOptimで解決策を管理しました。

#generate data for fit 
set.seed(2) 
c <- 1 
m <- 0.47 #slope of 2nd line 
d <- 4 
n<- 30 
sd <- 0.2 
b <- c-m*d 
x<- runif(n,0,10) 
y<- ifelse(x<=d,c+rnorm(n,0,sd),m*x+b+rnorm(n,0,sd)) #piecewise data for fit 
plot(x,y) 

#define slope value constriants 
l.slope.min <- 0 
l.slope.max <- 0.1 
r.slope.min <- 0.5 
r.slope.max <- 0.99 

# par[1] - transition point for pieces. 
# par[2] - slope of left piece. 
# par[3] - intercept of left piece. 
# par[4] - slope of right piece. 
# the intercept of the right piece is par[1]*(par[2]-par[4])+par[3] 
f <- function (x,par) {ifelse(x<par[1],par[2]*x+par[3],par[4]*x+(par[1]*(par[2]-par[4])+par[3]))} 
r2 <- function(data,par) {with(data, sum(((y-f(x,par))^2)))} 

#to set up the constraint matrix and vector follow instructions in constrOptim documentation 
ui <- cbind(c(0,0,0,0),c(-1,1,0,0),c(0,0,0,0),c(0,0,1,-1)) 
ci <- c(-l.slope.max,l.slope.min,r.slope.min,-r.slope.max) 
ci <- ci - 1e-6 #displace startinng point from boundry 
par <- c(5, 0.05,1.2, 0.8) #initial parameters 

(ui %*% par)-ci >=0 #confirm that starting values are in the feasible region 

result <-constrOptim(theta =par, f=r2 , grad= NULL, ui=ui, ci=ci,data = data.frame(x,y)) 

abline(a = result$par[3], b = result$par[2], col = "red") 
abline(a = result$par[1]*(result$par[2]-result$par[4])+result$par[3], b = result$par[4], col = "blue") 
関連する問題