ggplot置信帶地預測$適應並預測$ se.fit

我有以下變量：

prod：正整數

tenure：正數

cohort：因素

以下是這些規格的一些模擬數據。

set.seed(123) 
my_data <- data.frame(prod = rnbinom(10000, mu = 2.5, size = 1.5), 
         tenure = rexp(10000), 
         cohort = factor(sample(2011:2014, size = 10000, replace = TRUE, 
              prob = c(0.17, 0.49, 0.26, 0.08))))

我使用mgcv:gam符合以下模型：

library(mgcv) 
mod <- gam(prod ~ s(tenure, by = cohort) + cohort, data = my_data, family = nb())

將得到預測和他們的標準誤差：

preds <- predict(mod, se.fit = TRUE) 
my_data <- data.frame(my_data, 
         mu = exp(preds$fit), 
         low = exp(preds$fit - 1.96 * preds$se.fit), 
         high = exp(preds$fit + 1.96 * preds$se.fit))

這是相當容易使用package:ggplot2獲得平滑預測mu爲每個隊列（同時也迫使平滑器具有正值）：

library(magrittr) 
library(ggplot2) 
library(splines) 
my_plot <- 
    ggplot(my_data, aes(x = tenure, y = mu, color = cohort)) %>% 
    + geom_smooth(method = "glm", 
       formula = y ~ ns(x, 3), 
       family = "quasipoisson", 
       fill = NA)

但我想從GAM平滑信心帶。我如何添加這些？

沒有答案

刪除fill = NA。不。那些置信區間將是無限小的，因爲任職期間的預測在一個隊列中是完全相同的。
添加致電geom_ribbon(aes(x = tenure, ymin = low, ymax = high))。不。這給了我一個超級搖擺的，不平滑的信心樂隊。
使用package:ggvis！沒有package:ggvis的答案，請除非ggplot2沒有辦法做到這一點。我目前的繪圖框架是ggplot2，我現在堅持使用它，除非我必須切換才能做這個陰謀。

來源

2014-10-09 Brash Equilibrium

（如果你想運行的代碼，你應該提供數據。） – 2014-10-09 19:17:39

'geom_smooth'應該工作，如[的GGPLOT2示例]證明（http://docs.ggplot2.org/current /geom_smooth.html）。一個側面說明，爲什麼在你的ggplot2調用中使用'％>％'？我希望有一個加號。 – 2014-10-09 19:34:44

謝謝。正如你所看到的，我包括一個'+'運算符。我使用'％>％'，這樣我就可以在不同的行上繪製元素。我認爲這就是'ggvis'所做的（是的，我會在某一天切換到'ggvis'）。 – 2014-10-09 19:43:09

這對我有效。

require(ggplot2) 
require(mgcv) 

set.seed(123) 
my_data <- data.frame(prod = rnbinom(10000, mu = 2.5, size = 1.5), 
         tenure = rexp(10000), 
         cohort = factor(sample(2011:2014, size = 10000, replace = TRUE, 
              prob = c(0.17, 0.49, 0.26, 0.08)))) 
mod <- gam(prod ~ s(tenure, by = cohort) + cohort, data = my_data, family = nb()) 
preds <- predict(mod, se.fit = TRUE) 
my_data <- data.frame(my_data, 
         mu = exp(preds$fit), 
         low = exp(preds$fit - 1.96 * preds$se.fit), 
         high = exp(preds$fit + 1.96 * preds$se.fit)) 

ggplot(my_data, aes(x = tenure, y = prod, color = cohort)) + 
    geom_point() + 
    geom_smooth(aes(ymin = low, ymax = high, y = mu), stat = "identity")

enter image description here

來源

2014-10-10 06:40:32

ggplot置信帶地預測$適應並預測$ se.fit

回答

相關問題