Примерка лессовой линии в R

Question

У меня возникли проблемы с использованием loess() а также predict() в R. Я использовал следующий код для моделирования данных:

Overall=0.6
RR=1.1
Noise=0.05

x=seq(from=0.01, to=10, by=0.01) 
logRR=log(RR)
logBeta0=log(Overall)


linear.pred = logBeta0 + (logRR*x) + rnorm(length(x), 0, Noise*sqrt(x))
linear.pred.high = logBeta0 + (logRR*15) + rnorm(length(x), 0, Noise/5)

PoissonRate <- function (x) ifelse(x<=9, exp(linear.pred), exp(linear.pred.high))

xyplot(PoissonRate(x)~x) #the shape of the 'raw' data


loess_fit <- loess(x~PoissonRate(x))
lines(predict(loess_fit), col = "black")

Извиняюсь, но я не могу понять, как прикрепить картинку, чтобы показать, как это выглядит!

Последние две строки кода заканчиваются просто добавлением случайной черной линии наполовину от графика, хотя, когда я ранее использовал эту команду для разных (очень похожих) данных, казалось, что она работает нормально. Что мне не хватает?! Любая помощь была бы отличной, спасибо:)

Answer 1

Ты не должен звонить llines() (если это то, что вы имели в виду lines()) за пределами xyplot позвоните, по крайней мере, в моем понимании. ?llines имеет:

Description:

     These functions are intended to replace common low level
     traditional graphics functions, primarily for use in panel
     functions.

Следовательно, один из вариантов - сделать так, как он предлагает, и создать собственную функцию панели на лету. Вот пример:

set.seed(1)
dat <- data.frame(pr = PoissonRate(x), x = x)
loess_fit <- loess(pr ~ x, data = dat)

xyplot(PoissonRate(x) ~ x,
  panel = function(x, y, ...) {
    panel.xyplot(x, y, ...)
    llines(dat$x, predict(loess_fit), col.line = "red")
  })

который производит:

В общем, я бы, наверное, сделал это по-другому - я бы генерировал данные вне формулы. Я бы predict для новых местоположений данных распределены равномерно и в порядке в диапазоне x, Таким образом, даже если входные данные не в порядке возрастания в x вы получаете разумные прогнозы, которые вы можете использовать для линейного сюжета. Например

## following on from the above
pred <- with(dat, data.frame(x = seq(min(x), max(x), length = 100)))
pred <- transform(pred, pr = predict(loess_fit, newdata = x))

xyplot(PoissonRate(x) ~ x,
  panel = function(x, y, ...) {
    panel.xyplot(x, y, ...)
    with(pred, llines(x, pr, col.line = "red"))
  })