Постройте доверительный интервал квадратичного на полулогарифмическом участке - PullRequest
0 голосов
/ 06 июня 2018

При заданном фрейме данных Data в форме

     x          y
1  250 1.00000000
2  345 0.03567766
3  290 0.16654457
4  260 0.58363858
5  270 0.38754579
6  280 0.24713065
7  290 0.17142857
8  300 0.11709402
9  310 0.09047619
10 320 0.06439560
11 330 0.05098901

Я могу вывести и построить график соответствия для данных с помощью

library(ggplot2)
Data$x2<-Data$x^2
quadratic.model <- lm(log(Data$y) ~ Data$x + Data$x2)

fun_quad <- function(x){return(exp(
    quadratic.model$coef[[3]] * x ^ 2 + 
    quadratic.model$coef[[2]] * x    + 
    quadratic.model$coef[[1]]  
    ))}

chartObj <- ggplot() +

    stat_function(
        fun         = fun_quad,
        aes(color   = factor(0)),
        size        = 1.3, 
        linetype    = "dotdash"
        )+  

    geom_point(data = Data,
        aes(x = x, y = y, fill = factor(0)),
        color = "black", shape = 22, stroke = 0.7, size = 2.2) +

    coord_trans(y = 'log10', 
            limx = c(250,350), limy = c(.025,1))+ 

    theme_bw() + 

    guides(fill=F,color=F,linetype=F)

chartObj

, что делает

generic plot.

Я также попытался построить CI, используя confint и geom_ribbon.

ribbon.ymin <- function(x){return(exp(
    confint(quadratic.model)[[3]]*x^2 + 
    confint(quadratic.model)[[2]]*x   + 
    confint(quadratic.model)[[1]] 
    ))}

ribbon.ymax <- function(x){return(exp(
    confint(quadratic.model)[[6]]*x^2 + 
    confint(quadratic.model)[[5]]*x   + 
    confint(quadratic.model)[[4]]
    ))}


ribbonData <- as.data.frame(cbind(x = seq(250,350,.01)))
attach(ribbonData)
ribbonData$ymin     <- ribbon.ymin(x)
ribbonData$ymax     <- ribbon.ymax(x)
ribbonData$y        <- fun_quad(x)
detach(ribbonData)
head(ribbonData)

chartObj <- chartObj + 

    geom_ribbon( data = ribbonData,
            aes(x = x, y = 0:0,
                ymin = ymin, ymax = ymax,
                color = factor(0),fill = factor(0)),
            alpha = 0.3) 

однако, это выглядит так, как показано ниже, что снова кажется явно неправильным.

geom_ribbon plot

Итак, какЯ строю доверительный интервал, связанный с функцией, описанной quadratic.model?

Обновление

Я думаю, что я нашел почти то, что я ищу с использованием *Команда 1034 *, в частности, показана ниже, однако это все же оставляет желать лучшего, особенно неровности краев получаемой ленты.

Data$x2<-Data$x^2
quadratic.model <- lm(log(Data$y) ~ Data$x + Data$x2)

fun_quad <- function(x){return(exp(
    quadratic.model$coef[[3]] * x ^ 2 + 
    quadratic.model$coef[[2]] * x    + 
    quadratic.model$coef[[1]]  
    ))}

ribbonData<-predict(quadratic.model,data.frame(x=Data$x),interval="predict",level=.95)
# "predict" used over "confidence" in this example to show the rough edges better.
ribbonData<-as.data.frame(cbind(x=Data$x,fit=ribbonData[,1],lower=ribbonData[,2],upper=ribbonData[,3]))
ribbonData[,2:4]<-exp(ribbonData[,2:4])

chartObj <- ggplot() +

    geom_ribbon( data = ribbonData,
        aes(x = x, y = fit,
            ymin = lower, ymax = upper,
            color = factor(0),fill = factor(0)),
        alpha = 0.3) +

    stat_function(
        fun         = fun_quad,
        aes(color   = factor(0)),
        size        = 1.3, 
        linetype    = "dotdash"
        )+  

    geom_point(data = Data,
        aes(x = x, y = y, fill = factor(0)),
        color = "black", shape = 22, stroke = 0.7, size = 2.2) +

    coord_trans(y = 'log10', 
            limx = c(250,350), limy = c(.025,1))+ 

    theme_bw() + 

    guides(fill=F,color=F,linetype=F)

predict plot

Есть ли лучший способ представить информацию, представленную на графике выше?Чтобы сгладить неровные края ленты?

1 Ответ

0 голосов
/ 06 июня 2018

Это может «чувствовать себя явно неправильно», но оно показывает, что его спросили.Весь интервал не виден, потому что были установлены limx и limy:

ribbon <- function(x, level = 0.95) {
  data.frame(
    x,
    ymin = exp(
      confint(quadratic.model, level = level)[[3]] * x ^ 2 + 
        confint(quadratic.model, level = level)[[2]] * x + 
        confint(quadratic.model, level = level)[[1]] 
    ),
    ymax = exp(
      confint(quadratic.model, level = level)[[6]]*x^2 + 
        confint(quadratic.model, level = level)[[5]]*x   + 
        confint(quadratic.model, level = level)[[4]]
    )
  )
}

chartObj +
  coord_trans(y = 'log10') +
  geom_ribbon(data = ribbon(seq(250, 350, .01), level = 0.95),
              aes(x = x, ymin = ymin, ymax = ymax,
                  color = factor(0), fill = factor(0)),
              alpha = 0.3)

enter image description here

(Примечание: мой ответстрого о программировании с ggplot2 и ничего не говорит о статистической достоверности возведения в степень доверительного интервала).


Редактировать в ответ на обновленный вопрос OP (сгладить края ленты).

predict() больше очков:

quadratic.model <- lm(log(y) ~ x + x2, data = Data)

ribbonData <- data.frame(x = seq(250, 350, 0.01), x2 = seq(250, 350, 0.01) ^ 2)
ribbonData <- cbind(
  ribbonData,
  predict(quadratic.model, ribbonData,
          interval = "prediction", level = 0.95)
)
# "predict" used over "confidence" in this example to show the rough edges better.
ribbonData[, 3:5] <- exp(ribbonData[, 3:5])

ggplot() +

  geom_ribbon( data = ribbonData,
               aes(x = x, y = fit,
                   ymin = lwr, ymax = upr,
                   color = factor(0),fill = factor(0)),
               alpha = 0.3) +

  stat_function(
    fun         = fun_quad,
    aes(color   = factor(0)),
    size        = 1.3, 
    linetype    = "dotdash"
  ) +  

  geom_point(data = Data,
             aes(x = x, y = y, fill = factor(0)),
             color = "black", shape = 22, stroke = 0.7, size = 2.2) +

  coord_trans(y = 'log10', 
              limx = c(250, 350), limy = c(.025, 1)) + 

  theme_bw() + 

  guides(fill = F, color = F, linetype = F)

enter image description here

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...