2017-01-10 5 views
0

Ich versuche Daten an ein stückweise lineares Modell anzupassen. Das segmentierte Paket funktioniert wirklich gut, außer dass ich keine Möglichkeit finde, Constraints auf die Koeffizienten anzuwenden. Ich muss die Steigungen der linearen Segmente auf einen Bereich von Werten beschränken (z. B. 0-0,1 für Segment 1 und> 0,5 für Segment 2). Hier ist ein Beispiel ohne Einschränkungen. Vielen Dank.constrain Koeffizienten für segmentierte stückweise lineare Anpassung

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) 

Antwort

0

Dies ist vielleicht nicht der elegante Weg, dies zu tun, aber ich schaffte eine Lösung mit 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") 
Verwandte Themen