2013-11-01 22 views
19

kullanarak mgcv paketinden bir model takıyorum ve model sonucunu saklıyorum ve şu ana kadar plot(model) kullanarak sorunsuz bileşenlere bakıyordum. Kısa bir süre önce ggplot2 kullanmaya başladım ve çıkışını beğeniyorum. Bu yüzden merak ediyorum, bu grafikleri ggplot2 kullanarak çizmek mümkün mü? İşte Bir gamın pürüzsüz bileşenlerini ggplot2 ile çizmek mümkün mü?

bir örnek:

x1 = rnorm(1000) 
x2 = rnorm(1000) 
n = rpois(1000, exp(x1) + x2^2) 

model = gam(n ~ s(x1, k=10) + s(x2, k=20), family="poisson") 
plot(model, rug=FALSE, select=1) 
plot(model, rug=FALSE, select=2) 

Ve s(x1, k=10) ve s(x2, k=20) değil uyum içinde olan ilgiyi duyuyorum.

Kısmi cevap:

Ben plot.gam ve mgcv:::plot.mgcv.smooth derinliklerine kazılan ve pürüzsüz bileşenlerden tahmin etkileri ve standart hataları ayıklar kendi işlevi inşa etti. Tüm seçenekleri ve plot.gam durumlarını ele almıyor, bu yüzden sadece kısmi bir çözüm olarak düşünüyorum, ama benim için iyi çalışıyor.

EvaluateSmooths = function(model, select=NULL, x=NULL, n=100) { 
    if (is.null(select)) { 
    select = 1:length(model$smooth) 
    } 
    do.call(rbind, lapply(select, function(i) { 
    smooth = model$smooth[[i]] 
    data = model$model 

    if (is.null(x)) { 
     min = min(data[smooth$term]) 
     max = max(data[smooth$term]) 
     x = seq(min, max, length=n) 
    } 
    if (smooth$by == "NA") { 
     by.level = "NA" 
    } else { 
     by.level = smooth$by.level 
    } 
    range = data.frame(x=x, by=by.level) 
    names(range) = c(smooth$term, smooth$by) 

    mat = PredictMat(smooth, range) 
    par = smooth$first.para:smooth$last.para 

    y = mat %*% model$coefficients[par] 

    se = sqrt(rowSums(
     (mat %*% model$Vp[par, par, drop = FALSE]) * mat 
    )) 

    return(data.frame(
     label=smooth$label 
     , x.var=smooth$term 
     , x.val=x 
     , by.var=smooth$by 
     , by.val=by.level 
     , value = y 
     , se = se 
    )) 
    })) 
} 

Bu düz bileşenleri ile bir "erimiş" veri çerçevesini döndürür, böylece, yukarıda örnek ggplot kullanmak mümkündür: herkes bu sağlayan bir paket bilir

smooths = EvaluateSmooths(model) 

ggplot(smooths, aes(x.val, value)) + 
    geom_line() + 
    geom_line(aes(y=value + 2*se), linetype="dashed") + 
    geom_line(aes(y=value - 2*se), linetype="dashed") + 
    facet_grid(. ~ x.var) 

genel durum çok minnettar olurum.

+1

ggplot '=, geom_smooth' 'için' yani sadece yapılacak 'yöntem predict' kullanır –

+0

gam'' düzgün terimleri değil, uygunluğu gösterir. Yani bunun çözüm olduğunu düşünmüyorum. – unique2

+0

Veri kümesine bağlantı (sadece bir başlangıç ​​noktası olarak 'mgcv''den bir örnek ve çoğaltmaya çalıştığınız çizimden alıntı yapın) ve (muhtemelen) nasıl yapacağınızı gösterebiliriz. –

cevap

18

Visreg paketini, plyr paketi ile birlikte kullanabilirsiniz. visreg temel olarak, predict() öğesini kullanabileceğiniz herhangi bir modeli gösterir.

library(mgcv) 
library(visreg) 
library(plyr) 
library(ggplot2) 

# Estimating gam model: 
x1 = rnorm(1000) 
x2 = rnorm(1000) 
n = rpois(1000, exp(x1) + x2^2) 
model = gam(n ~ s(x1, k=10) + s(x2, k=20), family="poisson") 

# use plot = FALSE to get plot data from visreg without plotting 
plotdata <- visreg(model, type = "contrast", plot = FALSE) 

# The output from visreg is a list of the same length as the number of 'x' variables, 
# so we use ldply to pick the objects we want from the each list part and make a dataframe: 
smooths <- ldply(plotdata, function(part) 
    data.frame(Variable = part$meta$x, 
      x=part$fit[[part$meta$x]], 
      smooth=part$fit$visregFit, 
      lower=part$fit$visregLwr, 
      upper=part$fit$visregUpr)) 

# The ggplot: 
ggplot(smooths, aes(x, smooth)) + geom_line() + 
    geom_line(aes(y=lower), linetype="dashed") + 
    geom_line(aes(y=upper), linetype="dashed") + 
    facet_grid(. ~ Variable, scales = "free_x") 

Biz bir işlev içine her şeyi koymak ve modelde (res = DOĞRU) den artığı göstermek için bir seçenek ekleyebilirsiniz: ggplot with residuals Renkler http://colorbrewer2.org/ den toplanır

ggplot.model <- function(model, type="conditional", res=FALSE, 
         col.line="#7fc97f", col.point="#beaed4", size.line=1, size.point=1) { 
    require(visreg) 
    require(plyr) 
    plotdata <- visreg(model, type = type, plot = FALSE) 
    smooths <- ldply(plotdata, function(part) 
    data.frame(Variable = part$meta$x, 
      x=part$fit[[part$meta$x]], 
      smooth=part$fit$visregFit, 
      lower=part$fit$visregLwr, 
      upper=part$fit$visregUpr)) 
    residuals <- ldply(plotdata, function(part) 
    data.frame(Variable = part$meta$x, 
       x=part$res[[part$meta$x]], 
       y=part$res$visregRes)) 
    if (res) 
    ggplot(smooths, aes(x, smooth)) + geom_line(col=col.line, size=size.line) + 
     geom_line(aes(y=lower), linetype="dashed", col=col.line, size=size.line) + 
     geom_line(aes(y=upper), linetype="dashed", col=col.line, size=size.line) + 
     geom_point(data = residuals, aes(x, y), col=col.point, size=size.point) + 
     facet_grid(. ~ Variable, scales = "free_x") 
    else 
    ggplot(smooths, aes(x, smooth)) + geom_line(col=col.line, size=size.line) + 
     geom_line(aes(y=lower), linetype="dashed", col=col.line, size=size.line) + 
     geom_line(aes(y=upper), linetype="dashed", col=col.line, size=size.line) + 
     facet_grid(. ~ Variable, scales = "free_x") 
    } 

ggplot.model(model) 
ggplot.model(model, res=TRUE) 

ggplot without residuals .

+0

Artık bir geçici dosyaya çizim yapmak yerine, herhangi bir şey göstermeden çizim verilerini döndürmek için 'visreg 'öğesine' plot = FALSE' argümanını kullanabilirsiniz. Ancak, iade edilen nesnenin, varsaydığın şeyden değiştiğini düşünüyorum. – Spacedman

+0

Mesajların büyük olasılıkla bir güncellemeye ihtiyacı var. Yukarıdaki kodu çalıştırdığımda 'plyr :: ldply()' çağrısıyla yanlış bir şey varsa, 'smooths' nesnesini döndürür. –

+0

@ pat-s Teşekkürler, haklısınız. Gönderi şimdi güncellendi ve çalışmalı. –

2

Bilginize, visreg ile doğrudan çıkış gg nesne olabilmektedir: geom_smooth anlamak gibi

visreg(model, "x1", gg=TRUE) 

enter image description here

İlgili konular