разделение данных и подгонка в R

У меня есть следующие данные об осадках за 1951 год:

dat.1951=c(122,122,122,122,122,122,122,122,122,122,122,121,121,121,121,120,119,119,117,117,117,115,115,115,114,112,112,111,110,109,106,105,104,103,102,99,97,95,91,89,88,86,84,83,83,82,82,79,77,77,76,74,74,72,72,71,70,69,66,65,64,61,61,58,56,56,54,53,51,49,48,47,46,46,46,45,42,40,39,38,37,36,36,35,34,33,33,32,30,30,29,28,28,27,25,25,23,22,21,20,20,20,20,20,19,19,18,18,18,16,16,15,15,15,15,15,14,14,14,14,14,14,14,14,14,14,14,13,13,12,12,11,11,11,11,11,11,11,11,11,11,11,11,11,11,10,10,10,9,8,8,8,8,8,8,8,8,8,8,8,8,7,7,6,6,6,6,5,5,5,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,3,3,3,3,3,2,2,2,2,2,2,2,2,2,2,2,1,1)

Я хочу соответствовать этим данным. Я разбиваю эти данные на 2 области (голова и хвост). Тот, где баллов меньше 100 (голова), а остальные (ибо > 100) - хвост. Я мог бы поместить экспоненту в головную часть (см. код ниже). Для хвоста я хочу подобрать кривую и хочу построить обе части на одном графике вместе с данными. Кто-нибудь может помочь?

dat.1951<-dat.1951[dat.1951 > 0]
dat.1951.tail<-dat.1951[dat.1951 >= 100]
dat.1951.head<-dat.1951[dat.1951 < 100]
x.head<-seq(1,length(dat.1951.head))
log.data<-log(dat.1951.head)
idf.head<-data.frame(x.head,dat.1951.head)
idf.head$dat.1951.head<-log(idf.head$dat.1951.head)
model=lm(idf.head$dat.1951.head ~ idf.head$x.head,data=idf.head)
summary(model)

plot(dat.1951.head)
lines(idf.head$x.head,exp(fitted(model)),col="blue")

r data-fitting
person jkp    schedule 10.06.2013    source источник
comment
Я хочу соответствовать кривой. Что это значит?   -  person Roland    schedule 10.06.2013

Ответы (1)


arrow_upward
1
arrow_downward

Я не уверен, почему вы хотите (1) разбить данные на два региона, (2) исключить записи, где не было осадков, и (3) подогнать описанную вами модель. Вы можете проконсультироваться со статистиком по этим вопросам.

Однако, чтобы ответить на ваш вопрос, я придумал пример для второй модели и показываю соответствие обеих моделей на одном графике.

x <- seq(dat.1951)
sel <- dat.1951 >= 100
model1 <- lm(dat.1951[sel] ~ poly(x[sel], 2))
model2 <- lm(log(dat.1951[!sel]) ~ x[!sel])

plot(dat.1951, cex=1.5)
lines(x[sel], fitted(model1), col="blue", lwd=3)
lines(x[!sel], exp(fitted(model2)), col="navy", lwd=3)

Ради смеха я добавил третью модель, которая подходит ко всем данным, с помощью обобщенной аддитивной модели, используя функцию gam() из пакета mgcv.

library(mgcv)
model3 <- gam(dat.1951 ~ s(x))
lines(x, fitted(model3), col="orange", lwd=3, lty=2)
person Jean V. Adams    schedule 10.06.2013